function trapezoid(a,b,n) result(s)
    !
    ! 台形公式を用いた数値積分を行います. 
    !
    ! $ f(x) = x^2 - 2x + 2 $ の積分値を台形公式によって求めます.
    !
    ! * 積分区間は $ [a, b] $ とします.
    ! * 実際に計算する式は以下の通りです.
    !   ここで $ h $ は $ h = \frac{b-a}{n} $ です.
    !
    ! \[
    !    \int_{a}^{b} f(x) dx =
    !      \frac{h}{2} \left\{f(a) + f(b)\right\} + h\sum_{i=1}^{n-1} f(a+ih)
    ! \]
    !
    ! このサンプルプログラムは http://www.gfd-dennou.org/library/dcmodel
    ! にて入手できます. 
    !
    !

    real,    intent(in) :: a ! 積分区間の下限
    real,    intent(in) :: b ! 積分区間の上限
    integer, intent(in), optional :: n ! 積分区間の分割数
    real :: s  ! 積分値
    real :: h
    integer :: i, div = 10
    if (present(n)) div = n
    h = (b-a) / div
    s = (f(a)+f(b)) / 2.0
    do i = 1, n-1
      s = s + f(a + i*h)
    end do
    s = h*s
  end function trapezoid