第九章_動力學(xué)與振動講稿

1、第七章第七章 動力學(xué)與振動動力學(xué)與振動7.1軌跡軌跡7.2單自由度系統(tǒng)單自由度系統(tǒng)7.3多自由度系統(tǒng)多自由度系統(tǒng)7.1軌跡軌跡舉例說明舉例說明:重力場中有兩個物體重力場中有兩個物體,其中質(zhì)量為其中質(zhì)量為m2的物體固定的物體固定,而質(zhì)量為而質(zhì)量為m1的物體繞的物體繞m2做平面圓做平面圓周運動周運動.做圓周運動的做圓周運動的m1物體的軌道半徑用變物體的軌道半徑用變量量r表示表示,角度用變量角度用變量a表示表示. m1 r a m2 兩物體系統(tǒng)兩物體系統(tǒng)衛(wèi)星繞地球轉(zhuǎn)動時，衛(wèi)星繞地球轉(zhuǎn)動時，m2等于地球的質(zhì)量，等于地球的質(zhì)量，m1等于衛(wèi)星的質(zhì)量，等于衛(wèi)星的質(zhì)量，r為衛(wèi)星球心與地球為衛(wèi)星球心與地球球心間

2、的距離。其運動軌跡由下列方程組球心間的距離。其運動軌跡由下列方程組決定：決定：式中：式中： ，其中，其中t是時間變量，是時間變量，p為為物體在地球表面做圓周運動物體在地球表面做圓周運動的周期。在地球表面，r=6.373x106 m。0242222222ddaddrdadrrddardrdpt /用龍格用龍格庫塔法可以實現(xiàn)求解：庫塔法可以實現(xiàn)求解：引入新狀態(tài)變量：引入新狀態(tài)變量：帶入前面的微分方程組，可得四個一階微分方程。帶入前面的微分方程組，可得四個一階微分方程。ddaxaxddrxrx432114244321224122124xxxddxxddxxxxddxxddx0242222222dda

3、ddrdadrrddardrd建立函數(shù)文件建立函數(shù)文件orbit.mfunction xd=orbit(t,x)xd=x(2);x(1)*x(4)2-4.0*pi2/x(1)2; x(4);-2.0*x(2)*x(4)/x(1);三組初始條件三組初始條件(t=0)：組X1初始X2初始X3初始X4初始軌跡類型12001.5橢圓橢圓21002pi圓圓32004雙曲線雙曲線由初始條件建立執(zhí)行文件由初始條件建立執(zhí)行文件menu71.minitcond=2 0 0 1.5;1 0 0 2*pi;2 0 0 4;tspan=linspace(0,5,1000);options=odeset(RelTol,

4、1e-6,AbsTol,1e-6 1e-6 1e-6 1e-6);for i=1:3 t,x=ode45(orbit,tspan,initcond(i,:),options); polar(x(:,3),x(:,1); hold onendtext(0.5,-1.2,橢圓軌跡橢圓軌跡);text(-1.2,1,圓軌跡圓軌跡);text(1.75,2,雙曲線軌跡雙曲線軌跡);程序運行結(jié)果程序運行結(jié)果7.2單自由度系統(tǒng)單自由度系統(tǒng)7.2.1概述概述一一.力學(xué)模型力學(xué)模型 彈簧彈簧質(zhì)量質(zhì)量阻尼系統(tǒng)阻尼系統(tǒng)其中：振體質(zhì)量為其中：振體質(zhì)量為m，彈簧的線性系數(shù)為，彈簧的線性系數(shù)為k，非線，非線 性系數(shù)為性

5、系數(shù)為a，阻尼系數(shù)為，阻尼系數(shù)為c，外力，外力F（t）。）。mcK，aX（t）F（t）=X（0）kf（t）二二.運動微分方程運動微分方程用用x表示系統(tǒng)的位移，則運動微分方程為：表示系統(tǒng)的位移，則運動微分方程為：式中：式中： 固有頻率固有頻率 非線性系數(shù)非線性系數(shù) 阻尼因子阻尼因子引入新變量轉(zhuǎn)化狀態(tài)空間方程形式：引入新變量轉(zhuǎn)化狀態(tài)空間方程形式：)(20322fXxxddxdxdtnnmc2mkn2nmaddxxxx21)(203112221 fXxxxddxxddx7.2.2 線性系統(tǒng)的自由振動線性系統(tǒng)的自由振動一一.運動微分方程運動微分方程當(dāng)當(dāng) 時，得到線性振動系統(tǒng)的自由時，得到線性振動系統(tǒng)的

6、自由振動方程。振動方程。二二.MATLAB求解求解對應(yīng)的函數(shù)文件對應(yīng)的函數(shù)文件FreeOcillation.mfunction xdot=FreeOcillation(t,x,dummy,zeta)xdot=x(2);-2.0*zeta*x(2)-x(1);三種阻尼系數(shù)（三種阻尼系數(shù)（1）阻尼系數(shù)為）阻尼系數(shù)為0.1時是欠阻尼時是欠阻尼情況（情況（2）阻尼系數(shù)為）阻尼系數(shù)為1時是臨界阻尼情況時是臨界阻尼情況（3）阻尼系數(shù)為）阻尼系數(shù)為5時是過阻尼情況時是過阻尼情況0)( F0222xddxdxd 由初始條件（位移和速度均為由初始條件（位移和速度均為1時）建立執(zhí)行文件時）建立執(zhí)行文件menu72

7、.mzeta=0.1 1.0 5.0;tspan=linspace(0,40,400);%生成生成0-40的四百個線的四百個線性點性點lintype=-b -r -r;for i=1:3 t,x=ode45(FreeOcillation,tspan,1 1,zeta(i); subplot(2,1,1); plot(t,x(:,1),lintype(2*(i-1)+1:2*i); hold on subplot(2,1,2); plot(x(:,1),x(:,2),lintype(2*(i-1)+1:2*i); hold onendsubplot(2,1,1);xlabel(Time( tau

8、);ylabel(Displacement x( tau);title(Displacement as a function of( tau);axis(0 40 -2.0 2.0);text(2.7,-1.3,阻尼系數(shù)=0.1);text(3.6,-0.1,1.0);text(3.6,1.0,5.0);subplot(2,1,2);xlabel(Displacement);ylabel(Velocity);title(Phase portrait);axis(-2.0 2.0 -2.0 2.0);text(0.7,-1.25,阻尼系數(shù)=0.1);text(0.8,-0.65,1.0);tex

9、t(0.8,0.1,5.0);程序運行結(jié)果程序運行結(jié)果7.2.3 線性系統(tǒng)的強(qiáng)迫振動線性系統(tǒng)的強(qiáng)迫振動一一.運動微分方程運動微分方程二二.MATLAB求解求解若若對應(yīng)的函數(shù)文件對應(yīng)的函數(shù)文件ForceOcillation.mfunction xdot=ForceOcillation(t,x,dummy,zeta,Omega,x0)xdot=x(2);-2.0*zeta*x(2)-x(1)+x0*cos(Omega*t); )(2022 fXxddxdxd)cos()(0 XF為了獲得頻譜圖為了獲得頻譜圖,建立函數(shù)文件建立函數(shù)文件AmplitudeSpectrum.mfunctionf,ampl

10、itude=AmplitudeSpectrum(yy,Fs,Nstart,N);f=(Fs*(0:N-1)/N)*2.0*pi;amplitude=abs(fft(yy(Nstart:Nstart+N),N)/N;采樣速率采樣速率30/6000=0.005,則采樣頻率則采樣頻率1/0.005=200,這個頻率遠(yuǎn)遠(yuǎn)超出了必須達(dá)到的這個頻率遠(yuǎn)遠(yuǎn)超出了必須達(dá)到的采樣頻率采樣頻率,結(jié)果顯示截短頻譜圖結(jié)果顯示截短頻譜圖,需設(shè)置需設(shè)置Nstart=3200,N=211=2048。fft的應(yīng)用見的應(yīng)用見Help編制執(zhí)行文件編制執(zhí)行文件menu72f.mzeta=0.4;Omega=3.0;x0=50;tsp

11、an=linspace(0,30,6000);options=odeset(RelTol,1e-8,AbsTol,1e-8);lintype=-b; t,x=ode45(ForceOcillation,tspan,0 0,options,zeta,Omega,x0); subplot(2,1,1); plot(t,x(:,1); axis(0 30 -8 8); hold on subplot(2,1,2); function xdot=ForceOcillation(t,x,dummy,zeta,Omega,x0)xdot=x(2);-2.0*zeta*x(2)-x(1)+x0*cos(Om

12、ega*t);yy=x(:,1); N=2048;Nstart=3200;Fs=200;f,Amplitude=AmplitudeSpectrum(yy,Fs,Nstart,N); semilogy(f(1:40),2*Amplitude(1:40);xlabel(Frequency);ylabel(Amplitude);title(Response spectrum of a linear system); hold onsubplot(2,1,1);xlabel(Time( tau);ylabel(Displacement x( tau);title(Response of a linea

13、r system);hold on程序運行結(jié)果程序運行結(jié)果7.2.4 線性系統(tǒng)的頻率響應(yīng)、階躍響應(yīng)及線性系統(tǒng)的頻率響應(yīng)、階躍響應(yīng)及脈沖響應(yīng)脈沖響應(yīng)單自由度振動系統(tǒng)的強(qiáng)迫振動微分方程可為：單自由度振動系統(tǒng)的強(qiáng)迫振動微分方程可為：通過通過LAPLACE變換，得傳遞函數(shù)：變換，得傳遞函數(shù)：其中：其中： mknnmc2Logspace, rad2deg, loglog, semilogx see help一一.編制執(zhí)行文件編制執(zhí)行文件frequency72.m，求頻率響應(yīng)。，求頻率響應(yīng)。m = 1;zeta = 0.1:0.1:1;k = 1;wn = sqrt(k/m);10w = logspace

14、(-1,1,400);rad2deg = 180/pi;s = j*w;for cnt = 1:length(zeta)xfer(cnt,:)=(1/m) ./ (s.2 + 2*zeta(cnt)*wn*s + wn2);mag(cnt,:) = abs(xfer(cnt,:);phs(cnt,:) = angle(xfer(cnt,:)*rad2deg;endfor cnt = 1:length(zeta)figure(1)loglog(w,mag(cnt,:),k-)title(SDOF frequency response magnitudes for zeta = 0.2 to 1.

15、0 in steps of 0.2)xlabel(Frequency(rad/sec)ylabel(Magnitude)gridhold onendhold offfor cnt = 1:length(zeta)figure(2)semilogx(w,phs(cnt,:),k-)title(SDOF frequency response phases for zeta = 0.2 to 1.0 in steps of 0.2)xlabel(Frequency(rad/sec)ylabel(phase)gridhold onendhold off程序運行結(jié)果程序運行結(jié)果幅頻曲線幅頻曲線程序運行結(jié)果程序運行結(jié)果相頻曲線相頻曲線二二.求時頻響應(yīng)的基本函數(shù)命令求時頻響應(yīng)的基本函數(shù)命令可以通過上述命令求線性系統(tǒng)的波得圖、乃奎可以通過上述命令求線性系統(tǒng)的波得圖、乃奎斯特圖、斯特圖、階躍響應(yīng)、脈沖響應(yīng)、初始條件響應(yīng)、階躍響應(yīng)、脈沖響應(yīng)、初始條件響應(yīng)、輸入輸入u的響應(yīng)的響應(yīng)例：博得圖例：博得圖編制執(zhí)行文件編制執(zhí)行文件bode72.mm=1zeta=0.1k=1wn=sqrt(k/m)den=1 2*zeta*wn wn2num=1/msys=tf(num,den)bode(sys)Bode,nyquist see help例：求正弦輸入激勵響應(yīng)例：求正弦輸入