Jeremy Knight Math 639 4/5/11 Algorithm 11.2: Nonlinear Shooting with Newton's Method restart;
<Text-field style="Heading 1" size="10" layout="Heading 1"><Font size="10">NonlinearShooting(f,fy,fyp,a,b,alpha,beta,N,TOL,M)</Font></Text-field> To approximate the solution of the nolinear boundary-value problem LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYoLUkjbWlHRiQ2JlEjeSJGJy8lJXNpemVHUSMxMEYnLyUnaXRhbGljR1EldHJ1ZUYnLyUsbWF0aHZhcmlhbnRHUSdpdGFsaWNGJy1JI21vR0YkNi5RIj1GJ0YvL0Y2USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGQC8lKXN0cmV0Y2h5R0ZALyUqc3ltbWV0cmljR0ZALyUobGFyZ2VvcEdGQC8lLm1vdmFibGVsaW1pdHNHRkAvJSdhY2NlbnRHRkAvJSdsc3BhY2VHUSwwLjI3Nzc3NzhlbUYnLyUncnNwYWNlR0ZPLUYsNiZRImZGJ0YvRjJGNS1JKG1mZW5jZWRHRiQ2JS1GIzYqLUYsNiZRInhGJ0YvRjJGNS1GOTYuUSIsRidGL0Y8Rj4vRkJGNEZDRkVGR0ZJRksvRk5RJjAuMGVtRicvRlFRLDAuMzMzMzMzM2VtRictRiw2JlEieUYnRi9GMkY1RmduRl9vLUY5Ni5RIidGJ0YvRjxGPkZBRkNGRUZHRklGSy9GTlEsMC4xMTExMTExZW1GJy9GUUZcb0YvRjxGL0Y8Ri9GPA==, for LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYpLUkjbWlHRiQ2JlEiYUYnLyUlc2l6ZUdRIzEwRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkjbW9HRiQ2LlEmJmxlcTtGJ0YvL0Y2USdub3JtYWxGJy8lJmZlbmNlR1EmZmFsc2VGJy8lKnNlcGFyYXRvckdGQC8lKXN0cmV0Y2h5R0ZALyUqc3ltbWV0cmljR0ZALyUobGFyZ2VvcEdGQC8lLm1vdmFibGVsaW1pdHNHRkAvJSdhY2NlbnRHRkAvJSdsc3BhY2VHUSwwLjI3Nzc3NzhlbUYnLyUncnNwYWNlR0ZPLUYsNiZRInhGJ0YvRjJGNUY4LUYsNiZRImJGJ0YvRjJGNUYvRjw=, with LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYoLUkjbWlHRiQ2JlEieUYnLyUlc2l6ZUdRIzEwRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkobWZlbmNlZEdGJDYlLUYjNiUtRiw2JlEiYUYnRi9GMkY1Ri8vRjZRJ25vcm1hbEYnRi9GQC1JI21vR0YkNi5RIj1GJ0YvRkAvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkgvJSlzdHJldGNoeUdGSC8lKnN5bW1ldHJpY0dGSC8lKGxhcmdlb3BHRkgvJS5tb3ZhYmxlbGltaXRzR0ZILyUnYWNjZW50R0ZILyUnbHNwYWNlR1EsMC4yNzc3Nzc4ZW1GJy8lJ3JzcGFjZUdGVy1GLDYmUScmIzk0NTtGJ0YvL0YzRkhGQEYvRkA= and LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYoLUkjbWlHRiQ2JlEieUYnLyUlc2l6ZUdRIzEwRicvJSdpdGFsaWNHUSV0cnVlRicvJSxtYXRodmFyaWFudEdRJ2l0YWxpY0YnLUkobWZlbmNlZEdGJDYlLUYjNiUtRiw2JlEiYkYnRi9GMkY1Ri8vRjZRJ25vcm1hbEYnRi9GQC1JI21vR0YkNi5RIj1GJ0YvRkAvJSZmZW5jZUdRJmZhbHNlRicvJSpzZXBhcmF0b3JHRkgvJSlzdHJldGNoeUdGSC8lKnN5bW1ldHJpY0dGSC8lKGxhcmdlb3BHRkgvJS5tb3ZhYmxlbGltaXRzR0ZILyUnYWNjZW50R0ZILyUnbHNwYWNlR1EsMC4yNzc3Nzc4ZW1GJy8lJ3JzcGFjZUdGVy1GLDYmUScmIzk0NjtGJ0YvL0YzRkhGQEYvRkA= NonlinearShooting:=proc(f,fy,fyp,a,b,alpha,beta,N,TOL,M) #f = function f(x,y(x,t),y'(x,t)) #fy = function (df/dy)(x,y,y') #fyp = function (df/dy')(x,y,y') #a,b=endpoints; alpha,beta=boundary conditions #N>=2 is number of subintervals #TOL = tolerance #M = Maximum iterations local c,i,h,k,kp,w,u,x,TK,temp: #option trace; h:=(b-a)/N: c:=1: #counter TK:=(beta-alpha)/(b-a): w:=Array(1..2,1..N+1): #approximations array k:=Array(1..4,1..2): #Array for Runge-Kutta method kp:=Array(1..4,1..2): #k' Array for Runge-Kutta method u:=Array(1..2): while c<=M do w[1,1]:=alpha: w[2,1]:=TK: u[1]:=0: u[2]:=1: #Runge-Kutta Method loop for i from 2 to N+1 do x:=a+(i-2)*h: k[1,1]:=h*(w[2,i-1]): k[1,2]:=h*(f(Array([ x , w[1,i-1] , w[2,i-1] ]))): k[2,1]:=h*(w[2,i-1]+0.5*k[1,2]): k[2,2]:=h*(f(Array([ x+h/2.0 , w[1,i-1]+0.5*k[1,1] , w[2,i-1]+0.5*k[1,2]]))): k[3,1]:=h*(w[2,i-1]+0.5*k[2,2]): k[3,2]:=h*(f(Array([ x+h/2.0 , w[1,i-1]+0.5*k[2,1] , w[2,i-1]+0.5*k[2,2]]))): k[4,1]:=h*(w[2,i-1]+k[3,2]): k[4,2]:=h*(f(Array([x+h , w[1,i-1]+k[3,1] , w[2,i-1]+k[3,2] ]))): w[1,i]:=w[1,i-1]+(k[1,1]+2*k[2,1]+2*k[3,1]+k[4,1])/6.0: w[2,i]:=w[2,i-1]+(k[1,2]+2*k[2,2]+2*k[3,2]+k[4,2])/6.0: kp[1,1]:=h*(u[2]): temp:=Array([x , w[1,i-1] , w[2,i-1]]): kp[1,2]:=h*(fy(temp)*u[1] + fyp(temp)*u[2] ): kp[2,1]:=h*(u[2]+0.5*kp[1,2]): temp:=Array([ x+h/2.0 , w[1,i-1] , w[2,i-1] ]): kp[2,2]:=h*(fy(temp)*(u[1]+0.5*kp[1,1]) + fyp(temp)*(u[2]+0.5*kp[1,2]) ): kp[3,1]:=h*(u[2]+0.5*kp[2,2]): temp:=Array([ x+h/2.0 , w[1,i-1] , w[2,i-1] ]): kp[3,2]:=h*(fy(temp)*(u[1]+0.5*kp[2,1])+fyp(temp)*(u[2]+0.5*kp[2,2])): kp[4,1]:=h*(u[2]+kp[3,2]): temp:=Array([ x+h , w[1,i-1] , w[2,i-1] ]): kp[4,2]:=h*(fy(temp)*(u[1]+kp[3,1]) + fyp(temp)*(u[2]+kp[3,2])): u[1]:=u[1]+(kp[1,1]+2*kp[2,1]+2*kp[3,1]+kp[4,1])/6.0: u[2]:=u[2]+(kp[1,2]+2*kp[2,2]+2*kp[3,2]+kp[4,2])/6.0: end do: print (c, "TK=",TK); #check tolerance and print if reached. if abs(w[1,N+1]-beta)<=TOL then print("i, x, w[1,i], w[2,i]"): for i from 1 to N+1 do x:=a+(i-1)*h: print (i-1 ,x, w[1,i] , w[2,i]): end do: print("iterations=",c): break: end if: #Newton's Method is used to compute TK TK:=TK-(w[1,N+1]-beta)/u[1]: c:=c+1: end do: if c>M then print("Maximum number of iterations exceeded"): end if: end proc:
<Text-field style="Heading 1" size="10" layout="Heading 1"><Font size="10">Example pg 682</Font></Text-field> #define array w such that x=w[1] , y=w[2], y'=w[3] f:=(w::Array)->1.0/8.0*(32.0+2.0*w[1]^3-w[2]*w[3]); fy:=(w::Array)->-1.0/8.0*w[3]; fyp:=(w::Array)->-1.0/8.0*(w[2]); Zio2IydJIndHNiJJJkFycmF5RyUqcHJvdGVjdGVkR0YmNiRJKW9wZXJhdG9yR0YmSSZhcnJvd0dGJkYmKigkIiM1ISIiIiIiJCIjISlGL0YvLCgkIiQ/JEYvRjAqJiQiIz9GL0YwKiQmOSQ2I0YwIiIkRjBGMComJkY7NiMiIiNGMCZGOzYjRj1GMEYvRjBGJkYmRiY= Zio2IydJIndHNiJJJkFycmF5RyUqcHJvdGVjdGVkR0YmNiRJKW9wZXJhdG9yR0YmSSZhcnJvd0dGJkYmLCQqKCQiIzUhIiIiIiIkIiMhKUYwRjAmOSQ2IyIiJEYxRjBGJkYmRiY= Zio2IydJIndHNiJJJkFycmF5RyUqcHJvdGVjdGVkR0YmNiRJKW9wZXJhdG9yR0YmSSZhcnJvd0dGJkYmLCQqKCQiIzUhIiIiIiIkIiMhKUYwRjAmOSQ2IyIiI0YxRjBGJkYmRiY= NonlinearShooting(f,fy,fyp,1.0,3.0,17.0,43.0/3.0,20,10^(-5),10); NiUiIiJRJFRLPTYiJCErTkxMTDghIio= NiUiIiNRJFRLPTYiJCErPU5lPzshIik= NiUiIiRRJFRLPTYiJCErUE1mNjkhIik= NiUiIiVRJFRLPTYiJCErKTNaLFMiISIp NiUiIiZRJFRLPTYiJCErSC4tKzkhIik= UTVpLH54LH53WzEsaV0sfndbMixpXTYi NiYiIiEkIiM1ISIiJCIkcSJGJiQhK0guLSs5ISIp NiYiIiIkIisrKysrNiEiKiQiKzgmXGJkIiEiKSQhKyV5TUI1IkYp NiYiIiMkIisrKysrNyEiKiQiK0kqUXRaIiEiKSQhK0I8STYoKUYm NiYiIiQkIisrKysrOCEiKiQiK3FeeCpSIiEiKSQhK1lVaW5vRiY= NiYiIiUkIisrKysrOSEiKiQiK2NHJylROCEiKSQhKyJIN01PJkYm NiYiIiYkIisrKysrOiEiKiQiKyQqPW4iSCIhIikkIStGJ1E3NiVGJg== NiYiIickIisrKysrOyEiKiQiK01ZK2M3ISIpJCErKzE2XUlGJg== NiYiIigkIisrKysrPCEiKiQiKyhbIT1JNyEiKSQhK20jR2s4I0Ym NiYiIikkIisrKysrPSEiKiQiKy1CKkdAIiEiKSQhK2tfTlE4RiY= NiYiIiokIisrKysrPiEiKiQiKzIiM0o/IiEiKSQhKyNbZj9LJyEjNQ== NiYiIzUkIisrKysrPyEiKiQiKzxCKys3ISIpJCEnJWVQJyEjNQ== NiYiIzYkIisrKysrQCEiKiQiKypmMUg/IiEiKSQiKzRPRT1kISM1 NiYiIzckIisrKysrQSEiKiQiK0pURjY3ISIpJCIrLGk7JTQiRiY= NiYiIzgkIisrKysrQiEiKiQiKyo9YFlBIiEiKSQiKyIqR1F2OkYm NiYiIzkkIisrKysrQyEiKiQiK050bVU3ISIpJCIrTmA9QT9GJg== NiYiIzokIisrKysrRCEiKiQiK2kuK2w3ISIpJCIrKjNvKlJDRiY= NiYiIzskIisrKysrRSEiKiQiKzNaUSJIIiEiKSQiK2wnM0okR0Ym NiYiIzwkIisrKysrRiEiKiQiK11DZkA4ISIpJCIrVCMqPTBLRiY= NiYiIz0kIisrKysrRyEiKiQiK0QjR2FOIiEiKSQiK21SO2ZORiY= NiYiIz4kIisrKysrSCEiKiQiKzxPcyNSIiEiKSQiKyhmJ1soKlFGJg== NiYiIz8kIisrKysrSSEiKiQiK2xFTEw5ISIpJCIrcikzQUElRiY= NiRRLGl0ZXJhdGlvbnM9NiIiIiY=
LUklbXJvd0c2Iy9JK21vZHVsZW5hbWVHNiJJLFR5cGVzZXR0aW5nR0koX3N5c2xpYkdGJzYjLUkjbWlHRiQ2I1EhRic= JSFH
<Text-field style="Heading 1" layout="Heading 1">Exercises 11.2 pg 684</Text-field>
<Text-field style="Heading 2" layout="Heading 2">1.</Text-field> f:=(x::Array)->2*x[2]^3; fy:=(x::Array)->6*x[2]^2; fyp:=(x::Array)->0; Zio2IydJInhHNiJJJkFycmF5RyUqcHJvdGVjdGVkR0YmNiRJKW9wZXJhdG9yR0YmSSZhcnJvd0dGJkYmLCQqJCY5JDYjIiIjIiIkRjFGJkYmRiY= Zio2IydJInhHNiJJJkFycmF5RyUqcHJvdGVjdGVkR0YmNiRJKW9wZXJhdG9yR0YmSSZhcnJvd0dGJkYmLCQqJCY5JDYjIiIjRjEiIidGJkYmRiY= Zio2IydJInhHNiJJJkFycmF5RyUqcHJvdGVjdGVkR0YmNiRJKW9wZXJhdG9yR0YmSSZhcnJvd0dGJkYmIiIhRiZGJkYm NonlinearShooting(f,fy,fyp,-1.0,0.0,1/2.0,1.0/3.0,20,10^(-5),10); NiUiIiJRJFRLPTYiJCErbm1tbTshIzU= NiUiIiNRJFRLPTYiJCErJz1LVFsjISM1 NiUiIiRRJFRLPTYiJCErJW9XKipcIyEjNQ== UTVpLH54LH53WzEsaV0sfndbMixpXTYi NiYiIiEkISM1ISIiJCIrKysrK11GJSQhKyVvVyoqXCNGJQ== NiYiIiIkISsrKysrJiohIzUkIitvOjB5W0YmJCEraTBbekJGJg== NiYiIiMkISsrKysrISohIzUkIitXLiI+dyVGJiQhK256Xm5BRiY= NiYiIiQkISsrKysrJikhIzUkIismPXI2bCVGJiQhK0lgRmpARiY= NiYiIiUkISsrKysrISkhIzUkIit4ZFlYWEYmJCErdSllZzEjRiY= NiYiIiYkISsrKysrdiEjNSQiK1omZVdXJUYmJCErVDVEdj5GJg== NiYiIickISsrKysrcSEjNSQiKyxKJXlNJUYmJCEraDFJISo9RiY= NiYiIigkISsrKysrbCEjNSQiK0EiUmBEJUYmJCErcFhyNT1GJg== NiYiIikkISsrKysrZyEjNSQiK18nKm9tVEYmJCErSC8wTzxGJg== NiYiIiokISsrKysrYiEjNSQiKzUoZTszJUYmJCErdzAiZm0iRiY= NiYiIzUkISsrKysrXSEjNSQiKyQ+SCsrJUYmJCErNW4kKipmIkYm NiYiIzYkISsrKysrWCEjNSQiK0U1Z0BSRiYkISs0YCF5YCJGJg== NiYiIzckISsrKysrUyEjNSQiK0MmKj1ZUUYmJCErZU9BejlGJg== NiYiIzgkISsrKysrTiEjNSQiK1VSaXRQRiYkISslUkVSVSJGJg== NiYiIzkkISsrKysrSSEjNSQiKyM9WVBxJEYmJCErJFxzO1AiRiY= NiYiIzokISsrKysrRCEjNSQiK1YnNGtqJEYmJCErVUVDQThGJg== NiYiIzskISsrKysrPyEjNSQiKyk+eTlkJEYmJCErTHBWdjdGJg== NiYiIzwkISsrKysrOiEjNSQiK3JfIykzTkYmJCErUkcySjdGJg== NiYiIz0kISsrKysrNSEjNSQiKztJTFtNRiYkIStJTSkqKT0iRiY= NiYiIz4kISorKysrJiEjNSQiK3g6KikqUSRGJiQhKy1lLFw2RiY= NiYiIz8kIiIhRiUkIitAJSlSTEwhIzUkISsob0g1NiJGKA== NiRRLGl0ZXJhdGlvbnM9NiIiIiQ= #define array w such that x=w[1] , y=w[2], y'=w[3] f:=(w::Array)->(1.0/25.0)*(1+w[3]^2)^(1/2); fy:=(w::Array)->0; fyp:=(w::Array)->(1.0/25.0)*w[3]*(1+w[3]^2)^(-1/2); Zio2IydJIndHNiJJJkFycmF5RyUqcHJvdGVjdGVkR0YmNiRJKW9wZXJhdG9yR0YmSSZhcnJvd0dGJkYmKigkIiM1ISIiIiIiJCIkXSNGL0YvKiQsJkYwRjAqJCY5JDYjIiIkIiIjRjAjRjBGOkYwRiZGJkYm Zio2IydJIndHNiJJJkFycmF5RyUqcHJvdGVjdGVkR0YmNiRJKW9wZXJhdG9yR0YmSSZhcnJvd0dGJkYmIiIhRiZGJkYm Zio2IydJIndHNiJJJkFycmF5RyUqcHJvdGVjdGVkR0YmNiRJKW9wZXJhdG9yR0YmSSZhcnJvd0dGJkYmKiokIiM1ISIiIiIiJCIkXSNGL0YvJjkkNiMiIiRGMCokLCZGMEYwKiRGMyIiI0YwI0YvRjpGMEYmRiZGJg== NonlinearShooting(f,fy,fyp,0,20,10,10,20,.00001,20) Warning, inserted missing semicolon at end of statement NiUiIiJRJFRLPTYiIiIh NiUiIiNRJFRLPTYiJCErY0A2R1EhIzU= NiUiIiRRJFRLPTYiJCErZXhCM1QhIzU= NiUiIiVRJFRLPTYiJCErJil5XjJUISM1 NiUiIiZRJFRLPTYiJCErdEtfMlQhIzU= UTVpLH54LH53WzEsaV0sfndbMixpXTYi NiYiIiFGIyIjNSQhK3RLXzJUISM1 NiYiIiJGIyQiK2dDdzUnKiEiKiQhK2VhRXlPISM1 NiYiIiNGIyQiKzwyOmsjKiEiKiQhK2BPKltEJCEjNQ== NiYiIiRGIyQiKz4saGYqKSEiKiQhK3MudE9HISM1 NiYiIiVGIyQiK09MbCdwKSEiKiQhK21rNUJDISM1 NiYiIiZGIyQiKyJlZltaKSEiKiQhK1orTzg/ISM1 NiYiIidGIyQiK1BSKFFIKSEiKiQhK01hJG9nIiEjNQ== NiYiIihGIyQiKyl5MU06KSEiKiQhKyxAKUc/IiEjNQ== NiYiIilGIyQiK2RMQmAhKSEiKiQhK2lqYDMhKSEjNg== NiYiIipGIyQiK1hMPiQqeiEiKiQhK3pxMSxTISM2 NiYiIzVGIyQiK3YxPnR6ISIqJCEkZiQhIzY= NiYiIzZGIyQiK1FMPiQqeiEiKiQiK2hqMSxTISM2 NiYiIzdGIyQiK1ZMQmAhKSEiKiQiK1RjYDMhKSEjNg== NiYiIzhGIyQiK25uU2AiKSEiKiQiK0g/KUc/IiEjNQ== NiYiIzlGIyQiKzRSKFFIKSEiKiQiK2lgJG9nIiEjNQ== NiYiIzpGIyQiK1kmZltaKSEiKiQiK3UqZkwsIyEjNQ== NiYiIztGIyQiKyVIYG1wKSEiKiQiKyNSMUpVIyEjNQ== NiYiIzxGIyQiK3EraGYqKSEiKiQiKylISW4kRyEjNQ== NiYiIz1GIyQiK2cxOmsjKiEiKiQiK3lOKltEJCEjNQ== NiYiIz5GIyQiKyZSaTJoKiEiKiQiKyNRbCN5TyEjNQ== NiYiIz9GIyQiK0cqKioqKioqKiEiKiQiKyc+QnY1JSEjNQ== NiRRLGl0ZXJhdGlvbnM9NiIiIiY=
UTVpLH54LH53WzEsaV0sfndbMixpXTYi NiYiIiFGIyIjNSQhK3RLXzJUISM1 NiYiIiJGIyQiK2dDdzUnKiEiKiQhK2VhRXlPISM1 NiYiIiNGIyQiKzwyOmsjKiEiKiQhK2BPKltEJCEjNQ== NiYiIiRGIyQiKz4saGYqKSEiKiQhK3MudE9HISM1 NiYiIiVGIyQiK09MbCdwKSEiKiQhK21rNUJDISM1 NiYiIiZGIyQiKyJlZltaKSEiKiQhK1orTzg/ISM1 NiYiIidGIyQiK1BSKFFIKSEiKiQhK01hJG9nIiEjNQ== NiYiIihGIyQiKyl5MU06KSEiKiQhKyxAKUc/IiEjNQ== NiYiIilGIyQiK2RMQmAhKSEiKiQhK2lqYDMhKSEjNg== NiYiIipGIyQiK1hMPiQqeiEiKiQhK3pxMSxTISM2 NiYiIzVGIyQiK3YxPnR6ISIqJCEkZiQhIzY= NiYiIzZGIyQiK1FMPiQqeiEiKiQiK2hqMSxTISM2 NiYiIzdGIyQiK1ZMQmAhKSEiKiQiK1RjYDMhKSEjNg== NiYiIzhGIyQiK25uU2AiKSEiKiQiK0g/KUc/IiEjNQ== NiYiIzlGIyQiKzRSKFFIKSEiKiQiK2lgJG9nIiEjNQ== NiYiIzpGIyQiK1kmZltaKSEiKiQiK3UqZkwsIyEjNQ== NiYiIztGIyQiKyVIYG1wKSEiKiQiKyNSMUpVIyEjNQ== NiYiIzxGIyQiK3EraGYqKSEiKiQiKylISW4kRyEjNQ== NiYiIz1GIyQiK2cxOmsjKiEiKiQiK3lOKltEJCEjNQ== NiYiIz5GIyQiKyZSaTJoKiEiKiQiKyNRbCN5TyEjNQ== NiYiIz9GIyQiK0cqKioqKioqKiEiKiQiKyc+QnY1JSEjNQ== NiRRLGl0ZXJhdGlvbnM9NiIiIiY= JSFH