二维数值传热学C程序代码

维数值传热学代码档编制序号LLETD298-POI08]维数值传热学代码档编制序号LLETD298-POI08]#include<>#include<iostream>#include<>#include<fstream>/////////////////////////////////////////////////////////////////////////////////""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""""//个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个个/////////////////////////////////////////////////////////////////////////////////〃从此开始至下一个/*/标记带处之间的代码可以写在头文件里。////////////////////////////////////////////////////////////////////////////////计算参数宏定义开始////////////////////////////////////////////////////////////////////////////////#defineGRID_ROW82#defineGRID_COLUMN102〃网格横纵数目，边界处多出一列处理#defineDX#defineDY〃网格单元尺寸，m#defineDEN_COPPER8900#defineDEN_STEEL7850#defineDEN_CONCRETE2430//材料的密度，kg/(m”3)#defineC_COPPER390#defineC_STEEL460#defineC_CONCRETE970//材料比热容，J/(kg*K)#defineCOND_COPPER400#defineCOND_STEEL20#defineCOND_CONCRETE〃材料导热系数，W/(m*K)#defineCONV5//空气与各种材料的对流换热系数，W/(m"2*K)#defineINITL_TMPRT298//初始温度，K#defineT_C2K(x)((x)+273)//CtoK#defineT_K2C(x)((x)-273)//KtoC〃摄氏温度与绝对温度相互转换#defineT_AIR298;//空气温度,K#definet_END3600//计算时长，s#definet_STEP1//时间步长，s#defineCOND_COP_CON#defineCOND_COP_STL#defineCOND_CON_STL〃不同材料的两个相邻结点，导热系数取调和平均数，W/(m*K)////////////////////////////////////////////////////////////////////////////////计算参数宏定义结束////////////////////////////////////////////////////////////////////////////////enumProperty{copper,steel,concrete,iso834,adiabat,air};〃结点属性的枚举类型声明structGrid{doublemNodeT[GRID_ROW][GRID_COLUMN];//结点温度PropertymNodeProperty[GRID_ROW][GRID_COLUMN];//结点材料属性doublemNodeC[GRID_ROW][GRID_COLUMN];//结点比热容doublemNodeDen[GRID_ROW][GRID_COLUMN];//结点密度};//结构体Grid定义，此为计算的主要承载数据结构，将在计算中存储温度数据〃函数全局声明开始voidSetProperty(Grid&);//属性设置函数voidInitlTMPRT(Grid&);//温度初始化函数voidSetNodePara(Grid&);//结点物性参数设置函数voidCalculateTMPRT(Grid&);//温度计算函数，几乎全部数值计算的承担者voidwrite(int,Grid&);//写文件函数，每隔100秒输出温度分布文件到当前目录〃函数全局说明结束voidSetProperty(Grid&A){for(intj=0;j<GRID_COLUMN;++j)[0][j]=air;〃第一行为空气边界，计算区域外for(intj=0;j<GRID_COLUMN;++j)[GRID_ROW-1][j]=adiabat;//最后一行为绝热边界，计算区域外for(inti=0;i<GRID_ROW;++i)if((i>=0&&i<=40)||(i>=61&&i<GRID_ROW))[i][0]=air;else[i][0]=iso834;//第一列的空气边界和ISO834边界，计算区域外for(inti=0;i<GRID_ROW;++i)[i][GRID_COLUMN-1]=air;//最后一列，空气边界，计算区域外for(inti=1;i<GRID_ROW-1;++i)for(intj=1;j<GRID_COLUMN-1;++j){if((i>=1&&i<=44&&j>=1&&j<=20)||(i>=61&&i<GRID_ROW-2&&j>=1&&j<=(100-(80-i)*(49/19))))[i][j]=concrete;〃混凝土区域elseif((i-20)*(i-20)+(j-25)*(j-25)<=20*20)[i][j]=air;//圆孔空气区域elseif(i>=41&&i<=60&&j>=1&&j<=90&&((i-20)*(i-20)+(j-25)*(j-25)>=20*20))[i][j]=copper;//金属铜区域的一部分elseif(i>=1&&i<+40&&j>=21&&j<=90)[i][j]=copper;〃金属铜区域的另一部分else[i][j]=steel;〃剩下的即为金属铁的区域}}//属性设置函数，将Grid结构体A的所有网格赋予材料属性(枚举类型Property)voidInitlTMPRT(Grid&A){for(inti=0;i<GRID_ROW;++i)for(intj=0;j<GRID_COLUMN;++j)[i][j]=INITL_TMPRT;}//温度初始化函数，将Grid结构体A的全部结点温度设置为初始温度INITL_TMPRTvoidCalculateTMPRT(Grid&A){doubleQE=0,QW=0,QN=0,QS=0;//四个方向传来的热量for(inti=1;i<GRID_ROW-1;++i)for(intj=1;j<GRID_COLUMN-1;++j){switch[i][j]){{switch[i][j+1])//判断右侧属性，并计算QE{casecopper:QE=[i][j+1]-[i][j])*COND_COPPER;break;caseair:QE=[i][j+1]-[i][j])*CONV*DY;break;casesteel:QE=[i][j+1]-[i][j])*COND_COP_STL;break;caseconcrete:QE=[i][j+1]-[i][j])*COND_COP_CON;break;caseiso834:QE=2*[i][j+1]-[i][j])*COND_COPPER;break;caseadiabat:QE=0;break;default:std::cout<<i<<”,"<<j+1<<"结点属性未指定'n";getchar();}//std::cout<<QE<<"\n";switch[i][j-1])//判断左侧属性，并计算QW{casecopper:QW=[i][j-1]-[i][j])*COND_COPPER;break;caseair:QW=[i][j-1]-[i][j])*CONV*DY;break;casesteel:QW=[i][j-1]-[i][j])*COND_COP_STL;break;caseconcrete:QW=[i][j-1]-[i][j])*COND_COP_CON;break;caseiso834:QW=2*[i][j-1]-[i][j])*COND_COPPER;break;caseadiabat:QW=0;break;default:{std::cout<<i<<","<<j-1<<"结点属性未指定\n";getchar();}}switch[i-1][j])//判断上侧属性，并计算QN{casecopper:QN=[i-1][j]-[i][j])*COND_COPPER;break;caseair:QN=[i-1][j]-[i][j])*CONV*DX;break;casesteel:QN=[i-1][j]-[i][j])*COND_COP_STL;break;caseconcrete:QN=[i-1][j]-[i][j])*COND_COP_CON;break;caseiso834:QN=2*[i-1][j]-[i][j])*COND_COPPER;break;caseadiabat:QN=0;break;default:std::cout<<i-1<<”,"<<j<<"结点属性未指定\n”;getchar();}switch[i+1][j])//判断下侧属性，并计算QS{casecopper:QS=[i+1][j]-[i][j])*COND_COPPER;break;caseair:QS=[i+1][j]-[i][j])*CONV*DX;break;casesteel:QS=[i+1][j]-[i][j])*COND_COP_STL;break;caseconcrete:QS=[i+1][j]-[i][j])*COND_COP_CON;break;caseiso834:QS=2*[i+1][j]-[i][j])*COND_COPPER;break;caseadiabat:QS=0;break;default:std::cout<<i+1<<”,"<<j<<"结点属性未指定\n”;getchar();}break;}{switch[i][j+1])//判断右侧属性，并计算QE{casecopper:QE=[i][j+1]-[i][j])*COND_COP_STL;break;caseair:QE=[i][j+1]-[i][j])*CONV*DY;break;casesteel:QE=[i][j+1]-[i][j])*COND_STEEL;break;caseconcrete:QE=[i][j+1]-[i][j])*COND_CON_STL;break;caseiso834:QE=2*[i][j+1]-[i][j])*COND_STEEL;break;caseadiabat:QE=0;break;default:std::cout<<i<<”,"<<j+1<<"结点属性未指定\n”;getchar();}switch[i][j-1])//判断左侧属性，并计算QW{casecopper:QW=[i][j-1]-[i][j])*COND_COP_STL;break;caseair:QW=[i][j-1]-[i][j])*CONV*DY;break;casesteel:QW=[i][j-1]-[i][j])*COND_STEEL;break;caseconcrete:QW=[i][j-1]-[i][j])*COND_CON_STL;break;caseiso834:QW=2*[i][j-1]-[i][j])*COND_STEEL;break;caseadiabat:QW=0;break;default:std::cout<<i<<”,"<<j-1<<"结点属性未指定\n”;getchar();}switch[i-1][j])//判断上侧属性，并计算QN{casecopper:QN=[i-1][j]-[i][j])*COND_COP_STL;break;caseair:QN=[i-1][j]-[i][j])*CONV*DX;break;casesteel:QN=[i-1][j]-[i][j])*COND_STEEL;break;caseconcrete:QN=[i-1][j]-[i][j])*COND_CON_STL;break;caseiso834:QN=2*[i-1][j]-[i][j])*COND_STEEL;break;caseadiabat:QN=0;break;default:std::cout<<i-1<<”,"<<j<<"结点属性未指定\n”;getchar();}switch[i+1][j])//判断下侧属性，并计算QS{casecopper:QS=[i+1][j]-[i][j])*COND_COP_STL;break;caseair:QS=[i+1][j]-[i][j])*CONV*DX;break;casesteel:QS=[i+1][j]-[i][j])*COND_STEEL;break;caseconcrete:QS=[i+1][j]-[i][j])*COND_CON_STL;break;caseiso834:QS=2*[i+1][j]-[i][j])*COND_STEEL;break;caseadiabat:QS=0;break;default:std::cout<<i+1<<”,"<<j<<"结点属性未指定\n”;getchar();}break;}{switch[i][j+1])//判断右侧属性，并计算QE{casecopper:QE=[i][j+1]-[i][j])*COND_COP_CON;break;caseair:QE=[i][j+1]-[i][j])*CONV*DY;break;casesteel:QE=[i][j+1]-[i][j])*COND_CON_STL;break;caseconcrete:QE=[i][j+1]-[i][j])*COND_CONCRETE;break;caseiso834:break;caseadiabat:QE=0;break;default:std::cout<<i<<”,"<<j+1<<"结点属性未指定'n";getchar();}switch[i][j-1])//判断左侧属性，并计算QW{casecopper:break;caseair:QW=[i][j-1]-[i][j])*CONV*DY;break;casesteel:break;caseconcrete:QW=[i][j-1]-[i][j])*COND_CONCRETE;break;caseiso834:break;caseadiabat:QW=0;break;default:std::cout<<i<<”,"<<j-1<<"结点属性未指定\n”;getchar();}switch[i-1][j])//判断上侧属性，并计算QN{casecopper:QN=[i-1][j]-[i][j])*COND_COP_CON;break;caseair:QN=[i-1][j]-[i][j])*CONV*DX;break;casesteel:QN=[i-1][j]-[i][j])*COND_CON_STL;break;caseconcrete:QN=[i-1][j]-[i][j])*COND_CONCRETE;break;caseiso834:break;caseadiabat:QN=0;break;default:std::cout<<i-1<<”,"<<j<<"结点属性未指定\n”;getchar();}switch[i+1][j])//判断下侧属性，并计算QS{casecopper:QS=[i+1][j]-[i][j])*COND_COP_CON;break;caseair:QS=[i+1][j]-[i][j])*CONV*DX;break;casesteel:QS=[i+1][j]-[i][j])*COND_CON_STL;break;caseconcrete:QS=[i+1][j]-[i][j])*COND_CONCRETE;break;caseiso834:break;caseadiabat:QS=0;break;default:std::cout<<i+1<<”,"<<j<<"结点属性未指定\n”;getchar();}break;}{[i][j]=T_AIR;break;}default:break;}[i][j]=(QE+QN+QS+QW)/[i][j]*[i][j]*DX*DY);}}〃函数会判断每个结点相邻4个结点的属性，然后将四个结点传来的热量加和，〃最后根据能量守恒计算出中间结点的温度值(时间全显式)voidwrite(intt,Grid&A){charfilename[40];intx=1,y=1;sprintf_s(filename,"temperature%",t);std::ofstreamfout(filename);fout<<"X,Y,TEMPERATURE\n”;for(inti=GRID_ROW-2;i>0;--i){for(intj=1;j<GRID_COLUMN-1;++j){if(y>100)y=1;fout<<x<<"