#immersol)                                                                      
#declare                                                                        
   IMMERSOL Emitting and absorbing gray media (1D tests;                        
            solving for H1)                                                     
#cls                                                                            
TEXT(IMMERSOL 1D Radiative equilibrium  :203                                    
TITLE                                                                           
  **************************************************************                
   PHOTON USE                                                                   
   AUTOPLOT                                                                     
   FILE                                                                         
   PHIDA 3                                                                        
                                                                                
   CL;DA 1 TMP1;DA 1 T3;COL3 1;BLB4 2;REDR                                      
   MSG    Temperature and T3 profiles.                                          
   ENDUSE                                                                       
  **************************************************************                
  DISPLAY                                                                       
   This problem is similar to that of case 202; but here solution               
   for gas enthalpy replaces solution for temperature directly.                 
                                                                                
  ENDDIS                                                                        
#direct                                                                         
#temps                                                                          
MESG(  ENTER optical thickness Kr*Wg (default 1.):                              
READVDU(OPTHI,REAL,1.0)                                                         
WSL1 = 0.01;  WSL2= 0.01;  WGAP= 1.0;                                           
  *** CP111 is specific heat for 111; CP112 is for 112:                         
CP111= 473.;  CP112= 226.; HHOT= CP111*THOT; HCLD= CP112*TCLD                   
QRAD= SIGMA*(THOT**4-TCLD**4)                                                   
KROSS= OPTHI/WGAP;  SCATT= 0.0;  EMISS= KROSS-SCATT                             
MESG( Expected radiative heat flux (W/m**2):                                    
QRAD = QRAD/(0.75*OPTHI + 1./EMISH + 1./EMISC - 1.)                             
QRAD                                                                            
MESG( Expected gas temperatures near walls are:                                 
TGCL = ((SIGMA*TCLD**4 + QRAD*(1./EMISC-1.))/SIGMA)**0.25                       
TGHT = ((SIGMA*THOT**4 - QRAD*(1./EMISH-1.))/SIGMA)**0.25                       
TGCL                                                                            
TGHT                                                                            
SOLVE(H1); STORE(TMP1,PRPS,sph1,cp1)                                            
#geom1                                                                          
#radflux                                                                        
TERMS(H1,N,N,Y,N,Y,N);                                                          
INIADD= F;   FIINIT(PRPS)= 0.                                                   
FIINIT(TMP1) = 0.5*(THOT + TCLD)  ! arithmetic mean                             
FIINIT(H1)=FIINIT(TMP1) * 1000.0  ! temp * approx specific heat of air          
 TMP1=LINH                                                                      
 INIT(SOL1,PRPS,0.0,111.); INIT(SOL2,PRPS,0.0,112.)                             
                                                           
  SPEDAT(SET,EMISSIVITY,111,R,:EMISH:)                                          
  SPEDAT(SET,EMISSIVITY,112,R,:EMISC:)                                          
 COVAL(HOT,H1,1.0,HHOT); COVAL(COLD,H1,1.0,HCLD)                          
LSWEEP= 2000;  VARMAX(T3)= 1.1*THOT                 
IXMON=NX/2+1;IYMON=NY/2+1;IZMON=NZ/2+1;                                         
NXPRIN=1; NYPRIN=1; NZPRIN=1                                                    
RELAX(H1,LINRLX,0.1)                                                            
RELAX(T3,LINRLX,0.1)                                                            
OUTPUT(LTLS,Y,N,N,N,N,N); OUTPUT(WDIS,Y,N,N,N,N,N)                              
OUTPUT(WGAP,Y,N,N,N,N,N)                                                        
L($R193)                                                                        
                                                                                
CP1     = 1.005000E+03
FIINIT(CP1 ) =  CP1
VARMAX(T3)=THOT  ;VARMIN(T3)=TCLD                                               
VARMAX(TMP1)=THOT;VARMIN(TMP1)=TCLD                                             
fiinit(prps)=-1;prndtl(h1)=1.e10                                                
isg21=10  ! minimum number of sweeps, needed for selref to work                                                                                
TSTSWP=-1   
RESFAC=0.001
   Activating the following In-Form statements 
   (by moving them 2 spaces to the left) 
   enables the coefficients in the equations for TEM1 and T3 
   to be inspected.
    inform7begin                                                                  
  (stored var #3-1 is t3-tmp1)                                                    
  (stored var #rs1 is resi(h1))                                                   
  (stored var #rs3 is resi(t3))                                                   
  (stored var #ap1 is apco(h1))                                                   
  (stored var #ap3 is apco(t3))                                                   
  if(idir.eq.1) then                                                              
  (stored var #ae1 is aeco(h1))                                                   
  (stored var #ae3 is aeco(t3))                                                   
  endif                                                                           
  if(idir.eq.2) then                                                              
  (stored var #an1 is anco(h1))                                                   
  (stored var #an3 is anco(t3))                                                   
  endif                                                                           
  if(idir.eq.3) then                                                              
  (stored var #ah1 is ahco(h1))                                                   
  (stored var #ah3 is ahco(t3))                                                   
  endif                                                                           
    inform7end