2071 Lab 3 Summary
$Id: summary.txt,v 1.9 2016/12/22 21:03:19 mike Exp $
M. M. Sussman
{@(#) Mon Jan  4 19:21:39 2016 }

EXERCISE 1
1. copy stiff4_ode.m
2. copy stiff4_solution.m
3. make stiff55_ode.m, stiff50_solution.m, sim for 10000
3,4. plot with sine wave [0,2*pi] and direction field

EXERCISE 2
1. How many points make pleasing plot of stiff10000_solution?  40
2. solve: 40,000 steps is good for euler, rk3

EXERCISE 3
1-2. back_euler_lam.m + comments
3. nice, smooth solution with 40 steps
4. Implicit solution using 40 steps should be almost same line
   as explicit solution using 40000 steps.
5. numbers agree with those in lab.

EXERCISE 4
1. nsteps  back_euler_lam error   ratio
       40   1.2607e-04            2.8300   
       80   4.4547e-05            2.5269   
      160   1.7629e-05            2.3041   
      320   7.6513e-06            2.1647   
      640   3.5346e-06            2.0859 
     1280   1.6945e-06            2.0439
     2560   8.2903e-07

2. first order, p=1
3. nsteps  back_euler_lam error   euler
       40   1.2607e-04            3.8130e+33
       80   4.4547e-05            8.8607e+39
      160   1.7629e-05            3.6617e+08
      320   7.6513e-06           -5.3206e-06
      640   3.5346e-06           -2.9519e-06
     1280   1.6945e-06           -1.5488e-06
     2560   8.2903e-07           -7.9262e-07
4. same rate of convergence (ratios: 1.8025 1.9059 1.9540)

EXERCISE 5
1. change stiff10000_ode to include dfdy
2. change stiff55_ode
3,4 copy back_euler.m, newton4euler.m
5. add comments:
  back_euler.m
  % f_ode is the handle of a function whose signature is ???: BOTH or at least one
  % The following for loop performs the spatial stepping (on x): back_euler
  % The following statement computes the initial guess for Newton: back_euler
  newton4euler.m
  % f is a function whose signature is ???: BOTH or at least one
  % TOL = ??? (explain in words what TOL is used for): newton4euler
  % MAXITS = ??? (explain in words what MAXITS is used for): newton4euler
  % When F is a vector and D a matrix, the expression J/F maans ???
  % The Matlab function "eye" is used for ???: newton4euler
6. F=[-6    J=[ 4  4
      -3]      12  5];
7. variables are J, increment, F
8. comment on implementing Equation 5 (Newton update)
9. comment on implementing Equation 6 (def. F(Y) )
10. if length(Y)>1, is Y column vector.
11. F=yk+h*(LAMBDA*(-Y+sin(xkp1))-Y=0, yk=1,h=.1,Y=1,xkp1=2, 
    Y=0.909388038786895 = (1+10000*.1*sin(2))/(1+10000*.1)
12. newton4euler agrees exactly
13. solutions agree to all places, norm(inf) diff=1.11022e-16

EXERCISE 6
1. vanderpol_ode.m, fill in:
   yprime = [y(2); -a*(y(1)^2-1)*y(2)-y(1)+exp(-x)];
   df1dy2 = 1;
   df2dy1 = -2*a*y(1)*y(2)-1;
2. plot (req'd): back_euler is smooth, forward_euler is not.
3. plot (req'd): 640 points: two curves close.

EXERCISE 7
1-3. trapezoid.m, copied new lines
4-7. newton4trapezoid.m, comments, new method, MAXITS=30
8. 0.894104843976271
9. solve:
   Y = ( yk+lambda*h/2(-yk+sin(xk)+sin(xkp1)) )/( 1+LAMBDA*h/2 )
     = (1.0 +5.5/2*(-1+sin(1.9)+sin(2.0)))/(1+5.5/2) = 0.894104843976271
10. numSteps   trapezoid error         ratio
     10         5.2779e-03         217.98972
     20         2.4212e-05          -0.64668
     40        -3.7440e-05           4.00740
     80        -9.3427e-06           4.00185
     160       -2.3346e-06           4.00046
     320       -5.8358e-07

6. second order

EXERCISE 8
1. plots of 4 trapezoid vanderpol_ode solutions.  Mostly together
2. plots of 4 back_euler vanderpol_ode solutions.  Marching to convgce.
3. plots of 800-step trapezoid and 12800-step back_euler are similar.

EXERCISE 9
1. plot: sine wave and +/- jigglies around it.
2. plots 200,400,800,1600.  Only 1600 is required: "fur" initially.

EXERCISE 10
1. plot ode45('vanderpol_ode',[0,70],[0;0])
2. plots using ode45 and ode15s not required, ode15s less dense in flat areas.
3. no. steps=256,                  y_1 at 70=-1.128238369949555
4. no. steps with RelTol=1.e-8 is 2102.  y_1=-1.125147778949563

EXERCISE 11 (Extra Credit)
1. bdm2.m
2. stiff55: err=  -1.4746e-04 -3.7187e-05 -9.3219e-06 -2.3327e-06
            rat=  3.9654 3.9892 3.9962
3. plot: close
4. plot shows many bad overshoots.