P-1 radiation model

Appendix 1

P-1 radiation model: A stove with significant radiation

BY : Dr S V Zhubrin, CHAM Ltd

DATE : December 2001

FOR : Demonstration case for V3.4


The case presents the implementation in PHOENICS of P-1 radiation model. It is the simplest case of the more general P-N model based on the first order spherical harmonic expansion of the radiation intensity. The model is capable of taking into account the effect of anisotropic scattering in absorbing, emmiting and scattering media.

The model is easy to implement and solve with little computing efforts. It can also be applied to complex BFC geometries and is easily extendable to handle the particulate effects.

The model offers the advantages of simplicity, high computational efficiency and relatively good accuracy, if the optical thickness is not too small.

The model is applied to the reactive flow in a three-dimensional gas-fired stove.


The independent variables of the problem are the three components of cartesian coordinate system, namely X, Y and Z.

The main dependent (solved for) variables are:

Turbulence and combustion models

The K-epsilon model, KEMODL, closed by wall functions is used to calculate the distribution of turbulence energy and its dissipation rate from which the turbulence viscosity is deduced.

Combustion is treated as a single-step irreversible diffusion-controlled chemical reaction with a infinitely fast rate between fuel and oxidant. The gas composition and its enthalpy are related to the mixture fraction according to the Simple Chemical Reaction Scheme, SCRS, concept.

Properties and auxiliary relations

The gas density is computed from the local pressures, gas temperatures and local mixture molecular masses.

The specific enthalpies are related to gas temperatures, fuel mass fraction and the heat of combustion.

P-1 model of thermal radiation

Thermal radiation is modelled by the expanding the radiation intensity in terms of first order spherical harmonics. Assuming that only four terms representing the moments of the intensity are used, this leads to a incident radiation, CRAD or RI, W m-2,conservation equation:

div ( GradgradRI ) + a ( 4sT4 - RI) = 0

where s = 5.68 10-8, is Steffan-Boltzman constant, W m-2 K-4.

The exchange coefficient, Grad, is expressed by:

Grad = ( 3(a+s)- Cgs )-1


The symmetry factor is used to model anisotroping scattering by means of a linear-anisotropic scattering phase function. Cg ranges from -1 to +1 and represents the amount of radiation scattered in forward direction.

A positive value indicates that more radiant energy is scattered forward than backward with Cg=1 corresponding to complete forward scattering.

A negative value means that more radiant energy is scattered backward than forward with Cg=-1 standing for complete backward scattering.

A zero value of Cg defines the scattering that is equally likely in all directions, i.e. isotropic scattering

The volumetric source term for the mixture enthalpy, due to radiation, is then given by:

SH1,rad = a (RI - 4sT4)

For symmetry planes and perfectly-reflecting boundaries, the radiation boundary conditions are assumed to be zero flux type.

For the incident radiation equation, the following boundary sources per unit area are used at the walls:

SR, wall = 0.5 ew (4sTw4 - RI)(2 - ew)-1

where ew is the wall emmisivity.

The sources for incident radiation at the inlets and outlets are computed in the manner similar to the walls.

Often, it is safe to assume that the emmisivity of all flow inlets and outlets is unity (black body absorption). If the temperature outside the inlet or outlet considerably differs from that in the enclosure, the different temperatures should be used for radiation and convection fluxes at inlet and outlets.



The furnace chamber comprised a rectangular enclosure. The fuel, CH4, and oxidant, air, are introduced through separate streams. The fuel is introduced through the duct on the middle of the bottom wall, and the heated air is introduced through the inlet slots surrounding the fuel duct.

The gases are ignited on the entry and steady state combustion takes place. The combustion products with constant absorption coefficient of 0.4 and constant scattering coefficient of 0.1 flow, in a swirling motion, towards four exhaust openings located at the corners of the top wall.

The furnace is symmetrical, and only quarter of it is simulated and represented by computational domain as shown here.

Boundary conditions

The walls, inlets and outlets are considered as perfectly reflecting boundaries.

At all inlets, values are given of all dependent variables together with the prescribed flow rates.

At the outlet, the fixed exit pressure is set equal to zero and the computed pressures are relative to this pressure.

At the walls, standard "wall-functions" are used for the gas velocities, and the condition of zero flux is assumed for all other dependent variables.


The plots show the flow distribution, gas temperature and mixture composition within the furnace, as represented by the model in the absence of radiation and with radiation included.

Pictures are as follows :


All model settings have been made by PIL commands and PLANT settings of PHOENICS 3.4