====== 【数値情報処理a】第3回の答え ======
===== 課題 =====

$y=f(x)$ と $x$ 軸ではさまれた領域の面積をシンプソンの公式で数値積分するプログラムを作って提出しなさい。

  * 求める面積は下図の赤い部分である。
  * $y$ が負の領域も正の面積として求める必要がある。
  * 第2回課題の4つの解を必ず使うこと。

{{ report_integral.png?nolink |求める面積}}


===== 解答例 =====

$f(x)=-x^{4}-x^{3}+10x^{2}+x-2$ の場合

第2回課題の $f(x)=0$ の4つの解を使って積分範囲を3つに分割して積分を行う。
$f(x)<0$ の区間では積分の符号を変える。

<file rsplus report3.R[enable_line_numbers="true"]>
# シンプソンの公式

# 被積分関数の定義
f <- function(x) {
  return(-x^4 - x^3 + 10 * x^2 + x - 2)
}

# 区間数
n <- 1000

# シンプソンの公式による [a,b] の数値積分
simpson <- function(a, b) {
  # 幅
  h <- (b - a) / n

  s1 <- f(a + h * 0.5)
  s2 <- 0
  for (i in 1 : (n - 1)) {
    x <- a + i * h
    s1 <- s1 + f(x + h * 0.5)
    s2 <- s2 + f(x)
  }
  s <- (f(a) + 4 * s1 + 2 * s2 + f(b)) * h / 6
  return(s)
}

x1 <- -3.63423216189423481381
x2 <- -0.49319311896925860372
x3 <- 0.41073330471801877684
x4 <- 2.71669197614524948747

s1 <- simpson(x1, x2)
s2 <- simpson(x2, x3)
s3 <- simpson(x3, x4)
print(sprintf("S1 = %f", s1), quote = FALSE)
print(sprintf("S2 = %f", s2), quote = FALSE)
print(sprintf("S3 = %f", s3), quote = FALSE)
print(sprintf("S = %f", abs(s1) + abs(s2) + abs(s3)), quote = FALSE)
</file>

同じ処理を繰り返すために区間 $[a,b]$ の数値積分を行う ''simpson(a,b)'' という関数を定義している。
最後に負の面積を正にするため ''abs()'' を使用している。

----
==== 結果 ====

$f(x)=0$ の解やパラメーターの違いにより小数点以下の値が多少異なることがある。

$f(x)=-x^{4}-x^{3}+10x^{2}+x-2$ の場合

<code>
S = 87.251731
</code>

$f(x)=-x^{4}-2x^{3}+10x^{2}+x-2$ の場合

<code>
S = 138.341753
</code>

$f(x)=-x^{4}-3x^{3}+10x^{2}+x-2$ の場合

<code>
S = 246.731532
</code>

$f(x)=-x^{4}-4x^{3}+10x^{2}+x-2$ の場合

<code>
S = 447.533001
</code>

$f(x)=-x^{4}-5x^{3}+10x^{2}+x-2$ の場合

<code>
S = 795.972953
</code>

$f(x)=-x^{4}-6x^{3}+10x^{2}+x-2$ の場合

<code>
S = 1372.366483
</code>