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 trunk/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z – NEMO

source: trunk/NEMOGCM/NEMO/TOP_SRC/PISCES/P4Z/p4zagg.F90 @ 7698

Last change on this file since 7698 was 7698, checked in by mocavero, 7 years ago

update trunk with OpenMP parallelization
File size: 9.5 KB
Line 
1MODULE p4zagg
2   !!======================================================================
3   !!                         ***  MODULE p4zagg  ***
4   !! TOP :  PISCES  aggregation of particles
5   !!======================================================================
6   !! History :   1.0  !  2004     (O. Aumont) Original code
7   !!             2.0  !  2007-12  (C. Ethe, G. Madec)  F90
8   !!             3.4  !  2011-06  (O. Aumont, C. Ethe) Change aggregation formula
9   !!             3.5  !  2012-07  (O. Aumont) Introduce potential time-splitting
10   !!             3.6  !  2015-05  (O. Aumont) PISCES quota
11   !!----------------------------------------------------------------------
12   !!----------------------------------------------------------------------
13   !!   p4z_agg       :  Compute aggregation of particles
14   !!----------------------------------------------------------------------
15   USE oce_trc         !  shared variables between ocean and passive tracers
16   USE trc             !  passive tracers common variables
17   USE sms_pisces      !  PISCES Source Minus Sink variables
18   USE prtctl_trc      !  print control for debugging
19
20   IMPLICIT NONE
21   PRIVATE
22
23   PUBLIC   p4z_agg         ! called in p4zbio.F90
24
25   !!----------------------------------------------------------------------
26   !! NEMO/TOP 3.3 , NEMO Consortium (2010)
27   !! $Id: p4zagg.F90 3160 2011-11-20 14:27:18Z cetlod $
28   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
29   !!----------------------------------------------------------------------
30CONTAINS
31
32   SUBROUTINE p4z_agg ( kt, knt )
33      !!---------------------------------------------------------------------
34      !!                     ***  ROUTINE p4z_agg  ***
35      !!
36      !! ** Purpose :   Compute aggregation of particles
37      !!
38      !! ** Method  : - ???
39      !!---------------------------------------------------------------------
40      INTEGER, INTENT(in) :: kt, knt
41      INTEGER  ::   ji, jj, jk
42      REAL(wp) ::   zagg, zagg1, zagg2, zagg3, zagg4
43      REAL(wp) ::   zaggpoc1, zaggpoc2, zaggpoc3, zaggpoc4
44      REAL(wp) ::   zaggpoc , zaggfe, zaggdoc, zaggdoc2, zaggdoc3
45      REAL(wp) ::   zaggpon , zaggdon, zaggdon2, zaggdon3
46      REAL(wp) ::   zaggpop, zaggdop, zaggdop2, zaggdop3
47      REAL(wp) ::   zaggtmp, zfact, zmax
48      CHARACTER (len=25) :: charout
49      !!---------------------------------------------------------------------
50      !
51      IF( nn_timing == 1 )  CALL timing_start('p4z_agg')
52
53      !
54      !  Exchange between organic matter compartments due to coagulation/disaggregation
55      !  ---------------------------------------------------
56      IF( ln_p4z ) THEN
57         !
58!$OMP PARALLEL DO schedule(static) private(jk,jj,ji,zfact,zagg1,zagg2,zagg3,zagg4,zagg,zaggfe,zaggdoc,zaggdoc2,zaggdoc3)
59         DO jk = 1, jpkm1
60            DO jj = 1, jpj
61               DO ji = 1, jpi
62                  !
63                  zfact = xstep * xdiss(ji,jj,jk)
64                  !  Part I : Coagulation dependent on turbulence
65                  zagg1 = 25.9  * zfact * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jppoc)
66                  zagg2 = 4452. * zfact * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpgoc)
67
68                  ! Part II : Differential settling
69
70                  !  Aggregation of small into large particles
71                  zagg3 =  47.1 * xstep * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jpgoc)
72                  zagg4 =  3.3  * xstep * trb(ji,jj,jk,jppoc) * trb(ji,jj,jk,jppoc)
73
74                  zagg   = zagg1 + zagg2 + zagg3 + zagg4
75                  zaggfe = zagg * trb(ji,jj,jk,jpsfe) / ( trb(ji,jj,jk,jppoc) + rtrn )
76
77                  ! Aggregation of DOC to POC :
78                  ! 1st term is shear aggregation of DOC-DOC
79                  ! 2nd term is shear aggregation of DOC-POC
80                  ! 3rd term is differential settling of DOC-POC
81                  zaggdoc  = ( ( 0.369 * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4 * trb(ji,jj,jk,jppoc) ) * zfact       &
82                  &            + 2.4 * xstep * trb(ji,jj,jk,jppoc) ) * 0.3 * trb(ji,jj,jk,jpdoc)
83                  ! transfer of DOC to GOC :
84                  ! 1st term is shear aggregation
85                  ! 2nd term is differential settling
86                  zaggdoc2 = ( 3.53E3 * zfact + 0.1 * xstep ) * trb(ji,jj,jk,jpgoc) * 0.3 * trb(ji,jj,jk,jpdoc)
87                  ! tranfer of DOC to POC due to brownian motion
88                  zaggdoc3 =  114. * 0.3 * trb(ji,jj,jk,jpdoc) *xstep * 0.3 * trb(ji,jj,jk,jpdoc)
89
90                  !  Update the trends
91                  tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zagg + zaggdoc + zaggdoc3
92                  tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zagg + zaggdoc2
93                  tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zaggfe
94                  tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zaggfe
95                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zaggdoc - zaggdoc2 - zaggdoc3
96                  !
97                  conspoc(ji,jj,jk) = conspoc(ji,jj,jk) - zagg + zaggdoc + zaggdoc3
98                  prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zagg + zaggdoc2
99                  !
100               END DO
101            END DO
102         END DO
103      ELSE    ! ln_p5z
104        !
105!$OMP PARALLEL DO schedule(static) private(jk,jj,ji,zfact,zaggtmp,zaggfe,zaggpoc,zaggpoc1,zaggpoc2,zaggpoc3,zaggpoc4) &
106!$OMP& private(zaggpon,zaggpop,zaggdoc,zaggdon,zaggdop,zaggdoc2,zaggdon2,zaggdop2,zaggdoc3,zaggdon3,zaggdop3)
107         DO jk = 1, jpkm1
108            DO jj = 1, jpj
109               DO ji = 1, jpi
110                  !
111                  zfact = xstep * xdiss(ji,jj,jk)
112                  !  Part I : Coagulation dependent on turbulence
113                  zaggtmp = 25.9  * zfact * trb(ji,jj,jk,jppoc)
114                  zaggpoc1 = zaggtmp * trb(ji,jj,jk,jppoc)
115                  zaggtmp = 4452. * zfact * trb(ji,jj,jk,jpgoc)
116                  zaggpoc2 = zaggtmp * trb(ji,jj,jk,jppoc)
117
118                  ! Part II : Differential settling
119   
120                  !  Aggregation of small into large particles
121                  zaggtmp =  47.1 * xstep * trb(ji,jj,jk,jpgoc)
122                  zaggpoc3 = zaggtmp * trb(ji,jj,jk,jppoc)
123                  zaggtmp =  3.3  * xstep * trb(ji,jj,jk,jppoc)
124                  zaggpoc4 = zaggtmp * trb(ji,jj,jk,jppoc)
125
126                  zaggpoc   = zaggpoc1 + zaggpoc2 + zaggpoc3 + zaggpoc4
127                  zaggpon = zaggpoc * trb(ji,jj,jk,jppon) / ( trb(ji,jj,jk,jppoc) + rtrn)
128                  zaggpop = zaggpoc * trb(ji,jj,jk,jppop) / ( trb(ji,jj,jk,jppoc) + rtrn)
129                  zaggfe = zaggpoc * trb(ji,jj,jk,jpsfe) / ( trb(ji,jj,jk,jppoc)  + rtrn )
130
131                  ! Aggregation of DOC to POC :
132                  ! 1st term is shear aggregation of DOC-DOC
133                  ! 2nd term is shear aggregation of DOC-POC
134                  ! 3rd term is differential settling of DOC-POC
135                  zaggtmp = ( ( 0.369 * 0.3 * trb(ji,jj,jk,jpdoc) + 102.4 * trb(ji,jj,jk,jppoc) ) * zfact       &
136                  &            + 2.4 * xstep * trb(ji,jj,jk,jppoc) )
137                  zaggdoc  = zaggtmp * 0.3 * trb(ji,jj,jk,jpdoc)
138                  zaggdon  = zaggtmp * 0.3 * trb(ji,jj,jk,jpdon)
139                  zaggdop  = zaggtmp * 0.3 * trb(ji,jj,jk,jpdop)
140
141                  ! transfer of DOC to GOC :
142                  ! 1st term is shear aggregation
143                  ! 2nd term is differential settling
144                  zaggtmp = ( 3.53E3 * zfact + 0.1 * xstep ) * trb(ji,jj,jk,jpgoc)
145                  zaggdoc2 = zaggtmp * 0.3 * trb(ji,jj,jk,jpdoc)
146                  zaggdon2 = zaggtmp * 0.3 * trb(ji,jj,jk,jpdon)
147                  zaggdop2 = zaggtmp * 0.3 * trb(ji,jj,jk,jpdop)
148
149                  ! tranfer of DOC to POC due to brownian motion
150                  zaggtmp = ( 114. * 0.3 * trb(ji,jj,jk,jpdoc) ) * xstep
151                  zaggdoc3 =  zaggtmp * 0.3 * trb(ji,jj,jk,jpdoc)
152                  zaggdon3 =  zaggtmp * 0.3 * trb(ji,jj,jk,jpdon)
153                  zaggdop3 =  zaggtmp * 0.3 * trb(ji,jj,jk,jpdop)
154
155                  !  Update the trends
156                  tra(ji,jj,jk,jppoc) = tra(ji,jj,jk,jppoc) - zaggpoc + zaggdoc + zaggdoc3
157                  tra(ji,jj,jk,jppon) = tra(ji,jj,jk,jppon) - zaggpon + zaggdon + zaggdon3
158                  tra(ji,jj,jk,jppop) = tra(ji,jj,jk,jppop) - zaggpop + zaggdop + zaggdop3
159                  tra(ji,jj,jk,jpgoc) = tra(ji,jj,jk,jpgoc) + zaggpoc + zaggdoc2
160                  tra(ji,jj,jk,jpgon) = tra(ji,jj,jk,jpgon) + zaggpon + zaggdon2
161                  tra(ji,jj,jk,jpgop) = tra(ji,jj,jk,jpgop) + zaggpop + zaggdop2
162                  tra(ji,jj,jk,jpsfe) = tra(ji,jj,jk,jpsfe) - zaggfe
163                  tra(ji,jj,jk,jpbfe) = tra(ji,jj,jk,jpbfe) + zaggfe
164                  tra(ji,jj,jk,jpdoc) = tra(ji,jj,jk,jpdoc) - zaggdoc - zaggdoc2 - zaggdoc3
165                  tra(ji,jj,jk,jpdon) = tra(ji,jj,jk,jpdon) - zaggdon - zaggdon2 - zaggdon3
166                  tra(ji,jj,jk,jpdop) = tra(ji,jj,jk,jpdop) - zaggdop - zaggdop2 - zaggdop3
167                  !
168                  conspoc(ji,jj,jk) = conspoc(ji,jj,jk) - zaggpoc + zaggdoc + zaggdoc3
169                  prodgoc(ji,jj,jk) = prodgoc(ji,jj,jk) + zaggpoc + zaggdoc2
170                  !
171               END DO
172            END DO
173         END DO
174         !
175      ENDIF
176      !
177      IF(ln_ctl)   THEN  ! print mean trends (used for debugging)
178         WRITE(charout, FMT="('agg')")
179         CALL prt_ctl_trc_info(charout)
180         CALL prt_ctl_trc(tab4d=tra, mask=tmask, clinfo=ctrcnm)
181      ENDIF
182      !
183      IF( nn_timing == 1 )  CALL timing_stop('p4z_agg')
184      !
185   END SUBROUTINE p4z_agg
186
187   !!======================================================================
188END MODULE p4zagg
Note: See TracBrowser for help on using the repository browser.