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.
p4zsink_kriest.F90 in branches/dev_001_GM/NEMO/TOP_SRC/PISCES_SMS – NEMO

source: branches/dev_001_GM/NEMO/TOP_SRC/PISCES_SMS/p4zsink_kriest.F90 @ 775

Last change on this file since 775 was 775, checked in by gm, 16 years ago

dev_001_GM - PISCES in F90 : encapsulation of all p4z...F files in module F90 + doctor norme for local variables - compilation OK
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 11.1 KB
Line 
1MODULE p4zsink_kriest
2   !!======================================================================
3   !!                   ***  MODULE p4zsink_kriest  ***
4   !! TOP :   PISCES Compute vertical flux of particulate matter due to gravitational sinking
5   !!         Kriest parameterization
6   !!======================================================================
7   !! History :   1.0  !  2004     (O. Aumont) Original code
8   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
9   !!----------------------------------------------------------------------
10#if defined key_pisces  &&  defined key_kriest
11   !!----------------------------------------------------------------------
12   !!   'key_pisces'    and                                PISCES bio-model
13   !!   'key_kriest'                                          kriest option
14   !!----------------------------------------------------------------------
15   !!   p4z_sink_kriest :  Compute vertical flux of particulate matter due
16   !!                      to gravitational sinking (Kriest parameterization)
17   !!----------------------------------------------------------------------
18   USE oce_trc         !
19   USE trp_trc
20   USE sms
21   USE p4zsink2
22
23   IMPLICIT NONE
24   PRIVATE
25
26   PUBLIC   p4z_sink_kriest    ! called in p4zbio.F90
27
28   !!* Substitution
29#  include "domzgr_substitute.h90"
30   !!----------------------------------------------------------------------
31   !! NEMO/TOP 2.0 , LOCEAN-IPSL (2007)
32   !! $Header:$
33   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt)
34   !!----------------------------------------------------------------------
35
36CONTAINS
37
38   SUBROUTINE p4z_sink_kriest
39      !!---------------------------------------------------------------------
40      !!                ***  ROUTINE p4z_sink_kriest  ***
41      !!
42      !! ** Purpose :   Compute vertical flux of particulate matter due to
43      !!              gravitational sinking - Kriest parameterization
44      !!
45      !! ** Method  : - ???
46      !!---------------------------------------------------------------------
47      INTEGER  ::   ji, jj, jk
48      INTEGER  ::   iksed
49      REAL(wp) ::   zagg1, zagg2, zagg3, zagg4, zagg5, zaggsi, zaggsh
50      REAL(wp) ::   znum , zeps, zfm, zgm, zsm
51      REAL(wp) ::   zdiv , zdiv1, zdiv2, zdiv3, zdiv4, zdiv5
52      REAL(wp) ::   zval1, zval2, zval3, zval4
53      REAL(wp) ::   zstep
54#if defined key_trc_dia3d
55      REAL(wp) ::   zrfact2
56#endif
57      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   znum3d
58      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   sinking, sinking2
59      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   sinkfer
60      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   sinkcal, sinksil
61      !!---------------------------------------------------------------------
62
63       zstep=rfact2/rjjss      ! Time step duration for biology
64
65
66!     Initialisation of variables used to compute Sinking Speed
67!     ---------------------------------------------------------
68
69       znum3d(:,:,:) = 0.e0
70       iksed = 10
71       zval1 = 1. + xkr_zeta
72       zval2 = 1. + xkr_zeta + xkr_eta
73       zval3 = 1. + xkr_eta
74
75!     Computation of the vertical sinking speed : Kriest et Evans, 2000
76!     -----------------------------------------------------------------
77   
78      DO jk = 1, jpkm1
79         DO jj = 1, jpj
80            DO ji = 1, jpi
81               IF( tmask(ji,jj,jk) /= 0.e0 ) THEN
82                  znum = trn(ji,jj,jk,jppoc) / ( trn(ji,jj,jk,jpnum) + rtrn ) / xkr_massp
83! -------------- To avoid sinking speed over 50 m/day -------
84                  znum  = MIN( xnumm(jk), znum )
85                  znum  = MAX( 1.1      , znum )
86                  znum3d(ji,jj,jk) = znum
87!------------------------------------------------------------
88                  zeps  = ( zval1 * znum - 1. )/ ( znum - 1. )
89                  zfm   = xkr_frac**( 1. - zeps )
90                  zgm   = xkr_frac**( zval1 - zeps )
91                  zdiv  = MAX( 1.e-4, ABS( zeps - zval2 ) ) * SIGN( 1., ( zeps - zval2 ) )
92                  zdiv1 = zeps - zval3
93!!gmoptimisation possible below
94                  wsbio3(ji,jj,jk) = xkr_wsbio_min * ( zeps - zval1 ) / zdiv    &
95     &                             - xkr_wsbio_max *   zgm * xkr_eta  / zdiv
96                  wsbio4(ji,jj,jk) = xkr_wsbio_min *   ( zeps-1. )    / zdiv1   &
97     &                             - xkr_wsbio_max *   zfm * xkr_eta  / zdiv1
98                  IF( znum == 1.1)   wsbio3(ji,jj,jk) = wsbio4(ji,jj,jk)
99               ENDIF
100            END DO
101         END DO
102      END DO
103
104      wscal(:,:,:) = MAX( wsbio3(:,:,:), 50. )
105
106
107!   INITIALIZE TO ZERO ALL THE SINKING ARRAYS
108!   -----------------------------------------
109
110      sinking (:,:,:) = 0.e0
111      sinking2(:,:,:) = 0.e0
112      sinkcal (:,:,:) = 0.e0
113      sinkfer (:,:,:) = 0.e0
114      sinksil (:,:,:) = 0.e0
115
116!   Compute the sedimentation term using p4zsink2 for all
117!   the sinking particles
118!   -----------------------------------------------------
119
120      CALL p4z_sink2( wsbio3, sinking , jppoc )
121      CALL p4z_sink2( wsbio4, sinking2, jpnum )
122      CALL p4z_sink2( wsbio3, sinkfer , jpsfe )
123      CALL p4z_sink2( wscal , sinksil , jpdsi )
124      CALL p4z_sink2( wscal , sinkcal , jpcal )
125
126!  Exchange between organic matter compartments due to
127!  coagulation/disaggregation
128!  ---------------------------------------------------
129
130      zval1 = 1. + xkr_zeta
131      zval2 = 1. + xkr_eta
132      zval3 = 3. + xkr_eta
133      zval4 = 4. + xkr_eta
134
135      DO jk = 1,jpkm1
136         DO jj = 1,jpj
137            DO ji = 1,jpi
138               IF( tmask(ji,jj,jk) /= 0.e0 ) THEN
139
140                  znum = trn(ji,jj,jk,jppoc)/(trn(ji,jj,jk,jpnum)+rtrn) / xkr_massp
141! -------------- To avoid sinking speed over 50 m/day -------
142                  znum  = min(xnumm(jk),znum)
143                  znum  = MAX( 1.1,znum)
144!------------------------------------------------------------
145                  zeps  = ( zval1 * znum - 1.) / ( znum - 1.)
146                  zdiv  = MAX( 1.e-4, ABS( zeps - zval3) ) * SIGN( 1., zeps - zval3 )
147                  zdiv1 = MAX( 1.e-4, ABS( zeps - 4.   ) ) * SIGN( 1., zeps - 4.    )
148                  zdiv2 = zeps - 2.
149                  zdiv3 = zeps - 3.
150                  zdiv4 = zeps - zval2
151                  zdiv5 = 2.* zeps - zval4
152                  zfm   = xkr_frac**( 1.- zeps )
153                  zsm   = xkr_frac**xkr_eta
154
155!    Part I : Coagulation dependant on turbulence
156!    ----------------------------------------------
157
158                  zagg1 = ( 0.163 * trn(ji,jj,jk,jpnum)**2               &
159                     &            * 2.*( (zfm-1.)*(zfm*xkr_mass_max**3-xkr_mass_min**3)    &
160                     &            * (zeps-1)/zdiv1 + 3.*(zfm*xkr_mass_max-xkr_mass_min)    &
161                     &            * (zfm*xkr_mass_max**2-xkr_mass_min**2)                  &
162                     &            * (zeps-1.)**2/(zdiv2*zdiv3))            &
163# if defined key_off_degrad
164                     &                 *facvol(ji,jj,jk)       &
165# endif
166                     &    )
167
168                  zagg2 = (  2*0.163*trn(ji,jj,jk,jpnum)**2*zfm*                       &
169                     &                   ((xkr_mass_max**3+3.*(xkr_mass_max**2          &
170                     &                    *xkr_mass_min*(zeps-1.)/zdiv2                 &
171                     &                    +xkr_mass_max*xkr_mass_min**2*(zeps-1.)/zdiv3)    &
172                     &                    +xkr_mass_min**3*(zeps-1)/zdiv1)                  &
173                     &                    -zfm*xkr_mass_max**3*(1.+3.*((zeps-1.)/           &
174                     &                    (zeps-2.)+(zeps-1.)/zdiv3)+(zeps-1.)/zdiv1))      &
175#    if defined key_off_degrad
176                     &                 *facvol(ji,jj,jk)             &
177#    endif
178                     &    )
179
180                  zagg3 = (  0.163*trn(ji,jj,jk,jpnum)**2*zfm**2*8. * xkr_mass_max**3   &
181#    if defined key_off_degrad
182                     &                 *facvol(ji,jj,jk)             &
183#    endif
184                     &    )
185
186                  zaggsh = ( zagg1 + zagg2 + zagg3 ) * rfact2 * zdiss(ji,jj,jk) / 1000.
187
188!    Aggregation of small into large particles
189!    Part II : Differential settling
190!    ----------------------------------------------
191
192                  zagg4 = (  2.*3.141*0.125*trn(ji,jj,jk,jpnum)**2*                       &
193                     &                 xkr_wsbio_min*(zeps-1.)**2                         &
194                     &                 *(xkr_mass_min**2*((1.-zsm*zfm)/(zdiv3*zdiv4)      &
195                     &                 -(1.-zfm)/(zdiv*(zeps-1.)))-                       &
196                     &                 ((zfm*zfm*xkr_mass_max**2*zsm-xkr_mass_min**2)     &
197                     &                 *xkr_eta)/(zdiv*zdiv3*zdiv5) )                     &
198# if defined key_off_degrad
199                     &                 *facvol(ji,jj,jk)        &
200# endif
201                     &    )
202
203                  zagg5 = (  2.*3.141*0.125*trn(ji,jj,jk,jpnum)**2                         &
204                     &                 *(zeps-1.)*zfm*xkr_wsbio_min                        &
205                     &                 *(zsm*(xkr_mass_min**2-zfm*xkr_mass_max**2)         &
206                     &                 /zdiv3-(xkr_mass_min**2-zfm*zsm*xkr_mass_max**2)    &
207                     &                 /zdiv)                   &
208# if defined key_off_degrad
209                     &                 *facvol(ji,jj,jk)        &
210# endif
211                     &    )
212
213                  zaggsi = ( zagg4 + zagg5 ) * zstep / 10.
214
215                  xagg(ji,jj,jk) = 0.5 * xkr_stick * ( zaggsh + zaggsi )
216
217!     Aggregation of DOC to small particles
218!     --------------------------------------
219
220                  xaggdoc(ji,jj,jk) = (     0.4 * trn(ji,jj,jk,jpdoc)               &
221                     &                 + 1018.  * trn(ji,jj,jk,jppoc)  ) * zstep    &
222# if defined key_off_degrad
223                     &              * facvol(ji,jj,jk)                              &
224# endif
225                     &              * zdiss(ji,jj,jk) * trn(ji,jj,jk,jpdoc)
226
227               ENDIF
228            END DO
229         END DO
230      END DO
231
232#    if defined key_trc_dia3d
233      zrfact2 = 1.e3 * rfact2r
234      trc2d(:,:, 5)   = sinking (:,:,iksed+1) * zrfact2
235      trc2d(:,:, 6)   = sinking2(:,:,iksed+1) * zrfact2
236      trc2d(:,:, 7)   = sinkfer (:,:,iksed+1) * zrfact2
237      trc2d(:,:, 9)   = sinksil (:,:,iksed+1) * zrfact2
238      trc2d(:,:,10)   = sinkcal (:,:,iksed+1) * zrfact2
239      trc3d(:,:,:,12) = sinking (:,:,:)       * zrfact2
240      trc3d(:,:,:,13) = sinking2(:,:,:)       * zrfact2
241      trc3d(:,:,:,14) = sinksil (:,:,:)       * zrfact2
242      trc3d(:,:,:,15) = sinkcal (:,:,:)       * zrfact2
243      trc3d(:,:,:,16) = znum3d  (:,:,:)
244      trc3d(:,:,:,17) = wsbio3  (:,:,:)
245      trc3d(:,:,:,18) = wsbio4  (:,:,:)
246#    endif
247      !
248   END SUBROUTINE p4z_sink_kriest
249
250#else
251   !!======================================================================
252   !!  Dummy module :                                   No PISCES bio-model
253   !!======================================================================
254CONTAINS
255   SUBROUTINE p4z_sink_kriest                    ! Empty routine
256   END SUBROUTINE p4z_sink_kriest
257#endif 
258
259   !!======================================================================
260END MODULE  p4zsink_kriest
Note: See TracBrowser for help on using the repository browser.