zbjfsss.m

内燃机转子仿真

clear all
clc
close all
m=[2 0;0 1];
c=[0 0;0 0];
k=[6 -2;-2 4];
F=[0;10];
i=1;
f=F;
b=1.4;
n=length(F);
x10=zeros(n,1);
x20=zeros(n,1);
x30=inv(m)*(f-c*x20-k*x10)
dt=0.28;
t=6;
td=rem(t,dt);
tdd=t-td;
n1=tdd/dt;
xy=zeros(2*n+1,n1);
xy(1,i)=i;
xy(2:n+1,i)=x10;
xy(n+2:2*n+1,i)=x20;
  
a1=6/b^3/dt^2;
a2=6/b^2/dt;
a3=1-3/b;
a4=dt/2;
a5=dt^2/6;
a6=6/b^2/dt^2;
a7=6/b/dt;
a8=3/b/dt;
a9=b*dt/2;
k1=k+a2*m+a8*c;
f1=F;

%for i=2:n1
    %for j=1:60
  
     fa=f1+m*(a6*x10+a7*x20+2*x30)+c*(a8*x10+2*x20+a9*x30);
     x11=inv(k1)*fa
     x3=a1*(x11-x10)-a2*x20+a3*x30;
     x2=x20+a4*(x3+x30);
     x1=x10+dt*x20+(x3+2*x30)*a5;
     x10=x1
     x20=x2;
     x30=x3;

   
   
    %end
    
%end


 %xy(1,i)=i;
    %xy(2:n+1,i)=x10;
   % xy(n+2:2*n+1,i)=x20;
%figure(1);
%plot(xy(1,:),xy(2,:));
%figure(2);
%plot(xy(1,:),xy(3,:));