Changeset 7975


Ignore:
Timestamp:
2017-04-26T16:30:14+02:00 (3 years ago)
Author:
marc
Message:

Removed plankton processes from trcbio_medusa.F90 into extra routines

Location:
branches/UKMO/dev_r5518_medusa_chg_trc_bio_medusa/NEMOGCM/NEMO/TOP_SRC/MEDUSA
Files:
3 added
4 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/dev_r5518_medusa_chg_trc_bio_medusa/NEMOGCM/NEMO/TOP_SRC/MEDUSA/air_sea.F90

    r7958 r7975  
    6363      USE mocsy_wrapper,     ONLY: mocsy_interface 
    6464      USE oce,               ONLY: PCO2a_in_cpl 
    65       USE par_oce,           ONLY: jpim1, jpjm1, jpk 
     65      USE par_oce,           ONLY: jpim1, jpjm1 
    6666      USE sbc_oce,           ONLY: fr_i, lk_oasis, qsr, wndm 
    6767      USE sms_medusa,        ONLY: jdms, jdms_input, jdms_model,          & 
  • branches/UKMO/dev_r5518_medusa_chg_trc_bio_medusa/NEMOGCM/NEMO/TOP_SRC/MEDUSA/bio_medusa_mod.F90

    r7958 r7975  
    2020 
    2121   !! model state variables 
    22    REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::    zchn,zchd,zphn,zphd,zpds,zzmi 
    23    REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::    zzme,zdet,zdtc,zdin,zsil,zfer 
    24 # if defined key_roam 
    25    REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::    zdic, zalk, zoxy 
    26    REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::    ztmp, zsal 
     22   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zchn,zchd,zphn,zphd,zpds,zzmi 
     23   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zzme,zdet,zdtc,zdin,zsil,zfer 
     24# if defined key_roam 
     25   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zdic, zalk, zoxy 
     26   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: ztmp, zsal 
    2727# endif 
    2828# if defined key_mocsy 
    29    REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::    zpho 
    30 # endif 
    31  
     29   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: zpho 
     30# endif 
     31 
     32   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fthetan,faln,fchn1,fchn 
     33   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fjln,fprn,frn 
     34   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fthetad,fald,fchd1,fchd 
     35   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fjld,fprd,frd 
     36 
     37   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fjlim_pn, fjlim_pd 
     38   !! AXY (03/02/11): add in Liebig terms 
     39   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fpnlim, fpdlim 
     40   !! AXY (16/07/09): add in Eppley curve functionality 
     41   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fun_T,xvpnT,xvpdT 
     42 
     43   !! AXY (16/05/11): per Katya's prompting, add in new T-dependence 
     44   !!                 for phytoplankton growth only (i.e. no change 
     45   !!                 for remineralisation) 
     46   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fun_Q10 
    3247   !! AXY (01/03/10): add in mixed layer PP diagnostics 
    3348   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fprn_ml,fprd_ml 
     49   !! 
     50   !! nutrient limiting factors 
     51   !! N and Fe (renaming ffln to ffln2 to avoid conflict with 
     52   !! ffln in module sms_medusa - marc 25/4/17) 
     53   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fnln,ffln2 
     54   !! N, Fe and Si 
     55   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fnld,ffld,fsld,fsld2 
     56   !! 
     57   !! silicon cycle 
     58   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fsin,fnsi,fprds,fsdiss 
    3459 
    3560   !! Variable for iron-ligand system 
    3661   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: xFree 
     62 
     63   !! Microzooplankton grazing 
     64   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fmi1,fmi,fgmipn,fgmid,fgmidc 
     65   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: finmi,ficmi,fstarmi,fmith 
     66   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fmigrow,fmiexcr,fmiresp 
     67   !! 
     68   !! Mesozooplankton grazing 
     69   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fme1,fme,fgmepn,fgmepd 
     70   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fgmepds,fgmezmi,fgmed,fgmedc 
     71   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: finme,ficme,fstarme,fmeth 
     72   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fmegrow,fmeexcr,fmeresp 
     73   !! 
     74   !! mortality/Remineralisation (defunct parameter "fz" removed) 
     75   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fdpn,fdpd,fdpds,fdzmi,fdzme,fdd 
     76   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fdpn2,fdpd2,fdpds2,fdzmi2,fdzme2 
    3777 
    3878   !! Mortality/Remineralisation 
     
    4484 
    4585   !! Particle flux 
     86   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: fdep1,fcaco3 
     87 
    4688   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: ffastn,ffastsi,ffastfe 
    4789   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: ffastc,ffastca 
     
    67109   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: ftot_pn,ftot_pd 
    68110   REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: ftot_zmi,ftot_zme,ftot_det,ftot_dtc 
     111 
     112   !! diagnose fluxes (should only be used in 1D runs) 
     113   INTEGER                               :: idf, idfval 
    69114 
    70115   !! Nitrogen and silicon production and consumption 
     
    208253               ztmp(jpi,jpj),zsal(jpi,jpj),                           & 
    209254               zpho(jpi,jpj),                                         & 
     255               fthetan(jpi,jpj),faln(jpi,jpj),fchn1(jpi,jpj),         & 
     256               fchn(jpi,jpj),fjln(jpi,jpj),fprn(jpi,jpj),             & 
     257               frn(jpi,jpj),                                          & 
     258               fthetad(jpi,jpj),fald(jpi,jpj),fchd1(jpi,jpj),         & 
     259               fchd(jpi,jpj),fjld(jpi,jpj),fprd(jpi,jpj),             & 
     260               frd(jpi,jpj),                                          & 
     261               fjlim_pn(jpi,jpj), fjlim_pd(jpi,jpj),                  & 
     262               fpnlim(jpi,jpj), fpdlim(jpi,jpj),                      & 
     263               fun_T(jpi,jpj),xvpnT(jpi,jpj),xvpdT(jpi,jpj),          & 
     264               fun_Q10(jpi,jpj),                                      & 
    210265               fprn_ml(jpi,jpj),fprd_ml(jpi,jpj),                     & 
     266               fnln(jpi,jpj),ffln2(jpi,jpj),                          & 
     267               fnld(jpi,jpj),ffld(jpi,jpj),fsld(jpi,jpj),             & 
     268               fsld2(jpi,jpj),                                        & 
     269               fsin(jpi,jpj),fnsi(jpi,jpj),fprds(jpi,jpj),            & 
     270               fsdiss(jpi,jpj),                                       & 
    211271               xFree(jpi,jpj),                                        & 
     272               fmi1(jpi,jpj),fmi(jpi,jpj),fgmipn(jpi,jpj),            & 
     273               fgmid(jpi,jpj),fgmidc(jpi,jpj),                        & 
     274               finmi(jpi,jpj),ficmi(jpi,jpj),fstarmi(jpi,jpj),        & 
     275               fmith(jpi,jpj),fmigrow(jpi,jpj),fmiexcr(jpi,jpj),      & 
     276               fmiresp(jpi,jpj),                                      & 
     277               fme1(jpi,jpj),fme(jpi,jpj),fgmepn(jpi,jpj),            & 
     278               fgmepd(jpi,jpj),fgmepds(jpi,jpj),fgmezmi(jpi,jpj),     & 
     279               fgmed(jpi,jpj),fgmedc(jpi,jpj),                        & 
     280               finme(jpi,jpj),ficme(jpi,jpj),fstarme(jpi,jpj),        & 
     281               fmeth(jpi,jpj),fmegrow(jpi,jpj),fmeexcr(jpi,jpj),      & 
     282               fmeresp(jpi,jpj),                                      & 
     283               fdpn(jpi,jpj),fdpd(jpi,jpj),fdpds(jpi,jpj),            & 
     284               fdzmi(jpi,jpj),fdzme(jpi,jpj),fdd(jpi,jpj),            & 
     285               fdpn2(jpi,jpj),fdpd2(jpi,jpj),fdpds2(jpi,jpj),         & 
     286               fdzmi2(jpi,jpj),fdzme2(jpi,jpj),                       & 
    212287               fslownflux(jpi,jpj),fslowcflux(jpi,jpj),               & 
    213288               fregenfast(jpi,jpj),fregenfastsi(jpi,jpj),             & 
     
    215290               fregenfastc(jpi,jpj),                                  & 
    216291# endif 
     292               fdep1(jpi,jpj),fcaco3(jpi,jpj),                        & 
    217293               ffastn(jpi,jpj),ffastsi(jpi,jpj),ffastfe(jpi,jpj),     & 
    218294               ffastc(jpi,jpj),ffastca(jpi,jpj),                      & 
  • branches/UKMO/dev_r5518_medusa_chg_trc_bio_medusa/NEMOGCM/NEMO/TOP_SRC/MEDUSA/carb_chem.F90

    r7958 r7975  
    4040                                   f_schmidtco2, f_TALK, f_TDIC, f_xco2a, & 
    4141                                   zalk, zdic, zpho, zsal, zsil, ztmp  
    42       USE dom_oce,           ONLY: gdept_0, gdept_n, gdepw_n, gphit,      & 
    43                                    mbathy, tmask 
     42      USE dom_oce,           ONLY: gdept_0, gdept_n, gdepw_0, gdepw_n,    & 
     43                                   gphit, mbathy, tmask 
    4444      USE in_out_manager,    ONLY: lwp, numout 
    4545      USE mocsy_wrapper,     ONLY: mocsy_interface 
  • branches/UKMO/dev_r5518_medusa_chg_trc_bio_medusa/NEMOGCM/NEMO/TOP_SRC/MEDUSA/trcbio_medusa.F90

    r7958 r7975  
    100100      USE carb_chem_mod,              ONLY: carb_chem 
    101101      USE air_sea_mod,                ONLY: air_sea 
     102      USE plankton_mod,               ONLY: plankton 
    102103      USE bio_medusa_diag_slice_mod,  ONLY: bio_medusa_diag_slice 
    103104      USE bio_medusa_fin_mod,         ONLY: bio_medusa_fin 
     
    181182      !! 
    182183      !! primary production and chl related quantities       
    183       REAL(wp), DIMENSION(jpi,jpj) ::    fthetan,faln,fchn1,fchn,fjln,fprn,frn 
    184       REAL(wp), DIMENSION(jpi,jpj) ::    fthetad,fald,fchd1,fchd,fjld,fprd,frd 
     184!      REAL(wp), DIMENSION(jpi,jpj) ::    fthetan,faln,fchn1,fchn,fjln,fprn,frn 
     185!      REAL(wp), DIMENSION(jpi,jpj) ::    fthetad,fald,fchd1,fchd,fjld,fprd,frd 
    185186      !! AXY (23/11/16): add in light-only limitation term (normalised 0-1 range) 
    186       REAL(wp), DIMENSION(jpi,jpj) ::    fjlim_pn, fjlim_pd 
     187!      REAL(wp), DIMENSION(jpi,jpj) ::    fjlim_pn, fjlim_pd 
    187188      !! AXY (03/02/11): add in Liebig terms 
    188       REAL(wp), DIMENSION(jpi,jpj) ::    fpnlim, fpdlim 
     189!      REAL(wp), DIMENSION(jpi,jpj) ::    fpnlim, fpdlim 
    189190      !! AXY (16/07/09): add in Eppley curve functionality 
    190       REAL(wp), DIMENSION(jpi,jpj) ::    fun_T,xvpnT,xvpdT 
     191!      REAL(wp), DIMENSION(jpi,jpj) ::    fun_T,xvpnT,xvpdT 
    191192      INTEGER  ::    ieppley 
    192193      !! AXY (16/05/11): per Katya's prompting, add in new T-dependence 
    193194      !!                 for phytoplankton growth only (i.e. no change 
    194195      !!                 for remineralisation) 
    195       REAL(wp), DIMENSION(jpi,jpj) ::    fun_Q10 
     196!      REAL(wp), DIMENSION(jpi,jpj) ::    fun_Q10 
    196197      !! AXY (01/03/10): add in mixed layer PP diagnostics 
    197198!      REAL(wp), DIMENSION(jpi,jpj) ::    fprn_ml,fprd_ml 
    198199      !! 
    199200      !! nutrient limiting factors 
    200       REAL(wp), DIMENSION(jpi,jpj) ::    fnln,ffln            !! N and Fe 
    201       REAL(wp), DIMENSION(jpi,jpj) ::    fnld,ffld,fsld,fsld2 !! N, Fe and Si 
     201!      REAL(wp), DIMENSION(jpi,jpj) ::    fnln,ffln2            !! N and Fe 
     202!      REAL(wp), DIMENSION(jpi,jpj) ::    fnld,ffld,fsld,fsld2 !! N, Fe and Si 
    202203      !! 
    203204      !! silicon cycle 
    204       REAL(wp), DIMENSION(jpi,jpj) ::    fsin,fnsi,fprds,fsdiss 
     205!      REAL(wp), DIMENSION(jpi,jpj) ::    fsin,fnsi,fprds,fsdiss 
    205206      REAL(wp)                     ::    fsin1,fnsi1,fnsi2 
    206207      !! 
     
    223224      !! 
    224225      !! microzooplankton grazing 
    225       REAL(wp), DIMENSION(jpi,jpj) ::    fmi1,fmi,fgmipn,fgmid,fgmidc 
    226       REAL(wp), DIMENSION(jpi,jpj) ::    finmi,ficmi,fstarmi,fmith,fmigrow,fmiexcr,fmiresp 
     226!      REAL(wp), DIMENSION(jpi,jpj) ::    fmi1,fmi,fgmipn,fgmid,fgmidc 
     227!      REAL(wp), DIMENSION(jpi,jpj) ::    finmi,ficmi,fstarmi,fmith,fmigrow,fmiexcr,fmiresp 
    227228      !! 
    228229      !! mesozooplankton grazing 
    229       REAL(wp), DIMENSION(jpi,jpj) ::    fme1,fme,fgmepn,fgmepd,fgmepds,fgmezmi,fgmed,fgmedc 
    230       REAL(wp), DIMENSION(jpi,jpj) ::    finme,ficme,fstarme,fmeth,fmegrow,fmeexcr,fmeresp 
     230!      REAL(wp), DIMENSION(jpi,jpj) ::    fme1,fme,fgmepn,fgmepd,fgmepds,fgmezmi,fgmed,fgmedc 
     231!      REAL(wp), DIMENSION(jpi,jpj) ::    finme,ficme,fstarme,fmeth,fmegrow,fmeexcr,fmeresp 
    231232      !! 
    232233      !! mortality/Remineralisation (defunct parameter "fz" removed) 
    233       REAL(wp), DIMENSION(jpi,jpj) ::    fdpn,fdpd,fdpds,fdzmi,fdzme,fdd 
     234!      REAL(wp), DIMENSION(jpi,jpj) ::    fdpn,fdpd,fdpds,fdzmi,fdzme,fdd 
    234235# if defined key_roam 
    235236      REAL(wp), DIMENSION(jpi,jpj) ::    fddc 
    236237# endif 
    237       REAL(wp), DIMENSION(jpi,jpj) ::    fdpn2,fdpd2,fdpds2,fdzmi2,fdzme2 
     238!      REAL(wp), DIMENSION(jpi,jpj) ::    fdpn2,fdpd2,fdpds2,fdzmi2,fdzme2 
    238239      REAL(wp), DIMENSION(jpi,jpj) ::    fslown, fslowc 
    239240!      REAL(wp), DIMENSION(jpi,jpj) ::    fslownflux, fslowcflux 
     
    247248      !! 
    248249      !! particle flux 
    249       REAL(WP), DIMENSION(jpi,jpj) ::    fdep1,fcaco3 
     250!      REAL(WP), DIMENSION(jpi,jpj) ::    fdep1,fcaco3 
    250251      REAL(WP), DIMENSION(jpi,jpj) ::    ftempn,ftempsi,ftempfe,ftempc,ftempca 
    251252      REAL(wp), DIMENSION(jpi,jpj) ::    freminn,freminsi,freminfe,freminc,freminca 
     
    297298      !! 
    298299      !! diagnose fluxes (should only be used in 1D runs) 
    299       INTEGER  ::    idf, idfval 
     300!      INTEGER  ::    idf, idfval 
    300301      !! 
    301302      !! nitrogen and silicon production and consumption 
     
    844845# endif 
    845846 
     847         !!------------------------------------------------------------------ 
     848         !! Phytoplankton growth, zooplankton grazing and miscellaneous 
     849         !! plankton losses.  
     850         !!------------------------------------------------------------------ 
     851         CALL plankton( jk ) 
     852 
     853! Detritus process - marc 
     854 
    846855         DO jj = 2,jpjm1 
    847856         DO ji = 2,jpim1 
    848857            !! OPEN wet point IF..THEN loop 
    849858            if (tmask(ji,jj,jk) == 1) then 
    850  
    851                !!---------------------------------------------------------------------- 
    852                !! Chlorophyll calculations 
    853                !!---------------------------------------------------------------------- 
    854                !! 
    855                !! non-diatoms 
    856           if (zphn(ji,jj).GT.rsmall) then 
    857                   fthetan(ji,jj) = max(tiny(zchn(ji,jj)), (zchn(ji,jj) * xxi) / (zphn(ji,jj) + tiny(zphn(ji,jj)))) 
    858                   faln(ji,jj)    = xaln * fthetan(ji,jj) 
    859                else 
    860                   fthetan(ji,jj) = 0. 
    861                   faln(ji,jj)    = 0. 
    862                endif 
    863                !! 
    864                !! diatoms 
    865           if (zphd(ji,jj).GT.rsmall) then 
    866                   fthetad(ji,jj) = max(tiny(zchd(ji,jj)), (zchd(ji,jj) * xxi) / (zphd(ji,jj) + tiny(zphd(ji,jj)))) 
    867                   fald(ji,jj)    = xald * fthetad(ji,jj) 
    868                else 
    869                   fthetad(ji,jj) = 0. 
    870                   fald(ji,jj)    = 0. 
    871                endif 
    872  
    873 # if defined key_debug_medusa 
    874                !! report biological calculations 
    875                if (idf.eq.1.AND.idfval.eq.1) then 
    876                   IF (lwp) write (numout,*) '------------------------------' 
    877                   IF (lwp) write (numout,*) 'faln(',jk,') = ', faln(ji,jj) 
    878                   IF (lwp) write (numout,*) 'fald(',jk,') = ', fald(ji,jj) 
    879                endif 
    880 # endif 
    881  
    882                !!---------------------------------------------------------------------- 
    883                !! Phytoplankton light limitation 
    884                !!---------------------------------------------------------------------- 
    885                !! 
    886                !! It is assumed xpar is the depth-averaged (vertical layer) PAR  
    887                !! Light limitation (check self-shading) in W/m2 
    888                !! 
    889                !! Note that there is no temperature dependence in phytoplankton 
    890                !! growth rate or any other function.  
    891                !! In calculation of Chl/Phy ratio tiny(phyto) is introduced to avoid 
    892                !! NaNs in case of Phy==0.   
    893                !! 
    894                !! fthetad and fthetan are Chl:C ratio (gChl/gC) in diat and non-diat:  
    895                !! for 1:1 Chl:P ratio (mgChl/mmolN) theta=0.012 
    896                !! 
    897                !! AXY (16/07/09) 
    898                !! temperature for new Eppley style phytoplankton growth 
    899                fun_T(ji,jj)   = 1.066**(1.0 * tsn(ji,jj,jk,jp_tem)) 
    900                !! AXY (16/05/11): add in new Q10 (1.5, not 2.0) for 
    901                !phytoplankton 
    902                !!                 growth; remin. unaffected 
    903                fun_Q10(ji,jj) = jq10**((tsn(ji,jj,jk,jp_tem) - 0.0) / 10.0) 
    904                if (jphy.eq.1) then 
    905                   xvpnT(ji,jj) = xvpn * fun_T(ji,jj) 
    906                   xvpdT(ji,jj) = xvpd * fun_T(ji,jj) 
    907                elseif (jphy.eq.2) then 
    908                   xvpnT(ji,jj) = xvpn * fun_Q10(ji,jj) 
    909                   xvpdT(ji,jj) = xvpd * fun_Q10(ji,jj) 
    910                else 
    911                   xvpnT(ji,jj) = xvpn 
    912                   xvpdT(ji,jj) = xvpd 
    913                endif 
    914                !! 
    915                !! non-diatoms 
    916                fchn1(ji,jj)   = (xvpnT(ji,jj) * xvpnT(ji,jj)) + (faln(ji,jj) * faln(ji,jj) * xpar(ji,jj,jk) * xpar(ji,jj,jk)) 
    917                if (fchn1(ji,jj).GT.rsmall) then 
    918                   fchn(ji,jj)    = xvpnT(ji,jj) / (sqrt(fchn1(ji,jj)) + tiny(fchn1(ji,jj))) 
    919                else 
    920                   fchn(ji,jj)    = 0. 
    921                endif 
    922                fjln(ji,jj)    = fchn(ji,jj) * faln(ji,jj) * xpar(ji,jj,jk) !! non-diatom J term 
    923                fjlim_pn(ji,jj) = fjln(ji,jj) / xvpnT(ji,jj) 
    924                !! 
    925                !! diatoms 
    926                fchd1(ji,jj)   = (xvpdT(ji,jj) * xvpdT(ji,jj)) + (fald(ji,jj) * fald(ji,jj) * xpar(ji,jj,jk) * xpar(ji,jj,jk)) 
    927                if (fchd1(ji,jj).GT.rsmall) then 
    928                   fchd(ji,jj)    = xvpdT(ji,jj) / (sqrt(fchd1(ji,jj)) + tiny(fchd1(ji,jj))) 
    929                else 
    930                   fchd(ji,jj)    = 0. 
    931                endif 
    932                fjld(ji,jj)    = fchd(ji,jj) * fald(ji,jj) * xpar(ji,jj,jk) !! diatom J term 
    933                fjlim_pd(ji,jj) = fjld(ji,jj) / xvpdT(ji,jj) 
    934        
    935 # if defined key_debug_medusa 
    936                !! report phytoplankton light limitation 
    937                if (idf.eq.1.AND.idfval.eq.1) then 
    938                   IF (lwp) write (numout,*) '------------------------------' 
    939                   IF (lwp) write (numout,*) 'fchn(',jk,') = ', fchn(ji,jj) 
    940                   IF (lwp) write (numout,*) 'fchd(',jk,') = ', fchd(ji,jj) 
    941                   IF (lwp) write (numout,*) 'fjln(',jk,') = ', fjln(ji,jj) 
    942                   IF (lwp) write (numout,*) 'fjld(',jk,') = ', fjld(ji,jj) 
    943                endif 
    944 # endif 
    945  
    946                !!---------------------------------------------------------------------- 
    947                !! Phytoplankton nutrient limitation 
    948                !!---------------------------------------------------------------------- 
    949                !! 
    950                !! non-diatoms (N, Fe) 
    951                fnln(ji,jj) = zdin(ji,jj) / (zdin(ji,jj) + xnln) !! non-diatom Qn term 
    952                ffln(ji,jj) = zfer(ji,jj) / (zfer(ji,jj) + xfln) !! non-diatom Qf term 
    953                !! 
    954                !! diatoms (N, Si, Fe) 
    955                fnld(ji,jj) = zdin(ji,jj) / (zdin(ji,jj) + xnld) !! diatom Qn term 
    956                fsld(ji,jj) = zsil(ji,jj) / (zsil(ji,jj) + xsld) !! diatom Qs term 
    957                ffld(ji,jj) = zfer(ji,jj) / (zfer(ji,jj) + xfld) !! diatom Qf term 
    958  
    959 # if defined key_debug_medusa 
    960                !! report phytoplankton nutrient limitation 
    961                if (idf.eq.1.AND.idfval.eq.1) then 
    962                   IF (lwp) write (numout,*) '------------------------------' 
    963                   IF (lwp) write (numout,*) 'fnln(',jk,') = ', fnln(ji,jj) 
    964                   IF (lwp) write (numout,*) 'fnld(',jk,') = ', fnld(ji,jj) 
    965                   IF (lwp) write (numout,*) 'ffln(',jk,') = ', ffln(ji,jj) 
    966                   IF (lwp) write (numout,*) 'ffld(',jk,') = ', ffld(ji,jj) 
    967                   IF (lwp) write (numout,*) 'fsld(',jk,') = ', fsld(ji,jj) 
    968                endif 
    969 # endif 
    970  
    971                !!---------------------------------------------------------------------- 
    972                !! Primary production (non-diatoms) 
    973                !! (note: still needs multiplying by phytoplankton concentration) 
    974                !!---------------------------------------------------------------------- 
    975                !! 
    976                if (jliebig .eq. 0) then 
    977                   !! multiplicative nutrient limitation 
    978                   fpnlim(ji,jj) = fnln(ji,jj) * ffln(ji,jj) 
    979                elseif (jliebig .eq. 1) then 
    980                   !! Liebig Law (= most limiting) nutrient limitation 
    981                   fpnlim(ji,jj) = min(fnln(ji,jj), ffln(ji,jj)) 
    982                endif 
    983                fprn(ji,jj) = fjln(ji,jj) * fpnlim(ji,jj) 
    984  
    985                !!---------------------------------------------------------------------- 
    986                !! Primary production (diatoms) 
    987                !! (note: still needs multiplying by phytoplankton concentration) 
    988                !! 
    989                !! production here is split between nitrogen production and that of 
    990                !! silicon; depending upon the "intracellular" ratio of Si:N, model 
    991                !! diatoms will uptake nitrogen/silicon differentially; this borrows 
    992                !! from the diatom model of Mongin et al. (2006) 
    993                !!---------------------------------------------------------------------- 
    994                !! 
    995                if (jliebig .eq. 0) then 
    996                   !! multiplicative nutrient limitation 
    997                   fpdlim(ji,jj) = fnld(ji,jj) * ffld(ji,jj) 
    998                elseif (jliebig .eq. 1) then 
    999                   !! Liebig Law (= most limiting) nutrient limitation 
    1000                   fpdlim(ji,jj) = min(fnld(ji,jj), ffld(ji,jj)) 
    1001                endif 
    1002                !! 
    1003           if (zphd(ji,jj).GT.rsmall .AND. zpds(ji,jj).GT.rsmall) then 
    1004                   !! "intracellular" elemental ratios 
    1005                   ! fsin(ji,jj)  = zpds(ji,jj) / (zphd(ji,jj) + tiny(zphd(ji,jj))) 
    1006                   ! fnsi(ji,jj)  = zphd(ji,jj) / (zpds(ji,jj) + tiny(zpds(ji,jj))) 
    1007                   fsin(ji,jj) = 0.0 
    1008                   IF( zphd(ji,jj) .GT. rsmall) fsin(ji,jj)  = zpds(ji,jj) / zphd(ji,jj) 
    1009                   fnsi(ji,jj) = 0.0 
    1010                   IF( zpds(ji,jj) .GT. rsmall) fnsi(ji,jj)  = zphd(ji,jj) / zpds(ji,jj) 
    1011                   !! AXY (23/02/10): these next variables derive from Mongin et al. (2003) 
    1012                   fsin1 = 3.0 * xsin0 !! = 0.6 
    1013                   fnsi1 = 1.0 / fsin1 !! = 1.667 
    1014                   fnsi2 = 1.0 / xsin0 !! = 5.0 
    1015                   !! 
    1016                   !! conditionalities based on ratios 
    1017                   !! nitrogen (and iron and carbon) 
    1018                   if (fsin(ji,jj).le.xsin0) then 
    1019                      fprd(ji,jj)  = 0.0 
    1020                      fsld2(ji,jj) = 0.0 
    1021                   elseif (fsin(ji,jj).lt.fsin1) then 
    1022                      fprd(ji,jj)  = xuif * ((fsin(ji,jj) - xsin0) / (fsin(ji,jj) + tiny(fsin(ji,jj)))) * (fjld(ji,jj) * fpdlim(ji,jj)) 
    1023                      fsld2(ji,jj) = xuif * ((fsin(ji,jj) - xsin0) / (fsin(ji,jj) + tiny(fsin(ji,jj)))) 
    1024                   elseif (fsin(ji,jj).ge.fsin1) then 
    1025                      fprd(ji,jj)  = (fjld(ji,jj) * fpdlim(ji,jj)) 
    1026                      fsld2(ji,jj) = 1.0 
    1027                   endif 
    1028                   !! 
    1029                   !! silicon 
    1030                   if (fsin(ji,jj).lt.fnsi1) then 
    1031                      fprds(ji,jj) = (fjld(ji,jj) * fsld(ji,jj)) 
    1032                   elseif (fsin(ji,jj).lt.fnsi2) then 
    1033                      fprds(ji,jj) = xuif * ((fnsi(ji,jj) - xnsi0) / (fnsi(ji,jj) + tiny(fnsi(ji,jj)))) * (fjld(ji,jj) * fsld(ji,jj)) 
    1034                   else 
    1035                      fprds(ji,jj) = 0.0 
    1036                   endif      
    1037                else 
    1038                   fsin(ji,jj)  = 0.0 
    1039                   fnsi(ji,jj)  = 0.0 
    1040                   fprd(ji,jj)  = 0.0 
    1041                   fsld2(ji,jj) = 0.0 
    1042                   fprds(ji,jj) = 0.0 
    1043                endif 
    1044  
    1045 # if defined key_debug_medusa 
    1046                !! report phytoplankton growth (including diatom silicon submodel) 
    1047                if (idf.eq.1.AND.idfval.eq.1) then 
    1048                   IF (lwp) write (numout,*) '------------------------------' 
    1049                   IF (lwp) write (numout,*) 'fsin(',jk,')   = ', fsin(ji,jj) 
    1050                   IF (lwp) write (numout,*) 'fnsi(',jk,')   = ', fnsi(ji,jj) 
    1051                   IF (lwp) write (numout,*) 'fsld2(',jk,')  = ', fsld2(ji,jj) 
    1052                   IF (lwp) write (numout,*) 'fprn(',jk,')   = ', fprn(ji,jj) 
    1053                   IF (lwp) write (numout,*) 'fprd(',jk,')   = ', fprd(ji,jj) 
    1054                   IF (lwp) write (numout,*) 'fprds(',jk,')  = ', fprds(ji,jj) 
    1055                endif 
    1056 # endif 
    1057  
    1058                !!---------------------------------------------------------------------- 
    1059                !! Mixed layer primary production 
    1060                !! this block calculates the amount of primary production that occurs 
    1061                !! within the upper mixed layer; this allows the separate diagnosis 
    1062                !! of "sub-surface" primary production; it does assume that short- 
    1063                !! term variability in mixed layer depth doesn't mess with things 
    1064                !! though 
    1065                !!---------------------------------------------------------------------- 
    1066                !! 
    1067                if (fdep1(ji,jj).le.hmld(ji,jj)) then 
    1068                   !! this level is entirely in the mixed layer 
    1069                   fq0 = 1.0 
    1070                elseif (fsdepw(ji,jj,jk).ge.hmld(ji,jj)) then 
    1071                   !! this level is entirely below the mixed layer 
    1072                   fq0 = 0.0 
    1073                else 
    1074                   !! this level straddles the mixed layer 
    1075                   fq0 = (hmld(ji,jj) - fsdepw(ji,jj,jk)) / fse3t(ji,jj,jk) 
    1076                endif 
    1077                !! 
    1078                fprn_ml(ji,jj) = fprn_ml(ji,jj) + (fprn(ji,jj) * zphn(ji,jj) * fse3t(ji,jj,jk) * fq0) 
    1079                fprd_ml(ji,jj) = fprd_ml(ji,jj) + (fprd(ji,jj) * zphd(ji,jj) * fse3t(ji,jj,jk) * fq0) 
    1080                 
    1081                !!---------------------------------------------------------------------- 
    1082                !! Vertical Integral -- 
    1083                !!---------------------------------------------------------------------- 
    1084                ftot_pn(ji,jj)  = ftot_pn(ji,jj)  + (zphn(ji,jj) * fse3t(ji,jj,jk))   !! vertical integral non-diatom phytoplankton 
    1085                ftot_pd(ji,jj)  = ftot_pd(ji,jj)  + (zphd(ji,jj) * fse3t(ji,jj,jk))   !! vertical integral diatom phytoplankton 
    1086                ftot_zmi(ji,jj) = ftot_zmi(ji,jj) + (zzmi(ji,jj) * fse3t(ji,jj,jk))   !! vertical integral microzooplankton 
    1087                ftot_zme(ji,jj) = ftot_zme(ji,jj) + (zzme(ji,jj) * fse3t(ji,jj,jk))   !! vertical integral mesozooplankton 
    1088                ftot_det(ji,jj) = ftot_det(ji,jj) + (zdet(ji,jj) * fse3t(ji,jj,jk))   !! vertical integral slow detritus, nitrogen 
    1089                ftot_dtc(ji,jj) = ftot_dtc(ji,jj) + (zdtc(ji,jj) * fse3t(ji,jj,jk))   !! vertical integral slow detritus, carbon 
    1090                 
    1091                !!---------------------------------------------------------------------- 
    1092                !! More chlorophyll calculations 
    1093                !!---------------------------------------------------------------------- 
    1094                !! 
    1095                !! frn(ji,jj) = (xthetam / fthetan(ji,jj)) * (fprn(ji,jj) / (fthetan(ji,jj) * xpar(ji,jj,jk))) 
    1096                !! frd(ji,jj) = (xthetam / fthetad(ji,jj)) * (fprd(ji,jj) / (fthetad(ji,jj) * xpar(ji,jj,jk))) 
    1097                frn(ji,jj) = (xthetam * fchn(ji,jj) * fnln(ji,jj) * ffln(ji,jj)       ) / (fthetan(ji,jj) + tiny(fthetan(ji,jj))) 
    1098                !! AXY (12/05/09): there's potentially a problem here; fsld, silicic acid  
    1099                !!   limitation, is used in the following line to regulate chlorophyll  
    1100                !!   growth in a manner that is inconsistent with its use in the regulation  
    1101                !!   of biomass growth; the Mongin term term used in growth is more complex 
    1102                !!   than the simple multiplicative function used below 
    1103                !! frd(ji,jj) = (xthetam * fchd(ji,jj) * fnld(ji,jj) * ffld(ji,jj) * fsld(ji,jj)) / (fthetad(ji,jj) + tiny(fthetad(ji,jj))) 
    1104                !! AXY (12/05/09): this replacement line uses the new variable, fsld2, to 
    1105                !!   regulate chlorophyll growth 
    1106                frd(ji,jj) = (xthetamd * fchd(ji,jj) * fnld(ji,jj) * ffld(ji,jj) * fsld2(ji,jj)) / (fthetad(ji,jj) + tiny(fthetad(ji,jj))) 
    1107  
    1108 # if defined key_debug_medusa 
    1109                !! report chlorophyll calculations 
    1110                if (idf.eq.1.AND.idfval.eq.1) then 
    1111                   IF (lwp) write (numout,*) '------------------------------' 
    1112                   IF (lwp) write (numout,*) 'fthetan(',jk,') = ', fthetan(ji,jj) 
    1113                   IF (lwp) write (numout,*) 'fthetad(',jk,') = ', fthetad(ji,jj) 
    1114                   IF (lwp) write (numout,*) 'frn(',jk,')     = ', frn(ji,jj) 
    1115                   IF (lwp) write (numout,*) 'frd(',jk,')     = ', frd(ji,jj) 
    1116                endif 
    1117 # endif 
    1118  
    1119 ! MAYBE BUT A BREAK IN HERE, ZOOPLANKTON GRAZING - marc 20/4/17  
    1120 ! (plankton growth is 281 lines) 
    1121  
    1122                !!---------------------------------------------------------------------- 
    1123                !! Zooplankton Grazing  
    1124                !! this code supplements the base grazing model with one that 
    1125                !! considers the C:N ratio of grazed food and balances this against 
    1126                !! the requirements of zooplankton growth; this model is derived  
    1127                !! from that of Anderson & Pondaven (2003) 
    1128                !! 
    1129                !! the current version of the code assumes a fixed C:N ratio for 
    1130                !! detritus (in contrast to Anderson & Pondaven, 2003), though the 
    1131                !! full equations are retained for future extension 
    1132                !!---------------------------------------------------------------------- 
    1133                !! 
    1134                !!---------------------------------------------------------------------- 
    1135                !! Microzooplankton first 
    1136                !!---------------------------------------------------------------------- 
    1137                !! 
    1138                fmi1(ji,jj)    = (xkmi * xkmi) + (xpmipn * zphn(ji,jj) * zphn(ji,jj)) + (xpmid * zdet(ji,jj) * zdet(ji,jj)) 
    1139                fmi(ji,jj)     = xgmi * zzmi(ji,jj) / fmi1(ji,jj) 
    1140                fgmipn(ji,jj)  = fmi(ji,jj) * xpmipn * zphn(ji,jj) * zphn(ji,jj)   !! grazing on non-diatoms 
    1141                fgmid(ji,jj)   = fmi(ji,jj) * xpmid  * zdet(ji,jj) * zdet(ji,jj)   !! grazing on detrital nitrogen 
    1142 # if defined key_roam 
    1143                fgmidc(ji,jj)  = rsmall !acc 
    1144                IF ( zdet(ji,jj) .GT. rsmall ) fgmidc(ji,jj)  = (zdtc(ji,jj) / (zdet(ji,jj) + tiny(zdet(ji,jj)))) * fgmid(ji,jj)  !! grazing on detrital carbon 
    1145 # else 
    1146                !! AXY (26/11/08): implicit detrital carbon change 
    1147                fgmidc(ji,jj)  = xthetad * fgmid(ji,jj)              !! grazing on detrital carbon 
    1148 # endif 
    1149                !! 
    1150                !! which translates to these incoming N and C fluxes 
    1151                finmi(ji,jj)   = (1.0 - xphi) * (fgmipn(ji,jj) + fgmid(ji,jj)) 
    1152                ficmi(ji,jj)   = (1.0 - xphi) * ((xthetapn * fgmipn(ji,jj)) + fgmidc(ji,jj)) 
    1153                !! 
    1154                !! the ideal food C:N ratio for microzooplankton 
    1155                !! xbetan = 0.77; xthetaz = 5.625; xbetac = 0.64; xkc = 0.80 
    1156                fstarmi(ji,jj) = (xbetan * xthetazmi) / (xbetac * xkc) 
    1157                !! 
    1158                !! process these to determine proportioning of grazed N and C 
    1159                !! (since there is no explicit consideration of respiration, 
    1160                !! only growth and excretion are calculated here) 
    1161                fmith(ji,jj)   = (ficmi(ji,jj) / (finmi(ji,jj) + tiny(finmi(ji,jj)))) 
    1162                if (fmith(ji,jj).ge.fstarmi(ji,jj)) then 
    1163                   fmigrow(ji,jj) = xbetan * finmi(ji,jj) 
    1164                   fmiexcr(ji,jj) = 0.0 
    1165                else 
    1166                   fmigrow(ji,jj) = (xbetac * xkc * ficmi(ji,jj)) / xthetazmi 
    1167                   fmiexcr(ji,jj) = ficmi(ji,jj) * ((xbetan / (fmith(ji,jj) + tiny(fmith(ji,jj)))) - ((xbetac * xkc) / xthetazmi)) 
    1168                endif 
    1169 # if defined key_roam 
    1170                fmiresp(ji,jj) = (xbetac * ficmi(ji,jj)) - (xthetazmi * fmigrow(ji,jj)) 
    1171 # endif 
    1172  
    1173 # if defined key_debug_medusa 
    1174                !! report microzooplankton grazing 
    1175                if (idf.eq.1.AND.idfval.eq.1) then 
    1176                   IF (lwp) write (numout,*) '------------------------------' 
    1177                   IF (lwp) write (numout,*) 'fmi1(',jk,')    = ', fmi1(ji,jj) 
    1178                   IF (lwp) write (numout,*) 'fmi(',jk,')     = ', fmi(ji,jj) 
    1179                   IF (lwp) write (numout,*) 'fgmipn(',jk,')  = ', fgmipn(ji,jj) 
    1180                   IF (lwp) write (numout,*) 'fgmid(',jk,')   = ', fgmid(ji,jj) 
    1181                   IF (lwp) write (numout,*) 'fgmidc(',jk,')  = ', fgmidc(ji,jj) 
    1182                   IF (lwp) write (numout,*) 'finmi(',jk,')   = ', finmi(ji,jj) 
    1183                   IF (lwp) write (numout,*) 'ficmi(',jk,')   = ', ficmi(ji,jj) 
    1184                   IF (lwp) write (numout,*) 'fstarmi(',jk,') = ', fstarmi(ji,jj) 
    1185                   IF (lwp) write (numout,*) 'fmith(',jk,')   = ', fmith(ji,jj) 
    1186                   IF (lwp) write (numout,*) 'fmigrow(',jk,') = ', fmigrow(ji,jj) 
    1187                   IF (lwp) write (numout,*) 'fmiexcr(',jk,') = ', fmiexcr(ji,jj) 
    1188 #  if defined key_roam 
    1189                   IF (lwp) write (numout,*) 'fmiresp(',jk,') = ', fmiresp(ji,jj) 
    1190 #  endif 
    1191                endif 
    1192 # endif 
    1193  
    1194                !!---------------------------------------------------------------------- 
    1195                !! Mesozooplankton second 
    1196                !!---------------------------------------------------------------------- 
    1197                !! 
    1198                fme1(ji,jj)    = (xkme * xkme) + (xpmepn * zphn(ji,jj) * zphn(ji,jj)) + (xpmepd * zphd(ji,jj) * zphd(ji,jj)) + &  
    1199                          (xpmezmi * zzmi(ji,jj) * zzmi(ji,jj)) + (xpmed * zdet(ji,jj) * zdet(ji,jj)) 
    1200                fme(ji,jj)     = xgme * zzme(ji,jj) / fme1(ji,jj) 
    1201                fgmepn(ji,jj)  = fme(ji,jj) * xpmepn  * zphn(ji,jj) * zphn(ji,jj)  !! grazing on non-diatoms 
    1202                fgmepd(ji,jj)  = fme(ji,jj) * xpmepd  * zphd(ji,jj) * zphd(ji,jj)  !! grazing on diatoms 
    1203                fgmepds(ji,jj) = fsin(ji,jj) * fgmepd(ji,jj)                !! grazing on diatom silicon 
    1204                fgmezmi(ji,jj) = fme(ji,jj) * xpmezmi * zzmi(ji,jj) * zzmi(ji,jj)  !! grazing on microzooplankton 
    1205                fgmed(ji,jj)   = fme(ji,jj) * xpmed   * zdet(ji,jj) * zdet(ji,jj)  !! grazing on detrital nitrogen 
    1206 # if defined key_roam 
    1207                fgmedc(ji,jj)  = rsmall !acc 
    1208                IF ( zdet(ji,jj) .GT. rsmall ) fgmedc(ji,jj)  = (zdtc(ji,jj) / (zdet(ji,jj) + tiny(zdet(ji,jj)))) * fgmed(ji,jj)  !! grazing on detrital carbon 
    1209 # else 
    1210                !! AXY (26/11/08): implicit detrital carbon change 
    1211                fgmedc(ji,jj)  = xthetad * fgmed(ji,jj)              !! grazing on detrital carbon 
    1212 # endif 
    1213                !! 
    1214                !! which translates to these incoming N and C fluxes 
    1215                finme(ji,jj)   = (1.0 - xphi) * (fgmepn(ji,jj) + fgmepd(ji,jj) + fgmezmi(ji,jj) + fgmed(ji,jj)) 
    1216                ficme(ji,jj)   = (1.0 - xphi) * ((xthetapn * fgmepn(ji,jj)) + (xthetapd * fgmepd(ji,jj)) + & 
    1217                         (xthetazmi * fgmezmi(ji,jj)) + fgmedc(ji,jj)) 
    1218                !! 
    1219                !! the ideal food C:N ratio for mesozooplankton 
    1220                !! xbetan = 0.77; xthetaz = 5.625; xbetac = 0.64; xkc = 0.80 
    1221                fstarme(ji,jj) = (xbetan * xthetazme) / (xbetac * xkc) 
    1222                !! 
    1223                !! process these to determine proportioning of grazed N and C 
    1224                !! (since there is no explicit consideration of respiration, 
    1225                !! only growth and excretion are calculated here) 
    1226                fmeth(ji,jj)   = (ficme(ji,jj) / (finme(ji,jj) + tiny(finme(ji,jj)))) 
    1227                if (fmeth(ji,jj).ge.fstarme(ji,jj)) then 
    1228                   fmegrow(ji,jj) = xbetan * finme(ji,jj) 
    1229                   fmeexcr(ji,jj) = 0.0 
    1230                else 
    1231                   fmegrow(ji,jj) = (xbetac * xkc * ficme(ji,jj)) / xthetazme 
    1232                   fmeexcr(ji,jj) = ficme(ji,jj) * ((xbetan / (fmeth(ji,jj) + tiny(fmeth(ji,jj)))) - ((xbetac * xkc) / xthetazme)) 
    1233                endif 
    1234 # if defined key_roam 
    1235                fmeresp(ji,jj) = (xbetac * ficme(ji,jj)) - (xthetazme * fmegrow(ji,jj)) 
    1236 # endif 
    1237  
    1238 # if defined key_debug_medusa 
    1239                !! report mesozooplankton grazing 
    1240                if (idf.eq.1.AND.idfval.eq.1) then 
    1241                   IF (lwp) write (numout,*) '------------------------------' 
    1242                   IF (lwp) write (numout,*) 'fme1(',jk,')    = ', fme1(ji,jj) 
    1243                   IF (lwp) write (numout,*) 'fme(',jk,')     = ', fme(ji,jj) 
    1244                   IF (lwp) write (numout,*) 'fgmepn(',jk,')  = ', fgmepn(ji,jj) 
    1245                   IF (lwp) write (numout,*) 'fgmepd(',jk,')  = ', fgmepd(ji,jj) 
    1246                   IF (lwp) write (numout,*) 'fgmepds(',jk,') = ', fgmepds(ji,jj) 
    1247                   IF (lwp) write (numout,*) 'fgmezmi(',jk,') = ', fgmezmi(ji,jj) 
    1248                   IF (lwp) write (numout,*) 'fgmed(',jk,')   = ', fgmed(ji,jj) 
    1249                   IF (lwp) write (numout,*) 'fgmedc(',jk,')  = ', fgmedc(ji,jj) 
    1250                   IF (lwp) write (numout,*) 'finme(',jk,')   = ', finme(ji,jj) 
    1251                   IF (lwp) write (numout,*) 'ficme(',jk,')   = ', ficme(ji,jj) 
    1252                   IF (lwp) write (numout,*) 'fstarme(',jk,') = ', fstarme(ji,jj) 
    1253                   IF (lwp) write (numout,*) 'fmeth(',jk,')   = ', fmeth(ji,jj) 
    1254                   IF (lwp) write (numout,*) 'fmegrow(',jk,') = ', fmegrow(ji,jj) 
    1255                   IF (lwp) write (numout,*) 'fmeexcr(',jk,') = ', fmeexcr(ji,jj) 
    1256 #  if defined key_roam 
    1257                   IF (lwp) write (numout,*) 'fmeresp(',jk,') = ', fmeresp(ji,jj) 
    1258 #  endif 
    1259                endif 
    1260 # endif 
    1261  
    1262                fzmi_i(ji,jj)  = fzmi_i(ji,jj)  + fse3t(ji,jj,jk) * (  & 
    1263                   fgmipn(ji,jj) + fgmid(ji,jj) ) 
    1264                fzmi_o(ji,jj)  = fzmi_o(ji,jj)  + fse3t(ji,jj,jk) * (  & 
    1265                   fmigrow(ji,jj) + (xphi * (fgmipn(ji,jj) + fgmid(ji,jj))) + fmiexcr(ji,jj) + ((1.0 - xbetan) * finmi(ji,jj)) ) 
    1266                fzme_i(ji,jj)  = fzme_i(ji,jj)  + fse3t(ji,jj,jk) * (  & 
    1267                   fgmepn(ji,jj) + fgmepd(ji,jj) + fgmezmi(ji,jj) + fgmed(ji,jj) ) 
    1268                fzme_o(ji,jj)  = fzme_o(ji,jj)  + fse3t(ji,jj,jk) * (  & 
    1269                   fmegrow(ji,jj) + (xphi * (fgmepn(ji,jj) + fgmepd(ji,jj) + fgmezmi(ji,jj) + fgmed(ji,jj))) + fmeexcr(ji,jj) + ((1.0 - xbetan) * finme(ji,jj)) ) 
    1270  
    1271 ! MAYBE BUT A BREAK IN HERE, MISCELLANEOUS PLANKTON LOSSES - marc 20/4/17  
    1272 ! (zoo plankton grazing is 152 lines) 
    1273  
    1274                !!---------------------------------------------------------------------- 
    1275                !! Plankton metabolic losses 
    1276                !! Linear loss processes assumed to be metabolic in origin 
    1277                !!---------------------------------------------------------------------- 
    1278                !! 
    1279                fdpn2(ji,jj)  = xmetapn  * zphn(ji,jj) 
    1280                fdpd2(ji,jj)  = xmetapd  * zphd(ji,jj) 
    1281                fdpds2(ji,jj) = xmetapd  * zpds(ji,jj) 
    1282                fdzmi2(ji,jj) = xmetazmi * zzmi(ji,jj) 
    1283                fdzme2(ji,jj) = xmetazme * zzme(ji,jj) 
    1284  
    1285                !!---------------------------------------------------------------------- 
    1286                !! Plankton mortality losses 
    1287                !! EKP (26/02/09): phytoplankton hyperbolic mortality term introduced  
    1288                !! to improve performance in gyres 
    1289                !!---------------------------------------------------------------------- 
    1290                !! 
    1291                !! non-diatom phytoplankton 
    1292                if (jmpn.eq.1) fdpn(ji,jj) = xmpn * zphn(ji,jj)               !! linear 
    1293                if (jmpn.eq.2) fdpn(ji,jj) = xmpn * zphn(ji,jj) * zphn(ji,jj)        !! quadratic 
    1294                if (jmpn.eq.3) fdpn(ji,jj) = xmpn * zphn(ji,jj) * &           !! hyperbolic 
    1295                   (zphn(ji,jj) / (xkphn + zphn(ji,jj))) 
    1296                if (jmpn.eq.4) fdpn(ji,jj) = xmpn * zphn(ji,jj) * &           !! sigmoid 
    1297                   ((zphn(ji,jj) * zphn(ji,jj)) / (xkphn + (zphn(ji,jj) * zphn(ji,jj)))) 
    1298                !! 
    1299                !! diatom phytoplankton 
    1300                if (jmpd.eq.1) fdpd(ji,jj) = xmpd * zphd(ji,jj)               !! linear 
    1301                if (jmpd.eq.2) fdpd(ji,jj) = xmpd * zphd(ji,jj) * zphd(ji,jj)        !! quadratic 
    1302                if (jmpd.eq.3) fdpd(ji,jj) = xmpd * zphd(ji,jj) * &           !! hyperbolic 
    1303                   (zphd(ji,jj) / (xkphd + zphd(ji,jj))) 
    1304                if (jmpd.eq.4) fdpd(ji,jj) = xmpd * zphd(ji,jj) * &           !! sigmoid 
    1305                   ((zphd(ji,jj) * zphd(ji,jj)) / (xkphd + (zphd(ji,jj) * zphd(ji,jj)))) 
    1306                fdpds(ji,jj) = fdpd(ji,jj) * fsin(ji,jj) 
    1307                !! 
    1308                !! microzooplankton 
    1309                if (jmzmi.eq.1) fdzmi(ji,jj) = xmzmi * zzmi(ji,jj)            !! linear 
    1310                if (jmzmi.eq.2) fdzmi(ji,jj) = xmzmi * zzmi(ji,jj) * zzmi(ji,jj)     !! quadratic 
    1311                if (jmzmi.eq.3) fdzmi(ji,jj) = xmzmi * zzmi(ji,jj) * &        !! hyperbolic 
    1312                   (zzmi(ji,jj) / (xkzmi + zzmi(ji,jj))) 
    1313                if (jmzmi.eq.4) fdzmi(ji,jj) = xmzmi * zzmi(ji,jj) * &        !! sigmoid 
    1314                   ((zzmi(ji,jj) * zzmi(ji,jj)) / (xkzmi + (zzmi(ji,jj) * zzmi(ji,jj)))) 
    1315                !! 
    1316                !! mesozooplankton 
    1317                if (jmzme.eq.1) fdzme(ji,jj) = xmzme * zzme(ji,jj)            !! linear 
    1318                if (jmzme.eq.2) fdzme(ji,jj) = xmzme * zzme(ji,jj) * zzme(ji,jj)     !! quadratic 
    1319                if (jmzme.eq.3) fdzme(ji,jj) = xmzme * zzme(ji,jj) * &        !! hyperbolic 
    1320                   (zzme(ji,jj) / (xkzme + zzme(ji,jj))) 
    1321                if (jmzme.eq.4) fdzme(ji,jj) = xmzme * zzme(ji,jj) * &        !! sigmoid 
    1322                   ((zzme(ji,jj) * zzme(ji,jj)) / (xkzme + (zzme(ji,jj) * zzme(ji,jj)))) 
    1323  
    1324 ! MAYBE BUT A BREAK IN HERE, DETRITUS PROCESS- marc 20/4/17  
    1325 ! (misc plankton loses is 53 lines - maybe we can link with another section) 
    1326859 
    1327860               !!---------------------------------------------------------------------- 
     
    1400933! MAYBE BUT A BREAK IN HERE, IRON CHEMISTRY AND SCAVENGING - marc 20/4/17  
    1401934! (detritus processes is 74 lines) 
     935            ENDIF 
     936         ENDDO 
     937         ENDDO 
     938 
     939         DO jj = 2,jpjm1 
     940         DO ji = 2,jpim1 
     941            !! OPEN wet point IF..THEN loop 
     942            if (tmask(ji,jj,jk) == 1) then 
    1402943 
    1403944               !!---------------------------------------------------------------------- 
     
    17381279! MAYBE BUT A BREAK IN HERE, MISCELLANEOUS PROCESSES - marc 20/4/17  
    17391280! (iron chemistry and scavenging is 340 lines) 
     1281            ENDIF 
     1282         ENDDO 
     1283         ENDDO 
     1284 
     1285         DO jj = 2,jpjm1 
     1286         DO ji = 2,jpim1 
     1287            !! OPEN wet point IF..THEN loop 
     1288            if (tmask(ji,jj,jk) == 1) then 
    17401289 
    17411290               !!---------------------------------------------------------------------- 
     
    18171366! MAYBE BUT A BREAK IN HERE, FAST-SINKINIG DETRITUS - marc 20/4/17  
    18181367! (miscellaneous processes is 79 lines) 
     1368            ENDIF 
     1369         ENDDO 
     1370         ENDDO 
     1371 
     1372         DO jj = 2,jpjm1 
     1373         DO ji = 2,jpim1 
     1374            !! OPEN wet point IF..THEN loop 
     1375            if (tmask(ji,jj,jk) == 1) then 
    18191376 
    18201377               !!---------------------------------------------------------------------- 
     
    23941951! MAYBE BUT A BREAK IN HERE, BUSINESS AND UPDATING - marc 20/4/17  
    23951952! (fast sinking detritus is 576 lines) 
     1953            ENDIF 
     1954         ENDDO 
     1955         ENDDO 
     1956 
     1957         DO jj = 2,jpjm1 
     1958         DO ji = 2,jpim1 
     1959            !! OPEN wet point IF..THEN loop 
     1960            if (tmask(ji,jj,jk) == 1) then 
    23961961 
    23971962               !!====================================================================== 
     
    28402405! (this would make the previous section about 470 lines and the one below 
    28412406! about 700 lines) 
     2407            ENDIF 
     2408         ENDDO 
     2409         ENDDO 
     2410 
     2411         DO jj = 2,jpjm1 
     2412         DO ji = 2,jpim1 
     2413            !! OPEN wet point IF..THEN loop 
     2414            if (tmask(ji,jj,jk) == 1) then 
    28422415 
    28432416# if defined key_trc_diabio 
     
    29392512                  ENDIF 
    29402513                  IF( med_diag%PN_FELIM%dgsave ) THEN 
    2941                       ffln2d(ji,jj) = ffln2d(ji,jj) + (ffln(ji,jj)  * zphn(ji,jj) * fse3t(ji,jj,jk))  
     2514                      ffln2d(ji,jj) = ffln2d(ji,jj) + (ffln2(ji,jj)  * zphn(ji,jj) * fse3t(ji,jj,jk))  
    29422515                  ENDIF 
    29432516                  IF( med_diag%PD_JLIM%dgsave ) THEN 
     
    32262799                  ENDIF 
    32272800        IF ( med_diag%PNLIMFE3%dgsave ) THEN 
    3228                      pnlimfe3(ji,jj,jk) = ffln(ji,jj) 
     2801                     pnlimfe3(ji,jj,jk) = ffln2(ji,jj) 
    32292802                  ENDIF 
    32302803        IF ( med_diag%PDLIMJ3%dgsave  ) THEN 
     
    32812854                  trc2d(ji,jj,25) = trc2d(ji,jj,25) + (fjlim_pn(ji,jj) * zphn(ji,jj) * fse3t(ji,jj,jk)) !! non-diatom J  limitation term  
    32822855                  trc2d(ji,jj,26) = trc2d(ji,jj,26) + (fnln(ji,jj)  * zphn(ji,jj) * fse3t(ji,jj,jk))    !! non-diatom N  limitation term  
    3283                   trc2d(ji,jj,27) = trc2d(ji,jj,27) + (ffln(ji,jj)  * zphn(ji,jj) * fse3t(ji,jj,jk))    !! non-diatom Fe limitation term  
     2856                  trc2d(ji,jj,27) = trc2d(ji,jj,27) + (ffln2(ji,jj)  * zphn(ji,jj) * fse3t(ji,jj,jk))    !! non-diatom Fe limitation term  
    32842857                  trc2d(ji,jj,28) = trc2d(ji,jj,28) + (fjlim_pd(ji,jj) * zphd(ji,jj) * fse3t(ji,jj,jk)) !! diatom     J  limitation term  
    32852858                  trc2d(ji,jj,29) = trc2d(ji,jj,29) + (fnld(ji,jj)  * zphd(ji,jj) * fse3t(ji,jj,jk))    !! diatom     N  limitation term  
Note: See TracChangeset for help on using the changeset viewer.