programmer's documentation
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
Functions/Subroutines
codits.f90 File Reference

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

Functions/Subroutines

subroutine codits (idtvar, ivar, iconvp, idiffp, ireslp, ndircp, nitmap, imrgra, nswrsp, nswrgp, imligp, ircflp, ischcp, isstpp, iescap, imucpp, idftnp, iswdyp, imgrp, ncymxp, nitmfp, ipp, iwarnp, blencp, epsilp, epsrsp, epsrgp, climgp, extrap, relaxp, thetap, pvara, pvark, coefap, coefbp, cofafp, cofbfp, flumas, flumab, viscfm, viscbm, visccm, viscfs, viscbs, visccs, weighf, weighb, icvflb, icvfli, rovsdt, smbrp, pvar, dpvar, xcpp, eswork)
 

Detailed Description

This function solves an advection diffusion equation with source terms for one time step for the variable $ a $.

The equation reads:

\[ f_s^{imp}(a^{n+1}-a^n) + \divs \left( a^{n+1} \rho \vect{u} - \mu \grad a^{n+1} \right) = Rhs \]

This equation is rewritten as:

\[ f_s^{imp} \delta a + \divs \left( \delta a \rho \vect{u} - \mu \grad \delta a \right) = Rhs^1 \]

where $ \delta a = a^{n+1} - a^n$ and $ Rhs^1 = Rhs - \divs( a^n \rho \vect{u} - \mu \grad a^n)$

It is in fact solved with the following iterative process:

\[ f_s^{imp} \delta a^k + \divs \left(\delta a^k \rho \vect{u}-\mu\grad\delta a^k \right) = Rhs^k \]

where $Rhs^k=Rhs-f_s^{imp}(a^k-a^n) - \divs \left( a^k\rho\vect{u}-\mu\grad a^k \right)$

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

Function/Subroutine Documentation

subroutine codits ( integer  idtvar,
integer  ivar,
integer  iconvp,
integer  idiffp,
integer  ireslp,
integer  ndircp,
integer  nitmap,
integer  imrgra,
integer  nswrsp,
integer  nswrgp,
integer  imligp,
integer  ircflp,
integer  ischcp,
integer  isstpp,
integer  iescap,
integer  imucpp,
integer  idftnp,
integer  iswdyp,
integer  imgrp,
integer  ncymxp,
integer  nitmfp,
integer  ipp,
integer  iwarnp,
double precision  blencp,
double precision  epsilp,
double precision  epsrsp,
double precision  epsrgp,
double precision  climgp,
double precision  extrap,
double precision  relaxp,
double precision  thetap,
double precision, dimension(ncelet)  pvara,
double precision, dimension(ncelet)  pvark,
double precision, dimension(nfabor)  coefap,
double precision, dimension(nfabor)  coefbp,
double precision, dimension(nfabor)  cofafp,
double precision, dimension(nfabor)  cofbfp,
double precision, dimension(nfac)  flumas,
double precision, dimension(nfabor)  flumab,
double precision, dimension(nfac)  viscfm,
double precision, dimension(nfabor)  viscbm,
double precision, dimension(ncelet)  visccm,
double precision, dimension(nfac)  viscfs,
double precision, dimension(nfabor)  viscbs,
double precision, dimension(ncelet)  visccs,
double precision, dimension(2,nfac)  weighf,
double precision, dimension(nfabor)  weighb,
integer  icvflb,
integer, dimension(nfabor)  icvfli,
double precision, dimension(ncelet)  rovsdt,
double precision, dimension(ncelet)  smbrp,
double precision, dimension(ncelet)  pvar,
double precision, dimension(ncelet)  dpvar,
double precision, dimension(ncelet)  xcpp,
double precision, dimension(ncelet)  eswork 
)
Parameters
[in]idtvarindicateur du schema temporel
[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]imrgraindicator
  • 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]ischcpindicator
  • 1 centred
  • 0 2nd order
[in]isstppindicator
  • 1 without slope test
  • 0 with slope test
[in]iescapcompute the predictor indicator if 1
[in]imucppindicator
  • 0 do not multiply the convectiv term by Cp
  • 1 do multiply the convectiv term by Cp
[in]idftnpindicator
  • 0 the diffusivity is scalar
  • 1 the diffusivity is a diagonal tensor
  • 2 the diffusivity is a symmetric tensor
[in]iswdypindicator
  • 0 no dynamic relaxation
  • 1 dynamic relaxation depending on $ \delta \varia^k $
  • 2 dynamic relaxation depending on $ \delta \varia^k $ and $ \delta \varia^{k-1} $
[in]imgrpindicator
  • 0 no multi-grid
  • 1 otherwise
[in]ncymxpmax. number of multigrid cycles
[in]nitmfpnumber of equivalent iterations on fine mesh
[in]ippindex of the variable for post-processing
[in]iwarnpverbosity
[in]blencpfraction of upwinding
[in]epsilpprecision pour resol iter
[in]epsrsprelative precision for the iterative process
[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 $ a^n $
[in]pvarkvariable at the previous sub-iteration $ 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]visccmsymmetric cell tensor $ \tens{\mu}_\celli $
[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]visccssymmetric 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]rovsdt$ f_s^{imp} $
[in]smbrpRight hand side $ Rhs^k $
[in,out]pvarcurrent variable
[in,out]dpvarlast variable increment
[in]xcpparray of specific heat (Cp)
[out]esworkprediction-stage error estimator (if iescap > 0)