8.0
general documentation
cs_pressure_correction.c File Reference
#include "cs_defs.h"
#include <assert.h>
#include "bft_mem.h"
#include "bft_error.h"
#include "cs_ale.h"
#include "cs_array.h"
#include "cs_atmo.h"
#include "cs_balance.h"
#include "cs_blas.h"
#include "cs_boundary_conditions.h"
#include "cs_cf_thermo.h"
#include "cs_convection_diffusion.h"
#include "cs_divergence.h"
#include "cs_equation_iterative_solve.h"
#include "cs_face_viscosity.h"
#include "cs_field.h"
#include "cs_field_default.h"
#include "cs_field_operator.h"
#include "cs_field_pointer.h"
#include "cs_gradient.h"
#include "cs_halo.h"
#include "cs_lagr.h"
#include "cs_log.h"
#include "cs_matrix_building.h"
#include "cs_mesh_location.h"
#include "cs_parall.h"
#include "cs_parameters.h"
#include "cs_physical_constants.h"
#include "cs_physical_model.h"
#include "cs_porous_model.h"
#include "cs_prototypes.h"
#include "cs_sat_coupling.h"
#include "cs_sles_default.h"
#include "cs_time_step.h"
#include "cs_velocity_pressure.h"
#include "cs_volume_mass_injection.h"
#include "cs_vof.h"
#include "cs_post.h"
#include "cs_cdo_headers.h"
#include "cs_pressure_correction.h"
+ Include dependency graph for cs_pressure_correction.c:

Macros

#define CS_PRESSURE_CORRECTION_CDO_DBG   0
 

Functions

static void _pressure_correction_fv (int iterns, cs_lnum_t nfbpcd, cs_lnum_t ncmast, cs_lnum_t ifbpcd[nfbpcd], cs_lnum_t ltmast[], int isostd[], cs_real_t vel[restrict][3], cs_real_t da_uu[restrict][6], cs_real_t coefav[restrict][3], cs_real_t coefbv[restrict][3][3], cs_real_t coefa_dp[restrict], cs_real_t coefb_dp[restrict], cs_real_t spcond[restrict], cs_real_t svcond[restrict], cs_real_t frcxt[restrict][3], cs_real_t dfrcxt[restrict][3], cs_real_t i_visc[restrict], cs_real_t b_visc[restrict])
 Perform the pressure correction step of the Navier-Stokes equations for incompressible or slightly compressible flows. More...
 
static void _pressure_correction_cdo (cs_real_t vel[restrict][3], cs_real_t coefav[restrict][3], cs_real_t coefbv[restrict][3][3])
 Perform the pressure correction step of the Navier-Stokes equations for incompressible or slightly compressible flows. More...
 
static cs_pressure_correction_cdo_t_pressure_correction_cdo_create (void)
 Create and initialize a cs_pressure_correction_cdo_t structure. More...
 
void cs_pressure_correction_fv_activate (void)
 Activate the pressure increment solving with Legacy FV. More...
 
void cs_pressure_correction_cdo_activate (void)
 Activate the pressure increment solving with CDO. More...
 
void cs_pressure_correction_model_activate (void)
 Activate the pressure increment, either FV or CDO. More...
 
bool cs_pressure_correction_cdo_is_activated (void)
 Test if pressure solving with CDO is activated. More...
 
void cs_pressure_correction_cdo_init_setup (void)
 Start setting-up the pressure increment equation At this stage, numerical settings should be completely determined but connectivity and geometrical information is not yet available. More...
 
void cs_pressure_correction_cdo_finalize_setup (const cs_domain_t *domain)
 Finalize setting-up the pressure increment equation At this stage, numerical settings should be completely determined. More...
 
void cs_pressure_correction_cdo_destroy_all (void)
 Free the main structure related to the pressure correction. More...
 
void cs_pressure_correction (int iterns, cs_lnum_t nfbpcd, cs_lnum_t ncmast, cs_lnum_t ifbpcd[nfbpcd], cs_lnum_t ltmast[], int isostd[], cs_real_t vel[restrict][3], cs_real_t da_uu[restrict][6], cs_real_t coefav[restrict][3], cs_real_t coefbv[restrict][3][3], cs_real_t coefa_dp[restrict], cs_real_t coefb_dp[restrict], cs_real_t spcond[restrict], cs_real_t svcond[restrict], cs_real_t frcxt[restrict][3], cs_real_t dfrcxt[restrict][3], cs_real_t i_visc[restrict], cs_real_t b_visc[restrict])
 Perform the pressure correction step of the Navier-Stokes equations for incompressible or slightly compressible flows. More...
 

Variables

static const char _err_empty_prcdo []
 
static bool cs_pressure_correction_cdo_active = false
 
static cs_pressure_correction_cdo_tcs_pressure_correction_cdo = NULL
 

Detailed Description

Pressure correction step.

Macro Definition Documentation

◆ CS_PRESSURE_CORRECTION_CDO_DBG

#define CS_PRESSURE_CORRECTION_CDO_DBG   0

Function Documentation

◆ _pressure_correction_cdo()

static void _pressure_correction_cdo ( cs_real_t  vel[restrict][3],
cs_real_t  coefav[restrict][3],
cs_real_t  coefbv[restrict][3][3] 
)
static

Perform the pressure correction step of the Navier-Stokes equations for incompressible or slightly compressible flows.

\Remark:

  • CDO face-based scheme is used to solve the Poisson equation.

    \[ - DIV \cdot Hdg \cdot GRAD \left(\delta p \right) = - \divs \left( \rho \vect{\widetilde{u}}\right) \]

    The mass flux is then updated as follows:

    \[ \dot{m}^{n+1} = \dot{m}^{n} - Hodge \cdot GRAD \left(\delta p \right) \]

    Parameters
    [in]velvelocity
    [in]coefavboundary condition array for the variable (explicit part)
    [in]coefbvboundary condition array for the variable (implicit part)

◆ _pressure_correction_cdo_create()

static cs_pressure_correction_cdo_t* _pressure_correction_cdo_create ( void  )
static

Create and initialize a cs_pressure_correction_cdo_t structure.

Returns
a pointer to a newly allocated cs_pressure_correction_cdo_t structure

◆ _pressure_correction_fv()

static void _pressure_correction_fv ( int  iterns,
cs_lnum_t  nfbpcd,
cs_lnum_t  ncmast,
cs_lnum_t  ifbpcd[nfbpcd],
cs_lnum_t  ltmast[],
int  isostd[],
cs_real_t  vel[restrict][3],
cs_real_t  da_uu[restrict][6],
cs_real_t  coefav[restrict][3],
cs_real_t  coefbv[restrict][3][3],
cs_real_t  coefa_dp[restrict],
cs_real_t  coefb_dp[restrict],
cs_real_t  spcond[restrict],
cs_real_t  svcond[restrict],
cs_real_t  frcxt[restrict][3],
cs_real_t  dfrcxt[restrict][3],
cs_real_t  i_visc[restrict],
cs_real_t  b_visc[restrict] 
)
static

Perform the pressure correction step of the Navier-Stokes equations for incompressible or slightly compressible flows.

This function solves the following Poisson equation on the pressure:

\[ D \left( \Delta t, \delta p \right) = \divs \left( \rho \vect{\widetilde{u}}\right) - \Gamma^n + \dfrac{\rho^n - \rho^{n-1}}{\Delta t} \]

The mass flux is then updated as follows:

\[ \dot{m}^{n+1}_\ij = \dot{m}^{n}_\ij - \Delta t \grad_\fij \delta p \cdot \vect{S}_\ij \]

\Remark:

  • an iterative process is used to solve the Poisson equation.
  • if the arak coefficient is set to 1, the the Rhie & Chow filter is activated.

Please refer to the resopv section of the theory guide for more information.

Parameters
[in]iternsNavier-Stokes iteration number
[in]nfbpcdnumber of faces with condensation source term
[in]ncmastnumber of cells with condensation source terms
[in]ifbpcdindex of faces with condensation source term
[in]ltmastlist of cells with condensation source terms (1 to n numbering)
[in]isostdindicator of standard outlet and index of the reference outlet face
[in]velvelocity
[in,out]da_uuvelocity matrix
[in]coefavboundary condition array for the variable (explicit part)
[in]coefbvboundary condition array for the variable (implicit part)
[in]coefa_dpboundary conditions for the pressure increment
[in]coefb_dpboundary conditions for the pressure increment
[in]spcondvariable value associated to the condensation source term (for ivar=ipr, spcond is the flow rate $ \Gamma_{s,cond}^n $)
[in]svcondvariable value associated to the condensation source term (for ivar=ipr, svcond is the flow rate $ \Gamma_{v, cond}^n $)
[in]frcxtexternal forces making hydrostatic pressure
[in]dfrcxtvariation of the external forces composing the hydrostatic pressure
[in]i_viscvisc*surface/dist aux faces internes
[in]b_viscvisc*surface/dist aux faces de bord

◆ cs_pressure_correction()

void cs_pressure_correction ( int  iterns,
cs_lnum_t  nfbpcd,
cs_lnum_t  ncmast,
cs_lnum_t  ifbpcd[nfbpcd],
cs_lnum_t  ltmast[],
int  isostd[],
cs_real_t  vel[restrict][3],
cs_real_t  da_uu[restrict][6],
cs_real_t  coefav[restrict][3],
cs_real_t  coefbv[restrict][3][3],
cs_real_t  coefa_dp[restrict],
cs_real_t  coefb_dp[restrict],
cs_real_t  spcond[restrict],
cs_real_t  svcond[restrict],
cs_real_t  frcxt[restrict][3],
cs_real_t  dfrcxt[restrict][3],
cs_real_t  i_visc[restrict],
cs_real_t  b_visc[restrict] 
)

Perform the pressure correction step of the Navier-Stokes equations for incompressible or slightly compressible flows.

This function solves the following Poisson equation on the pressure:

\[ D \left( \Delta t, \delta p \right) = \divs \left( \rho \vect{\widetilde{u}}\right) - \Gamma^n + \dfrac{\rho^n - \rho^{n-1}}{\Delta t} \]

The mass flux is then updated as follows:

\[ \dot{m}^{n+1}_\ij = \dot{m}^{n}_\ij - \Delta t \grad_\fij \delta p \cdot \vect{S}_\ij \]

\Remark:

  • an iterative process is used to solve the Poisson equation.
  • if the arak coefficient is set to 1, the the Rhie & Chow filter is activated.

Please refer to the resopv section of the theory guide for more information.

Parameters
[in]iternsNavier-Stokes iteration number
[in]nfbpcdnumber of faces with condensation source term
[in]ncmastnumber of cells with condensation source terms
[in]ifbpcdindex of faces with condensation source term
[in]ltmastlist of cells with condensation source terms (1 to n numbering)
[in]isostdindicator of standard outlet and index of the reference outlet face
[in]velvelocity
[in,out]da_uuvelocity matrix
[in]coefavboundary condition array for the variable (explicit part)
[in]coefbvboundary condition array for the variable (implicit part)
[in]coefa_dpboundary conditions for the pressure increment
[in]coefb_dpboundary conditions for the pressure increment
[in]spcondvariable value associated to the condensation source term (for ivar=ipr, spcond is the flow rate $ \Gamma_{s,cond}^n $)
[in]svcondvariable value associated to the condensation source term (for ivar=ipr, svcond is the flow rate $ \Gamma_{v, cond}^n $)
[in]frcxtexternal forces making hydrostatic pressure
[in]dfrcxtvariation of the external forces composing the hydrostatic pressure
[in]i_viscvisc*surface/dist aux faces internes
[in]b_viscvisc*surface/dist aux faces de bord

◆ cs_pressure_correction_cdo_activate()

void cs_pressure_correction_cdo_activate ( void  )

Activate the pressure increment solving with CDO.

◆ cs_pressure_correction_cdo_destroy_all()

void cs_pressure_correction_cdo_destroy_all ( void  )

Free the main structure related to the pressure correction.

◆ cs_pressure_correction_cdo_finalize_setup()

void cs_pressure_correction_cdo_finalize_setup ( const cs_domain_t domain)

Finalize setting-up the pressure increment equation At this stage, numerical settings should be completely determined.

Parameters
[in]domainpointer to a cs_domaint_t structure

◆ cs_pressure_correction_cdo_init_setup()

void cs_pressure_correction_cdo_init_setup ( void  )

Start setting-up the pressure increment equation At this stage, numerical settings should be completely determined but connectivity and geometrical information is not yet available.

◆ cs_pressure_correction_cdo_is_activated()

bool cs_pressure_correction_cdo_is_activated ( void  )

Test if pressure solving with CDO is activated.

Returns
true if solving with CDO is requested, false otherwise

◆ cs_pressure_correction_fv_activate()

void cs_pressure_correction_fv_activate ( void  )

Activate the pressure increment solving with Legacy FV.

◆ cs_pressure_correction_model_activate()

void cs_pressure_correction_model_activate ( void  )

Activate the pressure increment, either FV or CDO.

Variable Documentation

◆ _err_empty_prcdo

const char _err_empty_prcdo[]
static
Initial value:
=
" Stop execution. The structure related to the pressure correction is"
" empty.\n Please check your settings.\n"

◆ cs_pressure_correction_cdo

cs_pressure_correction_cdo_t* cs_pressure_correction_cdo = NULL
static

◆ cs_pressure_correction_cdo_active

bool cs_pressure_correction_cdo_active = false
static