0

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

首先，我想绘制 y、y_L 和 y_R 与时间的关系（对于 Omega 的给定值），然后我想绘制不同数量的 Omega 的域（y、y_L 和 y_R 的最大值）。不幸的是，我不擅长 Octave。我在 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')
4

1 回答 1

0

你有一个奇怪的术语，而不是在矢量定义的末尾

... + 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]，因此必须注意表达式不会被解释为矩阵行。操作符号的两个位置上的空格切换回单表达式解释。

于 2020-10-06T09:47:18.440 回答