Blame view

src/extras/bio/hydrocarbon.F90 5.82 KB
84ed20a5   Dany Dumont   Construction d'un...
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
!$Id: nitrate.F90,v 1.11 2007-01-06 11:49:15 kbk Exp $
#include"cppdefs.h"
!-----------------------------------------------------------------------
!BOP
!
! !ROUTINE: The hydrocarbon equation \label{sec:hydrocarbon}
!
! !INTERFACE:
   subroutine hydrocarbon(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: avmolhc                       !CHG3
   use meanflow,     only: h,u,v,w,hcb,avh
!   use observations, only: dndx,dndy,n_adv              !CHG3
   use observations, only: w_adv_discr,w_adv_method
   use observations, only: hcprof,HCRelaxTau              !CHG3
   use airsea,       only: p_e
   use util,         only: Dirichlet,Neumann
   use util,         only: oneSided,zeroDivergence
   use bio_var,      only: bio_model,numc,ws

   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.*hcb(nlev)*p_e
   DiffSdw        = _ZERO_

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

!  compute total diffusivity
   do i=0,nlev
      avh(i)=nus(i)+avmolhc
   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
   if (bio_model.eq.7) then
      do i=1,nlev
         hcb(i) = cc(10,i)
      end do
   end if
#endif

!  do advection due to settling or rising
   call adv_center(nlev,dt,h,h,ws(10,:),AdvBcup,AdvBcdw,                    &
                          _ZERO_,_ZERO_,w_adv_discr,adv_mode,hcb)

!  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,hcb)
   end if

!  do diffusion step
   call diff_center(nlev,dt,cnpar,posconc,h,DiffBcup,DiffBcdw,          &
                    DiffSup,DiffSdw,avh,LSour,Qsour,HCRelaxTau,hcprof,hcb)

   return
   end subroutine hydrocarbon
!EOC

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