Blame view

src/extras/bio/ammonium.F90 5.8 KB
8473db4a   dumoda01   Les routines nitr...
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
!$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 ammonium(nlev,dt,cnpar,nus,gams,cc)
!
! !DESCRIPTION:
! This subroutine computes the balance of ammonium 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: avmola                       !CHG5
   use meanflow,     only: h,u,v,w,amm,avh              !CHG5
   use observations, only: dadx,dady,a_adv              !CHG5
   use observations, only: w_adv_discr,w_adv_method
   use observations, only: aprof,ARelaxTau              !CHG5
   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
   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.*amm(nlev)*p_e     !CHG5
   DiffSdw        = _ZERO_

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

!  compute total diffusivity
   do i=0,nlev
      avh(i)=nus(i)+avmola            !CHG5
   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 (a_adv) then
      do i=1,nlev
         Qsour(i) = Qsour(i) - u(i)*dadx(i) - v(i)*dady(i)        !CHG5
      end do
   end if

! redefinir amm apres un cyle bio
#ifdef BIO
!do i=1,nlev
!   !amm(i) = cc(6,i)     ! bio_model=2,4
!   amm(i) = cc(9,i)     ! bio_model=6
!end do
   if (bio_model.eq.2) then
      do i=1,nlev
         amm(i) = cc(6,i)
      end do
   else if (bio_model.eq.4) then
      do i=1,nlev
         amm(i) = cc(6,i)
      end do
   else if (bio_model.eq.6) then
      do i=1,nlev
         amm(i) = cc(9,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,amm)         !CHG5
   end if

!  do diffusion step
   call diff_center(nlev,dt,cnpar,posconc,h,DiffBcup,DiffBcdw,          &
                    DiffSup,DiffSdw,avh,LSour,Qsour,ARelaxTau,aprof,amm)   !CHG5

   return
   end subroutine ammonium
!EOC

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