My Project
programmer's 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.

Accessing variables and properties in Fortran:
  • 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).

  • Alternatively, if no global index is pointing to the searched field, variables and properties can be accessed using their field names:
Accessing variables and properties in C:
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.

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 just an allocatable array in Fortran while it is mapped as a field in C.

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(ivoidf), 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_s(ivarfl(ir11), cvar_r11) CS_F_(r11)->val Reynolds stress component $ R_{xx} $
call field_get_val_s(ivarfl(ir22), cvar_r22) CS_F_(r22)->val Reynolds stress component $ R_{yy} $
call field_get_val_s(ivarfl(ir33), cvar_r33) CS_F_(r33)->val Reynolds stress component $ R_{zz} $
call field_get_val_s(ivarfl(ir12), cvar_r12) CS_F_(r12)->val Reynolds stress component $ R_{xy} $
call field_get_val_s(ivarfl(ir23), cvar_r23) CS_F_(r23)->val Reynolds stress component $ R_{yz} $
call field_get_val_s(ivarfl(ir13), cvar_r13) CS_F_(r13)->val Reynolds stress component $ R_{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

Specific physics

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
call field_get_val_s(ivarfl(isca(itotwt)), cvar_totwt) 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 incs_field_pointer.h.

C code Description
->val Radiative luminance
->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
->val Radiative incident radiative flux density
->val Wall thermal conductivity
->val Wall thickness
->val Wall emissivity
->val Boundary radiative flux
->val Boundary radiative convective flux
->val Radiative exchange coefficient

Eulerian-Eulerian multiphase flows

Defined incs_field_pointer.h.

C code Description
->val Non-condensable gas mass fractions
->val Particles turbulent kinetic energy Q2
->val Covariance of the turbulent Q12
->val XX component of qfp
->val XY component of qfp
->val XZ component of qfp
->val YX component of qfp
->val YY component of qfp
->val YZ component of qfp
->val ZX component of qfp
->val ZY component of qfp
->val ZZ component of qfp
->val Interfacial mass transfer
->val Interfacial area
->val Droplets x2
->val Droplets Sauter mean diameter
->val Drag between phases
->val Added mass
->val Thermal diffusivity
->val Turbulent thermal diffusivity
->val dRho over dP
->val dRho over dH
->val Turbulent tau12 for particles
->val Particles lift
->val Particles turbulent dispersion
->val Particles drift velocity