求助一道matlab的题 ralston's 方法,a2=2/3a1=1/3, p1=q11=3/4y i+1=y i + [(1/3)k1+(2/3)k2]hk1=f(xi,yi)k2=f(xi+(3/4)h, yi+(3/4)k1h)

来源:学生作业帮助网 编辑:作业帮 时间:2024/05/02 01:24:10
求助一道matlab的题 ralston's 方法,a2=2/3a1=1/3, p1=q11=3/4y i+1=y i + [(1/3)k1+(2/3)k2]hk1=f(xi,yi)k2=f(xi+(3/4)h, yi+(3/4)k1h)

求助一道matlab的题 ralston's 方法,a2=2/3a1=1/3, p1=q11=3/4y i+1=y i + [(1/3)k1+(2/3)k2]hk1=f(xi,yi)k2=f(xi+(3/4)h, yi+(3/4)k1h)
求助一道matlab的题
ralston's 方法,
a2=2/3
a1=1/3, p1=q11=3/4
y i+1=y i + [(1/3)k1+(2/3)k2]h
k1=f(xi,yi)
k2=f(xi+(3/4)h, yi+(3/4)k1h)

求助一道matlab的题 ralston's 方法,a2=2/3a1=1/3, p1=q11=3/4y i+1=y i + [(1/3)k1+(2/3)k2]hk1=f(xi,yi)k2=f(xi+(3/4)h, yi+(3/4)k1h)
function [x,y]=Runge_kutta2(f,a,b,x0)
%2阶Runge_kutta解微分方程
%调用格式同ode45,
%f为微分方程函数,a b为积分区间,x0初值
%h为步长,默认为0.001
h=0.001;
xk=a:h:b;
n=(b-a)/h+1;
a=ones(length(x0),1);
y1(:,1)=x0';
for i=1:n
 x1=xk(i)*a;
k1=f(x1,y1(:,i));
k2=f(x1+(3/4)*h*a,y1(:,i)+(3/4)*k1*h);
y1(:,i+1)=y1(:,i)+((1/3)*k1+(2/3)*k2)*h;
y(i,:)=y1(:,i)';
end
x=xk;
--------------------------------------------------------------
odefun=@(t,x)[-10*x(1)*x(3)+x(2)
    10*x(1)*x(3)-x(2)
    -10*x(1)*x(3)+x(2)-2*x(3)];
[t,y]=ode45(odefun,[0 10],[50 0 40]);
[t1,y1]=Runge_kutta2(odefun,0,10,[50 0 40])
subplot(2,1,1),plot(t,y);legend('a-t','b-t','c-t');title('ode45')
subplot(2,1,2),plot(t1,y1);legend('a-t','b-t','c-t');title('Runge_kutta2')