matlab问题

zzz.dan

我的代码如下：
close all;
$ |$ ~- ^2 T+ K+ fclear all;
. h  N) Z% H3 E3 B+ j% H+ DL=1; EI=5; T=2; rho=2; M=6;
% q7 j3 m8 U: e1 F" j, }8 cnt=80000; nx=20; tf=15; Ttr=50;
dx=L/nx; dt=tf/nt;
x = linspace(0,L,nx+1); t_tr = linspace(0,tf,Ttr+1);
* u# d1 Y1 X: b6 Q0 Z2 _9 T# r0 xw=zeros(nx+1,nt+1);
) r' n6 i* h( E; L; dw3D_free=zeros(Ttr+1,nx+1);
d=zeros(nt+1,1);
1 A8 v! t; b( ~for j=1:nt+1
  O6 v- o$ M; ^2 V6 D$ l    d(j) = 1 + sin(2*pi*(j-1)*dt) + cos(3*pi*(j-1)*dt) ;
end
Br=90;Bl=-90;k=100;k1=1;k2=4;r=1;m=1;
w=zeros(nx+1,nt+1); w3D_model=zeros(Ttr+1,nx+1);
/ U9 B8 R% `+ Iu_model=zeros(nt+1,1);u3D_model=zeros(Ttr+1,1);
ua=zeros(nt+1,1); d_estimate=zeros(nt+1,1);
: P$ c$ ?+ B: N$ Rd1=zeros(nt+1,1);
dbar=max(abs(d));
d_estimate(1:2)= dbar;
for i=1:nx+1
  w(i,1)=(i-1)*dx;
8 h) S7 g6 I8 O4 u  w(i,2)=w(i,1);
5 i* Z' P9 a1 m9 d8 t8 }end
w(1,=0;w(2,=w(1,;
) y3 M, Z# v5 v+ P# l8 ^% vw3D_model(1,:)=w(:,1);w3D_model(2,:)=w(:,2);
for j=3:nt+1
4 k' H% R7 P, b    for i=3:nx-1
# C3 Z9 y3 p+ k% i  I8 v1 s% i) h        wxx=(w(i+1,j-1)-2*w(i,j-1)+w(i-1,j-1))/dx^2;
3 x3 m3 s& k& x4 a        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;
2 {. S2 t) @6 x0 O7 n! C* F    end

; Z' F/ ~" j1 I4 X    wxl=(w(nx+1,j-1)-w(nx,j-1))/dx;
+ F3 j& R% |4 s  P. m' ^* A/ M! H- a. @    wxxxl=(w(nx+1,j-1)-3*w(nx,j-1)+3*w(nx-1,j-1)-w(nx-2,j-1))/dx^3;
# T) |9 t$ s/ y# Z8 m' C7 b( B    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);
8 m2 ~! s; h1 R    wxtl=(w(nx+1,j-1)-w(nx,j-1)-w(nx+1,j-2)+w(nx,j-2))/(dx*dt);
: @0 q% m: W) _& {9 C    wtl=(w(nx+1,j-1)-w(nx+1,j-2))/dt;
/ V  K, C6 Z3 Q    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;
& ~  N, E, n5 H9 V* K4 ^" |    dv=(u_model(j)-u_model(j-1))/dt;
6 A% d  t# x# }* o4 D' |6 U# {    if  dv>0
# g6 L* L7 ]/ V/ s7 K        d0=-m*Br;
    else
" Z2 E& h2 `! i5 u1 h3 d6 T7 g        d0=-m*Bl;
; s! z' L3 s7 ~/ Q/ N    end
4 I' ^- p/ V# Q    d1(j)=d0+d(j);
' ~4 V# ?! g6 |8 V, ]    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;
* S  d; f4 w; D1 ^    w(nx,j)=(w(nx+1,j)+w(nx-1,j))/2;
    if mod(j-1,(nt/Ttr))==0
6 e4 C' R2 D7 }/ T        w3D_model(1+(j-1)*Ttr/nt,:)=w(:,j);
        u3D_model(1+(j-1)*Ttr/nt)=u_model(j);
( t# f4 c2 b3 W    end
end
' F7 P. ^$ A: @figure(1);
2 b6 \8 H# R% ysuRF(x,t_tr,w3D_model);
title('Displacement of beam with robust boundary control');
ylabel('Time (s)','Fontsize',12);
xlabel('x(m)','Fontsize',12);
zlabel('w(x,t)(m)','Fontsize',12);
7 Y5 ]) m0 U, W* i( n" M' yfigure(2);
6 [7 A. l% e4 q/ Rplot(t_tr,u3D_model,'r')
xlabel('time');
/ f. ^) f5 H5 U, r" K; l0 ~ylabel('u');
; j- A, k8 P" [8 F! r# s1 W
2 @0 V6 V# @. k! J/ }4 S
7 v+ U9 M' n8 O* D以上我的代码没出问题，但图像出不来,求助!
( W8 l' C3 F, F# ~4 r) _

小小鲁班

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
1 w  D8 \/ I$ \* B' s# E/ X" D修改意见可以减少迭代次数，或者改用sym符号计算试试

Uifhjvv

来学习一下

勇往直前11

zaiyiaaaa