6.3
general documentation
Examples of data settings for radiative transfers

Activation of the module

The module can be activated in the usppmo routine in cs_user_parameters.f90. The corresponding keyword is iirayo in the cs_glob_rad_transfer_params structure.

This member can take the values:

  • iirayo = 0: module desactivated.
  • iirayo = 1: the module is activated and the Discrete Ordinates Method is used.
  • iirayo = 2: the module is activated and the P1 model is used.

Radiation module specific parameters.

When the module is activated, its specific input parameters should be set in the cs_user_radiative_transfer_parameters function of the cs_user_radiative_transfer.c file.

Calculation options for the radiative transfer module.

Radiative transfer parameters may be defined using the cs_user_radiative_transfer_parameters function.

/* indicate whether the radiation variables should be
initialized (=0) or read from a restart file (=1) */
/* period of the radiation module */
/* Quadrature Sn (n(n+2) directions)
1: S4 (24 directions)
2: S6 (48 directions)
3: S8 (80 directions)
Quadrature Tn (8n^2 directions)
4: T2 (32 directions)
5: T4 (128 directions)
6: Tn (8*ndirec^2 directions)
*/
/* Number of directions, only for Tn quadrature */
/* Method used to calculate the radiative source term:
- 0: semi-analytic calculation (required with transparent media)
- 1: conservative calculation
- 2: semi-analytic calculation corrected
in order to be globally conservative
(If the medium is transparent, the choice has no effect) */
/* Verbosity level in the log concerning the calculation of
the wall temperatures (0, 1 or 2) */
/* Verbosity mode for the radiance (0, 1 or 2) */
/* Compute the absorption coefficient through Modak (if 1 or 2),
or do not use Modak (if 0).
Useful ONLY when gas or coal combustion is activated
- imodak = 1: ADF model with 8 wave length intervals
- imodak = 2: ADF model with 50 wave length intervals */
/* Compute the absorption coefficient via ADF model
Useful ONLY when coal combustion is activated
imoadf = 0: switch off the ADF model
imoadf = 1: switch on the ADF model (with 8 bands ADF08)
imoadf = 2: switch on the ADF model (with 50 bands ADF50) */
/* Compute the absorption coefficient through FSCK model (if 1)
Useful ONLY when coal combustion is activated
imfsck = 1: activated
imfsck = 0: not activated */
/* Activate 3D radiative models for atmospheric flows
atmo_model |= CS_RAD_ATMO_3D_INFRARED: Infrared
atmo_model |= CS_RAD_ATMO_3D_DIRECT_SOLAR: direct solar
atmo_model |= CS_RAD_ATMO_3D_DIFFUSE_SOLAR: diffuse solar
*/

Radiative transfer boundary conditions

Sketch of thermal flux in boundary walls

The radiative boundary condition is based on the calculation of a new wall temperature. This temperature is computed with a thermal flux balance:

\[{ Q_{conduction} = Q_{convection} + (Q_{rayt_{absorption}} - Q_{rayt_{emission}}}) \]

Therefore :

\[ \dfrac{xlamp}{epap} (T_{fluid} - T_{wall}) = h_{fluid} (T_{fluid} - T_{wall}) + epsp (Q_{incid} - \sigma * T_{wall}) \]

Note
In Code_Saturne the flux is positive when it is oriented from inside to outside.
Body Emissivity
polished steel 0.06
oxidized steel 0.80
steel rough 0.94
polished aluminium 0.04
oxidiezd aluminium (inside) 0.09
oxidized aluminium (wet air) 0.90
brick 0.93
concrete 0.93
paper 0.8 to 0.9
water 0.96

Boundary faces identification

Boundary faces may be identified using the getfbr function, or preferrably, through boundary zones, defined using the GUI or the cs_user_zones function..

Initialization and finalization

The following declaration and initialization block needs to be added for the following examples:

Remaining initialisation

ivar: number of the thermal variable

const cs_lnum_t ivart
= cs_field_get_key_int(fth, cs_field_key_id("variable_id")) - 1;

Min and Max values for the wall temperatures (clipping otherwise)

$ T_{min} $ and $T_{max} $ are given in Kelvin.

*tmin = 0.0;
*tmax = cs_math_big_r + tkelvi;

Assign boundary conditions to boundary wall

Zone definitions

For each boundary face face_id, a specific output (logging and postprocessing) zone id may be assigned. This allows realizing balance sheets by treating them separately for each zone. By default, the output zone id is set to the general (input) zone id associated to a face.

To access output zone ids (both for reading and modifying), use the cs_boundary_zone_face_zone_id function.

Wall characteristics

Warning
The unit of the temperature is the Kelvin

Mandatory data

  • isothp(ifac) boundary face type
    • itpimp -> Gray wall with fixed inside temperature
    • ipgrno -> Gray wall with fixed outside temperature
    • iprefl -> Reflecting wall with fixed outside temperature
    • ifgrno -> Gray wall with fixed conduction flux
    • ifrefl -> Reflecting wall with fixed conduction flux
  • tintp(ifac) inside wall temperature (Kelvin) initialize thwall at the first time step. If isothp = itpimp, the value of thwall is fixed to tintp In the other case, tintp is only for initialization.

Other data (depending of the isothp)

  • rcodcl = conduction flux
  • epsp = emissivity
  • xlamp = conductivity ( $W.m^{-1}.K^{-1}$)
  • epap = thickness ( $m$)
  • textp = outside temperature ( $K$)

Examples of boundary conditions

Here is a list of examples:

Gray or black wall with profil of fixed inside temperature

For wall boundary faces, selection criteria: color 1

zone = cs_boundary_zone_by_name("wall_1");
for (cs_lnum_t ilelt = 0; ilelt < zone->n_elts; ilelt++) {
cs_lnum_t face_id = zone->elt_ids[ilelt];
if (bc_type[face_id] == CS_SMOOTHWALL) {
/* logging zone number */
izfrdp[face_id] = 51;
/* Type of condition: gray or black wall with fixed inside temperature */
isothp[face_id] = cs_glob_rad_transfer_params->itpimp;
/* Emissivity */
epsp[face_id] = 0.1;
/* Fixed inside temperature */
tintp[face_id] = 200 + tkelvi;
}
}

Gray or black wall with fixed outside temperature \f$ T_{ext} \f$

For wall boundary faces, selection criteria: color 2

zone = cs_boundary_zone_by_name("wall_2");
for (cs_lnum_t ilelt = 0; ilelt < zone->n_elts; ilelt++) {
cs_lnum_t face_id = zone->elt_ids[ilelt];
if (bc_type[face_id] == CS_ROUGHWALL) {
/* logging zone number */
izfrdp[face_id] = 52;
/* Type of condition: gray or black wall with fixed
outside temperature TEXTP */
isothp[face_id] = cs_glob_rad_transfer_params->ipgrno;
/* Emissivity */
epsp[face_id] = 0.9;
/* Conductivity (W/m/K)*/
xlamp[face_id] = 3.0;
/* Thickness (m)*/
epap[face_id] = 0.1;
/* Fixed outside temperature: 473.15 K */
textp[face_id] = 200. + tkelvi;
/* Initial inside temperature: 473.15 K */
tintp[face_id] = 200. + tkelvi;
}
}

Reflecting wall (\f$ epsp = 0 \f$) with fixed outside temperature \f$ T_{ext} \f$

For wall boundary faces, selection criteria: color 3

zone = cs_boundary_zone_by_name("wall_3");
for (cs_lnum_t ilelt = 0; ilelt < zone->n_elts; ilelt++) {
cs_lnum_t face_id = zone->elt_ids[ilelt];
if (bc_type[face_id] == CS_SMOOTHWALL) {
/* log zone number */
izfrdp[face_id] = 53;
/* Type of condition: reflecting wall with fixed outside temperature TEXTP */
isothp[face_id] = cs_glob_rad_transfer_params->iprefl;
/* Conductivity (W/m/K) */
xlamp[face_id] = 3.0;
/* Thickness (m)*/
epap[face_id] = 0.10;
/* Fixed outside temperature: 473.15 K */
textp[face_id] = 200.0 + tkelvi;
/* Initial inside temperature: 473.15 K */
tintp[face_id] = 200.0 + tkelvi;
}
}

Gray or black wall and fixed conduction flux through the wall

For wall boundary faces which have the color 4:

\[ \begin{array}{rcl} \frac{\texttt{xlamp}}{\texttt{epap}} \cdot (T_{wall} - T_{ext}) &=& \text{fixed conduction flux in } W.m^{-2} \\ &=& \texttt{rodcl(ifac,ivar,3)} \end{array} \]

If the conduction flux is zero then the wall is adiabatic. The array $ \texttt{rcodcl(ifac,ivar,3)}$ has the value of the flux.
Flux density (< 0 if gain for the fluid)

  • For temperature $T$, in $ W.m^{-2}$:

\[ rcodcl(ifac,ivar,3)=C_p (viscls+\frac{visct}{\sigma})\cdot \grad{T}\cdot \vect{n} \]

  • For enthalpy $h$, in $ W.m^{-2} $:

    \[ RCODC(IFAC,IVAR,3)=(viscls+\frac{visct}{\sigma})\cdot \grad{H} \cdot \vect{n}\]

zone = cs_boundary_zone_by_name("wall_4");
for (cs_lnum_t ilelt = 0; ilelt < zone->n_elts; ilelt++) {
cs_lnum_t face_id = zone->elt_ids[ilelt];
if (bc_type[face_id] == CS_SMOOTHWALL) {
/* log zone number */
izfrdp[face_id] = 54;
/* Type of condition: gray or black wall with fixed conduction
flux through the wall */
isothp[face_id] = cs_glob_rad_transfer_params->ifgrno;
/* Emissivity */
epsp[face_id] = 0.9;
/* Conduction flux (W/m2) */
rcodcl[face_id + ivart * n_b_faces + 2 * nvar * n_b_faces] = 0.0;
/* Initial inside temperature: 473.15 K */
tintp[face_id] = 200.0 + tkelvi;
}
}

Reflecting wall and fixed conduction flux through the wall

For wall boundary faces which have the color 5:

\[ \frac{xlamp}{epap} \cdot (T_{wall} - T_{ext}) = \text{fixed conduction flux} \]

and $ epsp = 0 $

If the conduction flux is zero then the wall is adiabatic. Flux density (< 0 if gain for the fluid)

  • For temperatures $T$, in $ W.m^{-2} $:

    \[ rcodcl(ifac,ivar,3) = C_p (viscls+\frac{visct}{\sigma}) \cdot \grad{T}\cdot \vect{n} \]

  • For enthalpies $h$, in $ W.m^{-2} $:

    \[ rcodcl(ifac,ivar,3) = (viscls+\frac{visct}{\sigma}) \cdot \grad{H} \cdot \vect{n} \]

zone = cs_boundary_zone_by_name("wall_5");
for (cs_lnum_t ilelt = 0; ilelt < zone->n_elts; ilelt++) {
cs_lnum_t face_id = zone->elt_ids[ilelt];
if (bc_type[face_id] == CS_SMOOTHWALL) {
/* log zone number */
izfrdp[face_id] = 55;
/* Type of condition: reflecting wall with fixed conduction
flux through the wall */
isothp[face_id] = cs_glob_rad_transfer_params->ifrefl;
/* Conduction flux (W/m2)*/
rcodcl[face_id + ivart * n_b_faces + 2 * nvar * n_b_faces ] = 0.0;
/* Initial inside temperature: 473.15 K */
tintp[face_id] = 200.0 + tkelvi;
}
}

Absorption coefficient and net radiation flux

The absorption coefficient and the net radiation flux for the radiative module can be defined in cs_user_radiative_transfer.c through the cs_user_rad_transfer_absorption and Net radiation flux subroutines.

Absorption coefficient

The absorption coefficient is defined in cs_user_rad_transfer_absorption.

Computation of the absorption coefficient

/*
* Absorption coefficient of the medium (m-1)
*
* In the case of specific physics (gas/coal/fuel combustion, elec),
* Ck must not be defined here (it is determined automatically, possibly
* from the parametric file)
*
* In other cases, Ck must be defined (it is zero by default)
*/
for (cs_lnum_t cell_id = 0; cell_id < cs_glob_mesh->n_cells; cell_id++)
ck[cell_id] = 0.;
}

Net radiation flux

The net radiation flux is computed in Net radiation flux.

Local variables to be added

Initialization

At the end of the subroutine, if iok is different from zero, some faces have been forgotten and the calculation stops.

/* Initializations */

Computation of the net radiation flux

/* Net flux dendity for the boundary faces
* The provided examples are sufficient in most of cases.*/
/* If the boundary conditions given above have been modified
* it is necessary to change the way in which density is calculated from
* the net radiative flux consistently.*/
/* The rule is:
* the density of net flux is a balance between the emitting energy from a
* boundary face (and not the reflecting energy) and the absorbing radiative
* energy. Therefore if a wall heats the fluid by radiative transfer, the
* net flux is negative */
for (cs_lnum_t ifac = 0; ifac < cs_glob_mesh->n_b_faces; ifac++) {
/* Wall faces */
if ( bc_type[ifac] == CS_SMOOTHWALL
|| bc_type[ifac] == CS_ROUGHWALL)
net_flux[ifac] = eps[ifac] * (qincid[ifac] - stephn * pow(twall[ifac], 4));
/* Symmetry */
else if (bc_type[ifac] == CS_SYMMETRY)
net_flux[ifac] = 0.0;
/* Inlet/Outlet */
else if ( bc_type[ifac] == CS_INLET
|| bc_type[ifac] == CS_CONVECTIVE_INLET
|| bc_type[ifac] == CS_OUTLET
|| bc_type[ifac] == CS_FREE_INLET) {
net_flux[ifac] = qincid[ifac] - cs_math_pi * coefap[ifac];
net_flux[ifac] = 0.0;
}
/* Stop if there are forgotten faces */
else
(__FILE__, __LINE__, 0,
"In %s:\n"
" non-handled boundary faces for net flux calculation\n\n"
" Last face: %10ld; zone = %d; nature = %d\n",
__func__,
(long)ifac,
bc_type[ifac]);
}
cs_mesh_t::n_cells
cs_lnum_t n_cells
Definition: cs_mesh.h:96
cs_rad_transfer_params_t::ipgrno
int ipgrno
Definition: cs_rad_transfer.h:117
cs_rad_transfer_params_t::itpimp
int itpimp
Definition: cs_rad_transfer.h:116
cs_math_pi
const cs_real_t cs_math_pi
cs_boundary_zone_face_zone_id
const int * cs_boundary_zone_face_zone_id(void)
Return pointer to zone id associated with each boundary face.
Definition: cs_boundary_zone.c:828
cs_rad_transfer_params_t::ndirec
int ndirec
Definition: cs_rad_transfer.h:107
cs_boundary_zone_face_class_id
int * cs_boundary_zone_face_class_id(void)
Get pointer to optional boundary face class ids.
Definition: cs_boundary_zone.c:958
cs_rad_transfer_params_t::verbosity
int verbosity
Definition: cs_rad_transfer.h:99
cs_rad_transfer_params_t::ifgrno
int ifgrno
Definition: cs_rad_transfer.h:119
cs_zone_t
Definition: cs_zone.h:55
cs_rad_transfer_params_t::imodak
int imodak
Definition: cs_rad_transfer.h:100
dimens::nvar
integer, save nvar
number of solved variables (must be lower than nvarmx)
Definition: dimens.f90:42
cs_restart_present
int cs_restart_present(void)
Check if we have a restart directory.
Definition: cs_restart.c:2092
cs_glob_rad_transfer_params
cs_rad_transfer_params_t * cs_glob_rad_transfer_params
cs_rad_transfer_params_t::i_quadrature
int i_quadrature
Definition: cs_rad_transfer.h:106
CS_FREE_INLET
@ CS_FREE_INLET
Definition: cs_parameters.h:114
CS_RAD_ATMO_3D_DIRECT_SOLAR
@ CS_RAD_ATMO_3D_DIRECT_SOLAR
Definition: cs_rad_transfer.h:78
cs_real_t
double cs_real_t
Floating-point value.
Definition: cs_defs.h:304
cs_mesh_t::n_b_faces
cs_lnum_t n_b_faces
Definition: cs_mesh.h:98
cs_glob_physical_model_flag
int cs_glob_physical_model_flag[CS_N_PHYSICAL_MODEL_TYPES]
Definition: cs_physical_model.c:109
cs_rad_transfer_params_t::type
cs_rad_transfer_model_t type
Definition: cs_rad_transfer.h:95
cs_math_big_r
const cs_real_t cs_math_big_r
CS_INLET
@ CS_INLET
Definition: cs_parameters.h:102
bft_error
void bft_error(const char *const file_name, const int line_num, const int sys_error_code, const char *const format,...)
Calls the error handler (set by bft_error_handler_set() or default).
Definition: bft_error.c:193
cstphy::tkelvi
double precision tkelvi
Temperature in Kelvin correponding to 0 degrees Celsius (= +273,15)
Definition: cstphy.f90:44
cs_glob_mesh
cs_mesh_t * cs_glob_mesh
CS_OUTLET
@ CS_OUTLET
Definition: cs_parameters.h:103
cs_rad_transfer_params_t::imoadf
int imoadf
Definition: cs_rad_transfer.h:101
CS_CONVECTIVE_INLET
@ CS_CONVECTIVE_INLET
Definition: cs_parameters.h:116
cstphy::stephn
double precision stephn
Stephan constant for the radiative module in .
Definition: cstphy.f90:53
eps
@ eps
Definition: cs_field_pointer.h:71
cs_thermal_model_field
cs_field_t * cs_thermal_model_field(void)
Definition: cs_thermal_model.c:207
CS_RAD_TRANSFER_DOM
@ CS_RAD_TRANSFER_DOM
Definition: cs_rad_transfer.h:53
cs_field_get_key_int
int cs_field_get_key_int(const cs_field_t *f, int key_id)
Return a integer value for a given key associated with a field.
Definition: cs_field.c:2991
CS_SMOOTHWALL
@ CS_SMOOTHWALL
Definition: cs_parameters.h:105
cs_rad_transfer_params_t::ifrefl
int ifrefl
Definition: cs_rad_transfer.h:120
cs_rad_transfer_params_t::atmo_model
int atmo_model
Definition: cs_rad_transfer.h:124
cs_rad_transfer_params_t::imfsck
int imfsck
Definition: cs_rad_transfer.h:103
cs_rad_transfer_params_t::iimpar
int iimpar
Definition: cs_rad_transfer.h:98
cs_lnum_t
int cs_lnum_t
local mesh entity id
Definition: cs_defs.h:298
CS_RAD_ATMO_3D_INFRARED
@ CS_RAD_ATMO_3D_INFRARED
Definition: cs_rad_transfer.h:77
cs_rad_transfer_params_t::restart
int restart
Definition: cs_rad_transfer.h:111
cs_physical_constants_celsius_to_kelvin
const double cs_physical_constants_celsius_to_kelvin
Definition: cs_physical_constants.c:375
cs_boundary_zone_by_name
const cs_zone_t * cs_boundary_zone_by_name(const char *name)
Return a pointer to a boundary zone based on its name if present.
Definition: cs_boundary_zone.c:726
cs_zone_t::elt_ids
const cs_lnum_t * elt_ids
Definition: cs_zone.h:65
cs_field_key_id
int cs_field_key_id(const char *name)
Return an id associated with a given key name.
Definition: cs_field.c:2497
cs_zone_t::n_elts
cs_lnum_t n_elts
Definition: cs_zone.h:64
CS_RAD_ATMO_3D_DIFFUSE_SOLAR
@ CS_RAD_ATMO_3D_DIFFUSE_SOLAR
Definition: cs_rad_transfer.h:79
cs_rad_transfer_params_t::iprefl
int iprefl
Definition: cs_rad_transfer.h:118
CS_SYMMETRY
@ CS_SYMMETRY
Definition: cs_parameters.h:104
cs_rad_transfer_params_t::nfreqr
int nfreqr
Definition: cs_rad_transfer.h:112
cs_physical_constants_stephan
const double cs_physical_constants_stephan
Definition: cs_physical_constants.c:380
CS_PHYSICAL_MODEL_FLAG
@ CS_PHYSICAL_MODEL_FLAG
Definition: cs_physical_model.h:60
cs_field_t
Field descriptor.
Definition: cs_field.h:125
CS_RAD_TRANSFER_P1
@ CS_RAD_TRANSFER_P1
Definition: cs_rad_transfer.h:54
cs_rad_transfer_params_t::idiver
int idiver
Definition: cs_rad_transfer.h:105
CS_ROUGHWALL
@ CS_ROUGHWALL
Definition: cs_parameters.h:106