UDF Scalar Code: HT 1
For those interested in some code (see post 2) for scalar transport eqns under Fluent v5 here's some examples. . . I developed this for solving heat transfer between fluid and a packed bed solid (after a few weeks I moved to Fluent v4.5 which is better for this application).
Anyway this should allow you to solve for solid phase temperature as a scalar transport eqn as a special user porous media model assuming a packed bed of spheres. I've cut down the code so it may or may not compile and I don't provide any guarantees . . it may be more confusing than helpful!!! But if you play around I think you can get it to work  once compiled  use makefiles  don't think it will work under interpreted mode (but try??)  then link the source terms to the appropriate Fluent eqn. That's about it. Unfortunately, if it doesn't work/doesn't work for you I can't provide detailed responses to questions  sorry. A bit on how it works: it uses two main UDS's  ZETA = local porosity and SOLID_TEMP for local solid temperature. ZETA just stores the value, SOLID_TEMP can be solved for. Another six UDS are required for other debug type results that I found useful. So in total declare space for 8 UDS, but only solve for SOLID_TEMP  turn off solution of all others. (with mods you can update ZETA due to reaction for example) I would use the UDMI function in 5.4.8 for storing values at cells  but that wasn't available until recently. There are source terms for momentum  SRCE_Xmon etc.  these calculate the drag due to the solid spheres on the fluid phase using the Ergun eqn. Link these into the fluid domain sections in the Boundary Conditions panel. There are source terms for energy  fluid and solid  this just calculates the interphase heat transfer with a HTC correlation. Link these as source terms for fluid energy as solved by Fluent and the SOLID_TEMP UDS. If you solve steady probs there's really no point for this. Use the diffusivity routine to tell Fluent how to calculate the diffusivity for your UDS's. This only calculates for the SOLID_TEMP eqn where it uses alpha = k/rho*Cp = thermal diffusivity. You can use the init routine to help initialise the domain but will need mods depending on what you try out. Then set some b/cs  note to set also your SOLID_TEMP UDS b/c on each boundary  otherwise your solid temp will be 0 Kelvin which is not realistic!!! Greg Perkins ps: if you want to use a variable porosity use Fluent v4.5  its toooo hard to implement in Fluent v5.x 
UDF Scalar Code: HT 2
For those interested in some code (see post 2) for scalar transport eqns under Fluent v5 here's some examples. . . I developed this for solving heat transfer between fluid and a packed bed solid (after a few weeks I moved to Fluent v4.5 which is better for this application).
Anyway this should allow you to solve for solid phase temperature as a scalar transport eqn as a special user porous media model assuming a packed bed of spheres. I've cut down the code so it may or may not compile and I don't provide any guarantees . . it may be more confusing than helpful!!! But if you play around I think you can get it to work  once compiled  use makefiles  don't think it will work under interpreted mode (but try??)  then link the source terms to the appropriate Fluent eqn. That's about it. Unfortunately, if it doesn't work/doesn't work for you I can't provide detailed responses to questions  sorry. A bit on how it works: it uses two main UDS's  ZETA = local porosity and SOLID_TEMP for local solid temperature. ZETA just stores the value, SOLID_TEMP can be solved for. Another six UDS are required for other debug type results that I found useful. So in total declare space for 8 UDS, but only solve for SOLID_TEMP  turn off solution of all others. (with mods you can update ZETA due to reaction for example) I would use the UDMI function in 5.4.8 for storing values at cells  but that wasn't available until recently. There are source terms for momentum  SRCE_Xmon etc.  these calculate the drag due to the solid spheres on the fluid phase using the Ergun eqn. Link these into the fluid domain sections in the Boundary Conditions panel. There are source terms for energy  fluid and solid  this just calculates the interphase heat transfer with a HTC correlation. Link these as source terms for fluid energy as solved by Fluent and the SOLID_TEMP UDS. If you solve steady probs there's really no point for this. Use the diffusivity routine to tell Fluent how to calculate the diffusivity for your UDS's. This only calculates for the SOLID_TEMP eqn where it uses alpha = k/rho*Cp = thermal diffusivity. You can use the init routine to help initialise the domain but will need mods depending on what you try out. Then set some b/cs  note to set also your SOLID_TEMP UDS b/c on each boundary  otherwise your solid temp will be 0 Kelvin which is not realistic!!! Greg Perkins ps: if you want to use a variable porosity use Fluent v4.5  its toooo hard to implement in Fluent v5.x /************************************************** **********************/ /* Reaction Module Source Code */ /* */ /* By Greg Perkins */ /* CANCES, Australian Technology Park */ /* Ph: 02 9318 0004, Fax: 02 9319 2328 */ /* Email: perkinsg@cances.atp.com.au */ /* */ /* Started: 03012000 */ /* */ /* Revisions: */ /* */ /* 14022000 (i) Solid phase temperature and UDS (ii) General species structure, material properties and reactions */ /************************************************** **********************/ /*  These routines implement a hetrogeneous reaction model for the combustion and gasification reactions of carbon together with a porous media model */ #include "udf.h" #include "sg.h" /*  REACTION MODEL DEFINES  */ #define PI 3.141592654 #define KJCAL 4.184 #define UNIVERSAL_GAS_CONST 8.314 /*  [J/molK] */ #define SOLID_DENSITY(T) 1200 /*  Coal material density [kg/m3] */ #define SOLID_THERMAL_COND(T) 1.47/1.0*pow(T/273,0.5) /*  Hobbs et al. 1993, [W/mK] */ #define SOLID_CP(T) 3000.0*UNIVERSAL_GAS_CONST/12.0*(exp(1200/T)/pow((exp(1200/T)1)*T/1200,2.0)) /*  Hobbs et al. 1993, [J/kgK] */ #define alpha_s 0.5 /*  fraction of energy released in solid phase */ #define alpha_f 0.5 /*  fraction of energy released in fluid phase */ #define ZETA_CO 0.8 /*  Initial coal vol. frac. */ #define ZETA_VO 0.2 /*  Initial void vol. frac. */ #define R_CO 0.01 /*  Initial coal particle radius [m] */ #define ZETA_MAX 0.95 /*  Maximum porosity without checks */ #define DTHI_MAX 1.0e30 #define PERM_EXP_VOID 25.0 #define TEMP_KELVIN 0 #define TEMP_DEGC 1 #define TREF 273.0 /*  POROUS MEDIA MODEL DEFINES  */ /*  use a user defined scalar to track material properties of the domain */ enum { ZETA, SOLID_TEMP, SRCE_SOLID, SRCE_FLUID, DTEMP, DZETA, THI, RE, N_REQUIRED_UDS }; /*  Material Properties This routine returns the local properties of the material (C2,C2) for a given cell in the domain. 1. CHECK IMPLEMENTATION  */ void Material_Properties(cell_t c, Thread *t, real *C1, real *C2) { real zeta, zeta3, temp1; /*  use this routine to return the permeability, alpha and the coefficient 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) */ zeta = C_UDSI(c,t,ZETA); zeta3 = pow(zeta,3.0); temp1 = pow((1.0((zetaZETA_VO)/ZETA_CO)),(1.0/3.0)); temp1 = 1.0; /*  Dp fixed */ /*  Calculate local permeability = 1/C1 */ *C1 = (150.0*pow((1.0zeta),2.0))/(zeta3*4.0*pow(R_CO,2.0)*pow(temp1,2.0)); /*  Calculate local C2 coefficient */ *C2 = 1.75*(1.0zeta)/(R_CO*temp1*zeta3); } /*  X_Momentum_Source This routine returns the source term for the Xmomentum 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 C1, C2, constant1, constant2, Ux, source; /*  determine local properties */ Material_Properties(cell, thread, &C1, &C2); /*  determine constants 1,2 */ constant1 = C_MU_L(cell,thread) * C1; constant2 = 0.5 * C_R(cell,thread) * C2; /*  determine xvelocity */ Ux = C_U(cell,thread); source = (constant1*Ux + constant2 * fabs(Ux) * Ux); dS[eqn] = (constant1 + 2 * constant2 * fabs(Ux)); return source; } /*  Y_Momentum_Source This routine returns the source term for the Ymomentum 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 C1, C2, constant1, constant2, Uy, source; /*  determine local properties */ Material_Properties(cell, thread, &C1, &C2); /*  determine constants 1,2 */ constant1 = C_MU_L(cell,thread) * C1; constant2 = 0.5 * C_R(cell,thread) * C2; /*  determine yvelocity */ Uy = C_V(cell,thread); source = (constant1*Uy + constant2 * fabs(Uy) * Uy); dS[eqn] = (constant1 + 2 * constant2 * fabs(Uy)); return source; } /*  Z_Momentum_Source This routine returns the source term for the Zmomentum 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 C1, C2, constant1, constant2, Uz, source; /*  determine local properties */ Material_Properties(cell, thread, &C1, &C2); /*  determine constants 1,2 */ constant1 = C_MU_L(cell,thread) * C1; constant2 = 0.5 * C_R(cell,thread) * C2; /*  determine zvelocity */ Uz = C_W(cell,thread); source = (constant1*Uz + constant2 * fabs(Uz) * Uz); dS[eqn] = (constant1 + 2 * constant2 * fabs(Uz)); return source; } /*  Initialization UDF This routine initializes the solution in the computational domain.  */ DEFINE_INIT(INIT_Solution, domain) { cell_t c; face_t f; Thread *thread; real x[ND_ND]; /*  Set initial porosity to be ZETA_VO for y>0.5m */ thread_loop_c (thread, domain) { begin_c_loop (c, thread) { C_CENTROID(x,c,thread); if (x[1] > 0.0) C_UDSI(c, thread, ZETA) = ZETA_VO; else C_UDSI(c, thread, ZETA) = 1.0; } end_c_loop (c, thread) } } /*  UDS Diffusivity [m2/s] This routine returns the diffusivity for the user defined scalar transport equations. The only nonzero setting is for the solid phase temperature eqn (solid energy eqn). (RB002, p101)  */ DEFINE_DIFFUSIVITY(UDS_Diffusivity, cell, thread, i) { real source,Ts; if (i==SOLID_TEMP) { Ts = C_UDSI(cell,thread,SOLID_TEMP); source = SOLID_THERMAL_COND(Ts)/( SOLID_DENSITY(Ts)*SOLID_CP(Ts) ); } else source = 0.0; return source; } /*  GASSOLID Heat Transfer Coefficient [W/m2K] This routine calculates the gassolid heat transfer coefficient source: Wakao & Kaguei, (1982), (RB002, p103)  */ real Gas_Solid_HTC(cell_t cell, Thread *thread) { real Pr,Re,Us; /*  calculate superficial velocity */ Us = sqrt( pow(C_U(cell,thread),2.0) + pow(C_V(cell,thread),2.0) + pow(C_W(cell,thread),2.0) ); /*  calculate Prandtl number */ Pr = C_MU_EFF(cell,thread)*C_CP(cell,thread)/C_K_EFF(cell,thread); /*  calculate Reynolds no. */ Re = 2*R_CO*Us*C_R(cell,thread)/C_MU_EFF(cell,thread); return C_K_EFF(cell,thread)/(2*R_CO)*(2.0 + 1.1*pow(Pr,1.0/3.0)*pow(Re,0.6)); } /*  ENERGY SOURCE: Solid Phase [K/s] This routine calculates the solid phase energy source term. This is composed of two parts: the heat transfer between the solid and gas and the energy released in the solid phase, (RB002, p102)  */ DEFINE_SOURCE(SRCE_Solid_Energy,cell,thread,dS,eqn ) { real zeta, source, thi, dthi, Ts, HTC; /*  Current porosity */ zeta = C_UDSI(cell, thread, ZETA); /*  Specific Surface Area Term */ thi = (3*ZETA_CO/R_CO)*pow((1.0((zetaZETA_VO)/ZETA_CO)),2.0/3.0); if (zeta > ZETA_MAX) dthi = DTHI_MAX; else dthi = (2.0/R_CO)*pow((1.0((zetaZETA_VO)/ZETA_CO)),1.0/3.0); /*  Check that there is some solid in this cell*/ if (zeta>0.99) zeta=0.9999; /*  Solid Temperature & Heat Transfer Coeff */ Ts = C_UDSI(cell,thread,SOLID_TEMP); HTC = Gas_Solid_HTC(cell,thread); source = HTC*thi*(TsC_T(cell,thread)); source = source/( SOLID_DENSITY(Ts)*SOLID_CP(Ts)*(1zeta) ); dS[eqn] = HTC*thi/( SOLID_DENSITY(Ts)*SOLID_CP(Ts)*(1zeta) ); C_UDSI(cell,thread,SRCE_SOLID) = source; C_UDSI(cell,thread,RE) = HTC; return source; } /*  ENERGY SOURCE: Fluid Phase [J/m3s] This routine calculates the fluid phase energy source term. This is composed of two parts: the heat transfer between the solid and gas and the energy released in the fluid phase, (RB002, p102)  */ DEFINE_SOURCE(SRCE_Fluid_Energy,cell,thread,dS,eqn ) { real zeta, source, thi, dthi, Ts, HTC; /*  Current porosity */ zeta = C_UDSI(cell, thread, ZETA); /*  Specific Surface Area Term */ thi = (3*ZETA_CO/R_CO)*pow((1.0((zetaZETA_VO)/ZETA_CO)),2.0/3.0); if (zeta > ZETA_MAX) dthi = DTHI_MAX; else dthi = (2.0/R_CO)*pow((1.0((zetaZETA_VO)/ZETA_CO)),1.0/3.0); /*  Check that there is some fluid in this cell */ if (zeta<0.01) zeta=0.0001; /*  Solid Temperature & Heat Transfer Coeff */ Ts = C_UDSI(cell,thread,SOLID_TEMP); HTC = Gas_Solid_HTC(cell,thread); source = HTC*thi*(TsC_T(cell,thread))/zeta; dS[eqn] = HTC*thi/zeta; return source; } /*  REFERENCES Felder, R. M., and Rousseau, R. W., (1986), "Elementary principles of Chemical Processes",2nd Edition, John Wiley & Sons, New York. Hobbs, M., L., Radulovic, P. T., and Smoot, L. D., (1993), "Combustion and Gasification of Coals in Fixed Beds", Prog. Eenergy Combust. Sci., Vol. 19, pp.505586. Wakao, N., and Kaguei, S., (1982), "Heat and Mass Transfer in Packed Beds", Gordon & Breach Science Publishers, New York.  */ 
UDF Scalar Code: HT 2
Well, the code is actually at:
for Dos/Windows: www.perkinssoftware.com.au/Fluent/udf_scalar_dos.c for Unix: www.perkinssoftware.com.au/Fluent/udf_scalar_unix.c Have fun! Greg 
Re: UDF Scalar Code: HT 1
Greg,
a very interesting and helpful post. I just wanted to ask where did you find it particularly difficult to implement variable porosity in F5? I ask this because I've done this before (my PhD thesis, due to be handed in in a month, is based around this concept) and it wasn't as difficult. Then again, I've never used F4.5 since I've moved to Fluent straight from FIDAP, so maybe there is a big difference in implementation. Anyway, this is just my curiosity since maybe I could have developed one third of my thesis in a third of the time! Go figure :) Cheerio! AL  Alfonso Ferrandez  ferrandez@mecheng.leeds.ac.uk Computational BioFluid Mechanics Research Group The University of Leeds, Leeds, UK URL: http://www.mecheng.leeds.ac.uk/ferrandez/ 
Re: UDF Scalar Code: HT 1
Alfonso,
that's very intresting! How did you handle the change in porosity on the convective term in the various momentum and energy equations?????? When I looked into it, I worked out these terms need to be scaled by the porosity, which is a local cell property. Ie., similar to the equations for multiphase flow in Fluent v4.4/4.5 chapter 19. One option was to correct for this with a source term, but this seemed very difficult to implement. What did you do???? How???? Do you have the full set of equations you solve for in Fluent v5?? what type of UDFs did you write  different from mine??? Very interested to find out more, Regards Greg Perkins 
Re: UDF Scalar Code: HT 1
Hi Greg!
That was a quick reply if I ever saw one :) I did in fact write my own source term equations, since the porous media set in Fluent wasn't complex enough. This required a lot of literature searching on porous transport until eventually I came up with the proper source term functions. Unfortunately I can't send you the details as of yet, since my thesis has been set under a 2/3 year NDA (industry sponsors, go figure!) I had the same aversion to source term development as you, but after weeks of frustration, I made the bloody thing work! I suppose if I'd known about F4.5 I wouldn't have bothered :) Cheers! AL 
Re: UDF Scalar Code: HT 2
Greg,
I noticed your UDF deals with commbustion/gasification of carbon (char?) in a packed bed. Are you aware of the work of Prof. William Hallett at Univ. of Ottawa in Canada? http://by.genie.uottawa.ca/profs/hallett/packedbd.html Regards, Henrik 
Re: UDF Scalar Code: HT 2
Thanks Henrik,
yes I'm solving for underground coal gasification  in which part of the process can be represented like a packed bed surface gasifier.... Regards Greg 
Re: UDF Scalar Code: HT 1
Thanks Greg, very helpful code, do you mind if I ask you some questions regarding the scalar implementation? Regards Arturo

All times are GMT 4. The time now is 12:40. 