振动方程 数值计算 特征频率
对三自由度振动微分方程组采用线性加速度法进行数值计算，可以得到三个自由度的位移、速度及加速度的时程曲线，对三个自由度中某一自由度的加速度进行FFT变换时，从加速度的频响曲线中应该可以看到三个特征频率，却知道看两个，少了一个特征频率，请各位大侠指教。振动微分方程组的质量矩阵m=刚度矩阵为k=50*，理论计算结果得到的三个特征频率为0.5825Hz，1.5915Hz，2.1741Hz。数值计算程序如下。close all;
clc;
clear all;
m=; %质量矩阵
% c=; %阻尼矩阵
c=;
% k=50*;
k=50*; %刚度矩阵
x0='; %初位移
v0='; %初速度
delt=1/(12); %时间步长
fs=1/delt;
time=25; %仿真时间
n=time/delt; %循环次数
disp=zeros(n,3); %设定n行3列存储加速度矩阵
minv=inv(m+delt*c/2+delt^2*k/6);
i=1;
for t=0:delt:time
% f=';
f=';
if t==0
a0=inv(m)*(f-k*x0-c*v0); %计算初始加速度
else
a=minv*(f-c*(v0+delt/2*a0)-k*(x0+delt*v0+delt^2*a0/3));%计算加速度
v=v0+delt*(a0+a)/2; %计算速度
x=x0+delt*v0+delt^2/3*a0+delt^2/6*a; %计算位移
a0=a;v0=v;x0=x;i=i+1;
end
disp(i,:)=a0; %3个自由度的加速度矩阵
end
t=0:delt:time;
figure
% plot(t,disp(:,1),t,disp(:,2),t,disp(:,3));
plot(t,disp(:,1)); %三个自由度方向上加速度仿真曲线
grid on ;
xlabel('时间(s)'),title('3自由度时程曲线');
A=disp(:,1); %第一自由度方向的加速度向量
M=length(A);
n=0:M-1;
figure
plot((0:M-1)*fs/M,20*log10(abs(fft(disp(:,1)))));
grid;
仿真的图为
是自由响应？其它点的频谱呢？ 数值好像也不对啊 mni987654321 发表于 2017-6-5 09:08
是自由响应？其它点的频谱呢？
嗯，是自由响应，我这里把阻尼和外力都设为了0 Agoni 发表于 2017-6-5 10:48
数值好像也不对啊
大侠，哪里数值不对呀？ liuds123 发表于 2017-6-5 11:19
嗯，是自由响应，我这里把阻尼和外力都设为了0
阻尼为零的话，响应还能是自由响应吗？
把激励变为非零，阻尼要有，然后响应就正常了，上图是响应的fft,固有频率还是很清楚的 mni987654321 发表于 2017-6-6 09:05
把激励变为非零，阻尼要有，然后响应就正常了，上图是响应的fft,固有频率还是很清楚的
能不能把你改的程序发给我看看呀？ mni987654321 发表于 2017-6-6 09:05
把激励变为非零，阻尼要有，然后响应就正常了，上图是响应的fft,固有频率还是很清楚的
大侠，能不能看看你改过的程序呀？ liuds123 发表于 2017-6-6 09:27 static/image/common/back.gif
大侠，能不能看看你改过的程序呀？
很简单，把激振力改为非零，阻尼也改为非零，就行了，没有别的改动。就是改为冲击响应了。
mni987654321 发表于 2017-6-6 20:03
很简单，把激振力改为非零，阻尼也改为非零，就行了，没有别的改动。就是改为冲击响应了。
激振力是哪行代码？ Eminem 发表于 2017-6-7 14:36
激振力是哪行代码？
把箭头指的那行删掉，用它上面那行就可以了
mni987654321 发表于 2017-6-6 20:03
很简单，把激振力改为非零，阻尼也改为非零，就行了，没有别的改动。就是改为冲击响应了。
大侠，能不能把具体改的地方发给我学习一下，刚学，不怎么懂 你选定的自由度，在某阶模态上如果正好是“节点”，那么该自由度对应的响应中是不能反映这阶模态的~ mni987654321 发表于 2017-6-6 20:03
很简单，把激振力改为非零，阻尼也改为非零，就行了，没有别的改动。就是改为冲击响应了。
帅哥，冲击响应是怎么改的呀？在程序里改什么地方呀？能不能指导一下，一直不知道怎么弄
页:
[1]
2