一维非稳态常系数热传导方程(第一类边界条件)显式解法
内容来源:B站计算传热学大叔
链接: https://www.bilibili.com/video/BV1Wk4y1R73p.
C++代码:
/*
所解方程为一维常系数方程,边界条件为两侧都是第一类边界条件
*/
#include <iostream>
#include<iomanip>
#include <fstream>
using namespace std;
int main()
{
ofstream file;
file.open("1.txt", ios::out | ios::trunc);
file << fixed;
double length; //导热体的长度
double TL, TR; //左侧温度,右侧温度
double den, c, k;//物质密度,比热,导热系数
double s;//热源
double dt;//时间间隔
int n;//控制体的个数
int tstep, maxtstep;//时间步,最大时间步
//-----------------------------------------------------
double dx;
double T0[11], T[11];//从0-10共十一个元素,即左右边界两个体积为零的控制体和中间9个控制体,内部的元素对应11个温度的值
double ae0[11], aw0[11], ap0[11], ap1[11], b[11];//系数
int i;
length = 3;
TL = 3;
TR = 5;
den = 100;
c = 1000;
k = 10;
s = 10;
dt = 100;
n = 3; //当n值取得很大时,t值不变会引起非物理现象,即ap0为负数 如果n想取大于10的数需要改T0 T 以及系数
maxtstep = 3000;
//-------------------------------------------------
dx = length / n;
for (i = 2; i < n;i++)
{
ae0[i] = k / dx;
aw0[i] = k / dx;
ap0[i] = den * c * dx / dt - ae0[i] - aw0[i];//如果n值很大可能出现负数的情况,非物理现象,解决办法dt减小。
ap1[i] = ap0[i] + ae0[i] + aw0[i];
b[i] = s * dx;
}
i = 1;
{
ae0[i] = k / dx;
aw0[i] = k / (dx / 2);
ap0[i] = den * c * dx / dt - ae0[i] - aw0[i];
ap1[i] = ap0[i] + ae0[i] + aw0[i];
b[i] = s * dx;
}
i = n;
{
ae0[i] = k / (dx / 2);
aw0[i] = k / dx;
ap0[i] = den * c * dx / dt - ae0[i] - aw0[i];
ap1[i] = ap0[i] + ae0[i] + aw0[i];
b[i] = s * dx;
}
for (i = 1; i <= n; i++)
{
T0[i] = 3;
}
T0[0] = TL;
T0[n + 1] = TR;
for (tstep = 1; tstep <= maxtstep; tstep++)
{
cout << tstep << '\t';
for (i = 1; i <= n; i++)
{
T[i] = (ae0[i] * T0[i + 1] + aw0[i] * T0[i - 1] + ap0[i] * T0[i] + b[i]) / ap1[i];
cout << setprecision(5) << T[i] << '\t';
}
cout << '\n';
for (i = 1; i <= n; i++)
{
file << setprecision(5) << T[i] << '\t';
}
file << '\n';
for (i = 1; i <= n; i++)
{
T0[i] = T[i];
}
}
return 0;
}
特别感谢数值传热学大叔的网课,也欢迎大家找我交流讨论,我也在努力学习数值传热学。