source: branches/NERC/dev_r5518_GO6_CO2_cmip/NEMOGCM/NEMO/TOP_SRC/MEDUSA/iron_chem_scav.F90 @ 9309

Last change on this file since 9309 was 8441, checked in by frrh, 3 years ago

Commit changes relating to Met Office GMED ticket 339 for the modularisation of
of trcbio_medusa.F90.

Branch http://fcm3/projects/NEMO.xm/log/branches/NERC/dev_r5518_GO6_split_trcbiomedusa
from revisions 8394 to 8423 inclusive refer.

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