范德堡方程的龙格库塔算法求解.docx

文档介绍：精品文档，仅供学****与交流，如有侵权请联系网站删除
【精品文档】第 1 页
龙格库塔算法：
function [x,y]=lgkt(f,x0,xf,y0,h)
h=;
n=fix((xf-x0)/h);
x(1)=x0;
y(:,1)=y0;
for i=1:n-1
x(i+1)=x(i)+h;
k1=f(x(i),y(:,i));
k2=f(x(i)+h/2,y(:,i)+h*k1/2);
k3=f(x(i)+h/2,y(:,i)+h*k2/2);
k4=f(x(i)+h,y(:,i)+h*k3);
y(:,i+1)=y(:,i)+h*(k1+2*k2+2*k3+k4)/6;
end
函数化简：
function f=f(t,y)
f=[y(2);(1-y(1).^2).*y(2)-y(1)];
Command Window的输入：
x0=input(' x0=');
xf=input(' xf=');
while x0>=xf
disp('输入不合逻辑，请重新输入')
x0=input(' x0=');
xf=input(' xf=');
end
y0=[1 8];
h=;
[t,F1]=lgkt(***@f,x0,xf,y0,h)
plot(t,F1)
grid
title('范德堡的龙格库塔算法之求解')
legend('y原函数，y的一阶导')