Matematik och Datavetenskap, Chalmers Tekniska Högskola och Göteborgs Universitet
1. Read AMBS Ch 40.2. Write a program based on the Euler forward method for solving the initial value problem for a general system of ordinary differential equations:
u'(t)=f(t,u(t)), a< t <b,
u(a)=ua.
You must use the template my_ode.m. This program is VERY SIMILAR to my_exp.m. and my_trig.m that you have already written. Use this fact and copy stuff from your old programs.
2. You MUST NOT include any plot commands in the program my_ode. Plotting is a separate task that we may, or may not, want to do after we have solved the differential equation. It is convenient to put such plot commands in a script file called, for example, test1.m:
f='funk';
int=[0 1];
ua=[0;1];
h=1e-3;
[t,U]=my_ode(f,int,ua,h);
figure(1)
plot(t,U) % solution curves
figure(2)
plot(U(:,1),U(:,2)) % phase portrait
axis equal
To execute these commands you type test1 in the matlab command window. This is more convenient than to re-type the commands every time you want to run your program.
3. Solve the following equations by hand and then use them to test your program. (Note: you have already done this in earlier exercises. Now do it again.)
a. u'(t)=1/t, u(1)=0
b. u'=u, u(0)=1
c. u'=-u, u(0)=1
d. u''+u=0, u(0)=0, u'(0)=1
e. u''-u=0, u(0)=0, u'(0)=1
You can compare the solution formula and the numerical solution by plotting them in the same figure. For example, if the solution formula is u(t)=3*exp(t):
plot(t,U)
hold on
fplot('3*exp(x)',[0 1],'r')
NOTE: Matlabs built-in function exp.m is NOT AN EXACT SOLUTION, it is a numerical approximation just as the solution that you compute with my_ode and my_exp.
4. The Matlab programs ode23, ode45, ode15s do more or less the same thing as my_ode. Read about it: help ode23 or in helpdesk or in Pärt-Enander-Sjöberg. Use ode23 to solve some of the examples. Compare the numerical solutions obtained by my_ode and ode23.
[t,U]=ode23('funk',int,ua)
/stig