Code_Saturne
CFD tool
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Functions/Subroutines
coditv.f90 File Reference

This function solves an advection diffusion equation with source terms for one time step for the vector variable $ \vect{a} $. More...

Functions/Subroutines

subroutine coditv
 

Detailed Description

This function solves an advection diffusion equation with source terms for one time step for the vector variable $ \vect{a} $.

The equation reads:

\[ \tens{f_s}^{imp}(\vect{a}^{n+1}-\vect{a}^n) + \divv \left( \vect{a}^{n+1} \otimes \rho \vect {u} - \mu \gradt \vect{a}^{n+1}\right) = \vect{Rhs} \]

This equation is rewritten as:

\[ \tens{f_s}^{imp} \delta \vect{a} + \divv \left( \delta \vect{a} \otimes \rho \vect{u} - \mu \gradt \delta \vect{a} \right) = \vect{Rhs}^1 \]

where $ \delta \vect{a} = \vect{a}^{n+1} - \vect{a}^n$ and $ \vect{Rhs}^1 = \vect{Rhs} - \divv \left( \vect{a}^n \otimes \rho \vect{u} - \mu \gradt \vect{a}^n \right)$

It is in fact solved with the following iterative process:

\[ \tens{f_s}^{imp} \delta \vect{a}^k + \divv \left( \delta \vect{a}^k \otimes \rho \vect{u} - \mu \gradt \delta \vect{a}^k \right) = \vect{Rhs}^k \]

where $ \vect{Rhs}^k = \vect{Rhs} - \tens{f_s}^{imp} \left(\vect{a}^k-\vect{a}^n \right) - \divv \left( \vect{a}^k \otimes \rho \vect{u} - \mu \gradt \vect{a}^k \right)$

Be careful, it is forbidden to modify $ \tens{f_s}^{imp} $ here!

Function/Subroutine Documentation

subroutine coditv ( )
Parameters
[in]nvartotal number of variables
[in]nscaltotal number of scalars
[in]idtvarindicator of the temporal scheme
[in]ivarindex of the current variable
[in]iconvpindicator
  • 1 convection,
  • 0 sinon
[in]idiffpindicator
  • 1 diffusion,
  • 0 sinon
[in]ireslpindicator
  • 0 conjugate gradient
  • 1 jacobi
  • 2 bi-cgstab
[in]ndircpindicator (0 if the diagonal is stepped aside)
[in]nitmapmaximum number of iteration to solve the iterative process
[in]imrgraindicateur
  • 0 iterative gradient
  • 1 least square gradient
[in]nswrspnumber of reconstruction sweeps for the Right Hand Side
[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]ivisepindicator to take $ \divv \left(\mu \gradt \transpose{\vect{a}} \right) -2/3 \grad\left( \mu \dive \vect{a} \right)$
  • 1 take into account,
  • 0 otherwise
[in]ischcpindicator
  • 1 centred
  • 0 2nd order
[in]isstppindicator
  • 1 without slope test
  • 0 with slope test
[in]iescapcompute the predictor indicator if 1
[in]imgrpindicator
  • 0 no multi-grid
  • 1 otherwise
[in]ipp*index of the variable for post-processing
[in]iwarnpverbosity
[in]blencpfraction of upwinding
[in]epsilpprecision pour resol iter
[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]pvaravariable at the previous time step $ \vect{a}^n $
[in]pvarkvariable at the previous sub-iteration $ \vect{a}^k $. If you sub-iter on Navier-Stokes, then it allows to initialize by something else than pvara (usually pvar=pvara)
[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]flumasmass flux at interior faces
[in]flumabmass flux at boundary faces
[in]viscfm$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the matrix
[in]viscbm$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the matrix
[in]viscfs$ \mu_\fij \dfrac{S_\fij}{\ipf \jpf} $ at interior faces for the r.h.s.
[in]viscbs$ \mu_\fib \dfrac{S_\fib}{\ipf \centf} $ at border faces for the r.h.s.
[in]fimp$ \tens{f_s}^{imp} $
[in]smbrpRight hand side $ \vect{Rhs}^k $
[in,out]pvarcurrent variable
[out]esworkprediction-stage error estimator (if iescap > 0)

Here is the call graph for this function: