TEXT(RSTM_1D DEVELOPED PIPE FLOW :T602
TITLE
DISPLAY
The problem considered is 1d fully-developed turbulent flow
and heat and mass transfer in a pipe at a Reynolds number of
5.E4. The turbulence is simulated by use of the Reynolds stress
transport model (RSTM), and optionally by use of the k-eps model.
The energy equation may be solved by means of the H1 or TEM1
equation, and two scalar variables are solved which should
produce identical solutions to the energy equation. For heat and
scalar transport with the RSTM, one of the following models may
be employed: a simple gradient-diffusion model (IRSMSM=0); a
generalised gradient-diffusion model (IRSMSM=1); or a full
transport model (IRSMSM=2). The laminar Prandtl number is 3.0.
ENDDIS
For fully-smooth conditions pipe-flow data indicates that
the friction factor f = 0.021 and the Nusselt number Nu = 222.0.
The Petukhov correlation ( see below ) is used to estimate
the Nusselt number. The PHOENICS predictions yield the
following results:
k-eps model f = 0.021 Nu = 222
RSTM gradient-diffusion model f = 0.021 Nu = 218
RSTM generalised-diffusion model f = 0.021 Nu = 225
RSTM full transport model f = 0.021 Nu = 224.
These values agree closely with the data.
The PIL variable WALPRN has been set to T, thereby activating
printout of the local friction factor (sloc) and Stanton number
(Stloc) in the RESULT file. In order to convert these values to
f and Nu above, the following relations should be used:
f = 8.*sloc*(W1(NY)/WIN)**2
Nu = REY*PRNDTL(H1)*Stloc*W1(NY)*(TW-TEM1(NY))/[WIN*(TW-TB)]
where TB is the bulk temperature printed in the RESULT file. The
expected values are: sloc=4.923E-3 and Stloc=2.512E-3.
BOOLEAN(KEMOD,HEAT)
** CH1=H1 activates solution of the energy eqn via H1
=TEM1 '' '' '' '' '' '' '' TEM1
CHAR(CH1);CH1=TEM1
** KEMOD=T selects k-eps solution
** HEAT=T activates solution of the energy and scalar equations
KEMOD=F;HEAT=T;IRSMSM=2
REAL(REY,DIAM,WIN,US,DPDZ,FRIC,TKEIN,EPSIN,MIXL,DELT,DTF)
** REYS = Friction-velocity Reynolds number
** REYZ = Axial-velocity Reynolds number
DIAM=1.0;WIN=1.0; REY=5.E4; RHO1=1.0;ENUL=WIN*DIAM/REY
** estimate initial values of k and eps
FRIC=1./(1.82*LOG10(REY)-1.64)**2
US=WIN*(FRIC/8.)**0.5;DPDZ=0.5*RHO1*WIN*WIN*FRIC/DIAM
REY
US
FRIC
DPDZ
TKEIN=0.25*WIN*WIN*FRIC
MIXL=0.09*0.5*DIAM;EPSIN=0.1643*TKEIN**1.5/MIXL
** compute axial pressure drop given the Reynolds number
XULAST=0.1;CARTES=F;DELT=2.*40.*ENUL/US
NREGY=2; REGEXT(Y,0.5*DIAM);IREGY=1;GRDPWR(Y,29,0.5*DIAM-DELT,0.8)
IREGY=2;GRDPWR(Y,1,DELT,1.0);SOLVE(W1)
IF(HEAT) THEN
+ SOLVE(:CH1:);SOLVE(SC1,SC2)
ENDIF
IF(.NOT.KEMOD) THEN
+ PATCH(WAL1,NWALL,1,1,NY,NY,1,NZ,1,1)
ENDIF
** prescribed pressure-force for w1-equation
PATCH(PFOR,VOLUME,1,1,1,NY,1,NZ,1,1);COVAL(PFOR,W1,FIXFLU,DPDZ)
IF(KEMOD) THEN
+ TURMOD(KEMODL);STORE(ENUT,LEN1)
+ PATCH(WAL1,NWALL,1,1,NY,NY,1,NZ,1,1)
+ COVAL(WAL1,W1,LOGLAW,0.0);COVAL(WAL1,EP,LOGLAW,LOGLAW)
+ COVAL(WAL1,KE,LOGLAW,LOGLAW)
ELSE
+ STORE(V1,KE,DWDY,PVW,PW2,PV2,PU2,DVW)
+ STORE(PK,EPDK,FWAL,VWDK,U2DK,V2DK,W2DK)
+ PATCH(WAL1,NWALL,1,1,NY,NY,1,NZ,1,1);COVAL(WAL1,W1,1.0,0.0)
+ DTF=0.05;TURMOD(REYSTRS,DTF,WAL1)
+ PATCH(SMPLS,SOUTH,1,1,1,1,1,NZ,1,1);COVAL(SMPLS,VWRS,GRND1,0.0)
ENDIF
** deactivate convection for single-slab solution
TERMS(W1,P,N,P,P,P,P);TERMS(EP,P,N,P,P,P,P);DIFCUT=0
FIINIT(W1)=WIN;FIINIT(V1)=0.0;FIINIT(KE)=TKEIN;FIINIT(EP)=EPSIN
REAL(PI,AIN,FLOW);PI=3.14159265
AIN=PI*DIAM*DIAM/4.*XULAST/(2.*PI);FLOW=RHO1*WIN*AIN
FLOW
IF(KEMOD) THEN
+ TERMS(KE,P,N,P,P,P,P);IRSMSM=0
+ RELAX(U1,FALSDT,1.E2); RELAX(W1,FALSDT,1.E2)
+ RELAX(EP,FALSDT,10.0); RELAX(KE,FALSDT,10.0)
+ RESREF(KE)=1.E-12*FLOW*TKEIN
ELSE
+ TERMS(U2RS,P,N,P,P,P,P);TERMS(V2RS,P,N,P,P,P,P)
+ TERMS(W2RS,P,N,P,P,P,P);TERMS(VWRS,P,N,P,P,P,P)
+ FIINIT(W2RS)=2.*FIINIT(KE)/3.;FIINIT(V2RS)=FIINIT(W2RS)
+ FIINIT(U2RS)=FIINIT(W2RS);FIINIT(VWRS)=0.3*FIINIT(KE)
+ RELAX(W1,FALSDT,0.1); RELAX(EP,FALSDT,0.1)
+ RESREF(U2RS)=1.E-12*FLOW*FIINIT(U2RS)
+ RESREF(V2RS)=1.E-12*FLOW*FIINIT(V2RS)
+ RESREF(W2RS)=1.E-12*FLOW*FIINIT(W2RS)
+ RESREF(VWRS)=1.E-12*FLOW*FIINIT(VWRS)
ENDIF
RESREF(W1)=1.E-12*FLOW*WIN; RESREF(EP)=1.E-12*FLOW*EPSIN
IF(KEMOD) THEN
+ LSWEEP=20;NPLT=5
ELSE
+ LSWEEP=300;LITHYD=8;NPLT=20
ENDIF
NYPRIN=1;IYMON=NY-1;TSTSWP=-1;YPLS=T;WALPRN=T
** prescribe energy flow from slab and fluid specific heat
estimated from Dittus-Boelter Nu=0.023*Re**0.8*Pr**0.4
with (Tw-Tb)=5.0
IF(HEAT) THEN
+ REAL(NUSS,COND,CP,QIN,DTDZ,AWAL,TW,XR,DTSC)
+ SOLVE(:CH1:);SOLVE(SC1,SC2);DTSC=1.0
+ PRNDTL(H1)=3.0;PRNDTL(SC1)=PRNDTL(H1)
+ PRNDTL(SC2)=PRNDTL(H1)
+ NUSS=0.023*REY**0.8*PRNDTL(H1)**0.4
+ CP=1.0;COND=RHO1*CP*ENUL/PRNDTL(H1)
** AWAL is wall surface area per unit length
+ AWAL=0.5*DIAM*XULAST;QIN=NUSS*5.0*COND/DIAM
+ NUSS
+ QIN
** compute d(tbulk)/dz for input to single-slab
thermal solver
+ DTDZ=QIN*AWAL/(CP*FLOW)
IF(:CH1:.EQ.TEM1) THEN
DTDZ=CP*DTDZ
ENDIF
+ AWAL
** prescribe wall temperature
+ TW=10.;TERMS(:CH1:,N,N,P,P,P,P)
+ TERMS(SC1,N,N,P,P,P,P);TERMS(SC2,N,N,P,P,P,P)
+ FIINIT(:CH1:)=0.5*TW;FIINIT(SC1)=FIINIT(:CH1:)
+ FIINIT(SC2)=FIINIT(:CH1:);COVAL(WAL1,:CH1:,LOGLAW,TW)
+ COVAL(WAL1,SC1,LOGLAW,TW);COVAL(WAL1,SC2,LOGLAW,TW)
** temperature source/sink term for fully-developed flow
+ PATCH(FDFCWT,PHASEM,1,NX,1,NY,1,NZ,1,1)
+ COVAL(FDFCWT,:CH1:,DTDZ,TW)
+ COVAL(FDFCWT,SC1,DTDZ,TW);COVAL(FDFCWT,SC2,DTDZ,TW)
+ FDFSOL=T; RESREF(:CH1:)=1.E-12*QIN*AWAL*ZWLAST
+ RESREF(SC1)=RESREF(:CH1:); RESREF(SC2)=RESREF(SC1)
** compute expected Nusselt number from Petukhov
correlation and printout from satellite
+ XR=1.07+12.7*(PRNDTL(H1)**.666-1.)*(FRIC/8.)**0.5
+ NUSS=REY*PRNDTL(H1)*FRIC/(8.*XR)
+ NUSS
+ COND
+ DIAM
IF(IRSMSM.EQ.1) THEN
+ DTSC=0.5;STORE(DSDY)
+ RELAX(VTRS,LINRLX,0.2); RELAX(VSC1,LINRLX,0.2)
+ RELAX(VSC2,LINRLX,0.2)
+ OUTPUT(VTRS,Y,Y,Y,Y,Y,Y)
ENDIF
IF(IRSMSM.EQ.2) THEN
+ TERMS(VTRS,N,N,P,P,P,P);TERMS(VSC1,N,N,P,P,P,P)
+ TERMS(VSC2,N,N,P,P,P,P);STORE(DSDY)
+ COVAL(SMPLS,VTRS,GRND1,0.0);DTSC=0.5
+ COVAL(SMPLS,VSC1,GRND1,0.0);COVAL(SMPLS,VSC2,GRND1,0.0)
+ RESREF(VTRS)=RESREF(VWRS); RESREF(VSC1)=RESREF(VTRS)
+ RESREF(VSC2)=RESREF(VTRS)
ENDIF
RELAX(:CH1:,FALSDT,DTSC); RELAX(SC1,FALSDT,DTSC)
RELAX(SC2,FALSDT,DTSC)
IF(:CH1:.EQ.TEM1) THEN
PRNDTL(TEM1)=CONDFILE;
STORE(PRPS);FIINIT(PRPS)=36
** mat no. rho enul cp kond expan
** 1 air
CSG10='q1'
MATFLG=T;NMAT=1
36 1. 2.E-5 1.0 6.667E-6 0
ENDIF
ENDIF