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.
p4zsed.F90 in trunk/NEMO/TOP_SRC/PISCES – NEMO

source: trunk/NEMO/TOP_SRC/PISCES/p4zsed.F90 @ 1119

Last change on this file since 1119 was 1119, checked in by cetlod, 16 years ago

style of all top namelist has been modified ; update modules to take it into account, see ticket:196
  • Property svn:executable set to *
File size: 23.2 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   !!----------------------------------------------------------------------
9#if defined key_pisces
10   !!----------------------------------------------------------------------
11   !!   'key_pisces'                                       PISCES bio-model
12   !!----------------------------------------------------------------------
13   !!   p4z_sed        :  Compute loss of organic matter in the sediments
14   !!   p4z_sbc        :  Read and interpolate time-varying nutrients fluxes
15   !!   p4z_sed_init   :  Initialization of p4z_sed
16   !!----------------------------------------------------------------------
17   USE trc
18   USE oce_trc         !
19   USE sms_pisces
20   USE lib_mpp
21   USE prtctl_trc
22   USE p4zbio
23   USE p4zint
24   USE p4zopt
25   USE p4zsink
26   USE p4zrem
27   USE p4zlim
28   USE lbclnk
29   USE iom
30
31
32   IMPLICIT NONE
33   PRIVATE
34
35   PUBLIC   p4z_sed   
36
37   !! * Shared module variables
38   LOGICAL, PUBLIC ::    &
39     bdustfer  = .FALSE.      ,  &  !:
40     briver    = .FALSE.      ,  &  !:
41     bndepo    = .FALSE.      ,  &  !:
42     bsedinput = .FALSE.            !:
43
44   REAL(wp), PUBLIC ::   &
45     sedfeinput = 1.E-9_wp   ,  &  !:
46     dustsolub  = 0.014_wp         !:
47
48   !! * Module variables
49   INTEGER ::          &
50      numdust,          &  ! logical unit for surface fluxes data
51      nflx1 , nflx2,   &  !  first and second record used
52      nflx11, nflx12      ! ???
53   REAL(wp), DIMENSION(jpi,jpj,2) ::    &  !:
54     dustmo                                !: 2 consecutive set of dust fields
55   REAL(wp), DIMENSION(jpi,jpj)   ::    &
56     rivinp, cotdep, nitdep, dust
57   REAL(wp), DIMENSION(jpi,jpj,jpk)  ::   &
58     ironsed
59   REAL(wp) :: sumdepsi, rivalkinput, rivpo4input, nitdepinput
60
61   !!* Substitution
62#  include "domzgr_substitute.h90"
63   !!----------------------------------------------------------------------
64   !! NEMO/TOP 2.0 , LOCEAN-IPSL (2007)
65   !! $Header:$
66   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
67   !!----------------------------------------------------------------------
68
69CONTAINS
70
71   SUBROUTINE p4z_sed(kt, jnt)
72      !!---------------------------------------------------------------------
73      !!                     ***  ROUTINE p4z_sed  ***
74      !!
75      !! ** Purpose :   Compute loss of organic matter in the sediments. This
76      !!              is by no way a sediment model. The loss is simply
77      !!              computed to balance the inout from rivers and dust
78      !!
79      !! ** Method  : - ???
80      !!---------------------------------------------------------------------
81      INTEGER, INTENT(in) ::   kt, jnt ! ocean time step
82      INTEGER  ::   ji, jj, jk
83      INTEGER  ::   ikt
84      REAL(wp) ::   zsumsedsi, zsumsedpo4, zsumsedcal
85      REAL(wp) ::   zconctmp , zdenitot  , znitrpottot
86      REAL(wp) ::   zlim, zconctmp2, zstep, zfact
87      REAL(wp), DIMENSION(jpi,jpj)     ::   zsidep
88      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   znitrpot, zirondep
89      CHARACTER (len=25) :: charout
90      !!---------------------------------------------------------------------
91
92
93      IF( ( kt * jnt ) == nittrc000  )   CALL p4z_sed_init      ! Initialization (first time-step only)
94
95      IF( (jnt == 1) .and. (bdustfer) )  CALL p4z_sbc( kt )
96
97      zstep = rfact2 / rjjss      ! Time step duration for the biology
98
99      zirondep(:,:,:) = 0.e0          ! Initialisation of variables used to compute deposition
100      zsidep  (:,:)   = 0.e0
101
102      ! Iron and Si deposition at the surface
103      ! -------------------------------------
104
105      DO jj = 1, jpj
106         DO ji = 1, jpi
107            zirondep(ji,jj,1) = ( dustsolub * dust(ji,jj) / ( 55.85 * rmoss ) + 3.e-10 / raass )   &
108               &             * rfact2 / fse3t(ji,jj,1)
109            zsidep  (ji,jj)   = 8.8 * 0.075 * dust(ji,jj) * rfact2 / ( fse3t(ji,jj,1) * 28.1 * rmoss )
110         END DO
111      END DO
112
113      ! Iron solubilization of particles in the water column
114      ! ----------------------------------------------------
115
116      DO jk = 2, jpkm1
117         zirondep(:,:,jk) = dust(:,:) / ( 10. * 55.85 * rmoss ) * rfact2 * 1.e-4
118      END DO
119
120      ! Add the external input of nutrients, carbon and alkalinity
121      ! ----------------------------------------------------------
122
123      trn(:,:,1,jppo4) = trn(:,:,1,jppo4) + rivinp(:,:) * rfact2 
124      trn(:,:,1,jpno3) = trn(:,:,1,jpno3) + (rivinp(:,:) + nitdep(:,:)) * rfact2
125      trn(:,:,1,jpfer) = trn(:,:,1,jpfer) + rivinp(:,:) * 3.e-5 * rfact2
126      trn(:,:,1,jpsil) = trn(:,:,1,jpsil) + zsidep (:,:) + cotdep(:,:)   * rfact2 / 6.
127      trn(:,:,1,jpdic) = trn(:,:,1,jpdic) + rivinp(:,:) * 2.631 * rfact2
128      trn(:,:,1,jptal) = trn(:,:,1,jptal) + (cotdep(:,:) - rno3*(rivinp(:,:) +  nitdep(:,:) ) ) * rfact2
129
130
131      ! Add the external input of iron which is 3D distributed
132      ! (dust, river and sediment mobilization)
133      ! ------------------------------------------------------
134
135      DO jk = 1, jpkm1
136         trn(:,:,jk,jpfer) = trn(:,:,jk,jpfer)   &
137         &       + zirondep(:,:,jk) + ironsed(:,:,jk) * rfact2
138      END DO
139
140      ! Initialisation of variables used to compute Sinking Speed
141      ! ---------------------------------------------------------
142
143      zsumsedsi  = 0.e0
144      zsumsedpo4 = 0.e0
145      zsumsedcal = 0.e0
146
147      ! Loss of biogenic silicon, Caco3 organic carbon in the sediments.
148      ! First, the total loss is computed.
149      ! The factor for calcite comes from the alkalinity effect
150      ! -------------------------------------------------------------
151
152      DO jj = 1, jpj
153         DO ji = 1, jpi
154            ikt = MAX( mbathy(ji,jj)-1, 1 )
155            zfact = e1t(ji,jj) * e2t(ji,jj) / rjjss * tmask_i(ji,jj)
156
157# if defined key_kriest
158            zsumsedsi  = zsumsedsi  + zfact *  trn(ji,jj,ikt,jpdsi) * wscal (ji,jj,ikt)
159            zsumsedpo4 = zsumsedpo4 + zfact *  trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt)
160# else
161            zsumsedsi  = zsumsedsi  + zfact *  trn(ji,jj,ikt,jpdsi) * wsbio4(ji,jj,ikt)
162            zsumsedpo4 = zsumsedpo4 + zfact *( trn(ji,jj,ikt,jpgoc) * wsbio4(ji,jj,ikt)   &
163               &       + trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt) )
164# endif
165
166            zsumsedcal = zsumsedcal + zfact *  trn(ji,jj,ikt,jpcal) * wscal (ji,jj,ikt) * 2.e0
167
168         END DO
169      END DO
170
171      IF( lk_mpp ) THEN
172         CALL mpp_sum( zsumsedsi  )   ! sums over the global domain
173         CALL mpp_sum( zsumsedcal )   ! sums over the global domain
174         CALL mpp_sum( zsumsedpo4 )   ! sums over the global domain
175      ENDIF
176
177      ! Then this loss is scaled at each bottom grid cell for
178      ! equilibrating the total budget of silica in the ocean.
179      ! Thus, the amount of silica lost in the sediments equal
180      ! the supply at the surface (dust+rivers)
181      ! ------------------------------------------------------
182
183      DO jj = 1, jpj
184         DO ji = 1, jpi
185            ikt = MAX( mbathy(ji,jj) - 1, 1 )
186            zconctmp = trn(ji,jj,ikt,jpdsi) * zstep / fse3t(ji,jj,ikt)   &
187# if ! defined key_kriest
188     &             * wscal (ji,jj,ikt)
189# else
190     &             * wsbio4(ji,jj,ikt)
191# endif
192
193            trn(ji,jj,ikt,jpdsi) = trn(ji,jj,ikt,jpdsi) - zconctmp
194            trn(ji,jj,ikt,jpsil) = trn(ji,jj,ikt,jpsil) + zconctmp   &
195            &      * ( 1.- ( sumdepsi + rivalkinput / raass / 6. ) / zsumsedsi )
196         END DO
197      END DO
198
199      DO jj = 1, jpj
200         DO ji = 1, jpi
201            ikt = MAX( mbathy(ji,jj) - 1, 1 )
202            zconctmp = trn(ji,jj,ikt,jpcal) * wscal(ji,jj,ikt) * zstep / fse3t(ji,jj,ikt)
203            trn(ji,jj,ikt,jpcal) = trn(ji,jj,ikt,jpcal) - zconctmp
204            trn(ji,jj,ikt,jptal) = trn(ji,jj,ikt,jptal) + zconctmp   &
205               &   * ( 1.- ( rivalkinput / raass ) / zsumsedcal ) * 2.e0
206            trn(ji,jj,ikt,jpdic) = trn(ji,jj,ikt,jpdic) + zconctmp   &
207               &   * ( 1.- ( rivalkinput / raass ) / zsumsedcal )
208         END DO
209      END DO
210
211      DO jj = 1, jpj
212         DO ji = 1, jpi
213            ikt = MAX( mbathy(ji,jj) - 1, 1 )
214            zfact = zstep / fse3t(ji,jj,ikt)
215
216# if ! defined key_kriest
217            zconctmp  = trn(ji,jj,ikt,jpgoc)
218            zconctmp2 = trn(ji,jj,ikt,jppoc)
219            trn(ji,jj,ikt,jpgoc) = trn(ji,jj,ikt,jpgoc) - zconctmp  * wsbio4(ji,jj,ikt) * zfact
220            trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) - zconctmp2 * wsbio3(ji,jj,ikt) * zfact
221            trn(ji,jj,ikt,jpdoc) = trn(ji,jj,ikt,jpdoc)    &
222            &      + ( zconctmp  * wsbio4(ji,jj,ikt) + zconctmp2 * wsbio3(ji,jj,ikt) ) * zfact   &
223            &      * ( 1.- rivpo4input / (raass * zsumsedpo4 ) )
224            trn(ji,jj,ikt,jpbfe) = trn(ji,jj,ikt,jpbfe) - trn(ji,jj,ikt,jpbfe) * wsbio4(ji,jj,ikt) * zfact
225            trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) - trn(ji,jj,ikt,jpsfe) * wsbio3(ji,jj,ikt) * zfact
226
227# else
228
229            zconctmp  = trn(ji,jj,ikt,jpnum)
230            zconctmp2 = trn(ji,jj,ikt,jppoc)
231            trn(ji,jj,ikt,jpnum) = trn(ji,jj,ikt,jpnum)   &
232            &      - zconctmp  * wsbio4(ji,jj,ikt) * zfact
233            trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc)   &
234            &      - zconctmp2 * wsbio3(ji,jj,ikt) * zfact
235            trn(ji,jj,ikt,jpdoc) = trn(ji,jj,ikt,jpdoc)    &
236            &      + ( zconctmp2 * wsbio3(ji,jj,ikt) )   &
237            &      * zfact * ( 1.- rivpo4input / ( raass * zsumsedpo4 ) )
238            trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe)   &
239            &      - trn(ji,jj,ikt,jpsfe) * wsbio3(ji,jj,ikt) * zfact
240
241# endif
242         END DO
243      END DO
244
245      ! Nitrogen fixation (simple parameterization). The total gain
246      ! from nitrogen fixation is scaled to balance the loss by
247      ! denitrification
248      ! -------------------------------------------------------------
249
250      zdenitot = 0.e0
251      DO jk = 1, jpkm1
252         DO jj = 1,jpj
253            DO ji = 1,jpi
254               zdenitot = zdenitot + denitr(ji,jj,jk) * rdenit * e1t(ji,jj) * e2t(ji,jj)   &
255               &        *fse3t(ji,jj,jk) * tmask(ji,jj,jk) * tmask_i(ji,jj) * xnegtr(ji,jj,jk)
256            END DO
257         END DO
258      END DO
259
260      IF( lk_mpp )   CALL mpp_sum( zdenitot )      ! sum over the global domain
261
262      ! Potential nitrogen fication dependant on temperature and iron
263      ! -------------------------------------------------------------
264
265!CDIR NOVERRCHK
266      DO jk = 1, jpk
267!CDIR NOVERRCHK
268         DO jj = 1, jpj
269!CDIR NOVERRCHK
270            DO ji = 1, jpi
271               zlim = ( 1.- xnanono3(ji,jj,jk) - xnanonh4(ji,jj,jk) )
272               IF( zlim <= 0.2 )   zlim = 0.01
273               znitrpot(ji,jj,jk) = MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) / rjjss )   &
274# if defined key_off_degrad
275               &                  * facvol(ji,jj,jk)   &
276# endif
277               &                  * zlim * rfact2 * trn(ji,jj,jk,jpfer)   &
278               &                  / ( conc3 + trn(ji,jj,jk,jpfer) ) * ( 1.- EXP( -etot(ji,jj,jk) / 50.) )
279            END DO
280         END DO
281      END DO
282
283      znitrpottot = 0.e0
284      DO jk = 1, jpkm1
285         DO jj = 1, jpj
286            DO ji = 1, jpi
287               znitrpottot = znitrpottot + znitrpot(ji,jj,jk) * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,jk)   &
288                  &                                           * tmask(ji,jj,jk) * tmask_i(ji,jj) 
289            END DO
290         END DO
291      END DO
292
293      IF( lk_mpp )   CALL mpp_sum( znitrpottot )  ! sum over the global domain
294
295      ! Nitrogen change due to nitrogen fixation
296      ! ----------------------------------------
297
298      DO jk = 1, jpk
299         DO jj = 1, jpj
300            DO ji = 1, jpi
301# if ! defined key_c1d && ( defined key_orca_r4 || defined key_orca_r2 || defined key_orca_r05 || defined key_orca_r025 )
302!!             zfact = znitrpot(ji,jj,jk) * zdenitot / znitrpottot
303               zfact = znitrpot(ji,jj,jk) * 1.e-7
304# else
305               zfact = znitrpot(ji,jj,jk) * 1.e-7
306# endif
307               trn(ji,jj,jk,jpnh4) = trn(ji,jj,jk,jpnh4) + zfact
308               trn(ji,jj,jk,jpoxy) = trn(ji,jj,jk,jpoxy) + zfact   * o2nit
309               trn(ji,jj,jk,jppo4) = trn(ji,jj,jk,jppo4) + 30./ 46.* zfact
310            END DO
311         END DO
312      END DO
313
314# if defined key_trc_diaadd
315      DO jj = 1,jpj
316         DO ji = 1,jpi
317            trc2d(ji,jj,jp_pcs0_2d + 11) = zirondep(ji,jj,1) * 1.e+3 * rfact2r * fse3t(ji,jj,1)
318            trc2d(ji,jj,jp_pcs0_2d + 12) = znitrpot(ji,jj,1) * 1.e-7 * fse3t(ji,jj,1) * 1.e+3 / rfact2
319         END DO
320      END DO
321# endif
322      !
323       IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
324         WRITE(charout, FMT="('sed ')")
325         CALL prt_ctl_trc_info(charout)
326         CALL prt_ctl_trc(tab4d=trn, mask=tmask, clinfo=ctrcnm)
327       ENDIF
328
329   END SUBROUTINE p4z_sed
330
331   SUBROUTINE p4z_sbc(kt)
332
333      !!----------------------------------------------------------------------
334      !!                  ***  ROUTINE p4z_sbc  ***
335      !!
336      !! ** Purpose :   Read and interpolate the external sources of
337      !!                nutrients
338      !!
339      !! ** Method  :   Read the files and interpolate the appropriate variables
340      !!
341      !! ** input   :   external netcdf files
342      !!
343      !!----------------------------------------------------------------------
344      !! * arguments
345      INTEGER, INTENT( in  ) ::   kt   ! ocean time step
346
347      !! * Local declarations
348      INTEGER ::   &
349         imois, imois2,       &  ! temporary integers
350         i15  , iman             !    "          "
351      REAL(wp) ::   &
352         zxy                     !    "         "
353
354
355      !!---------------------------------------------------------------------
356
357      ! Initialization
358      ! --------------
359
360      i15 = INT( 2 * FLOAT( nday ) / ( FLOAT( nobis(nmonth) ) + 0.5 ) )
361      iman  = INT( raamo )
362      imois = nmonth + i15 - 1
363      IF( imois == 0 ) imois = iman
364      imois2 = nmonth
365
366      ! 1. first call kt=nit000
367      ! -----------------------
368
369      IF( kt == nit000 ) THEN
370         ! initializations
371         nflx1  = 0
372         nflx11 = 0
373         ! open the file
374         IF(lwp) THEN
375            WRITE(numout,*) ' '
376            WRITE(numout,*) ' **** Routine p4z_sbc'
377         ENDIF
378         CALL iom_open ( 'dust.orca.nc', numdust )
379      ENDIF
380
381
382     ! Read monthly file
383      ! ----------------
384
385      IF( kt == nit000 .OR. imois /= nflx1 ) THEN
386
387         ! Calendar computation
388
389         ! nflx1 number of the first file record used in the simulation
390         ! nflx2 number of the last  file record
391
392         nflx1 = imois
393         nflx2 = nflx1+1
394         nflx1 = MOD( nflx1, iman )
395         nflx2 = MOD( nflx2, iman )
396         IF( nflx1 == 0 )   nflx1 = iman
397         IF( nflx2 == 0 )   nflx2 = iman
398         IF(lwp) WRITE(numout,*) 'first record file used nflx1 ',nflx1
399         IF(lwp) WRITE(numout,*) 'last  record file used nflx2 ',nflx2
400
401         ! Read monthly fluxes data
402
403         ! humidity
404         CALL iom_get ( numdust, jpdom_data, 'dust', dustmo(:,:,1), nflx1 )
405         CALL iom_get ( numdust, jpdom_data, 'dust', dustmo(:,:,2), nflx2 )
406
407         IF(lwp .AND. nitend-nit000 <= 100 ) THEN
408            WRITE(numout,*)
409            WRITE(numout,*) ' read clio flx ok'
410            WRITE(numout,*)
411               WRITE(numout,*)
412               WRITE(numout,*) 'Clio month: ',nflx1,'  field: dust'
413               CALL prihre( dustmo(:,:,1),jpi,jpj,1,jpi,20,1,jpj,10,1e9,numout )
414         ENDIF
415
416      ENDIF
417
418     ! 3. at every time step interpolation of fluxes
419      ! ---------------------------------------------
420
421      zxy = FLOAT( nday ) / FLOAT( nobis(nflx1) ) + 0.5 - i15
422
423      dust(:,:) = ( (1.-zxy) * dustmo(:,:,1) + zxy * dustmo(:,:,2) )
424
425      IF( kt == nitend ) CALL iom_close (numdust)
426
427   END SUBROUTINE p4z_sbc
428
429
430   SUBROUTINE p4z_sed_init
431
432      !!----------------------------------------------------------------------
433      !!                  ***  ROUTINE p4z_sed_init  ***
434      !!
435      !! ** Purpose :   Initialization of the external sources of nutrients
436      !!
437      !! ** Method  :   Read the files and compute the budget
438      !!      called at the first timestep (nittrc000)
439      !!
440      !! ** input   :   external netcdf files
441      !!
442      !!----------------------------------------------------------------------
443
444      INTEGER ::   ji, jj, jk, jm
445      INTEGER , PARAMETER ::   jpmois = 12, jpan = 1
446      INTEGER :: numriv, numbath, numdep
447
448
449      REAL(wp) ::   zcoef
450      REAL(wp) ::   expide, denitide,zmaskt
451      REAL(wp) , DIMENSION (jpi,jpj)     ::   riverdoc, river, ndepo
452      REAL(wp) , DIMENSION (jpi,jpj,jpk) ::   cmask
453      REAL(wp), DIMENSION(jpi,jpj,12)    ::   zdustmo
454
455      NAMELIST/nampissed/ bdustfer, briver, bndepo, bsedinput, sedfeinput, dustsolub
456
457
458      REWIND( numnat )                     ! read numnat
459      READ  ( numnat, nampissed )
460
461      IF(lwp) THEN
462         WRITE(numout,*) ' '
463         WRITE(numout,*) ' Namelist : nampissed '
464         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~ '
465         WRITE(numout,*) '    Dust input from the atmosphere           bdustfer  = ', bdustfer
466         WRITE(numout,*) '    River input of nutrients                 briver    = ', briver
467         WRITE(numout,*) '    Atmospheric deposition of N              bndepo    = ', bndepo
468         WRITE(numout,*) '    Fe input from sediments                  bsedinput = ', bsedinput
469         WRITE(numout,*) '    Coastal release of Iron                  sedfeinput =', sedfeinput
470         WRITE(numout,*) '    Solubility of the dust                   dustsolub  =', dustsolub
471      ENDIF
472
473      ! Dust input from the atmosphere
474      ! ------------------------------
475      IF( bdustfer ) THEN
476         IF(lwp) WRITE(numout,*) '    Initialize dust input from atmosphere '
477         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ '
478         CALL iom_open ( 'dust.orca.nc', numdust )
479         DO jm = 1, jpmois
480            CALL iom_get( numdust, jpdom_data, 'dust', zdustmo(:,:,jm), jm )
481         END DO
482         CALL iom_close( numdust )
483      ELSE
484         zdustmo(:,:,:) = 0.e0
485         dust(:,:) = 0.0
486      ENDIF
487
488      ! Nutrient input from rivers
489      ! --------------------------
490      IF( briver ) THEN
491         IF(lwp) WRITE(numout,*) '    Initialize the nutrient input by rivers from river.orca.nc file'
492         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
493         CALL iom_open ( 'river.orca.nc', numriv )
494         CALL iom_get  ( numriv, jpdom_data, 'riverdic', river   (:,:), jpan )
495         CALL iom_get  ( numriv, jpdom_data, 'riverdoc', riverdoc(:,:), jpan )
496         CALL iom_close( numriv )
497      ELSE
498         river   (:,:) = 0.e0
499         riverdoc(:,:) = 0.e0
500      endif
501
502      ! Nutrient input from dust
503      ! ------------------------
504      IF( bndepo ) THEN
505         IF(lwp) WRITE(numout,*) '    Initialize the nutrient input by dust from ndeposition.orca.nc'
506         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
507         CALL iom_open ( 'ndeposition.orca.nc', numdep )
508         CALL iom_get  ( numdep, jpdom_data, 'ndep', ndepo(:,:), jpan )
509         CALL iom_close( numdep )
510      ELSE
511         ndepo(:,:) = 0.e0
512      ENDIF
513
514      ! Coastal and island masks
515      ! ------------------------
516      IF( bsedinput ) THEN     
517         IF(lwp) WRITE(numout,*) '    Computation of an island mask to enhance coastal supply of iron'
518         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
519         IF(lwp) WRITE(numout,*) '       from bathy.orca.nc file '
520         CALL iom_open ( 'bathy.orca.nc', numbath )
521         CALL iom_get  ( numbath, jpdom_data, 'bathy', cmask(:,:,:), jpan )
522         CALL iom_close( numbath )
523         !
524         DO jk = 1, 5
525            DO jj = 2, jpjm1
526               DO ji = 2, jpim1
527                  IF( tmask(ji,jj,jk) /= 0. ) THEN
528                     zmaskt = tmask(ji+1,jj,jk) * tmask(ji-1,jj,jk) * tmask(ji,jj+1,jk)    &
529                        &                       * tmask(ji,jj-1,jk) * tmask(ji,jj,jk+1)
530                     IF( zmaskt == 0. )   cmask(ji,jj,jk ) = 0.1
531                  ENDIF
532               END DO
533            END DO
534         END DO
535         DO jk = 1, jpk
536            DO jj = 1, jpj
537               DO ji = 1, jpi
538                  expide   = MIN( 8.,( fsdept(ji,jj,jk) / 500. )**(-1.5) )
539                  denitide = -0.9543 + 0.7662 * LOG( expide ) - 0.235 * LOG( expide )**2
540                  cmask(ji,jj,jk) = cmask(ji,jj,jk) * MIN( 1., EXP( denitide ) / 0.5 )
541               END DO
542            END DO
543         END DO
544      ELSE
545         cmask(:,:,:) = 0.e0
546      ENDIF
547
548      CALL lbc_lnk( cmask , 'T', 1. )      ! Lateral boundary conditions on cmask   (sign unchanged)
549
550
551      ! total atmospheric supply of Si
552      ! ------------------------------
553      sumdepsi = 0.e0
554      DO jm = 1, jpmois
555         DO jj = 2, jpjm1
556            DO ji = 2, jpim1
557               sumdepsi = sumdepsi + zdustmo(ji,jj,jm) / (12.*rmoss) * 8.8        &
558                  &     * 0.075/28.1 * e1t(ji,jj) * e2t(ji,jj) * tmask(ji,jj,1) * tmask_i(ji,jj)
559            END DO
560         END DO
561      END DO
562      IF( lk_mpp )   CALL mpp_sum( sumdepsi )  ! sum over the global domain
563
564      ! N/P and Si releases due to coastal rivers
565      ! -----------------------------------------
566      DO jj = 1, jpj
567         DO ji = 1, jpi
568            zcoef = raass * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1)
569            cotdep(ji,jj) =  river(ji,jj)                  *1E9 / ( 12. * zcoef + rtrn ) * tmask(ji,jj,1)
570            rivinp(ji,jj) = (river(ji,jj)+riverdoc(ji,jj)) *1E9 / ( 31.6* zcoef + rtrn ) * tmask(ji,jj,1)
571            nitdep(ji,jj) = 7.6 * ndepo(ji,jj)                  / ( 14E6*raass*fse3t(ji,jj,1) + rtrn ) * tmask(ji,jj,1)
572         END DO
573      END DO
574      ! Lateral boundary conditions on ( cotdep, rivinp, nitdep )   (sign unchanged)
575      CALL lbc_lnk( cotdep , 'T', 1. )  ;  CALL lbc_lnk( rivinp , 'T', 1. )  ;  CALL lbc_lnk( nitdep , 'T', 1. )
576
577      rivpo4input=0.e0
578      rivalkinput=0.e0
579      nitdepinput=0.e0
580      DO jj = 2 , jpjm1
581         DO ji = 2, jpim1
582            zcoef = e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1) * tmask(ji,jj,1) * tmask_i(ji,jj) * raass
583            rivpo4input = rivpo4input + rivinp(ji,jj) * zcoef
584            rivalkinput = rivalkinput + cotdep(ji,jj) * zcoef
585            nitdepinput = nitdepinput + nitdep(ji,jj) * zcoef
586         END DO
587     END DO
588      IF( lk_mpp ) THEN
589         CALL mpp_sum( rivpo4input )  ! sum over the global domain
590         CALL mpp_sum( rivalkinput )  ! sum over the global domain
591         CALL mpp_sum( nitdepinput )  ! sum over the global domain
592      ENDIF
593
594
595      ! Coastal supply of iron
596      ! -------------------------
597      DO jk = 1, jpkm1
598         ironsed(:,:,jk) = sedfeinput * cmask(:,:,jk) / ( fse3t(:,:,jk) * rjjss )
599      END DO
600      CALL lbc_lnk( ironsed , 'T', 1. )      ! Lateral boundary conditions on ( ironsed )   (sign unchanged)
601
602
603   END SUBROUTINE p4z_sed_init
604
605
606
607#else
608   !!======================================================================
609   !!  Dummy module :                                   No PISCES bio-model
610   !!======================================================================
611CONTAINS
612   SUBROUTINE p4z_sed                         ! Empty routine
613   END SUBROUTINE p4z_sed
614#endif 
615
616   !!======================================================================
617END MODULE  p4zsed
Note: See TracBrowser for help on using the repository browser.