| 查看: 352 | 回复: 1 | |||
[求助]
求助，帮我分析下程序哪里出了问题？求指教
|
|
function H=Ha(nx,nz,hm,epslong,elta,N,w2,Ps,d,L) rj=d/2; u2=2*pi*N/60; deltax=2*pi*rj/(nx-1) deltaz=L/(nz-1) delta=2*pi/(nx-1) h=zeros(nx+1,nz+1); a=zeros(nx+1,nz+1); b=zeros(nx+1,nz+1); c=zeros(nx+1,nz+1); d=zeros(nx+1,nz+1); e=zeros(nx+1,nz+1); f=zeros(nx+1,nz+1); ff=zeros(nx+1,nz+1); k=1; Ax=12*elta*u2*L/(hm^2*Ps); Az=12*elta*w2*L/(hm^2*Ps); for i=1:nx theta=delta*(i-1); h(i,  =hm*(1+epslong*cos(theta+delta/2));end S=0; T=0; ERR=1e-3; GAP=1; while GAP>ERR k=k+1 for j=2:1:nz-1 if i==1 a(i,j)=h(i,j)^3/deltax^2-3*h(i,j)^2*(h(i+1,j)-h(nx,j))/(4*deltax^2); b(i,j)=h(i,j)^3/deltax^2+3*h(i,j)^2*(h(i+1,j)-h(nx,j))/(4*deltax^2); c(i,j)=-2*h(i,j)^3/deltax^2-2*h(i,j)^3/deltaz^2; d(i,j)=h(i,j)^3/deltaz^2; e(i,j)=h(i,j)^3/deltaz^2; ff(i,j)=(a(i,j)*f(nx,j)+b(i,j)*f(i+1,j)+d(i,j)*f(i,j-1)+e(i,j)*f(i,j+1))/(-c(i,j)); elseif i==nx a(i,j)=h(i,j)^3/deltax^2-3*h(i,j)^2*(h(1,j)-h(i-1,j))/(4*deltax^2); b(i,j)=h(i,j)^3/deltax^2+3*h(i,j)^2*(h(1,j)-h(i-1,j))/(4*deltax^2); c(i,j)=-2*h(i,j)^3/deltax^2-2*h(i,j)^3/deltaz^2; d(i,j)=h(i,j)^3/deltaz^2; e(i,j)=h(i,j)^3/deltaz^2; pp(i,j)=(a(i,j)*f(i-1,j)+b(i,j)*f(1,j)+d(i,j)*f(i,j-1)+e(i,j)*f(i,j+1))/(-c(i,j)); else if pp(i,j)<0 pp(i,j)=0; else end end end for i=2:1:nx-1 if j==1 a(i,j)=h(i,j)^3/deltax^2-3*h(i,j)^2*(h(i+1,j)-h(i-1,j))/(4*deltax^2); b(i,j)=h(i,j)^3/deltax^2+3*h(i,j)^2*(h(i+1,j)-h(i-1,j))/(4*deltax^2); c(i,j)=-2*h(i,j)^3/deltax^2-2*h(i,j)^3/deltaz^2; d(i,j)=h(i,j)^3/deltaz^2; e(i,j)=h(i,j)^3/deltaz^2; ff(i,j)=(a(i,j)*f(i-1,j)+b(i,j)*f(i+1,j)+d(i,j)*Ps^2+e(i,j)*f(i,j+1))/(-c(i,j)); elseif j==nz a(i,j)=h(i,j)^3/deltax^2-3*h(i,j)^2*(h(i+1,j)-h(i-1,j))/(4*deltax^2); b(i,j)=h(i,j)^3/deltax^2+3*h(i,j)^2*(h(i+1,j)-h(i-1,j))/(4*deltax^2); c(i,j)=-2*h(i,j)^3/deltax^2-2*h(i,j)^3/deltaz^2; d(i,j)=h(i,j)^3/deltaz^2; e(i,j)=h(i,j)^3/deltaz^2; ff(i,j)=(a(i,j)*f(i-1,j)+b(i,j)*f(i+1,j)+d(i,j)*f(i,j-1)+e(i,j)*f(i,j+1))/(-c(i,j)); else if ff(i,j)<0 ff(i,j)=0; else end end end if i==1&&j==1 p(i,j)=Ps elseif i==nx&&j==1 p(i,j)=Ps elseif i==1&&j==nz a(i,j)=h(i,j)^3/deltax^2-3*h(i,j)^2*(h(1,j)-h(i-1,j))/(4*deltax^2); b(i,j)=h(i,j)^3/deltax^2+3*h(i,j)^2*(h(1,j)-h(i-1,j))/(4*deltax^2); c(i,j)=-2*h(i,j)^3/deltax^2-2*h(i,j)^3/deltaz^2; d(i,j)=h(i,j)^3/deltaz^2; e(i,j)=h(i,j)^3/deltaz^2; ff(i,j)=(a(i,j)*f(nx,j)+b(i,j)*f(i+1,j)+d(i,j)*f(i,j-1)+e(i,j)*f(i,nz))/(-c(i,j)); elseif i==nx&&j==nz a(i,j)=h(i,j)^3/deltax^2-3*h(i,j)^2*(h(1,j)-h(i-1,j))/(4*deltax^2); b(i,j)=h(i,j)^3/deltax^2+3*h(i,j)^2*(h(1,j)-h(i-1,j))/(4*deltax^2); c(i,j)=-2*h(i,j)^3/deltax^2-2*h(i,j)^3/deltaz^2; d(i,j)=h(i,j)^3/deltaz^2; e(i,j)=h(i,j)^3/deltaz^2; ff(i,j)=(a(i,j)*f(i-1,j)+b(i,j)*f(1,j)+d(i,j)*f(i,j-1)+e(i,j)*f(i,nz-1))/(-c(i,j)); else if ff(i,j)<0 ff(i,j)=0; else end end for i=2:1:nx for j=2:1:nz S=sum(sum(abs(ff-f))); T=sum(sum(abs(f))); end end GAP=S/T; f=ff end |
回复此楼» 猜你喜欢
寻找论文指导老师
已经有6人回复
-
南方科技大学机械与能源工程系-邓辉研究员招收2026级博士生
已经有0人回复
-
机械工程论文润色/翻译怎么收费?
已经有206人回复
-
南方科技大学机械与能源工程系-邓辉研究员招收2026级博士生
已经有0人回复
-
同一篇文章，用不同账号投稿对编辑决定是否送审有没有影响？
已经有3人回复
-
郑州大学田佳佳团队诚招2026年入学博士研究生
已经有0人回复
-
南方科技大学机械与能源工程系-邓辉研究员招收2026级博士生
已经有0人回复
-
材料科学基础
已经有0人回复
-
2025-2026国产空气检漏仪/气密测漏仪/气密性测试仪推荐生产厂家
已经有0人回复
-
南方科技大学机械与能源工程系-邓辉研究员招收2026级博士生
已经有0人回复
-
南方科技大学机械与能源工程系-邓辉研究员招收2026级博士生
已经有0人回复
2楼2016-03-19 00:24:51