newmark法和wilson_cita法为什么计算结果差很大
这是NewMark程序%NewMark Function
clc,clear
m1=341; % Quality of the former landing gear
m2=363; % Quality of the back landing gear
mc=24900; % Quality of the body
Jc=424378; % Moment of inertia of the body
Lz=5.805; % Distance between the former and the back tyre（m
a=0.91*Lz; % Distance of the formwe tyre spread
b=0.09*Lz; % Distance of the back tyre spread
k1=8.1123e5; % Rigidity coefficient of the former tyre
c1=3*1e4; % Damp coefficient of the former tyre
k3=1.6285e5; % Rigidity coefficient of the former Amortine shoring
c3=5*1e4; % Damp coefficient of the former Amortine shoring
k2=4.41046; % Rigidity coefficient of the back tyre
c2=3*1e4; % Damp coefficient of the back tyre
k4=7.3361e5; % Rigidity coefficent of the back Amortine shoring
c4=5*1e4; % Damp coefficient of the back Amortine shoring
g=9.81;
M =[(mc*b^2+Jc)/Lz^2(mc*a*b-Jc)/Lz^20 0
(mc*a*b-Jc)/Lz^2(mc*a^2+Jc)/Lz^20 0
0 0 m10
0 0 0 m2];
C =[c4 0-c4 0
0 c30-c3
-c40c2+c40
-c300 c1+c3];
K =[k4 0-k4 0
0 k30-k3
-k40k2+k40
-k300 k1+k3];
u=';
v=';
a=';
t(1)=0; %time
x(:,1)=u; %distance
x1(:,1)=v; %speed
x2(:,1)=a; %accelerrate
%newmar kcoefficient
gama=0.5;
dt=1e-2;
delta=0.25;
b0=1/(delta*dt^2);
b1=gama/(delta*dt);
b2=1/(delta*dt);
b3=1/(2*delta)-1;
b4=gama/delta-1;
b5=dt*(gama/(2*delta)-1);
b6=dt*(1-gama);
b7=gama*dt;
%等效刚度矩阵
K1=K+b0*M+b1*C;
t_max=10; %计算时间总长
i=1;
t(1)=0.01;
q=zeros(4,1);
while t(i)<t_max
Q(:,i)=';
q(:,i+1)=Q(:,i)+M*(b0*x(:,i)+b2*x1(:,i)+b3*x2(:,i))+C*(b1*x(:,i)+b4*x1(:,i)+b5*x2(:,i));
x(:,i+1)=inv(K1)*q(:,i+1);
x2(:,i+1)=b0*(x(:,i+1)-x(:,i))-b2*x1(:,i)-b3*x2(:,i);
x1(:,i+1)=b1*(x(:,i+1)-x(:,i))-b4*x1(:,i)-b5*x2(:,i);
i=i+1;
t(i)=t(i-1)+dt;
end
这是wilson_cita程序
%wilson_cita Function
clc,clear
m1=341; % Quality of the former landing gear
m2=363; % Quality of the back landing gear
mc=24900; % Quality of the body
Jc=424378; % Moment of inertia of the body
Lz=5.805; % Distance between the former and the back tyre（m
a=0.91*Lz; % Distance of the formwe tyre spread
b=0.09*Lz; % Distance of the back tyre spread
k1=8.1123e5; % Rigidity coefficient of the former tyre
c1=3*1e4; % Damp coefficient of the former tyre
k3=1.6285e5; % Rigidity coefficient of the former Amortine shoring
c3=5*1e4; % Damp coefficient of the former Amortine shoring
k2=4.41046; % Rigidity coefficient of the back tyre
c2=3*1e4; % Damp coefficient of the back tyre
k4=7.3361e5; % Rigidity coefficent of the back Amortine shoring
c4=5*1e4; % Damp coefficient of the back Amortine shoring
g=9.81;
M =[(mc*b^2+Jc)/Lz^2(mc*a*b-Jc)/Lz^20 0
(mc*a*b-Jc)/Lz^2(mc*a^2+Jc)/Lz^20 0
0 0 m10
0 0 0 m2];
C =[c4 0-c4 0
0 c30-c3
-c40c2+c40
-c300 c1+c3];
K =[k4 0-k4 0
0 k30-k3
-k40k2+k40
-k300 k1+k3];
u=';
v=';
a1=';
t(1)=0; %time
x(:,1)=u; %distance
x1(:,1)=v; %speed
x2(:,1)=a1; %accelerate
%cita coeffi
cita=1.42085;
dt=1e-2;
s=cita*dt;
%等效刚度矩阵
K1=K+6/s^2*M+3/s*C;
t_max=10; %full time
i=1;
t(1)=0.01;
while t(i)<t_max
Q=';
gt=M*(6/s^2*x(:,i)+6/s*x1(:,i)+2*x2(:,i))+C*(3/s*x(:,i)+2*x1(:,i)+s/2*x2(:,i))+Q;
us=inv(K1)*gt;
u2s=6/s^2*(us-x(:,i))-6/s*x1(:,i)-2*x2(:,i);
u1s=3/s*(us-x(:,i))-2*x1(:,i)-s/2*x2(:,i);
x2(:,i+1)=6*dt/s^3*(us-x(:,i))-6*dt/s^2*x1(:,i)+(1-3/cita)*x2(:,i);
x1(:,i+1)=x1(:,i)+dt/2*(u2s-x2(:,i));
x(:,i+1)=x(:,i)+dt*x1(:,i)+dt^2/6*(u2s+2*x2(:,i));
i=i+1;
t(i)=t(i-1)+dt;
end
两个程序计算结果相差很大，希望大家帮忙看看。
其中NewMark程序算出稳定解
而wilson_cita算出的竟然不稳定
各有优劣
那么你有没有wilson_cita的Matlab程序啊我觉的我的程序好像没有错 q(:,i+1)=Q(:,i)+M*(b0*x(:,i)+b2*x1(:,i)+b3*x2(:,i))+C*(b1*x(:,i)+b4*x1(:,i)+b5*x2(:,i));
好像这个语句有点问题，b0*x(:,i)应该是b6*x(:,i)吧 将你的程序改正如下。
x2(:,i+1)=6*dt/s^3*(us-x(:,i))-6*dt/s^2*x1(:,i)+(1-3/cita)*x2(:,i);
x1(:,i+1)=x1(:,i)+dt/2*(x2(:,i+1)+x2(:,i))
觉得
u1s=3/s*(us-x(:,i))-2*x1(:,i)-s/2*x2(:,i);
此行程序没有作用 威尔逊法:plot(x2(2,:))得到如下 图象,好象收敛
图象
回复 谢谢了，heyong 兄
页:
[1]