¨ç¼Æ¨D®Ú¡G¤û¹yªk¡]§t¼s¸q¤û¹yªk¡^
±q¤½¦¡¨Ó¤F¸Ñ¤û¹yªk
°²³] x1 ¤£¬O f(x) ªº®Ú¡]¦]¦¹ f(x1) ¤£µ¥©ó 0¡^¡A¦ý¤]°÷±µªñ¯u¥¿ªº®Ú x' ¤F¡C °²³]ÁÙ®t Δx (Δx = x' - x1)¡A«h¼g¤U f(x') ªº®õ°Ç®i¶}·|¬O
f(x1+Δx) = f(x1) + f'(x1)Δx + ...
©¿²¤°ª¦¸¶µ¡A¨Ã¥Î¤W f(x') = 0¡]®Úªº©w¸q¡^
0 = f(x1) + f'(x1)Δx
Δx = - f(x1) / f'(x1)
³oÓΔx ªºªí¥Ü¦¡¬O«Ü¦³¥Îªº¡A¥¦§i¶D§ÚÌ¡A¥Ø«e²q´úªº¦ì¸m x1 »P®Úªº¦ì¸m¤j¬ùÁÙ¦³ - f(x1) / f'(x1) ¨º»ò»·¡C
¦³¤F³oÓ¸ê°T¡A§ÚÌ·íµM´N¥i¥HÀò±o§ó±µªñ®Úªº²q´úÈ x1 +Δx¡A¦b¦¹¥i§â¥¦¥s°µ¬O x2¡C¤û¹yªk´N¦b³oÓ²³æ¦Ó¦³«Â¤Oªººtºâªk¤U¡A±q x1¡Bx2¡Bx3 ¤@ª½§ó·s¤U¥h¡A«Ü§Ö´N§ä¨ì¨ç¼Æ®Ú¤F¡C
¡]¦b¤£°Ñ¦Ò¥ô¦ó¸ê®Æªº±¡ªp¤U¡An¯à±À¾É¥X¤û¹yªkªº¤½¦¡¡^
±q¹Ï§Î¨Ó¤F¸Ñ¤û¹yªk
±q¹Ï¨Ó¬Ý¤ÏÂШDΔx = - f(x1) / f'(x1) ¨Ã¥B¨Ì§Ç§ä¥X x1¡Bx2¡Bx3 ¦b´X¦ó¤W§êºt¤°»ò¼Ëªº·N¸q¡C
§Ú̦P®É¤]n¤p¤ß¥H¤Uªºª¬ªp¡A
¤Ï¬M¦b¤½¦¡¤¤ªº¡A¥¿¬O¤À¥À¡]±×²v¡^¤£¯à¬°¹sªºn¨D¡C
¥t¥~¡AÁÙ¦³¦¬ÀĪº¹Lµ{¥d¦b¬Y¤@¹ïÂI¤W¡AµLªk¦A¶i¨B¡C
¨ç¦¡ rtnewt
FUNCTION rtnewt(funcd,x1,x2,xacc)
INTEGER JMAX
REAL rtnewt,x1,x2,xacc
EXTERNAL funcd
PARAMETER (JMAX=20) Set to maximum number of iterations.
Using the Newton-Raphson method, find the root of a function known to lie in the interval
[x1, x2]. The root rtnewt will be refined until its accuracy is known within ¡Óxacc. funcd
is a user-supplied subroutine that returns both the function value and the first derivative
of the function at the point x.
INTEGER j
REAL df,dx,f¨ä¤¤¡A¨Ï¥ÎªÌ¦Û¤vnµ¹ªº°Æµ{¦¡ funcd ¬O³o¼Ë³Q©I¥sªº
¤@Ó¯à¨Ï¥Î rtnewt ¨ç¦¡ªº¥Dµ{¦¡¡A¨ä¤¤¨Ã§tnÅý rtnewt ©I¥sªº¨Ï¥ÎªÌ´£¨Ñ¤§°Æµ{¦¡ funcd ¡]¨D¨ç¼Æ x3 ªº®Ú¡^¡A ½Ð¤j®a°Ñ¦Ò½m²ß¡G
program rtnewt_main
real rtnewt, x1, x2, xacc
external rtnewt, funcdwrite (*,*) 'Please type in interval x1 and x2 :'
read (*,*) x1, x2
write(*,*) 'Please typle in theaccuracy for root finding :'
read (*,*) xacc
write(*,*) 'The root is ', rtnewt(funcd,x1,x2,xacc)
write(*,*) 'within the accuracy of +/- ', xacc
end
subroutine funcd(x,f,df)
real x, f, df
f = x**3
df = 3*(x**2)
end
¨ç¦¡ rtsafe
¤û¹yªkÁöµM¾a±×²vªº¤Þ¾É¯à«Ü§Ö¦a§ä¨ì®Ú¡A¦ý¤]¦³«ez´XºØ¹³¬O¶]¥X½d³ò¡B¨«¤JµL½a°j°éµ¥ª¬ªp¡A¯Âºéªº¤û¹yªk¤£¤@©w¯à§ä¨ì¸Ñ¡Crtsafe µ²¦X¤F§Ö³t¤û¹yªk»Péwªº¤G¤Àªk¡A¦b¤û¹yªk¥¢®Äªº®ÉÔ§ï¥Î¤G¤Àªk¡A¥un¦³¥]³ò¨ì®Ú´N¤@©w·|§ä±o¨ì¡C
FUNCTION rtsafe(funcd,x1,x2,xacc)
INTEGER MAXIT
REAL rtsafe,x1,x2,xacc
EXTERNAL funcd
PARAMETER (MAXIT=100) Maximum allowed number of iterations.
Using a combination of Newton-Raphson and bisection, find the root of a function bracketed
between x1 and x2. The root, returned as the function value rtsafe, will be refined until
its accuracy is known within ¡Óxacc. funcd is a user-supplied subroutine which returns
both the function value and the first derivative of the function.
INTEGER j
REAL df,dx,dxold,f,fh,fl,temp,xh,xl¨Ï¥Î rtsafe »P¤W± rtnewt ¬O¤@¼Ò¤@¼Ëªº¡A½Ð¦Û¦æ½m²ß¡C
¦hºû«×ªº Newton-Raphson ¤èªk
¦hºû«×ªº®õ°Ç®i¶}
¡]¥H¤Uµe±¨ú¦Û½Ò¥»¡^
ª`·N¡]9.6.6¡^¬O¯x°}¤èµ{¦¡¡A»Ýn³z¹L¥ý«e¾Ç¹Lªº¸Ñ½u©Ê¥N¼Æ¤èµ{¦¡ªº¤èªk¨ÓÀò±o dx ¡C
°Æµ{¦¡ªº¨Ï¥Î
SUBROUTINE mnewt(ntrial,x,n,tolx,tolf)
INTEGER n,ntrial,NP
REAL tolf,tolx,x(n)
PARAMETER (NP=15) Up to NP variables.
C USES lubksb,ludcmp,usrfun
Given an initial guess x for a root in n dimensions, take ntrial Newton-Raphson steps to
improve the root. Stop if the root converges in either summed absolute variable increments
tolx or summed absolute function values tolf.
INTEGER i,k,indx(NP)
REAL d,errf,errx,fjac(NP,NP),fvec(NP),p(NP)¦p¦P«e±¦³»¡¨ì¡A³oùØ·|¥Î¨ì¸Ñ½u©Ê¥N¼Æ¤èµ{°ÝÃD¡C¦]¦¹·|¥Î¨ì lubksb ¤Î ludcmp ¡A¦ý³o¬O mnewt n¥s¥Îªº¡A§A¼gªº¥Dµ{¦¡¤£·|¥X²{¥s¥Î¥¦Ìªº°Ê§@¡A°O±o¦b½sĶ³sµ²®Én´£¨ì³o¨âӰƵ{¦¡«K¬O¡C
«nªº¬O¡A¬°¤FÅý mnewt ¥s¥Î¡A¨Ï¥ÎªÌ¥²¶·¦Û¦æ¼g¤@Ó usrfun ªº°Æµ{¦¡¡A¥¦¥]§t¤F¯àµ¹¥X F ¦V¶q¤Î J ¯x°}ªº¥\¯à¡C§A¦Û¤vn°Ê¤â¨D±o¨ç¼Æªº·L¤À«á¤½¦¡¡C
¥H¤W¬O usrfun ¦b mnewt ¤¤³Q¥s¥Îªº±¡¦æ¡A¥Ñ¨ä³Q©I¥sªº§Î¦¡¥i¨£¡A(1) usrfun À³«Å§i¬°°Æµ{¦¡¡A(2) ¨Ì·Ó x, n, NP, fvec, fjac ªº¶¶§Ç©w°Æµ{¦¡¤Þ¼Æªº¶¶§Ç¡]¤Þ¼Æ¦WºÙ³q±`·|»P©I¥s¥¦ªº¤W¼hµ{¦¡¥Îªº¤£¦P¡^¡A(3) «ö³o¨Ç¤Þ¼Æ¦b¤W¼h©I¥sµ{¦¡¤¤©Ò³Q«Å§iªº¼ÆÈ«¬ºA¡A°Æµ{¦¡¤¤¤]§@¤@¼Ëªº«Å§i¡C
½d¨Ò¥Dµ{¦¡¡]§t usrfun °Æµ{¦¡¡^¦p¤U¡G
program mnewt_main
implicit nonec integer n_max to define maximum n
integer n_max, i, i_again
parameter (n_max = 16)c integer for mnewt, n has to be smaller than NP=16
integer ntrial,n
parameter (n=2)
real x(n),tolx,tolfwrite(*,*) 'How many repeating steps for Newton methods ?'
read(*,*) ntrial
c write(*,*) 'How many variables is there in your problem ?'
c read(*,*) n
write(*,*) 'Please type in x tol and f tol for mnewt :'
read(*,*) tolx, tolf100 ¡@write(*,*) 'What is your inital guess on x(i) vector ?'
read(*,*) (x(i),i=1,n)
call mnewt(ntrial,x,n,tolx,tolf)
write(*,*) 'The root vector is :'
do i=1,n
write(*,*) x(i)
end do101 ¡@write(*,*) 'Find root again ? (1=Yes;0=No)'
read(*,*) i_again
if (i_again.eq.1) goto 100
if (i_again.eq.0) stop
goto 101end
subroutine usrfun (x,n,NP,fvec,fjac)
integer n,NP
real x(n),fvec(NP),fjac(NP,NP)fvec(1) = x(2) - x(1)
fvec(2) = x(2) + x(1)
fjac(1,1) = -1
fjac(2,1) = 1
fjac(1,2) = 1
fjac(2,2) = 1end
µù¡G¥H¤Wªº½d¨Ò¬O¨D F(x,y) ªº®Ú¡A¨ä¤¤ F1(x,y) = y - x¡BF2(x,y) = y + x¡C
Y Jacobian ¤£®e©ö°Ê¤â±À¾É±o¨ì¡A¤]¥i§ï¥Î¼ÆÈ·L¤Àªº fdjac ¦p¤U¡G
¦p¦ó¦b§A¦Û¤v¼gªº usrfun ¤¤¦A¥s¥Î fdjac¡A½Ð¬Ý²M·¡¨ä¨Ï¥Î¤èªk»¡©ú«á¥Î¤ß·Q¤@·Q¨Ã½T¹ê½m²ß¡C
¤û¹yªk»P¸H§Î (fractal)
«e±´¿¸gĵ§i¹L¤j®a¤û¹yªk¦³®É¤]·|¸I¨ì§ä¤£¨ì®Úªº±¡ªp¡A³o·íµM»P¨ç¼Æ§Îª¬ªº¯S®í©Ê¦³Ãö¡A¦Ó³oùئ³Ó¦³½ìªº¹ê§@¤j®a¥i¥H¸Õ¸Õ¬Ý¡A¥¦·|²£¥Í¸H§Î (http://en.wikipedia.org/wiki/Fractal) ªº¹Ï§Î¡A§Ų́ӬݥH¤U½ÆÅܨç¼Æªº®Ú
¨Ì·Ó¤û¹yªk¡A¬YÓ®Úªº²q´ú¦ì¸m¨ì¤UÓ®Úªº²q´ú¦ì¸m¶¡ªº¶ZÂ÷¡A¬O¦b¸ÓÂI¨útªº¨ç¼ÆȦA°£¥H¦b¨º¤@ÂI³Bªº±×²v¡C«ö¤W±¨ç¼Æªº©w¸q¡A§â Dz ¥[¨ì zj ¤W¨Ó±o¨ì zj+1 ¡A¥i¼g¤U¥H¤U¥i¥Î¨Ó¤@¨B¨B¡¥NªºÃö«Y¦¡
¦b½Æ¼Æ¥±¤Wªº¨C¤@ÓÂI¡A³£¥i®³¨Ó§@¬°²q´úªº°_©lÈ¡A¬Ý¥¦¬O§_·|¤@¨B¨B¹Gªñ¨ìº¡¨¬ì¨ç¼Æªº®Ú¡AY¬O¡A«h¸Ó°_©lȼХզâ¡AY¤£¬O«h¼Ð¶Â¦â¡C¦pªG§Ú̧â½Æ¼Æ¥±¤W©Ò¦³ªºÂI³£¼Ð©wÃC¦â¡A«h§ÚÌ·|¬Ý¨ì¸H§Î¡A³oºØ¹Ï¤£ºÞ¦p¦ó©ñ¤j¡A³£·|¦³«ÂÐÃþ¦üªºµ²ºc¤£Â_¥X²{¡C¤j®a¥i¥H¼g¤@Óµ{¦¡¡A¥Î pgplot µe¥X¨Ó¬Ý¬Ý¡C
§ó¶i¤@¨Bªº²Ó¸`¡A¥i°Ñ¨£§Ú¥t¥~½s¼gªº ¸H§Î ³æ¤¸¡C
ºô¸ô¤WÃö©ó¸H§Î¹Ï®w©Îµ{¦¡ªººô¯¸«Ü¦h¡A¦³¨Ç¬Æ¦Ü¿Ä¤J¤FÃÀ³N¤Æªº³B²z¡C
http://sprott.physics.wisc.edu/carlson/
http://www.fractalus.com/gumbycat/
²V¨PÃÀ´Y
http://www-chaos.umd.edu/gallery.htmlICM 2006 Benoit Mandelbrot ¸H§ÎÃÀ³NÄvÁÉ
http://www.fractalartcontests.com/2006/¥t¥~¡A¤@¥»«Ü¦nªº®Ñ¬O "The Beauty of Fractals"¡A¥Ñ H.-O. Petigen »P P.H. Richter ©ÒµÛ¡ASpringer-Verlag ¥Xª©¡C