Commit 74bf91b3e5f57a36b4081ff6adcf8d2a3f4d058a

Authored by Dany Dumont
1 parent 70a3afeb
Exists in master and in 1 other branch snow

Routine bio_npzd4.F90

Showing 1 changed file with 688 additions and 0 deletions   Show diff stats
src/extras/bio/bio_npzd4.F90 0 → 100644
... ... @@ -0,0 +1,688 @@
  1 +!$id: bio_npzd4.F90,v 1.11 2016-05-05 11:49:15 dd Exp $
  2 +#include"cppdefs.h"
  3 +!-----------------------------------------------------------------------
  4 +!BOP
  5 +!
  6 +! !MODULE: bio_fasham --- Fasham et al. biological model \label{sec:bio-fasham}
  7 +!
  8 +! !INTERFACE:
  9 + module bio_npzd4
  10 +!
  11 +! !DESCRIPTION:
  12 +! The model developed by \cite{Fashametal1990}
  13 +! uses nitrogen as 'currency' according to the evidence that in
  14 +! most cases nitrogen is the limiting macronutrient. It consists of
  15 +! seven state variables: phytoplankton, zooplankton, bacteria,
  16 +! particulate organic matter (detritus), dissolved organic matter
  17 +! and the nutrients nitrate and ammonium.
  18 +! The structure of the \cite{Fashametal1990} biogeochemical model
  19 +! is given in figure \ref{fig_fasham}.
  20 +! \begin{figure}
  21 +! \begin{center}
  22 +! \scalebox{0.5}{\includegraphics{figures/fasham_structure.eps}}
  23 +! \caption{Structure of the \cite{Fashametal1990} model with bacteria (bac),
  24 +! phytoplankton (phy), detritus (det), zooplankton (zoo), labile dissolved
  25 +! organic nitrogen (don), ammonium (amm) and nitrate (nit) as the seven
  26 +! state variables.
  27 +! The concentrations are in mmol N\,m$^{-3}$,
  28 +! all fluxes (green arrows) are conservative.
  29 +! }\label{fig_fasham}
  30 +! \end{center}
  31 +! \end{figure}
  32 +! A detailed mathematical description of all
  33 +! processes is given in section \ref{sec:bio-fasham-rhs}.
  34 +! The version of the \cite{Fashametal1990} model which is implemented includes
  35 +! slight modifications by \cite{KuehnRadach1997} and has been
  36 +! included into GOTM by \cite{Burchardetal05}.
  37 +
  38 +! !USES:
  39 +! default: all is private.
  40 + use bio_var
  41 + use output
  42 + use observations, only : aa,g2
  43 + private
  44 +!
  45 +! !PUBLIC MEMBER FUNCTIONS:
  46 + public init_bio_npzd4, init_var_npzd4, var_info_npzd4, &
  47 + light_npzd4, do_bio_npzd4, end_bio_npzd4
  48 +!
  49 +! !PRIVATE DATA MEMBERS:
  50 +!
  51 +! !REVISION HISTORY:!
  52 +! Original author(s): Hans Burchard & Karsten Bolding
  53 +!
  54 +! $Log: bio_fasham.F90,v $
  55 +! Revision 1.11 2007-01-06 11:49:15 kbk
  56 +! namelist file extension changed .inp --> .nml
  57 +!
  58 +! Revision 1.10 2006-10-26 13:12:46 kbk
  59 +! updated bio models to new ode_solver
  60 +!
  61 +! Revision 1.9 2005-12-02 20:57:27 hb
  62 +! Documentation updated and some bugs fixed
  63 +!
  64 +! Revision 1.8 2005-11-17 09:58:18 hb
  65 +! explicit argument for positive definite variables in diff_center()
  66 +!
  67 +! Revision 1.7 2005/09/12 14:48:33 kbk
  68 +! merged generic biological module support
  69 +!
  70 +! Revision 1.6.2.1 2005/07/05 20:25:35 hb
  71 +! added control over par calculation
  72 +!
  73 +! Revision 1.6 2004/08/09 11:53:39 hb
  74 +! bioshading now without detritus
  75 +!
  76 +! Revision 1.5 2004/08/02 08:34:36 hb
  77 +! updated init routines to reflect new internal bio interface
  78 +!
  79 +! Revision 1.4 2004/08/01 15:52:57 hb
  80 +! alpha now devided by seconds per day
  81 +!
  82 +! Revision 1.3 2004/07/30 09:22:20 hb
  83 +! use bio_var in specific bio models - simpliefied internal interface
  84 +!
  85 +! Revision 1.2 2004/07/28 11:34:29 hb
  86 +! Bioshade feedback may now be switched on or off, depending on bioshade_feedback set to .true. or .false. in bio.nml
  87 +!
  88 +! Revision 1.1 2004/06/29 08:03:16 hb
  89 +! Fasham et al. 1990 model implemented
  90 +!
  91 +! Revision 1.2 2003/10/16 15:42:16 kbk
  92 +! simple mussesl model implemented - filter only
  93 +!
  94 +! Revision 1.1 2003/07/23 12:27:31 hb
  95 +! more generic support for different bio models
  96 +!
  97 +! Revision 1.3 2003/04/05 07:01:41 kbk
  98 +! moved bioshade variable to meanflow - to compile properly
  99 +!
  100 +! Revision 1.2 2003/04/04 14:25:52 hb
  101 +! First iteration of four-compartment geobiochemical model implemented
  102 +!
  103 +! Revision 1.1 2003/04/01 17:01:00 hb
  104 +! Added infrastructure for geobiochemical model
  105 +!
  106 +! !LOCAL VARIABLES:
  107 +! from a namelist
  108 + REALTYPE :: p_initial= 0.056666666
  109 + REALTYPE :: z_initial= 0.05
  110 + REALTYPE :: b_initial= 0.001
  111 + REALTYPE :: d_initial= 0.416666666
  112 +! Nitrate and ammonium are initialized within the GOTM observation module
  113 +! REALTYPE :: n_initial= 8.3
  114 +! REALTYPE :: a_initial= 0.22
  115 + REALTYPE :: l_initial= 0.14
  116 + REALTYPE :: p0 = 0.0
  117 + REALTYPE :: z0 = 0.0
  118 + REALTYPE :: b0 = 0.0
  119 + REALTYPE :: vp = 1.5
  120 + REALTYPE :: alpha = 0.065
  121 + REALTYPE :: I_opt = 10.0
  122 + REALTYPE :: inib = 0.05
  123 + REALTYPE :: theta = 0.0
  124 + REALTYPE :: w_pmin = -0.06
  125 + REALTYPE :: w_pmax = -0.38
  126 + REALTYPE :: w_zmax = 100.0
  127 + REALTYPE :: k1 = 0.2
  128 + REALTYPE :: k2 = 0.8
  129 + REALTYPE :: mu1 = 0.05
  130 + REALTYPE :: k5 = 0.2
  131 + REALTYPE :: gamma = 0.05
  132 + REALTYPE :: w_p = -0.5
  133 + REALTYPE :: gmax = 1.0
  134 + REALTYPE :: k3 = 1.0
  135 + REALTYPE :: beta = 0.625
  136 + REALTYPE :: mu2 = 0.3
  137 + REALTYPE :: k6 = 0.2
  138 + REALTYPE :: delta = 0.1
  139 + REALTYPE :: epsi = 0.70
  140 + REALTYPE :: r1 = 0.60
  141 + REALTYPE :: r2 = 0.40
  142 + REALTYPE :: vb = 1.2
  143 + REALTYPE :: k4 = 0.5
  144 + REALTYPE :: mu3 = 0.15
  145 + REALTYPE :: eta = 0.0
  146 + REALTYPE :: mu4 = 0.02
  147 + REALTYPE :: mu5 = 0.02
  148 + REALTYPE :: w_d1 = -1.0
  149 + REALTYPE :: w_d2 = -10.0
  150 + REALTYPE :: w_d3 = -100.0
  151 + REALTYPE :: bertha = 0.05
  152 + REALTYPE :: parcrit = 0.02
  153 + REALTYPE :: pmin = 0.05
  154 + REALTYPE, public :: kc = 0.03
  155 + integer :: out_unit
  156 + integer, parameter :: p=1,z=2,b=3,d1=4,n=5,a=6,l=7,d2=8,d3=9,d4=10
  157 +!EOP
  158 +!-----------------------------------------------------------------------
  159 +
  160 + contains
  161 +
  162 +!-----------------------------------------------------------------------
  163 +!BOP
  164 +!
  165 +! !IROUTINE: Initialise the bio module
  166 +!
  167 +! !INTERFACE:
  168 + subroutine init_bio_npzd4(namlst,fname,unit)
  169 +!
  170 +! !DESCRIPTION:
  171 +! Here, the bio namelist {\tt bio\_fasham.nml} is read and
  172 +! various variables (rates and settling velocities)
  173 +! are transformed into SI units.
  174 +!
  175 +! !USES:
  176 + IMPLICIT NONE
  177 +!
  178 +! !INPUT PARAMETERS:
  179 + integer, intent(in) :: namlst
  180 + character(len=*), intent(in) :: fname
  181 + character(len=20) :: pfile
  182 + integer, intent(in) :: unit
  183 +!
  184 +! !REVISION HISTORY:
  185 +! Original author(s): Hans Burchard & Karsten Bolding
  186 +!
  187 +! !LOCAL VARIABLES:
  188 + namelist /bio_npzd4_nml/ numc, &
  189 + p_initial,z_initial,b_initial,d_initial, &
  190 + l_initial,p0,z0,b0,vp,alpha,k1,k2,mu1,k5, &
  191 + gamma,w_p,gmax,k3,beta,mu2,k6,delta,epsi,r1,r2, &
  192 + vb,k4,mu3,eta,mu4,mu5,w_d1,w_d2,w_d3,kc, &
  193 + I_opt,inib, &
  194 + theta,w_pmax,w_pmin, &
  195 + bertha,parcrit,w_zmax,pmin
  196 +!EOP
  197 +!-----------------------------------------------------------------------
  198 +!BOC
  199 + LEVEL2 'init_bio_npzd4'
  200 +
  201 + open(namlst,file=fname,action='read',status='old',err=98)
  202 + read(namlst,nml=bio_npzd4_nml,err=99)
  203 + close(namlst)
  204 +
  205 + numcc=numc
  206 +
  207 +! Print some parameter values in standard output
  208 +! and save them in a separate file [out_fn]_fasham.par
  209 + pfile = trim(out_fn) // '_npzd4.par'
  210 + open(10,status='unknown',action='write',file=pfile)
  211 + LEVEL3 'Biogeochemical parameters saved in ', pfile
  212 + write(*,900) ' vp = ',vp
  213 + write(10,901) vp
  214 + write(*,900) ' alpha = ',alpha
  215 + write(10,901) alpha
  216 + write(*,900) ' inib = ',inib
  217 + write(10,901) inib
  218 + write(*,900) ' k1 = ',k1
  219 + write(10,901) k1
  220 + write(*,900) ' k2 = ',k2
  221 + write(10,901) k2
  222 + write(*,900) ' w_p = ',w_p
  223 + write(10,901) w_p
  224 + write(*,900) ' gmax = ',gmax
  225 + write(10,901) gmax
  226 +
  227 +900 format (a,f8.5)
  228 +901 format (f8.5)
  229 +
  230 +! Conversion from day to second
  231 + vp = vp /secs_pr_day
  232 + vb = vb /secs_pr_day
  233 + mu1 = mu1 /secs_pr_day
  234 + mu2 = mu2 /secs_pr_day
  235 + mu3 = mu3 /secs_pr_day
  236 + mu4 = mu4 /secs_pr_day
  237 + mu5 = mu5 /secs_pr_day
  238 + gmax = gmax /secs_pr_day
  239 + w_p = w_p /secs_pr_day
  240 + w_pmin = w_pmin /secs_pr_day !DD
  241 + w_pmax = w_pmax /secs_pr_day !DD
  242 + w_zmax = w_zmax /secs_pr_day
  243 + theta = theta /secs_pr_day !DD
  244 + w_d1 = w_d1 /secs_pr_day
  245 + w_d2 = w_d2 /secs_pr_day
  246 + w_d3 = w_d3 /secs_pr_day
  247 + alpha= alpha/secs_pr_day
  248 + inib = inib /secs_pr_day !CHG1
  249 +
  250 + out_unit=unit
  251 +
  252 + LEVEL3 'NOCERA bio module initialised ...'
  253 +
  254 + return
  255 +
  256 +98 LEVEL2 'I could not open bio_npzd4.nml'
  257 + LEVEL2 'If thats not what you want you have to supply bio_npzd4.nml'
  258 + LEVEL2 'See the bio example on www.gotm.net for a working bio_npzd4.nml'
  259 + return
  260 +99 FATAL 'I could not read bio_npzd4.nml'
  261 + stop 'init_bio_npzd4'
  262 + end subroutine init_bio_npzd4
  263 +!EOC
  264 +
  265 +!-----------------------------------------------------------------------
  266 +!BOP
  267 +!
  268 +! !IROUTINE: Initialise the concentration variables
  269 +!
  270 +! !INTERFACE:
  271 + subroutine init_var_npzd4(nlev)
  272 +!
  273 +! !DESCRIPTION:
  274 +! Here, the the initial conditions are set and the settling velocities are
  275 +! transferred to all vertical levels. All concentrations are declared
  276 +! as non-negative variables, and it is defined which variables would be
  277 +! taken up by benthic filter feeders.
  278 +!
  279 +! !USES:
  280 + use observations, only: nprof,aprof !CHG3-5
  281 + use meanflow, only: nit,amm,T,S !CHG3-5
  282 +
  283 + IMPLICIT NONE
  284 +
  285 +!
  286 +! !INPUT PARAMETERS:
  287 + integer, intent(in) :: nlev
  288 +!
  289 +! !REVISION HISTORY:
  290 +! Original author(s): Hans Burchard & Karsten Bolding
  291 +
  292 +! !LOCAL VARIABLES:
  293 + integer :: i
  294 +!EOP
  295 +!-----------------------------------------------------------------------
  296 +!BOC
  297 + do i=1,nlev
  298 + cc(p ,i)=p_initial
  299 + cc(z ,i)=z_initial
  300 + cc(b ,i)=b_initial
  301 + cc(d1,i)=d_initial
  302 + cc(n ,i)=nprof(i) !CHG3
  303 + cc(a ,i)=aprof(i) !CHG5
  304 + cc(l ,i)=l_initial
  305 + cc(d2,i)=d_initial
  306 + cc(d3,i)=d_initial
  307 + cc(d4,i)=d_initial
  308 + end do
  309 +
  310 + do i=0,nlev
  311 + ws(z ,i)= _ZERO_
  312 + ws(b ,i)= _ZERO_
  313 + ws(n ,i)= _ZERO_
  314 + ws(a ,i)= _ZERO_
  315 + ws(l ,i)= _ZERO_
  316 + ws(p ,i)= w_p
  317 + ws(d1,i)= w_d1
  318 + ws(d2,i)= w_d2
  319 + ws(d3,i)= w_d3
  320 + end do
  321 +
  322 + posconc(p) = 1
  323 + posconc(z) = 1
  324 + posconc(b) = 1
  325 + posconc(d1) = 1
  326 + posconc(n) = 1
  327 + posconc(a) = 1
  328 + posconc(l) = 1
  329 + posconc(d2) = 1
  330 + posconc(d3) = 1
  331 +
  332 +#if 0
  333 + mussels_inhale(p) = .true.
  334 + mussels_inhale(z) = .true.
  335 + mussels_inhale(b) = .false.
  336 + mussels_inhale(d1) = .true.
  337 + mussels_inhale(n) = .false.
  338 + mussels_inhale(a) = .false.
  339 + mussels_inhale(l) = .true.
  340 + mussels_inhale(d2) = .true.
  341 + mussels_inhale(d3) = .true.
  342 +#endif
  343 +
  344 + LEVEL3 'NPZD4 variables initialised ...'
  345 +
  346 + return
  347 +
  348 + end subroutine init_var_npzd4
  349 +!EOC
  350 +
  351 +!-----------------------------------------------------------------------
  352 +!BOP
  353 +!
  354 +! !IROUTINE: Providing info on variables
  355 +!
  356 +! !INTERFACE:
  357 + subroutine var_info_npzd4()
  358 +!
  359 +! !DESCRIPTION:
  360 +! This subroutine provides information about the variable names as they
  361 +! will be used when storing data in NetCDF files.
  362 +!
  363 +! !USES:
  364 + IMPLICIT NONE
  365 +!
  366 +! !REVISION HISTORY:
  367 +! Original author(s): Hans Burchard & Karsten Bolding
  368 +!
  369 +! !LOCAL VARIABLES:
  370 +!EOP
  371 +!-----------------------------------------------------------------------
  372 +!BOC
  373 + var_names(1) = 'phy'
  374 + var_units(1) = 'mmol/m**3'
  375 + var_long(1) = 'phytoplankton'
  376 +
  377 + var_names(2) = 'zoo'
  378 + var_units(2) = 'mmol/m**3'
  379 + var_long(2) = 'zooplankton'
  380 +
  381 + var_names(3) = 'bac'
  382 + var_units(3) = 'mmol/m**3'
  383 + var_long(3) = 'bacteria'
  384 +
  385 + var_names(4) = 'dph'
  386 + var_units(4) = 'mmol/m**3'
  387 + var_long(4) = 'detritus'
  388 +
  389 + var_names(5) = 'nit'
  390 + var_units(5) = 'mmol/m**3'
  391 + var_long(5) = 'nitrate'
  392 +
  393 + var_names(6) = 'amm'
  394 + var_units(6) = 'mmol/m**3'
  395 + var_long(6) = 'ammonium'
  396 +
  397 + var_names(7) = 'ldn'
  398 + var_units(7) = 'mmol/m**3'
  399 + var_long(7) = 'labile_dissolved_organic_nitrogen'
  400 +
  401 + var_names(8) = 'fcp'
  402 + var_units(8) = 'mmol/m**3'
  403 + var_long(8) = 'fecal pellets'
  404 +
  405 + var_names(9) = 'dzo'
  406 + var_units(9) = 'mmol/m**3'
  407 + var_long(9) = 'dead zooplankton'
  408 +
  409 + var_names(10) = 'msn'
  410 + var_units(10) = 'mmol/m**3'
  411 + var_long(10) = 'marine snow'
  412 +
  413 + return
  414 + end subroutine var_info_npzd4
  415 +!EOC
  416 +
  417 +!-----------------------------------------------------------------------
  418 +!BOP
  419 +!
  420 +! !IROUTINE: Light properties for the Fasham model
  421 +!
  422 +! !INTERFACE:
  423 + subroutine light_npzd4(nlev,bioshade_feedback)
  424 +!
  425 +! !DESCRIPTION:
  426 +! Here, the photosynthetically available radiation is calculated
  427 +! by simply assuming that the short wave part of the total
  428 +! radiation is available for photosynthesis.
  429 +! The photosynthetically
  430 +! available radiation, $I_{PAR}$, follows from (\ref{light}).
  431 +! The user should make
  432 +! sure that this is consistent with the light class given in the
  433 +! {\tt extinct} namelist of the {\tt obs.nml} file.
  434 +! The self-shading effect is also calculated in this subroutine,
  435 +! which may be used to consider the effect of bio-turbidity also
  436 +! in the temperature equation (if {\tt bioshade\_feedback} is set
  437 +! to true in {\tt bio.nml}).
  438 +! For details, see section \ref{sec:do-bio}.
  439 +!
  440 +! !USES:
  441 + IMPLICIT NONE
  442 +!
  443 +! !INPUT PARAMETERS:
  444 + integer, intent(in) :: nlev
  445 + logical, intent(in) :: bioshade_feedback
  446 +!
  447 +! !REVISION HISTORY:
  448 +! Original author(s): Hans Burchard, Karsten Bolding
  449 +!
  450 +! !LOCAL VARIABLES:
  451 + integer :: i
  452 + REALTYPE :: zz,add
  453 +!EOP
  454 +!-----------------------------------------------------------------------
  455 +!BOC
  456 + zz = _ZERO_
  457 + add = _ZERO_
  458 + do i=nlev,1,-1
  459 + add=add+0.5*h(i)*(cc(p,i)+p0)
  460 + zz=zz+0.5*h(i)
  461 + par(i)=rad(nlev)*(1.-aa)*exp(-zz/g2)*exp(-kc*add)
  462 + add=add+0.5*h(i)*(cc(p,i)+p0)
  463 + zz=zz+0.5*h(i)
  464 + if (bioshade_feedback) bioshade_(i)=exp(-kc*add)
  465 + end do
  466 +
  467 +
  468 + return
  469 + end subroutine light_npzd4
  470 +!EOC
  471 +
  472 +!-----------------------------------------------------------------------
  473 +!BOP
  474 +!
  475 +! !IROUTINE: Right hand sides of geobiochemical model \label{sec:bio-fasham-rhs}
  476 +!
  477 +! !INTERFACE:
  478 + subroutine do_bio_npzd4(first,numc,nlev,cc,pp,dd)
  479 +!
  480 +! !DESCRIPTION:
  481 +!
  482 +! The \cite{Fashametal1990} model consisting of the $I=7$
  483 +! state variables phytoplankton, bacteria, detritus, zooplankton,
  484 +! nitrate, ammonium and dissolved organic nitrogen is described here
  485 +! in detail.
  486 +!
  487 +! Phytoplankton mortality and zooplankton grazing loss of phytoplankton:
  488 +! \begin{equation}\label{d13}
  489 +! d_{1,3} = \mu_1 \frac{c_1+c_{1}^{\min}}{K_5+c_1+c_{1}^{\min}}c_1+
  490 +! (1-\beta)\frac{g\rho_1 c_1^2}{K_3 \sum_{j=1}^3 \rho_jc_j
  491 +! + \sum_{j=1}^3 \rho_jc_j^2} (c_4+c_{4}^{\min}).
  492 +! \end{equation}
  493 +! Phytoplankton loss to LDON (labile dissolved organic nitrogen):
  494 +! \begin{equation}\label{d17}
  495 +! d_{1,7} = \gamma
  496 +! F(I_{PAR})\frac{\frac{c_5}{K_1}
  497 +! +\frac{c_6}{K_2}}{1+\frac{c_5}{K_1}+\frac{c_6}{K_2}}c_1,
  498 +! \end{equation}
  499 +! with
  500 +! \begin{equation}\label{FI}
  501 +! F(I_{PAR}) = \frac{V_p\alpha I_{PAR}(z)}{\left(V_p^2+\alpha^2(I_{PAR}(z))^2
  502 +! \right)^{1/2}}.
  503 +! \end{equation}
  504 +! With $I_{PAR}$ from (\ref{light}).
  505 +!
  506 +! Zooplankton grazing loss:
  507 +! \begin{equation}\label{di3}
  508 +! d_{2,3} = (1-\beta)\frac{g\rho_2 c_2^2}{K_3 \sum_{j=1}^3 \rho_jc_j
  509 +! + \sum_{j=1}^3 \rho_jc_j^2} (c_4+c_{4}^{\min}).
  510 +! \end{equation}
  511 +! Zooplankton grazing:
  512 +! \begin{equation}\label{di4}
  513 +! d_{i,4} = \beta\frac{g\rho_i c_i^2}{K_3 \sum_{j=1}^3 \rho_jc_j
  514 +! + \sum_{j=1}^3 \rho_jc_j^2} (c_4+c_{4}^{\min}), \quad i=1,\dots,3.
  515 +! \end{equation}
  516 +! Bacteria excretion rate:
  517 +! \begin{equation}\label{d26}
  518 +! d_{2,6} = \mu_3 c_2.
  519 +! \end{equation}
  520 +! Detritus breakdown rate:
  521 +! \begin{equation}\label{d37}
  522 +! d_{3,7} = \mu_4 c_3.
  523 +! \end{equation}
  524 +! Zooplankton losses to detritus, ammonium and LDON:
  525 +! \begin{equation}\label{d43}
  526 +! d_{4,3} = (1-\epsilon-\delta)\mu_2
  527 +! \frac{c_4+c_{4}^{\min}}{K_6+c_4+c_{4}^{\min}}c_4.
  528 +! \end{equation}
  529 +! \begin{equation}\label{d46}
  530 +! d_{4,6} = \epsilon\mu_2 \frac{c_4+c_{4}^{\min}}{K_6+c_4+c_{4}^{\min}}c_4.
  531 +! \end{equation}
  532 +! \begin{equation}\label{d47}
  533 +! d_{4,7} = \delta\mu_2 \frac{c_4+c_{4}^{\min}}{K_6+c_4+c_{4}^{\min}}c_4.
  534 +! \end{equation}
  535 +! Nitrate uptake by phytoplankton:
  536 +! \begin{equation}\label{d51}
  537 +! d_{5,1} = F(I_{PAR})\frac{\frac{c_5}{K_1}}{1+\frac{c_5}{K_1}
  538 +! +\frac{c_6}{K_2}}(c_1+c_{1}^{\min}).
  539 +! \end{equation}
  540 +! Ammonium uptake by phytoplankton:
  541 +! \begin{equation}\label{d61}
  542 +! d_{6,1} = F(I_{PAR})\frac{\frac{c_6}{K_2}}{1+\frac{c_5}{K_1}
  543 +! +\frac{c_6}{K_2}}(c_1+c_{1}^{\min}).
  544 +! \end{equation}
  545 +! Ammonium uptake by bacteria:
  546 +! \begin{equation}\label{d62}
  547 +! d_{6,2} = V_b \frac{\min(c_6,\eta c_7)}{K_4+\min(c_6,\eta c_7)+c_7}
  548 +! (c_2+c_{2}^{\min}).
  549 +! \end{equation}
  550 +! LDON uptake by bacteria:
  551 +! \begin{equation}\label{d72}
  552 +! d_{7,2} = V_b \frac{c_7}{K_4+\min(c_6,\eta c_7)+c_7} (c_2+c_{2}^{\min}).
  553 +! \end{equation}
  554 +!
  555 +! !USES:
  556 + IMPLICIT NONE
  557 +!
  558 +! !INPUT PARAMETERS:
  559 + logical, intent(in) :: first
  560 + integer, intent(in) :: numc,nlev
  561 + REALTYPE, intent(in) :: cc(1:numc,0:nlev)
  562 +!
  563 +! !OUTPUT PARAMETERS:
  564 + REALTYPE, intent(out) :: pp(1:numc,1:numc,0:nlev)
  565 + REALTYPE, intent(out) :: dd(1:numc,1:numc,0:nlev)
  566 +!
  567 +! !REVISION HISTORY:
  568 +! Original author(s): Hans Burchard, Karsten Bolding
  569 +!
  570 +! !LOCAL VARIABLES:
  571 + REALTYPE :: ff,fac,min67,q1,q2
  572 + REALTYPE :: Ps !CHG1
  573 + integer :: i,j,ci
  574 +!EOP
  575 +!-----------------------------------------------------------------------
  576 +!BOC
  577 +!KBK - is it necessary to initialise every time - expensive in a 3D model
  578 + pp = _ZERO_
  579 + dd = _ZERO_
  580 +
  581 + do ci=1,nlev
  582 +
  583 +!CHG1
  584 +! Smith (1936) - saturation (default)
  585 +! ff= vp*alpha*par(ci)/sqrt(vp**2+alpha**2*par(ci)**2)
  586 +! Blackman (1919)
  587 +! if (par(ci) .lt. vp/alpha) then
  588 +! ff=alpha*par(ci)
  589 +! else
  590 +! ff=vp
  591 +! endif
  592 +! Steele (1962) - inhibition
  593 +! ff= vp*((par(ci)/I_opt)*exp(1-(par(ci)/I_opt)))
  594 +! Parker (1974) - inhibition
  595 +! ff= vp*((par(ci)/I_opt)*exp(1-(par(ci)/I_opt)))**2
  596 +! Platt et al. (1980) - inhibition
  597 + Ps= vp/((alpha/(alpha+inib))*(alpha/(alpha+inib))**(inib/alpha))
  598 + ff= Ps*(1.-exp(-1.*alpha*par(ci)/Ps))*exp(-1.*inib*par(ci)/Ps)
  599 +! --------------------------------------------------------------------
  600 +
  601 +
  602 + q1=(cc(n,ci)/k1)/(1.+cc(n,ci)/k1+cc(a,ci)/k2)
  603 + q2=(cc(a,ci)/k2)/(1.+cc(n,ci)/k1+cc(a,ci)/k2)
  604 + fac=(cc(z,ci)+z0)/(k3*(r1*cc(p,ci)+r2*cc(b,ci))+ &
  605 + r1*cc(p,ci)**2+r2*cc(b,ci)**2)
  606 + min67=min(cc(a,ci),eta*cc(l,ci))
  607 +
  608 + ! Light and nutrient limitation factors
  609 + lumlim1(ci)=ff
  610 + nitlim1(ci)=q1
  611 + ammlim1(ci)=q2
  612 +
  613 + ! 1. We remove the grazing of fecal pellets by zooplankton since it
  614 + ! produces those
  615 + dd(p,d1,ci)=mu1*(cc(p,ci)+p0)/(k5+cc(p,ci)+p0)*cc(p,ci) &
  616 + +(1.-beta)*gmax*r1*cc(p,ci)**2*fac
  617 + dd(p,l,ci)=gamma*ff*(q1+q2)*cc(p,ci)
  618 + dd(b,d1,ci)=(1.-beta)*gmax*r2*cc(b,ci)**2*fac
  619 + dd(p,z,ci)=beta*gmax*r1*cc(p,ci)**2*fac
  620 + dd(b,z,ci)=beta*gmax*r2*cc(b,ci)**2*fac
  621 + dd(b,a,ci)=mu3*cc(b,ci)
  622 + dd(d1,l,ci)=mu4*cc(d1,ci)
  623 + dd(d2,l,ci)=mu4*cc(d2,ci)
  624 + dd(d3,l,ci)=mu4*cc(d3,ci)
  625 + dd(z,d2,ci)=mu5*gmax*(r1*cc(p,ci)**2+r2*cc(b,ci)**2)*fac
  626 + dd(z,d3,ci)=(1.-epsi-delta)*mu2*(cc(z,ci)+z0)/(k6+cc(z,ci)+z0)*cc(z,ci)
  627 + dd(z,a,ci)=epsi*mu2*(cc(z,ci)+z0)/(k6+cc(z,ci)+z0)*cc(z,ci)
  628 + dd(z,l,ci)=delta*mu2*(cc(z,ci)+z0)/(k6+cc(z,ci)+z0)*cc(z,ci)
  629 + dd(n,p,ci)=ff*q1*(cc(p,ci)+p0)
  630 + dd(a,p,ci)=ff*q2*(cc(p,ci)+p0)
  631 + dd(a,b,ci)=vb*min67/(k4+min67+cc(l,ci))*(cc(b,ci)+b0)
  632 + dd(l,b,ci)=vb*cc(l,ci)/(k4+min67+cc(l,ci))*(cc(b,ci)+b0)
  633 +
  634 + ppnet(ci) =(dd(a,p,ci)+dd(n,p,ci)-dd(p,d1,ci)-dd(p,z,ci)-dd(p,l,ci))*secs_pr_day
  635 +
  636 + !Diurnal migration of zooplancton as a function of light and phytoplancton
  637 + !(food) concentration
  638 + !ws(z,ci) = -1.0*w_zmax*tanh(bertha*(par(ci)-parcrit))
  639 + if (cc(p,ci) .lt. pmin) then
  640 + ws(z,ci) = -1.0*w_zmax*tanh(bertha*(par(ci)-parcrit))
  641 + else
  642 + ws(z,ci) = 0.0
  643 + end if
  644 +
  645 +
  646 + do i=1,numc
  647 + do j=1,numc
  648 + pp(i,j,ci)=dd(j,i,ci)
  649 + end do
  650 + end do
  651 + end do
  652 +
  653 + return
  654 + end subroutine do_bio_npzd4
  655 +!EOC
  656 +
  657 +!-----------------------------------------------------------------------
  658 +!BOP
  659 +!
  660 +! !IROUTINE: Finish the bio calculations
  661 +!
  662 +! !INTERFACE:
  663 + subroutine end_bio_npzd4
  664 +!
  665 +! !DESCRIPTION:
  666 +! Nothing done yet --- supplied for completeness.
  667 +!
  668 +! !USES:
  669 + IMPLICIT NONE
  670 +!
  671 +! !REVISION HISTORY:
  672 +! Original author(s): Hans Burchard & Karsten Bolding
  673 +!
  674 +!EOP
  675 +!-----------------------------------------------------------------------
  676 +!BOC
  677 +
  678 + return
  679 + end subroutine end_bio_npzd4
  680 +!EOC
  681 +
  682 +!-----------------------------------------------------------------------
  683 +
  684 + end module bio_npzd4
  685 +
  686 +!-----------------------------------------------------------------------
  687 +! Copyright by the GOTM-team under the GNU Public License - www.gnu.org
  688 +!-----------------------------------------------------------------------
... ...