马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。
您需要 登录 才可以下载或查看,没有账号?我要加入
x
我有一个偏微分方程边界条件是随时间变化的,请问怎么才能解,我试了好多方法!!谢谢指导.
- clear
- kongjian=1;
- shijian=50000;
- t=1000;
- h=0.02;
- n=(kongjian/h)+1;
- m=(shijian/t)+1;
- cs=2;
- k=0.001;
- %k=1;
- k1=0.000001
- u1=zeros(n,m);
- %a xiao c de kuo san xi shu
- a=0.0000001;
- r=1/(t*k/2+1);
- s=a*t/h^2;
- ug=20;
- lamud=5e-008;
- %bian jie tiao jian
- %u1(1,1:m)=20;
- % chu shi tiao jian
- %con=0
- u1(1:n,1:m)=1;
- con=4
- for j=2:m
- u1(1,j)=r*(a*t/h^2*(u1(2,j-1)-2*u1(1,j-1)+u1(1,j-1)+h*lamud/a*(ug-u1(1,j-1)))+(1-t*k/2)*u1(1,j-1)+t*k*cs)
- for i=2:n-1
- % if h*(i)<=(k1*1000*(j-1))^0.5
- if i<con-1
- u1(i,j)=r*(a*t/(h^2)*(u1(i+1,j-1)-2*u1(i,j-1)+u1(i-1,j-1))+(1-t*k/2)*u1(i,j-1)+t*k*cs)
- else
- if con-1==i
- u1(i,j)=r*(a*t/h^2*(u1(i,j-1)-2*u1(i,j-1)+u1(i,j-1)+h*lamud/a*(u1(i-1,j-1)-u1(i,j-1)))+(1-t*k/2)*u1(i,j-1)+t*k*cs)
- % continue
- else
- u1(i,j)=r*(a*t/(h^2)*(u1(i+1,j-1)-2*u1(i,j-1)+u1(i-1,j-1))+(1-t*k/2)*u1(i,j-1)+t*k*cs)
- %continue
- end
- end
- if h*(i)>(k1*t*(j-1))^0.5
- con=i
- break
- else
- %break
- continue
- end
- break
- u1(n,j)=u1(n-1,j)
- end
- u1(n,j)=u1(n-1,j)
- end
- w=linspace(0,shijian,m);
- w2=u1(1,:)
- w3=u1(2,:)
- w4=u1(4,:)
- w5=u1(6,:)
- w6=u1(8,:)
- w8=u1(10,:)
- %plot(w,w2,'b')
- ww=linspace(0,kongjian,n);
- ww2=u1(:,40)
- ww3=u1(:,10)
- ww4=u1(:,2)
- %plot(ww,ww2,'r')
- subplot(2,1,1),plot(w,w2,'b',w,w3,'r',w,w4,'y',w,w5,'c',w,w6,'k',w,w8,'g')
- subplot(2,1,2),plot(ww,ww2,'r',ww,ww3,'b',ww,ww4,'y')
- %fevel(pde,n-haha,m-hehe)
复制代码 |