汽車二自由度振動(dòng)模型計(jì)算

1、建立系統(tǒng)的動(dòng)力學(xué)方程建立系統(tǒng)的動(dòng)力學(xué)方程的方法：1、牛頓力學(xué)：牛頓第二定律；2、分析力學(xué)：拉格朗日方程。以x、為兩個(gè)變量建立二自由度系統(tǒng)動(dòng)力學(xué)方程；3、影響系數(shù)法-張量算法1、 根據(jù)牛頓第二定律， 2、 拉格朗日方程， -拉格朗日第二類方程，稱為拉氏函數(shù)，泛函；勢(shì)能函數(shù)：；動(dòng)能函數(shù)：；廣義坐標(biāo)對(duì)應(yīng)的非保守力： -保守系統(tǒng)的拉氏方程-利用上訴拉氏方程求解- ；將L=T-V代入拉氏方程可解的： 可見：與牛頓第二定律求得的結(jié)果一致。M+K=Q;對(duì)保守系統(tǒng)Q=0；M=()K=()q=()M+K=Q 和有明確的物理意義：彈性恢復(fù)力-Kq、慣性力-M與保守力Q平衡。如張量理論，可認(rèn)為是張量的坐標(biāo)，表示：使

2、系統(tǒng)僅產(chǎn)生沿坐標(biāo)的單位位移時(shí)，沿坐標(biāo)必須施加的外力，或者說Q的i分量在q的j分量上的影響量的投影。3、影響系數(shù)法-張量算法！下面由張量分量投影計(jì)算理論直接求解-！設(shè) m11 m12 * x” + k11 k12 * x m21 m22 ” + k21 k22 k11:僅當(dāng)x動(dòng)1單位時(shí)，x向的作用力為k1*1+k2*1=k1+k2;k12:僅當(dāng)動(dòng)1單位時(shí)，x向的作用力為-k1*(l1*1)+k2*(l2*1)=-k1*l1+k2*l2;k21:僅當(dāng)x動(dòng)1單位時(shí)，向的作用力為-(k1*1)*l1+(k2*1)*l2=-k1*l1+k2*l2;k22:僅當(dāng)動(dòng)1單位時(shí)，向的作用力為 k1*(l1*1)

3、*l1+k2*(l2*1)*l2=k1*l12+k2*l22;m11:僅當(dāng)x”動(dòng)1單位時(shí)，x”向的作用力為m*1;k12:僅當(dāng)”動(dòng)1單位時(shí)，x”向的作用力為0; %僅繞質(zhì)心轉(zhuǎn)動(dòng)時(shí)不影響x”向慣性力k21:僅當(dāng)x”動(dòng)1單位時(shí)，”向的作用力為0; %僅平動(dòng)時(shí)不影響”向慣性力？不過質(zhì)心時(shí)該怎么計(jì)算？若旋轉(zhuǎn)中心偏離質(zhì)心a,則變?yōu)閙a此時(shí)，k12: 僅當(dāng)”動(dòng)1單位時(shí)，x”向的作用力為m*(a*1);k21:僅當(dāng)x”動(dòng)1單位時(shí)，”向的作用力為m*1*a;k22:僅當(dāng)”動(dòng)1單位時(shí)，”向的作用力為 (J+m*a2)*1;k22:僅當(dāng)”動(dòng)1單位時(shí)，”向的作用力為 J*1;可見與上述結(jié)果一致。對(duì)建立的動(dòng)力學(xué)方程更

4、換坐標(biāo)求偏頻對(duì)上述系統(tǒng)建立前后輪縱向位移x1、x2的動(dòng)力學(xué)方程 :x1=x-l1*；x2=x+l2*。使用matlab的solve（x1=x-l1*；x2=x+l2*，x1, x2）可以直接解出: =(x1-x2)/(l1+l2);x=(x1*l2+l1*x2)/(l1+l2),代入前面創(chuàng)建的方程組： 消去x1和x2,可得到如下的方程：所以式中， 聯(lián)系系數(shù)，表示兩坐標(biāo)之間的聯(lián)系 偏頻，表示前后懸掛獨(dú)立振動(dòng)時(shí)的振動(dòng)頻率，即x1=0時(shí)的振動(dòng)頻率是w2,x2=0時(shí)的振動(dòng)頻率是w1,不同于系統(tǒng)的固有頻率（2自由度獨(dú)立時(shí)才相等）。 汽車?yán)@質(zhì)心軸的回轉(zhuǎn)半徑在汽車設(shè)計(jì)中，希望行車時(shí)一個(gè)懸掛的振動(dòng)不傳到另一個(gè)

5、懸掛上，為此，應(yīng)使車身質(zhì)量分布和前后輪位置滿足：質(zhì)量分配系數(shù)，這時(shí)對(duì)于一般質(zhì)量分配系數(shù)的耦合情況，可以用模態(tài)分析法求固有頻率及其通解：即可見，特征向量陣(模態(tài)矩陣)組成坐標(biāo)變換矩陣（由老基到新的主坐標(biāo)基的坐標(biāo)變換矩陣），xp=Tx為主坐標(biāo)，主振動(dòng)的坐標(biāo)，在該坐標(biāo)系上，各自由度獨(dú)立振動(dòng)（解耦）。分別是新的主坐標(biāo)基的兩個(gè)基矢量在老基下的投影坐標(biāo)。在新的主坐標(biāo)基下，Mp=TM為主質(zhì)量陣（主質(zhì)量組成的對(duì)角陣）；Kp=TK為主剛度陣（主剛度組成的對(duì)角陣）。這種矩陣變換的本質(zhì)是張量的坐標(biāo)變換。Mp+Kp=0 主坐標(biāo)方程組為解耦方程組。利用特征值分解找到系統(tǒng)的主坐標(biāo)基，通過坐標(biāo)變換進(jìn)行解耦、簡(jiǎn)化計(jì)算、然后再

6、變換回去，這是坐標(biāo)變換的意義所在。用matlab特征值分解法求平等與轉(zhuǎn)動(dòng)主模態(tài)（振型）%SH760小轎車空載主要參數(shù)m=1340;a=1.54;b=1.29;Ic=2395; %繞質(zhì)心的轉(zhuǎn)動(dòng)慣量k1=40000;k2=44000;M=m,0;0,Ic;K=k1+k2,-(k1*a-k2*b);-(k1*a-k2*b),k1*a2+k2*b2;eig_vec,eig_val = eig(inv(M)*K);omeg,w_order = sort(sqrt(diag(eig_val); %頻率mode_vec = eig_vec(:,w_order); %振型T=2.*pi./omeg; %周期m

7、ode_vec(:,1)=mode_vec(:,1)./mode_vec(1,1);mode_vec(:,2)=mode_vec(:,2)./mode_vec(1,2);subplot(2,1,1)plot(1;2,mode_vec(:,1)title(strcat(w1=,num2str(omeg(1);subplot(2,1,2)plot(1;2,mode_vec(:,2)title(strcat(w2=,num2str(omeg(2);因?yàn)閷?duì)特征值進(jìn)行了排序，所以w1Y=expm(A*t)*Y0 Y=x1;x2;x1;x2%expm(A)的意義是將坐標(biāo)先變換到主坐標(biāo)系，對(duì)對(duì)角值進(jìn)行exp

8、運(yùn)算后再變換到原坐標(biāo)系，如同張量坐標(biāo)變換help expmy0=x0;xd0; %四元變量的初始條件for i=1:round(tf/dt)+1 %設(shè)定計(jì)算點(diǎn)，作循環(huán)計(jì)算 tj(i)=dt*(i-1); y(:,i)=expm(A*tj(i)*y0; %循環(huán)計(jì)算矩陣指數(shù)endsubplot(2,2,4),plot(tj,y(1,:),y(2,:),grid可見，坐標(biāo)的選取對(duì)固有頻率沒有影響，但對(duì)振型有影響。W1_pianpin和w2_pianpin是前后的偏頻（假設(shè)質(zhì)量分配系數(shù)為1計(jì)算）。用matlab ode45()直接進(jìn)行仿真計(jì)算%SH760小轎車空載主要參數(shù)clear;m=1340;a=

9、1.54;b=1.29;l=a+b;Ic=2395; %繞質(zhì)心的轉(zhuǎn)動(dòng)慣量rou=sqrt(Ic/m);k1=40*1000;k2=44*1000;M=m*(b2+rou2)/l2,m*(a*b-rou2)/l2;m*(a*b-rou2)/l2,m*(a2+rou2)/l2;K=k1,0;0,k2;%用matlab ode45數(shù)值解-A=zeros(2,2),eye(2);-MK,zeros(2,2); %四階參數(shù)4X4矩陣X=AX-X=expm(A*t)*X0 X=x1;x2;x1;x2syms x1 x2 dx1 dx2df_sym=A*x1;x2;dx1;dx2;df_sym=subs(d

10、f_sym,x1,x2,dx1,dx2,x(1),x(2),x(3),x(4);n=length(df_sym);i=1;ss=; %先定義好很重要，否則再循環(huán)體中定義時(shí)，每一循環(huán)ss不累加。while in ss=strcat(ss,char(df_sym(i),;); i=i+1;endss=strcat(ss,char(df_sym(i);ss=strcat(ss,);f=inline(ss,t,x);t,x=ode45(f,0 10,1,0,0,0);%初始y=0,y=1%subplot(2,2,1)plot(t,x(:,1),x(:,2) %時(shí)間狀態(tài)系列用s-function進(jìn)行仿真

11、計(jì)算%sh760.mfunction sys,x0,str,ts=s_function(t,x,u,flag)switch flag,case 0, sys,x0,str,ts=mdlInitializeSizes;case 1, sys=mdlDerivatives(t,x,u);case 3, sys=mdlOutputs(t,x,u);case 2, 4, 9 sys = ;otherwise error(Unhandled flag = ,num2str(flag);endfunction sys,x0,str,ts=mdlInitializeSizessizes = simsizes

12、;sizes.NumContStates = 4;sizes.NumDiscStates = 0;sizes.NumOutputs = 2;sizes.NumInputs = 1;sizes.DirFeedthrough = 0;sizes.NumSampleTimes = 0;sys=simsizes(sizes);x0=1 0 0 0;str=;ts=;function sys=mdlDerivatives(t,x,u)m=1340;a=1.54;b=1.29;l=a+b;Ic=2395; %繞質(zhì)心的轉(zhuǎn)動(dòng)慣量rou=sqrt(Ic/m);k1=40*1000;k2=44*1000;M=m*(b2+rou2)/l2,m*(a*b-rou2)/l2;m*(a*b-rou2)/l2,m*(a2+rou2)/l2;K=k1,0;0,k2;A