一种用于求解瞬态热传导问题的时空数值方法与流程

本发明涉及陶瓷器件烧结、铸件凝固、炸药爆炸等含热工程领域温度场计算，具体涉及到一种用于求解瞬态热传导问题的时空数值方法，适用于二维、三维瞬态热传导工程问题。

背景技术：

近年来工程领域里涉及到温度场的计算越来越多，例如陶瓷器件烧结过程中温度场的数值模拟及仿真、弹道导弹弹头表面的温度场分布、炸药爆源的爆炸场热效应分析等工程领域问题都涉及到热力学的计算。知道温度场的分布对于这些工程问题尤为重要，有时即便是微小的温度误差都会造成工程灾难。因此，研究工程问题中温度场的分布尤为重要。

目前，通用的网格类数值方法主要有有限差分法、有限元法和边界元法。这些方法对于时间依赖性的偏微分方程，首先是通过时间差分等方法去除时间偏导数之后进行求解。但这些方法通常会消耗大量的时间，而且有时数学处理复杂。本发明提出的时空联合方法不需要对时间进行离散，而是把时间看作成一个空间维度，因此节省了大量时间，对研究含热工程问题具有重要的现实指导意见。

技术实现要素：

本发明的目的在于提供一种用于求解瞬态热传导问题的时空数值方法，通过把时间看作成一个空间维度，用于简化对瞬态热传导问题的处理过程。

为实现上述目的，本发明采取以下技术方案：

一种用于求解瞬态热传导问题的时空数值方法，所研究的瞬态热传导模型符合瞬态热传导控制方程：

式中，t(x,t)描述的是在空间坐标x＝{x1,...,xn}下的温度场函数，ρ(x)表示材料的密度，c(x)表示材料的比热容，k(x)表示材料的导热系数，f(x,t)表示热源，ω表示计算区域；t为时间变量。

具体包括以下步骤：

s1：将瞬态热传导几何区域分解为由若干个离散点组成的区域；

s2：对构成区域的离散点进行分类，分别归为内部点，边界点和初始点；

s3：对区域的离散点采用广义有限差分法，计算集合区域内每一个节点的温度值。

在将瞬态热传导几何区域分解时，时间将会被看作成一个空间维度，对由时间和空间构成的几何区域进行离散，区域中的离散点既可以是规则分布，也可以是不规则分布。

在对构成区域的离散点进行分类时，位于初始时刻的所有节点被归为初始点，之后任一时刻位于边界位置的节点被归为边界点，同理，任一时刻位于区域内部的节点被归为内部点。

在采用广义有限差分法进行计算时，初始点和边界点需要满足相对应的条

件：

初始点满足初始条件：

t(x,t0)＝t0(x),x∈ω

边界点满足第一边界条件(dirichlet边界条件)：

边界点满足第二边界条件(neumann边界条件)：

边界点满足第三边界条件(robin边界条件)：

式中γ1、γ2和γ3都是ω的一段边界，γ1∪γ2∪γ3＝γ，γ是ω的总边界，t0是初始时刻，α是表面的传热系数，ta是环境温度，表示法向导数，t0(x)代表初始时刻的温度，和分别代表边界上规定的温度和热流。

本发明的有益效果在于以下：

1、通过把时间看成一个空间维度，在广义有限差分法的基础上，发明了一种新的时空联合的数值方法，用于求解瞬态热传导工程问题；

2、本发明与传统的数值方法对比，具有求解过程简单，运算效率快等优势；

3、本发明可以运用到多种场合，不仅可以解决二维和三维单一材料的温度场计算，还可以解决二维和三维功能梯度材料的温度场计算。

附图说明

图1是本发明的流程示意图；

图2是模型的时空结构域；

图3是模型节点分布图；

图4是不同时间下x轴上节点温度的数值解和解析解对比图；

图5是不同的x轴上节点在0～1s时温度的数值解和解析解对比图。

具体实施方式

下面结合附图和实施例进一步阐述本发明。

本实施例研究的是二维瞬态热传导问题，将模型已知的密度、比热容和导热系数代入到控制方程里，则得到瞬态热传导控制方程：

具体包括以下步骤：

s1：将时间看成一个空间维度，得到时空结构域如图2所示。将该时空结构域分解为由若干个离散点组成的区域，如图3所示。

s2：对构成区域的离散点进行分类，图3的底面也就是位于初始时刻的所有节点被归为初始点，图3四个侧面也就是任一时刻位于边界位置的节点被归为边界点，最后剩下的所有节点被归为内部点。

s3：采用广义有限差分法进行计算时，初始点和边界点需要满足相对应的条件：

初始点满足初始条件：

t(x,y,t0)＝t0(x,y),(x,y)∈ω

边界点满足第一边界条件(dirichlet边界条件)：

边界点满足第二边界条件(neumann边界条件)：

边界点满足第三边界条件(robin边界条件)：

式中γ1∪γ2∪γ3＝γ，γ是ω的边界，α是表面的传热系数，ta是环境温度，表示法向导数，t0(x)代表初始时刻的温度，和分别代表边界上规定的温度和热流。

本实施例是一个二维正方形模型，顶部和底部是绝热的，左右两边有热流通过，左边热流为ql＝10e^x+t，右边热流为qr＝-10e^x+t，通过计算已经得到了解析解：

t＝10e^x+t

在本实例中取9261个节点，将运用本发明得到的温度场数值结果与温度场的解析解进行对比，对比结果如图4和图5所示。本发明得到的数值解和解析解吻合非常好。

虽然本发明实施例披露如上，但本发明并非限定于此。任何本领域技术人员，在不脱离本发明的精神和范围内，均可作各种更动与修改，因此本发明的保护范围应当以权利要求所限定的范围为准。