对于duffing方程 幅频响应 这一经典问题的求助
duffing方程 应该是很经典的问题了 试了很多方法 都不能很好得到其幅频响应 特指跳跃现象第一
用ezplot 直接或者扭曲的骨干曲线
第二
用 谐波平衡法求解当然用平均法 多尺度 都可以得到其关系式
但是再求解幅值时 利用的是盛金公式（见百度百科）
求出来的解跟实际情况有很大的差距
具体程序如下：
clc
clear all
n=0;
m=1;k=1;kn=0.0005;C=0.02;F=1;
for w=0.7:0.01:1.4
n=n+1;
a=9*kn^2/16;b=3*kn*(k-m*w^2)/2;c=(k-m*w^2)^2+C^2*w^2;d=-F^2;
X=b^2-3*a*c;
Y=b*c-9*a*d;
Z=c^2-3*b*d;
delt=Y^2-4*X*Z;
delt1(n)=Y^2-4*X*Z;
if X==0&&Y==0
A1(n)=-b/(3*a);
A2(n)=-c/(b);
A3(n)=-3*d/(c);
else if delt>0
W1=X*b+3*a*(-Y+abs(sqrt(Y^2-4*X*Z)))/2;
W2=X*b+3*a*(-Y-abs(sqrt(Y^2-4*X*Z)))/2;
%W1=X*b+3*a*(-Y+sqrt(Y^2-4*X*Z))/2;
% W2=X*b+3*a*(-Y-sqrt(Y^2-4*X*Z))/2;
A1(n)=(-b-(W1^(1/3)+W2^(1/3)))/(3*a);
A2(n)=(-b+1/2*(W1^(1/3)+W2^(1/3))+sqrt(3)/2*(W1^(1/3)-W2^(1/3))*i)/(3*a);
A3(n)=(-b+1/2*(W1^(1/3)+W2^(1/3))-sqrt(3)/2*(W1^(1/3)-W2^(1/3))*i)/(3*a);
else if delt==0
K=Y/X;
A1(n)=-b/a+K;
A2(n)=-K/2;
A3(n)=-K/2;
else if delt<0
T=(2*X*b-3*a*Y)/(2*sqrt(X^3));
sita=acos(T);
A1(n)=(-b-2*sqrt(X)*cos(sita/3))/(3*a);
A2(n)=(-b+sqrt(X)*(cos(sita/3)+sqrt(3)*sin(sita/3)))/(3*a);
A3(n)=(-b+sqrt(X)*(cos(sita/3)-sqrt(3)*sin(sita/3)))/(3*a);
end
end
end
end
end
A1=sqrt(A1);
A2=sqrt(A2);
A3=sqrt(A3);
w=0.7:0.01:1.4;
figure;plot(w,A1S,'+',w,A2S,'.',w,A3S,'*');
请高手指出其中的错误 谢谢
我记得以前直接matlab中的“solve”来求解三次方程，只保留正实根就可以了！ 还是觉得盛金公式太书面化。 lihaitao123 发表于 2014-4-18 15:36
我记得以前直接matlab中的“solve”来求解三次方程，只保留正实根就可以了！ 还是觉得盛金公式太书面化。
嗯 谢谢你 只留实根 选择自己需要的点 可以用伽辽金方法吗，我觉得这个方法很简单啊，为啥没人用 dum 发表于 2014-5-12 08:28
可以用伽辽金方法吗，我觉得这个方法很简单啊，为啥没人用
谢谢 用过最后我明白了 需要自己人为的去选择 实数根
页:
[1]