MATH661 Homework 10 - Multiple operators

This homework is a reading assignment meant as final examination preparation. The particular topic is numerical solution of ordinary differential equations, but the presentation is also a synthesis of the ideas encountered in the course. The scheduled course grade points are automatically awarded, but several suggestions are given for independent work that can be submitted for extra course credit up to the date prior to the final examination.

1Problem and algorithm statement

The problem considered here is finding approximations of y𝒴, y: a real-valued function with p continuous derivatives specified by an equality between the action of two operators

y=y,

with linear and nonlinear. The particular choice of operators =d/dt and =f corresponds to the ordinary differential equation

y'=f(t,y),

but the ideas underlying the approach are generic, applicable to other operator choices.

The function y: is approximated by a sample 𝒟={(ti,yi),i=0,,N}, with yiy(ti), ti=ih. Consider an approximation of =d/dt by linear combination of the sample values

(y)i1hk=0sakyi+k. (1)

Also introduce an approximation of y by linear combination of fi=f(t,yi)f(t,y(ti))

(y)ik=0sbkfi+k. (2)

Overall, the above approximations lead to the algorithm

1hk=0sakyi+k=k=0sbkfi+k. (3)

2Consistency

A first question that arises is whether the operator approximations converge. These are discrete approximations based upon sampling y(t) at ti=ih. The function y would be recovered in the limit h0, leading to the definition of consistency by

limh0(y)i=y'(ti),limh0(y)=f(ti,y(ti)),

stating that the exact operators are obtained.

Exercise 1. Carry out Taylor series expansion of the operator approximation (1) applied to the exact function y(t) to determine consistency conditions for coefficients ak.

Solution. Applied to the exact function, the operator approximation gives

y'(ti)(y(t))i=1hk=0saky(ti+k).

Taylor series expansion around ti gives

y(ti+k)=y(ti)+khy'(ti)+(kh)22!y''(ti)+.

Gathering terms,

(y(t))i=1h[(k=0sak)y(ti)+h(k=0skak)y'(ti)+h2(k=0sk2ak)y''(ti)++].

In order for limh0(y(t))i=y'(ti) the conditions

(k=0sak)=0,(k=0skak)=1,

must be met.

Exercise 2. Define a polynomial ρ(r)=(k=0sakrk). State the consistency conditions for as imposed values for ρ,ρ'.

Solution. The above consistency conditions correspond to

ρ(1)=0,ρ'(1)=1.

Exercise 3. Verify that the forward Euler, backward Euler, leapfrog schemes are consistent for f=0

Solution. The above schemes use operator series truncations

ddt1hΔ,ddt1h,ddt1hδ,

respectively. These all correspond to the same polynomial

ρ(r)=r-1,

after shifting indices. For example, in the leapfrog scheme

1h(yi+1/2-yi-1/2)=f(ti,yi)

shifting by 1.2 gives

1h(yi+1-yi)=f(ti+1/2,yi+1/2).

Since ρ(1)=0,ρ'(1)=1, the scheme is consistent.

Exercise 4. (1 Extra Credit point). Verify that the approximation in the first four Adams-Bashforth schemes are consistent.

Exercise 5. (1 Extra Credit point). Verify that the approximation in the first four Adams-Moulton schemes are consistent.

Exercise 6. Carry out Taylor series expansion of the operator approximation (2) applied to the exact function f(t,y(t)) to determine consistency conditions for coefficients bk.

Solution. Applied to the exact function, the operator approximation gives

f(ti,y(ti))(y(t))i=k=0sbkf(ti+k,y(ti+k)).

Taylor series expansion around ti gives

f(ti+k,y(ti+k))=f(ti,y(ti))+khft(ti,y(ti))+fy(ti,y(ti))(y(ti+k)-y(ti))+=

f(ti,y(ti))+khft(ti,y(ti))+fy(ti,y(ti))(y(ti)+khy'(ti)+-y(ti))+=

f(ti,y(ti))+kh[ft(ti,y(ti))+y'(ti)fy(ti,y(ti))]+𝒪(h2)

Gathering terms,

(y(t))i=[(k=0sbk)f(ti,y(ti))+𝒪(h)].

In order for limh0(y(t))i=f(ti,y(ti)) the condition

(k=0sbk)=1,

must be met.

Exercise 7. Define a polynomial σ(r)=(k=0sbkrk). State the consistency conditions for as imposed values for ρ,ρ'.

Solution. The above consistency condition corresponds to

σ(1)=1.

Exercise 8. State the consistency condition for (3) in terms of the ρ,σ polynomials

Solution. Gathering the above,

ρ(1)=0,ρ'(1)=1,σ(1)=1.

This can be slightly generalized to

ρ(1)=0,ρ'(1)=σ(1),

allowing a multiplicative factor in the evaluation of the derivative y'.

Exercise 9. Determine the consistency of the approximation for the forward, Euler, backward Euler, and leap frog schemes.

Solution. The σ polynomial for forward Euler is σ(r)=1, that for backward Euler is σ(r)=r. For leapfrog, rescale h2h, and obtain again σ(r)=r. All satisfy σ(1)=1.

Exercise 10. (1 Extra Credit point). Verify that the approximation in the first four Adams-Bashforth schemes are consistent.

Exercise 11. (1 Extra Credit point). Verify that the the approximation in the first four Adams-Moulton schemes are consistent.

Exercise 12. Apply the boundary locus method to determine the domain of linear stability for the forward, Euler, backward Euler, and leap frog schemes when f(t,y)=λy

Solution. The characteristic polynomial of the overall scheme is π(r,z)=ρ(r)-zσ(r), and the boundary locus method seeks roots of unit absolute value, r=eiθ.

Forward Euler

ρ(r)=r-1, σ(r)=1, π(r,z)=r-1-z=0z(θ)=eiθ-1, a circle centered at (-1,0)

n=90; h=2*π/n; θ=(0:n)*h; zFE=exp.(θ*im) .- 1;

Backward Euler

ρ(r)=r-1, σ(r)=r, π(r,z)=r-1-zr=0z(θ)=1-e-iθ, a circle centered at (1,0)

n=90; h=2*π/n; θ=(0:n)*h; zBE=1 .- exp.(θ*im);

Leapfrog

ρ(r)=r2-1, σ(r)=2r, π(r,z)=r2-1-zr=0z(θ)=eiθ-e-iθ=isinθ, a line segment on the imaginary axis

n=90; h=2*π/n; θ=(0:n)*h; zLF=im.*sin.(θ);
plot(real(zFE),imag(zFE),real(zBE),imag(zBE),real(zLF),imag(zLF));
axis("equal"); grid("on"); xlabel("Re(z)"); ylabel("Im(z)");
title("Boundary locii for forward Euler, backward Euler, leapfrog");
savefig("/home/student/courses/MATH661/images/H10Fig01.png");

Figure 1. Stability region for forward Euler is inside the circle, for backward Euler outside the circle, and for leapfrog within the slit on the imaginary axis.

Exercise 13. (2 Extra Credit points). Use the boundary locus method to find the linear stability region for the first four Adams-Bashforth schemes.

Exercise 14. (2 Extra Credit points). Use the boundary locus method to find the linear stability region for the first four Adams-Moulton schemes.