如何对模型施加白噪声激励

secondye

本人用 matlab写了一个如图所示的桁架的自由振动的程序。
现在想利用白噪声对其中一两个节点进行激励，得到各节点的反应时程数据，进而进行模态参数识别。
看到有朋友说用newmark法可以。
不知道各位有没有做过呢？

VibrationMaster

help randn

secondye

回复 VibrationMaster 的帖子

这个只是产生白噪声的吧？
如何利用产生的白噪声对模型进行激励呢？

VibrationMaster

secondye 发表于 2010-11-4 23:46

                               
登录/注册后可看大图

回复 VibrationMaster 的帖子

这个只是产生白噪声的吧？

把产生的白噪声当作已经离散的确定性激励计算.(注意时间步长, journal of sound &viibration上有一长文)
重复N次,取平均,计算方差

secondye

回复 VibrationMaster 的帖子

敢问 VibrationMaster，那篇文章问叫什么名字呢？
我下载下来读读
谢谢了

VibrationMaster

http://www.sciencedirect.com/science?_ob=ArticleURL&_udi=B6WM3-466TCF1-H&_user=1479053&_coverDate=07%2F18%2F2002&_alid=1528087241&_rdoc=7&_fmt=high&_orig=search&_origin=search&_zone=rslt_list_item&_cdi=6923&_sort=r&_st=13&_docanchor=&view=c&_ct=8&_acct=C000053031&_version=1&_urlVersion=0&_userid=1479053&md5=076abf6e516cc043f532a1aba12fb729&searchtype=a

secondye

回复 6 # VibrationMaster 的帖子

又有问题了，VibrationMaster，请教：
上次问您的，产生了白噪声作为激励荷载，然后我再用nermark-beta法计算了该桁架模型的各自由度的加速度响应。
我得到模型的刚度和质量直接K，M后，再构造阻尼矩阵。然后我就是利用 wave = wgn(20000,1,0);    产生了白噪声荷载，然后就直接利用Newmark-beta法了，得到响应信号，再做FFT，得到的FFT频谱图的峰值与对K,M做特征值分解eig(K,M)所得到的频率值能对得上。
但是有点疑问啊，就是白噪声激励是激励在哪个自由度上呢？还是全部自由度？
如果我选择激励其中某一个自由度，该怎么样选择呢？

这是对某一自由度的响应信号做FFT得到的频谱图：做eig(K,M)得到的频率值为：
10.67，32.60，42.359，63.052，93.495，114.548
从频谱图上得到的值为：
10.68，32.48，42.12，62.27，90.96，110.00
频谱图上得到的频率值比做特征值分析的要小，不知道为什么。

secondye

回复 6 # VibrationMaster 的帖子

我把代码贴上来，请指教一下，感谢：

Ndof=37;
[V,D] = eig(K,M);
d = sqrt(diag(D));
[ff,JJ]=sort(d/2/pi);
C = zeros( size( K ) ) ; % 无阻尼
% % 比例阻尼
% w1 = d(1);
% w2 = d(2);
% AA = 0.5*[1/w1 w1; 1/w2 w2];
% BB = [0.05 0.05]';     %% Assume that arbitary modal damping ratio is 5%
%
% CC = inv(AA)*BB;
%
% a1 = CC(1);
% a2 = CC(2);
%
% C  = a1*M + a2*K;
%%
% %------ Step 2. 施加地震荷载或白噪声激励   
    load quake;      % help quake (e,n,v)         
%    wave = 0.098*e; % 地震荷载
    wave = wgn(20000,1,0);    % Guass white noise,白噪声荷载
%%   
% %------ Step 3. Newmark-beta for time history response   
    % step3.1 初始计算--定义参数  
    gDeltaT = 0.001 ;
    gTimeEnd = 10000*gDeltaT  ;    % 计算时间为载荷通过所需时间的两倍
    timestep = floor(gTimeEnd/gDeltaT) ;
    gama = 0.5 ;
    beta = 0.25 ;
   
    % step3.2 Displacement,Velocity,Acceleration:定义位移，速度和加速度
    gDisp = zeros( Ndof, timestep ) ;
    gVelo = zeros( Ndof, timestep ) ;
    gAcce = zeros( Ndof, timestep ) ;
   
    % Initialization Conditions:初始条件
    gDisp(:,1) = zeros( Ndof, 1 ) ;
    gVelo(:,1) = ones( Ndof, 1 ) ;
   
    [N,N] = size( K ) ;
    alpha0 = 1/beta/gDeltaT^2 ;
    alpha1 = gama/beta/gDeltaT ;
    alpha2 = 1/beta/gDeltaT ;
    alpha3 = 1/2/beta - 1 ;
    alpha4 = gama/beta - 1 ;
    alpha5 = gDeltaT/2*(gama/beta-2) ;
    alpha6 = gDeltaT*(1-gama) ;
    alpha7 = gama*gDeltaT ;
    K1 = K + alpha0*M + alpha1*C;
    timestep = floor(gTimeEnd/gDeltaT) ;
   
    % step3.3 计算初始加速度   
    gAcce(:,1) = M\(-K*gDisp(:,1)-C*gVelo(:,1)) ;
    p01=-diag(M);
    % step3.4 迭代计算：对每一个时间步计算
    for i=2:1:timestep
        if mod(i,100) == 0
            fprintf( '当前时间步：%d\n', i ) ;
        end
        pi1 = p01*wave(i);
        f1 = pi1+M*(alpha0*gDisp(:,i-1)+alpha2*gVelo(:,i-1)+alpha3*gAcce(:,i-1)) ...
                  + C*(alpha1*gDisp(:,i-1)+alpha4*gVelo(:,i-1)+alpha5*gAcce(:,i-1)) ;
        
        gDisp(:,i) = inv(K1)*f1;
        gAcce(:,i) = alpha0*(gDisp(:,i)-gDisp(:,i-1)) - alpha2*gVelo(:,i-1) - alpha3*gAcce(:,i-1) ;
        gVelo(:,i) = gVelo(:,i-1) + alpha6*gAcce(:,i-1) + alpha7*gAcce(:,i) ;      
    end

VibrationMaster

1. 确切地说，激励作用在空间某位置，而非自由度
2. 振动是整体现象，一般来说作用在任何一点，系统各处都有响应
3. 只要激励合适，任何一点的响应都包含了整体参数，即频率和阻尼比的信息

dujunmin

首先，我觉得你用两种方法获得的结构模态是非常近似的，是在误差允许范围之内的。
关于你提出的问题，我觉得楼上回答得非常到位。只我想补充：从理论上讲，激励不论在那个方向，都能激起各点的响应，但是对一些特殊的结构，还是要考虑激励方向，比如一个梁，沿轴向激励，则弯曲模态就不能不激起来，因为这两个广义自由度上的传递函数为0；此外，在实际模态实验时，激振力的方向也很重要，要尽量选择你希望得到的模态的方向进行激励（如你希望得到x方向的弯曲模态，就激励x方向），因为其它方向上的响应可能会很小，采集的信号信噪比降低，从而导致不易被识别。如果你要得到所有的模态，那么可以多做几次实验，分别激励不同的方向，然后分别识别。

wang185396

你好  我想请教一下如何才获取白噪声激励， 下面是我在matlab里边的一个程序，得出来一个白噪声信号的时间变化图  和一个自功率谱。可是再往下 不知道该则么整了  请指点一下吧  谢谢
sf=input('sf=')
fi=input('fi=')
fa=input('fa=')
tl=input('tl=')
am=input('am=')


n=round(sf*tl)+1;
t=0:1/sf:(n-1)/sf;
nfft=2^nextpow2(n);
ni=round(fi*nfft/sf+1);
na=round(fa*nfft/sf+1);
r=zeros(1,nfft/2);
r(ni:na)=1;
a=2*pi*rand(1,nfft/2);
b=r.*exp(i*a)*nfft/2;
a=[b,b(nfft/2:-1:1)];
x=real(ifft(a,nfft));
y=am*x(1:n)/max(abs(x(1:n)));
nft=1024;
y1=psd(y,nft,sf,boxcar(nft),nft/2);
f=0: sf/nft: (nft/2-1) *sf/nft;
subplot(2,1,1);
plot(t,y);
xlabel('时间（s）');
ylabel('幅值');
grid on;

subplot(2,1,2);
plot(f(1:nft/4),y1(1:nft/4));
xlabel('频率(Hz)');
ylabel('幅值');
grid on;

secondye

回复 11 # wang185396 的帖子

这个是王济书上的程序吧，那个“y”就是噪声信号啊

欧阳中华

.
   secondye 你的识别研究的如何？有结果贴上来欣赏一下...

wang185396

回复 12 # secondye 的帖子

    但是我想要的是一个从20HZ-200Hz之间的  一个随频率变化的白噪声激励，   哎 搞不懂啊   还请赐教一下吧!

wang185396

回复 12 # secondye 的帖子

我就是想要一个 单位白噪声激励，激励频率为【20-200】   ，然后拿它去模拟  发动机的随机激励  ，求取车身的振动速度。我想这应该按照频响分析来计算吧，而不是瞬态响应。所以 应该获取一个按频率变化的白噪声激励吧。