Diffusion schemes -- grid centers (Source File: diff_center.F90)

INTERFACE:

    subroutine diff_center(N,dt,cnpar,posconc,h,Bcup,Bcdw, &
                           Yup,Ydw,nuY,Lsour,Qsour,Taur,Yobs,Y)

DESCRIPTION:

This subroutine solves the one-dimensional diffusion equation including source terms,

$$\dfrac{\partial {Y}}{\partial {t}} = \dfrac{\partial {}}{\partial... ... - \frac{1}{\tau_R}(Y-Y_{obs}) + Y L_{\text{sour}} + Q_{\text{sour}} \; , \quad$$ (231)

for al variables defined at the centers of the grid cells, and a diffusion coefficient $\nu_Y$ defined at the faces. Relaxation with time scale $\tau_R$ towards observed values $Y_{\text{obs}}$ is possible. $L_{\text{sour}}$ specifies a linear source term, and $Q_{\text{sour}}$ a constant source term. Central differences are used to discretize the problem as discussed in section 3.1.2. The diffusion term, the linear source term, and the linear part arising from the relaxation term are treated with an implicit method, whereas the constant source term is treated fully explicit.

The input parameters Bcup and Bcdw specify the type of the upper and lower boundary conditions, which can be either Dirichlet or Neumann-type. Bcup and Bcdw must have integer values corresponding to the parameters Dirichlet and Neumann defined in the module util, see section 8.2. Yup and Ydw are the values of the boundary conditions at the surface and the bottom. Depending on the values of Bcup and Bcdw, they represent either fluxes or prescribed values. The integer posconc indicates if a quantity is non-negative by definition (posconc=1, such as for concentrations) or not (posconc=0). For posconc=1 and negative boundary fluxes, the source term linearisation according to Patankar (1980) is applied.

Note that fluxes entering a boundary cell are counted positive by convention. The lower and upper position for prescribing these fluxes are located at the lowest und uppermost grid faces with index "0" and index "N", respectively. If values are prescribed, they are located at the centers with index "1" and index "N", respectively.

USES:

    use util,          only  : Dirichlet, Neumann
    use mtridiagonal
 
    IMPLICIT NONE

INPUT PARAMETERS:

 
    number of vertical layers
    integer,  intent(in)                :: N
 
    time step (s)
    REALTYPE, intent(in)                :: dt
 
    "implicitness" parameter
    REALTYPE, intent(in)                :: cnpar
 
    1: non-negative concentration, 0: else
    integer, intent(in)                 :: posconc
 
    layer thickness (m)
    REALTYPE, intent(in)                :: h(0:N)
 
    type of upper BC
    integer,  intent(in)                :: Bcup
 
    type of lower BC
    integer,  intent(in)                :: Bcdw
 
    value of upper BC
    REALTYPE, intent(in)                :: Yup
 
    value of lower BC
    REALTYPE, intent(in)                :: Ydw
 
    diffusivity of Y
    REALTYPE, intent(in)                :: nuY(0:N)
 
    linear source term
    (treated implicitly)
    REALTYPE, intent(in)                :: Lsour(0:N)
 
    constant source term
    (treated explicitly)
    REALTYPE, intent(in)                :: Qsour(0:N)
 
    relaxation time (s)
    REALTYPE, intent(in)                :: Taur(0:N)
 
    observed value of Y
    REALTYPE, intent(in)                :: Yobs(0:N)

INPUT/OUTPUT PARAMETERS:

    REALTYPE                            :: Y(0:N)

REVISION HISTORY:

    Original author(s): Lars Umlauf
    $Log: diff_center.F90,v $
    Revision 1.5  2010-09-17 12:53:52  jorn
    extensive code clean-up to ensure proper initialization and clean-up of all variables
    Revision 1.4  2005-11-17 09:58:20  hb
    explicit argument for positive definite variables in diff_center()
    Revision 1.3  2005/11/03 20:56:55  hb
    Source term linearisation now fully explicit again, reversion to old method
    Revision 1.2  2005/09/16 13:54:02  lars
    added missing IMPLICIT NONE
    Revision 1.1  2005/06/27 10:54:33  kbk
    new files needed