【Matlab代码】外力单自由度阻尼系统的数值解

function dx=odefun2(t,x)
dx=zeros(2,1);%存储x1,x2的导数
dx(1)=x(2);%第一个方程
dx(2)=-20*x(2)-100*x(1);%第二个方程

[t,x]=ode45('odefun2',[0,4],[1;0])

注：数值解仅仅是解的一个数据点列；数值解只能作定量分析，它能处理几乎所有的方程。

xt=dsolve('D2x+20*Dx+100*x=0','x(0)=1','Dx(0)=0')

给出数值解与解析解的对比图：

f=inline(xt)
plot(t,x(:,1),'r.',t,f(t),'b');
legend('数值解','解析解')

变量名=inline(’函数表达式’,‘变量名1’，‘变量名2’,…,‘变量名n’）

function dx=odefun3(t,x)
global c;
dx=zeros(2,1);%存储x1,x2的导数
dx(1)=x(2);%第一个方程
dx(2)=-20*c*x(2)-100*x(1)%第二个方程

注：通过全局变量c来接收外部阻尼系数的值。

function ode3(vc)
global c;%全局变量
hold on%hold住图形窗口
tspan=linspace(0,4,100);
for i=1:length(vc);
    c=vc(i);
    [t,x]=ode45('odefun3',tspan,[1,0]);
    text(t(10),x(10,1),['\leftarrow c',num2str(c)])
    plot(t,x(:,1))
end
hold off

vc=[0.2,0.4,0.6,0.8,1]
ode3(vc);

vc=[0.2,0.3,0.4,0.6,0.8,1,1.4,1.7];
ode3(vc);

微信“图像处理与模式识别研究所”关注我呦