用UDS求解换热方程（第三类边界条件，附代码）

正文共：2435 字 15 图 预计阅读时间： 7 分钟

1 引言

本案例将详细讲解如何用UDS求解换热问题，当然了，FLUENT打开能量模型自然就可以求解换热问题了。另外，本案例讨论的换热问题为非流固耦合换热，即计算域只有流体。本案例的宗旨在于进一步学习UDS求解变量的方法。

FLUENT UDS通用方程如下：

对流换热能量守恒方程如下：

边界条件：

2 问题分析

分别对纯导热问题和对流换热问题进行计算，边界条件取第三类。

2.1纯导热问题

以不可压缩二维稳态无源纯导热问题为例，如图1所示，一个1×0.1m2的方形平面空间，上下面对流换热边界，对流换热系数5W/m2℃，外部流体温度10℃，左边为恒温壁面（21℃），右边温度60℃。为方便问题分析，内部介质的导热系数取1W/m℃，此时只需要在Equations中关闭流动相关方程即可。不可压缩二维稳态无源纯导热方程如式2-1，与通用UDS方程对比可知，没有时间项，没有对流项，没有源项，扩散系数为λ/cp（对于第一类边界条件问题，扩散系数可以为1，即导热方程进行约分简化，也就是说，对于所有边界温度都已知的情况，导热系数对于温度分布没有影响）。因此，只需要在材料中设置扩散系数为λ/cp即可（可以用UDF编写）。

图1 问题描述

 采用FLUENT软件自带能量方程及边界条件直接进行计算，温度结果如图2所示，采用UDS计算的温度结果如图3所示，可以看出计算结果几乎完全一致。关于UDS第三类边界条件参考笔者之前关于第二类第三类边界条件转换成第一类的案例，采用DEFINE_PROFILE函数（代码见附录）定义边界条件（Specified Value）；

图2 FLUENT软件计算结果

图3 UDS计算结果

2.2 对流换热问题

同样采用以上的计算域，考虑对流换热。不可压缩二维稳态无源对流换热方程如式2-2，与通用UDS方程对比可知，没有时间项，没有源项，扩散系数为λ/cp（UDF编写），对流项为标准质量流率，也可以用DEFINE_UDS_FLUX函数自定义（代码见附录）

①采用标准k-e模型。    

②边界条件：入口速度2m/s，温度60℃，UDS值60（Specified Value）；

③上下壁面：wall，对流换热边界，对流换热系数5W/m2℃，外部流体温度10℃；参考笔者之前关于第二类第三类边界条件转换成第一类的案例，采用DEFINE_PROFILE函数（代码见附录）定义边界条件（Specified Value）；

④出口：outflow，UDS值0（Specified Flux）；（假定流动充分发展，流动方向上的梯度为0，详见数值传热学P232）

采用FLUENT软件自带边界条件直接进行计算，结果如图4所示，采用UDS计算的温度结果如图5所示。可以看出，温度分布几乎一致，但从中心线温度分布线看，两者存在一定误差，特别是越靠近出口，误差越大。笔者认为这是因为采用标准k-e湍流模型，而UDS的扩散系数是对NS方程而言的，湍流流动的扩散系数与NS方程下的扩散系数不一样。

图4 FLUENT软件自带能量方程及边界条件计算结果

图5 UDS计算结果

提示

    笔者尝试用DEFINE_PROFILE函数直接定义第三类UDS边界条件，即Specified Flux为所求UDS本身的函数关系，发现无法收敛，无法获得计算结果。因此才想到将第三类边界条件转换成第一类来处理。

参考文献：

[1]FLUENT帮助文档

[2]杨世铭,陶文铨.传热学(第4版)[M].北京:高等教育出版社

[3]陶文铨.数值传热学[M].西安:西安交通大学出版社

附录 UDF代码

#include "udf.h"
#include "prop.h"
#include "mem.h"
DEFINE_UDS_FLUX(MyFlux,f,t,i)
{
 return F_FLUX(f,t);
}
DEFINE_DIFFUSIVITY(MyDiff,c,t,i)
{
real lamda=C_K_L(c,t); //获取导热系数
real cp=C_CP(c,t); //获取比热
return lamda/cp;
}
DEFINE_PROFILE(UDSProfile3,thread,index)
{
real coef=5.0;  //对流换热系数
real tout=10;   //外部流体温度
real T;
face_t f;
cell_t c0;
Thread *t0;
cell_t c;
real thickness;
real thick;
real x[ND_ND];
real x0[ND_ND];
real NV_VEC(a);
real NV_VEC(A);
real cosphi;
begin_f_loop(f,thread)
{
t0=THREAD_T0(thread); //获取边界邻近网格的线
c0=F_C0(f,thread);  //获取边界邻近网格单元的ID
F_CENTROID(x,f,thread); //获取面心
C_CENTROID(x0,c0,t0);   //获取邻近单元重心
NV_VV(a, =, x0, -, x); //求面心和重心的向量(两心向量）
thick=NV_MAG(a);  //求面心和重心的距离（非法向）
F_AREA(A,f,thread);  //获取边界法向面积矢量
cosphi=fabs(NV_DOT(a,A))/NV_MAG(a)/NV_MAG(A); //求两心向量和法向的夹角余弦
thickness=thick*cosphi;  //求两心法向距离
T=(C_UDSI(c0,t0,0)+coef*thickness*tout/C_K_L(c0,t0))/(1+coef*thickness/C_K_L(c0,t0));
F_PROFILE(f,thread,index)=T;    
}  
end_f_loop(f,thread)
}