8.0
general documentation
How to access and manage variables and properties using the cs_field API?

Variables and properties can be accessed both in Fortran and in C using the cs_field API.

First, let us specify a few conventions used in code_saturne:

  • The dynamic variables designation includes velocity, pressure, and variables related to the turbulence model (k, ε, Rij, ω, ϕ, $ \overline{f} $, α, νt).
  • The designation “scalar” refers to (usually scalar) variables which are solution of an advection equation, apart from the dynamic variables listed above; for instance the temperature, scalars which may be passive or not, user-defined or model-defined.

    The mean value of the square of the fluctuations of a “scalar” is a “scalar”, too.

    Scalars may be divided into two groups: user-defined scalars and model-defined (sometimes referred to as “specific physics”) scalars.

  • Depending on the thermal model used; the solved thermal scalar variable may be the temperature, the enthalpy, or in the case of the compressible module, the total energy. When the temperature is not the solved thermal variable, it is usually available as a property.
  • When a user scalar represents the mean of the square of the fluctuations of another scalar (i.e. the average $ \overline{\varphi^\prime\varphi^\prime} $ for a fluctuating scalar ϕ). This can be made either via the interface or by declaring that scalar using cs_parameters_add_variable_variance in cs_user_parameters.c.
Accessing variables and properties in C:
  • The most general way of accessing a field is to use its name:

    cs_real_t *cvar_pr = cs_field_by_name("pressure")->val;
    then:
    cvar_pr[cell_id]

  • The main variables and properties can be accessed using the CS_F_ macro:

  • Indexed variables (such as user scalars) and indexed properties are accessed as:
    CS_FI_(name,ii-1)->val[cell_id].

  • The thermal scalar can be accessed using the cs_thermal_model_field function:
    field_t *tf = cs_thermal_model_field();
    then:
    tf->val[cell_id].
    if there is no thermal scalar, cs_thermal_model_field returns NULL.

  • In any case, all variables can be accessed using the function cs_field_by_name:
    cs_real_t *cvar_pr = cs_field_by_name("pressure")->val

  • A field's scalar_id key-word should be ≥ 0 if it is a scalar.

  • When a field is user-defined (rather than model-defined), its type flag should also match the CS_FIELD_USER mask:
    if (f->type & CS_FIELD_USER) ...
Remarks
Note that indexes in C begin at 0, while indexes in Fortran begin at 1. Thus, Fortran and C loop counters are related in the following as:
cell_id = iel-1.
Accessing variables and properties in Fortran:
  • The most general way of accessing a field is to use its name:

    call field_get_val_s_by_name("pressure", cvar_pr)
    pres = cvar_pr(iel)


  • Both variables and properties can be accessed via the cs_field API, using the global field indices and indirections arrays, as in the following examples:
    • For one-dimensional arrays:

      call field_get_val_s(ivarfl(ipr), cvar_pr)
      pres = cvar_pr(iel)
      ,

      call field_get_val_s(icp, cpro_cp)
      cp = cpro_cp(iel)


      The values of scalar variable can be accessed as follows:

      call field_get_val_s(ivarfl(isca(iscalt)), cvar_scalt)
      temp = cvar_scalt(iel)


    • For interleaved multidimensional arrays:

      call field_get_val_v(ivarfl(iu), cvar_vel)
      ux = cvar_vel(1,iel)


      ipr, iu are here variable indices and the array ivarfl allows to get the corresponding field indices.
      isca(iscalt) is also a variable index (iscalt is the scalar index of the thermal scalar).
      icp is directly a field index (there is no equivalent to the array ivarfl for field of type properties).

  • the solved variable number for scalar j, where 1 <= j <= nscal, is given by the iscal array. So to access the associated field, use:
    ivarfl(iscal(j))).

  • For model-defined scalars, the total scalar number for model scalar j, where 1 <= j <= nscapp, is given by the iscapp array. So to access the associated field, use:
    ivarfl(iscal(iscapp(j)))).

  • The thermal scalar can be accessed using the iscalt index in the list of the scalars. It's variable number is thus isca(iscalt) .
    if there is no thermal scalar, iscalt is equal to -1.

Cross-reference tables are available for the variables and properties of the standard solver and the specific physics features:

Variables

The Fortran variables indexes are defined in the files numvar.f90 (with the exception of ihm and iscal, which are respectively defined in ppincl.f90 and optcal.f90) and the C variables names are defined in cs_field_pointer.h.
Note that dt is accessible by its name (using the cs_field_by_name family of functions), but not through a permanent index..

Fortran code C code Description
call field_get_val_s(ivarfl(ipr), cvar_pr) CS_F_(p)->val Pressure
call field_get_val_v(ivarfl(iu), cvar_vel) CS_F_(vel)->val Velocity
call field_get_val_s(ivarfl(ivolf2), cvar_voidf) CS_F_(void_f)->val Void fraction for Volume of Fluid model
call field_get_val_s(ivarfl(ik ), cvar_k ) CS_F_(k)->val Turbulent kinetic energy $ k $
call field_get_val_s(ivarfl(iep ), cvar_eps) CS_F_(eps)->val Turbulent dissipation $ \varepsilon $
call field_get_val_v(ivarfl(irij), cvar_rij) CS_F_(rij)->val Reynolds stress tensor $ R $ (symmetric $ R_{ij} $ components, in the following order for each cell: xx, yy, zz, xy, yz, xz)
call field_get_val_s(ivarfl(iphi), cvar_phi) CS_F_(phi)->val $ \phi $ for $ \phi-f_b $ model
call field_get_val_s(ivarfl(ifb ), cvar_fb ) CS_F_(f_bar)->val $ f_b $ for $ \phi-f_b $ model
call field_get_val_s(ivarfl(ial ), cvar_al ) CS_F_(alp_bl)->val $ \alpha $ for $ Bl-v^2-k $
or EBRSM model
call field_get_val_s(ivarfl(iomg), cvar_omg) CS_F_(omg)->val $ \omega $ for $ k-\omega $ SST model
call field_get_val_s(ivarfl(inusa), cvar_nusa) CS_F_(nusa)->val $ \widetilde{\nu}_T $ for Spalart-Allmaras
call field_get_val_v(ivarfl(iuma), cvar_mesh_v) CS_F_(mesh_u)->val Mesh velocity
call field_get_val_s(ivarfl(isca(ihm)), cvar_hm) CS_F_(h)->val Enthalpy
call field_get_val_s(ivarfl(isca(iscalt)), cvar_scalt) CS_F_(t)->val Temperature

Properties

These properties are defined in the files numvar.f90 and cs_field_pointer.h.

Fortran code C code Description
dt CS_F_(dt)->val Local time step
call field_get_val_s(iviscl, cpro_viscl) CS_F_(mu)->val Molecular viscosity
call field_get_val_s(ivisct, cpro_visct) CS_F_(mu_t)->val Turbulent dynamic viscosity
call field_get_val_s(icp, cpro_cp) CS_F_(cp)->val Specific heat
call field_get_val_s(icrom, cpro_crom) CS_F_(rho)->val Density (at cells)
call field_get_val_s(ibrom, bpro_rho) CS_F_(rho_b)->val[face_id] Density (at boundary faces)
call field_get_val_s(ismago, cpro_smago) cs_real_t *cpro_smago = cs_field_by_name("smagorinsky_constant^2")->val Field id of the anisotropic turbulent viscosity
call field_get_val_s(icour, cpro_cour) cs_real_t *cpro_cour = cs_field_by_name("courant_number")->val Courant number
call field_get_val_s(ifour, cpro_four) cs_real_t *cpro_four = cs_field_by_name("fourier_number")->val Fourier number
call field_get_val_s(iprtot, cpro_prtot) cs_real_t *cpro_prtot = cs_field_by_name("total_pressure")->val Total pressure at cell centers
call field_get_val_s(ivisma, cpro_visma_s) cs_real_t *cpro_visma_s = cs_field_by_name("mesh_viscosity")->val Mesh velocity viscosity (scalar) for the ALE module
call field_get_val_v(ivisma, cpro_visma_s) cs_real_t *cpro_visma_v = cs_field_by_name("mesh_viscosity")->val Mesh velocity viscosity (vector) for the ALE module
call field_get_val_s(itsrho), cpro_tsrho ) cs_real_t *cpro_tsrho = cs_field_by_name("dila_st")->val Global dilatation source terms
call field_get_val_s(ibeta), cpro_beta ) cs_real_t *cpro_beta = cs_field_by_name("thermal_expansion")->val Thermal expansion coefficient
call field_get_val_s(ipori, cpro_ipori) CS_F_(poro)->val Porosity
call field_get_val_v(iporf, cpro_iporf) CS_F_(t_poro)->val Tensorial porosity
call field_get_val_v(iforbr, bpro_forbr) cs_real_t *bpro_forbr = cs_field_by_name("boundary_forces")->val Field id of the stresses at boundary
call field_get_val_s(iyplbr, bpro_yplus) cs_real_t *bpro_yplus = cs_field_by_name("yplus")->val Field id of $y^+$ at boundary
call field_get_val_v(idtten, dttens) cs_real_t *dttens = cs_field_by_name("dttens")->val Field id for the dttens tensor
call field_get_val_s(itempb, t_b) CS_F_(t_b)->val Boundary temperature

Some physical properties such as specific heat or diffusivity are often constant (depending on the model or user parameters). In that case, these properties are stored as a simple real numbers rather than in a field over all mesh cells.

  • This is the case for the specific heat Cp.
    • If Cp is constant, it is based on the reference value cs_glob_fluid_properties->cp0 (cp0 in Fortran).
      When this is the case, cs_glob_fluid_properties->icp (mapped as icp in Fortran) should remain set to 0.
      When icp is initialized to 1 (by the GUI, or in cs_user_parameters.c, it is automatically reset to the id of the cell-based property field referenced in the above table.
  • This is the same for the diffusivity K of each scalar.
    • When K is constant, its value is based on the field's reference diffusivity, accessible through the scalar field's diffusivity_ref keyword.

      When it is variable, the matching field can be specified and accessed using the base scalar field's diffusivity_id key (accessible using cs_field_key_id("diffusivity_id") in C, or kivisl in Fortran).
      For example, for a scalar field f:
      int k_f_id = cs_field_get_key_int(f, cs_field_key_id("diffusivity_id"));
      cs_field_t *kf = cs_field_by_id(k_f_id);

Reference physical values

Reference physical values represent either the fluid properties if they are constant, either simple mean values for the initialization if properties are variable and defined in cs_user_physical_properties. Reference values can be set in the GUI or in files cs_user_parameters.f90 and cs_user_parameters.c.

Fortran code C code Description
p0 cs_glob_fluid_properties->p0 Reference total pressure
ro0 cs_glob_fluid_properties->ro0 Reference density
viscl0 cs_glob_fluid_properties->viscl0 Reference molecular dynamic viscosity
cp0 cs_glob_fluid_properties->cp0 Reference specific heat

Specific physical models

Atmospheric

Defined in optcal.f90, atincl.f90, atvarp.f90 and cs_field_pointer.h.

Fortran code C code Description
call field_get_val_s(ivarfl(isca(iscalt)), cvar_scalt) CS_F_(pot_t)->val Potential temperature
- CS_F_(totwt)->val Total water content
call field_get_val_s(ivarfl(isca(intdrp)), cvar_intdrp) CS_F_(ntdrp)->val Total number of droplets
call field_get_val_s(ivarfl(isca(isca_chem(iesp))), cvar_sc) CS_FI_(chemistry,iesp-1)->val Chemistry species (indexed)

Coal combustion

Defined in ppincl.f90, ppcpfu.f90 and cs_field_pointer.h.

Fortran code C code Description
call field_get_val_s(isca(inp(iesp)), cvar_inpcl) CS_FI_(np,iesp-1)->val Particles per kg for coal class
call field_get_val_s(isca(ixch(iesp)), cvar_xchcl) CS_FI_(xch,iesp-1)->val Reactive coal mass fraction for coal class
call field_get_val_s(isca(ixck(iesp)), cvar_xckcl) CS_FI_(xck,iesp-1)->val Coke mass fraction for coal class
call field_get_val_s(isca(ixwt(iesp)), cvar_xwtcl) CS_FI_(xwt,iesp-1)->val Water mass fraction for coal class
call field_get_val_s(isca(ih2(iesp)), cvar_h2cl) CS_FI_(h2,iesp-1)->val Mass enthalpy for coal class (permeatic case)
call field_get_val_s(isca(if1m(iesp)), cvar_f1mcl) CS_FI_(f1m,iesp-1)->val Mean value light volatiles for coal class
call field_get_val_s(isca(if2m(iesp)), cvar_f2mcl) CS_FI_(f2m,iesp-1)->val Mean value heavy volatiles for coal class
call field_get_val_s(isca(if4m), cvar_f4m) CS_F_(f4m)->val Oxydant 2 mass fraction
call field_get_val_s(isca(if5m), cvar_f5m)) CS_F_(f5m)->val Oxydant 3 mass fraction
call field_get_val_s(isca(if6m), cvar_f6m)) CS_F_(f6m)->val Water from coal drying mass fraction
call field_get_val_s(isca(if7m), cvar_f7m)) CS_F_(f7m)->val Carbon from coal oxidyzed by O2 mass fraction
call field_get_val_s(isca(if8m), cvar_f8m)) CS_F_(f8m)->val Carbon from coal gasified by CO2 mass fraction
call field_get_val_s(isca(if9m), cvar_f9m)) CS_F_(f9m)->val Carbon from coal gasified by H2O mass fraction
call field_get_val_s(isca(ifvp2m), cvar_fvp2m) CS_F_(fvp2m)->val f1f2 variance
call field_get_val_s(isca(iyco2), cvar_yco2) CS_F_(yco2)->val CO2 fraction
call field_get_val_s(isca(iyhcn), cvar_yhnc) CS_F_(yhcn)->val HCN fraction
call field_get_val_s(isca(iyno), cvar, yno) CS_F_(yno)->val NO fraction
call field_get_val_s(isca(iynh3), cvar_ynh3) CS_F_(ynh3)->val NH3 enthalpy
call field_get_val_s(isca(ihox), cvar_hox) CS_F_(hox)->val Ox enthalpy

Compressible

Defined in ppincl.f90 and cs_field_pointer.h.

Fortran code C code Description
call field_get_val_s(isca(ienerg), cvar_energ) CS_F_(e_tot)->val Total energy
call field_get_val_s(isca(itempk), cvar_tempk) CS_F_(t_kelvin)->val Temperature, in Kelvin

Electric arcs

Defined in ppincl.f90 and cs_field_pointer.h.

Fortran code C code Description
call field_get_val_s_by_name("elec_pot_r", cvar_potr) CS_F_(potr)->val Electric potential, real part
call field_get_val_s_by_name("elec_pot_i", cvar_poti) CS_F_(poti)->val Electric potential, imaginary part
call field_get_val_v_by_name("vec_potential", cvar_potva) CS_F_(potva)->val Vector potential
call field_get_val_s_by_name("esl_fraction_01", cvar_ycoel_01) CS_FI_(ycoel,iesp-1)->val Constituent mass fraction

Gas combustion

Defined in ppincl.f90 and cs_field_pointer.h.

Fortran code C code Description
call field_get_val_s(isca(ifm), cvar_fm) CS_F_(fm)->val Mixture fraction
call field_get_val_s(isca(ifp2m), cvar_fp2m) CS_F_(fp2m)->val Mixture fraction variance
call field_get_val_s(isca(ifsm), cvar_fsm) CS_F_(fsm)->val Soot mass fraction
call field_get_val_s(isca(inpm), cvar_npm) CS_F_(npm)->val Soot precursor number
call field_get_val_s(isca(iygfm), cvar_ygfm) CS_F_(ygfm)->val Fresh gas fraction
call field_get_val_s(isca(iyfm), cvar_yfm) CS_F_(yfm)->val Mass fraction
call field_get_val_s(isca(iyfp2m), cvar_yfp2m) CS_F_(yfp2m)->val Mass fraction variance
call field_get_val_s(isca(icoyfp), cvar_coyfp) CS_F_(coyfp)->val Mass fraction covariance

Radiative transfer

Defined in cs_field_pointer.h.

C code Description
CS_F_(rad_lumin)->val Radiative luminance
CS_F_(rad_q)->val Radiative flux
CS_FI_(rad_ets,iesp-1)->val Radiative flux explicit source term
CS_FI_(rad_its,iesp-1)->val Radiative flux implicit source term
CS_FI_(rad_abs,iesp-1)->val Radiative absorption
CS_FI_(rad_emi,iesp-1)->val Radiative emission
CS_FI_(rad_cak,iesp-1)->val Radiative absorption coefficient
CS_F_(qinci)->val Radiative incident radiative flux density
CS_F_(xlam)->val Wall thermal conductivity
CS_F_(epa)->val Wall thickness
CS_F_(emissivity)->val Wall emissivity
CS_F_(fnet)->val Boundary radiative flux
CS_F_(fconv)->val Boundary radiative convective flux
CS_F_(hconv)->val Radiative exchange coefficient

Eulerian-Eulerian multiphase flows

Defined in cs_field_pointer.h.

C code Description
CS_FI_(yf_ncond,inc)->val Non-condensable gas mass fractions
CS_FI_(qp,ip)->val Particles turbulent kinetic energy Q2
CS_FI_(qfp,ip)->val Covariance of the turbulent Q12
CS_FI_(qfpxx,ip)->val XX component of qfp
CS_FI_(qfpxy,ip)->val XY component of qfp
CS_FI_(qfpxz,ip)->val XZ component of qfp
CS_FI_(qfpyx,ip)->val YX component of qfp
CS_FI_(qfpyy,ip)->val YY component of qfp
CS_FI_(qfpyz,ip)->val YZ component of qfp
CS_FI_(qfpzx,ip)->val ZX component of qfp
CS_FI_(qfpzy,ip)->val ZY component of qfp
CS_FI_(qfpzz,ip)->val ZZ component of qfp
CS_FI_(gamma,ip)->val Interfacial mass transfer
CS_FI_(ia,ip)->val Interfacial area
CS_FI_(x2,ip)->val Droplets x2
CS_FI_(d32,ip)->val Droplets Sauter mean diameter
CS_FI_(drag,ipcpl)->val Drag between phases
CS_FI_(ad_mass,ip)->val Added mass
CS_FI_(th_diff,ip)->val Thermal diffusivity
CS_FI_(th_diff_t,ip)->val Turbulent thermal diffusivity
CS_FI_(drho_dp,ip)->val dRho over dP
CS_FI_(drho_dh,ip)->val dRho over dH
CS_FI_(tau12_t,ip)->val Turbulent tau12 for particles
CS_FI_(lift,ip)->val Particles lift
CS_FI_(disp_t,ip)->val Particles turbulent dispersion
CS_FI_(drift_vel,ip)->val Particles drift velocity