New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
detritus_fast_sink.F90 in branches/UKMO/dev_r5518_medusa_chg_trc_bio_medusa/NEMOGCM/NEMO/TOP_SRC/MEDUSA – NEMO

source: branches/UKMO/dev_r5518_medusa_chg_trc_bio_medusa/NEMOGCM/NEMO/TOP_SRC/MEDUSA/detritus_fast_sink.F90 @ 7986

Last change on this file since 7986 was 7986, checked in by marc, 7 years ago

Pulled detritus processes out of trcbio_medusa.F90

File size: 45.4 KB
Line 
1MODULE detritus_fast_sink_mod
2   !!======================================================================
3   !!                         ***  MODULE detritus_fast_sink_mod  ***
4   !! Calculates fast-sinking detritus
5   !!======================================================================
6   !! History :
7   !!   -   ! 2017-04 (M. Stringer)        Code taken from trcbio_medusa.F90
8   !!----------------------------------------------------------------------
9#if defined key_medusa
10   !!----------------------------------------------------------------------
11   !!                                                   MEDUSA bio-model
12   !!----------------------------------------------------------------------
13
14   IMPLICIT NONE
15   PRIVATE
16     
17   PUBLIC   detritus_fast_sink        ! Called in detritus.F90
18
19   !!----------------------------------------------------------------------
20   !! NEMO/TOP 2.0 , LOCEAN-IPSL (2007)
21   !! $Id$
22   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
23   !!----------------------------------------------------------------------
24
25CONTAINS
26
27   SUBROUTINE detritus_fast_sink( jk, iball )
28      !!-------------------------------------------------------------------
29      !!                     ***  ROUTINE detritus_fast_sink  ***
30      !! This called from DETRITUS and calculates the fast-sinking detritus
31      !!-------------------------------------------------------------------
32      USE bio_medusa_mod,    ONLY: b0,                                     &
33                                   f_benout_c, f_benout_ca, f_benout_fe,   &
34                                   f_benout_lyso_ca, f_benout_n,           &
35                                   f_benout_si,                            &
36                                   f_fbenin_c, f_fbenin_ca, f_fbenin_fe,   &
37                                   f_fbenin_n, f_fbenin_si,                &
38                                   f_omcal, fcaco3,                        &
39                                   fccd, fccd_dep, fdep1,                  &
40                                   fdpd, fdpds, fdzme, fgmepds,            &
41                                   ffast2slowc, ffast2slown,               &
42                                   ffastc, ffastca, ffastfe, ffastn,       &
43                                   ffastsi,                                &
44                                   fifd_c, fifd_fe, fifd_n, fifd_si,       &
45                                   fofd_fe, fofd_n, fofd_si,               &
46                                   fprotf,                                 &
47                                   fregenfast, fregenfastsi,               &
48                                   freminc, freminca, freminfe,            &
49                                   freminn, freminsi,                      &
50                                   fsedc, fsedca, fsedn, fsedfe, fsedsi,   &
51                                   fslowc, fslown,                         &
52                                   ftempc, ftempca, ftempfe, ftempn,       &
53                                   ftempsi,                                &
54# if defined key_roam
55                                   fofd_c, fregenfastc,                    &
56# endif
57                                   idf, idfval
58      USE dom_oce,           ONLY: e3t_0, e3t_n, gdepw_0, gdepw_n,         &
59                                   gphit, mbathy, tmask
60      USE in_out_manager,    ONLY: lwp, numout
61      USE oce,               ONLY: tsn
62      USE par_kind,          ONLY: wp
63      USE par_oce,           ONLY: jpim1, jpjm1
64      USE sms_medusa,        ONLY: f2_ccd_cal, f3_omcal,                   &
65                                   jexport, jfdfate, jinorgben, jocalccd,  &
66                                   jorgben, jp_tem, jrratio,               &
67                                   ocal_ccd, xcaco3a, xcaco3b,             &
68                                   xfastc, xfastca, xfastsi,               &
69                                   xfdfrac1, xfdfrac2, xfdfrac3,           &
70                                   xmassc, xmassca, xmasssi,               &
71                                   xprotca, xprotsi, xrfn, xridg_r0,       &
72                                   xsedc, xsedca, xsedfe,xsedn, xsedsi,    &
73                                   xthetapd, xthetazme,                    &
74                                   zn_sed_c, zn_sed_ca, zn_sed_fe,         &
75                                   zn_sed_n, zn_sed_si
76
77   !!* Substitution
78#  include "domzgr_substitute.h90"
79
80      !! Level
81      INTEGER, INTENT( in ) :: jk
82      !! Fast detritus ballast scheme (0 = no; 1 = yes)
83      INTEGER, INTENT( in ) :: iball
84
85      !! Loop variables
86      INTEGER :: ji, jj
87
88      REAL(wp) ::    fb_val, fl_sst
89      !! temporary variables
90      REAL(wp) :: fq0,fq1,fq2,fq3,fq4,fq5,fq6,fq7,fq8
91
92      !!-------------------------------------------------------------------
93      !! Fast-sinking detritus terms
94      !! "local" variables declared so that conservation can be checked;
95      !! the calculated terms are added to the fast-sinking flux later on
96      !! only after the flux entering this level has experienced some
97      !! remineralisation
98      !! note: these fluxes need to be scaled by the level thickness
99      !!-------------------------------------------------------------------
100      DO jj = 2,jpjm1
101         DO ji = 2,jpim1
102            !! OPEN wet point IF..THEN loop
103            if (tmask(ji,jj,jk) == 1) then
104
105               !! nitrogen:   diatom and mesozooplankton mortality
106               ftempn(ji,jj)  = b0 * ((xfdfrac1 * fdpd(ji,jj))  +            &
107                                      (xfdfrac2 * fdzme(ji,jj)))
108               !!
109               !! silicon:    diatom mortality and grazed diatoms
110               ftempsi(ji,jj) = b0 * ((xfdfrac1 * fdpds(ji,jj)) +            &
111                                      (xfdfrac3 * fgmepds(ji,jj)))
112               !!
113               !! iron:       diatom and mesozooplankton mortality
114               ftempfe(ji,jj) = b0 * (((xfdfrac1 * fdpd(ji,jj)) +            &
115                                       (xfdfrac2 * fdzme(ji,jj))) * xrfn)
116               !!
117               !! carbon:     diatom and mesozooplankton mortality
118               ftempc(ji,jj)  = b0 * ((xfdfrac1 * xthetapd * fdpd(ji,jj)) +  & 
119                                      (xfdfrac2 * xthetazme * fdzme(ji,jj)))
120               !!
121            ENDIF
122         ENDDO
123      ENDDO
124
125# if defined key_roam
126      DO jj = 2,jpjm1
127         DO ji = 2,jpim1
128            if (tmask(ji,jj,jk) == 1) then
129               if (jrratio.eq.0) then
130                  !! CaCO3:      latitudinally-based fraction of total
131                  !!             primary production
132                  !!               0.10 at equator; 0.02 at pole
133                  fcaco3(ji,jj) = xcaco3a + ((xcaco3b - xcaco3a) *           &
134                                             ((90.0 - abs(gphit(ji,jj))) /   &
135                                               90.0))
136               elseif (jrratio.eq.1) then
137                  !! CaCO3:      Ridgwell et al. (2007) submodel, version 1
138                  !!             this uses SURFACE omega calcite to regulate
139                  !!             rain ratio
140                  if (f_omcal(ji,jj).ge.1.0) then
141                     fq1 = (f_omcal(ji,jj) - 1.0)**0.81
142                  else
143                     fq1 = 0.
144                  endif
145                  fcaco3(ji,jj) = xridg_r0 * fq1
146               elseif (jrratio.eq.2) then
147                  !! CaCO3:      Ridgwell et al. (2007) submodel, version 2
148                  !!             this uses FULL 3D omega calcite to regulate
149                  !!             rain ratio
150                  if (f3_omcal(ji,jj,jk).ge.1.0) then
151                     fq1 = (f3_omcal(ji,jj,jk) - 1.0)**0.81
152                  else
153                     fq1 = 0.
154                  endif
155                  fcaco3(ji,jj) = xridg_r0 * fq1
156               endif
157            ENDIF
158         ENDDO
159      ENDDO
160# else
161      DO jj = 2,jpjm1
162         DO ji = 2,jpim1
163            if (tmask(ji,jj,jk) == 1) then
164               !! CaCO3:      latitudinally-based fraction of total primary
165               !!              production
166               !!               0.10 at equator; 0.02 at pole
167               fcaco3(ji,jj) = xcaco3a + ((xcaco3b - xcaco3a) *              &
168                                          ((90.0 - abs(gphit(ji,jj))) / 90.0))
169            ENDIF
170         ENDDO
171      ENDDO
172# endif
173
174      DO jj = 2,jpjm1
175         DO ji = 2,jpim1
176            if (tmask(ji,jj,jk) == 1) then
177               !! AXY (09/03/09): convert CaCO3 production from function of
178               !! primary production into a function of fast-sinking material;
179               !! technically, this is what Dunne et al. (2007) do anyway; they
180               !! convert total primary production estimated from surface
181               !! chlorophyll to an export flux for which they apply conversion
182               !! factors to estimate the various elemental fractions (Si, Ca)
183               ftempca(ji,jj) = ftempc(ji,jj) * fcaco3(ji,jj)
184
185# if defined key_debug_medusa
186               !! integrate total fast detritus production
187               if (idf.eq.1) then
188                  fifd_n(ji,jj)  = fifd_n(ji,jj)  + (ftempn(ji,jj)  *         &
189                                                     fse3t(ji,jj,jk))
190                  fifd_si(ji,jj) = fifd_si(ji,jj) + (ftempsi(ji,jj) *         &
191                                                     fse3t(ji,jj,jk))
192                  fifd_fe(ji,jj) = fifd_fe(ji,jj) + (ftempfe(ji,jj) *         &
193                                                     fse3t(ji,jj,jk))
194#  if defined key_roam
195                  fifd_c(ji,jj)  = fifd_c(ji,jj)  + (ftempc(ji,jj)  *         &
196                                                     fse3t(ji,jj,jk))
197#  endif
198               endif
199
200               !! report quantities of fast-sinking detritus for each component
201               if (idf.eq.1.AND.idfval.eq.1) then
202                  IF (lwp) write (numout,*) '------------------------------'
203! These variables are not in this routine - marc 28/4/17
204!                  IF (lwp) write (numout,*) 'fdpd(',jk,')    = ', fdpd(ji,jj)
205!                  IF (lwp) write (numout,*) 'fdzme(',jk,')   = ', fdzme(ji,jj)
206                  IF (lwp) write (numout,*) 'ftempn(',jk,')  = ', ftempn(ji,jj)
207                  IF (lwp) write (numout,*) 'ftempsi(',jk,') = ', ftempsi(ji,jj)
208                  IF (lwp) write (numout,*) 'ftempfe(',jk,') = ', ftempfe(ji,jj)
209                  IF (lwp) write (numout,*) 'ftempc(',jk,')  = ', ftempc(ji,jj)
210                  IF (lwp) write (numout,*) 'ftempca(',jk,') = ', ftempca(ji,jj)
211                  IF (lwp) write (numout,*) 'flat(',jk,')    = ',             &
212                                            abs(gphit(ji,jj))
213                  IF (lwp) write (numout,*) 'fcaco3(',jk,')  = ', fcaco3(ji,jj)
214               endif
215# endif
216            ENDIF
217         ENDDO
218      ENDDO
219
220      !!----------------------------------------------------------
221      !! This version of MEDUSA offers a choice of three methods for
222      !! handling the remineralisation of fast detritus.  All three
223      !! do so in broadly the same way:
224      !!
225      !!   1.  Fast detritus is stored as a 2D array  [ ffastX  ]
226      !!   2.  Fast detritus is added level-by-level  [ ftempX  ]
227      !!   3.  Fast detritus is not remineralised in the top box
228      !!       [ freminX ]
229      !!   4.  Remaining fast detritus is remineralised in the
230      !!       bottom  [ fsedX   ] box
231      !!
232      !! The three remineralisation methods are:
233      !!   
234      !!   1.  Ballast model (i.e. that published in Yool et al.,
235      !!       2011)
236      !!  (1b. Ballast-sans-ballast model)
237      !!   2.  Martin et al. (1987)
238      !!   3.  Henson et al. (2011)
239      !!
240      !! The first of these couples C, N and Fe remineralisation to
241      !! the remineralisation of particulate Si and CaCO3, but the
242      !! latter two treat remineralisation of C, N, Fe, Si and CaCO3
243      !! completely separately.  At present a switch within the code
244      !! regulates which submodel is used, but this should be moved
245      !! to the namelist file.
246      !!
247      !! The ballast-sans-ballast submodel is an original development
248      !! feature of MEDUSA in which the ballast submodel's general
249      !! framework and parameterisation is used, but in which there
250      !! is no protection of organic material afforded by ballasting
251      !! minerals.  While similar, it is not the same as the Martin
252      !! et al. (1987) submodel.
253      !!
254      !! Since the three submodels behave the same in terms of
255      !! accumulating sinking material and remineralising it all at
256      !! the seafloor, these portions of the code below are common to
257      !! all three.
258      !!----------------------------------------------------------
259      if (jexport.eq.1) then
260         DO jj = 2,jpjm1
261            DO ji = 2,jpim1
262               if (tmask(ji,jj,jk) == 1) then
263                  !!=======================================================
264                  !! BALLAST SUBMODEL
265                  !!=======================================================
266                  !!
267                  !!-------------------------------------------------------
268                  !! Fast-sinking detritus fluxes, pt. 1: REMINERALISATION
269                  !! aside from explicitly modelled, slow-sinking detritus, the
270                  !! model includes an implicit representation of detrital
271                  !! particles that sink too quickly to be modelled with
272                  !! explicit state variables; this sinking flux is instead
273                  !! instantaneously remineralised down the water column using
274                  !! the version of Armstrong et al. (2002)'s ballast model
275                  !! used by Dunne et al. (2007); the version of this model
276                  !! here considers silicon and calcium carbonate ballast
277                  !! minerals; this section of the code redistributes the fast
278                  !! sinking material generated locally down the water column;
279                  !! this differs from Dunne et al. (2007) in that fast sinking
280                  !! material is distributed at *every* level below that it is
281                  !! generated, rather than at every level below some fixed
282                  !! depth; this scheme is also different in that sinking
283                  !! material generated in one level is aggregated with that
284                  !! generated by shallower levels; this should make the
285                  !! ballast model more self-consistent (famous last words)
286                  !!-------------------------------------------------------
287                  !!
288                  if (jk.eq.1) then
289                     !! this is the SURFACE OCEAN BOX (no remineralisation)
290                     !!
291                     freminc(ji,jj)  = 0.0
292                     freminn(ji,jj)  = 0.0
293                     freminfe(ji,jj) = 0.0
294                     freminsi(ji,jj) = 0.0
295                     freminca(ji,jj) = 0.0
296                  elseif (jk.le.mbathy(ji,jj)) then
297                     !! this is an OCEAN BOX (remineralise some material)
298                     !!
299                     !! set up CCD depth to be used depending on user choice
300                     if (jocalccd.eq.0) then
301                        !! use default CCD field
302                        fccd_dep(ji,jj) = ocal_ccd(ji,jj)
303                     elseif (jocalccd.eq.1) then
304                        !! use calculated CCD field
305                        fccd_dep(ji,jj) = f2_ccd_cal(ji,jj)
306                     endif
307                     !!
308                     !! === organic carbon ===
309                     !! how much organic C enters this box        (mol)
310                     fq0      = ffastc(ji,jj)
311                     if (iball.eq.1) then
312                        !! how much it weighs
313                        fq1      = (fq0 * xmassc)
314                        !! how much CaCO3 enters this box
315                        fq2      = (ffastca(ji,jj) * xmassca)
316                        !! how much  opal enters this box
317                        fq3      = (ffastsi(ji,jj) * xmasssi)
318                        !! total protected organic C
319                        fq4      = (fq2 * xprotca) + (fq3 * xprotsi)
320                        !! This next term is calculated for C but used for
321                        !! N and Fe as well
322                        !! It needs to be protected in case ALL C is protected
323                        if (fq4.lt.fq1) then
324                           !! protected fraction of total organic C (non-dim)
325                           fprotf(ji,jj)   = (fq4 / (fq1 + tiny(fq1)))
326                        else
327                           !! all organic C is protected (non-dim)
328                           fprotf(ji,jj)   = 1.0
329                        endif
330                        !! unprotected fraction of total organic C (non-dim)
331                        fq5      = (1.0 - fprotf(ji,jj))
332                        !! how much organic C is unprotected (mol)
333                        fq6      = (fq0 * fq5)
334                        !! how much unprotected C leaves this box (mol)
335                        fq7      = (fq6 * exp(-(fse3t(ji,jj,jk) / xfastc)))
336                        !! how much total C leaves this box (mol)
337                        fq8      = (fq7 + (fq0 * fprotf(ji,jj)))
338                        !! C remineralisation in this box (mol)
339                        freminc(ji,jj)  = (fq0 - fq8) / fse3t(ji,jj,jk)
340                        ffastc(ji,jj) = fq8
341# if defined key_debug_medusa
342                        !! report in/out/remin fluxes of carbon for this level
343                           if (idf.eq.1.AND.idfval.eq.1) then
344                              IF (lwp) write (numout,*)                       &
345                                       '------------------------------'
346                              IF (lwp) write (numout,*) 'totalC(',jk,')  = ', &
347                                       fq1
348                              IF (lwp) write (numout,*) 'prtctC(',jk,')  = ', &
349                                       fq4
350                              IF (lwp) write (numout,*) 'fprotf(',jk,')  = ', &
351                                       fprotf(ji,jj)
352                              IF (lwp) write (numout,*)                       &
353                                       '------------------------------'
354                              IF (lwp) write (numout,*) 'IN   C(',jk,')  = ', &
355                                       fq0
356                              IF (lwp) write (numout,*) 'LOST C(',jk,')  = ', &
357                                       freminc(ji,jj) * fse3t(ji,jj,jk)
358                              IF (lwp) write (numout,*) 'OUT  C(',jk,')  = ', &
359                                       fq8
360                              IF (lwp) write (numout,*) 'NEW  C(',jk,')  = ', &
361                                       ftempc(ji,jj) * fse3t(ji,jj,jk)
362                           endif
363# endif
364                        else
365                        !! how much organic C leaves this box (mol)
366                        fq1      = fq0 * exp(-(fse3t(ji,jj,jk) / xfastc))
367                        !! C remineralisation in this box (mol)
368                        freminc(ji,jj)  = (fq0 - fq1) / fse3t(ji,jj,jk)
369                        ffastc(ji,jj)  = fq1
370                     endif
371                     !!
372                     !! === organic nitrogen ===
373                     !! how much organic N enters this box (mol)
374                     fq0      = ffastn(ji,jj)
375                     if (iball.eq.1) then
376                        !! unprotected fraction of total organic N (non-dim)
377                        fq5      = (1.0 - fprotf(ji,jj))
378                        !! how much organic N is unprotected (mol)
379                        fq6      = (fq0 * fq5)
380                        !! how much unprotected N leaves this box (mol)
381                        fq7      = (fq6 * exp(-(fse3t(ji,jj,jk) / xfastc)))
382                        !! how much total N leaves this box (mol)
383                        fq8      = (fq7 + (fq0 * fprotf(ji,jj)))
384                        !! N remineralisation in this box (mol)
385                        freminn(ji,jj)  = (fq0 - fq8) / fse3t(ji,jj,jk)
386                        ffastn(ji,jj)  = fq8
387# if defined key_debug_medusa
388                        !! report in/out/remin fluxes of carbon for this level
389                        if (idf.eq.1.AND.idfval.eq.1) then
390                           IF (lwp) write (numout,*)                          &
391                                    '------------------------------'
392                           IF (lwp) write (numout,*) 'totalN(',jk,')  = ', fq1
393                           IF (lwp) write (numout,*) 'prtctN(',jk,')  = ', fq4
394                           IF (lwp) write (numout,*) 'fprotf(',jk,')  = ',    &
395                                    fprotf(ji,jj)
396                           IF (lwp) write (numout,*)                          &
397                                    '------------------------------'
398                           if (freminn(ji,jj) < 0.0) then
399                              IF (lwp) write (numout,*) '** FREMIN ERROR **'
400                           endif
401                           IF (lwp) write (numout,*) 'IN   N(',jk,')  = ', fq0
402                           IF (lwp) write (numout,*) 'LOST N(',jk,')  = ',    &
403                                    freminn(ji,jj) * fse3t(ji,jj,jk)
404                           IF (lwp) write (numout,*) 'OUT  N(',jk,')  = ', fq8
405                           IF (lwp) write (numout,*) 'NEW  N(',jk,')  = ',    &
406                                    ftempn(ji,jj) * fse3t(ji,jj,jk)
407                        endif
408# endif
409                     else
410                        !! how much organic N leaves this box (mol)
411                        fq1      = fq0 * exp(-(fse3t(ji,jj,jk) / xfastc))
412                        !! N remineralisation in this box (mol)
413                        freminn(ji,jj)  = (fq0 - fq1) / fse3t(ji,jj,jk)
414                        ffastn(ji,jj)  = fq1
415                     endif
416                     !!
417                     !! === organic iron ===
418                     !! how much organic Fe enters this box (mol)
419                     fq0      = ffastfe(ji,jj)
420                     if (iball.eq.1) then
421                        !! unprotected fraction of total organic Fe (non-dim)
422                        fq5      = (1.0 - fprotf(ji,jj))
423                        !! how much organic Fe is unprotected (mol)
424                        fq6      = (fq0 * fq5)
425                        !! how much unprotected Fe leaves this box (mol)
426                        fq7      = (fq6 * exp(-(fse3t(ji,jj,jk) / xfastc)))
427                        !! how much total Fe leaves this box (mol)
428                        fq8      = (fq7 + (fq0 * fprotf(ji,jj)))
429                        !! Fe remineralisation in this box (mol)
430                        freminfe(ji,jj) = (fq0 - fq8) / fse3t(ji,jj,jk)
431                        ffastfe(ji,jj) = fq8
432                     else
433                        !! how much total Fe leaves this box (mol)
434                        fq1      = fq0 * exp(-(fse3t(ji,jj,jk) / xfastc))
435                        !! Fe remineralisation in this box (mol)
436                        freminfe(ji,jj) = (fq0 - fq1) / fse3t(ji,jj,jk)
437                        ffastfe(ji,jj) = fq1
438                     endif
439                     !!
440                     !! === biogenic silicon ===
441                     !! how much  opal centers this box (mol)
442                     fq0      = ffastsi(ji,jj)
443                     !! how much  opal leaves this box (mol)
444                     fq1      = fq0 * exp(-(fse3t(ji,jj,jk) / xfastsi))
445                     !! Si remineralisation in this box (mol)
446                     freminsi(ji,jj) = (fq0 - fq1) / fse3t(ji,jj,jk)
447                     ffastsi(ji,jj) = fq1
448                     !!
449                     !! === biogenic calcium carbonate ===
450                     !! how much CaCO3 enters this box (mol)
451                     fq0      = ffastca(ji,jj)
452                     if (fsdepw(ji,jj,jk).le.fccd_dep(ji,jj)) then
453                        !! whole grid cell above CCD
454                        !! above lysocline, no Ca dissolves (mol)
455                        fq1      = fq0
456                        !! above lysocline, no Ca dissolves (mol)
457                        freminca(ji,jj) = 0.0
458                        !! which is the last level above the CCD?    (#)
459                        fccd(ji,jj) = real(jk)
460                     elseif (fsdepw(ji,jj,jk).ge.fccd_dep(ji,jj)) then
461                        !! whole grid cell below CCD
462                        !! how much CaCO3 leaves this box (mol)
463                        fq1      = fq0 * exp(-(fse3t(ji,jj,jk) / xfastca))
464                        !! Ca remineralisation in this box (mol)
465                        freminca(ji,jj) = (fq0 - fq1) / fse3t(ji,jj,jk)
466                     else
467                        !! partial grid cell below CCD
468                        !! amount of grid cell below CCD (m)
469                        fq2      = fdep1(ji,jj) - fccd_dep(ji,jj)
470                        !! how much CaCO3 leaves this box (mol)
471                        fq1      = fq0 * exp(-(fq2 / xfastca))
472                        !! Ca remineralisation in this box (mol)
473                        freminca(ji,jj) = (fq0 - fq1) / fse3t(ji,jj,jk)
474                     endif
475                     ffastca(ji,jj) = fq1 
476                  else
477                     !! this is BELOW THE LAST OCEAN BOX (do nothing)
478                     freminc(ji,jj)  = 0.0
479                     freminn(ji,jj)  = 0.0
480                     freminfe(ji,jj) = 0.0
481                     freminsi(ji,jj) = 0.0
482                     freminca(ji,jj) = 0.0             
483                  endif
484               ENDIF
485            ENDDO
486         ENDDO
487      elseif (jexport.eq.2.or.jexport.eq.3) then
488         DO jj = 2,jpjm1
489            DO ji = 2,jpim1
490               if (tmask(ji,jj,jk) == 1) then
491                  if (jexport.eq.2) then
492                     !!====================================================
493                     !! MARTIN ET AL. (1987) SUBMODEL
494                     !!====================================================
495                     !!
496                     !!----------------------------------------------------
497                     !! This submodel uses the classic Martin et al. (1987)
498                     !! curve to determine the attenuation of fast-sinking
499                     !! detritus down the water column.  All three organic
500                     !! elements, C, N and Fe, are handled identically, and
501                     !! their quantities in sinking particles attenuate
502                     !! according to a power relationship governed by
503                     !! parameter "b".  This is assigned a canonical value
504                     !! of -0.858.  Biogenic opal and calcium carbonate are
505                     !! attentuated using the same function as in the
506                     !! ballast submodel
507                     !!----------------------------------------------------
508                     !!
509                     fb_val = -0.858
510                  elseif (jexport.eq.3) then
511                     !!====================================================
512                     !! HENSON ET AL. (2011) SUBMODEL
513                     !!====================================================
514                     !!
515                     !!----------------------------------------------------
516                     !! This submodel reconfigures the Martin et al. (1987)
517                     !! curve by allowing the "b" value to vary
518                     !! geographically.  Its value is set, following Henson
519                     !! et al. (2011), as a function of local sea surface
520                     !! temperature:
521                     !!   b = -1.06 + (0.024 * SST)
522                     !! This means that remineralisation length scales are
523                     !! longer in warm, tropical areas and shorter in cold,
524                     !! polar areas.  This does seem back-to-front (i.e.
525                     !! one would expect GREATER remineralisation in warmer
526                     !! waters), but is an outcome of analysis of sediment
527                     !! trap data, and it may reflect details of ecosystem
528                     !! structure that pertain to particle production
529                     !! rather than simply Q10.
530                     !!----------------------------------------------------
531                     !!
532                     fl_sst = tsn(ji,jj,1,jp_tem)
533                     fb_val = -1.06 + (0.024 * fl_sst)
534                  endif
535                  !!   
536                  if (jk.eq.1) then
537                     !! this is the SURFACE OCEAN BOX (no remineralisation)
538                     !!
539                     freminc(ji,jj)  = 0.0
540                     freminn(ji,jj)  = 0.0
541                     freminfe(ji,jj) = 0.0
542                     freminsi(ji,jj) = 0.0
543                     freminca(ji,jj) = 0.0
544                  elseif (jk.le.mbathy(ji,jj)) then
545                     !! this is an OCEAN BOX (remineralise some material)
546                     !!
547                     !! === organic carbon ===
548                     !! how much organic C enters this box (mol)
549                     fq0      = ffastc(ji,jj)
550                     !! how much organic C leaves this box (mol)
551                     fq1      = fq0 * ((fdep1(ji,jj)/fsdepw(ji,jj,jk))**fb_val)
552                     !! C remineralisation in this box (mol)
553                     freminc(ji,jj)  = (fq0 - fq1) / fse3t(ji,jj,jk)
554                     ffastc(ji,jj)  = fq1
555                     !!
556                     !! === organic nitrogen ===
557                     !! how much organic N enters this box (mol)
558                     fq0      = ffastn(ji,jj)
559                     !! how much organic N leaves this box (mol)
560                     fq1      = fq0 * ((fdep1(ji,jj)/fsdepw(ji,jj,jk))**fb_val)
561                     !! N remineralisation in this box (mol)
562                     freminn(ji,jj)  = (fq0 - fq1) / fse3t(ji,jj,jk)
563                     ffastn(ji,jj)  = fq1
564                     !!
565                     !! === organic iron ===
566                     !! how much organic Fe enters this box (mol)
567                     fq0      = ffastfe(ji,jj)
568                     !! how much organic Fe leaves this box (mol)
569                     fq1      = fq0 * ((fdep1(ji,jj)/fsdepw(ji,jj,jk))**fb_val)
570                     !! Fe remineralisation in this box (mol)
571                     freminfe(ji,jj) = (fq0 - fq1) / fse3t(ji,jj,jk)
572                     ffastfe(ji,jj) = fq1
573                     !!
574                     !! === biogenic silicon ===
575                     !! how much  opal centers this box (mol)
576                     fq0      = ffastsi(ji,jj)
577                     !! how much  opal leaves this box (mol)
578                     fq1      = fq0 * exp(-(fse3t(ji,jj,jk) / xfastsi))
579                     !! Si remineralisation in this box (mol)
580                     freminsi(ji,jj) = (fq0 - fq1) / fse3t(ji,jj,jk)
581                     ffastsi(ji,jj) = fq1
582                     !!
583                     !! === biogenic calcium carbonate ===
584                     !! how much CaCO3 enters this box (mol)
585                     fq0      = ffastca(ji,jj)
586                     if (fsdepw(ji,jj,jk).le.ocal_ccd(ji,jj)) then
587                        !! whole grid cell above CCD
588                        !! above lysocline, no Ca dissolves (mol)
589                        fq1      = fq0
590                        !! above lysocline, no Ca dissolves (mol)
591                        freminca(ji,jj) = 0.0
592                        !! which is the last level above the CCD?    (#)
593                        fccd(ji,jj) = real(jk)
594                     elseif (fsdepw(ji,jj,jk).ge.ocal_ccd(ji,jj)) then
595                        !! whole grid cell below CCD
596                        !! how much CaCO3 leaves this box (mol)
597                        fq1      = fq0 * exp(-(fse3t(ji,jj,jk) / xfastca))
598                        !! Ca remineralisation in this box (mol)
599                        freminca(ji,jj) = (fq0 - fq1) / fse3t(ji,jj,jk)
600                     else
601                        !! partial grid cell below CCD
602                        !! amount of grid cell below CCD (m)
603                        fq2      = fdep1(ji,jj) - ocal_ccd(ji,jj)
604                        !! how much CaCO3 leaves this box (mol)
605                        fq1      = fq0 * exp(-(fq2 / xfastca))
606                        !! Ca remineralisation in this box (mol)
607                        freminca(ji,jj) = (fq0 - fq1) / fse3t(ji,jj,jk)
608                     endif
609                     ffastca(ji,jj) = fq1 
610                  else
611                     !! this is BELOW THE LAST OCEAN BOX (do nothing)
612                     freminc(ji,jj)  = 0.0
613                     freminn(ji,jj)  = 0.0
614                     freminfe(ji,jj) = 0.0
615                     freminsi(ji,jj) = 0.0
616                     freminca(ji,jj) = 0.0             
617                  endif
618               ENDIF
619            ENDDO
620         ENDDO
621      endif
622
623      DO jj = 2,jpjm1
624         DO ji = 2,jpim1
625            if (tmask(ji,jj,jk) == 1) then
626               !!----------------------------------------------------------
627               !! Fast-sinking detritus fluxes, pt. 2: UPDATE FAST FLUXES
628               !! here locally calculated additions to the fast-sinking
629               !! flux are added to the total fast-sinking flux; this is
630               !! done here such that material produced in a particular
631               !! layer is only remineralised below this layer
632               !!----------------------------------------------------------
633               !!
634               !! add sinking material generated in this layer to running
635               !! totals
636               !!
637               !! === organic carbon ===
638               !! (diatom and mesozooplankton mortality)
639               ffastc(ji,jj)  = ffastc(ji,jj)  + (ftempc(ji,jj)  *           &
640                                                  fse3t(ji,jj,jk))
641               !!
642               !! === organic nitrogen ===
643               !! (diatom and mesozooplankton mortality)
644               ffastn(ji,jj)  = ffastn(ji,jj)  + (ftempn(ji,jj)  *           &
645                                                  fse3t(ji,jj,jk))
646               !!
647               !! === organic iron ===
648               !! (diatom and mesozooplankton mortality)
649               ffastfe(ji,jj) = ffastfe(ji,jj) + (ftempfe(ji,jj) *          &
650                                                  fse3t(ji,jj,jk))
651               !!
652               !! === biogenic silicon ===
653               !! (diatom mortality and grazed diatoms)
654               ffastsi(ji,jj) = ffastsi(ji,jj) + (ftempsi(ji,jj) *          &
655                                                  fse3t(ji,jj,jk))
656               !!
657               !! === biogenic calcium carbonate ===
658               !! (latitudinally-based fraction of total primary production)
659               ffastca(ji,jj) = ffastca(ji,jj) + (ftempca(ji,jj) *          &
660                                                  fse3t(ji,jj,jk))
661            ENDIF
662         ENDDO
663      ENDDO
664
665      DO jj = 2,jpjm1
666         DO ji = 2,jpim1
667            if (tmask(ji,jj,jk) == 1) then
668               !!----------------------------------------------------------
669               !! Fast-sinking detritus fluxes, pt. 3: SEAFLOOR
670               !! remineralise all remaining fast-sinking detritus to dissolved
671               !! nutrients; the sedimentation fluxes calculated here allow the
672               !! separation of what's remineralised sinking through the final
673               !! ocean box from that which is added to the final box by the
674               !! remineralisation of material that reaches the seafloor (i.e.
675               !! the model assumes that *all* material that hits the seafloor
676               !! is remineralised and that none is permanently buried; hey,
677               !! this is a giant GCM model that can't be run for long enough
678               !! to deal with burial fluxes!)
679               !!
680               !! in a change to this process, in part so that MEDUSA behaves
681               !! a little more like ERSEM et al., fast-sinking detritus (N, Fe
682               !! and C) is converted to slow sinking detritus at the seafloor
683               !! instead of being remineralised; the rationale is that in
684               !! shallower shelf regions (... that are not fully mixed!) this
685               !! allows the detrital material to return slowly to dissolved
686               !! nutrient rather than instantaneously as now; the alternative
687               !! would be to explicitly handle seafloor organic material - a
688               !! headache I don't wish to experience at this point; note that
689               !! fast-sinking Si and Ca detritus is just remineralised as
690               !! per usual
691               !!
692               !! AXY (13/01/12)
693               !! in a further change to this process, again so that MEDUSA is
694               !! a little more like ERSEM et al., material that reaches the
695               !! seafloor can now be added to sediment pools and stored for
696               !! slow release; there are new 2D arrays for organic nitrogen,
697               !! iron and carbon and inorganic silicon and carbon that allow
698               !! fast and slow detritus that reaches the seafloor to be held
699               !! and released back to the water column more slowly; these
700               !! arrays are transferred via the tracer restart files between
701               !! repeat submissions of the model
702               !!----------------------------------------------------------
703               !!
704               ffast2slowc(ji,jj)  = 0.0
705               ffast2slown(ji,jj)  = 0.0
706! I don't think this is used - marc 10/4/17
707!               ffast2slowfe(ji,jj) = 0.0
708               !!
709               if (jk.eq.mbathy(ji,jj)) then
710                  !! this is the BOTTOM OCEAN BOX (remineralise everything)
711                  !!
712                  !! AXY (17/01/12): tweaked to include benthos pools
713                  !!
714                  !! === organic carbon ===
715                  if (jfdfate.eq.0 .and. jorgben.eq.0) then
716                     !! C remineralisation in this box (mol/m3)
717                     freminc(ji,jj)  = freminc(ji,jj) + (ffastc(ji,jj) /     &
718                                                         fse3t(ji,jj,jk))
719                  elseif (jfdfate.eq.1 .and. jorgben.eq.0) then
720                     !! fast C -> slow C (mol/m3)
721                     ffast2slowc(ji,jj) = ffastc(ji,jj) / fse3t(ji,jj,jk)
722                     fslowc(ji,jj)      = fslowc(ji,jj) + ffast2slowc(ji,jj)
723                  elseif (jfdfate.eq.0 .and. jorgben.eq.1) then
724                     !! fast C -> benthic C (mol/m2)
725                     f_fbenin_c(ji,jj)  = ffastc(ji,jj)
726                  endif
727                  !! record seafloor C (mol/m2)
728                  fsedc(ji,jj)   = ffastc(ji,jj)
729                  ffastc(ji,jj)  = 0.0
730                  !!
731                  !! === organic nitrogen ===
732                  if (jfdfate.eq.0 .and. jorgben.eq.0) then
733                     !! N remineralisation in this box (mol/m3)
734                     freminn(ji,jj)  = freminn(ji,jj) + (ffastn(ji,jj) /     &
735                                                         fse3t(ji,jj,jk))
736                  elseif (jfdfate.eq.1 .and. jorgben.eq.0) then
737                     !! fast N -> slow N (mol/m3)
738                     ffast2slown(ji,jj) = ffastn(ji,jj) / fse3t(ji,jj,jk)
739                     fslown(ji,jj)      = fslown(ji,jj) + ffast2slown(ji,jj)
740                  elseif (jfdfate.eq.0 .and. jorgben.eq.1) then
741                     !! fast N -> benthic N (mol/m2)
742                     f_fbenin_n(ji,jj)  = ffastn(ji,jj)
743                  endif
744                  !! record seafloor N (mol/m2)
745                  fsedn(ji,jj)   = ffastn(ji,jj)
746                  ffastn(ji,jj)  = 0.0
747                  !!
748                  !! === organic iron ===
749                  if (jfdfate.eq.0 .and. jorgben.eq.0) then
750                     !! Fe remineralisation in this box (mol/m3)
751                     freminfe(ji,jj) = freminfe(ji,jj) + (ffastfe(ji,jj) /   &
752                                                          fse3t(ji,jj,jk))
753! I don't think ffast2slowfe is used - marc 10/4/17
754!                  elseif (jfdfate.eq.1 .and. jorgben.eq.0) then
755!                     !! fast Fe -> slow Fe (mol/m3)
756!                     ffast2slowfe(ji,jj) = ffastn(ji,jj) / fse3t(ji,jj,jk)
757                  elseif (jfdfate.eq.0 .and. jorgben.eq.1) then
758                     !! fast Fe -> benthic Fe (mol/m2)
759                     f_fbenin_fe(ji,jj) = ffastfe(ji,jj)
760                  endif
761                  !! record seafloor Fe (mol/m2)
762                  fsedfe(ji,jj)  = ffastfe(ji,jj)
763                  ffastfe(ji,jj) = 0.0
764                  !!
765                  !! === biogenic silicon ===
766                  if (jinorgben.eq.0) then
767                     !! Si remineralisation in this box (mol/m3)
768                     freminsi(ji,jj) = freminsi(ji,jj) + (ffastsi(ji,jj) /   &
769                                                          fse3t(ji,jj,jk))
770                  elseif (jinorgben.eq.1) then
771                     !! fast Si -> benthic Si
772                     f_fbenin_si(ji,jj) = ffastsi(ji,jj)
773                  endif
774                  !! record seafloor Si (mol/m2)
775                  fsedsi(ji,jj)   = ffastsi(ji,jj)
776                  ffastsi(ji,jj) = 0.0
777                  !!
778                  !! === biogenic calcium carbonate ===
779                  if (jinorgben.eq.0) then
780                     !! Ca remineralisation in this box (mol/m3)
781                     freminca(ji,jj) = freminca(ji,jj) + (ffastca(ji,jj) /   &
782                                                          fse3t(ji,jj,jk))
783                  elseif (jinorgben.eq.1) then
784                     !! fast Ca -> benthic Ca (mol/m2)
785                     f_fbenin_ca(ji,jj) = ffastca(ji,jj)
786                  endif
787                  !! record seafloor Ca (mol/m2)
788                  fsedca(ji,jj)   = ffastca(ji,jj)
789                  ffastca(ji,jj) = 0.0
790               endif
791
792# if defined key_debug_medusa
793               if (idf.eq.1) then
794                  !!-------------------------------------------------------
795                  !! Integrate total fast detritus remineralisation
796                  !!-------------------------------------------------------
797                  !!
798                  fofd_n(ji,jj)  = fofd_n(ji,jj)  + (freminn(ji,jj)  *       &
799                                                     fse3t(ji,jj,jk))
800                  fofd_si(ji,jj) = fofd_si(ji,jj) + (freminsi(ji,jj) *       &
801                                                     fse3t(ji,jj,jk))
802                  fofd_fe(ji,jj) = fofd_fe(ji,jj) + (freminfe(ji,jj) *       &
803                                                     fse3t(ji,jj,jk))
804#  if defined key_roam
805                  fofd_c(ji,jj)  = fofd_c(ji,jj)  + (freminc(ji,jj)  *       &
806                                                     fse3t(ji,jj,jk))
807#  endif
808               endif
809# endif
810            ENDIF
811         ENDDO
812      ENDDO
813
814      DO jj = 2,jpjm1
815         DO ji = 2,jpim1
816            if (tmask(ji,jj,jk) == 1) then
817               !!----------------------------------------------------------
818               !! Sort out remineralisation tally of fast-sinking detritus
819               !!----------------------------------------------------------
820               !!
821               !! update fast-sinking regeneration arrays
822               fregenfast(ji,jj)   = fregenfast(ji,jj)   +                  &
823                                     (freminn(ji,jj)  * fse3t(ji,jj,jk))
824               fregenfastsi(ji,jj) = fregenfastsi(ji,jj) +                  &
825                                     (freminsi(ji,jj) * fse3t(ji,jj,jk))
826# if defined key_roam
827               fregenfastc(ji,jj)  = fregenfastc(ji,jj)  +                  &
828                                     (freminc(ji,jj)  * fse3t(ji,jj,jk))
829# endif
830            ENDIF
831         ENDDO
832      ENDDO
833
834      DO jj = 2,jpjm1
835         DO ji = 2,jpim1
836            if (tmask(ji,jj,jk) == 1) then
837               !!----------------------------------------------------------
838               !! Benthic remineralisation fluxes
839               !!----------------------------------------------------------
840               !!
841               if (jk.eq.mbathy(ji,jj)) then
842                  !!
843                  !! organic components
844                  if (jorgben.eq.1) then
845                     f_benout_n(ji,jj)  = xsedn  * zn_sed_n(ji,jj)
846                     f_benout_fe(ji,jj) = xsedfe * zn_sed_fe(ji,jj)
847                     f_benout_c(ji,jj)  = xsedc  * zn_sed_c(ji,jj)
848                  endif
849                  !!
850                  !! inorganic components
851                  if (jinorgben.eq.1) then
852                     f_benout_si(ji,jj) = xsedsi * zn_sed_si(ji,jj)
853                     f_benout_ca(ji,jj) = xsedca * zn_sed_ca(ji,jj)
854                     !!
855                     !! account for CaCO3 that dissolves when it shouldn't
856                     if ( fsdepw(ji,jj,jk) .le. fccd_dep(ji,jj) ) then
857                        f_benout_lyso_ca(ji,jj) = xsedca * zn_sed_ca(ji,jj)
858                     endif
859                  endif
860               endif
861               CALL flush(numout)
862
863            ENDIF
864         ENDDO
865      ENDDO
866
867   END SUBROUTINE detritus_fast_sink
868
869#else
870   !!======================================================================
871   !!  Dummy module :                                   No MEDUSA bio-model
872   !!======================================================================
873CONTAINS
874   SUBROUTINE detritus_fast_sink( )                    ! Empty routine
875      WRITE(*,*) 'detritus_fast_sink: You should not have seen this print! error?'
876   END SUBROUTINE detritus_fast_sink
877#endif 
878
879   !!======================================================================
880END MODULE detritus_fast_sink_mod
Note: See TracBrowser for help on using the repository browser.