【正文】 ????????mmmM2000000 ，???????????????KKKKKKKK2202302 算出 K 的逆陣及系統(tǒng)的動(dòng)力矩陣為 ????????????22111111kK，???????????? ?5214212111kmMKA 若 ? ?TX 1111 ? ，第一次迭代后得到 ????????????87411 kmAXY ， ? ?????????????0 0 0 0 0 8 5 7 0 0 5 0 0 0 0 1312 YYX 重復(fù)上述步驟，各次的迭代結(jié)果列于表 41。為防止迭代過程中迭代列陣的元素變得過大或過小，每次迭代后需要使列陣歸一化，例如使它最后一個(gè)元素成為 1。 記 1X 為初始迭代列陣，由展開定理， 1X 可以表示為 111 ?aX ? nnaa ?? ??? ?22 （ 310） 對(duì)上式左乘矩陣 A，由式（ 39）得知第一次迭代后所得的列陣為 nnnaaaAXX ?????? ????? ?22211112 = ???????? ??? nnnaaa ????????1212111?（ 311） 如果特征值 1? 不是特征方程的重根，那么上式中的1n1312 ?????? 、 ? 都小于 1，因此比起其他主振型 1? 在 2X 內(nèi)占的比重相對(duì)地比在 1X 中占的比重大，換句話說，用矩陣 A 迭代計(jì)算一次后，擴(kuò)大了迭代列陣中第一階主振型的優(yōu)勢(shì)。 將 i???? ?? 和i 帶入 公式 中，得 iii ??? ?A （ 39） 若將上式左端看作新列陣，上式表示：對(duì)于精確的主振型。應(yīng)為 ??C 一般粘性阻尼矩陣，不一定滿足矩陣可對(duì)角化條件，為了把系統(tǒng)方程 （ 32） 轉(zhuǎn)化為一般特征值問題公式，需引入恒等式： ? ? ? ? ? ? ? ?( ) 0p M p M X?? （ 36） 將此 式與 （ 32） 合并： ? ? ? ? ? ? ? ?()p A B Y F ??? （ 37） 其中 遼寧工程技術(shù)大學(xué)畢業(yè)設(shè)計(jì)（論文） 19 ? ? ? ? ? ?? ? ? ?0 MAMC??? ????, ? ? ? ? ? ?? ? ? ?00MB K???? ????, ? ? ? ?? ?pXY X?????????, ? ? ? ?? ?0F F????? ??????? 令 ??F? =??0 , 則 （ 37） 的特征值滿足下列方程： ? ? ? ? 0p A B?? （ 38） 對(duì)于 N 自由度系統(tǒng)，此方程有 2N 個(gè)復(fù)共軛對(duì)出現(xiàn)的特征根： i i ii i ijj? ? ?? ? ???? ??? ?? ?? 其中 i? 阻尼因子； i? 為阻尼固有頻率。 式 （ 35）的分 母，叫做系統(tǒng)的特征方程。把這個(gè)時(shí)域矩陣方程變換到拉氏域（變數(shù)為 p ），并假定初始位移和初始速度為零，則得： ? ? ? ? ? ? ? ? ? ?2( ) ( ) ( )p M p C K X p F p? ? ? （ 32） 或 ? ?? ? ? ?( ) ( ) ( )Z p X p F p? （ 33） 式中 ? ?()Zp ：動(dòng)剛度矩陣。) 0 0 . 2 0 . 4 0 . 6 0 . 8 1 1 . 2 1 . 4 1 . 6 1 . 8 2 0 . 2 0 . 1 5 0 . 1 0 . 0 500 . 0 50 . 10 . 1 5tx(t)E x 3 . 1 1 本章小 結(jié) 基于 MATLAB 對(duì)單自由度自由振動(dòng)繪制振動(dòng)圖像，進(jìn)行粘性阻尼，強(qiáng)迫振動(dòng)振幅放大因子繪圖進(jìn)行數(shù)據(jù)分析，使振動(dòng)數(shù)據(jù)更加明顯。 title(39。x(t)39。)。 xlabel(39。 x(i)=x0_dot*sin(wn*t(i))/wn+(x0f_0/(wn^2w^2))*cos(wn*t(i)) +f_0/(wn^2w^2)*cos(w*t(i))。 f_0=F0/m。 遼寧工程技術(shù)大學(xué)畢業(yè)設(shè)計(jì)（論文） 17 x0=。 m=5。 系統(tǒng)全解形式如下： tftfxtxtx nnnnn ???????? c o sc o s)(s i n)( 22 022 000 ?????? ? 式中， sr a dsr a dmkmFf n /30,/20,20510000 ?????? ?? 利用 MATLAB 繪制解曲線上式的程序如下： % F0=100。 王超：基于 MATLAB 的振動(dòng)系統(tǒng)編程分析 16 0 0 . 5 1 1 . 5 2 2 . 5 3 3 . 5 4 4 . 5 500 . 511 . 522 . 53頻率比振幅放大因子0 . 0 50 . 1 00 . 1 50 . 2 50 . 3 7 50 . 5 01 . 01 . 1 4 1 算例利用 MATLAB，繪制彈簧 質(zhì)量系統(tǒng)在簡(jiǎn)諧力作用下的響應(yīng)曲線。 beta=sqrt((1+(2*kesai*lamda).^2)./((1lamda.^2).^2+(2*kesai*lamda).^2))。 plot(lamda,beta) hold on end axis([0 5 0 3])。 0 0 . 5 1 1 . 5 2 2 . 5 3 3 . 5 4 4 . 5 500 . 511 . 522 . 53頻率比振幅放大因子1 . 00 . 50 . 3 7 500 . 0 50 . 1 50 . 2 5 偏心質(zhì)量引起的強(qiáng)迫振動(dòng)振幅放大因子 2222)2()1( ??? ? ???meMB (227) 程序如下： for kesai=[,] lamda=0::。 beta=1./(sqrt((1lamda.^2).^2+(2*kesai*lamda).^2))。上述結(jié)果也可以由直接設(shè) tBx ?sin? 并代入下列方程而得到： tPkxxm ?sin0???? （ 224） 為了具體討論影響穩(wěn)定響應(yīng)的振幅和相位差的各種因素，記 kPB 00? （ 225） 0B 實(shí)際是質(zhì)量塊在激振力幅靜作用下的最大位移。 無阻尼系統(tǒng)對(duì)簡(jiǎn)諧激勵(lì)的穩(wěn)態(tài)響應(yīng)可以從式（ 226）得出。將式（ 221）代入（ 220），解得 ????? nn imPB 21220 ??? （ 222） 記 ? 為頻率比，它定義為 n???? （ 223） 則式（ 222）可以寫成 ??????? ii BeekpikPB ?? ??????? 222020 )2()1( 121 1 （ 224） 式中 2220 )2()1(1 ??? ??? kpB （ 225） 21 12 ???? ?? ?tg （ 226） 將式（ 224）代入（ 221），得到復(fù)數(shù)形式的特解為 )( ???? tiBex （ 227） 比較方程（ 217）與（ 218），可知（ 219）中的位移 x 是（ 220）中復(fù)數(shù) x 的虛部，因此（ 225）的虛部就是方程（ 212）的特解，即有 )sin( ?? ?? tBx （ 228） 其中 B 為振幅， ? 為相位差。 將方程（ 215）的兩端同 除以質(zhì)量 m ，并且令 22 nmc ??? （ 218） 其中 ? 為相對(duì)阻尼系數(shù)， n? 為相應(yīng)的無阻尼系統(tǒng)的固有頻率，則方程（ 215）成為 xm? P(t) kx xc? c k x m P(t) 王超：基于 MATLAB 的振動(dòng)系統(tǒng)編程分析 12 tmpxxx nn ???? s in2 02 ??? ??? （ 219） 上述方程特解可以通過 )sin( ?? ?? tBx 或者 tBtAx ?? s inc o s ?? 來求得，這里介紹用復(fù)數(shù)方法求式（ 219）的特解。從圖的受力分析，得到運(yùn)動(dòng)微分方程為： tpkxxcxm ?s in0??? ??? (216) 由常微分方程理論知道，方程 ()的通解 x 由相應(yīng)的齊次 方程的通解 hx 和非齊次方程的任意特解 px 兩部分組成，即 )()()( txtxtx ph ?? (217) 當(dāng)欠阻尼時(shí)，式中 )(txh 為有阻尼自由振動(dòng)，它的特點(diǎn)是振動(dòng)頻率為阻尼固有頻率，振幅按指數(shù)規(guī)律衰減，稱為瞬態(tài)振動(dòng)或瞬態(tài)響應(yīng)； )(txp 是一種持續(xù)的等幅振動(dòng)，它是由于簡(jiǎn)諧激勵(lì)振力的持續(xù)作用而產(chǎn)生的，稱之為穩(wěn)態(tài)強(qiáng)迫振動(dòng)或穩(wěn)態(tài)振動(dòng)，在間隔充分長(zhǎng)時(shí)間考慮的振動(dòng)就是這種穩(wěn)態(tài)振動(dòng)，而在剛受到外界激勵(lì)時(shí)，系統(tǒng)的響應(yīng)則是上述兩種振動(dòng)之和。圖所示的彈簧質(zhì)量系統(tǒng)中，質(zhì)量塊上作用有簡(jiǎn)諧激振力 tPtP ?sin)( 0? (215) 遼寧工程技術(shù)大學(xué)畢業(yè)設(shè)計(jì)（論文） 11 其中 0P 為激振力幅， ? 為激振頻率。)。 title(39。delta_s_t39。)。gtext(39。 hold on。w_n39。 end plot(t,w)。w(i)=(g/t(i))^。可以利用下列 MATLAB 程序畫出 st? 在 0~ 范圍內(nèi) n? 和 n? 的變換曲線： % g=。 在一定的 ? 之下，欠阻尼系統(tǒng)能夠更快地達(dá)到穩(wěn)態(tài)值；而過阻尼系統(tǒng)反應(yīng)遲飩，動(dòng)作緩慢，所以系統(tǒng)通常設(shè)計(jì)成欠阻尼系統(tǒng)， ? 取值為 2 遼寧工程技術(shù)大學(xué)畢業(yè)設(shè)計(jì)（論文） 9 0 2 4 6 8 10 12 14 16 18 0 . 200 . 20 . 40 . 60 . 811 . 21 . 4t ( s )x(t)a = 1 . 0a = 2 . 0a = 0 . 5 圖 22 算例繪制無阻尼單自由度系統(tǒng)的固有頻率和周期隨靜變形的變化曲線。 1?? 時(shí)，振蕩系統(tǒng)等同于兩個(gè)一階系統(tǒng)串聯(lián)。)。 figure(1),plot(t,y,39。wd=w0*sqrt(a*a1)。t=0::18。)。 figure(1),plot(t,y,39。x1=wd。 w0=1。hold on a=。r39。 wd=w0.*sqrt(1a*a)。 k=1。t=0::18。程序如下： 王超：基于 MATLAB 的振動(dòng)系統(tǒng)編程分析 8 clear,format pact。以靜平衡位置為原點(diǎn)建立如圖坐標(biāo)，由牛頓定律得運(yùn)動(dòng)方程為 [13]: 0??? kxxcxm ??? (21) 令 mkm n ?? 2,2 ? 其中 n 稱為衰減系數(shù)，單位為 s1 ； n? 是相應(yīng)的無阻尼時(shí)的固有頻率，式 (21)可以寫為 ： 02 2 ??? xxnx n???? (22) 如果進(jìn)一步令 nn??? (23) 其中無量綱的 ? 稱為相對(duì)阻尼系數(shù)，則式 (22)可寫為： 02 2 ??? xxx nn ??? ??? (24) 為了求解，令 x xc? kx x m m c m 遼寧工程技術(shù)大學(xué)畢業(yè)設(shè)計(jì)（論文） 7 stex? (26) 代入 (24)后得到特征方程： 02 22 ??? nnss ???