source: LMDZ6/branches/Amaury_dev/libf/dyn3d_common/advzp.f90 @ 5195

Last change on this file since 5195 was 5159, checked in by abarral, 5 months ago

Put dimensions.h and paramet.h into modules

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 13.1 KB
RevLine 
[524]1! $Header$
[5099]2
[5136]3SUBROUTINE ADVZP(LIMIT, DTZ, W, SM, S0, SSX, SY, SZ &
4        , SSXX, SSXY, SSXZ, SYY, SYZ, SZZ, ntra)
5  USE lmdz_comgeom
[524]6
[5159]7  USE lmdz_dimensions, ONLY: iim, jjm, llm, ndm
8  USE lmdz_paramet
[5105]9  IMPLICIT NONE
[524]10
[5105]11  !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
12  !                                                                 C
13  !  second-order moments (SOM) advection of tracer in Z direction  C
14  !                                                                 C
15  !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
16  !                                                                 C
17  !  Source : Pascal Simon ( Meteo, CNRM )                          C
18  !  Adaptation : A.A. (LGGE)                                       C
19  !  Derniere Modif : 19/11/95 LAST                                 C
20  !                                                                 C
21  !  sont les arguments d'entree pour le s-pg                       C
22  !                                                                 C
23  !  argument de sortie du s-pg                                     C
24  !                                                                 C
25  !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
26  !CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
[5159]27
[5105]28  ! Rem : Probleme aux poles il faut reecrire ce cas specifique
29  !    Attention au sens de l'indexation
30  !
[524]31
[5159]32
[5105]33  !  parametres principaux du modele
34  !
[5159]35
36
37
[5105]38  !  Arguments :
39  !  ----------
40  !  dty : frequence fictive d'appel du transport
41  !  parbu,pbarv : flux de masse en x et y en Pa.m2.s-1
[5159]42
[5136]43  INTEGER :: lon, lat, niv
44  INTEGER :: i, j, jv, k, kp, l, lp
45  INTEGER :: ntra
46  ! PARAMETER (ntra = 1)
[5159]47
[5136]48  REAL :: dtz
49  REAL :: w (iip1, jjp1, llm)
[5159]50
[5105]51  !  moments: SM  total mass in each grid box
52  !       S0  mass of tracer in each grid box
53  !       Si  1rst order moment in i direction
[5159]54
[5136]55  REAL :: SM(iip1, jjp1, llm) &
56          , S0(iip1, jjp1, llm, ntra)
57  REAL :: SSX(iip1, jjp1, llm, ntra) &
58          , SY(iip1, jjp1, llm, ntra) &
59          , SZ(iip1, jjp1, llm, ntra) &
60          , SSXX(iip1, jjp1, llm, ntra) &
61          , SSXY(iip1, jjp1, llm, ntra) &
62          , SSXZ(iip1, jjp1, llm, ntra) &
63          , SYY(iip1, jjp1, llm, ntra) &
64          , SYZ(iip1, jjp1, llm, ntra) &
65          , SZZ(iip1, jjp1, llm, ntra)
[5159]66
[5105]67  !  Local :
68  !  -------
[5159]69
[5105]70  !  mass fluxes across the boundaries (UGRI,VGRI,WGRI)
71  !  mass fluxes in kg
72  !  declaration :
[5159]73
[5136]74  REAL :: WGRI(iip1, jjp1, 0:llm)
[524]75
[5105]76  ! Rem : UGRI et VGRI ne sont pas utilises dans
77  !  cette SUBROUTINE ( advection en z uniquement )
78  !  Rem 2 :le dimensionnement de VGRI depend de celui de pbarv
[5136]79  ! attention a celui de WGRI
[5159]80
[5105]81  !  the moments F are similarly defined and used as temporary
82  !  storage for portions of the grid boxes in transit
[5159]83
[5105]84  !  the moments Fij are used as temporary storage for
85  !  portions of the grid boxes in transit at the current level
[5159]86
[5105]87  !  work arrays
[5159]88
89
[5136]90  REAL :: F0(iim, llm, ntra), FM(iim, llm)
91  REAL :: FX(iim, llm, ntra), FY(iim, llm, ntra)
92  REAL :: FZ(iim, llm, ntra)
93  REAL :: FXX(iim, llm, ntra), FXY(iim, llm, ntra)
94  REAL :: FXZ(iim, llm, ntra), FYY(iim, llm, ntra)
95  REAL :: FYZ(iim, llm, ntra), FZZ(iim, llm, ntra)
[5105]96  REAL :: S00(ntra)
97  REAL :: SM0             ! Just temporal variable
[5159]98
[5105]99  !  work arrays
[5159]100
[5136]101  REAL :: ALF(iim), ALF1(iim)
102  REAL :: ALFQ(iim), ALF1Q(iim)
103  REAL :: ALF2(iim), ALF3(iim)
[5105]104  REAL :: ALF4(iim)
105  REAL :: TEMPTM          ! Just temporal variable
[5136]106  REAL :: SLPMAX, S1MAX, S1NEW, S2NEW
[5159]107
[5136]108  REAL :: sqi, sqf
[5105]109  LOGICAL :: LIMIT
[524]110
[5105]111  lon = iim         ! rem : Il est possible qu'un pbl. arrive ici
112  lat = jjp1        ! a cause des dim. differentes entre les
113  niv = llm         !       tab. S et VGRI
[524]114
[5105]115  !-----------------------------------------------------------------
116  ! *** Test : diag de la qtite totale de traceur dans
[5136]117  ! l'atmosphere avant l'advection en Y
[5159]118
[5105]119  sqi = 0.
120  sqf = 0.
[5159]121
[5136]122  DO l = 1, llm
123    DO j = 1, jjp1
124      DO i = 1, iim
125        sqi = sqi + S0(i, j, l, ntra)
126      END DO
127    END DO
[5105]128  END DO
[5136]129  PRINT*, '---------- DIAG DANS ADVZP - ENTREE --------'
130  PRINT*, 'sqi=', sqi
[524]131
[5105]132  !-----------------------------------------------------------------
133  !  Interface : adaptation nouveau modele
134  !  -------------------------------------
[5159]135
[5105]136  !  Conversion des flux de masses en kg
137
[5136]138  DO l = 1, llm
139    DO j = 1, jjp1
140      DO i = 1, iip1
141        wgri (i, j, llm + 1 - l) = w (i, j, l)
142      END DO
143    END DO
[5105]144  END DO
[5158]145  DO j = 1, jjp1
146    DO i = 1, iip1
[5136]147      wgri(i, j, 0) = 0.
148    enddo
[5105]149  enddo
[5159]150
[5105]151  !AA rem : Je ne suis pas sur du signe
152  !AA       Je ne suis pas sur pour le 0:llm
[5159]153
[5105]154  !-----------------------------------------------------------------
155  !---------------------- START HERE -------------------------------
[5159]156
[5105]157  !  boucle sur les latitudes
[5159]158
[5136]159  DO K = 1, LAT
[5159]160
[5136]161    !  place limits on appropriate moments before transport
162    !  (if flux-limiting is to be applied)
[5159]163
[5136]164    IF(.NOT.LIMIT) GO TO 101
[5159]165
[5136]166    DO JV = 1, NTRA
167      DO L = 1, NIV
168        DO I = 1, LON
169          IF(S0(I, K, L, JV)>0.) THEN
170            SLPMAX = S0(I, K, L, JV)
171            S1MAX = 1.5 * SLPMAX
172            S1NEW = AMIN1(S1MAX, AMAX1(-S1MAX, SZ(I, K, L, JV)))
173            S2NEW = AMIN1(2. * SLPMAX - ABS(S1NEW) / 3., &
174                    AMAX1(ABS(S1NEW) - SLPMAX, SZZ(I, K, L, JV)))
175            SZ (I, K, L, JV) = S1NEW
176            SZZ(I, K, L, JV) = S2NEW
177            SSXZ(I, K, L, JV) = AMIN1(SLPMAX, AMAX1(-SLPMAX, SSXZ(I, K, L, JV)))
178            SYZ(I, K, L, JV) = AMIN1(SLPMAX, AMAX1(-SLPMAX, SYZ(I, K, L, JV)))
179          ELSE
180            SZ (I, K, L, JV) = 0.
181            SZZ(I, K, L, JV) = 0.
182            SSXZ(I, K, L, JV) = 0.
183            SYZ(I, K, L, JV) = 0.
184          ENDIF
185        END DO
186      END DO
187    END DO
[5159]188
[5136]189    101   CONTINUE
[5159]190
[5136]191    !  boucle sur les niveaux intercouches de 1 a NIV-1
192    !   (flux nul au sommet L=0 et a la base L=NIV)
[5159]193
[5136]194    !  calculate flux and moments between adjacent boxes
195    ! (flux from LP to L if WGRI(L).lt.0, from L to LP if WGRI(L).gt.0)
196    !  1- create temporary moments/masses for partial boxes in transit
197    !  2- reajusts moments remaining in the box
[5159]198
[5136]199    DO L = 1, NIV - 1
200      LP = L + 1
[5159]201
[5136]202      DO I = 1, LON
[5159]203
[5136]204        IF(WGRI(I, K, L)<0.) THEN
205          FM(I, L) = -WGRI(I, K, L) * DTZ
206          ALF(I) = FM(I, L) / SM(I, K, LP)
207          SM(I, K, LP) = SM(I, K, LP) - FM(I, L)
[5105]208        ELSE
[5136]209          FM(I, L) = WGRI(I, K, L) * DTZ
210          ALF(I) = FM(I, L) / SM(I, K, L)
211          SM(I, K, L) = SM(I, K, L) - FM(I, L)
[5105]212        ENDIF
[5159]213
[5136]214        ALFQ (I) = ALF(I) * ALF(I)
215        ALF1 (I) = 1. - ALF(I)
216        ALF1Q(I) = ALF1(I) * ALF1(I)
217        ALF2 (I) = ALF1(I) - ALF(I)
218        ALF3 (I) = ALF(I) * ALFQ(I)
219        ALF4 (I) = ALF1(I) * ALF1Q(I)
[5159]220
[5136]221      END DO
[5159]222
[5136]223      DO JV = 1, NTRA
224        DO I = 1, LON
[5159]225
[5136]226          IF(WGRI(I, K, L)<0.) THEN
[5159]227
[5136]228            F0 (I, L, JV) = ALF (I) * (S0(I, K, LP, JV) - ALF1(I) * &
229                    (SZ(I, K, LP, JV) - ALF2(I) * SZZ(I, K, LP, JV)))
230            FZ (I, L, JV) = ALFQ(I) * (SZ(I, K, LP, JV) - 3. * ALF1(I) * SZZ(I, K, LP, JV))
231            FZZ(I, L, JV) = ALF3(I) * SZZ(I, K, LP, JV)
232            FXZ(I, L, JV) = ALFQ(I) * SSXZ(I, K, LP, JV)
233            FYZ(I, L, JV) = ALFQ(I) * SYZ(I, K, LP, JV)
234            FX (I, L, JV) = ALF (I) * (SSX(I, K, LP, JV) - ALF1(I) * SSXZ(I, K, LP, JV))
235            FY (I, L, JV) = ALF (I) * (SY(I, K, LP, JV) - ALF1(I) * SYZ(I, K, LP, JV))
236            FXX(I, L, JV) = ALF (I) * SSXX(I, K, LP, JV)
237            FXY(I, L, JV) = ALF (I) * SSXY(I, K, LP, JV)
238            FYY(I, L, JV) = ALF (I) * SYY(I, K, LP, JV)
[5159]239
[5136]240            S0 (I, K, LP, JV) = S0 (I, K, LP, JV) - F0 (I, L, JV)
241            SZ (I, K, LP, JV) = ALF1Q(I) &
242                    * (SZ(I, K, LP, JV) + 3. * ALF(I) * SZZ(I, K, LP, JV))
243            SZZ(I, K, LP, JV) = ALF4 (I) * SZZ(I, K, LP, JV)
244            SSXZ(I, K, LP, JV) = ALF1Q(I) * SSXZ(I, K, LP, JV)
245            SYZ(I, K, LP, JV) = ALF1Q(I) * SYZ(I, K, LP, JV)
246            SSX (I, K, LP, JV) = SSX (I, K, LP, JV) - FX (I, L, JV)
247            SY (I, K, LP, JV) = SY (I, K, LP, JV) - FY (I, L, JV)
248            SSXX(I, K, LP, JV) = SSXX(I, K, LP, JV) - FXX(I, L, JV)
249            SSXY(I, K, LP, JV) = SSXY(I, K, LP, JV) - FXY(I, L, JV)
250            SYY(I, K, LP, JV) = SYY(I, K, LP, JV) - FYY(I, L, JV)
[5159]251
[5136]252          ELSE
[5159]253
[5136]254            F0 (I, L, JV) = ALF (I) * (S0(I, K, L, JV) &
255                    + ALF1(I) * (SZ(I, K, L, JV) + ALF2(I) * SZZ(I, K, L, JV)))
256            FZ (I, L, JV) = ALFQ(I) * (SZ(I, K, L, JV) + 3. * ALF1(I) * SZZ(I, K, L, JV))
257            FZZ(I, L, JV) = ALF3(I) * SZZ(I, K, L, JV)
258            FXZ(I, L, JV) = ALFQ(I) * SSXZ(I, K, L, JV)
259            FYZ(I, L, JV) = ALFQ(I) * SYZ(I, K, L, JV)
260            FX (I, L, JV) = ALF (I) * (SSX(I, K, L, JV) + ALF1(I) * SSXZ(I, K, L, JV))
261            FY (I, L, JV) = ALF (I) * (SY(I, K, L, JV) + ALF1(I) * SYZ(I, K, L, JV))
262            FXX(I, L, JV) = ALF (I) * SSXX(I, K, L, JV)
263            FXY(I, L, JV) = ALF (I) * SSXY(I, K, L, JV)
264            FYY(I, L, JV) = ALF (I) * SYY(I, K, L, JV)
[5159]265
[5136]266            S0 (I, K, L, JV) = S0 (I, K, L, JV) - F0(I, L, JV)
267            SZ (I, K, L, JV) = ALF1Q(I) * (SZ(I, K, L, JV) - 3. * ALF(I) * SZZ(I, K, L, JV))
268            SZZ(I, K, L, JV) = ALF4 (I) * SZZ(I, K, L, JV)
269            SSXZ(I, K, L, JV) = ALF1Q(I) * SSXZ(I, K, L, JV)
270            SYZ(I, K, L, JV) = ALF1Q(I) * SYZ(I, K, L, JV)
271            SSX (I, K, L, JV) = SSX (I, K, L, JV) - FX (I, L, JV)
272            SY (I, K, L, JV) = SY (I, K, L, JV) - FY (I, L, JV)
273            SSXX(I, K, L, JV) = SSXX(I, K, L, JV) - FXX(I, L, JV)
274            SSXY(I, K, L, JV) = SSXY(I, K, L, JV) - FXY(I, L, JV)
275            SYY(I, K, L, JV) = SYY(I, K, L, JV) - FYY(I, L, JV)
[5159]276
[5136]277          ENDIF
[5159]278
[5136]279        END DO
280      END DO
[5159]281
[5136]282    END DO
[5159]283
[5136]284    !  puts the temporary moments Fi into appropriate neighboring boxes
[5159]285
[5136]286    DO L = 1, NIV - 1
287      LP = L + 1
[5159]288
[5136]289      DO I = 1, LON
[5159]290
[5136]291        IF(WGRI(I, K, L)<0.) THEN
292          SM(I, K, L) = SM(I, K, L) + FM(I, L)
293          ALF(I) = FM(I, L) / SM(I, K, L)
294        ELSE
295          SM(I, K, LP) = SM(I, K, LP) + FM(I, L)
296          ALF(I) = FM(I, L) / SM(I, K, LP)
297        ENDIF
[5159]298
[5136]299        ALF1(I) = 1. - ALF(I)
300        ALFQ(I) = ALF(I) * ALF(I)
301        ALF1Q(I) = ALF1(I) * ALF1(I)
302        ALF2(I) = ALF(I) * ALF1(I)
303        ALF3(I) = ALF1(I) - ALF(I)
[5159]304
[5136]305      END DO
[5159]306
[5136]307      DO JV = 1, NTRA
308        DO I = 1, LON
[5159]309
[5136]310          IF(WGRI(I, K, L)<0.) THEN
[5159]311
[5136]312            TEMPTM = -ALF(I) * S0(I, K, L, JV) + ALF1(I) * F0(I, L, JV)
313            S0 (I, K, L, JV) = S0(I, K, L, JV) + F0(I, L, JV)
314            SZZ(I, K, L, JV) = ALFQ(I) * FZZ(I, L, JV) + ALF1Q(I) * SZZ(I, K, L, JV) &
315                    + 5. * (ALF2(I) * (FZ(I, L, JV) - SZ(I, K, L, JV)) + ALF3(I) * TEMPTM)
316            SZ (I, K, L, JV) = ALF (I) * FZ (I, L, JV) + ALF1 (I) * SZ (I, K, L, JV) &
317                    + 3. * TEMPTM
318            SSXZ(I, K, L, JV) = ALF (I) * FXZ(I, L, JV) + ALF1 (I) * SSXZ(I, K, L, JV) &
319                    + 3. * (ALF1(I) * FX (I, L, JV) - ALF  (I) * SSX (I, K, L, JV))
320            SYZ(I, K, L, JV) = ALF (I) * FYZ(I, L, JV) + ALF1 (I) * SYZ(I, K, L, JV) &
321                    + 3. * (ALF1(I) * FY (I, L, JV) - ALF  (I) * SY (I, K, L, JV))
322            SSX (I, K, L, JV) = SSX (I, K, L, JV) + FX (I, L, JV)
323            SY (I, K, L, JV) = SY (I, K, L, JV) + FY (I, L, JV)
324            SSXX(I, K, L, JV) = SSXX(I, K, L, JV) + FXX(I, L, JV)
325            SSXY(I, K, L, JV) = SSXY(I, K, L, JV) + FXY(I, L, JV)
326            SYY(I, K, L, JV) = SYY(I, K, L, JV) + FYY(I, L, JV)
[5159]327
[5136]328          ELSE
[5159]329
[5136]330            TEMPTM = ALF(I) * S0(I, K, LP, JV) - ALF1(I) * F0(I, L, JV)
331            S0 (I, K, LP, JV) = S0(I, K, LP, JV) + F0(I, L, JV)
332            SZZ(I, K, LP, JV) = ALFQ(I) * FZZ(I, L, JV) + ALF1Q(I) * SZZ(I, K, LP, JV) &
333                    + 5. * (ALF2(I) * (SZ(I, K, LP, JV) - FZ(I, L, JV)) - ALF3(I) * TEMPTM)
334            SZ (I, K, LP, JV) = ALF (I) * FZ(I, L, JV) + ALF1(I) * SZ(I, K, LP, JV) &
335                    + 3. * TEMPTM
336            SSXZ(I, K, LP, JV) = ALF(I) * FXZ(I, L, JV) + ALF1(I) * SSXZ(I, K, LP, JV) &
337                    + 3. * (ALF(I) * SSX(I, K, LP, JV) - ALF1(I) * FX(I, L, JV))
338            SYZ(I, K, LP, JV) = ALF(I) * FYZ(I, L, JV) + ALF1(I) * SYZ(I, K, LP, JV) &
339                    + 3. * (ALF(I) * SY(I, K, LP, JV) - ALF1(I) * FY(I, L, JV))
340            SSX (I, K, LP, JV) = SSX (I, K, LP, JV) + FX (I, L, JV)
341            SY (I, K, LP, JV) = SY (I, K, LP, JV) + FY (I, L, JV)
342            SSXX(I, K, LP, JV) = SSXX(I, K, LP, JV) + FXX(I, L, JV)
343            SSXY(I, K, LP, JV) = SSXY(I, K, LP, JV) + FXY(I, L, JV)
344            SYY(I, K, LP, JV) = SYY(I, K, LP, JV) + FYY(I, L, JV)
[5159]345
[5136]346          ENDIF
[5159]347
[5136]348        END DO
349      END DO
[5159]350
[5136]351    END DO
[5159]352
[5136]353    !  fin de la boucle principale sur les latitudes
[5159]354
[5105]355  END DO
[5159]356
[5136]357  DO l = 1, llm
358    DO j = 1, jjp1
359      SM(iip1, j, l) = SM(1, j, l)
360      S0(iip1, j, l, ntra) = S0(1, j, l, ntra)
361      SSX(iip1, j, l, ntra) = SSX(1, j, l, ntra)
362      SY(iip1, j, l, ntra) = SY(1, j, l, ntra)
363      SZ(iip1, j, l, ntra) = SZ(1, j, l, ntra)
364    ENDDO
[5105]365  ENDDO
[5136]366  ! C-------------------------------------------------------------
367  ! *** Test : diag de la qqtite totale de tarceur
368  ! dans l'atmosphere avant l'advection en z
369  DO l = 1, llm
370    DO j = 1, jjp1
371      DO i = 1, iim
372        sqf = sqf + S0(i, j, l, ntra)
373      ENDDO
374    ENDDO
[5105]375  ENDDO
[5136]376  PRINT*, '-------- DIAG DANS ADVZ - SORTIE ---------'
377  PRINT*, 'sqf=', sqf
[5105]378
379  RETURN
380END SUBROUTINE ADVZP
Note: See TracBrowser for help on using the repository browser.