2 二分法

関数$y=f(x)$が区間$[a,b]$で連続であるとき、 $f(a)f(b)<0$ならば、その区間に解が少なくとも一つ存在する。いま、その区間に 解が一つしか存在しない場合を考えることにする。すると、 $f(a)f(b)<0$のとき、$a$と$b$の中点 $c={{a+b}\over{2}}$に対し、

1. $f(a)f(c)<0$ならば、解の存在区間は$[a,c]$に狭まるし、逆に
2. $f(b)f(c)<0$ならば、解の存在区間は$[c,b]$に狭まる。

そこで、狭まった解の存在区間を更に二分して、同様の判定を行ない、解の存在区間 を更に狭める。以下、同様に、この操作を繰り返せば、解を得る事が出来る。 収束判定条件としては、

1. $\vert b-a\vert<\varepsilon$
2. $\vert f(c)\vert<\varepsilon$

がある。

この方法では、収束解を得るのに必要な反復計算回数$n$を、初期値$a$、$b$を用いて、

$${{\vert b-a\vert}\over{2^n}}<\varepsilon$$ (14)

から、予め知る事が出来る。したがって、$n$回の反復計算の後も、収束解が得られな い場合には、プログラムに誤りがあると判断できる。

以下の外部関数rtbisは、関数$f(x)$をfuncとして、$f(x)=0$となる解を 二分法で解くプログラム例である 。

      function rtbis(func,x1,x2,xacc)
      implicit NONE
      Integer JMAX
      Real rtbis,x1,x2,xacc,func
      external func
      Parameter (JMAX=40)      !Maximum allowed number of bisections.
!	Using bisection, find the root of a function 'func' known to lie 
!	betwee x1 and x2. The root, returned as 'rtbis', will be refined  
!	until its accuracy is +-'xacc'.
      Integer j
      Real dx,f,fmid,xmid
      fmid=func(x2)
      f=func(x1)
      if(f*fmid.ge.0.) pause 'root must be bracketed in rtbis'
      if(f.lt.0.)then          !Orient the search so that f>0 lies at x+dx.
        rtbis=x1
        dx=x2-x1
      else
        rtbis=x2
        dx=x1-x2
      endif
      do j=1,JMAX              !Bisection loop.
        dx=dx*.5
        xmid=rtbis+dx
        fmid=func(xmid)
        if(fmid.le.0.)rtbis=xmid
        if(abs(dx).lt.xacc .or. fmid.eq.0.) return
      end do
      pause 'too many bisections in rtbis'
      end