Math 2280-2, Project #1, Problem 5
We first solve the crossbow example 2.3.3 using improved Euler's method with n=100
with(LinearAlgebra): with(plots): Digits:=16:
t0:=0.; tn:=10.; # first and last points in the interval
v0:=49.; # initial velocity
n:=100; # number of steps
h:=(tn-t0)/n; # step size
f:=(t,v)-> -(0.0011)*v*abs(v) -9.8; # slope function (rhs in DE dy/dx = f(x,y))
v:=v0: t:=t0:
ts1:=Vector(n+1): vs1:=Vector(n+1): ts1[1]:=t0: vs1[1]:=v0:
printf("%15s, %15s\134n","t","v"):
for i from 1 to n do
k1:= f(t,v): # left hand slope
k2:= f(t+h,v+h*k1): # approximation to right hand slope
k:=(k1+k2)/2: # averaged slope
v:= v + h*k: # improved Euler update
t:= t+h: # increase x
ts1[i+1]:=t: vs1[i+1]:=v:
if frac(i/10)=0 then
printf("%15.1f, %15.5f\134n",t,v):
end if:
od: # end for i loop
JCIiIUYj
JCIjNSIiIQ==
JCIjXCIiIQ==
IiQrIg==
JCIxKysrKysrKzUhIzs=
Zio2JEkidEc2IkkidkdGJUYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlLCYqKCQiIzYhIiUiIiI5JUYvLUkkYWJzRyUqcHJvdGVjdGVkRzYjRjBGLyEiIiQiIykqRjVGNUYlRiVGJQ==
t, v
1.0, 37.15469
2.0, 26.24277
3.0, 15.94529
4.0, 6.00406
5.0, -3.80200
6.0, -13.51045
7.0, -22.93557
8.0, -31.89841
9.0, -40.25568
10.0, -47.90661
And now with 200 points
t0:=0.; tn:=10.; # first and last points in the interval
v0:=49.; # initial velocity
n:=200; # number of steps
h:=(tn-t0)/n; # step size
f:=(t,v)-> -(0.0011)*v*abs(v) -9.8; # slope function (rhs in DE dy/dx = f(x,y))
v:=v0: t:=t0:
ts2:=Vector(n+1): vs2:=Vector(n+1): ts2[1]:=t0: vs2[1]:=v0:
printf("%15s, %15s\134n","t","v"):
for i from 1 to n do
k1:= f(t,v): # left hand slope
k2:= f(t+h,v+h*k1): # approximation to right hand slope
k:=(k1+k2)/2: # averaged slope
v:= v + h*k: # improved Euler update
t:= t+h: # increase x
ts2[i+1]:=t: vs2[i+1]:=v:
if frac(i/20)=0 then
printf("%15.1f, %15.5f\134n",t,v):
end if:
od: # end for i loop
JCIiIUYj
JCIjNSIiIQ==
JCIjXCIiIQ==
IiQrIw==
JCIxKysrKysrK10hIzw=
Zio2JEkidEc2IkkidkdGJUYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlLCYqKCQiIzYhIiUiIiI5JUYvLUkkYWJzRyUqcHJvdGVjdGVkRzYjRjBGLyEiIiQiIykqRjVGNUYlRiVGJQ==
t, v
1.0, 37.15474
2.0, 26.24292
3.0, 15.94555
4.0, 6.00444
5.0, -3.80159
6.0, -13.51019
7.0, -22.93545
8.0, -31.89845
9.0, -40.25588
10.0, -47.90696
We see that the two approximations agree with each other to two digits.
Here is a plot of both along with the true solution (from the file crossbow.mw that we used in class). All are indistinguishable!
rho2:=0.0011: g:=9.8:
C1:=arctan(v0*sqrt(rho2/g)): C2:=arctanh(v0*sqrt(rho2/g)):
tmax1:=C1/sqrt(rho2*g): tmax2:=C2/sqrt(rho2*g):
v3:=t->piecewise(t<tmax1,sqrt(g/rho2)*tan(C1-t*sqrt(rho2*g)),t>=tmax1,sqrt(g/rho2)*tanh(C2-(tmax2+t-tmax1)*sqrt(rho2*g)));
Zio2I0kidEc2IkYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlLUkqcGllY2V3aXNlRyUqcHJvdGVjdGVkRzYmMjkkSSZ0bWF4MUdGJSomLUklc3FydEc2JEYrSShfc3lzbGliR0YlNiMqJkkiZ0dGJSIiIkklcmhvMkdGJSEiIkY4LUkkdGFuR0YzNiMsJkkjQzFHRiVGOComRi5GOC1GMjYjKiZGOUY4RjdGOEY4RjpGODFGL0YuKiZGMUY4LUkldGFuaEdGMzYjLCZJI0MyR0YlRjgqJiwoSSZ0bWF4MkdGJUY4Ri5GOEYvRjpGOEZBRjhGOkY4RiVGJUYl
unassign('t'):
p1:=pointplot(Matrix([ts1,vs1]),symbol=diamond):
p2:=pointplot(Matrix([ts2,vs2]),symbol=cross):
p3:=plot(v3(t),t=0..10):
display({p1,p2,p3});
NiktJSdDVVJWRVNHNiQ3UzckJCIiISEiIiQiMS0rKysrKytcISM5NyQkIjFtbW07YXJ6QCEjOyQiMEpCNj9CPmolISM4NyQkIjFMTCRlOXVpMiUhIzskIjFXdi0hR3pNUyUhIzk3JCQiMW1tbSJ6XyI0aSEjOyQiMUYkPidcYF9eVCEjOTckJCIxbG1tVCZwaE4pISM7JCIxJEdoNFw1RyFSISM5NyQkIjFMTGUqPSlIXDUhIzokIjEiMzM1IikqeWZPISM5NyQkIjFtbSJ6LzN1QyIhIzokIjE9aG4uKSo9UU0hIzk3JCQiLURKJFJEWCIhIzYkIjFmQjQ2az83SyEjOTckJCIxbW0ielInb2s7ISM6JCIxTEVuU0MiPilIISM5NyQkIi1EMUo6dz0hIzYkIjEoKWVPaGFfYkYhIzk3JCQiMUxMTDNFbiQ0IyEjOiQiMXlIXikzbGNfIyEjOTckJCIxbW07L1JFJkcjISM6JCIxXGJlQjRdREIhIzk3JCQiK0QuJjRdIyEiKiQiMSZRIjRNMlwtQCEjOTckJCIrdkJfPEYhIiokIjFFJkdSTSN6ISk9ISM5NyQkIit2J0hpI0ghIiokIjFheFQxJT4hcDshIzk3JCQiMW1tInoqZXY6SiEjOiQiMXgjUl4pZTZ5OSEjOTckJCIxTExMMzQ3VEwhIzokIjFcaUFicWlfNyEjOTckJCIxTExMTFkuS04hIzokIjF4KSlIJD02RjEiISM5NyQkIjEqKipcN283VHYkISM6JCIxKjNnJlFsWkclKSEjOjckJCIxTExMJFEqb11SISM6JCIxYiJSJTR3IyoqWychIzo3JCQiMSoqKlw3PWxqOyUhIzokIjFgWCl5L0IjcFYhIzo3JCQiMSoqKlxQYVI8UCUhIzokIjFgYl4hUURSTiMhIzo3JCQiMUxMJGU5RWdlJSEjOiQiMUwlejFuJUhNRCEjOzckJCIxTExlUiIzR3klISM6JCExZl8+aillW24iISM6NyQkIjFtbTsvVDEmKlwhIzokITFjcSdIP2pKdiQhIzo3JCQiMW1tInpSUWJAJiEjOiQhMXM4ZV1zMDNmISM6NyQkIjEqKipcKD0+WTJhISM6JCExS0s1Y3YmKnl4ISM6NyQkIjFtbTt6WHU5YyEjOiQhMExfaGtaRXoqISM5NyQkIjEqKioqKipceSkpR2UhIzokITFSOE5XV04nPSIhIzk3JCQiMSoqKipcaV9RUWchIzokITEoKmUpPSFmJHlRIiEjOTckJCIxKioqXDd5JTNUaSEjOiQhMVRvIVt4ZTplIiEjOTckJCIxKioqKlxQIVtoWSchIzokITFOLSc0WGBdeiIhIzk3JCQiMUtMTCRReCRvbSEjOiQhMTomKXlbRkcmKT4hIzk3JCQiMSoqKioqXFArVilvISM6JCExaT9ZWkthJz0jISM5NyQkIjFtbSJ6cGUqenEhIzokIS9OS15SN25CISM3NyQkIjEqKioqKlwjXCdRSCghIzokITF4dDYseFdpRCEjOTckJCIxS0xlOVM4JlwoISM6JCExRz1UdXM0V0YhIzk3JCQiMSoqKlxpPz1icSghIzokITE4VT11bmlKSCEjOTckJCIxS0xMM3M/NnohIzokITAlZV11KkdENiQhIzg3JCQiLURKWGFFIikhIzYkITB4V3pZLyMqSCQhIzg3JCQiMW1tbW0qUlJMKSEjOiQhMSY+Z1AyZ2laJCEjOTckJCIxa207YTwuWSYpISM6JCExV1w3XEpXYU8hIzk3JCQiMUxMZTl0T2MoKSEjOiQhMGgiKkc9dyJHUSEjODckJCIqJlFrXCopISIpJCExQz1Eejs4JilSISM5NyQkIjFKTCQzZGc2PCohIzokITEocFl6KVJ4aFQhIzk3JCQiMWttbW14R3AkKiEjOiQhMWs+VT1JdztWISM5NyQkIjEpKipcN29LMGUqISM6JCExKj4oekgxJCl5VyEjOTckJCIxKSoqXCg9NXMjeSohIzokITFFZTpAVCMzaiUhIzk3JCQiJCsiISIiJCExZF0/OHBxIXolISM5LSUmQ09MT1JHNiYlJFJHQkckIiM1ISIiJCIiISEiIiQiIiEhIiItJSdQT0lOVFNHNmZ3NyQkIiIhISIiJCIkIVwhIiI3JCQiIiYhIiMkIjFQTHF6NSd6JFshIzk3JCQiIiIhIiIkIjFKRCpwVl5peCUhIzk3JCQiIzohIiMkIjExcVhkXidbciUhIzk3JCQiIiMhIiIkIjExWD1Ca3pgWSEjOTckJCIjRCEiIyQiMShHYC5eUklmJSEjOTckJCIiJCEiIiQiMSJbKSlveSllS1ghIzk3JCQiI04hIiMkIjE2MD01KFFDWiUhIzk3JCQiIiUhIiIkIjFITidbI1FlN1chIzk3JCQiI1ghIiMkIjEuPmxnKD1JTiUhIzk3JCQiIiYhIiIkIjAoR0JJI1FQSCUhIzg3JCQiI2IhIiMkIjFHTkBGcXRNVSEjOTckJCIiJyEiIiQiMSdSP1QtNWc8JSEjOTckJCIjbCEiIyQiMWkpXC88X3Y2JSEjOTckJCIiKCEiIiQiMCd5dyFcZSRmUyEjODckJCIjdiEiIyQiMTFzJkczQzkrJSEjOTckJCIiKSEiIiQiMWw1ImY2V1AlUiEjOTckJCIjJikhIiMkIjFuNXhHUUonKVEhIzk3JCQiIiohIiIkIjFjJzMjRyZHIkhRISM5NyQkIiMmKiEiIyQiMSJIW3NlJD1zUCEjOTckJCIjNSEiIiQiMXI2YFZXWjpQISM5NyQkIiQwIiEiIyQiMS0jPXhmJyoqZU8hIzk3JCQiIzYhIiIkIjFpNCQ+aFhGZyQhIzk3JCQiJDoiISIjJCIxaSVIIzNycllOISM5NyQkIiM3ISIiJCIxW2I3bm4hNFwkISM5NyQkIiREIiEiIyQiMT0zOEUuSk5NISM5NyQkIiM4ISIiJCIxbXpMeU4jKnpMISM5NyQkIiROIiEiIyQiMTRYLnJCdUNMISM5NyQkIiM5ISIiJCIxd3lNL0V3cEshIzk3JCQiJFgiISIjJCIxdDMjKkgtKVxAJCEjOTckJCIjOiEiIiQiMXZsaFw3UmdKISM5NyQkIiRiIiEiIyQiMUE8RDk8KmY1JCEjOTckJCIjOyEiIiQiMUF4TkF4eF5JISM5NyQkIiRsIiEiIyQiLyVvKj1hdSgqSCEjNzckJCIjPCEiIiQiMVpPViUqNCpRJUghIzk3JCQiJHYiISIjJCIxayJlS281LSpHISM5NyQkIiM9ISIiJCIxLlgnSHcrbiRHISM5NyQkIiQmPSEiIyQiMVsnekhiZEx5IyEjOTckJCIjPiEiIiQiMSJSV05UeCx0IyEjOTckJCIkJj4hIiMkIjF6WmtXbjp4RSEjOTckJCIjPyEiIiQiMVdcJ1wpPkhDRSEjOTckJCIkMCMhIiMkIjFBMSczaHo6ZCMhIzk3JCQiI0AhIiIkIjBzX2A4OyE+RCEjODckJCIkOiMhIiMkIjF3LTkyIilmbUMhIzk3JCQiI0EhIiIkIi9BajRASzlDISM3NyQkIiREIyEiIyQiMSc9KCoqZlo9aUIhIzk3JCQiI0IhIiIkIjE3N0IzRj01QiEjOTckJCIkTiMhIiMkIjFiP0NPRUplQSEjOTckJCIjQyEiIiQiMWwtJnBEcmw/IyEjOTckJCIkWCMhIiMkIjFhZ1M4YCZcOiMhIzk3JCQiI0QhIiIkIjEmZUB6ZGhNNSMhIzk3JCQiJGIjISIjJCIxTCI9NyZvM18/ISM5NyQkIiNFISIiJCIwaihlaHojMysjISM4NyQkIiRsIyEiIyQiMTAib1N3Im9cPiEjOTckJCIjRiEiIiQiMVZEalJeaykqPSEjOTckJCIkdiMhIiMkIjEsMlIlKlxyWj0hIzk3JCQiI0chIiIkIjFKTSEpZSMpKW96IiEjOTckJCIkJkchIiMkIjF0ITRwKT07WTwhIzk3JCQiI0ghIiIkIjF6OWNiR2AmcCIhIzk3JCQiJCZIISIjJCIwRXpPOykqXGsiISM4NyQkIiNJISIiJCIxUWVdSltiJWYiISM5NyQkIiQwJCEiIyQiMS4qeikqKikqPlc6ISM5NyQkIiNKISIiJCIxKillXkgvJFJcIiEjOTckJCIkOiQhIiMkIjF5U0grTnVWOSEjOTckJCIjSyEiIiQiMWkqZjBAT09SIiEjOTckJCIkRCQhIiMkIjEjXEtrbjBPTSIhIzk3JCQiI0whIiIkIjFjMzdKIVxPSCIhIzk3JCQiJE4kISIjJCIxTiVcVVVqUEMiISM5NyQkIiNNISIiJCIxJWUhPkBnJVI+IiEjOTckJCIkWCQhIiMkIjEwXFMtUz5XNiEjOTckJCIjTiEiIiQiMXcnKXlpWF0lNCIhIzk3JCQiJGIkISIjJCIxNCJHNSJcKFsvIiEjOTckJCIjTyEiIiQiMTE2ZypvQUkmKiohIzo3JCQiJGwkISIjJCIxdkJSNChReVgqISM6NyQkIiNQISIiJCIxY181SydwSicqKSEjOjckJCIkdiQhIiMkIjE3KnlHLikpKm8lKSEjOjckJCIjUSEiIiQiMSwheUxmbV8oeiEjOjckJCIkJlEhIiMkIjF3IXpjNHk+WyghIzo3JCQiI1IhIiIkIjFHUD87YTQqKXAhIzo3JCQiJCZSISIjJCIxJGYmcD46ZidcJyEjOjckJCIjUyEiIiQiMVQkUUtYUlcrJyEjOjckJCIkMCUhIiMkIjExPUpTQmg3YiEjOjckJCIjVCEiIiQiMFkkKVtQJDNAXSEjOTckJCIkOiUhIiMkIjFZQ2A6ZSMpSFghIzo3JCQiI1UhIiIkIjFBTWd6SCIpUVMhIzo3JCQiJEQlISIjJCIxJVsoUVAjPSFbTiEjOjckJCIjViEiIiQiMWh1SzFdVGRJISM6NyQkIiROJSEiIyQiMScqW0NYbihwYyMhIzo3JCQiI1chIiIkIjFOZ2VbcG53PyEjOjckJCIkWCUhIiMkIjElPSVvUyIqWydlIiEjOjckJCIjWCEiIiQiMW9kLnFvUSc0IiEjOjckJCIkYiUhIiMkIjFFWyhlLlBNMSchIzs3JCQiI1khIiIkIjFpKj0uQEtMOyIhIzs3JCQiJGwlISIjJCExSzYkUTZMbXQkISM7NyQkIiNaISIiJCExWXlPJWYqUU8nKSEjOzckJCIkdiUhIiMkITF1cXBoIW9OTiIhIzo3JCQiI1shIiIkITFrd1dBVVVWPSEjOjckJCIkJlshIiMkITFaSXBoNT1MQiEjOjckJCIjXCEiIiQhMWYlUTNBN0cjRyEjOjckJCIkJlwhIiMkITE+K09vOEg3TCEjOjckJCIjXSEiIiQhMFtpYT8jZixRISM5NyQkIiQwJiEiIyQhMT0/c3IlKW8hSCUhIzo3JCQiI14hIiIkITEiWzYyJlJielohIzo3JCQiJDomISIjJCExX2JndkM7b18hIzo3JCQiI18hIiIkITEkR2paJHpbY2QhIzo3JCQiJEQmISIjJCExJFsmR3hVXVdpISM6NyQkIiNgISIiJCExWldhPWI9S24hIzo3JCQiJE4mISIjJCExbE5vWGRdPnMhIzo3JCQiI2EhIiIkITFPVjNCIlJrcSghIzo3JCQiJFgmISIjJCExbSQpeigqKWZIPikhIzo3JCQiI2IhIiIkITAlPiVcU1UheicpISM5NyQkIiRiJiEiIyQhMTszMHQ1bWsiKiEjOjckJCIjYyEiIiQhMU5HUkgvelwnKiEjOjckJCIkbCYhIiMkITE4T181MFc4NSEjOTckJCIjZCEiIiQhMWRxNSUpeiU9MSIhIzk3JCQiJHYmISIjJCExQiF5IlwqKT41NiEjOTckJCIjZSEiIiQhMFp3SiE0XGU2ISM4NyQkIiQmZSEiIyQhMTAyK2M4czE3ISM5NyQkIiNmISIiJCExTCczMCR5KVtEIiEjOTckJCIkJmYhIiMkITFcSyxqeSlISSIhIzk3JCQiI2chIiIkITF2QHguIT41TiIhIzk3JCQiJDAnISIjJCExUTVaPCl5KilSIiEjOTckJCIjaCEiIiQhMSRcK1ApWydvVyIhIzk3JCQiJDonISIjJCExL25VKHp1WVwiISM5NyQkIiNpISIiJCExJzRiJXBoU1U6ISM5NyQkIiREJyEiIyQhMSxzKW9pYytmIiEjOTckJCIjaiEiIiQhMTsvZDhRaVA7ISM5NyQkIiROJyEiIyQhMTYwYSFSMF5vIiEjOTckJCIjayEiIiQhMTxuWE8hKlxLPCEjOTckJCIkWCchIiMkITFSI0ghW0MhKXo8ISM5NyQkIiNsISIiJCExTCFSL004cSM9ISM5NyQkIiRiJyEiIyQhMTMoXHhXSFQoPSEjOTckJCIjbSEiIiQhMTM5SkImWzYjPiEjOTckJCIkbCchIiMkITFIbDpTJG8hbz4hIzk3JCQiI24hIiIkITFeclEicCcpWywjISM5NyQkIiR2JyEiIyQhMWRTYiFSLDsxIyEjOTckJCIjbyEiIiQhMThzLXMtQDNAISM5NyQkIiQmbyEiIyQhMTB3TiI+Nlo6IyEjOTckJCIjcCEiIiQhMEtJYy0tNj8jISM4NyQkIiQmcCEiIyQhMXAoM1FuIVFaQSEjOTckJCIjcSEiIiQhMWwqcHEwWE5IIyEjOTckJCIkMCghIiMkITFgMjg+SmZSQiEjOTckJCIjciEiIiQhMUFcY0VHXyZRIyEjOTckJCIkOighIiMkITE5NjZwQExKQyEjOTckJCIjcyEiIiQhMWg5KCpmIj5xWiMhIzk3JCQiJEQoISIjJCExKikzNU89ZUFEISM5NyQkIiN0ISIiJCExQHFbZSM9IW9EISM5NyQkIiROKCEiIyQhMVIxVTdsSzhFISM5NyQkIiN1ISIiJCExVm12MlpdZUUhIzk3JCQiJFgoISIjJCExJFtsIno0Yi5GISM5NyQkIiN2ISIiJCExPF1QJ1tqJVtGISM5NyQkIiRiKCEiIyQhMHgtV1RTS3ojISM4NyQkIiN3ISIiJCExcyZ5UCgqenkkRyEjOTckJCIkbCghIiMkITFkdnYrL1EjKUchIzk3JCQiI3ghIiIkITE8Tl9kKlJuI0ghIzk3JCQiJHYoISIjJCExJmYjUUtwJjQoSCEjOTckJCIjeSEiIiQhMUV0JSpSJ0hdLCQhIzk3JCQiJCZ5ISIjJCExSzVKQGsmKmVJISM5NyQkIiN6ISIiJCExd0RAV2N0LUohIzk3JCQiJCZ6ISIjJCExLDo+LmRPWUohIzk3JCQiIyEpISIiJCExcE50Pl0lKSo9JCEjOTckJCIkMCkhIiMkITB0MUMvc0pCJCEjODckJCIjIikhIiIkITFjdilwQ1hqRiQhIzk3JCQiJDopISIjJCExeiJ5bDhqJD5MISM5NyQkIiMjKSEiIiQhMXNQclRVQWlMISM5NyQkIiREKSEiIyQhMUIxWT9yI1xTJCEjOTckJCIjJCkhIiIkITFmW11lLlpaTSEjOTckJCIkTikhIiMkITFzPUJwRCYpKlskISM5NyQkIiMlKSEiIiQhMTtseiRScz9gJCEjOTckJCIkWCkhIiMkITFMVD0sJkdUZCQhIzk3JCQiIyYpISIiJCExI1xpIyllPmdoJCEjOTckJCIkYikhIiMkITExWiMpelZ1ZE8hIzk3JCQiIycpISIiJCExNktpRztJKnAkISM5NyQkIiRsKSEiIyQhMSYqW1I6LHBTUCEjOTckJCIjKCkhIiIkITFodSgpWyczPnkkISM5NyQkIiR2KSEiIyQhMTlyImUxY0gjUSEjOTckJCIjKSkhIiIkITF1dyc0QkpRJ1EhIzk3JCQiJCYpKSEiIyQhMTY3M1BJYC9SISM5NyQkIiMqKSEiIiQhMTAuKltTZ10lUiEjOTckJCIkJiopISIjJCExTz8zJEc3YSlSISM5NyQkIiMhKiEiIiQhMTxRRVt3ZURTISM5NyQkIiQwKiEiIyQhMSI9QFxdJmVsUyEjOTckJCIjIiohIiIkITFOd08mKVtTMFQhIzk3JCQiJDoqISIjJCExLW5vXFsvWFQhIzk3JCQiIyMqISIiJCExdWltJlsvWD0lISM5NyQkIiREKiEiIyQhMSlSRCgzSHlCVSEjOTckJCIjJCohIiIkITE5TyQ9RXpHRSUhIzk3JCQiJE4qISIjJCEwKkhVOkZ6LFYhIzg3JCQiIyUqISIiJCExbUhIbkNfU1YhIzk3JCQiJFgqISIjJCExYSIzRHVuIXpWISM5NyQkIiMmKiEiIiQhMUIjKkckekZ1VCUhIzk3JCQiJGIqISIjJCEuKCopKSo9Z2JXISM2NyQkIiMnKiEiIiQhMCUqNGJPKmUkXCUhIzg3JCQiJGwqISIjJCExKDMmNEUmKlFKWCEjOTckJCIjKCohIiIkITEib3MtdSwhcFghIzk3JCQiJHYqISIjJCExWVg8JVJEa2clISM5NyQkIiMpKiEiIiQhMTBrSCsqZk9rJSEjOTckJCIkJikqISIjJCExbVVNKHAvMm8lISM5NyQkIiMqKiEiIiQhMUZdMl0jZnZyJSEjOTckJCIkJioqISIjJCExVDk4XElBYVohIzk3JCQiJCsiISIiJCExJWZyM2gmcCF6JSEjOS0lJ1NZTUJPTEc2IyUmQ1JPU1NHLSUnUE9JTlRTRzZicTckJCIiISEiIiQiJCFcISIiNyQkIiIiISIiJCIxJHBPQjFeaXglISM5NyQkIiIjISIiJCIxPU0hb2gmemBZISM5NyQkIiIkISIiJCIxIkh3YlwoZUtYISM5NyQkIiIlISIiJCIxJzQhKiopKj5lN1chIzk3JCQiIiYhIiIkIjFFdVpCZXQkSCUhIzk3JCQiIichIiIkIjEjKlxBJCpwK3dUISM5NyQkIiIoISIiJCIxNCE0YHphJGZTISM5NyQkIiIpISIiJCIxUCIpMz0oUlAlUiEjOTckJCIiKiEiIiQiMSx4ayNSQiJIUSEjOTckJCIjNSEiIiQiMW1INVAmb2FyJCEjOTckJCIjNiEiIiQiMXJAZS4qUUZnJCEjOTckJCIjNyEiIiQiMSIpPnNGIyoqM1wkISM5NyQkIiM4ISIiJCIxOCE+MD06KnpMISM5NyQkIiM5ISIiJCIxbDdSQUx2cEshIzk3JCQiIzohIiIkIjFrRUpmNVFnSiEjOTckJCIjOyEiIiQiMUBFKjNnbTwwJCEjOTckJCIjPCEiIiQiMS8wRz8qeVElSCEjOTckJCIjPSEiIiQiMV9CeTp4b09HISM5NyQkIiM+ISIiJCIwLCUzdUw7SUYhIzg3JCQiIz8hIiIkIjEuNylcJHBGQ0UhIzk3JCQiI0AhIiIkIjAnSF9kKys+RCEjODckJCIjQSEiIiQiMTEtWygpXEk5QyEjOTckJCIjQiEiIiQiMSNRY2dfayxKIyEjOTckJCIjQyEiIiQiMTk7dyoqPmIxQSEjOTckJCIjRCEiIiQiMEN3OEJUTTUjISM4NyQkIiNFISIiJCIxN04kPl4xMysjISM5NyQkIiNGISIiJCIxX1xsdERpKSo9ISM5NyQkIiNHISIiJCIxa194algnb3oiISM5NyQkIiNIISIiJCIwJm8+PiEzYnAiISM4NyQkIiNJISIiJCIxcG8iPiUpR1hmIiEjOTckJCIjSiEiIiQiMXVDPHZLIVJcIiEjOTckJCIjSyEiIiQiMTpCRCEpeWckUiIhIzk3JCQiI0whIiIkIjFyWSNSXj5PSCIhIzk3JCQiI00hIiIkIjFeWVYxYCJSPiIhIzk3JCQiI04hIiIkIjEvWy1TRVolNCIhIzk3JCQiI08hIiIkIjFoQk56N3BfKiohIzo3JCQiI1AhIiIkIjEneTMlUWYjRycqKSEjOjckJCIjUSEiIiQiLy0mb141XCh6ISM4NyQkIiNSISIiJCIxMFtwZW9zKSlwISM6NyQkIiNTISIiJCIxJW9Id0plUysnISM6NyQkIiNUISIiJCIxYz0hXGMqbz9dISM6NyQkIiNVISIiJCIwV14nKlIxJVFTISM5NyQkIiNWISIiJCIxeSRIOGMmKnAwJCEjOjckJCIjVyEiIiQiMWgoW1FhV2kyIyEjOjckJCIjWCEiIiQiMUJWNzY5JWY0IiEjOjckJCIjWSEiIiQiMSNcIWZPaHVlNiEjOzckJCIjWiEiIiQhMUIiUltiXTNrKSEjOzckJCIjWyEiIiQhMSRlWTUocCZRJT0hIzo3JCQiI1whIiIkITEiMydbPTpCQkchIzo3JCQiI10hIiIkITFkWyEqW3oqPiFRISM6NyQkIiNeISIiJCExbzYyUGclKnpaISM6NyQkIiNfISIiJCExcyN6O0Vtb3YmISM6NyQkIiNgISIiJCEwQ0NWKSpcRHQnISM5NyQkIiNhISIiJCExd3cnb2kqeTF4ISM6NyQkIiNiISIiJCExVmMobyUpeSR6JykhIzo3JCQiI2MhIiIkITE8aCIzcjcsbCohIzo3JCQiI2QhIiIkITAtcmt5eT0xIiEjODckJCIjZSEiIiQhMURSPnEtX2U2ISM5NyQkIiNmISIiJCExej16X2QiXEQiISM5NyQkIiNnISIiJCExa2Y4c2EvXjghIzk3JCQiI2ghIiIkITE0Wj0qKSkqKW9XIiEjOTckJCIjaSEiIiQhMTghW0xxSENhIiEjOTckJCIjaiEiIiQhMVVDVm5la1A7ISM5NyQkIiNrISIiJCExS1g5LSc+RHQiISM5NyQkIiNsISIiJCExPSs4NUMuRj0hIzk3JCQiI20hIiIkITF3I0cmKjNtNiM+ISM5NyQkIiNuISIiJCExc0UtWkYhXCwjISM5NyQkIiNvISIiJCExUipvLiJbQTNAISM5NyQkIiNwISIiJCExMl9RU102LEEhIzk3JCQiI3EhIiIkITEneilRVWxiJEgjISM5NyQkIiNyISIiJCEwVm1xeEtiUSMhIzg3JCQiI3MhIiIkITFjZ3dxdi14QyEjOTckJCIjdCEiIiQhMUQrP0ReLW9EISM5NyQkIiN1ISIiJCExVTRiRSteZUUhIzk3JCQiI3YhIiIkITFXSShSRG4lW0YhIzk3JCQiI3chIiIkITFCTlsoPSMpeSRHISM5NyQkIiN4ISIiJCExTChSX2hTbiNIISM5NyQkIiN5ISIiJCExWygqPlMoR10sJCEjOTckJCIjeiEiIiQhMEBsaD1MRjUkISM4NyQkIiMhKSEiIiQhMSUpbz8uNSUpKj0kISM5NyQkIiMiKSEiIiQhMVxaWnMnUmpGJCEjOTckJCIjIykhIiIkITFcZVI1ckBpTCEjOTckJCIjJCkhIiIkITFvV0VzO1laTSEjOTckJCIjJSkhIiIkITEqW0lfOmk/YCQhIzk3JCQiIyYpISIiJCExW0lwK3krO08hIzk3JCQiIycpISIiJCExJnkyaEgpRypwJCEjOTckJCIjKCkhIiIkITFLckB3UCo9eSQhIzk3JCQiIykpISIiJCExX1RyQlsiUSdRISM5NyQkIiMqKSEiIiQhMXQ5TnBDL1hSISM5NyQkIiMhKiEiIiQhL3ldIj5vYi0lISM3NyQkIiMiKiEiIiQhMSVSRF8iUlEwVCEjOTckJCIjIyohIiIkITFCRnQ1P1slPSUhIzk3JCQiIyQqISIiJCExZHNZIkhiR0UlISM5NyQkIiMlKiEiIiQhMWxTaDZxXFNWISM5NyQkIiMmKiEiIiQhMF53SiczUzxXISM4NyQkIiMnKiEiIiQhMSkpZWZzNGMkXCUhIzk3JCQiIygqISIiJCExZXQkcCo9KCpvWCEjOTckJCIjKSohIiIkITFUKmYnPidHT2slISM5NyQkIiMqKiEiIiQhMUEoKSlmYUV2ciUhIzk3JCQiJCsiISIiJCExLHghelxoMXolISM5LSUnU1lNQk9MRzYjJShESUFNT05ERy0lJVZJRVdHNiQ7JCIiISEiIiQiJCsiISIiJShERUZBVUxURy0lKkFYRVNTVFlMRUc2IyUnTk9STUFMRy0lKFNDQUxJTkdHNiMlLlVOQ09OU1RSQUlORURHLSUlUk9PVEc2Jy0lKUJPVU5EU19YRzYjJCIkIVIhIiItJSlCT1VORFNfWUc2IyQiIyEqISIiLSUtQk9VTkRTX1dJRFRIRzYjJCIlIVwkISIiLSUuQk9VTkRTX0hFSUdIVEc2IyQiJTVRISIiLSUpQ0hJTERSRU5HNiI=
(a) to obtain the value where v=0, we run improved Euler's method and select the time for which v is closest to 0. Since all computations are already stored in a vector here, we can compute the minimum in Maple (as done below). Normally (in higher dimensional problems) you do not want to store the value of your solution at all points, so the numerical method's loop could be modified to find the time for which v is minimum (with no additional memory cost).
# find the index of the minimum value in the vector vs2
res:=rtable_scanblock( abs(vs2), [rtable_dims(vs2)],(val,ind,res) -> `if`(val<res[2],[ind,val],res),[[1],abs(vs2[1])]);
printf("approximation = %4.2f seconds\134n",ts2[res[1][1]]);
NyQ3IyIjJCokIjFpKj0uQEtMOyIhIzs=
approximation = 4.60 seconds
(b) To estimate the impact velocity after 9.41 seconds, we could use improved Euler with tn=9.41, to get:
t0:=0.; tn:=9.41; # first and last points in the interval
v0:=49.; # initial velocity
n:=200; # number of steps
h:=(tn-t0)/n; # step size
f:=(t,v)-> -(0.0011)*v*abs(v) -9.8; # slope function (rhs in DE dy/dx = f(x,y))
v:=v0: t:=t0:
for i from 1 to n do
k1:= f(t,v): # left hand slope
k2:= f(t+h,v+h*k1): # approximation to right hand slope
k:=(k1+k2)/2: # averaged slope
v:= v + h*k: # improved Euler update
t:= t+h: # increase x
od: # end for i loop
printf("v( %4.2f ) approx. equals to %4.2f m/s\134n",t,v);
JCIiIUYj
JCIkVCohIiM=
JCIjXCIiIQ==
IiQrIw==
JCIxKysrKysrMFohIzw=
Zio2JEkidEc2IkkidkdGJUYlNiRJKW9wZXJhdG9yR0YlSSZhcnJvd0dGJUYlLCYqKCQiIzYhIiUiIiI5JUYvLUkkYWJzRyUqcHJvdGVjdGVkRzYjRjBGLyEiIiQiIykqRjVGNUYlRiVGJQ==
v( 9.41 ) approx. equals to -43.48 m/s