8.0
general documentation
cs_divergence.h File Reference
#include "cs_base.h"
#include "cs_halo.h"
#include "cs_mesh.h"
#include "cs_mesh_quantities.h"
+ Include dependency graph for cs_divergence.h:

Go to the source code of this file.

Functions

void inimav (const int *const f_id, const int *const itypfl, const int *const iflmb0, const int *const init, const int *const inc, const int *const imrgra, const int *const nswrgu, const int *const imligu, const int *const iwarnu, const cs_real_t *const epsrgu, const cs_real_t *const climgu, const cs_real_t rom[], const cs_real_t romb[], const cs_real_3_t vel[], const cs_real_3_t coefav[], const cs_real_33_t coefbv[], cs_real_t i_massflux[], cs_real_t b_massflux[])
 
void divmas (const int *const init, const cs_real_t i_massflux[], const cs_real_t b_massflux[], cs_real_t diverg[])
 
void divmat (const int *const init, const cs_real_3_t i_massflux[], const cs_real_3_t b_massflux[], cs_real_3_t diverg[])
 
void projts (const int *const init, const int *const nswrgu, const cs_real_3_t frcxt[], const cs_real_t cofbfp[], cs_real_t i_massflux[], cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t viselx[], const cs_real_t visely[], const cs_real_t viselz[])
 
void projtv (const int *const init, const int *const nswrgu, const int *const ircflp, const cs_real_3_t frcxt[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t viscel[], const cs_real_2_t weighf[], cs_real_t i_massflux[], cs_real_t b_massflux[])
 
void divrij (const int *const f_id, const int *const itypfl, const int *const iflmb0, const int *const init, const int *const inc, const int *const imrgra, const int *const nswrgu, const int *const imligu, const int *const iwarnu, const cs_real_t *const epsrgu, const cs_real_t *const climgu, const cs_real_t rom[], const cs_real_t romb[], const cs_real_6_t tensorvel[], const cs_real_6_t coefav[], const cs_real_66_t coefbv[], cs_real_3_t i_massflux[], cs_real_3_t b_massflux[])
 
void cs_mass_flux (const cs_mesh_t *m, const cs_mesh_quantities_t *fvq, int f_id, int itypfl, int iflmb0, int init, int inc, int imrgra, int nswrgu, int imligu, int iwarnu, double epsrgu, double climgu, const cs_real_t rom[], const cs_real_t romb[], const cs_real_3_t vel[], const cs_real_3_t coefav[], const cs_real_33_t coefbv[], cs_real_t *restrict i_massflux, cs_real_t *restrict b_massflux)
 Add $ \rho \vect{u} \cdot \vect{s}_\ij$ to the mass flux $ \dot{m}_\ij $. More...
 
void cs_divergence (const cs_mesh_t *m, int init, const cs_real_t i_massflux[], const cs_real_t b_massflux[], cs_real_t *restrict diverg)
 Add the integrated mass flux on the cells. More...
 
void cs_tensor_divergence (const cs_mesh_t *m, int init, const cs_real_3_t i_massflux[], const cs_real_3_t b_massflux[], cs_real_3_t *restrict diverg)
 Add the integrated mass flux on the cells for a tensor variable. More...
 
void cs_ext_force_flux (const cs_mesh_t *m, cs_mesh_quantities_t *fvq, int init, int nswrgu, const cs_real_3_t frcxt[], const cs_real_t cofbfp[], cs_real_t *restrict i_massflux, cs_real_t *restrict b_massflux, const cs_real_t i_visc[], const cs_real_t b_visc[], const cs_real_t viselx[], const cs_real_t visely[], const cs_real_t viselz[])
 Project the external source terms to the faces in coherence with cs_face_diffusion_scalar for the improved hydrostatic pressure algorithm (iphydr=1). More...
 
void cs_ext_force_anisotropic_flux (const cs_mesh_t *m, cs_mesh_quantities_t *fvq, int init, int nswrgp, int ircflp, const cs_real_3_t frcxt[], const cs_real_t cofbfp[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_6_t viscel[], const cs_real_2_t weighf[], cs_real_t *restrict i_massflux, cs_real_t *restrict b_massflux)
 Project the external source terms to the faces in coherence with cs_face_anisotropic_diffusion_scalar for the improved hydrostatic pressure algorithm (iphydr=1). More...
 
void cs_tensor_face_flux (const cs_mesh_t *m, cs_mesh_quantities_t *fvq, int f_id, int itypfl, int iflmb0, int init, int inc, int imrgra, int nswrgu, int imligu, int iwarnu, double epsrgu, double climgu, const cs_real_t c_rho[], const cs_real_t b_rho[], const cs_real_6_t c_var[], const cs_real_6_t coefav[], const cs_real_66_t coefbv[], cs_real_3_t *restrict i_massflux, cs_real_3_t *restrict b_massflux)
 Add $ \rho \tens{r} \vect{s}_\ij$ to a flux. More...
 

Function Documentation

◆ cs_divergence()

void cs_divergence ( const cs_mesh_t m,
int  init,
const cs_real_t  i_massflux[],
const cs_real_t  b_massflux[],
cs_real_t *restrict  diverg 
)

Add the integrated mass flux on the cells.

\[ \dot{m}_i = \dot{m}_i + \sum_{\fij \in \Facei{\celli}} \dot{m}_\ij \]

Parameters
[in]mpointer to mesh
[in]initindicator
  • 1 initialize the divergence to 0
  • 0 otherwise
[in]i_massfluxmass flux at interior faces
[in]b_massfluxmass flux at boundary faces
[in,out]divergmass flux divergence

◆ cs_ext_force_anisotropic_flux()

void cs_ext_force_anisotropic_flux ( const cs_mesh_t m,
cs_mesh_quantities_t fvq,
int  init,
int  nswrgp,
int  ircflp,
const cs_real_3_t  frcxt[],
const cs_real_t  cofbfp[],
const cs_real_t  i_visc[],
const cs_real_t  b_visc[],
cs_real_6_t  viscel[],
const cs_real_2_t  weighf[],
cs_real_t *restrict  i_massflux,
cs_real_t *restrict  b_massflux 
)

Project the external source terms to the faces in coherence with cs_face_anisotropic_diffusion_scalar for the improved hydrostatic pressure algorithm (iphydr=1).

Parameters
[in]mpointer to mesh
[in]fvqpointer to finite volume quantities
[in]initindicator
  • 1 initialize the mass flux to 0
  • 0 otherwise
[in]nswrgpnumber of reconstruction sweeps for the gradients
[in]ircflpindicator
  • 1 flux reconstruction,
  • 0 otherwise
[in]frcxtbody force creating the hydrostatic pressure
[in]cofbfpboundary condition array for the diffusion of the variable (implicit part)
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the r.h.s.
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the r.h.s.
[in]viscelsymmetric cell tensor $ \tens{\mu}_\celli $
[in]weighfinternal face weight between cells i j in case of tensor diffusion
[in,out]i_massfluxmass flux at interior faces
[in,out]b_massfluxmass flux at boundary faces

◆ cs_ext_force_flux()

void cs_ext_force_flux ( const cs_mesh_t m,
cs_mesh_quantities_t fvq,
int  init,
int  nswrgu,
const cs_real_3_t  frcxt[],
const cs_real_t  cofbfp[],
cs_real_t *restrict  i_massflux,
cs_real_t *restrict  b_massflux,
const cs_real_t  i_visc[],
const cs_real_t  b_visc[],
const cs_real_t  viselx[],
const cs_real_t  visely[],
const cs_real_t  viselz[] 
)

Project the external source terms to the faces in coherence with cs_face_diffusion_scalar for the improved hydrostatic pressure algorithm (iphydr=1).

Parameters
[in]mpointer to mesh
[in]fvqpointer to finite volume quantities
[in]initindicator
  • 1 initialize the mass flux to 0
  • 0 otherwise
[in]nswrgunumber of reconstruction sweeps for the gradients
[in]frcxtbody force creating the hydrostatic pressure
[in]cofbfpboundary condition array for the diffusion of the variable (implicit part)
[in,out]i_massfluxmass flux at interior faces
[in,out]b_massfluxmass flux at boundary faces
[in]i_visc$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the r.h.s.
[in]b_visc$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the r.h.s.
[in]viselxviscosity by cell, dir x
[in]viselyviscosity by cell, dir y
[in]viselzviscosity by cell, dir z

◆ cs_mass_flux()

void cs_mass_flux ( const cs_mesh_t m,
const cs_mesh_quantities_t fvq,
int  f_id,
int  itypfl,
int  iflmb0,
int  init,
int  inc,
int  imrgra,
int  nswrgu,
int  imligu,
int  iwarnu,
double  epsrgu,
double  climgu,
const cs_real_t  rom[],
const cs_real_t  romb[],
const cs_real_3_t  vel[],
const cs_real_3_t  coefav[],
const cs_real_33_t  coefbv[],
cs_real_t *restrict  i_massflux,
cs_real_t *restrict  b_massflux 
)

Add $ \rho \vect{u} \cdot \vect{s}_\ij$ to the mass flux $ \dot{m}_\ij $.

For the reconstruction, $ \gradt \left(\rho \vect{u} \right) $ is computed with the following approximated boundary conditions:

  • $ \vect{a}_{\rho u} = \rho_\fib \vect{a}_u $
  • $ \tens{b}_{\rho u} = \tens{b}_u $

For the mass flux at the boundary we have:

\[ \dot{m}_\ib = \left[ \rho_\fib \vect{a}_u + \rho_\fib \tens{b}_u \vect{u} + \tens{b}_u \left(\gradt \vect{u} \cdot \vect{\centi \centip}\right)\right] \cdot \vect{s}_\ij \]

The last equation uses some approximations detailed in the theory guide.

Parameters
[in]mpointer to mesh
[in]fvqpointer to finite volume quantities
[in]f_idfield id (or -1)
[in]itypflindicator (take rho into account or not)
  • 1 compute $ \rho\vect{u}\cdot\vect{s} $
  • 0 compute $ \vect{u}\cdot\vect{s} $
[in]iflmb0the mass flux is set to 0 on walls and symmetries if = 1
[in]initthe mass flux is initialized to 0 if > 0
[in]incindicator
  • 0 solve an increment
  • 1 otherwise
[in]imrgraindicator
  • 0 iterative gradient
  • 1 least square gradient
[in]nswrgunumber of sweeps for the reconstruction of the gradients
[in]imliguclipping gradient method
  • < 0 no clipping
  • = 0 thanks to neighbooring gradients
  • = 1 thanks to the mean gradient
[in]iwarnuverbosity
[in]epsrgurelative precision for the gradient reconstruction
[in]climguclipping coefficient for the computation of the gradient
[in]romcell density
[in]rombdensity at boundary faces
[in]velvector variable
[in]coefavboundary condition array for the variable (explicit part - vector array )
[in]coefbvboundary condition array for the variable (implicit part - 3x3 tensor array)
[in,out]i_massfluxmass flux at interior faces $ \dot{m}_\fij $
[in,out]b_massfluxmass flux at boundary faces $ \dot{m}_\fib $

For the reconstruction, $ \gradt \left(\rho \vect{u} \right) $ is computed with the following approximated boundary conditions:

  • $ \vect{a}_{\rho u} = \rho_\fib \vect{a}_u $
  • $ \tens{b}_{\rho u} = \tens{b}_u $

For the mass flux at the boundary we have:

\[ \dot{m}_\ib = \left[ \rho_\fib \vect{a}_u + \rho_\fib \tens{b}_u \vect{u} + \tens{b}_u \left(\gradt \vect{u} \cdot \vect{\centi \centip}\right)\right] \cdot \vect{s}_\ij \]

The last equation uses some approximations detailed in the theory guide.

Parameters
[in]mpointer to mesh
[in]fvqpointer to finite volume quantities
[in]f_idfield id (or -1)
[in]itypflindicator (take rho into account or not)
  • 1 compute $ \rho\vect{u}\cdot\vect{s} $
  • 0 compute $ \vect{u}\cdot\vect{s} $
[in]iflmb0the mass flux is set to 0 on symmetries if = 1
[in]initthe mass flux is initialized to 0 if > 0
[in]incindicator
  • 0 solve an increment
  • 1 otherwise
[in]imrgraindicator
  • 0 iterative gradient
  • 1 least square gradient
[in]nswrgunumber of sweeps for the reconstruction of the gradients
[in]imliguclipping gradient method
  • < 0 no clipping
  • = 0 thanks to neighbooring gradients
  • = 1 thanks to the mean gradient
[in]iwarnuverbosity
[in]epsrgurelative precision for the gradient reconstruction
[in]climguclipping coefficient for the computation of the gradient
[in]romcell density
[in]rombdensity at boundary faces
[in]velvector variable
[in]coefavboundary condition array for the variable (explicit part - vector array )
[in]coefbvboundary condition array for the variable (implicit part - 3x3 tensor array)
[in,out]i_massfluxmass flux at interior faces $ \dot{m}_\fij $
[in,out]b_massfluxmass flux at boundary faces $ \dot{m}_\fib $

◆ cs_tensor_divergence()

void cs_tensor_divergence ( const cs_mesh_t m,
int  init,
const cs_real_3_t  i_massflux[],
const cs_real_3_t  b_massflux[],
cs_real_3_t *restrict  diverg 
)

Add the integrated mass flux on the cells for a tensor variable.

\[ \dot{m}_i = \dot{m}_i + \sum_{\fij \in \Facei{\celli}} \dot{m}_\ij \]

Parameters
[in]mpointer to mesh
[in]initindicator
  • 1 initialize the divergence to 0
  • 0 otherwise
[in]i_massfluxmass flux vector at interior faces
[in]b_massfluxmass flux vector at boundary faces
[in,out]divergmass flux divergence vector

◆ cs_tensor_face_flux()

void cs_tensor_face_flux ( const cs_mesh_t m,
cs_mesh_quantities_t fvq,
int  f_id,
int  itypfl,
int  iflmb0,
int  init,
int  inc,
int  imrgra,
int  nswrgu,
int  imligu,
int  iwarnu,
double  epsrgu,
double  climgu,
const cs_real_t  c_rho[],
const cs_real_t  b_rho[],
const cs_real_6_t  c_var[],
const cs_real_6_t  coefav[],
const cs_real_66_t  coefbv[],
cs_real_3_t *restrict  i_massflux,
cs_real_3_t *restrict  b_massflux 
)

Add $ \rho \tens{r} \vect{s}_\ij$ to a flux.

Parameters
[in]mpointer to mesh
[in]fvqpointer to finite volume quantities
[in]f_idfield id (or -1)
[in]itypflindicator (take rho into account or not)
  • 1 compute $ \rho\vect{u}\cdot\vect{s} $
  • 0 compute $ \vect{u}\cdot\vect{s} $
[in]iflmb0the mass flux is set to 0 on walls and symmetries if = 1
[in]initthe mass flux is initialized to 0 if > 0
[in]incindicator
  • 0 solve an increment
  • 1 otherwise
[in]imrgraindicator
  • 0 iterative gradient
  • 1 least square gradient
[in]nswrgunumber of sweeps for the reconstruction of the gradients
[in]imliguclipping gradient method
  • < 0 no clipping
  • = 0 thanks to neighbooring gradients
  • = 1 thanks to the mean gradient
[in]iwarnuverbosity
[in]epsrgurelative precision for the gradient reconstruction
[in]climguclipping coefficient for the computation of the gradient
[in]c_rhocell density
[in]b_rhodensity at boundary faces
[in]c_varvariable
[in]coefavboundary condition array for the variable (explicit part - symmetric tensor array)
[in]coefbvboundary condition array for the variable (implicit part - 6x6 symmetric tensor array)
[in,out]i_massfluxmass flux at interior faces $ \dot{m}_\fij $
[in,out]b_massfluxmass flux at boundary faces $ \dot{m}_\fib $
[in]mpointer to mesh
[in]fvqpointer to finite volume quantities
[in]f_idfield id (or -1)
[in]itypflindicator (take rho into account or not)
  • 1 compute $ \rho\vect{u}\cdot\vect{s} $
  • 0 compute $ \vect{u}\cdot\vect{s} $
[in]iflmb0the mass flux is set to 0 on symmetries if = 1
[in]initthe mass flux is initialized to 0 if > 0
[in]incindicator
  • 0 solve an increment
  • 1 otherwise
[in]imrgraindicator
  • 0 iterative gradient
  • 1 least square gradient
[in]nswrgunumber of sweeps for the reconstruction of the gradients
[in]imliguclipping gradient method
  • < 0 no clipping
  • = 0 thanks to neighbooring gradients
  • = 1 thanks to the mean gradient
[in]iwarnuverbosity
[in]epsrgurelative precision for the gradient reconstruction
[in]climguclipping coefficient for the computation of the gradient
[in]c_rhocell density
[in]b_rhodensity at boundary faces
[in]c_varvariable
[in]coefavboundary condition array for the variable (explicit part - symmetric tensor array)
[in]coefbvboundary condition array for the variable (implicit part - 6x6 symmetric tensor array)
[in,out]i_massfluxmass flux at interior faces $ \dot{m}_\fij $
[in,out]b_massfluxmass flux at boundary faces $ \dot{m}_\fib $

◆ divmas()

void divmas ( const int *const  init,
const cs_real_t  i_massflux[],
const cs_real_t  b_massflux[],
cs_real_t  diverg[] 
)

◆ divmat()

void divmat ( const int *const  init,
const cs_real_3_t  i_massflux[],
const cs_real_3_t  b_massflux[],
cs_real_3_t  diverg[] 
)

◆ divrij()

void divrij ( const int *const  f_id,
const int *const  itypfl,
const int *const  iflmb0,
const int *const  init,
const int *const  inc,
const int *const  imrgra,
const int *const  nswrgu,
const int *const  imligu,
const int *const  iwarnu,
const cs_real_t *const  epsrgu,
const cs_real_t *const  climgu,
const cs_real_t  rom[],
const cs_real_t  romb[],
const cs_real_6_t  tensorvel[],
const cs_real_6_t  coefav[],
const cs_real_66_t  coefbv[],
cs_real_3_t  i_massflux[],
cs_real_3_t  b_massflux[] 
)

◆ inimav()

void inimav ( const int *const  f_id,
const int *const  itypfl,
const int *const  iflmb0,
const int *const  init,
const int *const  inc,
const int *const  imrgra,
const int *const  nswrgu,
const int *const  imligu,
const int *const  iwarnu,
const cs_real_t *const  epsrgu,
const cs_real_t *const  climgu,
const cs_real_t  rom[],
const cs_real_t  romb[],
const cs_real_3_t  vel[],
const cs_real_3_t  coefav[],
const cs_real_33_t  coefbv[],
cs_real_t  i_massflux[],
cs_real_t  b_massflux[] 
)

◆ projts()

void projts ( const int *const  init,
const int *const  nswrgu,
const cs_real_3_t  frcxt[],
const cs_real_t  cofbfp[],
cs_real_t  i_massflux[],
cs_real_t  b_massflux[],
const cs_real_t  i_visc[],
const cs_real_t  b_visc[],
const cs_real_t  viselx[],
const cs_real_t  visely[],
const cs_real_t  viselz[] 
)

◆ projtv()

void projtv ( const int *const  init,
const int *const  nswrgu,
const int *const  ircflp,
const cs_real_3_t  frcxt[],
const cs_real_t  cofbfp[],
const cs_real_t  i_visc[],
const cs_real_t  b_visc[],
cs_real_6_t  viscel[],
const cs_real_2_t  weighf[],
cs_real_t  i_massflux[],
cs_real_t  b_massflux[] 
)