nitrate.F90 6.03 KB
!$Id: nitrate.F90,v 1.11 2007-01-06 11:49:15 kbk Exp $
#include"cppdefs.h"
!-----------------------------------------------------------------------
!BOP
!
! !ROUTINE: The nitrate equation \label{sec:nitrate}
!
! !INTERFACE:
   subroutine nitrate(nlev,dt,cnpar,nus,gams,cc)
!
! !DESCRIPTION:
! This subroutine computes the balance of nitrate in the form
!  \begin{equation}
!   \label{SEq}
!    \dot{S}
!    = {\cal D}_S
!    - \frac{1}{\tau^S_R}(S-S_{obs})
!    \comma
!  \end{equation}
!  where $\dot{S}$ denotes the material derivative of the nitrate $N$, and
!  ${\cal D}_N$ is the sum of the turbulent and viscous transport
!  terms modelled according to
!  \begin{equation}
!   \label{DN}
!    {\cal D}_N
!    = \frstder{z}
!     \left(
!        \left( \nu^N_t + \nu^N \right) \partder{N}{z} - \tilde{\Gamma}_N
!        \right)
!    \point
!  \end{equation}
!  In this equation, $\nu^N_t$ and $\nu^N$ are the turbulent and
!  molecular diffusivities of nitrate, respectively,
!  and $\tilde{\Gamma}_N$
!  denotes the non-local flux of nitrate, see
!  \sect{sec:turbulenceIntro}. In the current version of GOTM,
!  we set $\nu^N_t = \nu^\N$ for simplicity.
!
!  Horizontal advection is optionally
!  included  (see {\tt obs.nml}) by means of prescribed
!  horizontal gradients $\partial_xN$ and $\partial_yN$ and
!  calculated horizontal mean velocities $U$ and $V$.
!  Relaxation with the time scale $\tau^N_R$
!  towards a precribed (changing in time)
!  profile $N_{obs}$ is possible.

!  Inner sources or sinks are not considered.
!  The surface freshwater flux is given by means of the precipitation
!  - evaporation data read in as $P-E$ through the {\tt airsea.nml} namelist:
!  \begin{equation}
!     \label{S_sbc}
!    {\cal D}_S =  S (P-E),
!    \qquad \mbox{at } z=\zeta,
!  \end{equation}
!  with $P-E$ given as a velocity (note that ${\cal D}_S$ is the flux in the
!  direction of $z$, and thus positive for a \emph{loss} of nitrate) .
!  Diffusion is numerically treated implicitly,
!  see equations (\ref{sigmafirst})-(\ref{sigmalast}).
!  The tri-diagonal matrix is solved then by a simplified Gauss elimination.
!  Vertical advection is included, and it must be non-conservative,
!  which is ensured by setting the local variable {\tt adv\_mode=0},
!  see section \ref{sec:advectionMean} on page \pageref{sec:advectionMean}.
!
! !USES:

   use meanflow,     only: avmoln                       !CHG3
   use meanflow,     only: h,u,v,w,nit,avh              !CHG3
   use observations, only: dndx,dndy,n_adv              !CHG3
   use observations, only: w_adv_discr,w_adv_method
   use observations, only: nprof,NRelaxTau              !CHG3
   use airsea,       only: p_e
   use util,         only: Dirichlet,Neumann
   use util,         only: oneSided,zeroDivergence
   use bio_var,      only: bio_model,numc

   IMPLICIT NONE
!

! !INPUT PARAMETERS:


!  number of vertical layers
   integer, intent(in)                 :: nlev

!   time step (s)
   REALTYPE, intent(in)                :: dt

!  numerical "implicitness" parameter
   REALTYPE, intent(in)                :: cnpar

!  diffusivity of nitrate (m^2/s) (!CHG3 same as salinity)
   REALTYPE, intent(in)                :: nus(0:nlev)

!  non-local salinity flux (psu m/s)
   REALTYPE, intent(in)                :: gams(0:nlev)

!  nitrate concentration after bio loop
!DD Not independent of the bio model, only works with 7 compartments model (Fasham)
!   REALTYPE, intent(in), optional      :: cc(1:10,0:nlev)   !CHG3
   REALTYPE, intent(in), optional      :: cc(1:numc,0:nlev)   !CHG3
!
! !REVISION HISTORY:
!  Original author(s): Dany Dumont (dany_dumont@ete.inrs.ca)
!
!  $Log: nitrate.F90,v $
!  Creation 1.0  2008-06-11 14:27:00  dd
!  based on salinity.F90
!
!EOP
!
! !LOCAL VARIABLES:
   integer                   :: adv_mode=0
   integer                   :: posconc=1
   integer                   :: i
   integer                   :: DiffBcup,DiffBcdw
   integer                   :: AdvBcup,AdvBcdw
   REALTYPE                  :: DiffSup,DiffSdw
   REALTYPE                  :: AdvSup,AdvSdw
   REALTYPE                  :: Lsour(0:nlev)
   REALTYPE                  :: Qsour(0:nlev)
!
!-----------------------------------------------------------------------
!BOC
!
!  set boundary conditions
   DiffBcup       = Neumann
   DiffBcdw       = Neumann
   DiffSup        = -1.*nit(nlev)*p_e     !CHG3
   DiffSdw        = _ZERO_

   AdvBcup       = zeroDivergence
   AdvBcdw       = oneSided
   AdvSup        = _ZERO_
   AdvSdw        = _ZERO_

!  compute total diffusivity
   do i=0,nlev
      avh(i)=nus(i)+avmoln            !CHG3
   end do

!  add contributions to source term
   Lsour=_ZERO_
   Qsour=_ZERO_

   do i=1,nlev
!     from non-local turbulence
#ifdef NONLOCAL
      Qsour(i) = Qsour(i) - ( gams(i) - gams(i-1) )/h(i)
#endif
   end do

!  ... and from lateral advection
   if (n_adv) then
      do i=1,nlev
         Qsour(i) = Qsour(i) - u(i)*dndx(i) - v(i)*dndy(i)        !CHG3
      end do
   end if

! redefinir nit apres un cyle bio
#ifdef BIO
!do i=1,nlev
!   !nit(i) = cc(7,i)    ! bio_model=2,4
!   nit(i) =  cc(1,i)   ! bio_model=1,6
!end do
   if (bio_model.eq.1) then
      do i=1,nlev
         nit(i) = cc(1,i)
      end do
   else if (bio_model.eq.2) then
      do i=1,nlev
         nit(i) = cc(7,i)
      end do
   else if (bio_model.eq.4) then
      do i=1,nlev
         nit(i) = cc(7,i)
      end do
  else if (bio_model.eq.6) then
      do i=1,nlev
         nit(i) = cc(1,i)
      end do
   end if
#endif

!  do advection step
   if (w_adv_method .ne. 0) then
      call adv_center(nlev,dt,h,h,w,AdvBcup,AdvBcdw,                    &
                          AdvSup,AdvSdw,w_adv_discr,adv_mode,nit) !CHG3
   end if

!  do diffusion step
   call diff_center(nlev,dt,cnpar,posconc,h,DiffBcup,DiffBcdw,          &
                    DiffSup,DiffSdw,avh,LSour,Qsour,NRelaxTau,nprof,nit)     !CHG3

   return
   end subroutine nitrate
!EOC

!-----------------------------------------------------------------------
! Copyright by the GOTM-team under the GNU Public License - www.gnu.org
!-----------------------------------------------------------------------