h=0.5;
f=inline('5*exp(2*x).*sin(x)-2*y+2*z','x','y','z');
x=0:0.5:20;
y1(x)=ones(1,length(x));
y1(1)=-2;
y2=y1;
y3=y1;
z1=ones(1,length(x));
z1(1)=-3;
z2=z1;
z3=z1;
for k=1:length(x)-1
z1(k+1)=z1(k)+h*f(x(k),y1(k),z1(k));
y1(k+1)=y1(k)+h*z1(k);
z2(k+1)=z2(k)+h/2*(f(x(k),y2(k),z2(k))+f(x(k+1),y1(k+1),z1(k+1)));
y2(k+1)=y2(k)+h/2*(f(x(k),y2(k),z2(k))+f(x(k+1),y1(k+1),z1(k+1)));
K1=f(x(k),y3(k),z3(k));
K2=f(x(k)+h/2,y3(k)+h/2*K1,z3(k)+h/2*K1);
K3=f(x(k)+h/2,y3(k)+h/2*K2,z3(k)+h/2*K2);
K4=f(x(k)+h,y3(k)+h*K3,z3+h*K3);
z3(k+1)=z3(k)+h*(K1+2*K2+2*K3+K4)/6;
K1=z3(k)+h*K1;
K2=z3(k)+h/2*K2;
K3=z3(k)+h/2*K3;
K4=z3(k)+h*K4;
y3(k+1)=y3(k)+h*(K1+2*K2+2*K3+K4)/6;
end
plot(x1,y1,'r*',x2,y2,'k.',x3,y3,'+')
legend('歐拉','改進歐拉','四階經典龍格庫塔');
hold on
下標索引必須為正整數型別或邏輯型別。
這個到底要怎么做呀 小白啥也不會TT TT
轉載請註明出處，本文鏈接：https://www.uj5u.com/houduan/24084.html
標籤：新手樂園
上一篇：C++矩陣相乘后結果特別離譜