Eugene Johnson - Lab 10
% This function first calls RK4 method to get first three
% values of t, this is a jump start because the way the AB/AM
% methods work.
function XX = ode_solver(Xinit,h,f,tfinal);
XX = zeros(5,(tfinal/h)+1);
XX(:,1) = Xinit;
% getting next three from rk4 method
for i = 1:3,
XX(:,i+1) = rk4(XX(:,i),h,f);
% Calling ab_am function to get predictor/corrector.
for i = 4:(tfinal/h),
XX(:,i+1) = ab_am(XX(:,(i-3):i),h,f);
function M = make_movie(XX)
% make_movie takes as input the XX matrix which is the output of
% ode_solver. This function should pick about twenty equally-spaced
% columns of XX and make a movie of them.
% clear the figure window and initialize M to be big enough to hold
% fifty frames. If MATLAB whines about running out of memory at this
% point, make your figure window smaller.
% You really don't get the true picture because the laws of physics
% the physics law is Kepler's law and this causes the shuttle to
% zip around the earth so fast that it looks like it passes through
% the earth. I tryed to make it look more reallistic. I got to where
% I wass take 200 pictures and it didn't change the movie so I
% decided that it was best to write a little comment and leave it
% at 50 clips. So here we go.
clf reset; M = moviein(50);
% xpos and ypos should end up being 1-by-50 (or smaller) row vectors
% picked at equal intervals from the second and third row of XX
% ------ finish these lines ------ %
xpos = XX(2,1:(6200/50):6201);
ypos = XX(3,1:(6200/50):6201);
% determine how to fix the axis limits so the movie doesn't wiggle around
% as it plays
% loop through the values of xpos and ypos. Plot the earth as a white 'o',
% the moon as a white 'x', and the shuttle as a white '+'.
u = .0121285;
for i = 1:length(xpos),
axis([min(xpos)-.5 max(xpos)+.5 min(ypos) max(ypos)]); axis equal
M(:,i) = getframe;
% This function is used to jump start everything. It gets the first
% three values of t
function Z = rk4(X,h,f);
f1 = func(X,f);
f2 = func(X+h/2*f1,f);
f3 = func(X+h/2*f2,f);
f4 = func(X+h*f3,f);
% getting the values of Z for the whole thing to meet
Z = X + h/6*(f1+2*f2+2*f3+f4);
% This function uses the AB/AM methods (predictor/corrector) to
% evaluate the values of the differential equations.
% implementing the predictor/corrector methods
function Z = ab_am(XX,h,f);
% Temp is a temperary variable which holds the predicted value
% Getting predicted value.
TEMP = XX(:,4) + h/24*(55*func(XX(:,4),f) - 59*func(XX(:,3),f) + ...
37*func(XX(:,2),f) - 9*func(XX(:,1),f));
% Correcting the predicted value
Z = XX(:,4) + h/24*(9*func(TEMP,f) + 19*func(XX(:,4),f) - ...
5*func(XX(:,3),f) + func(XX(:,2),f));