[SOLVED] Amath 481 / 581 homework 1 – initial value problems

1. Consider the ODE
dy(t)
dt = −3y(t) sin t, y(t = 0) = π
√
2
,
which has the exact solution y(t) = πe3(cos t−1)/
√
2 (you can verify that). Implement the methods
forward Euler and Heun’s for this ODE to test the error as a function of ∆t . In particular:(a) Solve the ODE numerically using the forward Euler method:
y(tn+1) = y(tn) + ∆tf(tn, y(tn))
with t = [0 : ∆t : 5], where ∆t = 2−2
, 2
−3
, 2
−4
, . . . , 2
−8
.For each of these ∆t values calculate the
error E = mean(abs(ytrue −ynum)) of the numerical method. Plot log(∆t) on the x axis and log(E)
on the y axis. Using polyfit, find the slope of the best fit line through this data. This is the order
of the forward Euler method.ANSWER: Save your last numerical solution (∆t = 2−8
) as a column vector in A1. Save the
error values in a row vector with seven components in A2. Save the slope of the line in A3.(b) Solve the ODE numerically using Heun’s method:
y(tn+1) = y(tn) + ∆t
2
[f(tn, y(tn)) + f(tn + ∆t, y(tn) + ∆tf(tn, y(tn)))]
with t = [0 : ∆t : 5], where ∆t = 2−2
, 2
−3
, 2
−4
, . . . , 2
−8
.For each of these ∆t values, calculate
the error E = mean(abs(ytrue − ynum)) of the numerical method. Plot log(∆t) on the x axis and
log(E) on the y axis. Using polyfit, find the slope of the best fit line through this data. This is
the order of the Heun’s method.ANSWER: Save your last numerical solution (∆t = 2−8
) as a column vector in A4. Save the
error values in a row vector with seven components in A5. Save the slope of the line in A6.2. Consider the van der Pol oscillator
d
2y(t)
dt2
+ ǫ[y
2
(t) − 1]dy(t)
dt + y(t) = 0
with ǫ being a parameter.(a) With ǫ = 0.1, solve the equation for t = [0 : 0.5 : 32] using ode45. The initial conditions are
y(t = 0) = √
3 and dy(t = 0)/dt = 1. Repeat this for ǫ = 1 and ǫ = 20.ANSWER: Save the solutions y(t) for different ǫ as a matrix of 3 columns in A7.(b) Using the time span t = [0, 32] (the step size for displaying the result is not specified), solve
the van der Pol’s equation with ode45. Use ǫ = 1 and the initial conditions y(t = 0) = 2 and
1
dy(t = 0)/dt = π
2
.Below is an example on how to control the error tolerance TOL in ode45:
TOL = 1e-4;
options = odeset(‘AbsTol’,TOL,‘RelTol’,TOL);
[T,Y] = ode45(‘rhs’,tspan,y0,options);Using the diff and mean commands on the vector T shown above, calculate the average stepsize t needed to solve the problem for each of the following tolerance values: 10−4
, 10−5
, . . . , 10−10
.
Plot log(∆t) on the x axis and log(TOL) on the y axis. Using polyfit, find the slope of the best fit
line through this data. This is the order of the local truncation error of ode45. Repeat this with
ode23 and ode113.ANSWER: The slopes should be saved in variables A8 – A10 for ode45, ode23, and ode113 respectively.3. To explore interaction between neurons, implement two Fitzhugh neurons coupled via linear
coupling:
dv1
dt = −v
3
1 + (1 + a1)v
2
1 − a1v1 − w1 + I + d12v2
dw1
dt = bv1 − cw1
dv2
dt = −v
3
2 + (1 + a2)v
2
2 − a2v2 − w2 + I + d21v1
dw2
dt = bv2 − cw2
with parameters a1 = 0.05, a2 = 0.25 , b = c = 0.01 and, I = 0.1. Start the simulations with the
initial condition of (v1(0), v2(0)) = (0.1, 0.1) and (w1(0), w2(0)) = (0, 0) and use the ode15s solver.Set the interaction parameters such that d12 is negative and d21 is positive. What do you observe
from the different graphical representations of the solutions?
ANSWERS: Set the interaction parameters to 5 different values
(d12, d21): (0,0), (0,0.2), (-0.1,0.2), (-0.3,0.2), (-0.5,0.2).For each interaction value solve the system
for t = [0 : 0.5 : 100] and save the computed solution, (v1, v2, w1, w2), in a 201 × 4 matrix. Write
out the 5 different solutions, each of which corresponds to an interaction parameter, in A11 – A15.