matlab问题

zzz.dan

我的代码如下：
1 B. _  U/ h- j5 vclose all;
: v& l, t4 i9 V  T! V) i5 ~+ `clear all;
& N  M" X( W6 p7 L2 Y+ ML=1; EI=5; T=2; rho=2; M=6;
1 L4 n- H# C) Y9 L; i5 xnt=80000; nx=20; tf=15; Ttr=50;
5 h" O8 Q# G. ]0 t$ K. w$ M1 Rdx=L/nx; dt=tf/nt;
! R0 p# y: r$ _5 ~, R! tx = linspace(0,L,nx+1); t_tr = linspace(0,tf,Ttr+1);
6 H# g& `3 T2 B: m9 X+ Vw=zeros(nx+1,nt+1);
w3D_free=zeros(Ttr+1,nx+1);
d=zeros(nt+1,1);
for j=1:nt+1
    d(j) = 1 + sin(2*pi*(j-1)*dt) + cos(3*pi*(j-1)*dt) ;
/ b- d- o& n1 h2 ?end
- C. I6 O: o3 D& K5 ZBr=90;Bl=-90;k=100;k1=1;k2=4;r=1;m=1;
0 i3 n" N  B- D' M# ~) M8 {& |w=zeros(nx+1,nt+1); w3D_model=zeros(Ttr+1,nx+1);
u_model=zeros(nt+1,1);u3D_model=zeros(Ttr+1,1);
, n( J; c7 F2 _; _ua=zeros(nt+1,1); d_estimate=zeros(nt+1,1);
: q2 x' B# g3 I& Sd1=zeros(nt+1,1);
6 v0 i. k) z0 Odbar=max(abs(d));
+ h6 x- U1 u% Z' _d_estimate(1:2)= dbar;
  f: z9 I2 j9 E' X, {for i=1:nx+1
% E2 W  `% G2 S$ c8 W( H  w(i,1)=(i-1)*dx;
" D7 {& m* p, r; \  w(i,2)=w(i,1);
end
w(1,=0;w(2,=w(1,;
/ Y" W, g7 V- k- ~w3D_model(1,:)=w(:,1);w3D_model(2,:)=w(:,2);
  C7 _. v9 T7 c6 G. Ifor j=3:nt+1
0 M1 n% P$ m" R' B& l    for i=3:nx-1
" ]& o6 ~4 w5 D8 F6 I1 L/ O        wxx=(w(i+1,j-1)-2*w(i,j-1)+w(i-1,j-1))/dx^2;
        wxxxx=(w(i+2,j-1)-4*w(i+1,j-1)+6*w(i,j-1)-4*w(i-1,j-1)+w(i-2,j-1))/dx^4;
        w(i,j)=2*w(i,j-1)-w(i,j-2)+(T*wxx-EI*wxxxx)*dt^2/rho;
& i8 w* v& `) @6 p, J) e; C4 `    end
" d! l- Q2 k* @! S7 @" U
' t9 }% [8 n1 W    wxl=(w(nx+1,j-1)-w(nx,j-1))/dx;
    wxxxl=(w(nx+1,j-1)-3*w(nx,j-1)+3*w(nx-1,j-1)-w(nx-2,j-1))/dx^3;
) ^. p- x3 H& U- \4 k' M    wxxxtl=(w(nx+1,j-1)-3*w(nx,j-1)+3*w(nx-1,j-1)-w(nx-2,j-1)-w(nx+1,j-2)+3*w(nx,j-2)-3*w(nx-1,j-2)+w(nx-2,j-2))/(dt*dx^3);
    wxtl=(w(nx+1,j-1)-w(nx,j-1)-w(nx+1,j-2)+w(nx,j-2))/(dx*dt);
    wtl=(w(nx+1,j-1)-w(nx+1,j-2))/dt;
    ua(j)=wtl+k1*wxl-k2*wxxxl;
    d_estimate(j)= (abs(ua(j))+d_estimate(j-1))/(1+r*dt);
    u_model(j)= (-EI*wxxxl+T*wxl-k1*M*wxtl+k2*M*wxxxtl-k*ua(j)-sign(ua(j))*d_estimate(j))/m;
    dv=(u_model(j)-u_model(j-1))/dt;
, S4 Q# m, F7 k7 w  o3 f    if  dv>0
        d0=-m*Br;
( N6 F' Q8 F/ R4 k    else
8 V6 N& m2 g! U: _/ N        d0=-m*Bl;
  W& J& v, u/ u+ |' [, X; i    end
& ]0 |& \6 [  V: F: w    d1(j)=d0+d(j);
3 L) g/ j- o7 t    w(nx+1,j)=2*w(nx+1,j-1) - w(nx+1,j-2) + (m*u_model(j-1)+d1(j-1)+EI*wxxxl-T*wxl)*dt^2/M;
- `9 ]- X. d  q8 Z    w(nx,j)=(w(nx+1,j)+w(nx-1,j))/2;
    if mod(j-1,(nt/Ttr))==0
9 S0 \1 r' f$ s$ G- V        w3D_model(1+(j-1)*Ttr/nt,:)=w(:,j);
        u3D_model(1+(j-1)*Ttr/nt)=u_model(j);
5 X! F7 {! \& Y6 e" Q7 R" \    end
3 z. Y9 f5 U4 ]; f) ?$ H+ @end
figure(1);
! p- O$ w; c* r0 m% A( ?) VsuRF(x,t_tr,w3D_model);
' [. R' O; s3 ^8 a4 @title('Displacement of beam with robust boundary control');
ylabel('Time (s)','Fontsize',12);
8 [1 X- f$ B8 S' C8 d  Lxlabel('x(m)','Fontsize',12);
zlabel('w(x,t)(m)','Fontsize',12);
. h) ]7 b  I* O8 c7 kfigure(2);
plot(t_tr,u3D_model,'r')
xlabel('time');
ylabel('u');

0 v" C* }1 f4 b1 n
7 `+ `& Y7 J% o* m+ J4 T  a以上我的代码没出问题，但图像出不来,求助!
/ N6 u: ~7 ]+ V" T' Q/ ~

小小鲁班

w(nx+1,j)=2*w(nx+1,j-1) - w(nx+1,j-2) + (m*u_model(j-1)+d1(j-1)+EI*wxxxl-T*wxl)*dt^2/M;
这行代码近似看成w(nx+1,j)=2*w(nx+1,j-1)，每次迭代乘以二，当迭代超过1000多次就会出现inf，进而出现NAN
修改意见可以减少迭代次数，或者改用sym符号计算试试

Uifhjvv

来学习一下

勇往直前11

zaiyiaaaa