CN108763688B - 一种静电悬浮液态合金传热的数值计算方法 - Google Patents

Abstract

本发明涉及一种静电悬浮液态合金传热的数值计算方法，涉及一种在静电悬浮条件下规则和非规则具有旋转对称结构形态液态合金传热的计算方法。利用有限体积方法数值迭代求解非规则形态液态合金内部的传热方程，得到液态合金内部的温度场及随时间的变化情况。本发明解决了目前实验测量静电悬浮液态合金滴温度困难且误差较大的问题，且传热计算模型仅得到合金液滴质点的温度及变化情况的局限，采用二维静电悬浮液态合金传热数值计算方法可获取(1)静电悬浮条件下规则及非规则具有旋转对称结构形态的合金液滴传热；(2)液态合金内部的温度分布、温度梯度及冷却速率等丰富准确的传热信息。

Description

一种静电悬浮液态合金传热的数值计算方法

技术领域

本发明属于传热计算、静电悬浮技术领域，涉及一种静电悬浮液态合金的传热数值计算方法，特别涉及一种在静电悬浮条件下规则和非规则形态合金液滴的传热计算。

背景技术

静电悬浮状态下的液态合金由于避免了与容器壁接触，可获得较大的过冷度，因而有助于深过冷液态合金的凝固机理与热物理性质研究。而实时观测及调控静电悬浮状态下液态合金的温度及变化趋势是进行上述研究的前提，因此准确获得液滴温度分布及其变化规律显得尤为重要。

近年来悬浮状态下液态合金温度测量的实验手段主要有：接触式热电偶测温和非接触式红外测温。前者把温度信号转换成热电动势信号，通过电气仪表转换成被测物质的温度，主要缺点是破坏了实验样品的悬浮无容器状态，且测量过程需要一定的热平衡反应时间，测量结果会出现一定的热延迟；后者具有无接触的优势，但受物体发射率、测量距离等外界因素影响较大，且该方法是将液态合金看作质点仅测得合金表面的温度，无法获取液态合金内部的温度场。而液态合金传热计算主要依据牛顿传热模型，将液态合金看作温度质点通过对流及辐射向外界散热，计算得到液态合金的温度变化情况，同样不能获得内部温度场信息。B.Q.Li等人模拟了静电悬浮和电磁悬浮条件下的液滴内部的流动及传热，但关注更多的是合金液滴在电场/磁场作用下内部的流动及表面变形。见文献：B.Q.Li,andS.P.Song,Thermal and Fluid Flow Aspects of Electromagnetic and ElectrostaticLevitation-A Comparative Modeling Study，Microgravity Science and Technology，1998,11(4):134-143.以及Huo Y,and Li B Q,Surface Deformation and Convection inElectrostatically-Positioned Droplets of Immiscible Liquids UnderMicrogravity,Journal of Heat Transfer,2006，128(6):520-529.

目前已有的实验手段及主要计算模型无法获得液态合金内部的温度场分布及变化情况。随着静电悬浮合金液滴的增大或冷却速率的提高，液滴内部温度分布并不均匀，传统模型所假设的均匀温度场不能客观地反映液滴内部温度分布情况。而且在重力、静电力、表面张力作用下，直径较大的合金液滴形态不再保持标准的球形，其形状的改变也会一定程度上影响合金液滴内部的温度场。二维静电悬浮液态合金传热计算方法针对在静电悬浮条件下标准球状及非规则具有旋转对称结构形态的合金液滴，专注于合金液滴内部的传热，采用有限体积方法，获取合金内部的温度场分布及变化、温度梯度、冷却速率等信息。

发明内容

要解决的技术问题

为了避免现有技术的不足之处，本发明提出一种静电悬浮液态合金传热的数值计算方法，针对目前实验测定及主要传热计算，仅获取规则合金液滴质点温度的局限，采用二维静电悬浮液态合金传热数值计算方法获得规则及非规则具有旋转对称结构形态的合金液滴内部具体温度分布及变化情况，提供更为丰富的传热信息。

技术方案

一种静电悬浮液态合金传热的数值计算方法，其特征在于：在静电悬浮条件下规则和非规则具有旋转对称结构形态合金液滴的传热计算步骤如下：

步骤1、传热控制方程和边界热交换条件：

1、传热控制方程：

其中：T₀(x,y,z；t)为温度场函数，ρ、c、λ分别为液滴密度、比热、导热系数；

2、边界热交换条件为辐射散热：

其中：T_l，T_g分别为合金液滴和外界环境温度，n为静电悬浮合金液滴的单位外法向向量，ε为静电悬浮合金液滴热辐射系数，σ为Stefan-Boltzmann常数；

步骤2：将传热控制方程通过坐标变换从三维直角坐标系变换到二维柱坐标系下：

三维直角坐标系的三维传热控制方程为：

变换到二维柱标系下的传热控制方程为：

其中温度场函数由T₀(x,y,z；t)变为T(r,z；t)；

步骤3：采用有限体积方法处理规则及非规则形态静电悬浮合金液滴的传热控制方程，包括求解区域的网格划分及相应的对偶单元剖分、传热控制方程在控制单元内的积分：

1、将静电悬浮合金液滴几何形态的正视图作为求解区域，静电悬浮合金液滴的直径为0.1～10mm，z_m/r_m＝1～1.57，其中z_m和r_m分别为液滴垂直方向和水平方向径长的最大值；

2、对求解区域进行三角网格划分以得到网格节点，确定每个网格节点的邻接点及其所在的三角形单元从而进行重心对偶剖分，以确定每个节点的控制单元；

3、将二维柱标系下的传热控制方程

用矢量形式表示为：

4、对矢量形式表示的传热控制方程在任一节点控制单元内进行积分：

节点控制单元内积分：

其中：

表示节点的控制单元；

通过高斯散度定理，将在控制单元内的体积分转化为在控制单元边界上积分：

其中：

表示节点的控制单元的边界；

步骤4：对积分方程

进行离散处理，包括积分方程左端时间项和右端扩散项的离散：

1、对左端时间项的离散，采用向前差分：

其中

和

分别为节点P₀的r坐标值和控制单元的面积，

分别表示节点

P₀的第n，n+1时刻的温度，τ为时间步长；

2、对右端扩散项的离散，分为内点和界点两种情形：

内点：

其中c_i为节点P_i处的温度T_Pi的系数：

式中

分别为节点P_i的r和z方向的坐标值，S_i为节点P₀控制单元中第i个三角形单元的面积；

界点时需代入边界条件：

其中节点系数c_i为：

为节点P₀边界相关系数：

步骤5、迭代计算静电悬浮合金液滴的瞬态传热：构造固定迭代矩阵A及需每步更新的右端项矩阵bⁿ，确定合适的网格参数及时间步长，采用隐格式迭代求解ATⁿ⁺¹＝bⁿ，直到到达瞬态传热终止时刻，获得液态合金内部温度场及随时间变化情况：

1、采用隐格式：

对于内点，由式(8)和(9)得：

得内点离散方程：

其中：

对于界点，由(8)和(10)得：

得界点离散方程：

其中：

为第n时刻边界相关系数

的值；

2、根据内点离散方程(12)和界点离散方程(14)的左端项构造固定迭代矩阵A：

对节点P₀及邻接点P₁，P₂，…，P_m，其编号为i₀,i₁,i₂,…i_m，那么矩阵A的第i₀行元素为：

矩阵A的第i₀行其他元素都为零；

3、根据内点离散方程(12)和界点离散方程(14)的右端项，构造右端项(bⁿ)_np×1，需每步更新：

若节点P₀为内点则右端项为

若节点P₀为边界点则右端项为

4、通过高斯消元法直接求解或者高斯-赛德尔迭代法求解线性方程组ATⁿ⁺¹＝bⁿ，由第n时刻悬浮液滴内部的温度场得到n+1时刻的温度场。

有益效果

本发明提出的一种静电悬浮液态合金传热的数值计算方法，涉及一种在静电悬浮条件下规则和非规则具有旋转对称结构形态液态合金传热的计算方法。利用有限体积方法数值迭代求解非规则形态液态合金内部的传热方程，得到液态合金内部的温度场及随时间的变化情况。包括以下计算流程：

(1)确定静电悬浮液态合金的传热计算模型：处于真空中的静电悬浮合金液滴内部采用热传导方程，表面采用辐射换热条件。

(2)传热控制方程的坐标变换：将传热控制方程通过坐标变换从三维直角坐标系转换到二维柱坐标系。

(3)确定数值计算方法：采用有限体积方法，包括求解区域的网格划分及相应的对偶单元剖分、传热控制方程在控制单元内的积分及离散处理。

(4)迭代计算瞬态传热：构造固定迭代矩阵A及需每步更新的右端项矩阵bⁿ，确定合适的网格参数及时间步长，采用隐格式迭代求解ATⁿ⁺¹＝bⁿ。

针对静电悬浮条件下直径为2080μm的Ni-Ti合金液滴，利用本发明的计算方法，给出了具体的实施方式，并计算得到了合金内部的温度场及随时间变化情况。

本发明的积极效果：解决了目前实验测量静电悬浮液态合金滴温度困难且误差较大的问题，且传热计算模型仅得到合金液滴质点的温度及变化情况的局限，采用二维静电悬浮液态合金传热数值计算方法可获取(1)静电悬浮条件下规则及非规则具有旋转对称结构形态的合金液滴传热；(2)液态合金内部的温度分布、温度梯度及冷却速率等丰富准确的传热信息。

附图说明

图1是二维静电悬浮液态合金传热计算流程示意图

图2是静电悬浮合金液滴形态图：(a)实验实际形态；(b)计算模拟形态。

图3求解区域的三角网格划分及相应的对偶剖分：(a)求解区域的三角网格划分；(b)内点对偶剖分；(c)界点对偶剖分。

图4(a)不同时刻液态合金中心的温度变化情况；(b)液态合金内部最大温差δT的变化情况。

图5是不同时刻液态合金内部的温度场：(a)t＝0.4s；b)t＝2.0s；(c)t＝8.0s；

具体实施方式

现结合实施例、附图对本发明作进一步描述：

1.确定液滴几何形态，并将其正视图作为计算求解区域：

通过实验得到实际形态或者计算得到模拟形态，如图2.

静电悬浮合金液滴直径0.1～10mm，z_m/r_m＝1～1.57，其中z_m和r_m分别为液滴垂直方向和水平方向径长的最大值。

2.确立传热控制方程及边界热交换条件：

(1)液滴传热的控制方程为：

其中T(x,y,z；t)为关于位置(x,y,z)和时间t的温度场函数，u为液滴内部流速，ρ、c、λ分别为液滴密度、比热、导热系数。

因为处于真空条件下静电悬浮的液滴内部流速很小，所以可忽略上式左端第二项内部对流。

(2)边界热交换条件为辐射散热：

其中n为静电悬浮合金液滴的单位外法向向量，T_l，T_g分别为合金液滴和外界环境温度，ε为热辐射系数，σ为Stefan-Boltzmann常数。

3.传热控制方程的坐标变换：

直角坐标系的三维传热控制方程为

其中温度场函数为T₀(x,y,z；t)

变换到柱标系下的传热控制方程为

其中温度场函数变为

由于静电悬浮条件下的合金液滴形态具有选择轴对称结构，与角度参量

无关，因此传热控制方程在柱标系可表示为：

其中温度场函数变为T(r,z；t)

两边同乘以r

用矢量形式可表示为

4.确定数值求解方法：

对于非规则形态的液滴，本专利采用有限体积方法数值求解传热控制方程。

5.网格划分及相应的对偶剖分：

对非规则区域进行三角网格划分以确定网格节点，如图3(a)。

确定网格节点的邻接点及所在的三角形单元，进行重心对偶剖分，以确定节点的控制单元如图3(b)和(c)。

6.传热控制方程在控制单元内积分及离散

(1)在节点P₀控制单元

内积分

其中

表示节点的控制单元。

(2)通过高斯散度定理，将在控制单元内的积分转化为在控制单元边界上积分

其中

表示节点的控制单元的边界。

(3)对积分方程(10)左端时间项的离散，采用向前差分：

其中

和

分别为节点P₀的r坐标值和控制单元的面积，

分别表示节点P₀的第n，n+1时刻的温度，τ为时间步长。

(4)对积分方程(10)右端扩散项的离散

①三角形单元内温度函数可插值表示为

T(r,z)＝L_iT_i+L_jT_j+L_kT_k (26)

式中T_i，T_j，T_k分别为节点i，j，k的温度，L_i，L_j，L_k为面积坐标系数

S为三角形单元面积，r_i，r_j，r_k，z_i，z_j，z_k，分别为三角形顶点坐标。

②梯度算子的离散

梯度算子的离散

③对积分方程(10)右端扩散项的离散，分为内点和界点两种情形：

对于内点求节点系数c_i

式中M_i为P₀P_i边的中线，Q_i为第i个三角形单元的重心，如图3(b)。

节点系数c_i：

而对于界点，需代入边界条件，求节点系数c_i,及边界相关系数

式中节点系数c_i为

边界相关系数

为

7.离散方程的迭代格式

本发明采用隐格式，稳定性好

对于内点，由式(11)和(16)得

整理后得内点离散方程

其中

对于界点，由式(11)和(18)得

整理后得界点离散方程

其中

为第n时刻边界相关系数

的值。

8.根据内点离散方程(22)和界点离散方程(24)的左端项，构造固定迭代矩阵A。

对节点P₀及邻接点P₁，P₂，…，P_m，其编号为i₀,i₁,i₂,…i_m，那么矩阵A的第i₀行：

矩阵A的第i₀行其他元素都为零。

那么对于np个节点构成的迭代矩阵为A_np×np。

9.根据内点离散方程(22)和界点离散方程(24)的右端项，构造右端项(bⁿ)_np×1，需每步更新。

若节点P₀为内点则右端项为

若节点P₀为边界点则右端项为

三角网格单元数nt：200～3000，空间步长dx:5×10^-3～10^-1mm，时间步长τ：1～10^{- 4}s。

10.求解线性方程组ATⁿ⁺¹＝bⁿ，由第n时刻悬浮液滴内部的温度场得到n+1时刻的温度场。

通过高斯消元法直接求解或者高斯-赛德尔迭代法求解。

传热计算结果

选择了半径为1040μm的Ni₅₅Ti₄₅液滴作为计算对象，初始温度设置为1630K，通过辐射散热，进行冷却。悬浮液滴形变量z_m/r_m＝1.02。

计算实例具体参数设置：

计算结果：如图4(a)所示，液滴中心温度随冷却时间变化曲线，可以看出由于液滴与外界温差的减小，冷却速率慢慢降低；图4(b)得到了不同时刻液滴内部的最大温差变化情况；图5给出了三个不同时刻液滴内部的温度场。

Claims (1)

1.一种静电悬浮液态合金传热的数值计算方法，其特征在于：在静电悬浮条件下规则和非规则具有旋转对称结构形态合金液滴的传热计算步骤如下：

步骤1:传热控制方程和边界热交换条件：

①传热控制方程：

其中：T₀(x,y,z；t)为温度场函数，ρ、c、λ分别为液滴密度、比热、导热系数；

②边界热交换条件为辐射散热：

其中：T_l，T_g分别为合金液滴和外界环境温度，n为静电悬浮合金液滴的单位外法向向量，ε为静电悬浮合金液滴热辐射系数，σ为Stefan-Boltzmann常数；

步骤2：将传热控制方程通过坐标变换从三维直角坐标系变换到二维柱坐标系下：

三维直角坐标系的三维传热控制方程为：

变换到二维柱标系下的传热控制方程为：

其中温度场函数由T₀(x,y,z；t)变为T(r,z；t)；

步骤3：采用有限体积方法处理规则及非规则形态静电悬浮合金液滴的传热控制方程，包括求解区域的网格划分及相应的对偶单元剖分、传热控制方程在控制单元内的积分：

①将静电悬浮合金液滴几何形态的正视图作为求解区域，静电悬浮合金液滴的直径为0.1～10mm，z_m/r_m＝1～1.57，其中z_m和r_m分别为液滴垂直方向和水平方向径长的最大值；

②对求解区域进行三角网格划分以得到网格节点，确定每个网格节点的邻接点及其所在的三角形单元从而进行重心对偶剖分，以确定每个节点的控制单元；

③将二维柱标系下的传热控制方程

用矢量形式表示为：

④对矢量形式表示的传热控制方程在任一节点控制单元内进行积分：

节点控制单元内积分：

其中：

表示节点的控制单元；

通过高斯散度定理，将在控制单元内的体积分转化为在控制单元边界上积分：

其中：

表示节点的控制单元的边界；

步骤4：对积分方程

进行离散处理，包括积分方程左端时间项和右端扩散项的离散：

①对左端时间项的离散，采用向前差分：

其中

和

分别为节点P₀的r坐标值和控制单元的面积，

分别表示节点P₀的第n+1 ，n时刻的温度，τ为时间步长；

②对右端扩散项的离散，分为内点和界点两种情形：

内点：

其中c_i为节点P_i处的温度

的系数：

式中

分别为节点P_i的r和z方向的坐标值，S_i为节点P₀控制单元中第i个三角形单元的面积；

界点时需代入边界条件：

其中节点系数c_i为：

为节点P₀边界相关系数：

步骤5:迭代计算静电悬浮合金液滴的瞬态传热：构造固定迭代矩阵A及需每步更新的右端项矩阵bⁿ，确定合适的网格参数及时间步长，采用隐格式迭代求解ATⁿ⁺¹＝bⁿ，直到到达瞬态传热终止时刻，获得液态合金内部温度场及随时间变化情况：

①采用隐格式：

对于内点，由式(8)和(9)得：

得内点离散方程：

其中：

对于界点，由(8)和(10)得：

得界点离散方程：

其中：

为第n时刻边界相关系数

的值；

②根据内点离散方程(12)和界点离散方程(14)的左端项构造固定迭代矩阵A：

对节点P₀及邻接点P₁，P₂，…，P_m，其编号为i₀,i₁,i₂,…i_m，那么矩阵A的第i₀行元素为：

矩阵A的第i₀行其他元素都为零；

③根据内点离散方程(12)和界点离散方程(14)的右端项，构造右端项(bⁿ)_np×1，需每步更新：

若节点P₀为内点则右端项为

若节点P₀为边界点则右端项为

④通过高斯消元法直接求解或者高斯-赛德尔迭代法求解线性方程组ATⁿ⁺¹＝bⁿ，由第n时刻悬浮液滴内部的温度场得到n+1时刻的温度场。

Images