PHOTON USE AUTOPLOT file phi 5 cl msg LAMINAR PIPE FLOW WITH HEAT TRANSFER msg constant wall heat flux msg Reynolds number = 10 Prandtl number = 1.0 msg Velocity (W1) profile msg Blue line --- PHOENICS solution msg crosses --- analytical solution da 1 w1;da 1 w1a col3 1;blb4 2 msg pressto continue pause cl msg constant wall heat flux msg Reynolds number = 10 Prandtl number = 1.0 msg Temperature (H1) profile msg Blue line --- PHOENICS solution msg crosses --- analytical solution da 1 h1;da 1 h1a col3 1;blb4 2 msg press to continue pause msg press to end pause end END_USE DISPLAY TEXT(1D Laminar Pipe Flow And Heat Trans TITLE DISPLAY The case considered is laminar flow and heat transfer of an incompressible, constant-property fluid in the hydrodynamically and thermally developed region of a of a circular tube. The heat transfer problem can be solved either with a constant wall heat flux ( QFCON=T ) or a constant wall temperature ( QFCON=F ). The Prandtl number is set to 1.0, but it should be noted that the Nusselt number is a function only of the thermal boundary condition, as follows: Nu = 3.656 for a constant surface temperature and Nu = 4.364 for a constant heat flux QIN. For the latter case, an analytical solution exists for the radial temperature profile: T = Tw-GAM*(0.1875+0.0625*(R/RO)**4-0.25*(R/RO)**2) where RO is the tube radius and GAM=4.*QIN/(K*R0) wherein K is the thermal conductivity of the fluid. ENDDIS For this case the velocity profile is parabolic: W = 2.*WIN*(1.-(R/RO)**2) where WIN is the bulk velocity and RO is the tube radius. The H1 equation is used as a vehicle for solving the temperature, so that the H1 equation is divided through by the fluid specific heat. For hydrodynamically and thermally developed flow, the following source/sink term appears in the temperature equation: -rho1*w1*dT/dz. For a constant wall heat flux, dT/dz=dTb/dz= dTw/dx where Tb is the bulk temperature and Tw is the wall temperature. For a constant wall temperature, dT/dz is not constant, rather it is equal to [(Tw-T)/(Tw-Tb)]*dTb/dz. For both cases the user must precribe dTb/dz, which may be evaluated from an overall energy balance for the slab. For the case of a uniform wall temperature, this means that the user must presently specify the energy flow out of the slab. GROUP 1. Run title and other preliminaries ** QFCON=T sets constant heat-flux boundary condition =F " " wall-temperature " " BOOLEAN(QFCON);QFCON=T REAL(REY,RIN,DIN,WIN,AIN,DPDZ,FLOWIN) REAL(QIN,DTDZ,CP,COND,AWAL,TC,TB,TW,BETA,GAM) REY=10.;RIN=0.1;DIN=2.*RIN;WIN=1.0 GROUP 2. Transience; time-step specification CARTES=F GROUP 3. X-direction grid specification AIN=0.5*RIN*RIN*XULAST GROUP 4. Y-direction grid specification NY=30;GRDPWR(Y,NY,RIN,1.0) GROUP 5. Z-direction grid specification ** AWAL is wall surface area per unit length AWAL=RIN*XULAST GROUP 7. Variables stored, solved & named ** W1A used to store analytical solution for W1 SOLVE(W1,H1);STORE(W1A) GROUP 8. Terms (in differential equations) & devices ** deactivate convection and built-in H1 sources TERMS(W1,N,N,P,P,P,P);TERMS(H1,N,N,P,P,P,P) GROUP 9. Properties of the medium (or media) PRNDTL(H1)=1.0;RHO1=1.0;ENUT=0.;ENUL=WIN*DIN/REY FLOWIN=RHO1*WIN*AIN ** compute pressure drop for printout from satellite Q1 DPDZ=8.*RHO1*ENUL*WIN/RIN/RIN DPDZ ** define input parameters for thermally-developed flow IF(QFCON) THEN ** specify wall flux, fluid specific heat & datum temperature at pipe axis + QIN=1.0;CP=1.0;TC=0.0 ** compute d(tbulk)/dz for input to single-slab thermal solver + DTDZ=QIN*AWAL/(CP*FLOWIN) ** H1A is a store for analytical temperature solution with constant heat-flux boundary condition + STORE(H1A) + COND=RHO1*CP*ENUL/PRNDTL(H1) + BETA=QIN*2.*RIN/COND + TB=TC+7.*BETA/48.;TW=TB+11.*BETA/48. + GAM=4.*QIN*RIN/COND ELSE + QIN=1.0;CP=1.0 + DTDZ=QIN*AWAL/(CP*FLOWIN) + TW=10. ENDIF GROUP 11. Initialization of variable or porosity fields FIINIT(W1)=WIN;FIINIT(H1)=0.4*DTDZ*ZWLAST ** compute analytical solutions REAL(WA,TA,GR,GR2);INTEGER(JJM1) DO JJ=1,NY +PATCH(IN:JJ:,INIVAL,1,NX,JJ,JJ,1,NZ,1,1) +GR=0.5*YFRAC(JJ) IF(JJ.NE.1) THEN +JJM1=JJ-1 +GR=YFRAC(JJM1)+0.5*(YFRAC(JJ)-YFRAC(JJM1)) ENDIF +GR=GR*YVLAST;GR2=(GR/RIN)**2 +WA=2.*WIN*(1.-GR2) +INIT(IN:JJ:,W1A,ZERO,WA) IF(QFCON) THEN +TA=TW-GAM*(0.1875+0.0625*(GR2**2)-0.25*GR2) +INIT(IN:JJ:,H1A,ZERO,TA) ENDIF ENDDO GROUP 13. Boundary conditions and special sources PATCH(WALL,NWALL,1,NX,NY,NY,1,NZ,1,1);COVAL(WALL,W1,1.0,0.0) ** activate pressure-drop calculation in single-slab solver FDFSOL=T;USOURC=T PATCH(FDFW1DP,VOLUME,1,NX,1,NY,1,NZ,1,1) COVAL(FDFW1DP,W1,FLOWIN,GRND1) ** temperature equation sources & boundary conditions IF(QFCON) THEN ** constant heat-flux boundary condition + PATCH(FDFCHF,PHASEM,1,NX,1,NY,1,NZ,1,1) + COVAL(FDFCHF,H1,DTDZ,GRND1);COVAL(WALL,H1,FIXFLU,QIN) ** set centre-line temperature as datum + PATCH(TFIX,CELL,1,NX,1,1,1,NZ,1,1) + COVAL(TFIX,H1,FIXP,TC) ELSE ** constant wall-temperature boundary condition + PATCH(FDFCWT,PHASEM,1,NX,1,NY,1,NZ,1,1) + COVAL(FDFCWT,H1,DTDZ,TW);COVAL(WALL,H1,1.0/PRNDTL(H1),TW) ENDIF GROUP 15. Termination of sweeps LSWEEP=8;LITHYD=2;LITC=5;ISWC1=4 GROUP 16. Termination of iterations RESREF(W1)=1.E-12*DPDZ*ZWLAST*AIN RESREF(H1)=1.E-12*QIN*ZWLAST*AWAL GROUP 17. Under-relaxation devices REAL(DTF);DTF=50.*(YVLAST/NY)**2/ENUL RELAX(W1,FALSDT,DTF) GROUP 19. Data communicated by satellite to GROUND COND=RHO1*CP*ENUL/PRNDTL(H1) COND DIN GROUP 21. Print-out of variables GROUP 22. Spot-value print-out IYMON=NY;TSTSWP=-1 GROUP 23. Field print-out and plot control NPLT=1;NYPRIN=1;NZPRIN=1 GROUP 24. Dumps for restarts