/************************************************** **********************/ /* Porous Media Source Code */ /* */ /* By Greg Perkins */ /* CANCES, Australian Technology Park */ /* Ph: 02 9318 0004, Fax: 02 9319 2328 */ /* Email: perkinsg@cances.atp.com.au */ /* */ /* Started: 10-11-99 */ /* */ /* Revised: 22-11-99 */ /************************************************** **********************/
/* ---- These routines implement a porous media model (same as standard
Fluent model) that can vary spatially and temporily.
To use in Fluent compile and link source library as described
in the Fluent UDF manual. Then use the udf in the source terms
for momentum equations */
#include "udf.h" #include "sg.h"
/* ---- use a user defined scalar to track material properties of the domain */ enum { PERM, N_REQUIRED_UDS };
/* --- misc defines */ #define PERM_EXP_VOID 35.0 #define PERM_EXP_MATERIAL -20.0 #define PERM_ADJ_FALSE 0 #define PERM_ADJ_TRUE 1 #define UDS_PERM(alpha) log10(alpha) /* --- permeability function for scalar-0 */ #define WEIGHT 1.0e20 /* --- weighting co-eff in linearized eqn */ #define WEIGHT2 1.0e20 #define INLET_VELOCITY 10.0 /* --- inlet velocity in m/s */
/* ---- define criteria for various possible configurations */ #define MAT_CHANNEL(x,y,z) ((y<=0.25) || (y>=0.75)) #define MAT_STEP(x,y,z) ( ((x<=1.0)&&(y<=0.5)) || ((x>1.0)&&(y<=0.25)) ) #define MAT_ONECHANNEL(x,y,z) (y<=0.5) #define MAT_CHANNEL2(x,y,z) ( ((x<=1.0)&&(y<=0.5)) || ((x>1.0)&&(y<=0.25)) || (y>=0.95) )
/* ------------------------------------------------------------------------
Material Properties
This routine returns the local properties of the material (alpha,C2) for
a given location in the domain (x,y,z).
------------------------------------------------------------------------ */ void Material_Properties(real x, real y, real z, real *alpha, real *C2) {
/* --- use this routine to return the permeability, alpha and the
co-efficient C2, for each location x,y,z in the domain. No
Modifications needed for 2D. For time dependent porous media
use the RP_Get_Real("flow time") function to find out t (secs) */
if MAT_STEP(x,y,z)
{
*alpha = pow(10.0,PERM_EXP_MATERIAL);
*C2 = 0.0;
}
else
{
*alpha = pow(10.0,PERM_EXP_VOID);
*C2 = 0.0;
} }
/* ------------------------------------------------------------------------
X_Momentum_Source
This routine returns the source term for the X-momentum term for each
control volume in the domain. The local properties are obtained by
calling Material_Properties
------------------------------------------------------------------------ */ DEFINE_SOURCE(SRCE_Xmom,cell,thread,dS,eqn) {
real x[ND_ND];
real alpha, C2, constant1, constant2, Ux, source;
/* --- determine x,y,z co-ordinates of cell */
C_CENTROID(x,cell,thread);
/* --- determine local properties */
Material_Properties(x[0], x[1], x[2], &alpha, &C2);
/* --- determine constants 1,2 */
constant1 = C_MU_L(cell,thread)/alpha;
constant2 = 0.5 * C_R(cell,thread) * C2;
/* --- determine x-velocity */
Ux = C_U(cell,thread);
source = -(constant1*Ux + constant2 * fabs(Ux) * Ux);
dS[eqn] = -(constant1 + 2 * constant2 * fabs(Ux)); /* XXX CHECK */
return source; }
/* ------------------------------------------------------------------------
Y_Momentum_Source
This routine returns the source term for the Y-momentum term for each
control volume in the domain. The local properties are obtained by
calling Material_Properties
------------------------------------------------------------------------ */ DEFINE_SOURCE(SRCE_Ymom, cell, thread, dS, eqn) {
real x[ND_ND];
real alpha, C2, constant1, constant2, Uy, source;
/* --- determine x,y,z co-ordinates of cell */
C_CENTROID(x,cell,thread);
/* --- determine local properties */
Material_Properties(x[0], x[1], x[2], &alpha, &C2);
/* --- determine constants 1,2 */
constant1 = C_MU_L(cell,thread)/alpha;
constant2 = 0.5 * C_R(cell,thread) * C2;
/* --- determine y-velocity */
Uy = C_V(cell,thread);
source = -(constant1*Uy + constant2 * fabs(Uy) * Uy);
dS[eqn] = -(constant1 + 2 * constant2 * fabs(Uy)); /* XXX CHECK */
return source; }
/* ------------------------------------------------------------------------
Z_Momentum_Source
This routine returns the source term for the Z-momentum term for each
control volume in the domain. The local properties are obtained by
calling Material_Properties
------------------------------------------------------------------------ */ DEFINE_SOURCE(SRCE_Zmom, cell, thread, dS, eqn) {
real x[ND_ND];
real alpha, C2, constant1, constant2, Uz, source;
/* --- determine x,y,z co-ordinates of cell */
C_CENTROID(x,cell,thread);
/* --- determine local properties */
Material_Properties(x[0], x[1], x[2], &alpha, &C2);
/* --- determine constants 1,2 */
constant1 = C_MU_L(cell,thread)/alpha;
constant2 = 0.5 * C_R(cell,thread) * C2;
/* --- determine z-velocity */
Uz = C_W(cell,thread);
source = -(constant1*Uz + constant2 * fabs(Uz) * Uz);
dS[eqn] = -(constant1 + 2 * constant2 * fabs(Uz)); /* XXX CHECK */
return source; }
