programmer's documentation
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
Macros | Functions
cs_convection_diffusion_balance.c File Reference
#include "cs_defs.h"
#include <assert.h>
#include <errno.h>
#include <stdio.h>
#include <stdarg.h>
#include <string.h>
#include <math.h>
#include <float.h>
#include "bft_error.h"
#include "bft_mem.h"
#include "bft_printf.h"
#include "cs_blas.h"
#include "cs_halo.h"
#include "cs_halo_perio.h"
#include "cs_log.h"
#include "cs_mesh.h"
#include "cs_field.h"
#include "cs_gradient.h"
#include "cs_gradient_perio.h"
#include "cs_ext_neighborhood.h"
#include "cs_mesh_quantities.h"
#include "cs_parameters.h"
#include "cs_prototypes.h"
#include "cs_timer.h"
#include "cs_convection_diffusion_balance.h"
Include dependency graph for cs_convection_diffusion_balance.c:

Macros

#define THR_MIN   128
 

Functions

void bilsc2 (const cs_int_t *const idtvar, const cs_int_t *const f_id, const cs_int_t *const iconvp, const cs_int_t *const idiffp, const cs_int_t *const nswrgp, const cs_int_t *const imligp, const cs_int_t *const ircflp, const cs_int_t *const ischcp, const cs_int_t *const isstpp, const cs_int_t *const icvflb, const cs_int_t *const inc, const cs_int_t *const imrgra, const cs_int_t *const iccocg, const cs_int_t *const ifaccp, const cs_int_t *const iwarnp, const cs_real_t *const blencp, const cs_real_t *const epsrgp, const cs_real_t *const climgp, const cs_real_t *const extrap, const cs_real_t *const relaxp, const cs_real_t *const thetap, cs_real_t pvar[], const cs_real_t pvara[], const cs_int_t bc_type[], const cs_int_t icvfli[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_t rhs[])
 
void cs_convection_diffusion_scalar (int idtvar, int f_id, int iconvp, int idiffp, int nswrgp, int imligp, int ircflp, int ischcp, int isstpp, int icvflb, int inc, int imrgra, int iccocg, int ifaccp, int iwarnp, double blencp, double epsrgp, double climgp, double extrap, double relaxp, double thetap, cs_real_t *restrict pvar, const cs_real_t *restrict pvara, const cs_int_t bc_type[], const cs_int_t icvfli[], const cs_real_t coefap[], const cs_real_t coefbp[], const cs_real_t cofafp[], const cs_real_t cofbfp[], const cs_real_t i_massflux[], const cs_real_t b_massflux[], const cs_real_t i_visc[], const cs_real_t b_visc[], cs_real_t *restrict rhs)
 This function adds the explicit part of the convection/diffusion terms of a standard transport equation of a scalar field $ \varia $. More...
 

Macro Definition Documentation

#define THR_MIN   128

Function Documentation

void bilsc2 ( const cs_int_t *const  idtvar,
const cs_int_t *const  f_id,
const cs_int_t *const  iconvp,
const cs_int_t *const  idiffp,
const cs_int_t *const  nswrgp,
const cs_int_t *const  imligp,
const cs_int_t *const  ircflp,
const cs_int_t *const  ischcp,
const cs_int_t *const  isstpp,
const cs_int_t *const  icvflb,
const cs_int_t *const  inc,
const cs_int_t *const  imrgra,
const cs_int_t *const  iccocg,
const cs_int_t *const  ifaccp,
const cs_int_t *const  iwarnp,
const cs_real_t *const  blencp,
const cs_real_t *const  epsrgp,
const cs_real_t *const  climgp,
const cs_real_t *const  extrap,
const cs_real_t *const  relaxp,
const cs_real_t *const  thetap,
cs_real_t  pvar[],
const cs_real_t  pvara[],
const cs_int_t  bc_type[],
const cs_int_t  icvfli[],
const cs_real_t  coefap[],
const cs_real_t  coefbp[],
const cs_real_t  cofafp[],
const cs_real_t  cofbfp[],
const cs_real_t  i_massflux[],
const cs_real_t  b_massflux[],
const cs_real_t  i_visc[],
const cs_real_t  b_visc[],
cs_real_t  rhs[] 
)

(end ignore by Doxygen)

void cs_convection_diffusion_scalar ( int  idtvar,
int  f_id,
int  iconvp,
int  idiffp,
int  nswrgp,
int  imligp,
int  ircflp,
int  ischcp,
int  isstpp,
int  icvflb,
int  inc,
int  imrgra,
int  iccocg,
int  ifaccp,
int  iwarnp,
double  blencp,
double  epsrgp,
double  climgp,
double  extrap,
double  relaxp,
double  thetap,
cs_real_t *restrict  pvar,
const cs_real_t *restrict  pvara,
const cs_int_t  bc_type[],
const cs_int_t  icvfli[],
const cs_real_t  coefap[],
const cs_real_t  coefbp[],
const cs_real_t  cofafp[],
const cs_real_t  cofbfp[],
const cs_real_t  i_massflux[],
const cs_real_t  b_massflux[],
const cs_real_t  i_visc[],
const cs_real_t  b_visc[],
cs_real_t *restrict  rhs 
)

This function adds the explicit part of the convection/diffusion terms of a standard transport equation of a scalar field $ \varia $.

More precisely, the right hand side $ Rhs $ is updated as follows:

\[ Rhs = Rhs - \sum_{\fij \in \Facei{\celli}} \left( \dot{m}_\ij \left( \varia_\fij - \varia_\celli \right) - \mu_\fij \gradv_\fij \varia \cdot \vect{S}_\ij \right) \]

Warning:

  • $ Rhs $ has already been initialized before calling bilsc2!
  • mind the sign minus

Options:

  • blencp = 0: upwind scheme for the advection
  • blencp = 1: no upwind scheme except in the slope test
  • ischcp = 0: second order
  • ischcp = 1: centred
Parameters
[in]idtvarindicator of the temporal scheme
[in]f_idfield id (or -1)
[in]iconvpindicator
  • 1 convection,
  • 0 sinon
[in]idiffpindicator
  • 1 diffusion,
  • 0 sinon
[in]nswrgpnumber of reconstruction sweeps for the gradients
[in]imligpclipping gradient method
  • < 0 no clipping
  • = 0 thank to neighbooring gradients
  • = 1 thank to the mean gradient
[in]ircflpindicator
  • 1 flux reconstruction,
  • 0 otherwise
[in]ischcpindicator
  • 1 centred
  • 0 2nd order
[in]isstppindicator
  • 1 without slope test
  • 0 with slope test
[in]icvflbglobal indicator of boundary convection flux
  • 0 upwind scheme at all boundary faces
  • 1 imposed flux at some boundary faces
[in]incindicator
  • 0 when solving an increment
  • 1 otherwise
[in]imrgraindicator
  • 0 iterative gradient
  • 1 least square gradient
[in]iccocgindicator
  • 1 re-compute cocg matrix (for iterativ gradients)
  • 0 otherwise
[in]ifaccpindicator
  • 1 coupling activated
  • 0 coupling not activated
[in]iwarnpverbosity
[in]blencpfraction of upwinding
[in]epsrgprelative precision for the gradient reconstruction
[in]climgpclipping coeffecient for the computation of the gradient
[in]extrapcoefficient for extrapolation of the gradient
[in]relaxpcoefficient of relaxation
[in]thetapweightening coefficient for the theta-schema,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centred scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]pvarsolved variable (current time step)
[in]pvarasolved variable (previous time step)
[in]bc_typeboundary condition type
[in]icvfliboundary face indicator array of convection flux
  • 0 upwind scheme
  • 1 imposed flux
[in]coefapboundary condition array for the variable (Explicit part)
[in]coefbpboundary condition array for the variable (Impplicit part)
[in]cofafpboundary condition array for the diffusion of the variable (Explicit part)
[in]cofbfpboundary condition array for the diffusion of the variable (Implicit part)
[in]i_massfluxmass flux at interior faces
[in]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,out]smbrpright hand side $ \vect{Rhs} $