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.
p4zfechem.F90 in branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z – NEMO

source: branches/CNRS/dev_r6270_PISCES_QUOTA/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zfechem.F90 @ 8003

Last change on this file since 8003 was 8003, checked in by aumont, 7 years ago

modification in the code to remove unnecessary parts such as kriest and non iomput options

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