matlab 汽车振动,汽车振动matlab编程问题

clear all

m=1235;Jx=472;Jy=1731;v=801000/60/60;

csp=4e-7;

mt1=49;mt2=49;mt3=49;mt4=49;

M=[m 0 0 0 0 0 0;0 Jx 0 0 0 0 0;0 0 Jy 0 0 0 0;0 0 0 mt1 0 0 0;0 0 0 0 mt2 0 0;0 0 0 0 0 mt3 0;0 0 0 0 0 0 mt4]

k1=22741;k2=22741;k3=26144;k4=26144;

kt1=302342;kt2=302342;kt3=492982;kt4=492982;

c1=1228;c2=1228;c3=1210;c4=1210

ct1=430;ct2=430;ct3=430;ct4=430

l=2.655;%轴距

l1=1.1135;l2=1.5415;l3=0.7075;l4=0.7075;%l1,l2 前轴距与后轴距，l3,l4：1/2 轮距

C=[c1+c2+c3+c4 -l1*c1+l2*c2-l1*c3+l2*c4 -l3*c1-l3*c2+l4*c3+l4*c4 -c1 -c2 -c3 -c4;...

-l1*c1+l2*c2-l1*c3+l2*c4 l1^2*c1+l2^2*c2+l1^2*c3+l2^2*c4 l1*l3*c1-l2*l3*c2-l1*l4*c3+l2*l4*c4 l1*c1 -l2*c2 l1*c3 -l2*c4;...

-l3*c1-l3*c2+l4*c3+l4*c4 l1*l3*c1-l2*l3*c2-l1*l4*c3+l2*l4*c4 l3^2*c1+l3^2*c2+l4^2*c3+l4^2*c4 l3*c1 l3*c2 -l4*c3 -l4*c4;...

-c1 l1*c1 l3*c1 c1+ct1 0 0 0;...

-c2 -l2*c2 l3*c2 0 c2+ct2 0 0;...

-c3 l1*c3 -l4*c3 0 0 c3+ct3 0;...

-c4 -l2*c4 -l4*c4 0 0 0 c4+ct4];

K=[k1+k2+k3+k4 -l1*k1+l2*k2-l1*k3+l2*k4 -l3*k1-l3*k2+l4*k3+l4*k4 -k1 -k2 -k3 -k4;...

-l1*k1+l2*k2-l1*k3+l2*k4 l1^2*k1+l2^2*k2+l1^2*k3+l2^2*k4 l1*l3*k1-l2*l3*k2-l1*l4*k3+l2*l4*k4 l1*k1 -l2*k2 l1*k3 -l2*k4;...

-l3*k1-l3*k2+l4*k3+l4*k4 l1*l3*k1-l2*l3*k2-l1*l4*k3+l2*l4*k4 l3^2*k1+l3^2*k2+l4^2*k3+l4^2*k4 l3*k1 l3*k2 -l4*k3 -l4*k4;...

-k1 l1*k1 l3*k1 k1+kt1 0 0 0;...

-k2 -l2*k2 l3*k2 0 k2+kt2 0 0;...

-k3 l1*k3 -l4*k3 0 0 k3+kt3 0;...

-k4 -l2*k4 -l4*k4 0 0 0 k4+kt4];

Ct=[0 0 0 0;0 0 0 0;0 0 0 0;ct1 0 0 0;0 ct2 0 0;0 0 ct3 0;0 0 0 ct4];

Kt=[0 0 0 0;0 0 0 0;0 0 0 0;kt1 0 0 0;0 kt2 0 0;0 0 kt3 0;0 0 0 kt4];

for n=1:1000

f=0.1*n;

Kq=Kt+2*pi*f*j*Ct

H=inv(K-(2*pi*f)^2*M+(2*pi*f)*j*C)*Kq;%(:,:,n)7X4 频率响应函数

[color=Red]sq=csp*v/f/f;%路面输入谱矩阵

tao=1/v;

tt=2*pi*f*tao*j

Sq=sq*[1 exp(-tt) 0.18 0.18*exp(-tt);exp(tt) 1 0.18*exp(tt)  0.18;...

0.18 0.18*exp(-tt) 1 exp(-tt);0.18*exp(tt) 0.18 exp(tt) 1];[/color]

%第一种求悬架动挠度响应谱

Hfd1=H(1,:)-l1*H(2,:)-l3*H(3,:)-H(4,:);%1x4

Hfd2=H(1,:)+l2*H(2,:)-l3*H(3,:)-H(5,:);

Hfd3=H(1,:)-l1*H(2,:)+l4*H(3,:)-H(6,:);

Hfd4=H(1,:)+l2*H(2,:)+l4*H(3,:)-H(7,:);

Sfd_q1=conj(Hfd1)*Sq*Hfd1.';

Sfd_q2=conj(Hfd2)*Sq*Hfd2.';

Sfd_q3=conj(Hfd3)*Sq*Hfd3.';

Sfd_q4=conj(Hfd4)*Sq*Hfd4.';

Sfd1(n)=abs(sqrt(Sfd_q1));

Sfd2(n)=abs(sqrt(Sfd_q2));

Sfd3(n)=abs(sqrt(Sfd_q3));

Sfd4(n)=abs(sqrt(Sfd_q4));

%教材上的扩大倍数求法

Sz1=conj(H)*Sq*H.';%conj(H):H 的共轭矩阵

Sz=sqrt(Sz1)

%Szq(n)=abs(Sz(1,1))*(2*pi*f)^4 %车身质心加速度功率谱

Sd1=Sz(1,1)-l1*Sz(2,2)-l3*Sz(3,3)-Sz(4,4);

Sd2=Sz(1,1)+l2*Sz(2,2)-l3*Sz(3,3)-Sz(4,4);%Sz(5,5)

Sd3=Sz(1,1)-l1*Sz(2,2)+l4*Sz(3,3)-Sz(4,4);%Sz(6,6)

Sd4=Sz(1,1)+l2*Sz(2,2)+l4*Sz(3,3)-Sz(4,4);%Sz(7,7)

Sf1(n)=abs(sqrt((k1*Sd1)^2+(c1*2*pi*f*Sd1)^2));%

Sf2(n)=abs(sqrt((k2*Sd2)^2+(c2*2*pi*f*Sd2)^2));

Sf3(n)=abs(sqrt((k3*Sd3)^2+(c3*2*pi*f*Sd3)^2));

Sf4(n)=abs(sqrt((k4*Sd4)^2+(c4*2*pi*f*Sd4)^2))

end

f1=0.1:0.1:100;

loglog(f1,Sfd1,f1,Sfd2,f1,Sfd3,f1,Sfd4)

figure

loglog(f1,Sf1,f1,Sf2,f1,Sf3,f1,Sf4)