programmer's documentation
Functions/Subroutines
bilscts.f90 File Reference

Wrapper to the function which adds the explicit part of the convection/diffusion terms of a transport equation of a tensor field $ \tens{\varia} $. More...

Functions/Subroutines

subroutine bilscts (idtvar, ivar, iconvp, idiffp, nswrgp, imligp, ircflp, ischcp, isstpp, inc, imrgra, iwarnp, idftnp, blencp, epsrgp, climgp, relaxp, thetap, pvar, pvara, coefa, coefb, cofaf, cofbf, flumas, flumab, viscf, viscb, viscce, weighf, weighb, icvflb, icvfli, smbr)
 

Detailed Description

Wrapper to the function which adds the explicit part of the convection/diffusion terms of a transport equation of a tensor field $ \tens{\varia} $.

More precisely, the right hand side $ \vect{Rhs} $ is updated as follows:

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

Warning:

Options for the diffusive scheme:

Options for the convective scheme:

Function/Subroutine Documentation

subroutine bilscts ( integer  idtvar,
integer  ivar,
integer  iconvp,
integer  idiffp,
integer  nswrgp,
integer  imligp,
integer  ircflp,
integer  ischcp,
integer  isstpp,
integer  inc,
integer  imrgra,
integer  iwarnp,
integer  idftnp,
double precision  blencp,
double precision  epsrgp,
double precision  climgp,
double precision  relaxp,
double precision  thetap,
double precision, dimension (6 ,ncelet)  pvar,
double precision, dimension (6 ,ncelet)  pvara,
double precision, dimension(6 ,nfabor)  coefa,
double precision, dimension(6,6,nfabor)  coefb,
double precision, dimension(6 ,nfabor)  cofaf,
double precision, dimension(6,6,nfabor)  cofbf,
double precision, dimension(nfac)  flumas,
double precision, dimension(nfabor)  flumab,
double precision, dimension (*)  viscf,
double precision, dimension (nfabor)  viscb,
double precision, dimension(6,ncelet)  viscce,
double precision, dimension(2,nfac)  weighf,
double precision, dimension(nfabor)  weighb,
integer  icvflb,
integer, dimension(nfabor)  icvfli,
double precision, dimension(6,ncelet)  smbr 
)
Parameters
[in]idtvarindicator of the temporal scheme
[in]ivarindex of the current variable
[in]iconvpindicator
  • 1 convection,
  • 0 otherwise
[in]idiffpindicator
  • 1 diffusion,
  • 0 otherwise
[in]nswrgpnumber of reconstruction sweeps for the gradients
[in]imligpclipping gradient method
  • < 0 no clipping
  • = 0 by neighboring gradients
  • = 1 by the mean gradient
[in]ircflpindicator
  • 1 flux reconstruction,
  • 0 otherwise
[in]ischcpindicator
  • 1 centered
  • 0 2nd order
[in]isstppindicator
  • 1 without slope test
  • 0 with slope test
[in]incindicator
  • 0 when solving an increment
  • 1 otherwise
[in]imrgraindicator
  • 0 iterative gradient
  • 1 least squares gradient
[in]iwarnpverbosity
[in]idftnpindicator
  • 1 scalar diffusivity
  • 6 symmetric tensor diffusivity
[in]blencpfraction of upwinding
[in]epsrgprelative precision for the gradient reconstruction
[in]climgpclipping coefficient for the computation of the gradient
[in]relaxpcoefficient of relaxation
[in]thetapweighting coefficient for the theta-schema,
  • thetap = 0: explicit scheme
  • thetap = 0.5: time-centered scheme (mix between Crank-Nicolson and Adams-Bashforth)
  • thetap = 1: implicit scheme
[in]pvarsolved velocity (current time step)
[in]pvarasolved velocity (previous time step)
[in]coefaboundary condition array for the variable (Explicit part)
[in]coefbboundary condition array for the variable (Impplicit part)
[in]cofafboundary condition array for the diffusion of the variable (Explicit part)
[in]cofbfboundary condition array for the diffusion of the variable (Implicit part)
[in]flumasmass flux at interior faces
[in]flumabmass flux at boundary faces
[in]viscf$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the r.h.s.
[in]viscb$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at boundary faces for the r.h.s.
[in]visccesymmetric cell tensor $ \tens{\mu}_\celli $
[in]weighfinternal face weight between cells i j in case of tensor diffusion
[in]weighbboundary face weight for cells i in case of tensor diffusion
[in]icvflbglobal indicator of boundary convection flux
  • 0 upwind scheme at all boundary faces
  • 1 imposed flux at some boundary faces
[in]icvfliboundary face indicator array of convection flux
  • 0 upwind scheme
  • 1 imposed flux
[in,out]smbrright hand side $ \vect{Rhs} $