2010年十月写的旧代码。
第一类边界条件是给定边界温度。
第二类是对流边界。
区域都是如下形状的:
--------------------------------
| |
| ----------------------
| |
| |
| |
| |
--------
用C++纯属蛋疼。
第一类边界条件:
/*
* 等温边界
*/
#include<iostream>
#include<cmath>
using namespace std;
const double out_temp=30.0;//外边界温度
const double in_temp=0.0;//内边界温度
const double accuracy=0.00000001;//精度
const double lambda=0.53;//导热系数
const int width=15;//上部点阵宽度
const int width_bottom=4;//下部点阵宽度
const int height=4;//上部点阵高度
const int height_bottom=7;//下部点阵高度
/*
//另一组参数
const int width=16;
const int width_bottom=6;
const int height=6;
const int height_bottom=6;
*/
//总共点数
const int num_of_points=width*height+width_bottom*height_bottom;
//单个点
class point
{
public:
double temp;//温度
int up;//数组下标
int down;
int left;
int right;
point(){
temp=0.0;//初始化成0摄氏度
up=down=right=left=0;
}
};
ostream & operator<<(ostream & src_stream,const point & src){
src_stream<<"temp="<<src.temp;
src_stream<<" up="<<src.up<<" down="<<src.down;
src_stream<<" left="<<src.left<<" right="<<src.right;
return src_stream;
}
void print_grid(point * points){//输出网格各点的温度
cout<<endl;
for(int position=0;position>=0;position=points[position].down){
for(int tmp=position;tmp>=0;tmp=points[tmp].right){
cout.width(10);
cout<<points[tmp].temp;
}
cout<<endl<<endl<<endl<<endl;
}
}
double calc_direction(point * points,int direction,double &opposit_weight){
//计算给定方向邻点的温度和反方向的权重。
switch (direction) {
case -1:
return out_temp;
break;
case -2:
return in_temp;
break;
case -3:
opposit_weight*=2;
return 0.0;
break;
default:
return points[direction].temp;
}
}
void compu_point(point * points,int now){
//根据周围四个点算出指定点的温度
double left_temp;
double right_temp;
double up_temp;
double down_temp;
double up_wight=0.25;//上方权重,默认0.25
double down_wight=0.25;
double left_wight=0.25;
double right_wight=0.25;
left_temp=calc_direction(points,points[now].left,right_wight);
right_temp=calc_direction(points,points[now].right,left_wight);
up_temp=calc_direction(points,points[now].up,down_wight);
down_temp=calc_direction(points,points[now].down,up_wight);
points[now].temp=left_temp*left_wight+right_temp*right_wight+up_temp*up_wight+down_temp*down_wight;
}
void rec_walk(point * points){
// "左=>右"里嵌"上=>下"遍历