博大广阔 发表于 2011-10-30 21:53

非线性系统的相平面的绘制

%要求的一阶方程组的形式
%dx=f(t,x,y);    %位移的导数,x代表位移,y代表速度
%dy=g(t,x,y);    %速度的导数
%在末尾出设置
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%二阶常微分方程组的求解——画相位图
function Runge_Kutta()
clc
clear all

%算法参数预设
h=0.1;Tmax=0.6;                                    %步长为h,计算时间为20s,n用于计数
xold=0;yold=2;                                     %设置初始值
%四阶龙格算法
for t=0:h:Tmax
t=0;
m1=f(t,xold,yold)*h;k1=g(t,xold,yold)*h;         %将h放在这里减少计算舍去误差
m2=f(t+0.5*h,xold+0.5*m1,yold+0.5*k1)*h;k2=g(t+0.5*h,xold+0.5*m1,yold+0.5*k1)*h;
m3=f(t+0.5*h,xold+0.5*m2,yold+0.5*k2)*h;k3=g(t+0.5*h,xold+0.5*m2,yold+0.5*k2)*h;
   m4=f(t+h,xold+m3,yold+k3)*h;k4=g(t+h,xold+m3,yold+k3)*h;

%保存
if t==0; n=1;end
TT(n,:)=t;XX(n,:)=xold;YY(n,:)=yold;n=n+1;
%更新
xold=xold+(m1+2*m2+2*m3+m4)/6;
yold=yold+(k1+2*k2+2*k3+k4)/6;
t=t+h;
end
%计算结果输出
%plot(TT,XX,'b');grid on;
plot(XX,YY,'r');grid on
end


%要求的位移的导数表达式dx=f(t,x,y);
function dx=f(t,x,y)
dx=2*x+4*y;               %不同的系统,需要修改
end
%要求的速度的导数表达式dy=g(t,x,y);
function dy=g(t,x,y)
%t是时间,x是位移,或者速度的导数,y是速度,或者速度的导数
dy=-x+6*y;                     %不同的系统,需要修改
end




dongyuanzhixing 发表于 2011-10-31 14:17

除了有点乱 ,其他的还 好哦

博大广阔 发表于 2011-10-31 18:07

仅三个函数。。。不乱吧
页: [1]
查看完整版本: 非线性系统的相平面的绘制