Math242 Lab2 September 2009 M242lab2.mws
| > |
with(plots): # plots package |
House keeping item for L'Hopital's Rules
| > |
limit( x*exp(1/x) - x, x=infinity); # please note the name of infinity |
 |
(1) |
| > |
limit( (1 - cos(x))/x^2, x=0); |
 |
(2) |
| > |
limit( (1 + 2/n)^n, n = infinity); |
 |
(3) |
Try this #37 in 4.4
 |
(4) |
From this answer above, how would you approximate the value of cos(0.1) 'by hand'? The value of cos(0.1) is
 |
(5) |
Simple introduction to 'do' loop and the command 'seq'
Suppose we would like to generate a set of numbers obtained from the Newton's method, say in solving for the first positive root of tan(x) = x, i.e. f(x) = tan(x) - x = 0.
Note: After you have plotted the graph, insert after y = -10..10 the option ' discont = true ' . Separate the option by a comma.
| > |
plot( [ tan(x), x], x=0..3*Pi/2, y=-10..10); |
We see that the first positive root is somewhere larger than 4.
| > |
F:=x -> (tan(x) - x)/(sec(x)^2 - 1): |
| > |
for k from 1 to 10 do a[k]:= a[k-1] - F(a[k-1]): od: #make sure you use colon |
Command sequence
| > |
tan(a[10]) - a[10]: # check |
This example shows that the sequence approaches the limit slowly. Replace the colons by semi-colons.
| > |
seq((1+ 1/n)^n, n=3..20): evalf(%): evalf(exp(1)): |
Simplier Example of seq - Fibonacci sequence.
| > |
a:= 'a': # clear the value of a, if needed |
| > |
for k from 3 to 20 do a[k]:= a[k-2] + a[k-1]: od: |
| > |
evalf(a[20]/a[19]); fsolve(x^2-x-1 = 0, x); #the larger root is the golden mean. |
 |
 |
(7) |
--------------------------------- Pendulum problem -------------------------------------
| > |
G:=(a,b)-> int(1/sqrt(a^2*cos(t)^2 + b^2*sin(t)^2), t=0..Pi/2); |
 |
(8) |
| > |
G(3,4); evalf(%); #You should try other values of a and b. |
 |
(9) |
 |
(9) |
AM and GM sequences
| > |
p[0]:= 3.: q[0]:= 4.: # Pick initial values you chose above. |
| > |
for k from 1 to 6 do p[k]:= (p[k-1]+q[k-1])/2: q[k]:= sqrt(p[k-1]*q[k-1]): od: |
| > |
seq([p[k],q[k]], k=1..6): |
Part A of your lab exercise is to write a simple do loop using the Babylonian algorithm to compute the square root of each of p[k-1]*q[k-1] . See the first hand out.
-------------- Example of solving a first order differential equations ------------------
The first part may look difficult. Actually it is in Math 302 syllabus. But don't worry here.
| > |
deq:= diff(y(t), t$2) + 2*sin(y(t)) = 0; |
 |
(10) |
| > |
sol:= dsolve({deq, y(0)=Pi/3, D(y)(0)=0}, y(t), numeric); |
 |
(11) |
| > |
p1:= odeplot(sol, t=0..10): |
As seen in class, by using the energy method, this second order equation is equivalent to solving for
| > |
Diff(y(t), t)^2 - 2*cos(y(t)) = C; |
 |
(12) |
The constant C on the right hand side is calculated from the initial conditions: y(0) = Pi/3, y'(0) = 0.
Using the first difference scheme as discussed, we have z[k+1] = z[k] + f(z[k])*h .
| > |
z[0]:= evalf(Pi/3): h:=0.01: |
| > |
for k from 0 to 200 do z[k+1]:= z[k]- sqrt(4*cos(z[k]) - 2)*h : od: |
| > |
pendpts:= seq([k*h, z[k]],k=0..200): |
| > |
p2:= plot( {pendpts}, style=point, symbol= point, color=blue): |
We'll display the two graphs, the true graph and the approximate one together.
In the graph above, try to increase the number of points. However, you'll see the blue dots never could get over the 'hump'. Give a thought on this aspect of setting up an algorithm to enable the blue graph go over the hump and still follows the original one closely.
Second lab exercise.
Do the following problem on population growth: Let P(t) represent the population of a specie at time t. Under some general assumption, we let P' = 1/2*P*(10 - P). We assume the initial population P(0) = 3 units say.
| > |
eq2:= diff(P(t),t) = 1/2*P(t)*(10 - P(t)); |
 |
(13) |
| > |
sol2:= dsolve({eq2, P(0)= 4}, P(t)); |
 |
(14) |
| > |
plot(rhs(sol2), t=0..4):# new command here 'rhs' |
Part B. Now solve for the solution by first difference scheme by choosing h = 0.01 again. Imitate the commands above to compute the approximate solution by the first order finite difference scheme. Iterate 200 times, i.e. over the time interval [0, 2]. Display the true graph sol2 together with the approximate solution together. Your answer should look like this.
Please start on a new work sheet and work out the two parts A and B nicely. Delete any errors and tidy things up.
