{VERSION 2 3 "IBM INTEL NT" "2.3" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 }{CSTYLE "2D Math" -1 2 "Times" 0 1 0 0 0 0 0 0 2 0 0 0 0 0 0 }{CSTYLE "2D Comment" 2 18 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 } {CSTYLE "TXT CMD" -1 256 "MS Sans Serif" 0 0 128 0 0 1 0 0 0 0 0 0 0 0 0 }{CSTYLE "Bookmark" 20 257 "" 0 0 0 128 0 1 1 0 0 0 0 1 0 0 0 } {CSTYLE "word" -1 258 "" 0 0 128 0 0 1 0 0 0 0 0 0 0 0 0 }{CSTYLE "boo kmark" -1 259 "" 0 0 0 128 0 1 0 0 0 0 0 0 0 0 0 }{CSTYLE "" -1 260 " " 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 }{CSTYLE "" -1 261 "" 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 }{CSTYLE "" -1 262 "" 0 1 0 0 0 0 0 0 1 0 0 0 0 0 0 } {CSTYLE "" -1 263 "" 0 1 0 0 0 0 0 1 1 0 0 0 0 0 0 }{PSTYLE "Normal" -1 0 1 {CSTYLE "" -1 -1 "" 0 0 0 128 0 0 0 0 0 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Text Output" -1 2 1 {CSTYLE "" -1 -1 "Courier" 1 10 0 0 255 1 0 0 0 0 0 1 3 0 0 }1 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Heading 1" 0 3 1 {CSTYLE "" -1 -1 "" 1 18 255 255 0 0 0 1 0 0 0 0 0 0 0 }1 0 0 0 6 6 0 0 0 0 0 0 -1 0 }{PSTYLE "Heading 2 " 3 4 1 {CSTYLE "" -1 -1 "" 1 14 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 0 0 -1 4 4 0 0 0 0 0 0 -1 0 }{PSTYLE "" 2 6 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 2 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "newpage " -1 256 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 1 0 -1 0 }{PSTYLE "vfill" -1 257 1 {CSTYLE "" -1 -1 " " 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 } {PSTYLE "Definition" -1 258 1 {CSTYLE "" -1 -1 "" 0 0 0 64 128 1 0 0 0 0 0 0 0 0 0 }0 0 0 -1 4 4 0 0 0 0 0 0 -1 0 }{PSTYLE "Theorem" -1 259 1 {CSTYLE "" -1 -1 "" 0 0 219 36 36 1 0 0 0 0 0 0 0 0 0 }0 0 0 -1 4 4 0 0 0 0 0 0 -1 0 }{PSTYLE "Problem" -1 260 1 {CSTYLE "" -1 -1 "" 0 0 0 0 255 1 0 0 0 0 0 0 0 0 0 }0 0 0 -1 4 4 3 4 0 0 0 0 -1 0 } {PSTYLE "dblnorm" -1 261 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 2 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "code" 2 262 1 {CSTYLE "" -1 -1 "Comic Sans MS" 0 0 128 0 128 1 0 1 0 0 0 0 3 0 0 }0 0 0 -1 -1 -1 3 12 0 0 0 0 -1 0 }{PSTYLE "asis" 0 263 1 {CSTYLE "" -1 -1 "Arial Narrow" 1 12 128 64 0 1 0 0 0 0 0 1 3 0 0 }0 0 0 -1 -1 -1 3 6 0 0 0 0 0 0 }{PSTYLE "subproblem" 0 264 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "dia gram" -1 265 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "dblnorm.mws" -1 266 1 {CSTYLE " " -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }0 2 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "Item" 0 267 1 {CSTYLE "" -1 -1 "Lucida Sans" 1 12 0 0 0 1 0 0 0 0 0 0 0 0 0 }0 0 0 -1 6 -1 3 0 0 0 0 0 -1 0 }{PSTYLE "" 4 268 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 269 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 270 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 271 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 4 272 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 273 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 274 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 275 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }{PSTYLE "" 0 276 1 {CSTYLE "" -1 -1 "" 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 }3 0 0 -1 -1 -1 0 0 0 0 0 0 -1 0 }} {SECT 0 {EXCHG {PARA 268 "" 0 "" {TEXT 262 53 "Use of Taylor's theorem in numerical differentiation." }}{PARA 0 "" 0 "" {TEXT -1 45 "The der ivative of f at x = a is defined as " }{XPPEDIT 18 0 "f*`'`(a)=limit ((f(a+h)-f(a))/h,h=0" "/*&%\"fG\"\"\"-%\"'G6#%\"aGF%-%&limitG6$*&,&-F$ 6#,&F)F%%\"hGF%F%-F$6#F)!\"\"F%F2F5/F2\"\"!" }{TEXT -1 262 ". Numeri cally, however, when h becomes so small that f(a+h) - f(a) = 0 in th e floating point arithmetic we are using, the ratio drops to 0 and so \+ we must seek other estimates to the derivative. One way is to use Ta ylor's theorem to expand f(a+h) and f(a-h)" }}{PARA 269 "" 0 "" {XPPEDIT 18 0 "f(a+h) = f(a) + f*`'`(a)*h +f*`''`(a)/2*h^2 + `...`" "/ -%\"fG6#,&%\"aG\"\"\"%\"hGF(,*-F$6#F'F(*(F$F(-%\"'G6#F'F(F)F(F(**F$F(- %#''G6#F'F(\"\"#!\"\"F)\"\"#F(%$...GF(" }{TEXT -1 1 " " }}{PARA 270 " " 0 "" {XPPEDIT 18 0 "f(a-h) = f(a) - f*`'`(a)*h +f*`''`(a)/2*h^2 - `. ..` " "/-%\"fG6#,&%\"aG\"\"\"%\"hG!\"\",*-F$6#F'F(*(F$F(-%\"'G6#F'F(F) F(F***F$F(-%#''G6#F'F(\"\"#F*F)\"\"#F(%$...GF*" }}{PARA 0 "" 0 "" {TEXT -1 30 "Subtracting, and solving for " }{XPPEDIT 18 0 "f*`'`(a) " "*&%\"fG\"\"\"-%\"'G6#%\"aGF$" }{TEXT -1 8 ", we get" }}{PARA 271 " " 0 "" {TEXT -1 6 "(1) " }{XPPEDIT 18 0 "f*`'`(a)= (f(a+h)-f(a-h))/( 2*h) - f*`'''`(a)/3!*h^2 - `...`" "/*&%\"fG\"\"\"-%\"'G6#%\"aGF%,(*&,& -F$6#,&F)F%%\"hGF%F%-F$6#,&F)F%F0!\"\"F4F%*&\"\"#F%F0F%F4F%**F$F%-%$'' 'G6#F)F%-%*factorialG6#\"\"$F4F0\"\"#F4%$...GF4" }}{PARA 0 "" 0 "" {TEXT -1 64 "What this indicates (but doesn't quite prove) is that the ratio " }{XPPEDIT 18 0 "(f(a+h)-f(a-h))/(2*h)" "*&,&-%\"fG6#,&%\"aG\" \"\"%\"hGF)F)-F%6#,&F(F)F*!\"\"F.F)*&\"\"#F)F*F)F." }{TEXT -1 32 " con verges to f'(a) with Order " }{XPPEDIT 18 0 "h^2" "*$%\"hG\"\"#" } {TEXT -1 4 ". " }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 260 "" 0 "" {TEXT 260 8 "Exercise" }{TEXT -1 117 ": We can compare the convergence of this approximation to f'(a) with the convergence of the standard \+ approximation " }{XPPEDIT 18 0 "(f(a+h)-f(a))/h" "*&,&-%\"fG6#,&%\"aG \"\"\"%\"hGF)F)-F%6#F(!\"\"F)F*F-" }{TEXT -1 489 " . Modify the progr am ndiv below (from a previous worksheet) to include this ratio. T est it on some functions at some points and make a (preliminary, of co urse) judgement about which approximation is better. Inputs to ndiv: f is the function, a is the point where we want the derivative f'(a) , inc is the initial value for h (it can be negative), and n is the nu mber of terms in the sequence of approximations to the derivative f'(a ), each time cutting the the value for h in half." }}}{EXCHG {PARA 0 " > " 0 "" {MPLTEXT 1 0 289 "ndiv := proc(f,a,inc,n)\n local h,i,aprx,df ;\n df := D(f);\n h := evalf(inc);\nprintf(`%12.7s %12.7s %12.7s\\n`, `arg `,\n`aprox `,`error `);\n for i from 1 to n do\naprx := e valf((f(a+h)-f(a))/h);\nprintf(`%12.7e %12.7e %12.7e\\n`, h, aprx,eva lf(df(a))-aprx);\n h := .5*h;\n od; \n NULL:\nend:" }}}{EXCHG {PARA 0 "" 0 "" {TEXT 263 10 "Solution: " }{TEXT -1 74 " First make a copy of ndiv by coloring and pasting it below the original." }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 239 "We can a line computing the symmetric approxim ation to f'(a) aprx2 and a column for the error in aprx2. While we' re at it, we can add a variable dfa to store the value of f'(a) rather than computing it over and over like we did in ndiv." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 405 "ndiv := proc(f,a,inc,n)\n local h, i,aprx,aprx2,dfa;\n dfa := evalf( D(f)(a));\n h := evalf(inc);\nprintf (`%8.7s %12.7s %12.7s %12.7s %12.7s\\n`, `arg `,\n` aprx `,` error1`,` aprx2` ,`error2`);\n for i from 1 to n do\naprx := evalf((f( a+h)-f(a))/h);\naprx2 := evalf((f(a+h)-f(a-h))/(2*h));\nprintf(`%12.7e % 12.7e % 12.7e % 12.7e % 12.7e\\n`, h, aprx,dfa-aprx,aprx2,dfa-aprx2 );\n h := .5*h;\n od; \n NULL:\nend:" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 117 "Now to test it. Try the sin function at a = .2 with an initia l value for h = inc of .1 and 20 terms in the sequence." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 19 "ndiv(sin,.2,.1,20);" }}{PARA 6 "" 1 "" {TEXT -1 1548 " arg aprx error1 ap rx2 error2\n1.0000000e-01 9.6850875e-01 1.1557818e-02 9.784 3395e-01 1.6326273e-03\n5.0000000e-02 9.7469257e-01 5.3740078e-03 \+ 9.7965826e-01 4.0830980e-04\n2.5000000e-02 9.7748125e-01 2.5853258e -03 9.7996449e-01 1.0208780e-04\n1.2500000e-02 9.7879938e-01 1.267 1938e-03 9.8004105e-01 2.5525800e-05\n6.2500000e-03 9.7943936e-01 \+ 6.2721780e-04 9.8006020e-01 6.3778000e-06\n3.1250000e-03 9.7975456e -01 3.1201780e-04 9.8006497e-01 1.6018000e-06\n1.5625000e-03 9.799 1097e-01 1.5560180e-04 9.8006617e-01 4.0180000e-07\n7.8125000e-04 \+ 9.7998886e-01 7.7713800e-05 9.8006649e-01 8.1800000e-08\n3.9062500e -04 9.8002764e-01 3.8929800e-05 9.8006656e-01 1.7800000e-08\n1.953 1250e-04 9.8004736e-01 1.9217800e-05 9.8006656e-01 1.7800000e-08\n 9.7656250e-05 9.8005708e-01 9.4898000e-06 9.8006630e-01 2.7380000e -07\n4.8828125e-05 9.8006016e-01 6.4178000e-06 9.8006528e-01 1.297 8000e-06\n2.4414062e-05 9.8006630e-01 2.7380000e-07 9.8006835e-01 - 1.7742000e-06\n1.2207031e-05 9.8006630e-01 2.7380000e-07 9.8006630e -01 2.7380000e-07\n6.1035156e-06 9.8005811e-01 8.4658000e-06 9.800 5811e-01 8.4658000e-06\n3.0517578e-06 9.8009087e-01 -2.4302000e-05 \+ 9.8009088e-01 -2.4302200e-05\n1.5258789e-06 9.8009087e-01 -2.4301700e -05 9.8009087e-01 -2.4301700e-05\n7.6293945e-07 9.8002534e-01 4.123 4300e-05 9.8002534e-01 4.1234300e-05\n3.8146972e-07 9.8015641e-01 - 8.9837600e-05 9.8015641e-01 -8.9837700e-05\n1.9073486e-07 9.7989427e -01 1.7230640e-04 9.7989427e-01 1.7230630e-04" }}}{EXCHG {PARA 0 " " 0 "" {TEXT -1 232 "Upon inspection, we see that the error in the get s smaller faster than the error in approximation 2. Both approximatio ns become unreliable when h gets near 10*machine epsilon (which is .5/ 10^9 for the defualt value of Digits = 10)." }}}{EXCHG {PARA 272 "" 0 "" {TEXT 261 31 "A better approximation to f'(a)" }}{PARA 0 "" 0 "" {TEXT -1 142 "The idea used in getting the approximation (1) to f'(a) above can be used to get a second approximation to f'(a) which conve rges with Order " }{XPPEDIT 18 0 "h^4" "*$%\"hG\"\"%" }{TEXT -1 29 ". \+ Name the approximation " }}{PARA 273 "" 0 "" {TEXT -1 4 " " } {XPPEDIT 18 0 "g[1] (h)= (f(a+h)-f(a-h))/(2*h)" "/-&%\"gG6#\"\"\"6#%\" hG*&,&-%\"fG6#,&%\"aG\"\"\"F)F1F1-F-6#,&F0F1F)!\"\"F5F1*&\"\"#F1F)F1F5 " }{TEXT -1 1 "." }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 37 "From equation (1) above, we can write" }}{PARA 274 "" 0 "" {TEXT -1 3 " " } {XPPEDIT 18 0 "g[1](h) = f*`'`(a) + a[1]*h^2 + a[2]*h^4 + a[3]*h^6+ ` ...`" "/-&%\"gG6#\"\"\"6#%\"hG,,*&%\"fG\"\"\"-%\"'G6#%\"aGF-F-*&&F16# \"\"\"F-*$F)\"\"#F-F-*&&F16#\"\"#F-*$F)\"\"%F-F-*&&F16#\"\"$F-*$F)\"\" 'F-F-%$...GF-" }}{PARA 0 "" 0 "" {TEXT -1 30 "The idea is to get rid o f the " }{XPPEDIT 18 0 "h^2" "*$%\"hG\"\"#" }{TEXT -1 33 " term. So e valuate g[1] at h/2." }}{PARA 275 "" 0 "" {TEXT -1 3 " " }{XPPEDIT 18 0 "g[1](h/2) = f*`'`(a) + a[1]*(h/2)^2 + a[2]*(h/2)^4 + a[3]*(h/2)^ 6+ `...`" "/-&%\"gG6#\"\"\"6#*&%\"hG\"\"\"\"\"#!\"\",,*&%\"fGF+-%\"'G6 #%\"aGF+F+*&&F46#\"\"\"F+*$*&F*F+\"\"#F-\"\"#F+F+*&&F46#\"\"#F+*$*&F*F +\"\"#F-\"\"%F+F+*&&F46#\"\"$F+*$*&F*F+\"\"#F-\"\"'F+F+%$...GF+" }} {PARA 0 "" 0 "" {TEXT -1 90 "Multiply this equation by 4, subtract the equation above and divide both sides by 3 to get" }}{PARA 276 "" 0 " " {TEXT -1 3 " " }{XPPEDIT 18 0 "(4*g[1](h/2) - g[1](h))/3= f*`'`(a) + b[1]*h^4 + b[2]*h^6 + `...`" "/*&,&*&\"\"%\"\"\"-&%\"gG6#\"\"\"6#* &%\"hGF'\"\"#!\"\"F'F'-&F*6#\"\"\"6#F/F1F'\"\"$F1,**&%\"fGF'-%\"'G6#% \"aGF'F'*&&%\"bG6#\"\"\"F'*$F/\"\"%F'F'*&&FA6#\"\"#F'*$F/\"\"'F'F'%$.. .GF'" }}{PARA 0 "" 0 "" {TEXT -1 59 "This approximation appears to con verge to f'(a) with Order " }{XPPEDIT 18 0 "h^4" "*$%\"hG\"\"%" } {TEXT -1 1 "." }}{PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 260 "" 0 "" {TEXT -1 97 "Exercise: Modify your modification of ndiff and test thi s approximation to f'(a). Is it better?" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "" 0 "" {TEXT -1 0 "" }}{PARA 260 "" 0 "" {TEXT -1 39 "Exercise: Call this ne w approximation " }{XPPEDIT 18 0 "g[2](h)" "-&%\"gG6#\"\"#6#%\"hG" } {TEXT -1 53 ". Use the same idea to get yet another approximation " } {XPPEDIT 18 0 "g[3](h)" "-&%\"gG6#\"\"$6#%\"hG" }{TEXT -1 37 " which c onverges to f'(a) with Order " }{XPPEDIT 18 0 "h^6" "*$%\"hG\"\"'" } {TEXT -1 1 "." }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}} {MARK "7 0 0" 232 }{VIEWOPTS 1 1 0 3 4 1802 }