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.
p4zagg.F90 in NEMO/trunk/src/TOP/PISCES/P4Z – NEMO

source: NEMO/trunk/src/TOP/PISCES/P4Z/p4zagg.F90 @ 15459

Last change on this file since 15459 was 15459, checked in by cetlod, 3 years ago

Various bug fixes and more comments in PISCES routines ; sette test OK in debug mode, nn_hls=1/2, with tiling ; run.stat unchanged ; of course tracer.stat different
  • Property svn:keywords set to Id
File size: 10.2 KB
RevLine 
[7162]1MODULE p4zagg
2   !!======================================================================
3   !!                         ***  MODULE p4zagg  ***
[15459]4   !! TOP :  PISCES  aggregation of particles (DOC, POC, GOC)
5   !!        This module is the same for both PISCES and PISCES-QUOTA
[7162]6   !!======================================================================
7   !! History :   1.0  !  2004     (O. Aumont) Original code
8   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
9   !!             3.4  !  2011-06  (O. Aumont, C. Ethe) Change aggregation formula
10   !!             3.5  !  2012-07  (O. Aumont) Introduce potential time-splitting
[7176]11   !!             3.6  !  2015-05  (O. Aumont) PISCES quota
[7162]12   !!----------------------------------------------------------------------
[9124]13
[7162]14   !!----------------------------------------------------------------------
15   !!   p4z_agg       :  Compute aggregation of particles
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
[13286]20   USE prtctl          !  print control for debugging
[7162]21
22   IMPLICIT NONE
23   PRIVATE
24
25   PUBLIC   p4z_agg         ! called in p4zbio.F90
26
[12377]27   !! * Substitutions
28#  include "do_loop_substitute.h90"
[7162]29   !!----------------------------------------------------------------------
[10067]30   !! NEMO/TOP 4.0 , NEMO Consortium (2018)
[10069]31   !! $Id$
[10068]32   !! Software governed by the CeCILL license (see ./LICENSE)
[7162]33   !!----------------------------------------------------------------------
34CONTAINS
35
[12377]36   SUBROUTINE p4z_agg ( kt, knt, Kbb, Krhs )
[7162]37      !!---------------------------------------------------------------------
38      !!                     ***  ROUTINE p4z_agg  ***
39      !!
[15459]40      !! ** Purpose :   Compute aggregation of particle. Aggregation by
41      !!                brownian motion, differential settling and shear
42      !!                are considered.
[7162]43      !!
[15459]44      !! ** Method  : - Aggregation rates are computed assuming a fixed and
45      !!                constant size spectrum in the different particulate
46      !!                pools. The coagulation rates have been computed
47      !!                externally using dedicated programs (O. Aumont). They
48      !!                are hard-coded because they can't be changed
49      !!                independently of each other.
[7162]50      !!---------------------------------------------------------------------
[9124]51      INTEGER, INTENT(in) ::   kt, knt   !
[12377]52      INTEGER, INTENT(in) ::   Kbb, Krhs ! time level indices
[9124]53      !
[7162]54      INTEGER  ::   ji, jj, jk
[7176]55      REAL(wp) ::   zagg, zagg1, zagg2, zagg3, zagg4
56      REAL(wp) ::   zaggpoc1, zaggpoc2, zaggpoc3, zaggpoc4
57      REAL(wp) ::   zaggpoc , zaggfe, zaggdoc, zaggdoc2, zaggdoc3
58      REAL(wp) ::   zaggpon , zaggdon, zaggdon2, zaggdon3
59      REAL(wp) ::   zaggpop, zaggdop, zaggdop2, zaggdop3
60      REAL(wp) ::   zaggtmp, zfact, zmax
[7162]61      CHARACTER (len=25) :: charout
62      !!---------------------------------------------------------------------
63      !
[9124]64      IF( ln_timing )   CALL timing_start('p4z_agg')
[7162]65      !
[15459]66      !  Exchange between organic matter compartments due to
67      !  coagulation/disaggregation
[7162]68      !  ---------------------------------------------------
[15459]69
70      ! PISCES part
[7176]71      IF( ln_p4z ) THEN
72         !
[15090]73         DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 )
[12377]74            !
75            zfact = xstep * xdiss(ji,jj,jk)
[15459]76            ! Part I : Coagulation dependent on turbulence
77            !  The stickiness has been assumed to be 0.1
[12377]78            !  Part I : Coagulation dependent on turbulence
[15459]79            zagg1 = 12.5  * zfact * tr(ji,jj,jk,jppoc,Kbb) * tr(ji,jj,jk,jppoc,Kbb)
80            zagg2 = 169.7 * zfact * tr(ji,jj,jk,jppoc,Kbb) * tr(ji,jj,jk,jpgoc,Kbb)
[7162]81
[12377]82            ! Part II : Differential settling
[15459]83            ! Aggregation of small into large particles
84            ! The stickiness has been assumed to be 0.1
85            zagg3 =  8.63  * xstep * tr(ji,jj,jk,jppoc,Kbb) * tr(ji,jj,jk,jppoc,Kbb)
86            zagg4 =  132.8 * xstep * tr(ji,jj,jk,jppoc,Kbb) * tr(ji,jj,jk,jpgoc,Kbb)
[7162]87
[12377]88            zagg   = zagg1 + zagg2 + zagg3 + zagg4
89            zaggfe = zagg * tr(ji,jj,jk,jpsfe,Kbb) / ( tr(ji,jj,jk,jppoc,Kbb) + rtrn )
[7162]90
[12377]91            ! Aggregation of DOC to POC :
92            ! 1st term is shear aggregation of DOC-DOC
93            ! 2nd term is shear aggregation of DOC-POC
94            ! 3rd term is differential settling of DOC-POC
[15459]95            ! 1/3 of DOC is supposed to experience aggregation (HMW)
96            zaggdoc  = ( ( 12.0 * 0.3 * tr(ji,jj,jk,jpdoc,Kbb) + 9.05 * tr(ji,jj,jk,jppoc,Kbb) ) * zfact       &
97            &            + 2.49 * xstep * tr(ji,jj,jk,jppoc,Kbb) ) * 0.3 * tr(ji,jj,jk,jpdoc,Kbb)
[12377]98            ! transfer of DOC to GOC :
99            ! 1st term is shear aggregation
[15459]100            ! 1/3 of DOC is supposed to experience aggregation (HMW)
101            zaggdoc2 = ( 1.94 * zfact + 1.37 * xstep ) * tr(ji,jj,jk,jpgoc,Kbb) * 0.3 * tr(ji,jj,jk,jpdoc,Kbb)
[12377]102            ! tranfer of DOC to POC due to brownian motion
[15459]103            ! The temperature dependency has been omitted.
104            zaggdoc3 =  ( 127.8 * 0.3 * tr(ji,jj,jk,jpdoc,Kbb) + 725.7 * tr(ji,jj,jk,jppoc,Kbb) ) * xstep * 0.3 * tr(ji,jj,jk,jpdoc,Kbb)
[7162]105
[12377]106            !  Update the trends
107            tr(ji,jj,jk,jppoc,Krhs) = tr(ji,jj,jk,jppoc,Krhs) - zagg + zaggdoc + zaggdoc3
108            tr(ji,jj,jk,jpgoc,Krhs) = tr(ji,jj,jk,jpgoc,Krhs) + zagg + zaggdoc2
109            tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) - zaggfe
110            tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + zaggfe
111            tr(ji,jj,jk,jpdoc,Krhs) = tr(ji,jj,jk,jpdoc,Krhs) - zaggdoc - zaggdoc2 - zaggdoc3
112            !
113            conspoc(ji,jj,jk) = conspoc(ji,jj,jk) - zagg + zaggdoc + zaggdoc3
114            prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zagg + zaggdoc2
115            !
116         END_3D
[7176]117      ELSE    ! ln_p5z
[15459]118        ! PISCES-QUOTA part
[7176]119        !
[15090]120         DO_3D( nn_hls, nn_hls, nn_hls, nn_hls, 1, jpkm1 )
[12377]121            !
122            zfact = xstep * xdiss(ji,jj,jk)
123            !  Part I : Coagulation dependent on turbulence
[15459]124            ! The stickiness has been assumed to be 0.1
[12377]125            zaggtmp = 25.9  * zfact * tr(ji,jj,jk,jppoc,Kbb)
126            zaggpoc1 = zaggtmp * tr(ji,jj,jk,jppoc,Kbb)
127            zaggtmp = 4452. * zfact * tr(ji,jj,jk,jpgoc,Kbb)
128            zaggpoc2 = zaggtmp * tr(ji,jj,jk,jppoc,Kbb)
[15459]129                 
[12377]130            ! Part II : Differential settling
[15459]131            ! The stickiness has been assumed to be 0.1
132   
[12377]133            !  Aggregation of small into large particles
134            zaggtmp =  47.1 * xstep * tr(ji,jj,jk,jpgoc,Kbb)
135            zaggpoc3 = zaggtmp * tr(ji,jj,jk,jppoc,Kbb)
136            zaggtmp =  3.3  * xstep * tr(ji,jj,jk,jppoc,Kbb)
137            zaggpoc4 = zaggtmp * tr(ji,jj,jk,jppoc,Kbb)
[7176]138
[15459]139            zaggpoc = zaggpoc1 + zaggpoc2 + zaggpoc3 + zaggpoc4
[12377]140            zaggpon = zaggpoc * tr(ji,jj,jk,jppon,Kbb) / ( tr(ji,jj,jk,jppoc,Kbb) + rtrn)
141            zaggpop = zaggpoc * tr(ji,jj,jk,jppop,Kbb) / ( tr(ji,jj,jk,jppoc,Kbb) + rtrn)
[15459]142            zaggfe  = zaggpoc * tr(ji,jj,jk,jpsfe,Kbb) / ( tr(ji,jj,jk,jppoc,Kbb)  + rtrn )
[7176]143
[12377]144            ! Aggregation of DOC to POC :
145            ! 1st term is shear aggregation of DOC-DOC
146            ! 2nd term is shear aggregation of DOC-POC
147            ! 3rd term is differential settling of DOC-POC
[15459]148            ! 1/3 of DOC is supposed to experience aggregation (HMW)
149            zaggtmp = ( ( 0.37 * 0.3 * tr(ji,jj,jk,jpdoc,Kbb) + 20.5 * tr(ji,jj,jk,jppoc,Kbb) ) * zfact       &
150            &            + 0.15 * xstep * tr(ji,jj,jk,jppoc,Kbb) )
[12377]151            zaggdoc  = zaggtmp * 0.3 * tr(ji,jj,jk,jpdoc,Kbb)
152            zaggdon  = zaggtmp * 0.3 * tr(ji,jj,jk,jpdon,Kbb)
153            zaggdop  = zaggtmp * 0.3 * tr(ji,jj,jk,jpdop,Kbb)
[7176]154
[12377]155            ! transfer of DOC to GOC :
156            ! 1st term is shear aggregation
157            ! 2nd term is differential settling
[15459]158            ! 1/3 of DOC is supposed to experience aggregation (HMW)
159            zaggtmp = 655.4 * zfact * tr(ji,jj,jk,jpgoc,Kbb)
[12377]160            zaggdoc2 = zaggtmp * 0.3 * tr(ji,jj,jk,jpdoc,Kbb)
161            zaggdon2 = zaggtmp * 0.3 * tr(ji,jj,jk,jpdon,Kbb)
162            zaggdop2 = zaggtmp * 0.3 * tr(ji,jj,jk,jpdop,Kbb)
[7176]163
[12377]164            ! tranfer of DOC to POC due to brownian motion
[15459]165            ! 1/3 of DOC is supposed to experience aggregation (HMW)
166            zaggtmp = ( 260.2 * 0.3 * tr(ji,jj,jk,jpdoc,Kbb) +  418.5 * tr(ji,jj,jk,jppoc,Kbb) ) * xstep
[12377]167            zaggdoc3 =  zaggtmp * 0.3 * tr(ji,jj,jk,jpdoc,Kbb)
168            zaggdon3 =  zaggtmp * 0.3 * tr(ji,jj,jk,jpdon,Kbb)
169            zaggdop3 =  zaggtmp * 0.3 * tr(ji,jj,jk,jpdop,Kbb)
170
[15459]171
[12377]172            !  Update the trends
173            tr(ji,jj,jk,jppoc,Krhs) = tr(ji,jj,jk,jppoc,Krhs) - zaggpoc + zaggdoc + zaggdoc3
174            tr(ji,jj,jk,jppon,Krhs) = tr(ji,jj,jk,jppon,Krhs) - zaggpon + zaggdon + zaggdon3
175            tr(ji,jj,jk,jppop,Krhs) = tr(ji,jj,jk,jppop,Krhs) - zaggpop + zaggdop + zaggdop3
176            tr(ji,jj,jk,jpgoc,Krhs) = tr(ji,jj,jk,jpgoc,Krhs) + zaggpoc + zaggdoc2
177            tr(ji,jj,jk,jpgon,Krhs) = tr(ji,jj,jk,jpgon,Krhs) + zaggpon + zaggdon2
178            tr(ji,jj,jk,jpgop,Krhs) = tr(ji,jj,jk,jpgop,Krhs) + zaggpop + zaggdop2
179            tr(ji,jj,jk,jpsfe,Krhs) = tr(ji,jj,jk,jpsfe,Krhs) - zaggfe
180            tr(ji,jj,jk,jpbfe,Krhs) = tr(ji,jj,jk,jpbfe,Krhs) + zaggfe
181            tr(ji,jj,jk,jpdoc,Krhs) = tr(ji,jj,jk,jpdoc,Krhs) - zaggdoc - zaggdoc2 - zaggdoc3
182            tr(ji,jj,jk,jpdon,Krhs) = tr(ji,jj,jk,jpdon,Krhs) - zaggdon - zaggdon2 - zaggdon3
183            tr(ji,jj,jk,jpdop,Krhs) = tr(ji,jj,jk,jpdop,Krhs) - zaggdop - zaggdop2 - zaggdop3
184            !
185            conspoc(ji,jj,jk) = conspoc(ji,jj,jk) - zaggpoc + zaggdoc + zaggdoc3
186            prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zaggpoc + zaggdoc2
187            !
188         END_3D
[7176]189         !
190      ENDIF
[7162]191      !
[12377]192      IF(sn_cfctl%l_prttrc)   THEN  ! print mean trends (used for debugging)
[7162]193         WRITE(charout, FMT="('agg')")
[13286]194         CALL prt_ctl_info( charout, cdcomp = 'top' )
195         CALL prt_ctl(tab4d_1=tr(:,:,:,:,Krhs), mask1=tmask, clinfo=ctrcnm)
[7162]196      ENDIF
197      !
[9124]198      IF( ln_timing )   CALL timing_stop('p4z_agg')
[7162]199      !
200   END SUBROUTINE p4z_agg
201
202   !!======================================================================
203END MODULE p4zagg
Note: See TracBrowser for help on using the repository browser.