38 PURE SUBROUTINE func(x,fx,plist)
41 REAL,
INTENT(IN),
DIMENSION(:),
OPTIONAL :: plist
42 REAL,
INTENT(OUT) :: fx
46 PURE SUBROUTINE funcd(x,fx,dfx,plist)
49 REAL,
INTENT(IN),
DIMENSION(:),
OPTIONAL :: plist
50 REAL,
INTENT(OUT) :: fx,dfx
55 INTEGER,
PARAMETER :: num_tests = 13
56 INTEGER,
PARAMETER :: num_methods = 8
57 CHARACTER(LEN=16),
PARAMETER,
DIMENSION(NUM_METHODS) :: method_name = (/ &
66 REAL,
PARAMETER,
DIMENSION(2,NUM_TESTS) :: bounds = reshape((/ 0.0, 1.2, &
80 DOUBLE PRECISION,
PARAMETER,
DIMENSION(NUM_TESTS) :: ref_roots = &
81 (/ 3.9942229171096819451d-01, 8.0413309750366432374d-01, 9.0340766319186021294d-01, &
82 7.7014241346192677110d-02, 2.5920449372984746773d-01, 5.3674166257799978186d-01, &
83 4.4754176206055907112d-01, 1.1111111111111111111d-01, 5.0000003403025908310d-01, &
84 6.7980892150470050192d-01, -3.5938136638046273022d-01, 1.6487212707001281468d-00, &
85 1.0000000000000000000d-00 &
88 REAL :: root, xm, dx_rel, dx_acc, plist(1)
89 INTEGER :: i,k, iter, error
90 CHARACTER(LEN=64) :: tap_message
91 LOGICAL :: verbose_results = .false.
94 WRITE (*,
"(A)")
"====================================================================="
95 WRITE (*,
"(A)")
" << Testing root finding algorithms >> "
96 WRITE (*,
"(A)")
"====================================================================="
99 tap_plan((num_methods-1)*num_tests)
108 xm = 0.5*(bounds(1,k)+bounds(2,k))
109 print
'(A,I2,A,2(ES10.2))',
"Test #",k,
" Search Interval: ",bounds(1,k),bounds(2,k)
110 print
'(A)',
"---------------------------------------------------------------------"
111 IF (verbose_results) print
'(A20,ES27.19)',
"reference result", ref_roots(k)
134 dx_rel = abs(1.0-root/ref_roots(k))
135 IF (verbose_results)
CALL printresults(method_name(i),root,ref_roots(k),error,iter)
144 WRITE (tap_message,
'(A16,A11,ES9.2,A3,ES9.2,A7,I4)') method_name(i),
": dx_rel = ", &
145 dx_rel,
" < ", dx_acc,
" iter = ", iter
146 tap_check(dx_rel.LE.dx_acc,tap_message)
150 print
'(A)',
"---------------------------------------------------------------------"
160 CHARACTER(LEN=*),
INTENT(IN) :: method
161 REAL,
INTENT(IN) :: root
162 DOUBLE PRECISION,
INTENT(IN) :: ref_root
163 INTEGER,
INTENT(IN) :: error,iter
167 print
'(A20,ES23.15,ES9.1,I4)',trim(method),root,abs(1.0-root/ref_root),iter
191PURE SUBROUTINE func(x,fx,plist)
194 REAL,
INTENT(IN) :: x
195 REAL,
INTENT(IN),
DIMENSION(:),
OPTIONAL :: plist
196 REAL,
INTENT(OUT) :: fx
200 fx = x*x*(x*x/3. + sqrt(2.)*sin(x)) - sqrt(3.)/18.
206 fx = 2*(x*exp(-9.)-exp(-9*x)) + 1.
210 fx = (x-1.)*exp(-9*x) + x**9
212 fx = x*x + sin(x/9.) - 0.25
216 fx = tan(x) - x - 0.0463025
218 fx = x*(x+sin(sqrt(75.)*x)) - 0.2
222 fx = log(x) + 0.5*x*x/exp(1.) -2*x/sqrt(exp(1.)) + 1.
229PURE SUBROUTINE funcd(x,fx,dfx,plist)
232 REAL,
INTENT(IN) :: x
233 REAL,
INTENT(IN),
DIMENSION(:),
OPTIONAL :: plist
234 REAL,
INTENT(OUT) :: fx,dfx
238 fx = x*x*(x*x/3. + sqrt(2.)*sin(x)) - sqrt(3.)/18.
239 dfx= 2*x*(sqrt(2.)*sin(x)+x*x/3.)+x*x*(sqrt(2.)*cos(x)+(2*x)/3.)
245 dfx= 35*35*x**34 - 1.0
247 fx = 2*(x*exp(-9.)-exp(-9*x)) + 1.
248 dfx= 2*(9*exp(-9*x)+exp(-9.))
251 dfx= 2*x + 9*(1.-x)**8
253 fx = (x-1.)*exp(-9*x) + x**9
254 dfx= (10.-9*x)*exp(-9*x) + 9*x**8
256 fx = x*x + sin(x/9.) - 0.25
257 dfx= 2*x + cos(x/9.)/9.
262 fx = tan(x) - x - 0.0463025
263 dfx= 1./cos(x)**2 - 1.
265 fx = x*(x+sin(sqrt(75.)*x)) - 0.2
266 dfx= sin(sqrt(75.)*x) + x*(2. + sqrt(75.)*cos(sqrt(75.)*x))
271 fx = log(x) + 0.5*x*x/exp(1.) - 2*x/sqrt(exp(1.)) + 1.
272 dfx= 1./x + x/exp(1.) - 2./sqrt(exp(1.))
pure subroutine funcd(y, fy, dfy, plist)
subroutine, public getroot_brentdekker(func, x1, x2, root, error, plist, xm, iterations)
subroutine, public getroot_bisection(func, x1, x2, root, error, plist, xm, iterations)
subroutine, public getroot_regulafalsi(func, x1, x2, root, error, plist, xm, iterations)
subroutine, public getroot_king(func, x1, x2, root, error, plist, xm, iterations)
subroutine, public getroot_ridder(func, x1, x2, root, error, plist, xm, iterations)
subroutine, public getroot_andersonbjoerk(func, x1, x2, root, error, plist, xm, iterations)
subroutine, public getroot_newton(funcd, x1, x2, root, error, plist, xm, iterations)
subroutine, public getroot_pegasus(func, x1, x2, root, error, plist, xm, iterations)
character(len=64) function, public geterrormessage(error)
real, parameter, public default_accuracy
pure subroutine func(x, fx, plist)
subroutine printresults(method, root, ref_root, error, iter)