Skip to content

Instantly share code, notes, and snippets.

@andyreagan

andyreagan/gist:7905737

Last active December 31, 2015 00:09
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save andyreagan/7905737 to your computer and use it in GitHub Desktop.
Save andyreagan/7905737 to your computer and use it in GitHub Desktop.
Download ZIP
just took the derivative of rk4 by hand...
Raw
gistfile1.matlab
L = L*(eye(dim) + h/6*(feval(J,t,y,params)+2*feval(J,t+h/2,y+h/2*s1,params)*(eye(dim)+h/2*feval(J,t,y,params))+2*feval(J,t+h/2,y+h/2*s2,params)*(eye(dim)+h/2*feval(J,t+h/2,y+h/2*s1,params))*(eye(dim)+h/2*feval(J,t,y,params))+feval(J,t+h,y+h*s3,params)*(eye(dim)+h*feval(J,t+h/2,y+h/2*s2,params))*(eye(dim)+h/2*feval(J,t+h/2,y+h/2*s1,params))*(eye(dim)+h/2*feval(J,t,y,params))));
%% BETTER
L1 = feval(J,t,y,params);
L2 = feval(J,t+h/2,y+h/2*s1,params)*(eye(dim)+h/2*L1);
L3 = feval(J,t+h/2,y+h/2*s2,params)*(eye(dim)+h/2*L2);
L4 = feval(J,t+h,y+h*s3,params)*(eye(dim)+h*L3);
L = L*(eye(dim) + h/6*(L1 + 2*L2 + 2*L3 + L4));
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment