program main !==================================================================== ! A single root of an equation f(x)=0 in [x1,x2] interval ! Method: Bisectional !==================================================================== implicit none double precision f,x1,x2,root,eps integer flag external f x1 = -10.0 x2 = 10.0 eps = 1.0e-6 write (*,100) call bisection(f,x1,x2,eps,root,flag) ! print solutions if(flag == 0) then write(*,104) stop end if write(*,101) flag, eps if(flag > 0) then write(*,102) Root else write(*,103) Root end if ! test solution write(*,105) f(Root) 100 format(' A single root of f(x) ',/, & ' Method - Bisectional') 101 format(' iterations = ',i3,/,& ' tolerance = ',1pe12.5) 102 format(' root = ',1pe12.5) 103 format(' Singularity = ',1pe12.5) 104 format(' no roots for Bisectional method') 105 format(' f(root) = ',1pe12.5) end program main Subroutine bisection(f,x1,x2,eps,Root,flag) !============================================================ ! Solutions of equation f(x)=0 on [x1,x2] interval ! Method: Bisectional (closed domain) (a single root) ! Alex G. January 2010 !------------------------------------------------------------ ! input ... ! f - function - evaluates f(x) for any x in [x1,x2] ! x1 - left endpoint of initial interval ! x2 - right endpoint of initial interval ! eps - desired uncertainity of the root as |b-a|0 - a single root found, flag=number of iterations ! 0 - no solutions for the bisectional method ! <0 - not a root but singularity, flag=number of iterations ! ! Comments: Function f(x) has to change sign between x1 and x2 ! Max number of iterations - 200 (accuracy (b-a)/2**200) !==================================================================== implicit none double precision f, x1, x2, eps, Root double precision a, b, c integer i, flag integer, parameter:: iter=200 !* check the bisection condition if(f(x1)*f(x2)>0.0) then flag = 0 return end if !* initialize calculations a=x1 b=x2 !* Iterative refining the solution do i=1,iter c=(b+a)/2.0 if(f(a)*f(c).le.0.0) then b = c else a=c end if ! condition(s) to stop iterations) if(abs(b-a)<= eps) exit end do Root=(b+a)/2.0 !* check if it is a root or singularity if (abs(f(Root)) < 1.0) then flag=i else flag = -i end if end subroutine bisection !----------------------------------------------------- Function f(x) implicit none double precision f, x f = x + cos(x) end function f