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 branches/2011/dev_r2787_LOCEAN3_TRA_TRP/NEMOGCM/NEMO/TOP_SRC/PISCES – NEMO

source: branches/2011/dev_r2787_LOCEAN3_TRA_TRP/NEMOGCM/NEMO/TOP_SRC/PISCES/p4zsed.F90 @ 2819

Last change on this file since 2819 was 2819, checked in by cetlod, 13 years ago

Improvment of branch dev_r2787_LOCEAN3_TRA_TRP

  • Property svn:keywords set to Id
File size: 22.8 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 oce_trc         !  shared variables between ocean and passive tracers
18   USE trc             !  passive tracers common variables
19   USE sms_pisces      !  PISCES Source Minus Sink variables
20   USE p4zsink         !  vertical flux of particulate matter due to sinking
21   USE p4zopt          !  optical model
22   USE p4zlim          !  Co-limitations of differents nutrients
23   USE p4zrem          !  Remineralisation of organic matter
24   USE prtctl_trc      !  print control for debugging
25   USE iom             !  I/O manager
26
27   IMPLICIT NONE
28   PRIVATE
29
30   PUBLIC   p4z_sed   
31   PUBLIC   p4z_sed_init   
32   PUBLIC   p4z_sed_alloc
33
34   !! * Shared module variables
35   LOGICAL, PUBLIC :: ln_dustfer  = .FALSE.    !: boolean for dust input from the atmosphere
36   LOGICAL, PUBLIC :: ln_river    = .FALSE.    !: boolean for river input of nutrients
37   LOGICAL, PUBLIC :: ln_ndepo    = .FALSE.    !: boolean for atmospheric deposition of N
38   LOGICAL, PUBLIC :: ln_sedinput = .FALSE.    !: boolean for Fe input from sediments
39
40   REAL(wp), PUBLIC :: sedfeinput = 1.E-9_wp   !: Coastal release of Iron
41   REAL(wp), PUBLIC :: dustsolub  = 0.014_wp   !: Solubility of the dust
42
43   !! * Module variables
44   REAL(wp) :: ryyss                  !: number of seconds per year
45   REAL(wp) :: ryyss1                 !: inverse of ryyss
46   REAL(wp) :: rmtss                  !: number of seconds per month
47   REAL(wp) :: rday1                  !: inverse of rday
48
49   INTEGER , PARAMETER :: jpmth = 12  !: number of months per year
50   INTEGER , PARAMETER :: jpyr  = 1   !: one year
51
52   INTEGER ::  numdust                !: logical unit for surface fluxes data
53   INTEGER ::  nflx1 , nflx2          !: first and second record used
54   INTEGER ::  nflx11, nflx12
55
56   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: dustmo    !: set of dust fields
57   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:) :: dust      !: dust fields
58   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:) :: rivinp, cotdep    !: river input fields
59   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:) :: nitdep    !: atmospheric N deposition
60   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ironsed   !: Coastal supply of iron
61
62   REAL(wp) :: sumdepsi, rivalkinput, rivpo4input, nitdepinput
63
64   !!* Substitution
65#  include "top_substitute.h90"
66   !!----------------------------------------------------------------------
67   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
68   !! $Header:$
69   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt)
70   !!----------------------------------------------------------------------
71
72CONTAINS
73
74
75   SUBROUTINE p4z_sed( kt, jnt )
76      !!---------------------------------------------------------------------
77      !!                     ***  ROUTINE p4z_sed  ***
78      !!
79      !! ** Purpose :   Compute loss of organic matter in the sediments. This
80      !!              is by no way a sediment model. The loss is simply
81      !!              computed to balance the inout from rivers and dust
82      !!
83      !! ** Method  : - ???
84      !!---------------------------------------------------------------------
85      USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released
86      USE wrk_nemo, ONLY: zsidep => wrk_2d_1, zwork => wrk_2d_2, zwork1 => wrk_2d_3
87      USE wrk_nemo, ONLY: znitrpot => wrk_3d_2, zirondep => wrk_3d_3
88      !
89      INTEGER, INTENT(in) ::   kt, jnt ! ocean time step
90      INTEGER  ::   ji, jj, jk, ikt
91#if ! defined key_sed
92      REAL(wp) ::   zsumsedsi, zsumsedpo4, zsumsedcal
93      REAL(wp) ::   zrivalk, zrivsil, zrivpo4
94#endif
95      REAL(wp) ::   zdenitot, znitrpottot, zlim, zfact
96      REAL(wp) ::   zwsbio3, zwsbio4, zwscal
97      CHARACTER (len=25) :: charout
98      !!---------------------------------------------------------------------
99
100      IF( ( wrk_in_use(2, 1,2,3) ) .OR. ( wrk_in_use(3, 2,3) ) ) THEN
101         CALL ctl_stop('p4z_sed: requested workspace arrays unavailable')  ;  RETURN
102      END IF
103
104      IF( jnt == 1  .AND.  ln_dustfer  )  CALL p4z_sbc( kt )
105
106      ! Iron and Si deposition at the surface
107      ! -------------------------------------
108
109      DO jj = 1, jpj
110         DO ji = 1, jpi
111            zirondep(ji,jj,1) = ( dustsolub * dust(ji,jj) / ( 55.85 * rmtss ) + 3.e-10 * ryyss1 )   &
112               &             * rfact2 / fse3t(ji,jj,1)
113            zsidep  (ji,jj)   = 8.8 * 0.075 * dust(ji,jj) * rfact2 / ( fse3t(ji,jj,1) * 28.1 * rmtss )
114         END DO
115      END DO
116
117      ! Iron solubilization of particles in the water column
118      ! ----------------------------------------------------
119
120      DO jk = 2, jpkm1
121         zirondep(:,:,jk) = dust(:,:) / ( 10. * 55.85 * rmtss ) * rfact2 * 1.e-4
122      END DO
123
124      ! Add the external input of nutrients, carbon and alkalinity
125      ! ----------------------------------------------------------
126
127      trn(:,:,1,jppo4) = trn(:,:,1,jppo4) + rivinp(:,:) * rfact2 
128      trn(:,:,1,jpno3) = trn(:,:,1,jpno3) + (rivinp(:,:) + nitdep(:,:)) * rfact2
129      trn(:,:,1,jpfer) = trn(:,:,1,jpfer) + rivinp(:,:) * 3.e-5 * rfact2
130      trn(:,:,1,jpsil) = trn(:,:,1,jpsil) + zsidep (:,:) + cotdep(:,:)   * rfact2 / 6.
131      trn(:,:,1,jpdic) = trn(:,:,1,jpdic) + rivinp(:,:) * 2.631 * rfact2
132      trn(:,:,1,jptal) = trn(:,:,1,jptal) + (cotdep(:,:) - rno3*(rivinp(:,:) +  nitdep(:,:) ) ) * rfact2
133
134
135      ! Add the external input of iron which is 3D distributed
136      ! (dust, river and sediment mobilization)
137      ! ------------------------------------------------------
138
139      DO jk = 1, jpkm1
140         trn(:,:,jk,jpfer) = trn(:,:,jk,jpfer) + zirondep(:,:,jk) + ironsed(:,:,jk) * rfact2
141      END DO
142
143
144#if ! defined key_sed
145      ! Loss of biogenic silicon, Caco3 organic carbon in the sediments.
146      ! First, the total loss is computed.
147      ! The factor for calcite comes from the alkalinity effect
148      ! -------------------------------------------------------------
149      DO jj = 1, jpj
150         DO ji = 1, jpi
151            ikt = mbkt(ji,jj) 
152# if defined key_kriest
153            zwork (ji,jj) = trn(ji,jj,ikt,jpdsi) * wscal (ji,jj,ikt)
154            zwork1(ji,jj) = trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt)
155# else
156            zwork (ji,jj) = trn(ji,jj,ikt,jpdsi) * wsbio4(ji,jj,ikt)
157            zwork1(ji,jj) = trn(ji,jj,ikt,jpgoc) * wsbio4(ji,jj,ikt) + trn(ji,jj,ikt,jppoc) * wsbio3(ji,jj,ikt) 
158# endif
159         END DO
160      END DO
161      zsumsedsi  = glob_sum( zwork (:,:) * e1e2t(:,:) ) * rday1
162      zsumsedpo4 = glob_sum( zwork1(:,:) * e1e2t(:,:) ) * rday1
163      DO jj = 1, jpj
164         DO ji = 1, jpi
165            ikt = mbkt(ji,jj) 
166            zwork (ji,jj) = trn(ji,jj,ikt,jpcal) * wscal (ji,jj,ikt)
167         END DO
168      END DO
169      zsumsedcal = glob_sum( zwork (:,:) * e1e2t(:,:) ) * 2.0 * rday1
170#endif
171
172      ! Then this loss is scaled at each bottom grid cell for
173      ! equilibrating the total budget of silica in the ocean.
174      ! Thus, the amount of silica lost in the sediments equal
175      ! the supply at the surface (dust+rivers)
176      ! ------------------------------------------------------
177
178      DO jj = 1, jpj
179         DO ji = 1, jpi
180            ikt = mbkt(ji,jj)
181            zfact = xstep / fse3t(ji,jj,ikt)
182            zwsbio3 = 1._wp - zfact * wsbio3(ji,jj,ikt)
183            zwsbio4 = 1._wp - zfact * wsbio4(ji,jj,ikt)
184            zwscal  = 1._wp - zfact * wscal (ji,jj,ikt)
185            !
186# if defined key_kriest
187            trn(ji,jj,ikt,jpdsi) = trn(ji,jj,ikt,jpdsi) * zwsbio4
188            trn(ji,jj,ikt,jpnum) = trn(ji,jj,ikt,jpnum) * zwsbio4
189            trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) * zwsbio3
190            trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) * zwsbio3
191# else
192            trn(ji,jj,ikt,jpdsi) = trn(ji,jj,ikt,jpdsi) * zwscal 
193            trn(ji,jj,ikt,jpgoc) = trn(ji,jj,ikt,jpgoc) * zwsbio4
194            trn(ji,jj,ikt,jppoc) = trn(ji,jj,ikt,jppoc) * zwsbio3
195            trn(ji,jj,ikt,jpbfe) = trn(ji,jj,ikt,jpbfe) * zwsbio4
196            trn(ji,jj,ikt,jpsfe) = trn(ji,jj,ikt,jpsfe) * zwsbio3
197# endif
198            trn(ji,jj,ikt,jpcal) = trn(ji,jj,ikt,jpcal) * zwscal
199         END DO
200      END DO
201
202#if ! defined key_sed
203      zrivsil =  1._wp - ( sumdepsi + rivalkinput * ryyss1 / 6. ) / zsumsedsi 
204      zrivalk =  1._wp - ( rivalkinput * ryyss1 ) / zsumsedcal 
205      zrivpo4 =  1._wp - ( rivpo4input * ryyss1 ) / zsumsedpo4 
206      DO jj = 1, jpj
207         DO ji = 1, jpi
208            ikt = mbkt(ji,jj)
209            zfact = xstep / fse3t(ji,jj,ikt)
210            zwsbio3 = zfact * wsbio3(ji,jj,ikt)
211            zwsbio4 = zfact * wsbio4(ji,jj,ikt)
212            zwscal  = zfact * wscal (ji,jj,ikt)
213            trn(ji,jj,ikt,jptal) =  trn(ji,jj,ikt,jptal) + trn(ji,jj,ikt,jpcal) * zwscal  * zrivalk * 2.0
214            trn(ji,jj,ikt,jpdic) =  trn(ji,jj,ikt,jpdic) + trn(ji,jj,ikt,jpcal) * zwscal  * zrivalk
215# if defined key_kriest
216            trn(ji,jj,ikt,jpsil) =  trn(ji,jj,ikt,jpsil) + trn(ji,jj,ikt,jpdsi) * zwsbio4 * zrivsil 
217            trn(ji,jj,ikt,jpdoc) =  trn(ji,jj,ikt,jpdoc) + trn(ji,jj,ikt,jppoc) * zwsbio3 * zrivpo4 
218# else
219            trn(ji,jj,ikt,jpsil) =  trn(ji,jj,ikt,jpsil) + trn(ji,jj,ikt,jpdsi) * zwscal  * zrivsil 
220            trn(ji,jj,ikt,jpdoc) =  trn(ji,jj,ikt,jpdoc)   &
221            &                     + ( trn(ji,jj,ikt,jppoc) * zwsbio3 + trn(ji,jj,ikt,jpgoc) * zwsbio4 ) * zrivpo4
222# endif
223         END DO
224      END DO
225# endif
226
227      ! Nitrogen fixation (simple parameterization). The total gain
228      ! from nitrogen fixation is scaled to balance the loss by
229      ! denitrification
230      ! -------------------------------------------------------------
231
232      zdenitot = glob_sum( denitr(:,:,:)  * cvol(:,:,:) * xnegtr(:,:,:) ) * rdenit
233
234      ! Potential nitrogen fixation dependant on temperature and iron
235      ! -------------------------------------------------------------
236
237!CDIR NOVERRCHK
238      DO jk = 1, jpk
239!CDIR NOVERRCHK
240         DO jj = 1, jpj
241!CDIR NOVERRCHK
242            DO ji = 1, jpi
243               zlim = ( 1.- xnanono3(ji,jj,jk) - xnanonh4(ji,jj,jk) )
244               IF( zlim <= 0.2 )   zlim = 0.01
245               znitrpot(ji,jj,jk) = MAX( 0.e0, ( 0.6 * tgfunc(ji,jj,jk) - 2.15 ) * rday1 )   &
246                  &                  * zlim * rfact2 * trn(ji,jj,jk,jpfer)   &
247# if defined key_degrad
248                  &                  * facvol(ji,jj,jk)
249# endif
250                  &                  / ( conc3 + trn(ji,jj,jk,jpfer) ) * ( 1.- EXP( -etot(ji,jj,jk) / 50.) )
251            END DO
252         END DO
253      END DO
254
255      znitrpottot = glob_sum( znitrpot(:,:,:) * cvol(:,:,:) )
256
257      ! Nitrogen change due to nitrogen fixation
258      ! ----------------------------------------
259
260      DO jk = 1, jpk
261         DO jj = 1, jpj
262            DO ji = 1, jpi
263               zfact = znitrpot(ji,jj,jk) * 1.e-7
264               trn(ji,jj,jk,jpnh4) = trn(ji,jj,jk,jpnh4) + zfact
265               trn(ji,jj,jk,jpoxy) = trn(ji,jj,jk,jpoxy) + zfact   * o2nit
266               trn(ji,jj,jk,jppo4) = trn(ji,jj,jk,jppo4) + 30./ 46.* zfact
267            END DO
268         END DO
269      END DO
270
271      IF( ln_diatrc ) THEN
272         !
273         zfact = 1.e+3 * rfact2r
274         zwork (:,:)  =  ( zirondep(:,:,1) + ironsed(:,:,1) * rfact2 ) * zfact * fse3t(:,:,1) * tmask(:,:,1) 
275         zwork1(:,:)  =    znitrpot(:,:,1) * 1.e-7                     * zfact * fse3t(:,:,1) * tmask(:,:,1)
276         IF( jnt == nrdttrc ) THEN
277            CALL iom_put( "Irondep", zwork  )  ! surface downward net flux of iron
278            CALL iom_put( "Nfix"   , zwork1 )  ! nitrogen fixation at surface
279         ENDIF
280# if ! defined key_iomput
281         trc2d(:,:,jp_pcs0_2d + 11) = zirondep(:,:,1)         * zfact * fse3t(:,:,1) * tmask(:,:,1)
282         trc2d(:,:,jp_pcs0_2d + 12) = znitrpot(:,:,1) * 1.e-7 * zfact * fse3t(:,:,1) * tmask(:,:,1)
283# endif
284         !
285      ENDIF
286      !
287       IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
288         WRITE(charout, FMT="('sed ')")
289         CALL prt_ctl_trc_info(charout)
290         CALL prt_ctl_trc(tab4d=trn, mask=tmask, clinfo=ctrcnm)
291       ENDIF
292
293      IF( ( wrk_not_released(2, 1,2,3) ) .OR. ( wrk_not_released(3, 2,3) ) )   &
294        &         CALL ctl_stop('p4z_sed: failed to release workspace arrays')
295
296   END SUBROUTINE p4z_sed
297
298   SUBROUTINE p4z_sbc( kt )
299
300      !!----------------------------------------------------------------------
301      !!                  ***  ROUTINE p4z_sbc  ***
302      !!
303      !! ** Purpose :   Read and interpolate the external sources of
304      !!                nutrients
305      !!
306      !! ** Method  :   Read the files and interpolate the appropriate variables
307      !!
308      !! ** input   :   external netcdf files
309      !!
310      !!----------------------------------------------------------------------
311      !! * arguments
312      INTEGER, INTENT( in  ) ::   kt   ! ocean time step
313
314      !! * Local declarations
315      INTEGER :: imois, i15, iman 
316      REAL(wp) :: zxy
317
318      !!---------------------------------------------------------------------
319
320      ! Initialization
321      ! --------------
322
323      i15 = nday / 16
324      iman  = INT( raamo )
325      imois = nmonth + i15 - 1
326      IF( imois == 0 ) imois = iman
327
328      ! Calendar computation
329      IF( kt == nit000 .OR. imois /= nflx1 ) THEN
330
331         IF( kt == nit000 )  nflx1  = 0
332
333         ! nflx1 number of the first file record used in the simulation
334         ! nflx2 number of the last  file record
335
336         nflx1 = imois
337         nflx2 = nflx1 + 1
338         nflx1 = MOD( nflx1, iman )
339         nflx2 = MOD( nflx2, iman )
340         IF( nflx1 == 0 )   nflx1 = iman
341         IF( nflx2 == 0 )   nflx2 = iman
342         IF(lwp) WRITE(numout,*) 
343         IF(lwp) WRITE(numout,*) ' p4z_sbc : first record file used nflx1 ',nflx1
344         IF(lwp) WRITE(numout,*) ' p4z_sbc : last  record file used nflx2 ',nflx2
345
346      ENDIF
347
348      ! 3. at every time step interpolation of fluxes
349      ! ---------------------------------------------
350
351      zxy = FLOAT( nday + 15 - 30 * i15 ) / 30
352      dust(:,:) = ( (1.-zxy) * dustmo(:,:,nflx1) + zxy * dustmo(:,:,nflx2) )
353
354   END SUBROUTINE p4z_sbc
355
356   SUBROUTINE p4z_sed_init
357
358      !!----------------------------------------------------------------------
359      !!                  ***  ROUTINE p4z_sed_init  ***
360      !!
361      !! ** Purpose :   Initialization of the external sources of nutrients
362      !!
363      !! ** Method  :   Read the files and compute the budget
364      !!      called at the first timestep (nit000)
365      !!
366      !! ** input   :   external netcdf files
367      !!
368      !!----------------------------------------------------------------------
369      USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released
370      USE wrk_nemo, ONLY: zriverdoc => wrk_2d_1, zriver => wrk_2d_2, zndepo => wrk_2d_3
371      USE wrk_nemo, ONLY: zcmask => wrk_3d_2
372      !
373      INTEGER :: ji, jj, jk, jm
374      INTEGER :: numriv, numbath, numdep
375      REAL(wp) ::   zcoef
376      REAL(wp) ::   expide, denitide,zmaskt
377      !
378      NAMELIST/nampissed/ ln_dustfer, ln_river, ln_ndepo, ln_sedinput, sedfeinput, dustsolub
379      !!----------------------------------------------------------------------
380
381      IF( ( wrk_in_use(2, 1,2,3) ) .OR. ( wrk_in_use(3, 2) ) ) THEN
382         CALL ctl_stop('p4z_sed_init: requested workspace arrays unavailable')  ;  RETURN
383      END IF
384
385      REWIND( numnatp )                     ! read numnatp
386      READ  ( numnatp, nampissed )
387
388      IF(lwp) THEN
389         WRITE(numout,*) ' '
390         WRITE(numout,*) ' Namelist : nampissed '
391         WRITE(numout,*) ' ~~~~~~~~~~~~~~~~~ '
392         WRITE(numout,*) '    Dust input from the atmosphere           ln_dustfer  = ', ln_dustfer
393         WRITE(numout,*) '    River input of nutrients                 ln_river    = ', ln_river
394         WRITE(numout,*) '    Atmospheric deposition of N              ln_ndepo    = ', ln_ndepo
395         WRITE(numout,*) '    Fe input from sediments                  ln_sedinput = ', ln_sedinput
396         WRITE(numout,*) '    Coastal release of Iron                  sedfeinput  =', sedfeinput
397         WRITE(numout,*) '    Solubility of the dust                   dustsolub   =', dustsolub
398      ENDIF
399
400      ! Dust input from the atmosphere
401      ! ------------------------------
402      IF( ln_dustfer ) THEN
403         IF(lwp) WRITE(numout,*) '    Initialize dust input from atmosphere '
404         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ '
405         CALL iom_open ( 'dust.orca.nc', numdust )
406         DO jm = 1, jpmth
407            CALL iom_get( numdust, jpdom_data, 'dust', dustmo(:,:,jm), jm )
408         END DO
409         CALL iom_close( numdust )
410      ELSE
411         dustmo(:,:,:) = 0.e0
412         dust(:,:) = 0.0
413      ENDIF
414
415      ! Nutrient input from rivers
416      ! --------------------------
417      IF( ln_river ) THEN
418         IF(lwp) WRITE(numout,*) '    Initialize the nutrient input by rivers from river.orca.nc file'
419         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
420         CALL iom_open ( 'river.orca.nc', numriv )
421         CALL iom_get  ( numriv, jpdom_data, 'riverdic', zriver   (:,:), jpyr )
422         CALL iom_get  ( numriv, jpdom_data, 'riverdoc', zriverdoc(:,:), jpyr )
423         CALL iom_close( numriv )
424      ELSE
425         zriver   (:,:) = 0.e0
426         zriverdoc(:,:) = 0.e0
427      endif
428
429      ! Nutrient input from dust
430      ! ------------------------
431      IF( ln_ndepo ) THEN
432         IF(lwp) WRITE(numout,*) '    Initialize the nutrient input by dust from ndeposition.orca.nc'
433         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
434         CALL iom_open ( 'ndeposition.orca.nc', numdep )
435         CALL iom_get  ( numdep, jpdom_data, 'ndep', zndepo(:,:), jpyr )
436         CALL iom_close( numdep )
437      ELSE
438         zndepo(:,:) = 0.e0
439      ENDIF
440
441      ! Coastal and island masks
442      ! ------------------------
443      IF( ln_sedinput ) THEN     
444         IF(lwp) WRITE(numout,*) '    Computation of an island mask to enhance coastal supply of iron'
445         IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~'
446         IF(lwp) WRITE(numout,*) '       from bathy.orca.nc file '
447         CALL iom_open ( 'bathy.orca.nc', numbath )
448         CALL iom_get  ( numbath, jpdom_data, 'bathy', zcmask(:,:,:), jpyr )
449         CALL iom_close( numbath )
450         !
451         DO jk = 1, 5
452            DO jj = 2, jpjm1
453               DO ji = fs_2, fs_jpim1
454                  IF( tmask(ji,jj,jk) /= 0. ) THEN
455                     zmaskt = tmask(ji+1,jj,jk) * tmask(ji-1,jj,jk) * tmask(ji,jj+1,jk)    &
456                        &                       * tmask(ji,jj-1,jk) * tmask(ji,jj,jk+1)
457                     IF( zmaskt == 0. )   zcmask(ji,jj,jk ) = MAX( 0.1, zcmask(ji,jj,jk) ) 
458                  ENDIF
459               END DO
460            END DO
461         END DO
462         DO jk = 1, jpk
463            DO jj = 1, jpj
464               DO ji = 1, jpi
465                  expide   = MIN( 8.,( fsdept(ji,jj,jk) / 500. )**(-1.5) )
466                  denitide = -0.9543 + 0.7662 * LOG( expide ) - 0.235 * LOG( expide )**2
467                  zcmask(ji,jj,jk) = zcmask(ji,jj,jk) * MIN( 1., EXP( denitide ) / 0.5 )
468               END DO
469            END DO
470         END DO
471      ELSE
472         zcmask(:,:,:) = 0.e0
473      ENDIF
474
475      CALL lbc_lnk( zcmask , 'T', 1. )      ! Lateral boundary conditions on zcmask   (sign unchanged)
476
477
478      !                                    ! Number of seconds per year and per month
479      ryyss  = nyear_len(1) * rday
480      rmtss  = ryyss / raamo
481      rday1  = 1. / rday
482      ryyss1 = 1. / ryyss
483      !                                    ! ocean surface cell
484
485      ! total atmospheric supply of Si
486      ! ------------------------------
487      sumdepsi = 0.e0
488      DO jm = 1, jpmth
489         zcoef = 1. / ( 12. * rmtss ) * 8.8 * 0.075 / 28.1       
490         sumdepsi = sumdepsi + glob_sum( dustmo(:,:,jm) * e1e2t(:,:) ) * zcoef
491      ENDDO
492
493      ! N/P and Si releases due to coastal rivers
494      ! -----------------------------------------
495      DO jj = 1, jpj
496         DO ji = 1, jpi
497            zcoef = ryyss * e1e2t(ji,jj)  * fse3t(ji,jj,1) * tmask(ji,jj,1) 
498            cotdep(ji,jj) =  zriver(ji,jj)                  *1E9 / ( 12. * zcoef + rtrn )
499            rivinp(ji,jj) = (zriver(ji,jj)+zriverdoc(ji,jj)) *1E9 / ( 31.6* zcoef + rtrn )
500            nitdep(ji,jj) = 7.6 * zndepo(ji,jj)                  / ( 14E6*ryyss*fse3t(ji,jj,1) + rtrn )
501         END DO
502      END DO
503      ! Lateral boundary conditions on ( cotdep, rivinp, nitdep )   (sign unchanged)
504      CALL lbc_lnk( cotdep , 'T', 1. )  ;  CALL lbc_lnk( rivinp , 'T', 1. )  ;  CALL lbc_lnk( nitdep , 'T', 1. )
505
506      rivpo4input = glob_sum( rivinp(:,:) * cvol(:,:,1) ) * ryyss
507      rivalkinput = glob_sum( cotdep(:,:) * cvol(:,:,1) ) * ryyss
508      nitdepinput = glob_sum( nitdep(:,:) * cvol(:,:,1) ) * ryyss
509
510
511      ! Coastal supply of iron
512      ! -------------------------
513      DO jk = 1, jpkm1
514         ironsed(:,:,jk) = sedfeinput * zcmask(:,:,jk) / ( fse3t(:,:,jk) * rday )
515      END DO
516      CALL lbc_lnk( ironsed , 'T', 1. )      ! Lateral boundary conditions on ( ironsed )   (sign unchanged)
517
518      IF( ( wrk_not_released(2, 1,2,3) ) .OR. ( wrk_not_released(3, 2) ) )   &
519        &         CALL ctl_stop('p4z_sed_init: failed to release workspace arrays')
520
521   END SUBROUTINE p4z_sed_init
522
523   INTEGER FUNCTION p4z_sed_alloc()
524      !!----------------------------------------------------------------------
525      !!                     ***  ROUTINE p4z_sed_alloc  ***
526      !!----------------------------------------------------------------------
527
528      ALLOCATE( dustmo(jpi,jpj,jpmth), dust(jpi,jpj)       ,     &
529        &       rivinp(jpi,jpj)      , cotdep(jpi,jpj)     ,     &
530        &       nitdep(jpi,jpj)      , ironsed(jpi,jpj,jpk), STAT=p4z_sed_alloc ) 
531
532      IF( p4z_sed_alloc /= 0 ) CALL ctl_warn('p4z_sed_alloc : failed to allocate arrays.')
533
534   END FUNCTION p4z_sed_alloc
535#else
536   !!======================================================================
537   !!  Dummy module :                                   No PISCES bio-model
538   !!======================================================================
539CONTAINS
540   SUBROUTINE p4z_sed                         ! Empty routine
541   END SUBROUTINE p4z_sed
542#endif 
543
544   !!======================================================================
545END MODULE  p4zsed
Note: See TracBrowser for help on using the repository browser.