# solutions.txt Thu Jan 14 191538 2010 1 PROGRAM roots by tfl42712

VIEWS: 16 PAGES: 2

• pg 1
```									solutions.txt            Thu Jan 14 19:15:38 2010         1                solutions.txt       Thu Jan 14 19:15:38 2010        2
PROGRAM roots                                                              ============

!                                                                             1   -1.000000    1.000000     0.000000   -1.000000
!    ******* Main program for root finding                                    2    0.000000    1.000000     0.500000   -0.500000
!                                                                             3    0.500000    1.000000     0.750000    0.031250
4    0.500000    0.750000     0.625000   -0.277344
USE main_mod !*** Use the required module                                 5    0.625000    0.750000     0.687500   -0.135254
6    0.687500    0.750000     0.718750   -0.055237
7    0.718750    0.750000     0.734375   -0.012825
!REAL :: a=-20,b=0     !** Equation One                                   8    0.734375    0.750000     0.742188    0.009002
!REAL :: a=-1,b=1      !** Equation Two                                   9    0.734375    0.742188     0.738281   -0.001964
REAL :: a=1,b=2       !** Equation Three                                 10    0.738281    0.742188     0.740234    0.003506
!REAL :: x=12                                                            11    0.738281    0.740234     0.739258    0.000768
12    0.738281    0.739258     0.738770   -0.000599
13    0.738770    0.739258     0.739014    0.000084
REAL, PARAMETER :: tol=0.0001                                            14    0.738770    0.739014     0.738892   -0.000258
Successful convergence
INTEGER :: i                                                            The root is estimated as :    0.7389526
LOGICAL, DIMENSION(2) :: test

!test=newt_raph(x,tol) !** Return root in x                            Equation Two
test=bisection(a,b,tol) !** Return root as (a+b)/2                     ============

1    1.000000    2.000000     1.500000    2.375000
IF (test(1)) THEN                                                         2    1.000000    1.500000     1.250000   -1.796875
PRINT*,"Successful convergence"                                        3    1.250000    1.500000     1.375000    0.162109
PRINT*, "The root is estimated as : ",(a+b)/2 !** Bisection            4    1.250000    1.375000     1.312500   -0.848389
!PRINT*, "The root is estimated as : ",x !** Newt-Raph                 5    1.312500    1.375000     1.343750   -0.350983
6    1.343750    1.375000     1.359375   -0.096409
ELSE                                                                      7    1.359375    1.375000     1.367188    0.032356
PRINT*,"Unsuccessful Convergence"                                      8    1.359375    1.367188     1.363281   -0.032150
ENDIF                                                                     9    1.363281    1.367188     1.365234    0.000072
10    1.363281    1.365234     1.364258   -0.016047
IF (test(2)) PRINT*, "Max iters exceeded : Unsuccessful Convergence"     11    1.364258    1.365234     1.364746   -0.007989
12    1.364746    1.365234     1.364990   -0.003960
END PROGRAM roots                                                            13    1.364990    1.365234     1.365112   -0.001944
Successful convergence
The root is estimated as :     1.365173
*******************************************
*******************************************
Newton-Raphson
==================        RESULTS     =================                    ==============

Equation One                                                                  1    6.083333   35.006947
============                                                                  2    3.206050    8.278759
3    1.914935    1.666978
1 -20.000000     0.000000        -10.000000    0.767387                    4    1.479679    0.189449
2 -20.000000 -10.000000          -15.000000    0.414978                    5    1.415662    0.004098
3 -20.000000 -15.000000          -17.500000    0.188819                    6    1.414214    0.000002
4 -20.000000 -17.500000          -18.750000    0.068040                    7    1.414214    0.000000
5 -20.000000 -18.750000          -19.375000    0.006490                  Successful convergence
6 -20.000000 -19.375000          -19.687500   -0.024472                  The root is estimated as :     1.414214
7 -19.687500 -19.375000          -19.531250   -0.008979
8 -19.531250 -19.375000          -19.453125   -0.001241                 *******************************************
9 -19.453125 -19.375000          -19.414062    0.002626                 *******************************************
10 -19.453125 -19.414062          -19.433594    0.000693
11 -19.453125 -19.433594          -19.443359   -0.000274                 MODULE main_mod
12 -19.443359 -19.433594          -19.438477    0.000209
13 -19.443359 -19.438477          -19.440918   -0.000032                   IMPLICIT NONE
14 -19.440918 -19.438477          -19.439697    0.000089
15 -19.440918 -19.439697          -19.440308    0.000028                   CONTAINS
16 -19.440918 -19.440308          -19.440613   -0.000002
17 -19.440613 -19.440308          -19.440460    0.000013                   ! ************************************************
Successful convergence
The root is estimated as :         -19.44054                                FUNCTION newt_raph(x,tol)
!**
Equation Two                                                                 !** To find the root of a function f(x)
solutions.txt       Thu Jan 14 19:15:38 2010        3                   solutions.txt         Thu Jan 14 19:15:38 2010      4
!**
!*** Function return declaration
!*** Dummy declaration                                                    REAL :: f
REAL, INTENT(INOUT) :: x
REAL, INTENT(IN) :: tol                                                   !f=COS(x/20)+SIN(x/20)**3 !** Equation One
!f=2*x**3-x**2+x-1         !** Equation Two
!*** Local variable declarations                                          f=x**3+4*x**2-10          !** Equation Three
REAL :: old_root                                                          !f=x**2-2                  !** Newton-Raph
INTEGER, SAVE :: iters=0
END FUNCTION f
!*** Function return declaration
LOGICAL, DIMENSION(2) :: newt_raph                                     ! ************************************************

newt_raph=.FALSE. !** Assume failiure to start                         FUNCTION df(x)
!**
DO !*** Start looping.                                                 !** Returns the value of the derivative of the function "f" at "x"
old_root=x                                                           !**
!** EXIT loop if either flag is true
x=x-f(x)/df(x)                           !* Calc. new root              !*** Dummy declaration
newt_raph=end_loop(x,old_root,tol,iters) !* Test if to exit loop        REAL, INTENT(IN) :: x
PRINT ’(I4,2F12.6)’,iters,x,f(x) !** Output results
IF (newt_raph(1) .OR. newt_raph(2)) EXIT                                !*** Function return declaration
END DO                                                                    REAL :: df

END FUNCTION newt_raph                                                      df=2*x

! ************************************************                       END FUNCTION df

FUNCTION bisection(a,b,tol)                                              ! ************************************************

!*** Dummy declarations                                                FUNCTION end_loop(num1,num2,tol,iters)
REAL, INTENT(INOUT) :: a,b                                             !*** Returns a logical type. If .TRUE. then the tolerance has been met or
REAL, INTENT(IN) :: tol                                                !*** the maximum number of iterations has been exceeded.

!*** Function return declaration                                         !*** Dummy variables
LOGICAL, DIMENSION(2) :: bisection                                       REAL, INTENT(IN) :: num1,num2
REAL, INTENT(IN) :: tol
INTEGER, INTENT(OUT) :: iters
!*** Local variable declarations
REAL :: p                                                                !*** Function return declaration
INTEGER, SAVE :: iters=0                                                 LOGICAL, DIMENSION(2) :: end_loop

bisection=.FALSE.                                                        end_loop(1)=(ABS(num1-num2)/2 < tol)
DO                                                                       end_loop(2)=(iters>1000)
bisection=end_loop(a,b,tol,iters)
IF (bisection(1) .OR. bisection(2)) EXIT                              iters=iters+1
p=(a+b)/2
PRINT ’(I4,4F12.6)’,iters,a,b,p,f(p) !** Output results             END FUNCTION end_loop

IF (f(a)*f(p) > 0) THEN     !** IF f(a) & f(p) same sign             ! ***********************************************************************
a=p
ELSE
b=p                                                              END MODULE main_mod
ENDIF
END DO                                                               *******************************************
*******************************************
END FUNCTION bisection

! ************************************************

FUNCTION f(x)
!**
!** Returns the value of the function "f" at "x"
!**

!*** Dummy declaration
REAL, INTENT(IN) :: x

```
To top