```
PHOTON USE
AUTOPLOT
file
PHI 5

cl;d 1 h1;d 1 ha;col3 1;blb4 2;redr
msg    temperature profile; press  to 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
TITLE
DISPLAY
The problem considered is that of radiative heat transfer
in a 1d plane-parallel slab in radiative equilibrium. The
slab is of thickness L and the gray medium may absorb,
emit and scatter radiation. At x=L the wall is a diffuse
emitter and reflecter kept at a fixed temperature. At the x=0
ENDDIS

Since the energy transfer is by pure radiation, the energy
equation is given by:

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 = E   ( radiative equilibrium )

E  = Ew + Qrad*[ 1./emw - 0.5 + 0.5*(a+s)*L*(1 - x/L) ]

where Ew is the emissive power at the wall and emw is the
emissivity of the wall. For the case considered below, Qrad
is taken to be Qrad = 2.*Ew.

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       wall temperature at x=L         {  K   }

GSIGMA=5.6697E-8;SCAT=0.5;ABSORB=0.5;EMIW=1.0;TWAL=1000.0
GROUP 3,4,5. X,Y,Z-direction grid specification
GRDPWR(X,50,1.0,1.0)
GROUP 7. Variables stored, solved & named
CP1=1.0
MESG( Enter required energy variable ? (TEM1 or H1)
IF(:CH1:.EQ.) THEN
ENDIF
IF(:CH1:.EQ.TEM1) THEN
+ MESG( TEM1 solution selected
ELSE
+ MESG( H1 solution selected
+ TMP1=LINH;TMP1B=1.0/CP1
ENDIF
GROUP 8. Terms (in differential equations) & devices
** Deactive conduction & any built-in sources
TERMS(:CH1:,N,N,N,N,P,P)
GROUP 11. Initialization of variable or porosity fields
** analytical solution
BET=0.5*(ABSORB+SCAT)*XULAST
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
+RAN=EG;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(WALLRA,WEST,1,1,1,NY,1,NZ,1,1)
PATCH(WALLRB,EAST,NX,NX,1,NY,1,NZ,1,1)
GROUP 15. Termination of sweeps
LSWEEP=50
GROUP 16. Termination of iterations