台形則（数値積分）

台形則とは

• 数値積分の解法の一つ
• 関数$f(x)$において、微小区間$[x_0,x_1]$内での関数値を一次方程式で近似する

算法

微小区間$[x_0,x_1]$の区間幅$x_1-x_0=h$とし、各$y$座標を$f(x_0)=y_0,f(x_1)=y_1$とする。
区間内の関数値を一次方程式で近似すると、この区間内の積分値は
$$ \int_{x_0}^{x_1}f(x)dx=\frac{h}{2}(y_0+y_1) $$
となるので、区間$[a,b]$内の積分は各微小区間内での積分値の和で表せるから、
$$ \int_a^bf(x)dx=\frac{h}{2}(y_0+y_1)+\frac{h}{2}(y_1+y_2)+\cdots+\frac{h}{2}(y_{n-1}+y_n) $$$$ =\frac{h}{2}(y_0+2y_1+2y_2+\cdots+2y_{n-1}+y_n)$$

サンプルコード

$f(x)=\sqrt{1-x^2}$において、区間$[0,1]$の定積分の値を求めるプログラム。
分割数は4とする。
解析解は$\pi/4$です（単位円の$1/4$）。

trapezoidal_rule.c

#include<stdio.h>
#include<math.h>

double f (double x) {
  return sqrt(1-x*x);
}

/* 台形則（区間[a,b]をn分割） */
double trapezoidal_rule (double a, double b, int n) {
  double h;
  int i;
  double value=0;

  h = (b - a) / n;   // 区間幅の計算
  for (i = 0; i <= n; i++) {
    if (i == 0 || i == n) value += f(a + i*h);
    else value += 2 * f(a + i*h);
  }
  value = value*h/2;
  return value;
}

int main (void) {
  printf("Analytical solution: %f\n", M_PI/4);
  printf("Numerical solution : %f\n", trapezoidal_rule(0, 1, 4));
  return 0;
}

実行結果

Analytical solution: 0.785398
Numerical solution : 0.748927

特徴

• 全区間の誤差
  $$ E=\frac{h^2(b-a)}{12}f^{\prime\prime}(\xi) \ \ , \ \ (a<\xi<b) $$
• 分割数を倍にすると誤差は$1/4$に減少
• でも、区間数を増やしすぎると丸め誤差が増加する