rkf45.m

来自「自己编写的变步长龙格库塔方法解强非线性微分方程」· M 代码 · 共 11 行

function [y1,y2]=rkf45(t,y,h)
k1=myfun(t,y);
k2=myfun(t+h/4,y+h/4*k1);
k3=myfun(t+3/8*h,y+3/32*k1+9/32*k2);
k4=myfun(t+12/13*h,y+1932/2197*k1-7200/2197*k2+7296/2197*k3);
k5=myfun(t+h,y+439/216*k1-8*k2+3680/513*k3-845/4104*k4);
k6=myfun(t+h/2,y-8/27*k1+2*k2-3544/2565*k3+1859/4104*k4-11/40*k5);
y1=y+h*(25/216*k1+1408/2565*k3+2197/4101*k4-1/5*k5);
y2=y+h*(16/135*k1+6656/12825*k3+28561/56430*k4-9/50*k5+2/55*k6);