计算科学 - 确定气缸内稳态温度的有限体积法 - 吾爱随笔录

我试图在 Matlab 中使用有限体积法 (FVM) 对圆柱体的稳态行为进行建模，该方法在 Matlab 中受到各种边界条件的影响。首先，我将圆柱体视为轴对称的，因此我只确定 rZ 平面中的温度分布。我设置了一个 2D 网格，使 r 坐标位于“x”轴上，Z 坐标位于“y”轴上。由于圆柱体是轴对称的，因此 r = 0 处的热通量为 0（绝缘）。我还在 Z=0 的底面上应用绝缘边界条件。侧面有对流，对流系数恒定，空气温度为 310 K。圆柱顶部，z=Z，有 1000 W/m^2 的热通量，来自 310K 空气的对流，还有表面对环境辐射，周围环境为 298 K。气缸的初始温度为 298 K（初始猜测值）。我遇到的问题是温度分布与 Comsol 的结果相匹配，网格为 40x40，热导率为 1。当我将电池数量增加到 50x50、60x60 等时，或者如果我将热导率更改为其他值大于 1 则温度曲线仍然是应有的温度，但所有温度都明显低于预期。我在 Matlab 中的控制循环如下：或者，如果我将热导率更改为 1 以外的值，那么温度曲线仍然是应该的，但所有温度都明显低于预期。我在 Matlab 中的控制循环如下：或者，如果我将热导率更改为 1 以外的值，那么温度曲线仍然是应该的，但所有温度都明显低于预期。我在 Matlab 中的控制循环如下：

使用初始值和问题常数求解顶部和侧面的初始表面温度。
使用表面温度和问题几何来确定网格中每个单元的 aE、aN、aS、aW、aP 和 b。
在两个方向上逐行使用 TDMA 来确定温度曲线。
计算新温度曲线与前一个温度曲线之间的误差。
由于使用先前迭代值的非线性辐射项，更新表面温度和 aP 和 b。
再次逐行执行 TDMA 并检查收敛。如果误差不小于 0.001，则重复该过程。

对于边界条件，我对顶面和侧面使用了能量平衡。顶面包括表面对环境的辐射，我计算出的表面温度为

T_{b} = \frac{α G + h T_{f} + \frac{k T_{p}}{δ z} + ϵ σ (T_{s u r}^{4} + 3 T_{b}^{* 4})}{h + 4 ϵ σ T_{b}^{* 3} + \frac{k}{δ z}}

$T_{b}=\frac{\alpha G+hT_{f}+\frac{kT_{p}}{\delta z}+\epsilon \sigma(T_{sur}^{4}+3T_{b}^{*4})}{h+4 \epsilon \sigma T_{b}^{*3}+\frac{k}{\delta z}}$

我也对通过对流加热的一侧的表面温度做同样的事情。由于顶面温度取决于其先前的迭代值，因此我对顶面使用初始猜测值 298。然后我根据它和我为第一次迭代提供的方程计算初始表面温度，

问题在于表面温度取决于热导率和网格间距。当我使用许多像 80x80 这样的电池时，第一次 TDMA 迭代的顶部表面温度几乎是 298，这很难相信，因为它被 1000 W/m^2 热通量加热，并且它也被 310 K 空气温度加热. 这导致第一次迭代的误差非常小，因此收敛只会产生 325 的最大温度，而它应该是 450。但是如果我使用像 40x40 的网格，那么第一次迭代的表面温度要高得多，并且收敛导致最终的温度曲线看起来几乎与预期的完全一样。

同样的事情发生在热导率上，如果我改变这个值，那么初始表面温度会发生巨大变化，这会导致收敛值远低于预期。

似乎这一切都与初始表面温度有关，但我不确定如何处理它。任何帮助都会令人惊叹，因为我花了很多个晚上直到凌晨 5 点试图弄清楚如何解决这个问题。我的联系电子邮件在我的个人资料中，如果您想查看我的 matlab 代码，或者您想更多地与我讨论我的问题。