24小时热门版块排行榜

返回列表

我原是我以为

新虫 (小有名气)

应助: 0 (幼儿园)
金币: 343.8
散金: 34
红花: 2
帖子: 163
在线: 84.1小时
虫号: 4254249
注册: 2015-11-29
性别: GG
专业: 工程热物理与能源利用

[求助] 求陶老师数值传热学4.8节长方形截面对流换热C程序，或者帮忙看看程序哪里不对

求大神帮忙看看loop2，求Nu过程，哪里不对，fRe和书上差不多
#include <math.h>
#include <stdio.h>
#define NT 70
void main()
{
/*定义变量*/
int N,M,i,j,Iter;
double W[NT][NT],W0[NT][NT],sita[NT][NT],sita0[NT][NT],T[NT][NT],De,Re,Nu;
double dltx,dlty,dx,dy,Wm,D,x,y,a,b;
double Eps,DT,DTmax,SD;
double ap,aw,ae,as,an,JF,b1,b2,JA,ap1;
FILE*fp;

/*读入参数*/
printf("请输入x方向节点数N y方向节点数M\n"

;
scanf("%d%d",&N,&M);
printf("请输入界面长度a 界面宽度b\n"

;
scanf("%lf%lf",&a,&b);

/*已知参数*/
SD=0;
JF=0;
dltx=a;
dlty=b;
De=2*dltx*dlty/(dltx+dlty);
D=a;/*特征尺寸*/
N;
M;
Eps=1.e-6;/*计算精度*/

/*基本参数*/
dx=dltx/(N-1);
dy=dlty/(M-1);
ae=aw=dy/dx;
an=as=dx/dy;

/*速度边界条件处理*/
/*固体表面速度为0值*/
for(j=0;j<M;j++)
  {
W[0][j]=0;
W[N-1][j]=0;
}
for(i=1;i<N-1;i++)
{
W[0]=0;
W[M-1]=0;
}
/*速度场赋初值*/
for(j=1;j<M-1;j++)
  for(i=1;i<N-1;i++)
{
W[j]=0.00008;
W0[j]=W[j];/*对上一次迭代值赋值*/
}
Iter=0;
/*Gauss-Seidel迭代计算速度*/
  loop1:
for(j=1;j<M-1;j++)
{
for(i=1;i<N-1;i++)
  {
  ap=ae+aw+an+as;
b1=(dx*dy)/(D*D);
W[j]=ae/ap*W[i+1][j]+aw/ap*W[i-1][j]+an/ap*W[j+1]+as/ap*W[j-1]+b1/ap;
}
}
/*计算两次迭代最大误差*/
Iter=Iter+1;
DTmax=0.0;
for(j=1;j<M-1;j++)
  for(i=1;i<N-1;i++)
{
DT=fabs(W[j]-W0[j]);
W0[j]=W[j];
if(DT>DTmax)
DTmax=DT;
}
if(DTmax>Eps)
goto loop1;

for(j=1;j<M-1;j++)
for(i=1;i<N-1;i++)
{
SD+=W[j]*dx*dy;
}
Wm=SD/(a*b);
printf("无量纲速度平均值:\n";
printf("Wm=%5.6f\n",Wm);

/*边界条件处理*/
/*固体表面无量纲为0,取T=JA*sita进行迭代求JA*/
for(j=0;j<M;j++)
  {
sita[0][j]=0;
sita[N-1][j]=0;
}
for(i=1;i<N-1;i++)
{
sita[0]=0;
sita[M-1]=0;
}
/*sita场赋初值*/
for(j=1;j<M-1;j++)
  for(i=1;i<N-1;i++)
{
sita[j]=0.003;
sita0[j]=sita[j];/*对上一次迭代值赋值*/
}

Iter=0;
/*Gauss-Seidel迭代计算*/
loop2:
for(j=1;j<M-1;j++)
{
for(i=1;i<N-1;i++)
{
JF+=W[j]*sita[j];
}
}
JA=a*b/JF/Wm;
for(j=1;j<M-1;j++)
for(i=1;i<N-1;i++)
{
b2=dx*dy*JA*W[j]/(D*D*Wm);
ap1=ap-b2;
sita[j]=ae/ap1*sita[i+1][j]+aw/ap1*sita[i-1][j]+an/ap1*sita[j+1]+as/ap1*sita[j-1];
}

/*计算两次迭代最大误差*/
Iter=Iter+1;
DTmax=0.0;
for(j=1;j<M-1;j++)
  for(i=1;i<N-1;i++)
{
DT=fabs(sita[j]-sita0[j]);
sita0[j]=sita[j];
if(DT>DTmax)
DTmax=DT;
}
if(DTmax>Eps)
goto loop2;

/*Re计算*/
printf("Re数为:\n";
Re=2*(De/D)*(De/D)/Wm;
printf("%5.6f\n",Re);

/*Nu计算*/
printf("Nu数为:\n";
Nu=0.25*(De/D)*(De/D)*JA;
printf("%5.6f\n",Nu);

for(j=1;j<M-1;j++)
  for(i=1;i<N-1;i++)
  {
   T[j]=JA*sita[j];
  }

  fp=fopen("result.dat","w";
  printf("无量纲温度:\n";
  for(j=M-1;j>=0;j--)
  {
   for(i=0;i<N;i++)
   {
   fprintf(fp,"%5.6f",T[j]);
   }
   printf("\n";
  }
  fclose(fp);
}

发自小木虫IOS客户端

回复此楼

» 猜你喜欢

【修正版】单环用磁约束低速电子实验方案（简化版）已经有0人回复
桂林理工大学物理学专业招收调剂，还有三个名额！！！已经有22人回复
物理学I论文润色/翻译怎么收费? 已经有74人回复
考博自荐已经有4人回复
山东大学第二批博士研究生招生已经有0人回复
中国科学院东莞材料科学与技术研究所-2026年博士招生-吴昊研究员-磁学与自旋电子学已经有0人回复
《电磁学》教材推荐已经有1人回复
【急招】合肥工大核聚变材料计算方向2026级工程博士生已经有4人回复
大豆异黄酮分离已经有0人回复
湖南大学材料学院急招2026年博士生，临时增加一名博士联培指标已经有10人回复
天津理工大学晶体材料全国重点实验室刘红军教授课题组招收博士生1-2名已经有0人回复

1楼 2017-04-12 19:48:44

已阅回复此楼关注TA 给TA发消息送TA红花 TA的回帖

相关版块跳转我要订阅楼主我原是我以为的主题更新

返回列表