Changeset 3198 for LMDZ6/trunk/libf/phylmd/flott_gwd_rando_m.F90

Timestamp:

Feb 12, 2018, 1:24:03 AM (6 years ago)

Author:

fhourdin

Message:

Retour vers l'insensibilite au decoupage en sous domaine.
Les routines gwd_rando incluait le calcul de niveaux de reference
sur la base d'un profile pris au milieu du domaine (en klon/2).
Rempace par un test en presnivs.

Une autre intercation entre routines concernant la tke a fait apparaitre
que la tke n'était pas passee correctement au niveau klev+1 au moment
du regroupement des mailles sous les sous surface.

Ces changements garantissent la convergence numerique si
addtkeoro=0
iflag_pbl<12
et
ok_gwd_rando=n
La convergence n'est pas garantie pour les dernieres versions des physiq.def
mais les differences devraient etre mineures.

File:

• LMDZ6/trunk/libf/phylmd/flott_gwd_rando_m.F90 (modified) (8 diffs)

LMDZ6/trunk/libf/phylmd/flott_gwd_rando_m.F90

@@ -18 +18 @@
       use dimphy, only: klon, klev
       use assert_m, only: assert
+      USE ioipsl_getin_p_mod, ONLY : getin_p
+      USE vertical_layers_mod, ONLY : presnivs
+
       include "YOMCST.h"
       include "clesphys.h"
@@ -103 +106 @@
     REAL PH(KLON, KLEV + 1) ! Pressure at 1/2 levels
     REAL PSEC ! Security to avoid division by 0 pressure
-    REAL PHM1(KLON, KLEV + 1) ! 1/Press at 1/2 levels
     REAL BV(KLON, KLEV + 1) ! Brunt Vaisala freq. (BVF) at 1/2 levels
     REAL BVSEC ! Security to avoid negative BVF
+    REAL RAN_NUM_1,RAN_NUM_2,RAN_NUM_3
+
+    REAL, DIMENSION(klev+1) ::HREF
+
+    LOGICAL, SAVE :: gwd_reproductibilite_mpiomp=.true.
+    LOGICAL, SAVE :: firstcall = .TRUE.
+  !$OMP THREADPRIVATE(firstcall,gwd_reproductibilite_mpiomp)
+
+    CHARACTER (LEN=20) :: modname='flott_gwd_rando'
+    CHARACTER (LEN=80) :: abort_message
+
+
+
+  IF (firstcall) THEN
+    ! Cle introduite pour resoudre un probleme de non reproductibilite
+    ! Le but est de pouvoir tester de revenir a la version precedenete
+    ! A eliminer rapidement
+    CALL getin_p('gwd_reproductibilite_mpiomp',gwd_reproductibilite_mpiomp)
+    IF (NW+3*NA>=KLEV) THEN
+       abort_message = 'NW+3*NA>=KLEV Probleme pour generation des ondes'
+       CALL abort_physic (modname,abort_message,1)
+    ENDIF
+    firstcall=.false.
+  ENDIF
+
 
     !-----------------------------------------------------------------
@@ -156 +183 @@
     ZOISEC = 1.E-6 ! Security FOR 0 INTRINSIC FREQ
 
+IF (1==0) THEN
     !ONLINE
         call assert(klon == (/size(pp, 1), size(tt, 1), size(uu, 1), &
@@ -167 +195 @@
           "FLOTT_GWD_RANDO klev")
     !END ONLINE
+ENDIF
 
     IF(DELTAT < DTIME)THEN
@@ -183 +212 @@
     DO LL = 2, KLEV
        PH(:, LL) = EXP((LOG(PP(:, LL)) + LOG(PP(:, LL - 1))) / 2.)
-       PHM1(:, LL) = 1. / PH(:, LL) 
-    end DO
-
+    end DO
     PH(:, KLEV + 1) = 0. 
-    PHM1(:, KLEV + 1) = 1. / PSEC
     PH(:, 1) = 2. * PP(:, 1) - PH(:, 2)
 
     ! Launching altitude
+
+    !Pour revenir a la version non reproductible en changeant le nombre de process
+    IF (gwd_reproductibilite_mpiomp) THEN
+       ! Reprend la formule qui calcule PH en fonction de PP=play
+       DO LL = 2, KLEV
+          HREF(LL) = EXP((LOG(presnivs(LL)) + LOG(presnivs(LL - 1))) / 2.)
+       end DO
+       HREF(KLEV + 1) = 0.
+       HREF(1) = 2. * presnivs(1) - HREF(2)
+    ELSE
+       HREF(1:KLEV)=PH(KLON/2,1:KLEV)
+    ENDIF
 
     LAUNCH=0
     LTROP =0
     DO LL = 1, KLEV
-       IF (PH(KLON / 2, LL) / PH(KLON / 2, 1) > XLAUNCH) LAUNCH = LL
+       IF (HREF(LL) / HREF(1) > XLAUNCH) LAUNCH = LL
     ENDDO
     DO LL = 1, KLEV
-       IF (PH(KLON / 2, LL) / PH(KLON / 2, 1) > XTROP) LTROP = LL
+       IF (HREF(LL) / HREF(1) > XTROP) LTROP = LL
     ENDDO
+    !LAUNCH=22 ; LTROP=33
+!   print*,'LAUNCH=',LAUNCH,'LTROP=',LTROP
 
     ! Log pressure vert. coordinate
@@ -245 +285 @@
     ! waves characteristics in an almost stochastic way
 
-    JW = 0
-    DO JP = 1, NP
-       DO JK = 1, NK
-          DO JO = 1, NO
-             JW = JW + 1
+    DO JW = 1, NW
              ! Angle
              DO II = 1, KLON
                 ! Angle (0 or PI so far)
-                ZP(JW, II) = (SIGN(1., 0.5 - MOD(TT(II, JW) * 10., 1.)) + 1.) &
+                RAN_NUM_1=MOD(TT(II, JW) * 10., 1.)
+                RAN_NUM_2= MOD(TT(II, JW) * 100., 1.)
+                ZP(JW, II) = (SIGN(1., 0.5 - RAN_NUM_1) + 1.) &
                      * RPI / 2.
                 ! Horizontal wavenumber amplitude
-                ZK(JW, II) = KMIN + (KMAX - KMIN) * MOD(TT(II, JW) * 100., 1.)
+                ZK(JW, II) = KMIN + (KMAX - KMIN) *RAN_NUM_2
                 ! Horizontal phase speed
                 CPHA = 0.
                 DO JJ = 1, NA
+                    RAN_NUM_3=MOD(TT(II, JW+3*JJ)**2, 1.)
                     CPHA = CPHA + &
-         CMAX*2.*(MOD(TT(II, JW+3*JJ)**2, 1.)-0.5)*SQRT(3.)/SQRT(NA*1.)
+                    CMAX*2.*(RAN_NUM_3 -0.5)*SQRT(3.)/SQRT(NA*1.)
                 END DO
                 IF (CPHA.LT.0.)  THEN
@@ -276 +315 @@
                 RUW0(JW, II) = RUWMAX
              ENDDO
-          end DO
-       end DO
-    end DO
+    ENDDO
 
     ! 4. COMPUTE THE FLUXES
@@ -417 +454 @@
     ENDDO
 
+
   END SUBROUTINE FLOTT_GWD_RANDO