PHOTON USE AUTOPLOT file PHI 5 cl;d 1 h1;d 1 ha;col3 1;blb4 2;redr msg temperature profile; pressto continue pause cl;d 1 radx;d 1 ra;col3 1;blb4 2;redr msg radiation-flux profile; press e to end pause;end ENDUSE GROUP 1. Run title and other preliminaries TEXT(Radiative Transfer + Heat Source TITLE DISPLAY The problem considered is that of pure radiative heat transfer in a 1d plane-parallel slab containing a uniformly distributed heat source. The medium may absorb, emit and scatter radiation and the boundaries of the slab are diffuse emitters and reflecters kept at the same uniform temperature. The medium is gray and because the wall temperatures are equal, symmetry can be exploited. ENDDIS Since the energy transfer is by pure radiation the energy equation is given by: -d/dx (Qrad) + Qvol = 0 where Qrad is the radiative heat flux and Qvol is the uniform volumetric heat generation rate in the medium. The equation for the composite radiative heat flux is given by: d/dx ( 1/(a+s) d/dx (Rx) ) + a (E - Rx) = 0 where a is the absorption coefficient, s is the scattering coefficient, Rx is the composite radiation flux defined as the average of the +ve and -ve radiation fluxes, and E is the black-body emissive power. It may be noted that the radiative heat flux is given by: d/dx (Qrad) = 2a (E - Rx) The black-body emissive power E=sig*T**4 where sig is the Stefan-Boltzmann constant and T is the temperature of the medium. The problem is the determination of the temperature and composite radiative-flux distributions, as given by the following analytical solutions: Rx = Ew + 0.5*Qvol*L* [ 2/emw-1+0.5*(a+s)*L*{1-(x/L)**2} ] E = Ew + Qvol*L*[ 1./(2*a*L) + 1./emw - 0.5 + 0.25*(a+s)*L*{1 - (x/L)**2} ] where Ew is the emissive power at the wall, L is the slab width from symmetry plane to wall, and emw is the emissivity of the wall. For the case considered below, Qvol is taken to be Qvol = Ew/L/(1./emw - 0.5). The locally-defined parameters are as follows: GSIGMA Stefan-Boltzmann constant {W/m**2/K**4} SCAT Scattering coefficient {1/m} ABSORB Absorption coefficient {1/m} EMIW emissivity of the wall TWAL hot wall temperature { K } QVOL volumetric internal heat source {W/m**3} CHAR(CH1);REAL(GSIGMA,SCAT,ABSORB,EMIW,TWAL,QVOL) GSIGMA=5.6697E-8;SCAT=0.5;ABSORB=0.5;EMIW=1.0;TWAL=1000.0 GROUPS 3,4,5. X,Y,Z-direction grid specification REAL(LENGTH);LENGTH=1.0;GRDPWR(X,50,LENGTH,1.0) GROUP 7. Variables stored, solved & named GROUP 7. Variables stored, solved & named CP1=1.0 MESG( Enter required energy variable ? (TEM1 or H1) READVDU(CH1,CHAR,H1) IF(:CH1:.EQ.TEM1) THEN + MESG( TEM1 solution selected ELSE + MESG( H1 solution selected + TMP1=LINH;TMP1B=1.0/CP1 ENDIF RADIAT(FLUX,ABSORB,SCAT,:CH1:);STORE(EMPO) CP1=1.0 MESG( Enter required energy variable ? (TEM1 or H1) IF(:CH1:.EQ.' ') THEN READVDU(CH1,CHAR,H1) ENDIF IF(:CH1:.EQ.TEM1) THEN + MESG( TEM1 solution selected ELSE + MESG( H1 solution selected + TMP1=LINH;TMP1B=1.0/CP1 ENDIF RADIAT(FLUX,ABSORB,SCAT,:CH1:);STORE(EMPO) GROUP 8. Terms (in differential equations) & devices ** Deactive built-in sources, convection and conduction TERMS(:CH1:,N,N,N,N,P,P) GROUP 11. Initialization of variable or porosity fields FIINIT(RADX)=GSIGMA*TWAL**4;FIINIT(:CH1:)=TWAL ** Define Qvol = Ew/L/ (1./emw - 0.5) QVOL=FIINIT(RADX)/XULAST/(1./EMIW-0.5) ** analytical solution REAL(EW,EG,QRADA,ALF,BET,GAM,ACON,TA,RAN,GX,GXRAT) STORE(HA,RA);INTEGER(JJM1) EW=GSIGMA*TWAL**4;QRADA=2.*EW;ALF=2./EMIW-1.0;GAM=0.5*QVOL*XULAST BET=0.5*(ABSORB+SCAT)*XULAST;ACON=1.0/(2.*ABSORB*XULAST)+0.5*ALF DO JJ=1,NX +PATCH(IN:JJ:,INIVAL,JJ,JJ,1,NY,1,NZ,1,1) +GX=0.5*XFRAC(JJ) IF(JJ.NE.1) THEN +JJM1=JJ-1;GX=XFRAC(JJM1)+0.5*(XFRAC(JJ)-XFRAC(JJM1)) ENDIF +GX=GX*XULAST;GXRAT=(GX/XULAST)**2 +EG=EW+QVOL*XULAST*(ACON+0.5*BET*(1.-GXRAT)) +RAN=EW+GAM*(ALF+BET*(1.0-GXRAT));TA=(EG/GSIGMA)**0.25 +INIT(IN:JJ:,RA,ZERO,RAN);INIT(IN:JJ:,HA,ZERO,TA) ENDDO GROUP 13. Boundary conditions and special sources ** Net radiation flux from wall PATCH(WALLR,EAST,NX,NX,1,NY,1,NZ,1,1) COVAL(WALLR,RADX,EMIW/(2.0-EMIW),GSIGMA*TWAL**4) ** uniformly-distributed volumetric heat source PATCH(QHEAT,VOLUME,1,NX,1,NY,1,NZ,1,1) COVAL(QHEAT,:CH1:,FIXFLU,QVOL) GROUP 15. Termination of sweeps LSWEEP=50 GROUP 16. Termination of iterations RESREF(:CH1:)=1.E-12*QVOL*LENGTH;RESREF(RADX)=0.5*RESREF(:CH1:) GROUP 17. Under-relaxation devices RELAX(:CH1:,FALSDT,100./QVOL) GROUP 22. Spot-value print-out IXMON=NX/2;NPLT=5;NXPRIN=NX/10 GROUP 23. Field print-out and plot control OUTPUT(:CH1:,Y,N,N,Y,Y,Y) PATCH(XWISE,PROFIL,1,NX,1,NY,1,NZ,1,1) PLOT(XWISE,:CH1:,0.0,0.0);PLOT(XWISE,RADX:,0.0,0.0) GROUP 24. Dumps for restarts TSTSWP=-1