如何绘制带有非线性微分方程的振动系统的频率响应？

问题描述

我有一个用于3自由度振动系统的非线性微分方程组。 system of differential equations

首先，我想针对时间绘制y，y_L和y_R（对于给定的Omega值），然后我针对各种Omega绘制域（y，y_L和y_R的最大值）。 不幸的是，我不擅长八度音阶。我已经在Octave中编写了以下代码（基于其中一个用户给出的示例），但最终出现此错误：“匿名函数体必须是单个表达式”。

如果有人能帮助我，我将不胜感激。

代码如下：

Me = 4000;
me = 20;
c = 2000;
c1 = 700;
c2 = 700;
k = 20000;
k1 = 250000;
k2 = 20000;
a0 = 0.01;
om = 25;
mu1 = (c+2*c2)/(Me);
mu2 = (c2)/(Me);
mu3 = (c1+c2)/(me);
mu4 = (c2)/(me);
w12 = (2*k2)/(Me);
w22 = (k1+k2)/(me);
a1 = (k2)/(me);
a2 = (k)/(Me);
F0 = (k1*a0)/(Me);
couplode = @(t,y) [y(2); mu4*y(4) - mu3*y(2) - w22*y(1) + a1*y(3) + F0*cos(om*t); y(4); mu2*(y(2)+y(6)) - mu1*y(4) - w12*y(3) + 0.5*w12*(y(1)+y(5)) + a2((y(3)).^3; y(6); mu4*y(4) - mu3*y(6) - w22*y(5) + a1*y(3) + F0*cos(om*t)];

[t,y] = ode45(couplode,[0 0.49*pi],[1;1;1;1;1;1]*1E-8); 

figure(1)
plot(t,y)
grid
str = {'$$ \dot{y_L} $$','$$ y_L $$','$$ \dot{y} $$','$$ y $$','$$ \dot{y_R} $$','$$ y_R $$'};
legend(str,'Interpreter','latex','Location','NW')

解决方法

您在向量定义的末尾有一个奇怪的术语

... + a2((y(3)).^3

您当然是想

... + a2*y(3).^3

通过将其分成几行，您可以获得更好的可见性和更容易的调试

couplode = @(t,y) [ y(2); 
                    mu4*y(4)-mu3*y(2)-w22*y(1)+a1*y(3)+F0*cos(om*t); 
                    y(4); 
                    mu2*(y(2)+y(6)) - mu1*y(4) - w12*y(3) + 0.5*w12*(y(1)+y(5)) + a2*y(3).^3; 
                    y(6); 
                    mu4*y(4)-mu3*y(6)-w22*y(5)+a1*y(3)+F0*cos(om*t)];

至少以这种形式，空格或无空格都没有区别。通常，在matlab / octave中，[a +b -c]与[a,+b,-c]相同，因此必须注意不要将表达式解释为矩阵行。操作符两个位置上的空格切换回单表达式解释。