source: NEMO/branches/2019/dev_r12072_TOP-01_ENHANCE-11_CEthe/src/TOP/PISCES/P4Z/p4zsed.F90 @ 12110

Last change on this file since 12110 was 12110, checked in by cetlod, 8 months ago

merge dev_r11219_TOP-01_cethe_PISCES_LBC onto dev_r12072_TOP-01_ENHANCE-11_CEthe

  • Property svn:keywords set to Id
File size: 21.0 KB
Line 
1MODULE p4zsed
2   !!======================================================================
3   !!                         ***  MODULE p4sed  ***
4   !! TOP :   PISCES Compute loss of organic matter in the sediments
5   !!======================================================================
6   !! History :   1.0  !  2004-03 (O. Aumont) Original code
7   !!             2.0  !  2007-12 (C. Ethe, G. Madec)  F90
8   !!             3.4  !  2011-06 (C. Ethe) USE of fldread
9   !!             3.5  !  2012-07 (O. Aumont) improvment of river input of nutrients
10   !!----------------------------------------------------------------------
11   !!   p4z_sed        :  Compute loss of organic matter in the sediments
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 p4zlim          !  Co-limitations of differents nutrients
17   USE p4zint          !  interpolation and computation of various fields
18   USE sed             !  Sediment module
19   USE iom             !  I/O manager
20   USE prtctl_trc      !  print control for debugging
21
22   IMPLICIT NONE
23   PRIVATE
24
25   PUBLIC   p4z_sed 
26   PUBLIC   p4z_sed_init
27   PUBLIC   p4z_sed_alloc
28 
29   REAL(wp), PUBLIC ::   nitrfix      !: Nitrogen fixation rate
30   REAL(wp), PUBLIC ::   diazolight   !: Nitrogen fixation sensitivty to light
31   REAL(wp), PUBLIC ::   concfediaz   !: Fe half-saturation Cste for diazotrophs
32
33   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: nitrpot    !: Nitrogen fixation
34   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:  ) :: sdenit     !: Nitrate reduction in the sediments
35   !
36   REAL(wp), SAVE :: r1_rday         
37   REAL(wp), SAVE :: sedsilfrac, sedcalfrac
38
39   !!----------------------------------------------------------------------
40   !! NEMO/TOP 4.0 , NEMO Consortium (2018)
41   !! $Id$
42   !! Software governed by the CeCILL license (see ./LICENSE)
43   !!----------------------------------------------------------------------
44CONTAINS
45
46   SUBROUTINE p4z_sed( kt, knt )
47      !!---------------------------------------------------------------------
48      !!                     ***  ROUTINE p4z_sed  ***
49      !!
50      !! ** Purpose :   Compute loss of organic matter in the sediments. This
51      !!              is by no way a sediment model. The loss is simply
52      !!              computed to balance the inout from rivers and dust
53      !!
54      !! ** Method  : - ???
55      !!---------------------------------------------------------------------
56      !
57      INTEGER, INTENT(in) ::   kt, knt ! ocean time step
58      INTEGER  ::  ji, jj, jk, ikt
59      REAL(wp) ::  zrivalk, zrivsil, zrivno3
60      REAL(wp) ::  zlim, zfact, zfactcal
61      REAL(wp) ::  zo2, zno3, zflx, zpdenit, z1pdenit, zolimit
62      REAL(wp) ::  zsiloss, zcaloss, zws3, zws4, zwsc, zdep
63      REAL(wp) ::  zwstpoc, zwstpon, zwstpop
64      REAL(wp) ::  ztrfer, ztrpo4s, ztrdp, zwdust, zmudia, ztemp
65      REAL(wp) ::  xdiano3, xdianh4
66      !
67      CHARACTER (len=25) :: charout
68      REAL(wp), DIMENSION(jpi,jpj    ) :: zdenit2d, zbureff, zwork
69      REAL(wp), DIMENSION(jpi,jpj    ) :: zwsbio3, zwsbio4
70      REAL(wp), DIMENSION(jpi,jpj    ) :: zsedcal, zsedsi, zsedc
71      REAL(wp), DIMENSION(jpi,jpj,jpk) :: zsoufer, zlight
72      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: ztrpo4, ztrdop, zirondep, zpdep
73      !!---------------------------------------------------------------------
74      !
75      IF( ln_timing )  CALL timing_start('p4z_sed')
76      !
77
78      ! Allocate temporary workspace
79      ALLOCATE( ztrpo4(jpi,jpj,jpk) )
80      IF( ln_p5z )    ALLOCATE( ztrdop(jpi,jpj,jpk) )
81
82      zdenit2d(:,:) = 0.e0
83      zbureff (:,:) = 0.e0
84      zwork   (:,:) = 0.e0
85      zsedsi  (:,:) = 0.e0
86      zsedcal (:,:) = 0.e0
87      zsedc   (:,:) = 0.e0
88
89      IF( .NOT.lk_sed ) THEN
90         ! OA: Warning, the following part is necessary to avoid CFL problems above the sediments
91         ! --------------------------------------------------------------------
92         DO jj = 1, jpj
93            DO ji = 1, jpi
94               ikt  = mbkt(ji,jj)
95               zdep = e3t_n(ji,jj,ikt) / xstep
96               zwsbio4(ji,jj) = MIN( 0.99 * zdep, wsbio4(ji,jj,ikt) )
97               zwsbio3(ji,jj) = MIN( 0.99 * zdep, wsbio3(ji,jj,ikt) )
98            END DO
99         END DO
100
101         ! Computation of the sediment denitrification proportion: The metamodel from midlleburg (2006) is being used
102         ! Computation of the fraction of organic matter that is permanently buried from Dunne's model
103         ! -------------------------------------------------------
104         DO jj = 1, jpj
105            DO ji = 1, jpi
106              IF( tmask(ji,jj,1) == 1 ) THEN
107                 ikt = mbkt(ji,jj)
108                 zflx = (  trb(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj)   &
109                   &     + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) )  * 1E3 * 1E6 / 1E4
110                 zflx  = LOG10( MAX( 1E-3, zflx ) )
111                 zo2   = LOG10( MAX( 10. , trb(ji,jj,ikt,jpoxy) * 1E6 ) )
112                 zno3  = LOG10( MAX( 1.  , trb(ji,jj,ikt,jpno3) * 1E6 * rno3 ) )
113                 zdep  = LOG10( gdepw_n(ji,jj,ikt+1) )
114                 zdenit2d(ji,jj) = -2.2567 - 1.185 * zflx - 0.221 * zflx**2 - 0.3995 * zno3 * zo2 + 1.25 * zno3    &
115                   &                + 0.4721 * zo2 - 0.0996 * zdep + 0.4256 * zflx * zo2
116                 zdenit2d(ji,jj) = 10.0**( zdenit2d(ji,jj) )
117                   !
118                 zflx = (  trb(ji,jj,ikt,jpgoc) * zwsbio4(ji,jj)   &
119                   &     + trb(ji,jj,ikt,jppoc) * zwsbio3(ji,jj) ) * 1E6
120                 zbureff(ji,jj) = 0.013 + 0.53 * zflx**2 / ( 7.0 + zflx )**2
121              ENDIF
122            END DO
123         END DO 
124         !
125      ENDIF
126
127      ! This loss is scaled at each bottom grid cell for equilibrating the total budget of silica in the ocean.
128      ! Thus, the amount of silica lost in the sediments equal the supply at the surface (dust+rivers)
129      ! ------------------------------------------------------
130      IF( .NOT.lk_sed )  zrivsil = 1._wp - sedsilfrac
131
132      DO jj = 1, jpj
133         DO ji = 1, jpi
134            ikt  = mbkt(ji,jj)
135            zdep = xstep / e3t_n(ji,jj,ikt) 
136            zwsc = zwsbio4(ji,jj) * zdep
137            zsiloss = trb(ji,jj,ikt,jpgsi) * zwsc
138            zcaloss = trb(ji,jj,ikt,jpcal) * zwsc
139            !
140            tra(ji,jj,ikt,jpgsi) = tra(ji,jj,ikt,jpgsi) - zsiloss
141            tra(ji,jj,ikt,jpcal) = tra(ji,jj,ikt,jpcal) - zcaloss
142         END DO
143      END DO
144      !
145      IF( .NOT.lk_sed ) THEN
146         DO jj = 1, jpj
147            DO ji = 1, jpi
148               ikt  = mbkt(ji,jj)
149               zdep = xstep / e3t_n(ji,jj,ikt) 
150               zwsc = zwsbio4(ji,jj) * zdep
151               zsiloss = trb(ji,jj,ikt,jpgsi) * zwsc
152               zcaloss = trb(ji,jj,ikt,jpcal) * zwsc
153               tra(ji,jj,ikt,jpsil) = tra(ji,jj,ikt,jpsil) + zsiloss * zrivsil 
154               !
155               zfactcal = MIN( excess(ji,jj,ikt), 0.2 )
156               zfactcal = MIN( 1., 1.3 * ( 0.2 - zfactcal ) / ( 0.4 - zfactcal ) )
157               zrivalk  = sedcalfrac * zfactcal
158               tra(ji,jj,ikt,jptal) =  tra(ji,jj,ikt,jptal) + zcaloss * zrivalk * 2.0
159               tra(ji,jj,ikt,jpdic) =  tra(ji,jj,ikt,jpdic) + zcaloss * zrivalk
160               zsedcal(ji,jj) = (1.0 - zrivalk) * zcaloss * e3t_n(ji,jj,ikt) 
161               zsedsi (ji,jj) = (1.0 - zrivsil) * zsiloss * e3t_n(ji,jj,ikt) 
162            END DO
163         END DO
164      ENDIF
165      !
166      DO jj = 1, jpj
167         DO ji = 1, jpi
168            ikt  = mbkt(ji,jj)
169            zdep = xstep / e3t_n(ji,jj,ikt) 
170            zws4 = zwsbio4(ji,jj) * zdep
171            zws3 = zwsbio3(ji,jj) * zdep
172            tra(ji,jj,ikt,jpgoc) = tra(ji,jj,ikt,jpgoc) - trb(ji,jj,ikt,jpgoc) * zws4 
173            tra(ji,jj,ikt,jppoc) = tra(ji,jj,ikt,jppoc) - trb(ji,jj,ikt,jppoc) * zws3
174            tra(ji,jj,ikt,jpbfe) = tra(ji,jj,ikt,jpbfe) - trb(ji,jj,ikt,jpbfe) * zws4
175            tra(ji,jj,ikt,jpsfe) = tra(ji,jj,ikt,jpsfe) - trb(ji,jj,ikt,jpsfe) * zws3
176         END DO
177      END DO
178      !
179      IF( ln_p5z ) THEN
180         DO jj = 1, jpj
181            DO ji = 1, jpi
182               ikt  = mbkt(ji,jj)
183               zdep = xstep / e3t_n(ji,jj,ikt) 
184               zws4 = zwsbio4(ji,jj) * zdep
185               zws3 = zwsbio3(ji,jj) * zdep
186               tra(ji,jj,ikt,jpgon) = tra(ji,jj,ikt,jpgon) - trb(ji,jj,ikt,jpgon) * zws4
187               tra(ji,jj,ikt,jppon) = tra(ji,jj,ikt,jppon) - trb(ji,jj,ikt,jppon) * zws3
188               tra(ji,jj,ikt,jpgop) = tra(ji,jj,ikt,jpgop) - trb(ji,jj,ikt,jpgop) * zws4
189               tra(ji,jj,ikt,jppop) = tra(ji,jj,ikt,jppop) - trb(ji,jj,ikt,jppop) * zws3
190            END DO
191         END DO
192      ENDIF
193
194      IF( .NOT.lk_sed ) THEN
195         ! The 0.5 factor in zpdenit is to avoid negative NO3 concentration after
196         ! denitrification in the sediments. Not very clever, but simpliest option.
197         DO jj = 1, jpj
198            DO ji = 1, jpi
199               ikt  = mbkt(ji,jj)
200               zdep = xstep / e3t_n(ji,jj,ikt) 
201               zws4 = zwsbio4(ji,jj) * zdep
202               zws3 = zwsbio3(ji,jj) * zdep
203               zrivno3 = 1. - zbureff(ji,jj)
204               zwstpoc = trb(ji,jj,ikt,jpgoc) * zws4 + trb(ji,jj,ikt,jppoc) * zws3
205               zpdenit  = MIN( 0.5 * ( trb(ji,jj,ikt,jpno3) - rtrn ) / rdenit, zdenit2d(ji,jj) * zwstpoc * zrivno3 )
206               z1pdenit = zwstpoc * zrivno3 - zpdenit
207               zolimit = MIN( ( trb(ji,jj,ikt,jpoxy) - rtrn ) / o2ut, z1pdenit * ( 1.- nitrfac(ji,jj,ikt) ) )
208               tra(ji,jj,ikt,jpdoc) = tra(ji,jj,ikt,jpdoc) + z1pdenit - zolimit
209               tra(ji,jj,ikt,jppo4) = tra(ji,jj,ikt,jppo4) + zpdenit + zolimit
210               tra(ji,jj,ikt,jpnh4) = tra(ji,jj,ikt,jpnh4) + zpdenit + zolimit
211               tra(ji,jj,ikt,jpno3) = tra(ji,jj,ikt,jpno3) - rdenit * zpdenit
212               tra(ji,jj,ikt,jpoxy) = tra(ji,jj,ikt,jpoxy) - zolimit * o2ut
213               tra(ji,jj,ikt,jptal) = tra(ji,jj,ikt,jptal) + rno3 * (zolimit + (1.+rdenit) * zpdenit )
214               tra(ji,jj,ikt,jpdic) = tra(ji,jj,ikt,jpdic) + zpdenit + zolimit 
215               sdenit(ji,jj) = rdenit * zpdenit * e3t_n(ji,jj,ikt)
216               zsedc(ji,jj)   = (1. - zrivno3) * zwstpoc * e3t_n(ji,jj,ikt)
217               IF( ln_p5z ) THEN
218                  zwstpop              = trb(ji,jj,ikt,jpgop) * zws4 + trb(ji,jj,ikt,jppop) * zws3
219                  zwstpon              = trb(ji,jj,ikt,jpgon) * zws4 + trb(ji,jj,ikt,jppon) * zws3
220                  tra(ji,jj,ikt,jpdon) = tra(ji,jj,ikt,jpdon) + ( z1pdenit - zolimit ) * zwstpon / (zwstpoc + rtrn)
221                  tra(ji,jj,ikt,jpdop) = tra(ji,jj,ikt,jpdop) + ( z1pdenit - zolimit ) * zwstpop / (zwstpoc + rtrn)
222               ENDIF
223            END DO
224         END DO
225       ENDIF
226
227
228      ! Nitrogen fixation process
229      ! Small source iron from particulate inorganic iron
230      !-----------------------------------
231      DO jk = 1, jpkm1
232         zlight (:,:,jk) =  ( 1.- EXP( -etot_ndcy(:,:,jk) / diazolight ) ) * ( 1. - fr_i(:,:) ) 
233         zsoufer(:,:,jk) = zlight(:,:,jk) * 2E-11 / ( 2E-11 + biron(:,:,jk) )
234      ENDDO
235      IF( ln_p4z ) THEN
236         DO jk = 1, jpkm1
237            DO jj = 1, jpj
238               DO ji = 1, jpi
239                  !                      ! Potential nitrogen fixation dependant on temperature and iron
240                  ztemp = tsn(ji,jj,jk,jp_tem)
241                  zmudia = MAX( 0.,-0.001096*ztemp**2 + 0.057*ztemp -0.637 ) * 7.625
242                  !       Potential nitrogen fixation dependant on temperature and iron
243                  xdianh4 = trb(ji,jj,jk,jpnh4) / ( concnnh4 + trb(ji,jj,jk,jpnh4) )
244                  xdiano3 = trb(ji,jj,jk,jpno3) / ( concnno3 + trb(ji,jj,jk,jpno3) ) * (1. - xdianh4)
245                  zlim = ( 1.- xdiano3 - xdianh4 )
246                  IF( zlim <= 0.1 )   zlim = 0.01
247                  zfact = zlim * rfact2
248                  ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) )
249                  ztrpo4(ji,jj,jk) = trb(ji,jj,jk,jppo4) / ( 1E-6 + trb(ji,jj,jk,jppo4) )
250                  ztrdp = ztrpo4(ji,jj,jk)
251                  nitrpot(ji,jj,jk) =  zmudia * r1_rday * zfact * MIN( ztrfer, ztrdp ) * zlight(ji,jj,jk)
252               END DO
253            END DO
254         END DO
255      ELSE       ! p5z
256         DO jk = 1, jpkm1
257            DO jj = 1, jpj
258               DO ji = 1, jpi
259                  !                      ! Potential nitrogen fixation dependant on temperature and iron
260                  ztemp = tsn(ji,jj,jk,jp_tem)
261                  zmudia = MAX( 0.,-0.001096*ztemp**2 + 0.057*ztemp -0.637 ) * 7.625
262                  !       Potential nitrogen fixation dependant on temperature and iron
263                  xdianh4 = trb(ji,jj,jk,jpnh4) / ( concnnh4 + trb(ji,jj,jk,jpnh4) )
264                  xdiano3 = trb(ji,jj,jk,jpno3) / ( concnno3 + trb(ji,jj,jk,jpno3) ) * (1. - xdianh4)
265                  zlim = ( 1.- xdiano3 - xdianh4 )
266                  IF( zlim <= 0.1 )   zlim = 0.01
267                  zfact = zlim * rfact2
268                  ztrfer = biron(ji,jj,jk) / ( concfediaz + biron(ji,jj,jk) )
269                  ztrpo4(ji,jj,jk) = trb(ji,jj,jk,jppo4) / ( 1E-6 + trb(ji,jj,jk,jppo4) )
270                  ztrdop(ji,jj,jk) = trb(ji,jj,jk,jpdop) / ( 1E-6 + trb(ji,jj,jk,jpdop) ) * (1. - ztrpo4(ji,jj,jk))
271                  ztrdp = ztrpo4(ji,jj,jk) + ztrdop(ji,jj,jk)
272                  nitrpot(ji,jj,jk) =  zmudia * r1_rday * zfact * MIN( ztrfer, ztrdp ) * zlight(ji,jj,jk)
273               END DO
274            END DO
275         END DO
276      ENDIF
277
278      ! Nitrogen change due to nitrogen fixation
279      ! ----------------------------------------
280      IF( ln_p4z ) THEN
281         DO jk = 1, jpkm1
282            DO jj = 1, jpj
283               DO ji = 1, jpi
284                  zfact = nitrpot(ji,jj,jk) * nitrfix
285                  tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zfact / 3.0
286                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * zfact / 3.0
287                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) - zfact * 2.0 / 3.0
288                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zfact * 1.0 / 3.0
289                  tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zfact * 1.0 / 3.0 * 2.0 / 3.0
290                  tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zfact * 1.0 / 3.0 * 1.0 / 3.0
291                  tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + ( o2ut + o2nit ) * zfact * 2.0 / 3.0 + o2nit * zfact / 3.0
292                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - 30E-6 * zfact * 1.0 / 3.0
293                  tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + 30E-6 * zfact * 1.0 / 3.0 * 2.0 / 3.0
294                  tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + 30E-6 * zfact * 1.0 / 3.0 * 1.0 / 3.0
295                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * rfact2 / rday
296                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) + concdnh4 / ( concdnh4 + trb(ji,jj,jk,jppo4) ) &
297                  &                     * 0.001 * trb(ji,jj,jk,jpdoc) * xstep
298              END DO
299            END DO
300         END DO
301      ELSE    ! p5z
302         DO jk = 1, jpkm1
303            DO jj = 1, jpj
304               DO ji = 1, jpi
305                  zfact = nitrpot(ji,jj,jk) * nitrfix
306                  tra(ji,jj,jk,jpnh4) = tra(ji,jj,jk,jpnh4) + zfact / 3.0
307                  tra(ji,jj,jk,jptal) = tra(ji,jj,jk,jptal) + rno3 * zfact / 3.0
308                  tra(ji,jj,jk,jppo4) = tra(ji,jj,jk,jppo4) - 16.0 / 46.0 * zfact * ( 1.0 - 1.0 / 3.0 ) &
309                  &                     * ztrpo4(ji,jj,jk) / (ztrpo4(ji,jj,jk) + ztrdop(ji,jj,jk) + rtrn)
310                  tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) + zfact * 1.0 / 3.0
311                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) + zfact * 1.0 / 3.0
312                  tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) + 16.0 / 46.0 * zfact / 3.0  &
313                  &                     - 16.0 / 46.0 * zfact * ztrdop(ji,jj,jk)   &
314                  &                     / (ztrpo4(ji,jj,jk) + ztrdop(ji,jj,jk) + rtrn)
315                  tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) + zfact * 1.0 / 3.0 * 2.0 / 3.0
316                  tra(ji,jj,jk,jppon) = tra(ji,jj,jk,jppon) + zfact * 1.0 / 3.0 * 2.0 /3.0
317                  tra(ji,jj,jk,jppop) = tra(ji,jj,jk,jppop) + 16.0 / 46.0 * zfact * 1.0 / 3.0 * 2.0 /3.0
318                  tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zfact * 1.0 / 3.0 * 1.0 / 3.0
319                  tra(ji,jj,jk,jpgon) = tra(ji,jj,jk,jpgon) + zfact * 1.0 / 3.0 * 1.0 /3.0
320                  tra(ji,jj,jk,jpgop) = tra(ji,jj,jk,jpgop) + 16.0 / 46.0 * zfact * 1.0 / 3.0 * 1.0 /3.0
321                  tra(ji,jj,jk,jpoxy) = tra(ji,jj,jk,jpoxy) + ( o2ut + o2nit ) * zfact * 2.0 / 3.0 + o2nit * zfact / 3.0
322                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) - 30E-6 * zfact * 1.0 / 3.0 
323                  tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) + 30E-6 * zfact * 1.0 / 3.0 * 2.0 / 3.0
324                  tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + 30E-6 * zfact * 1.0 / 3.0 * 1.0 / 3.0
325                  tra(ji,jj,jk,jpfer) = tra(ji,jj,jk,jpfer) + 0.002 * 4E-10 * zsoufer(ji,jj,jk) * rfact2 / rday
326              END DO
327            END DO
328         END DO
329         !
330      ENDIF
331
332      IF( lk_iomput ) THEN
333         IF( knt == nrdttrc ) THEN
334            zfact = 1.e+3 * rfact2r !  conversion from molC/l/kt  to molN/m3/s
335            IF( iom_use("Nfix"   ) ) CALL iom_put( "Nfix", nitrpot(:,:,:) * nitrfix * rno3 * zfact * tmask(:,:,:) )  ! nitrogen fixation
336            IF( iom_use("INTNFIX") ) THEN   ! nitrogen fixation rate in ocean ( vertically integrated )
337               zwork(:,:) = 0.
338               DO jk = 1, jpkm1
339                 zwork(:,:) = zwork(:,:) + nitrpot(:,:,jk) * nitrfix * rno3 * zfact * e3t_n(:,:,jk) * tmask(:,:,jk)
340               ENDDO
341               CALL iom_put( "INTNFIX" , zwork ) 
342            ENDIF
343            IF( iom_use("SedCal" ) ) CALL iom_put( "SedCal", zsedcal(:,:) * zfact )
344            IF( iom_use("SedSi" ) )  CALL iom_put( "SedSi",  zsedsi (:,:) * zfact )
345            IF( iom_use("SedC" ) )   CALL iom_put( "SedC",   zsedc  (:,:) * zfact )
346            IF( iom_use("Sdenit" ) ) CALL iom_put( "Sdenit", sdenit (:,:) * zfact * rno3 )
347         ENDIF
348      ENDIF
349      !
350      IF(ln_ctl) THEN  ! print mean trends (USEd for debugging)
351         WRITE(charout, fmt="('sed ')")
352         CALL prt_ctl_trc_info(charout)
353         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
354      ENDIF
355      !
356      IF( ln_p5z )    DEALLOCATE( ztrpo4, ztrdop )
357      !
358      IF( ln_timing )  CALL timing_stop('p4z_sed')
359      !
360   END SUBROUTINE p4z_sed
361
362   SUBROUTINE p4z_sed_init
363      !!----------------------------------------------------------------------
364      !!                  ***  routine p4z_sed_init  ***
365      !!
366      !! ** purpose :   initialization of some parameters
367      !!
368      !!----------------------------------------------------------------------
369      !!----------------------------------------------------------------------
370      INTEGER  :: ji, jj, jk, jm
371      INTEGER  :: ios                 ! Local integer output status for namelist read
372      !
373      !!
374      NAMELIST/nampissed/ nitrfix, diazolight, concfediaz
375      !!----------------------------------------------------------------------
376      !
377      IF(lwp) THEN
378         WRITE(numout,*)
379         WRITE(numout,*) 'p4z_sed_init : initialization of sediment mobilisation '
380         WRITE(numout,*) '~~~~~~~~~~~~ '
381      ENDIF
382      !                            !* set file information
383      REWIND( numnatp_ref )              ! Namelist nampissed in reference namelist : PISCES sediment
384      READ  ( numnatp_ref, nampissed, IOSTAT = ios, ERR = 901)
385901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nampissed in reference namelist', lwp )
386      REWIND( numnatp_cfg )              ! Namelist nampissbc in configuration namelist : PISCES sediment
387      READ  ( numnatp_cfg, nampissed, IOSTAT = ios, ERR = 902 )
388902   IF( ios >  0 )   CALL ctl_nam ( ios , 'nampissed in configuration namelist', lwp )
389      IF(lwm) WRITE ( numonp, nampissed )
390
391      IF(lwp) THEN
392         WRITE(numout,*) '   Namelist : nampissed '
393         WRITE(numout,*) '      nitrogen fixation rate                       nitrfix = ', nitrfix
394         WRITE(numout,*) '      nitrogen fixation sensitivty to light    diazolight  = ', diazolight
395         WRITE(numout,*) '      Fe half-saturation cste for diazotrophs  concfediaz  = ', concfediaz
396      ENDIF
397      !
398      r1_rday  = 1. / rday
399      !
400      sedsilfrac = 0.03     ! percentage of silica loss in the sediments
401      sedcalfrac = 0.6      ! percentage of calcite loss in the sediments
402      !
403      lk_sed = ln_sediment .AND. ln_sed_2way 
404      !
405   END SUBROUTINE p4z_sed_init
406
407   INTEGER FUNCTION p4z_sed_alloc()
408      !!----------------------------------------------------------------------
409      !!                     ***  ROUTINE p4z_sed_alloc  ***
410      !!----------------------------------------------------------------------
411      ALLOCATE( nitrpot(jpi,jpj,jpk), sdenit(jpi,jpj), STAT=p4z_sed_alloc )
412      !
413      IF( p4z_sed_alloc /= 0 )   CALL ctl_stop( 'STOP', 'p4z_sed_alloc: failed to allocate arrays' )
414      !
415   END FUNCTION p4z_sed_alloc
416
417   !!======================================================================
418END MODULE p4zsed
Note: See TracBrowser for help on using the repository browser.