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.
iron_chem_scav.F90 in branches/UKMO/dev_r5518_GO6_fix_key_comp/NEMOGCM/NEMO/TOP_SRC/MEDUSA – NEMO

source: branches/UKMO/dev_r5518_GO6_fix_key_comp/NEMOGCM/NEMO/TOP_SRC/MEDUSA/iron_chem_scav.F90 @ 9991

Last change on this file since 9991 was 9991, checked in by frrh, 6 years ago

Fixes to allow MEDUSA to compile with C1D without
the need for multiple (apparently) unrelated CPP keys
merely to satisfy spurious code interdependencies.

File size: 25.1 KB
Line 
1MODULE iron_chem_scav_mod
2   !!======================================================================
3   !!                         ***  MODULE iron_chem_scav_mod  ***
4   !! Calculate the iron chemistry and scavenging.
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   iron_chem_scav        ! Called in trcbio_medusa.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 iron_chem_scav( jk )
28      !!-------------------------------------------------------------------
29      !!                     ***  ROUTINE iron_chem_scav  ***
30      !! This called from TRC_BIO_MEDUSA and
31      !!  -
32      !!-------------------------------------------------------------------
33      USE bio_medusa_mod,    ONLY: ffastc, ffastca, ffastsi,              &
34                                   ffetop, ffebot, ffescav, xfree,        & 
35                                   zdet, zfer, zphd, zphn, zzme, zzmi,    &
36                                   idf, idfval                         
37      USE dom_oce,           ONLY: e3t_0, gdepw_0, mbathy, tmask
38      USE par_kind,          ONLY: wp
39      USE in_out_manager,    ONLY: lwp, numout
40      USE par_oce,           ONLY: jpi, jpim1, jpj, jpjm1
41      USE sms_medusa,        ONLY: i0500, jiron, xfe_sed, xfe_sol,        &
42                                   xfe_mass,                              &
43                                   xk_FeL, xk_sc_Fe, xLgT,                &
44                                   xmassc, xmassca, xmasssi,              &
45                                   xthetad, xthetapd, xthetapn,           &
46                                   xthetazme, xthetazmi,                  &
47                                   zirondep
48
49   !!* Substitution
50#  include "domzgr_substitute.h90"
51
52      !! Level
53      INTEGER, INTENT( in ) :: jk
54
55      !! iron cycle; includes parameters for Parekh et al. (2005) iron scheme
56      !! state variables for iron-ligand system
57      REAL(wp), DIMENSION(jpi,jpj) :: xFeT, xFeF, xFeL
58      REAL(wp) :: xLgF
59      !! iron-ligand parameters
60      REAL(wp) :: xb_coef_tmp, xb2M4ac
61      !! max Fe' parameters
62      REAL(wp) :: xmaxFeF,fdeltaFe
63      !!
64      !! local parameters for Moore et al. (2004) alternative scavenging
65      !! scheme
66      REAL(wp) :: fbase_scav,fscal_sink,fscal_part,fscal_scav
67      !!
68      !! local parameters for Moore et al. (2008) alternative scavenging
69      !! scheme
70      REAL(wp) :: fscal_csink,fscal_sisink,fscal_casink
71      !!
72      !! local parameters for Galbraith et al. (2010) alternative
73      !! scavenging scheme.
74      !! organic portion of scavenging
75      REAL(wp) :: xCscav1, xCscav2, xk_org, xORGscav
76      !! inorganic portion of scavenging
77      REAL(wp) :: xk_inorg, xINORGscav
78
79      INTEGER :: ji, jj
80
81      !!------------------------------------------------------------------
82      !! Iron chemistry and fractionation
83      !! following the Parekh et al. (2004) scheme adopted by the Met.
84      !! Office, Medusa models total iron but considers "free" and
85      !! ligand-bound forms for the purposes of scavenging (only "free"
86      !! iron can be scavenged
87      !!------------------------------------------------------------------
88      DO jj = 2,jpjm1
89         DO ji = 2,jpim1
90            !! OPEN wet point IF..THEN loop
91            if (tmask(ji,jj,jk) == 1) then
92               !!
93               !! total iron concentration (mmol Fe / m3 -> umol Fe / m3)
94               xFeT(ji,jj) = zfer(ji,jj) * 1.e3
95               !!
96               !! calculate fractionation (based on Diat-HadOCC; in turn
97               !! based on Parekh et al., 2004)
98               xb_coef_tmp = xk_FeL * (xLgT - xFeT(ji,jj)) - 1.0
99               xb2M4ac     = max(((xb_coef_tmp * xb_coef_tmp) +              &
100                                  (4.0 * xk_FeL * xLgT)), 0.0)
101               !!
102               !! "free" ligand concentration
103               xLgF        = 0.5 * (xb_coef_tmp + (xb2M4ac**0.5)) / xk_FeL
104               !!
105               !! ligand-bound iron concentration
106               xFeL(ji,jj) = xLgT - xLgF
107               !!
108               !! "free" iron concentration (and convert to mmol Fe / m3)
109               xFeF(ji,jj) = (xFeT(ji,jj) - xFeL(ji,jj)) * 1.e-3
110               xFree(ji,jj)= xFeF(ji,jj) / (zfer(ji,jj) + tiny(zfer(ji,jj)))
111            ENDIF
112         ENDDO
113      ENDDO
114
115
116      !!
117      !! scavenging of iron (multiple schemes); I'm only really
118      !! happy with the first one at the moment - the others
119      !! involve assumptions (sometimes guessed at by me) that
120      !! are potentially questionable
121      !!
122      if (jiron.eq.1) then
123         !!------------------------------------------------------
124         !! Scheme 1: Dutkiewicz et al. (2005)
125         !! This scheme includes a single scavenging term based
126         !! solely on a fixed rate and the availablility of
127         !! "free" iron
128         !!------------------------------------------------------
129         DO jj = 2,jpjm1
130            DO ji = 2,jpim1
131               IF (tmask(ji,jj,jk) == 1) THEN
132                  !! = mmol/m3/d
133                  ffescav(ji,jj)     = xk_sc_Fe * xFeF(ji,jj)
134                  !!
135                  !!------------------------------------------------------
136                  !!
137                  !! Mick's code contains a further (optional) implicit
138                  !! "scavenging" of iron that sets an upper bound on
139                  !! "free" iron concentration, and essentially caps the
140                  !! concentration of total iron as xFeL + "free" iron;
141                  !! since the former is constrained by a fixed total
142                  !! ligand concentration (= 1.0 umol/m3), and the latter
143                  !! isn't allowed above this upper bound, total iron is
144                  !! constrained to a maximum of ...
145                  !!
146                  !!    xFeL(ji,jj) + min(xFeF(ji,jj), 0.3 umol/m3) = 1.0 + 0.3
147                  !!                                  = 1.3 umol / m3
148                  !!
149                  !! In Mick's code, the actual value of total iron is
150                  !! reset to this sum (i.e. TFe = FeL + Fe'; but
151                  !! Fe' <= 0.3 umol/m3); this isn't our favoured approach
152                  !! to tracer updating here (not least because of the
153                  !! leapfrog), so here the amount scavenged is augmented
154                  !! by an additional amount that serves to drag total
155                  !! iron back towards that expected from this limitation
156                  !! on iron concentration ...
157                  !!
158                  !! = umol/m3
159                  xmaxFeF     = min((xFeF(ji,jj) * 1.e3), 0.3)
160                  !!
161                  !! Here, the difference between current total Fe and
162                  !! (FeL + Fe') is calculated and added to the scavenging
163                  !! flux already calculated above ...
164                  !!
165                  !! = mmol/m3
166                  fdeltaFe    = (xFeT(ji,jj) - (xFeL(ji,jj) + xmaxFeF)) * 1.e-3
167                  !!
168                  !! This assumes that the "excess" iron is dissipated
169                  !! with a time-scale of 1 day; seems reasonable to me
170                  !! ... (famous last words)
171                  !!
172                  !! = mmol/m3/d
173                  ffescav(ji,jj)     = ffescav(ji,jj) + fdeltaFe
174                  !!
175# if defined key_deep_fe_fix
176                  !! AXY (17/01/13)
177                  !! stop scavenging for iron concentrations below
178                  !! 0.5 umol / m3 at depths greater than 1000 m; this
179                  !! aims to end MEDUSA's continual loss of iron at depth
180                  !! without impacting things at the surface too much; the
181                  !! justification for this is that it appears to be what
182                  !! Mick Follows et al. do in their work (as evidenced by
183                  !! the iron initial condition they supplied me with); to
184                  !! be honest, it looks like Follow et al. do this at
185                  !! shallower depths than 1000 m, but I'll stick with this
186                  !! for now; I suspect that this seemingly arbitrary
187                  !! approach effectively "parameterises" the
188                  !! particle-based scavenging rates that other models use
189                  !! (i.e. at depth there are no sinking particles, so
190                  !! scavenging stops); it might be fun justifying this in
191                  !! a paper though!
192                  !!
193                  if ((fsdepw(ji,jj,jk).gt.1000.) .and.                       &
194                       (xFeT(ji,jj).lt.0.5)) then
195                     ffescav(ji,jj) = 0.
196                  endif
197# endif
198               ENDIF
199            ENDDO
200         ENDDO
201      elseif (jiron.eq.2) then
202         !!------------------------------------------------------
203         !! Scheme 2: Moore et al. (2004)
204         !! This scheme includes a single scavenging term that
205         !! accounts for both suspended and sinking particles in
206         !! the water column; this term scavenges total iron rather
207         !! than "free" iron
208         !!------------------------------------------------------
209         DO jj = 2,jpjm1
210            DO ji = 2,jpim1
211               IF (tmask(ji,jj,jk) == 1) THEN
212                  !!
213                  !! total iron concentration (mmol Fe / m3 -> umol Fe / m3)
214                  xFeT(ji,jj) = zfer(ji,jj) * 1.e3
215                  !!
216                  !! this has a base scavenging rate (12% / y) which is
217                  !! modified by local particle concentration and sinking
218                  !! flux (and dust - but I'm ignoring that here for now)
219                  !! and which is accelerated when Fe concentration gets
220                  !! 0.6 nM (= 0.6 umol/m3 = 0.0006 mmol/m3), and decreased
221                  !! as concentrations below 0.4 nM (= 0.4 umol/m3 =
222                  !! 0.0004 mmol/m3)
223                  !!
224                  !! base scavenging rate (0.12 / y)
225                  fbase_scav = 0.12 / 365.25
226                  !!
227                  !! calculate sinking particle part of scaling factor
228                  !! this takes local fast sinking carbon (mmol C / m2 / d)
229                  !! and gets it into nmol C / cm3 / s ("rdt" below is the
230                  !! number of seconds in a model timestep)
231                  !!
232                  !! fscal_sink = ffastc(ji,jj) * 1.e2 / (86400.)
233                  !!
234                  !! ... actually, re-reading Moore et al.'s equations, it
235                  !! looks like he uses his sinking flux directly, without
236                  !! scaling it by time-step or anything, so I'll copy this
237                  !! here ...
238                  !!
239                  fscal_sink = ffastc(ji,jj) * 1.e2
240                  !!
241                  !! calculate particle part of scaling factor
242                  !! this totals up the carbon in suspended particles
243                  !! (Pn, Pd, Zmi, Zme, D),
244                  !! which comes out in mmol C / m3 (= nmol C / cm3), and
245                  !! then multiplies it by a magic factor, 0.002, to get it
246                  !! into nmol C / cm2 / s
247                  !!
248                  fscal_part = ( (xthetapn * zphn(ji,jj)) +                  &
249                                 (xthetapd * zphd(ji,jj)) +                  &
250                                 (xthetazmi * zzmi(ji,jj)) +                 &
251                                 (xthetazme * zzme(ji,jj)) +                 &
252                                 (xthetad * zdet(ji,jj)) ) * 0.002
253                  !!
254                  !! calculate scaling factor for base scavenging rate
255                  !! this uses the (now correctly scaled) sinking flux and
256                  !! standing
257                  !! particle concentration, divides through by some sort
258                  !! of reference value (= 0.0066 nmol C / cm2 / s) and
259                  !! then uses this, or not if its too high, to rescale the
260                  !! base scavenging rate
261                  !!
262                  fscal_scav = fbase_scav *                                  &
263                               min(((fscal_sink + fscal_part) / 0.0066), 4.0)
264                  !!
265                  !! the resulting scavenging rate is then scaled further
266                  !! according to the local iron concentration (i.e.
267                  !! diminished in low iron regions; enhanced in high iron
268                  !! regions; less alone in intermediate iron regions)
269                  !!
270                  if (xFeT(ji,jj).lt.0.4) then
271                     !!
272                     !! low iron region
273                     !!
274                     fscal_scav = fscal_scav * (xFeT(ji,jj) / 0.4)
275                     !!
276                  elseif (xFeT(ji,jj).gt.0.6) then
277                     !!
278                     !! high iron region
279                     !!
280                     fscal_scav = fscal_scav + ((xFeT(ji,jj) / 0.6) *        &
281                                                (6.0 / 1.4))
282                     !!
283                  else
284                     !!
285                     !! intermediate iron region: do nothing
286                     !!
287                  endif
288                  !!
289                  !! apply the calculated scavenging rate ...
290                  !!
291                  ffescav(ji,jj) = fscal_scav * zfer(ji,jj)
292                  !!
293               ENDIF
294            ENDDO
295         ENDDO
296      elseif (jiron.eq.3) then
297         !!------------------------------------------------------
298         !! Scheme 3: Moore et al. (2008)
299         !! This scheme includes a single scavenging term that
300         !! accounts for sinking particles in the water column,
301         !! and includes organic C, biogenic opal, calcium
302         !! carbonate and dust in this (though the latter is
303         !! ignored here until I work out what units the incoming
304         !! "dust" flux is in); this term scavenges total iron
305         !! rather than "free" iron
306         !!------------------------------------------------------
307         DO jj = 2,jpjm1
308            DO ji = 2,jpim1
309               IF (tmask(ji,jj,jk) == 1) THEN
310                  !!
311                  !! total iron concentration (mmol Fe / m3 -> umol Fe / m3)
312                  xFeT(ji,jj) = zfer(ji,jj) * 1.e3
313                  !!
314                  !! this has a base scavenging rate which is modified by
315                  !! local particle sinking flux (including dust - but I'm
316                  !! ignoring that here for now) and which is accelerated
317                  !! when Fe concentration is > 0.6 nM (= 0.6 umol/m3 =
318                  !! 0.0006 mmol/m3), and decreased as concentrations <
319                  !! 0.5 nM (= 0.5 umol/m3 = 0.0005 mmol/m3)
320                  !!
321                  !! base scavenging rate (Fe_b in paper; units may be
322                  !! wrong there)
323                  fbase_scav = 0.00384 ! (ng)^-1 cm
324                  !!
325                  !! calculate sinking particle part of scaling factor;
326                  !! this converts mmol / m2 / d fluxes of organic carbon,
327                  !! silicon and calcium carbonate into ng / cm2 / s
328                  !! fluxes; it is assumed here that the mass conversions
329                  !! simply consider the mass of the main element
330                  !! (C, Si and Ca) and *not* the mass of the molecules
331                  !! that they are part of; Moore et al. (2008) is unclear
332                  !! on the conversion that should be used
333                  !!
334                  !! milli -> nano; mol -> gram; /m2 -> /cm2; /d -> /s
335                  !! ng C  / cm2 / s
336                  fscal_csink  = (ffastc(ji,jj)  * 1.e6 * xmassc  *          &
337                                  1.e-4 / 86400.)
338                  !! ng Si / cm2 / s
339                  fscal_sisink = (ffastsi(ji,jj) * 1.e6 * xmasssi *          &
340                                  1.e-4 / 86400.)
341                  !! ng Ca / cm2 / s
342                  fscal_casink = (ffastca(ji,jj) * 1.e6 * xmassca *          &
343                                  1.e-4 / 86400.)
344                  !!
345                  !! sum up these sinking fluxes and convert to ng / cm
346                  !! by dividing through by a sinking rate of
347                  !! 100 m / d = 1.157 cm / s
348                  !! ng / cm
349                  fscal_sink   = ((fscal_csink * 6.) + fscal_sisink +        &
350                                  fscal_casink) / (100. * 1.e3 / 86400)
351                  !!
352                  !! now calculate the scavenging rate based upon the base
353                  !! rate and this particle flux scaling; according to the
354                  !! published units, the result actually has *no* units,
355                  !! but as it must be expressed per unit time for it to
356                  !! make any sense, I'm assuming a missing "per second"
357                  !! / s
358                  fscal_scav = fbase_scav * fscal_sink
359                  !!
360                  !! the resulting scavenging rate is then scaled further
361                  !! according to the local iron concentration (i.e.
362                  !! diminished in low iron regions; enhanced in high iron
363                  !! regions; less alone in intermediate iron regions)
364                  !!
365                  if (xFeT(ji,jj).lt.0.5) then
366                     !!
367                     !! low iron region (0.5 instead of the 0.4 in Moore
368                     !! et al., 2004)
369                     !!
370                     fscal_scav = fscal_scav * (xFeT(ji,jj) / 0.5)
371                     !!
372                  elseif (xFeT(ji,jj).gt.0.6) then
373                     !!
374                     !! high iron region (functional form different in
375                     !! Moore et al., 2004)
376                     !!
377                     fscal_scav = fscal_scav + ((xFeT(ji,jj) - 0.6) * 0.00904)
378                     !!
379                  else
380                     !!
381                     !! intermediate iron region: do nothing
382                     !!
383                  endif
384                  !!
385                  !! apply the calculated scavenging rate ...
386                  !!
387                  ffescav(ji,jj) = fscal_scav * zfer(ji,jj)
388               ENDIF
389            ENDDO
390         ENDDO
391      elseif (jiron.eq.4) then
392         !!------------------------------------------------------
393         !! Scheme 4: Galbraith et al. (2010)
394         !! This scheme includes two scavenging terms, one for
395         !! organic, particle-based scavenging, and another for
396         !! inorganic scavenging; both terms scavenge "free" iron
397         !! only
398         !!------------------------------------------------------
399         DO jj = 2,jpjm1
400            DO ji = 2,jpim1
401               IF (tmask(ji,jj,jk) == 1) THEN
402                  !!
403                  !! Galbraith et al. (2010) present a more straightforward
404                  !! outline of the scheme in Parekh et al. (2005) ...
405                  !!
406                  !! sinking particulate carbon available for scavenging
407                  !! this assumes a sinking rate of 100 m / d (Moore &
408                  !! Braucher, 2008),
409                  xCscav1    = (ffastc(ji,jj) * xmassc) / 100. ! = mg C / m3
410                  !!
411                  !! scale by Honeyman et al. (1981) exponent coefficient
412                  !! multiply by 1.e-3 to express C flux in g C rather than
413                  !! mg C
414                  xCscav2    = (xCscav1 * 1.e-3)**0.58
415                  !!
416                  !! multiply by Galbraith et al. (2010) scavenging rate
417                  xk_org     = 0.5 ! ((g C m/3)^-1) / d
418                  xORGscav   = xk_org * xCscav2 * xFeF(ji,jj)
419                  !!
420                  !! Galbraith et al. (2010) also include an inorganic bit ...
421                  !!
422                  !! this occurs at a fixed rate, again based on the
423                  !! availability of "free" iron
424                  !!
425                  !! k_inorg = 1000 d**-1 nmol Fe**-0.5 kg**-0.5
426                  !!
427                  !! to implement this here, scale xFeF by 1026 to put in
428                  !! units of umol/m3 which approximately equal nmol/kg
429                  !!
430                  xk_inorg   = 1000. ! ((nmol Fe / kg)^1.5)
431                  xINORGscav = (xk_inorg * (xFeF(ji,jj) * 1026.)**1.5) * 1.e-3
432                  !!
433                  !! sum these two terms together
434                  ffescav(ji,jj) = xORGscav + xINORGscav
435               ENDIF
436            ENDDO
437         ENDDO
438      else
439         !!------------------------------------------------------
440         !! No Scheme: you coward!
441         !! This scheme puts its head in the sand and eskews any
442         !! decision about how iron is removed from the ocean;
443         !! prepare to get deluged in iron you fool!
444         !!------------------------------------------------------
445         DO jj = 2,jpjm1
446            DO ji = 2,jpim1
447               IF (tmask(ji,jj,jk) == 1) THEN
448                  ffescav(ji,jj) = 0.
449               ENDIF
450            ENDDO
451         ENDDO
452      endif
453
454      !!---------------------------------------------------------
455      !! Other iron cycle processes
456      !!---------------------------------------------------------
457      !!
458      !! aeolian iron deposition
459      !! zirondep      is in mmol-Fe / m2 / day
460      !! ffetop(ji,jj) is in mmol-dissolved-Fe / m3 / day
461      if (jk == 1) then
462         DO jj = 2,jpjm1
463            DO ji = 2,jpim1
464               IF (tmask(ji,jj,jk) == 1) THEN
465                  ffetop(ji,jj)  = zirondep(ji,jj) * xfe_sol / fse3t(ji,jj,jk) 
466               ENDIF
467            ENDDO
468         ENDDO
469      else
470         DO jj = 2,jpjm1
471            DO ji = 2,jpim1
472               IF (tmask(ji,jj,jk) == 1) THEN
473                  ffetop(ji,jj)  = 0.0
474                ENDIF
475            ENDDO
476         ENDDO
477      endif
478      !!
479      !! seafloor iron addition
480      DO jj = 2,jpjm1
481         DO ji = 2,jpim1
482            IF (tmask(ji,jj,jk) == 1) THEN
483               !! AXY (10/07/12): amended to only apply sedimentary flux up
484               !! to ~500 m down
485               !! if (jk.eq.(mbathy(ji,jj)-1).AND.jk.lt.i1100) then
486               if ((jk.eq.mbathy(ji,jj)).AND.jk.le.i0500) then
487                  !! Moore et al. (2004) cite a coastal California value of
488                  !! 5 umol/m2/d, but adopt a global value of 2 umol/m2/d
489                  !! for all areas < 1100 m; here we use this latter value
490                  !! but apply it everywhere
491                  !! AXY (21/07/09): actually, let's just apply it below
492                  !! 1100 m (levels 1-37)
493                  ffebot(ji,jj)  = (xfe_sed / fse3t(ji,jj,jk))
494               else
495                  ffebot(ji,jj)  = 0.0
496               endif
497            ENDIF
498         ENDDO
499      ENDDO
500
501      !! AXY (16/12/09): remove iron addition/removal processes
502      !! For the purposes of the quarter degree run, the iron
503      !! cycle is being pegged to the initial condition supplied
504      !! by Mick Follows via restoration with a 30 day period;
505      !! iron addition at the seafloor is still permitted with
506      !! the idea that this extra iron will be removed by the
507      !! restoration away from the source
508      !! ffescav(ji,jj) = 0.0
509      !! ffetop(ji,jj)  = 0.0
510      !! ffebot(ji,jj)  = 0.0
511
512# if defined key_debug_medusa
513      !! report miscellaneous calculations
514      !! report miscellaneous calculations
515      if (idf.eq.1.AND.idfval.eq.1) then
516         DO jj = 2,jpjm1
517            DO ji = 2,jpim1
518               IF (tmask(ji,jj,jk) == 1) THEN
519                  IF (lwp) write (numout,*) '------------------------------'
520                  IF (lwp) write (numout,*) 'xfe_sol  = ', xfe_sol
521                  IF (lwp) write (numout,*) 'xfe_mass = ', xfe_mass
522                  IF (lwp) write (numout,*) 'ffetop(',jk,')  = ', ffetop(ji,jj)
523                  IF (lwp) write (numout,*) 'ffebot(',jk,')  = ', ffebot(ji,jj)
524                  IF (lwp) write (numout,*) 'xFree(',jk,')   = ', xFree(ji,jj)
525                  IF (lwp) write (numout,*) 'ffescav(',jk,') = ', ffescav(ji,jj)
526               ENDIF
527            ENDDO
528         ENDDO
529      endif
530# endif
531
532   END SUBROUTINE iron_chem_scav
533
534#else
535   !!======================================================================
536   !!  Dummy module :                                   No MEDUSA bio-model
537   !!======================================================================
538CONTAINS
539   SUBROUTINE iron_chem_scav( )                    ! Empty routine
540      WRITE(*,*) 'iron_chem_scav: You should not have seen this print! error?'
541   END SUBROUTINE iron_chem_scav
542#endif 
543
544   !!======================================================================
545END MODULE iron_chem_scav_mod
Note: See TracBrowser for help on using the repository browser.