
matlab求常微分方程的,数值逼近。
2.对于一阶微分方程初值问题,用Euler法,改进Euler法,二阶R-K方法求解,并作图比较。复制的时候没有方程...
2. 对于一阶微分方程初值问题 ,用Euler法,改进Euler法,二阶R-K方法求解,并作图比较。
复制的时候没有方程 展开
复制的时候没有方程 展开
1个回答
展开全部
f=inline('cos(x)+sin(y)','x','y'); %微分方程的右边项
dx=0.05; %x方向步长
xleft=pi/2; %区域的左边界
xright=3*pi/2; %区域的右边界
xx=xleft:dx:xright; %一系列离散的点
n=length(xx); %点的个数
y0=0;
%%(1)欧拉法
Euler=y0;
for i=2:n
Euler(i)=Euler(i-1)+dx*f(xx(i-1),Euler(i-1));
end
%%(2)改进欧拉法
MEuler=y0;
for i=2:n
MEuler(i)=MEuler(i-1)+dx/2*(f(xx(i-1),MEuler(i-1))+f(xx(i),MEuler(i-1)+dx*f(xx(i-1),MEuler(i-1))));
end
%%(3)龙格库塔法
RK=y0;
for i=2:n
k1=f(xx(i-1),RK(i-1));
k2=f(xx(i-1)+dx/2,RK(i-1)+k1*dx/2);
k3=f(xx(i-1)+dx/2,RK(i-1)+k2*dx/2);
k4=f(xx(i-1)+dx,RK(i-1)+k3*dx);
RK(i)=RK(i-1)+dx*(k1+2*k2+2*k3+k4)/6;
end
%%Euler和MEuler,RK法结果比较
plot(xx,Euler,xx,MEuler,xx,RK)
hold on
legend('Euler','MEuler','Runge-Kutta')
dx=0.05; %x方向步长
xleft=pi/2; %区域的左边界
xright=3*pi/2; %区域的右边界
xx=xleft:dx:xright; %一系列离散的点
n=length(xx); %点的个数
y0=0;
%%(1)欧拉法
Euler=y0;
for i=2:n
Euler(i)=Euler(i-1)+dx*f(xx(i-1),Euler(i-1));
end
%%(2)改进欧拉法
MEuler=y0;
for i=2:n
MEuler(i)=MEuler(i-1)+dx/2*(f(xx(i-1),MEuler(i-1))+f(xx(i),MEuler(i-1)+dx*f(xx(i-1),MEuler(i-1))));
end
%%(3)龙格库塔法
RK=y0;
for i=2:n
k1=f(xx(i-1),RK(i-1));
k2=f(xx(i-1)+dx/2,RK(i-1)+k1*dx/2);
k3=f(xx(i-1)+dx/2,RK(i-1)+k2*dx/2);
k4=f(xx(i-1)+dx,RK(i-1)+k3*dx);
RK(i)=RK(i-1)+dx*(k1+2*k2+2*k3+k4)/6;
end
%%Euler和MEuler,RK法结果比较
plot(xx,Euler,xx,MEuler,xx,RK)
hold on
legend('Euler','MEuler','Runge-Kutta')
推荐律师服务:
若未解决您的问题,请您详细描述您的问题,通过百度律临进行免费专业咨询