source: NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/PISCES/P4Z/p4zfechem.F90 @ 10368

Last change on this file since 10368 was 10368, checked in by smasson, 2 years ago

dev_r10164_HPC09_ESIWACE_PREP_MERGE: merge with trunk@10365, see #2133

  • Property svn:keywords set to Id
File size: 23.0 KB
Line 
1MODULE p4zfechem
2   !!======================================================================
3   !!                         ***  MODULE p4zfechem  ***
4   !! TOP :   PISCES Compute iron chemistry and scavenging
5   !!======================================================================
6   !! History :   3.5  !  2012-07 (O. Aumont, A. Tagliabue, C. Ethe) Original code
7   !!             3.6  !  2015-05  (O. Aumont) PISCES quota
8   !!----------------------------------------------------------------------
9   !!   p4z_fechem       : Compute remineralization/scavenging of iron
10   !!   p4z_fechem_init  : Initialisation of parameters for remineralisation
11   !!   p4z_fechem_alloc : Allocate remineralisation variables
12   !!----------------------------------------------------------------------
13   USE oce_trc         ! shared variables between ocean and passive tracers
14   USE trc             ! passive tracers common variables
15   USE sms_pisces      ! PISCES Source Minus Sink variables
16   USE p4zche          ! chemical model
17   USE p4zsbc          ! Boundary conditions from sediments
18   USE prtctl_trc      ! print control for debugging
19   USE iom             ! I/O manager
20
21   IMPLICIT NONE
22   PRIVATE
23
24   PUBLIC   p4z_fechem        ! called in p4zbio.F90
25   PUBLIC   p4z_fechem_init   ! called in trcsms_pisces.F90
26
27   LOGICAL          ::   ln_fechem    !: boolean for complex iron chemistry following Tagliabue and voelker
28   LOGICAL          ::   ln_ligvar    !: boolean for variable ligand concentration following Tagliabue and voelker
29   REAL(wp), PUBLIC ::   xlam1        !: scavenging rate of Iron
30   REAL(wp), PUBLIC ::   xlamdust     !: scavenging rate of Iron by dust
31   REAL(wp), PUBLIC ::   ligand       !: ligand concentration in the ocean
32   REAL(wp), PUBLIC ::   kfep         !: rate constant for nanoparticle formation
33
34   REAL(wp) :: kl1, kl2, kb1, kb2, ks, kpr, spd, con, kth      !!gm  <<<== DOCTOR names SVP !!!
35
36   !!----------------------------------------------------------------------
37   !! NEMO/TOP 4.0 , NEMO Consortium (2018)
38   !! $Id$
39   !! Software governed by the CeCILL license (see ./LICENSE)
40   !!----------------------------------------------------------------------
41CONTAINS
42
43   SUBROUTINE p4z_fechem( kt, knt )
44      !!---------------------------------------------------------------------
45      !!                     ***  ROUTINE p4z_fechem  ***
46      !!
47      !! ** Purpose :   Compute remineralization/scavenging of iron
48      !!
49      !! ** Method  :   2 different chemistry models are available for iron
50      !!                (1) The simple chemistry model of Aumont and Bopp (2006)
51      !!                    based on one ligand and one inorganic form
52      !!                (2) The complex chemistry model of Tagliabue and
53      !!                    Voelker (2009) based on 2 ligands, 2 inorganic forms
54      !!                    and one particulate form (ln_fechem)
55      !!---------------------------------------------------------------------
56      INTEGER, INTENT(in) ::   kt, knt   ! ocean time step
57      !
58      INTEGER  ::   ji, jj, jk, jic, jn
59      REAL(wp) ::   zdep, zlam1a, zlam1b, zlamfac
60      REAL(wp) ::   zkeq, zfeequi, zfesatur, zfecoll, fe3sol
61      REAL(wp) ::   zdenom1, zscave, zaggdfea, zaggdfeb, zcoag
62      REAL(wp) ::   ztrc, zdust
63      REAL(wp) ::   zdenom2
64      REAL(wp) ::   zzFeL1, zzFeL2, zzFe2, zzFeP, zzFe3, zzstrn2
65      REAL(wp) ::   zrum, zcodel, zargu, zlight
66      REAL(wp) ::   zkox, zkph1, zkph2, zph, zionic, ztligand
67      REAL(wp) ::   za, zb, zc, zkappa1, zkappa2, za0, za1, za2
68      REAL(wp) ::   zxs, zfunc, zp, zq, zd, zr, zphi, zfff, zp3, zq2
69      REAL(wp) ::   ztfe, zoxy, zhplus, zxlam
70      REAL(wp) ::   zaggliga, zaggligb
71      REAL(wp) ::   dissol, zligco
72      REAL(wp) :: zrfact2
73      CHARACTER (len=25) :: charout
74      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zTL1, zFe3, ztotlig, precip, zFeL1
75      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zcoll3d, zscav3d, zlcoll3d
76      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zFeL2, zTL2, zFe2, zFeP
77      REAL(wp), ALLOCATABLE, DIMENSION(:,:  ) ::   zstrn, zstrn2
78      !!---------------------------------------------------------------------
79      !
80      IF( ln_timing )   CALL timing_start('p4z_fechem')
81      !
82      zFe3 (:,:,:) = 0.
83      zFeL1(:,:,:) = 0.
84      zTL1 (:,:,:) = 0.
85      IF( ln_fechem ) THEN
86         ALLOCATE( zstrn(jpi,jpj), zstrn2(jpi,jpj) )
87         ALLOCATE( zFe2(jpi,jpj,jpk), zFeL2(jpi,jpj,jpk), zTL2(jpi,jpj,jpk), zFeP(jpi,jpj,jpk) )
88         zFe2 (:,:,:) = 0.
89         zFeL2(:,:,:) = 0.
90         zTL2 (:,:,:) = 0.
91         zFeP (:,:,:) = 0.
92      ENDIF
93
94      ! Total ligand concentration : Ligands can be chosen to be constant or variable
95      ! Parameterization from Tagliabue and Voelker (2011)
96      ! -------------------------------------------------
97      IF( ln_ligvar ) THEN
98         ztotlig(:,:,:) =  0.09 * trb(:,:,:,jpdoc) * 1E6 + ligand * 1E9
99         ztotlig(:,:,:) =  MIN( ztotlig(:,:,:), 10. )
100      ELSE
101        IF( ln_ligand ) THEN  ;   ztotlig(:,:,:) = trb(:,:,:,jplgw) * 1E9
102        ELSE                  ;   ztotlig(:,:,:) = ligand * 1E9
103        ENDIF
104      ENDIF
105
106      IF( ln_fechem ) THEN
107         ! compute the day length depending on latitude and the day
108         zrum = REAL( nday_year - 80, wp ) / REAL( nyear_len(1), wp )
109         zcodel = ASIN(  SIN( zrum * rpi * 2._wp ) * SIN( rad * 23.5_wp )  )
110
111         ! day length in hours
112         zstrn(:,:) = 0.
113         DO jj = 1, jpj
114            DO ji = 1, jpi
115               zargu = TAN( zcodel ) * TAN( gphit(ji,jj) * rad )
116               zargu = MAX( -1., MIN(  1., zargu ) )
117               zstrn(ji,jj) = MAX( 0.0, 24. - 2. * ACOS( zargu ) / rad / 15. )
118            END DO
119         END DO
120
121         ! Maximum light intensity
122         zstrn2(:,:) = zstrn(:,:) / 24.
123         WHERE( zstrn(:,:) < 1.e0 ) zstrn(:,:) = 24.
124         zstrn(:,:) = 24. / zstrn(:,:)
125
126         ! ------------------------------------------------------------
127         ! NEW FE CHEMISTRY ROUTINE from Tagliabue and Volker (2009)
128         ! This model is based on two ligands, Fe2+, Fe3+ and Fep
129         ! Chemistry is supposed to be fast enough to be at equilibrium
130         ! ------------------------------------------------------------
131         DO jn = 1, 2
132            DO jk = 1, jpkm1
133               DO jj = 1, jpj
134                  DO ji = 1, jpi
135                     zlight  = etot(ji,jj,jk) * zstrn(ji,jj) * REAL( 2-jn, wp )
136                     zzstrn2 = zstrn2(ji,jj) * REAL( 2-jn, wp ) + (1. - zstrn2(ji,jj) ) * REAL( jn-1, wp )
137                     ! Calculate ligand concentrations : assume 2/3rd of excess goes to
138                     ! strong ligands (L1) and 1/3rd to weak ligands (L2)
139                     ztligand       = ztotlig(ji,jj,jk) - ligand * 1E9
140                     zTL1(ji,jj,jk) =                0.000001 + 0.67 * ztligand
141                     zTL2(ji,jj,jk) = ligand * 1E9 - 0.000001 + 0.33 * ztligand
142                     ! ionic strength from Millero et al. 1987
143                     zph    = -LOG10( MAX( hi(ji,jj,jk), rtrn) )
144                     zoxy   = trb(ji,jj,jk,jpoxy)
145                     ! Fe2+ oxydation rate from Santana-Casiano et al. (2005)
146                     zkox   = 35.407 - 6.7109 * zph + 0.5342 * zph * zph - 5362.6 / ( tempis(ji,jj,jk) + 273.15 )  &
147                     &        - 0.04406 * SQRT( salinprac(ji,jj,jk) ) - 0.002847 * salinprac(ji,jj,jk)
148                     zkox   = ( 10.** zkox ) * spd
149                     zkox   = zkox * MAX( 1.e-6, zoxy) / ( chemo2(ji,jj,jk) + rtrn )
150                     ! PHOTOREDUCTION of complexed iron : Tagliabue and Arrigo (2006)
151                     zkph2 = MAX( 0., 15. * zlight / ( zlight + 2. ) ) * (1. - fr_i(ji,jj))
152                     zkph1 = zkph2 / 5.
153                     ! pass the dfe concentration from PISCES
154                     ztfe  = trb(ji,jj,jk,jpfer) * 1e9
155                     ! ----------------------------------------------------------
156                     ! ANALYTICAL SOLUTION OF ROOTS OF THE FE3+ EQUATION
157                     ! As shown in Tagliabue and Voelker (2009), Fe3+ is the root of a 3rd order polynom.
158                     ! ----------------------------------------------------------
159                     ! calculate some parameters
160                     za = 1.0 + ks / kpr
161                     zb = 1.0 + zkph2 / ( zkox )
162                     zc = 1.0 + ( zkph1 + kth ) / ( zkox )
163                     zkappa1 = ( kb1 + zkph1 + kth ) / kl1
164                     zkappa2 = ( kb2 + zkph2 ) / kl2
165                     za2 = zTL2(ji,jj,jk) * zb / za + zTL2(ji,jj,jk) * zc / za + zkappa1 + zkappa2 - ztfe / za
166                     za1 = zkappa1 * zTL2(ji,jj,jk) * zb / za + zkappa2 * zTL1(ji,jj,jk) * zc / za &
167                      & + zkappa1 * zkappa2 - ( zkappa1 + zkappa2 ) * ztfe / za
168                     za0 = -zkappa1 * zkappa2 * ztfe / za
169                     zp  = za1 - za2 * za2 / 3.
170                     zq  = za2 * za2 * za2 * 2. / 27. - za2 * za1 / 3. + za0
171                     zp3 = zp / 3.
172                     zq2 = zq / 2.
173                     zd  = zp3 * zp3 * zp3 + zq2 * zq2
174                     zr  = zq / ABS( zq ) * SQRT( ABS( zp ) / 3. )
175                     ! compute the roots
176                     IF( zp > 0.) THEN
177                        ! zphi = ASINH( zq / ( 2. * zr * zr * zr ) )
178                        zphi =  zq / ( 2. * zr * zr * zr ) 
179                        zphi = LOG( zphi + SQRT( zphi * zphi + 1 ) )  ! asinh(x) = log(x + sqrt(x^2+1))
180                        zxs  = -2. * zr * SINH( zphi / 3. ) - za1 / 3.
181                     ELSE
182                        IF( zd > 0. ) THEN
183                           zfff = MAX( 1., zq / ( 2. * zr * zr * zr ) )
184                           ! zphi = ACOSH( zfff )
185                           zphi = LOG( zfff + SQRT( zfff * zfff - 1 ) )  ! acosh(x) = log(x + sqrt(x^2-1))
186                           zxs = -2. * zr * COSH( zphi / 3. ) - za1 / 3.
187                        ELSE
188                           zfff = MIN( 1., zq / ( 2. * zr * zr * zr ) )
189                           zphi = ACOS( zfff )
190                           DO jic = 1, 3
191                              zfunc = -2 * zr * COS( zphi / 3. + 2. * REAL( jic - 1, wp ) * rpi / 3. ) - za2 / 3.
192                              IF( zfunc > 0. .AND. zfunc <= ztfe)  zxs = zfunc
193                           END DO
194                        ENDIF
195                     ENDIF
196                     ! solve for the other Fe species
197                     zzFe3  = MAX( 0., zxs )
198                     zzFep  = MAX( 0., ( ks * zzFe3 / kpr ) )
199                     zzFeL1 = MAX( 0., ( zzFe3 * zTL1(ji,jj,jk) ) / ( zkappa1 + zzFe3 ) )
200                     zzFeL2 = (ztfe - za * zzFe3 - zc * zzFeL1 ) / zb
201                     zzFe2  = MAX( 0., ( ( ( zkph1 + kth ) * zzFeL1 + zkph2 * zzFeL2 ) / zkox ) )
202                     zzFep  = ztfe - zzFe3 - zzFe2 - zzFeL1 - zzFeL2
203                     zFe3(ji,jj,jk)  = zFe3(ji,jj,jk)  + zzFe3 * zzstrn2
204                     zFe2(ji,jj,jk)  = zFe2(ji,jj,jk)  + zzFe2 * zzstrn2
205                     zFeL2(ji,jj,jk) = zFeL2(ji,jj,jk) + zzFeL2 * zzstrn2
206                     zFeL1(ji,jj,jk) = zFeL1(ji,jj,jk) + zzFeL1 * zzstrn2
207                     zFeP(ji,jj,jk)  = zFeP(ji,jj,jk)  + zzFeP * zzstrn2
208                  END DO
209               END DO
210            END DO
211         END DO
212      ELSE
213         ! ------------------------------------------------------------
214         ! OLD FE CHEMISTRY ROUTINE from Aumont and Bopp (2006)
215         ! This model is based on one ligand and Fe'
216         ! Chemistry is supposed to be fast enough to be at equilibrium
217         ! ------------------------------------------------------------
218         DO jk = 1, jpkm1
219            DO jj = 1, jpj
220               DO ji = 1, jpi
221                  zTL1(ji,jj,jk)  = ztotlig(ji,jj,jk)
222                  zkeq            = fekeq(ji,jj,jk)
223                  zfesatur        = zTL1(ji,jj,jk) * 1E-9
224                  ztfe            = trb(ji,jj,jk,jpfer) 
225                  ! Fe' is the root of a 2nd order polynom
226                  zFe3 (ji,jj,jk) = ( -( 1. + zfesatur * zkeq - zkeq * ztfe )               &
227                     &              + SQRT( ( 1. + zfesatur * zkeq - zkeq * ztfe )**2       &
228                     &              + 4. * ztfe * zkeq) ) / ( 2. * zkeq )
229                  zFe3 (ji,jj,jk) = zFe3(ji,jj,jk) * 1E9
230                  zFeL1(ji,jj,jk) = MAX( 0., trb(ji,jj,jk,jpfer) * 1E9 - zFe3(ji,jj,jk) )
231              END DO
232            END DO
233         END DO
234         !
235      ENDIF
236
237      zdust = 0.         ! if no dust available
238      DO jk = 1, jpkm1
239         DO jj = 1, jpj
240            DO ji = 1, jpi
241               ! Scavenging rate of iron. This scavenging rate depends on the load of particles of sea water.
242               ! This parameterization assumes a simple second order kinetics (k[Particles][Fe]).
243               ! Scavenging onto dust is also included as evidenced from the DUNE experiments.
244               ! --------------------------------------------------------------------------------------
245               zhplus  = max( rtrn, hi(ji,jj,jk) )
246               fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  &
247               &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     &
248               &         + fesol(ji,jj,jk,5) / zhplus )
249               IF( ln_fechem ) THEN
250                  zfeequi = ( zFe3(ji,jj,jk) + zFe2(ji,jj,jk) + zFeP(ji,jj,jk) ) * 1E-9
251                  zfecoll = ( 0.3 * zFeL1(ji,jj,jk) + 0.5 * zFeL2(ji,jj,jk) ) * 1E-9
252                  precip(ji,jj,jk) = 0.0
253               ELSE
254                  zfeequi = zFe3(ji,jj,jk) * 1E-9
255                  zhplus  = max( rtrn, hi(ji,jj,jk) )
256                  fe3sol  = fesol(ji,jj,jk,1) * ( zhplus**3 + fesol(ji,jj,jk,2) * zhplus**2  &
257                  &         + fesol(ji,jj,jk,3) * zhplus + fesol(ji,jj,jk,4)     &
258                  &         + fesol(ji,jj,jk,5) / zhplus )
259                  zfecoll = 0.5 * zFeL1(ji,jj,jk) * 1E-9
260                  ! precipitation of Fe3+, creation of nanoparticles
261                  precip(ji,jj,jk) = MAX( 0., ( zFe3(ji,jj,jk) * 1E-9 - fe3sol ) ) * kfep * xstep
262               ENDIF
263               !
264               ztrc   = ( trb(ji,jj,jk,jppoc) + trb(ji,jj,jk,jpgoc) + trb(ji,jj,jk,jpcal) + trb(ji,jj,jk,jpgsi) ) * 1.e6 
265               IF( ln_dust )  zdust  = dust(ji,jj) / ( wdust / rday ) * tmask(ji,jj,jk) &
266               &  * EXP( -gdept_n(ji,jj,jk) / 540. )
267               IF (ln_ligand) THEN
268                  zxlam  = xlam1 * MAX( 1.E-3, EXP(-2 * etot(ji,jj,jk) / 10. ) * (1. - EXP(-2 * trb(ji,jj,jk,jpoxy) / 100.E-6 ) ))
269               ELSE
270                  zxlam  = xlam1 * 1.0
271               ENDIF
272               zlam1b = 3.e-5 + xlamdust * zdust + zxlam * ztrc
273               zscave = zfeequi * zlam1b * xstep
274
275               ! Compute the different ratios for scavenging of iron
276               ! to later allocate scavenged iron to the different organic pools
277               ! ---------------------------------------------------------
278               zdenom1 = zxlam * trb(ji,jj,jk,jppoc) / zlam1b
279               zdenom2 = zxlam * trb(ji,jj,jk,jpgoc) / zlam1b
280
281               !  Increased scavenging for very high iron concentrations found near the coasts
282               !  due to increased lithogenic particles and let say it is unknown processes (precipitation, ...)
283               !  -----------------------------------------------------------
284               zlamfac = MAX( 0.e0, ( gphit(ji,jj) + 55.) / 30. )
285               zlamfac = MIN( 1.  , zlamfac )
286               zdep    = MIN( 1., 1000. / gdept_n(ji,jj,jk) )
287               zcoag   = 1E-4 * ( 1. - zlamfac ) * zdep * xstep * trb(ji,jj,jk,jpfer)
288
289               !  Compute the coagulation of colloidal iron. This parameterization
290               !  could be thought as an equivalent of colloidal pumping.
291               !  It requires certainly some more work as it is very poorly constrained.
292               !  ----------------------------------------------------------------
293               zlam1a   = ( 0.369  * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4  * trb(ji,jj,jk,jppoc) ) * xdiss(ji,jj,jk)    &
294                   &      + ( 114.   * 0.3 * trb(ji,jj,jk,jpdoc) )
295               zaggdfea = zlam1a * xstep * zfecoll
296               !
297               zlam1b   = 3.53E3 * trb(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk)
298               zaggdfeb = zlam1b * xstep * zfecoll
299               !
300               tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - zscave - zaggdfea - zaggdfeb &
301               &                     - zcoag - precip(ji,jj,jk)
302               tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + zscave * zdenom1 + zaggdfea
303               tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zscave * zdenom2 + zaggdfeb
304               zscav3d(ji,jj,jk)   = zscave
305               zcoll3d(ji,jj,jk)   = zaggdfea + zaggdfeb
306               !
307            END DO
308         END DO
309      END DO
310      !
311      !  Define the bioavailable fraction of iron
312      !  ----------------------------------------
313      IF( ln_fechem ) THEN  ;  biron(:,:,:) = MAX( 0., trb(:,:,:,jpfer) - zFeP(:,:,:) * 1E-9 )
314      ELSE                  ;  biron(:,:,:) = trb(:,:,:,jpfer) 
315      ENDIF
316      !
317      IF( ln_ligand ) THEN
318         !
319         DO jk = 1, jpkm1
320            DO jj = 1, jpj
321               DO ji = 1, jpi
322                  zlam1a   = ( 0.369  * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4  * trb(ji,jj,jk,jppoc) ) * xdiss(ji,jj,jk)    &
323                      &    + ( 114.   * 0.3 * trb(ji,jj,jk,jpdoc) )
324                  !
325                  zlam1b   = 3.53E3 *   trb(ji,jj,jk,jpgoc) * xdiss(ji,jj,jk)
326                  zligco   = 0.5 * trn(ji,jj,jk,jplgw)
327                  zaggliga = zlam1a * xstep * zligco
328                  zaggligb = zlam1b * xstep * zligco
329                  tra(ji,jj,jk,jpfep) = tra(ji,jj,jk,jpfep) + precip(ji,jj,jk)
330                  tra(ji,jj,jk,jplgw) = tra(ji,jj,jk,jplgw) - zaggliga - zaggligb
331                  zlcoll3d(ji,jj,jk)  = zaggliga + zaggligb
332               END DO
333            END DO
334         END DO
335         !
336         IF( .NOT.ln_fechem) THEN
337            plig(:,:,:) =  MAX( 0., ( ( zFeL1(:,:,:) * 1E-9 ) / ( trb(:,:,:,jpfer) +rtrn ) ) )
338         ENDIF
339         !
340      ENDIF
341      !  Output of some diagnostics variables
342      !     ---------------------------------
343      IF( lk_iomput ) THEN
344         IF( knt == nrdttrc ) THEN
345         zrfact2 = 1.e3 * rfact2r  ! conversion from mol/L/timestep into mol/m3/s
346         IF( iom_use("Fe3")    )  CALL iom_put("Fe3"    , zFe3   (:,:,:)       * tmask(:,:,:) )   ! Fe3+
347         IF( iom_use("FeL1")   )  CALL iom_put("FeL1"   , zFeL1  (:,:,:)       * tmask(:,:,:) )   ! FeL1
348         IF( iom_use("TL1")    )  CALL iom_put("TL1"    , zTL1   (:,:,:)       * tmask(:,:,:) )   ! TL1
349         IF( iom_use("Totlig") )  CALL iom_put("Totlig" , ztotlig(:,:,:)       * tmask(:,:,:) )   ! TL
350         IF( iom_use("Biron")  )  CALL iom_put("Biron"  , biron  (:,:,:)  * 1e9 * tmask(:,:,:) )   ! biron
351         IF( iom_use("FESCAV") )  CALL iom_put("FESCAV" , zscav3d(:,:,:)  * 1e9 * tmask(:,:,:) * zrfact2 )
352         IF( iom_use("FECOLL") )  CALL iom_put("FECOLL" , zcoll3d(:,:,:)  * 1e9 * tmask(:,:,:) * zrfact2 )
353         IF( iom_use("LGWCOLL"))  CALL iom_put("LGWCOLL", zlcoll3d(:,:,:) * 1e9 * tmask(:,:,:) * zrfact2 )
354         IF( ln_fechem ) THEN
355            IF( iom_use("Fe2")  ) CALL iom_put("Fe2"    , zFe2   (:,:,:)       * tmask(:,:,:) )   ! Fe2+
356            IF( iom_use("FeL2") ) CALL iom_put("FeL2"   , zFeL2  (:,:,:)       * tmask(:,:,:) )   ! FeL2
357            IF( iom_use("FeP")  ) CALL iom_put("FeP"    , zFeP   (:,:,:)       * tmask(:,:,:) )   ! FeP
358            IF( iom_use("TL2")  ) CALL iom_put("TL2"    , zTL2   (:,:,:)       * tmask(:,:,:) )   ! TL2
359         ENDIF
360         ENDIF
361      ENDIF
362
363      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
364         WRITE(charout, FMT="('fechem')")
365         CALL prt_ctl_trc_info(charout)
366         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
367      ENDIF
368      !
369      IF( ln_fechem )  THEN
370         DEALLOCATE( zstrn, zstrn2 )
371         DEALLOCATE( zFe2, zFeL2, zTL2, zFeP )
372      ENDIF
373      !
374      IF( ln_timing )   CALL timing_stop('p4z_fechem')
375      !
376   END SUBROUTINE p4z_fechem
377
378
379   SUBROUTINE p4z_fechem_init
380      !!----------------------------------------------------------------------
381      !!                  ***  ROUTINE p4z_fechem_init  ***
382      !!
383      !! ** Purpose :   Initialization of iron chemistry parameters
384      !!
385      !! ** Method  :   Read the nampisfer namelist and check the parameters
386      !!      called at the first timestep
387      !!
388      !! ** input   :   Namelist nampisfer
389      !!
390      !!----------------------------------------------------------------------
391      INTEGER ::   ios   ! Local integer
392      !!
393      NAMELIST/nampisfer/ ln_fechem, ln_ligvar, xlam1, xlamdust, ligand, kfep 
394      !!----------------------------------------------------------------------
395      !
396      IF(lwp) THEN
397         WRITE(numout,*)
398         WRITE(numout,*) 'p4z_rem_init : Initialization of iron chemistry parameters'
399         WRITE(numout,*) '~~~~~~~~~~~~'
400      ENDIF
401      !
402      REWIND( numnatp_ref )            ! Namelist nampisfer in reference namelist : Pisces iron chemistry
403      READ  ( numnatp_ref, nampisfer, IOSTAT = ios, ERR = 901)
404901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nampisfer in reference namelist', lwp )
405      REWIND( numnatp_cfg )            ! Namelist nampisfer in configuration namelist : Pisces iron chemistry
406      READ  ( numnatp_cfg, nampisfer, IOSTAT = ios, ERR = 902 )
407902   IF( ios >  0 )   CALL ctl_nam ( ios , 'nampisfer in configuration namelist', lwp )
408      IF(lwm) WRITE( numonp, nampisfer )
409
410      IF(lwp) THEN                     ! control print
411         WRITE(numout,*) '   Namelist : nampisfer'
412         WRITE(numout,*) '      enable complex iron chemistry scheme      ln_fechem    =', ln_fechem
413         WRITE(numout,*) '      variable concentration of ligand          ln_ligvar    =', ln_ligvar
414         WRITE(numout,*) '      scavenging rate of Iron                   xlam1        =', xlam1
415         WRITE(numout,*) '      scavenging rate of Iron by dust           xlamdust     =', xlamdust
416         WRITE(numout,*) '      ligand concentration in the ocean         ligand       =', ligand
417         WRITE(numout,*) '      rate constant for nanoparticle formation  kfep         =', kfep
418      ENDIF
419      !
420      IF (ln_ligand .AND. ln_fechem) CALL ctl_stop( 'STOP', 'p4z_fechem_init: ln_ligand and ln_fechem are incompatible')
421      !
422      IF( ln_fechem ) THEN             ! set some constants used by the complexe chemistry scheme
423         !
424         spd = 3600. * 24.
425         con = 1.E9
426         ! LIGAND KINETICS (values from Witter et al. 2000)
427         ! Weak (L2) ligands
428         ! Phaeophytin
429         kl2 = 12.2E5  * spd / con
430         kb2 = 12.3E-6 * spd
431         ! Strong (L1) ligands
432         ! Saccharides
433         ! kl1 = 12.2E5  * spd / con
434         ! kb1 = 12.3E-6 * spd
435         ! DFOB-
436         kl1 = 19.6e5  * spd / con
437         kb1 = 1.5e-6  * spd
438         ! pcp and remin of Fe3p
439         ks  = 0.075
440         kpr = 0.05
441         ! thermal reduction of Fe3
442         kth = 0.0048 * 24.
443         !
444      ENDIF
445      !
446   END SUBROUTINE p4z_fechem_init
447   
448   !!======================================================================
449END MODULE p4zfechem
Note: See TracBrowser for help on using the repository browser.