中科大FLUENT讲稿 第七章 自定义函数

Define?Models?User-Defined Scalar…

关于用户自定义函数的具体说明可以参考相应章节。

/**************************************************************/ /* Implementation of the P1 model using user-defined scalars */

/**************************************************************/ #include \#include \#include \

/* Define which user-defined scalars to use. */ enum { P1,

N_REQUIRED_UDS， };

static real abs_coeff = 1.0; /* absorption coefficient */ static real scat_coeff = 0.0; /* scattering coefficient */

static real las_coeff = 0.0; /* linear-anisotropic scattering coefficient */ static real epsilon_w = 1.0; /* wall emissivity */

DEFINE_ADJUST(p1_adjust, domain) {

/* Make sure there are enough user defined-scalars. */ if (n_uds < N_REQUIRED_UDS)

Internal_Error(\}

DEFINE_SOURCE(energy_source, c, t, dS, eqn) {

dS[eqn] = -16.*abs_coeff*SIGMA_SBC*pow(C_T(c,t),3.);

return -abs_coeff*(4.*SIGMA_SBC*pow(C_T(c,t),4.) - C_UDSI(c,t,P1)); }

DEFINE_SOURCE(p1_source, c, t, dS, eqn) {

dS[eqn] = -abs_coeff;

return abs_coeff*(4.*SIGMA_SBC*pow(C_T(c,t),4.) - C_UDSI(c,t,P1)); }

DEFINE_DIFFUSIVITY(p1_diffusivity, c, t, i) {

return 1./(3.*abs_coeff + (3. - las_coeff)*scat_coeff); }

DEFINE_PROFILE(p1_bc, thread, position) {

face_t f;

real A[ND_ND],At;

real dG[ND_ND],dr0[ND_ND],es[ND_ND],ds,A_by_es; real aterm,alpha0,beta0,gamma0,Gsource,Ibw; real Ew = epsilon_w/(2.*(2. - epsilon_w)); Thread *t0=thread->t0;

/* Do nothing if areas aren't computed yet or not next to fluid. */ if (!Data_Valid_P() || !FLUID_THREAD_P(t0)) return; begin_f_loop (f,thread) {

cell_t c0 = F_C0(f,thread);

BOUNDARY_FACE_GEOMETRY(f,thread,A,ds,es,A_by_es,dr0); At = NV_MAG(A);

if (NULLP(T_STORAGE_R_NV(t0,SV_UDSI_G(P1)))) Gsource = 0.; /* if gradient not stored yet */ else

BOUNDARY_SECONDARY_GRADIENT_SOURCE(Gsource,SV_UDSI_G(P1), dG,es,A_by_es,1.);

gamma0 = C_UDSI_DIFF(c0,t0,P1); alpha0 = A_by_es/ds; beta0 = Gsource/alpha0; aterm = alpha0*gamma0/At;

Ibw = SIGMA_SBC*pow(WALL_TEMP_OUTER(f,thread),4.)/M_PI; /* Specify the radiative heat flux. */

F_PROFILE(f,thread,position) =aterm*Ew/(Ew + aterm)*(4.*M_PI*Ibw - C_UDSI(c0,t0,P1)

+ beta0);

}

end_f_loop (f,thread) }

函数DEFINE_ADJUST用于检验用户自定义标量输运方程是否已经激活。由于只有一个用户自定义输运方程，引用时为C_UDSI(cell,thread,0)，所以在前面枚举类型定义中使P1为零，这样就可以使用C_UDSI(cell,thread,P1)引用用户自定义的第一个标量输运方程求解的标量。

加入用户自定义输运方程之后，UDFs的函数p1_diffusivity可以在Materials面板选择，p1_source和energy_source可以在Fluid面板选择，p1_bc可以在相应的边界条件面板中选择。

P1_source 中宏SIGMA_SB给出Stefan-Boltzmann常数，σ＝5.672×10－8

W/m2K4；p1_bc中宏FLUID_THREAD_P(t0)检测当前的网格线是否属于流体区域； BOUNDARY_FACE_GEOMETRY(f,thread,A,ds,es,A_by_es,dr0)计算网格和各面(face area,cell,face)的质心的坐标，以及网格和各面质心之间的距离；宏

NULLP(T_STORAGE_R_NV(t0,SV_UDSI_G(P1)))用于检查用户自定义标量输运方程计算的标量梯度是否已经存在。

7.5.4 离散相模型

本节包含FLUENT离散相模型的三个UDFs，由于在每一步计算中都要使用这些UDFs，所以采用Compiled型以节约时间。本节的三个实例为：

1．沿颗粒轨道的熔化指数（Melting Index） 2．带电颗粒的磁场力 3．颗粒拉力系数

7.5.4.1 沿颗粒轨道的熔化指数（Melting Index）

本例使用宏DPM_SCALAR_UPDATE定义的函数，计算沿颗粒轨道的熔化指数（Melting Index）：

meltingindex??t10?dt

宏DEFINE_INIT函数初始化颗粒变量，宏DPM_OUTPUT输出特定面上的熔化指数（Melting Index），宏NULL_P，((p) = =NULL)，检测参数是否为空值。

/***********************************************************************/ /* UDF for computing the melting index along a particle trajectory */

/***********************************************************************/ #include \#include \#include \static real viscosity_0;

DEFINE_INIT(melt_setup, domain) {

/* if memory for the particle variable titles has not been allocated yet, do it now */ if (NULLP(user_particle_vars)) Init_User_Particle_Vars(); /* now set the name and label */

strcpy(user_particle_vars[0].name,\strcpy(user_particle_vars[0].label,\}

/* update the user scalar variables */

DEFINE_DPM_SCALAR_UPDATE(melting_index, cell, thread, initialize, p) {

cphase_state_t *c = &(p->cphase); if (initialize) {

/* this is the initialization call, set:

* p->user[0] contains the melting index, initialize to 0

* viscosity_0 contains the viscosity at the start of a time step */ p->user[0] = 0.;

viscosity_0 = c->mu;

} else {

/* use a trapezoidal rule to integrate the melting index */ p->user[0] += P_DT(p) * .5 * (1/viscosity_0 + 1/c->mu); /* save current fluid viscosity for start of next step */ viscosity_0 = c->mu; }

}

/* write melting index when sorting particles at surfaces */

DEFINE_DPM_OUTPUT(melting_output, header, fp, p, thread, plane) {

char name[100]; if (header) {

if (NNULLP(thread))

cxprintf(fp,\else

cxprintf(fp,\cxprintf(fp,\\

\} else {

sprintf(name,\cxprintf(fp,

\\p->state.pos[0], p->state.pos[1], p->state.pos[2], p->state.V[0], p->state.V[1], p->state.V[2],

p->state.diam, p->state.temp, p->flow_rate, p->state.time, p->user[0], name); } }

7.5.4.2 带电颗粒的磁场力

本例计算带电颗粒的磁场力。带电颗粒从上游进入层流流动的流体内，然后在t＝tstart时进入下游的磁场区域，做近似圆周运动（并不完全是圆周运动，因为颗粒所受的磁场力大小是变化的）。