一维非稳态常系数热传导方程（第一类边界条件）

一维非稳态常系数热传导方程（第一类边界条件）显式解法

内容来源：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;
}

特别感谢数值传热学大叔的网课，也欢迎大家找我交流讨论，我也在努力学习数值传热学。