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

VIEWS: 16 PAGES: 2

									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