碰撞问题

xiejinling

有没有做碰撞问题的，我下面作了一个单自由度的碰撞的matlab程序，但是运行不出结果，不知道问题出在那里，请指教。谢谢！
%方程
%d2x/dt2+k*dx/dt+x=cos(w*t)
%边界条件 x(t+)=x(t-),dx/dt(t+)=-r* dx/dt(t-)


clear;
global w;
a=0.01;
pi=3.1415;
r=0.6;
w=2.5;
fid= fopen('xjl.txt','w');
%fid1= fopen('xjl1.txt','w');

x=1.6;
v=0;
t=0;
t1=0.05
n=1;

while j<500
t=t+j*t1
A=(1-w*w)/((1-w*w)^2+(2*w*a)^2);
B=2*w*a/((1-w*w)^2+(2*w*a)^2);
w1=sqrt(1-a^2);
c=2*n*pi/w;
d=exp(-c*a);
e=-r*d;
w2=w1*c;
f=d*sin(w2)/(d-1);
g=w1+e*a*sin(w2)-e*w1*cos(w2);
c1=f*(e*g*w-g*w)/(g-a*f*cos(w2)-e*f*w1*sin(w2));
c2=(e*g*w-g*w)/(g-a*f*cos(w2)-e*f*w1*sin(w2));

x=exp(-a*t)*(c1*cos(w1*t)+c2*sin(w1*t))+A*sin(w*t)+B*cos(w*t);
v=exp(-a*t)*((c2*w1-a*c1)*cos(w1*t)-(a*c2+c1*w1)*sin(w1*t))+A*w*cos(w*t)-B*w*sin(w*t);
  
if x<=0
      x=x;
      v=-r*v;
  end
  
  if  j>200
          fprintf(fid,'%12.8f %12.8f %12.8f \n',t,x,v);
  end
end

fclose(fid);
load xjl.txt;
figure
plot(xjl(:,2),xjl(:,3));

xiejinling

希望高手指点，邮箱：xiejl19840228@163.com