From 78c74d696de734b192e6ecba57d0596d852f3c24 Mon Sep 17 00:00:00 2001
From: Quentin Rodier <quentin.rodier@meteo.fr>
Date: Fri, 29 Apr 2022 23:03:39 +0200
Subject: [PATCH] Quentin 29/04/2022: Expand in physical points :
 turb_ver_sv_flux

---
 src/common/turb/mode_tridiag.F90          |  98 +++++++++--------
 src/common/turb/mode_turb_ver_sv_flux.F90 | 126 ++++++++++++----------
 2 files changed, 120 insertions(+), 104 deletions(-)

diff --git a/src/common/turb/mode_tridiag.F90 b/src/common/turb/mode_tridiag.F90
index eba6cbdb0..a3fc17662 100644
--- a/src/common/turb/mode_tridiag.F90
+++ b/src/common/turb/mode_tridiag.F90
@@ -139,7 +139,7 @@ REAL, DIMENSION(D%NIT,D%NJT,D%NKT)  :: ZY ,ZGAM
 REAL, DIMENSION(D%NIT,D%NJT)         :: ZBET
                                          ! 2D work array
 INTEGER                              :: JI,JJ,JK            ! loop counter
-INTEGER                              :: IKB,IKE       ! inner vertical limits
+INTEGER                              :: IKB,IKE,IIB,IIE,IJB,IJE       ! inner vertical limits
 INTEGER                              :: IKT           ! array size in k direction
 INTEGER                              :: IKTB,IKTE     ! start, end of k loops in physical domain 
 
@@ -157,27 +157,31 @@ IKTB=D%NKTB
 IKTE=D%NKTE
 IKB=D%NKB
 IKE=D%NKE
+IIE=D%NIEC
+IIB=D%NIBC
+IJE=D%NJEC
+IJB=D%NJBC
 !
-!$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-ZY(:,:,IKB) = PVARM(:,:,IKB)  + PTSTEP*PSOURCE(:,:,IKB) -   &
-  PEXPL / PRHODJ(:,:,IKB) * PA(:,:,IKB+D%NKL) * (PVARM(:,:,IKB+D%NKL) - PVARM(:,:,IKB))
-!$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+!$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+ZY(IIB:IIE,IJB:IJE,IKB) = PVARM(IIB:IIE,IJB:IJE,IKB)  + PTSTEP*PSOURCE(IIB:IIE,IJB:IJE,IKB) -   &
+  PEXPL / PRHODJ(IIB:IIE,IJB:IJE,IKB) * PA(IIB:IIE,IJB:IJE,IKB+D%NKL) * (PVARM(IIB:IIE,IJB:IJE,IKB+D%NKL) - PVARM(IIB:IIE,IJB:IJE,IKB))
+!$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
 !
 DO JK=IKTB+1,IKTE-1
-  !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-  ZY(:,:,JK)= PVARM(:,:,JK)  + PTSTEP*PSOURCE(:,:,JK) -               &
-      PEXPL / PRHODJ(:,:,JK) *                                           &
-                             ( PVARM(:,:,JK-D%NKL)*PA(:,:,JK)                &
-                              -PVARM(:,:,JK)*(PA(:,:,JK)+PA(:,:,JK+D%NKL))   &
-                              +PVARM(:,:,JK+D%NKL)*PA(:,:,JK+D%NKL)              &
+  !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+  ZY(IIB:IIE,IJB:IJE,JK)= PVARM(IIB:IIE,IJB:IJE,JK)  + PTSTEP*PSOURCE(IIB:IIE,IJB:IJE,JK) -               &
+      PEXPL / PRHODJ(IIB:IIE,IJB:IJE,JK) *                                           &
+                             ( PVARM(IIB:IIE,IJB:IJE,JK-D%NKL)*PA(IIB:IIE,IJB:IJE,JK)                &
+                              -PVARM(IIB:IIE,IJB:IJE,JK)*(PA(IIB:IIE,IJB:IJE,JK)+PA(IIB:IIE,IJB:IJE,JK+D%NKL))   &
+                              +PVARM(IIB:IIE,IJB:IJE,JK+D%NKL)*PA(IIB:IIE,IJB:IJE,JK+D%NKL)              &
                              ) 
-  !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+  !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
 END DO
 ! 
-!$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-ZY(:,:,IKE)= PVARM(:,:,IKE) + PTSTEP*PSOURCE(:,:,IKE) +               &
-  PEXPL / PRHODJ(:,:,IKE) * PA(:,:,IKE) * (PVARM(:,:,IKE)-PVARM(:,:,IKE-D%NKL))
-!$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+!$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+ZY(IIB:IIE,IJB:IJE,IKE)= PVARM(IIB:IIE,IJB:IJE,IKE) + PTSTEP*PSOURCE(IIB:IIE,IJB:IJE,IKE) +               &
+  PEXPL / PRHODJ(IIB:IIE,IJB:IJE,IKE) * PA(IIB:IIE,IJB:IJE,IKE) * (PVARM(IIB:IIE,IJB:IJE,IKE)-PVARM(IIB:IIE,IJB:IJE,IKE-D%NKL))
+!$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
 !
 !
 !*       2.  INVERSION OF THE TRIDIAGONAL SYSTEM
@@ -188,53 +192,53 @@ IF ( PIMPL > 1.E-10 ) THEN
   !
   !  going up
   !
-  !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-  ZBET(:,:) = 1. - PIMPL * PA(:,:,IKB+D%NKL) / PRHODJ(:,:,IKB)  ! bet = b(ikb)
-  PVARP(:,:,IKB) = ZY(:,:,IKB) / ZBET(:,:)
-  !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)               
+  !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+  ZBET(IIB:IIE,IJB:IJE) = 1. - PIMPL * PA(IIB:IIE,IJB:IJE,IKB+D%NKL) / PRHODJ(IIB:IIE,IJB:IJE,IKB)  ! bet = b(ikb)
+  PVARP(IIB:IIE,IJB:IJE,IKB) = ZY(IIB:IIE,IJB:IJE,IKB) / ZBET(IIB:IIE,IJB:IJE)
+  !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)               
   !
   DO JK = IKB+D%NKL,IKE-D%NKL,D%NKL
-    !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-    ZGAM(:,:,JK) = PIMPL * PA(:,:,JK) / PRHODJ(:,:,JK-D%NKL) / ZBET(:,:)  
+    !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+    ZGAM(IIB:IIE,IJB:IJE,JK) = PIMPL * PA(IIB:IIE,IJB:IJE,JK) / PRHODJ(IIB:IIE,IJB:IJE,JK-D%NKL) / ZBET(IIB:IIE,IJB:IJE)  
                                                     ! gam(k) = c(k-1) / bet
-    ZBET(:,:)    = 1. - PIMPL * (  PA(:,:,JK) * (1. + ZGAM(:,:,JK))  &
-                                 + PA(:,:,JK+D%NKL)                      &
-                                ) / PRHODJ(:,:,JK)  
+    ZBET(IIB:IIE,IJB:IJE)    = 1. - PIMPL * (  PA(IIB:IIE,IJB:IJE,JK) * (1. + ZGAM(IIB:IIE,IJB:IJE,JK))  &
+                                 + PA(IIB:IIE,IJB:IJE,JK+D%NKL)                      &
+                                ) / PRHODJ(IIB:IIE,IJB:IJE,JK)  
                                                     ! bet = b(k) - a(k)* gam(k)  
-    PVARP(:,:,JK)= ( ZY(:,:,JK) - PIMPL * PA(:,:,JK) / PRHODJ(:,:,JK) &
-                    * PVARP(:,:,JK-D%NKL)                                 &
-                   ) / ZBET(:,:)
+    PVARP(IIB:IIE,IJB:IJE,JK)= ( ZY(IIB:IIE,IJB:IJE,JK) - PIMPL * PA(IIB:IIE,IJB:IJE,JK) / PRHODJ(IIB:IIE,IJB:IJE,JK) &
+                    * PVARP(IIB:IIE,IJB:IJE,JK-D%NKL)                                 &
+                   ) / ZBET(IIB:IIE,IJB:IJE)
                                         ! res(k) = (y(k) -a(k)*res(k-1))/ bet 
-    !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+    !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
   END DO
-  !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+  !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
   ! special treatment for the last level
-  ZGAM(:,:,IKE) = PIMPL * PA(:,:,IKE) / PRHODJ(:,:,IKE-D%NKL) / ZBET(:,:) 
+  ZGAM(IIB:IIE,IJB:IJE,IKE) = PIMPL * PA(IIB:IIE,IJB:IJE,IKE) / PRHODJ(IIB:IIE,IJB:IJE,IKE-D%NKL) / ZBET(IIB:IIE,IJB:IJE) 
                                                     ! gam(k) = c(k-1) / bet
-  ZBET(:,:)    = 1. - PIMPL * (  PA(:,:,IKE) * (1. + ZGAM(:,:,IKE))  &
-                              ) / PRHODJ(:,:,IKE)  
+  ZBET(IIB:IIE,IJB:IJE)    = 1. - PIMPL * (  PA(IIB:IIE,IJB:IJE,IKE) * (1. + ZGAM(IIB:IIE,IJB:IJE,IKE))  &
+                              ) / PRHODJ(IIB:IIE,IJB:IJE,IKE)  
                                                     ! bet = b(k) - a(k)* gam(k)  
-  PVARP(:,:,IKE)= ( ZY(:,:,IKE) - PIMPL * PA(:,:,IKE) / PRHODJ(:,:,IKE) &
-                                * PVARP(:,:,IKE-D%NKL)                      &
-                 ) / ZBET(:,:)
+  PVARP(IIB:IIE,IJB:IJE,IKE)= ( ZY(IIB:IIE,IJB:IJE,IKE) - PIMPL * PA(IIB:IIE,IJB:IJE,IKE) / PRHODJ(IIB:IIE,IJB:IJE,IKE) &
+                                * PVARP(IIB:IIE,IJB:IJE,IKE-D%NKL)                      &
+                 ) / ZBET(IIB:IIE,IJB:IJE)
                                        ! res(k) = (y(k) -a(k)*res(k-1))/ bet 
   !
   !  going down
   !
-  !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+  !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
   DO JK = IKE-D%NKL,IKB,-1*D%NKL
-  !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-    PVARP(:,:,JK) = PVARP(:,:,JK) - ZGAM(:,:,JK+D%NKL) * PVARP(:,:,JK+D%NKL) 
-  !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+  !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+    PVARP(IIB:IIE,IJB:IJE,JK) = PVARP(IIB:IIE,IJB:IJE,JK) - ZGAM(IIB:IIE,IJB:IJE,JK+D%NKL) * PVARP(IIB:IIE,IJB:IJE,JK+D%NKL) 
+  !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
   END DO
 !
 ELSE
 ! 
-  !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+  !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
   DO JK=IKTB,IKTE
-    PVARP(:,:,JK) = ZY(:,:,JK)
+    PVARP(IIB:IIE,IJB:IJE,JK) = ZY(IIB:IIE,IJB:IJE,JK)
   END DO
-  !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+  !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
 !
 END IF 
 !
@@ -242,10 +246,10 @@ END IF
 !*       3.  FILL THE UPPER AND LOWER EXTERNAL VALUES
 !            ----------------------------------------
 !
-!$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-PVARP(:,:,D%NKA)=PVARP(:,:,IKB)
-PVARP(:,:,D%NKU)=PVARP(:,:,IKE)
-!$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+!$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+PVARP(IIB:IIE,IJB:IJE,D%NKA)=PVARP(IIB:IIE,IJB:IJE,IKB)
+PVARP(IIB:IIE,IJB:IJE,D%NKU)=PVARP(IIB:IIE,IJB:IJE,IKE)
+!$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
 !
 !-------------------------------------------------------------------------------
 !
diff --git a/src/common/turb/mode_turb_ver_sv_flux.F90 b/src/common/turb/mode_turb_ver_sv_flux.F90
index e960864b4..73c08e7d9 100644
--- a/src/common/turb/mode_turb_ver_sv_flux.F90
+++ b/src/common/turb/mode_turb_ver_sv_flux.F90
@@ -225,6 +225,7 @@ USE MODI_GRADIENT_U
 USE MODI_GRADIENT_V
 USE MODI_GRADIENT_W
 USE MODI_GRADIENT_M
+USE SHUMAN_PHY , ONLY : DZM_PHY, MZM_PHY, MZF_PHY
 USE MODI_SHUMAN , ONLY : DZM, MZM, MZF
 USE MODE_TRIDIAG, ONLY: TRIDIAG
 USE MODE_EMOIST, ONLY: EMOIST
@@ -292,10 +293,10 @@ REAL, DIMENSION(D%NIT,D%NJT,D%NKT)  ::  &
        ZFLXZ,  &   ! vertical flux of the treated variable
        ZSOURCE,  & ! source of evolution for the treated variable
        ZKEFF,    & ! effectif diffusion coeff = LT * SQRT( TKE )
-       ZWORK1,ZWORK2! working var. for shuman operators (array syntax)
-INTEGER             :: IKB,IKE      ! I index values for the Beginning and End
-                                    ! mass points of the domain in the 3 direct.
+       ZWORK1,ZWORK2,&
+       ZWORK3,ZWORK4! working var. for shuman operators (array syntax)
 INTEGER             :: IKT          ! array size in k direction
+INTEGER             :: IIE,IIB,IJE,IJB,IKB,IKE      ! index value for the mass points of the domain 
 INTEGER             :: IKTB,IKTE    ! start, end of k loops in physical domain 
 INTEGER             :: JSV          ! loop counters
 INTEGER             :: JI,JJ,JK           ! loop
@@ -318,14 +319,21 @@ IKT=D%NKT
 IKTB=D%NKTB          
 IKTE=D%NKTE
 IKB=D%NKB
-IKE=D%NKE             
+IKE=D%NKE
+IIE=D%NIE
+IIB=D%NIB
+IJE=D%NJE
+IJB=D%NJB            
 !
 IF (OHARAT) THEN
-  !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)
-  ZKEFF(:,:,:) =  PLM(:,:,:) * SQRT(PTKEM(:,:,:)) + 50.*MFMOIST(:,:,:)
-  !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)
+  !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
+  ZKEFF(IIB:IIE,IJB:IJE,IKB:IKE) =  PLM(IIB:IIE,IJB:IJE,IKB:IKE) * SQRT(PTKEM(IIB:IIE,IJB:IJE,IKB:IKE)) + 50.*MFMOIST(IIB:IIE,IJB:IJE,IKB:IKE)
+  !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
 ELSE
-  ZKEFF(:,:,:) = MZM(PLM(:,:,:) * SQRT(PTKEM(:,:,:)), D%NKA, D%NKU, D%NKL)
+  !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
+  ZWORK1(IIB:IIE,IJB:IJE,1:D%NKT) = PLM(IIB:IIE,IJB:IJE,1:D%NKT)*SQRT(PTKEM(IIB:IIE,IJB:IJE,1:D%NKT))
+  !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
+  CALL MZM_PHY(D,ZWORK1,ZKEFF)
 ENDIF
 !
 IF(OBLOWSNOW) THEN
@@ -339,23 +347,23 @@ ENDIF
 !*       8.   SOURCES OF PASSIVE SCALAR VARIABLES
 !             -----------------------------------
 !
-ZWORK1=MZM(PRHODJ, D%NKA, D%NKU, D%NKL)
+CALL MZM_PHY(D,PRHODJ,ZWORK1)
 DO JSV=1,KSV
 !
   IF (ONOMIXLG .AND. JSV >= KSV_LGBEG .AND. JSV<= KSV_LGEND) CYCLE
 !
 ! Preparation of the arguments for TRIDIAG 
     IF (OHARAT) THEN
-      !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)  
-      ZA(:,:,:) = -PTSTEP * ZKEFF(:,:,:) * ZWORK1(:,:,:) / PDZZ(:,:,:)**2
-      !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)
+      !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)  
+      ZA(IIB:IIE,IJB:IJE,IKB:IKE) = -PTSTEP * ZKEFF(IIB:IIE,IJB:IJE,IKB:IKE) * ZWORK1(IIB:IIE,IJB:IJE,IKB:IKE) / PDZZ(IIB:IIE,IJB:IJE,IKB:IKE)**2
+      !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
     ELSE
-      !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)
-      ZA(:,:,:) = -PTSTEP*ZCSV*PPSI_SV(:,:,:,JSV) *   &
-                   ZKEFF(:,:,:) * ZWORK1(:,:,:) / PDZZ(:,:,:)**2
-      !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)
+      !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
+      ZA(IIB:IIE,IJB:IJE,IKB:IKE) = -PTSTEP*ZCSV*PPSI_SV(IIB:IIE,IJB:IJE,IKB:IKE,JSV) *   &
+                   ZKEFF(IIB:IIE,IJB:IJE,IKB:IKE) * ZWORK1(IIB:IIE,IJB:IJE,IKB:IKE) / PDZZ(IIB:IIE,IJB:IJE,IKB:IKE)**2
+      !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
     ENDIF
-  ZSOURCE(:,:,:) = 0.
+  ZSOURCE(IIB:IIE,IJB:IJE,IKB:IKE) = 0.
 !
 ! Compute the sources for the JSVth scalar variable
 
@@ -364,70 +372,74 @@ DO JSV=1,KSV
 !* in 1DIM case, the part of energy released in horizontal flux
 ! is taken into account in the vertical part
   IF (HTURBDIM=='3DIM') THEN
-    !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-    ZSOURCE(:,:,IKB) = (PIMPL*PSFSVP(:,:,JSV) + PEXPL*PSFSVM(:,:,JSV)) / &
-                       PDZZ(:,:,IKB) * PDIRCOSZW(:,:)                    &
-                     * 0.5 * (1. + PRHODJ(:,:,D%NKA) / PRHODJ(:,:,IKB))
-    !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+    !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+    ZSOURCE(IIB:IIE,IJB:IJE,IKB) = (PIMPL*PSFSVP(IIB:IIE,IJB:IJE,JSV) + PEXPL*PSFSVM(IIB:IIE,IJB:IJE,JSV)) / &
+                       PDZZ(IIB:IIE,IJB:IJE,IKB) * PDIRCOSZW(IIB:IIE,IJB:IJE)                    &
+                     * 0.5 * (1. + PRHODJ(IIB:IIE,IJB:IJE,D%NKA) / PRHODJ(IIB:IIE,IJB:IJE,IKB))
+    !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
   ELSE
-    !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-    ZSOURCE(:,:,IKB) = (PIMPL*PSFSVP(:,:,JSV) + PEXPL*PSFSVM(:,:,JSV)) / &
-                       PDZZ(:,:,IKB) / PDIRCOSZW(:,:)                    &
-                     * 0.5 * (1. + PRHODJ(:,:,D%NKA) / PRHODJ(:,:,IKB))
-    !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+    !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+    ZSOURCE(IIB:IIE,IJB:IJE,IKB) = (PIMPL*PSFSVP(IIB:IIE,IJB:IJE,JSV) + PEXPL*PSFSVM(IIB:IIE,IJB:IJE,JSV)) / &
+                       PDZZ(IIB:IIE,IJB:IJE,IKB) / PDIRCOSZW(IIB:IIE,IJB:IJE)                    &
+                     * 0.5 * (1. + PRHODJ(IIB:IIE,IJB:IJE,D%NKA) / PRHODJ(IIB:IIE,IJB:IJE,IKB))
+    !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
   END IF
-  ZSOURCE(:,:,IKTB+1:IKTE-1) = 0.
-  ZSOURCE(:,:,IKE) = 0.
+  ZSOURCE(IIB:IIE,IJB:IJE,IKTB+1:IKTE-1) = 0.
+  ZSOURCE(IIB:IIE,IJB:IJE,IKE) = 0.
 !
 ! Obtention of the split JSV scalar variable at t+ deltat  
   CALL TRIDIAG(D,PSVM(:,:,:,JSV),ZA,PTSTEP,PEXPL,PIMPL,PRHODJ,ZSOURCE,ZRES)
 !
 !  Compute the equivalent tendency for the JSV scalar variable
-  !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)
-  PRSVS(:,:,:,JSV)= PRSVS(:,:,:,JSV)+    &
-                    PRHODJ(:,:,:)*(ZRES(:,:,:)-PSVM(:,:,:,JSV))/PTSTEP
-  !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)
+  !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
+  PRSVS(IIB:IIE,IJB:IJE,IKB:IKE,JSV)= PRSVS(IIB:IIE,IJB:IJE,IKB:IKE,JSV)+    &
+                    PRHODJ(IIB:IIE,IJB:IJE,IKB:IKE)*(ZRES(IIB:IIE,IJB:IJE,IKB:IKE)-PSVM(IIB:IIE,IJB:IJE,IKB:IKE,JSV))/PTSTEP
+  !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
 !
   IF ( (OTURB_FLX .AND. TPFILE%LOPENED) .OR. OLES_CALL ) THEN
     ! Diagnostic of the cartesian vertical flux
     !
-    ZWORK1 = MZM(PLM*SQRT(PTKEM), D%NKA, D%NKU, D%NKL) 
-    ZWORK2 = DZM(PIMPL*ZRES(:,:,:) + PEXPL*PSVM(:,:,:,JSV), D%NKA, D%NKU, D%NKL)
-    !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)
-    ZFLXZ(:,:,:) = -ZCSV * PPSI_SV(:,:,:,JSV) * ZWORK1(:,:,:) / PDZZ(:,:,:) * &
-                  ZWORK2(:,:,:)
-    !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT,JK=1:D%NKT)
+    !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
+    ZWORK1(IIB:IIE,IJB:IJE,1:D%NKT) = PLM(IIB:IIE,IJB:IJE,1:D%NKT)*SQRT(PTKEM(IIB:IIE,IJB:IJE,1:D%NKT))
+    ZWORK2(IIB:IIE,IJB:IJE,1:D%NKT) = PIMPL*ZRES(IIB:IIE,IJB:IJE,1:D%NKT) + PEXPL*PSVM(IIB:IIE,IJB:IJE,1:D%NKT,JSV)
+    !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
+    CALL MZM_PHY(D,ZWORK1,ZWORK3)
+    CALL DZM_PHY(D,ZWORK2,ZWORK4)
+    !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
+    ZFLXZ(IIB:IIE,IJB:IJE,IKB:IKE) = -ZCSV * PPSI_SV(IIB:IIE,IJB:IJE,IKB:IKE,JSV) * ZWORK3(IIB:IIE,IJB:IJE,IKB:IKE) / PDZZ(IIB:IIE,IJB:IJE,IKB:IKE) * &
+                  ZWORK4(IIB:IIE,IJB:IJE,IKB:IKE)
+    !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE,JK=1:D%NKT)
     ! surface flux
     !* in 3DIM case, a part of the flux goes vertically, and another goes horizontally
     ! (in presence of slopes)
     !* in 1DIM case, the part of energy released in horizontal flux
     ! is taken into account in the vertical part
     IF (HTURBDIM=='3DIM') THEN
-      !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-      ZFLXZ(:,:,IKB) = (PIMPL*PSFSVP(:,:,JSV) + PEXPL*PSFSVM(:,:,JSV))  &
-                       * PDIRCOSZW(:,:)  
-      !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+      !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+      ZFLXZ(IIB:IIE,IJB:IJE,IKB) = (PIMPL*PSFSVP(IIB:IIE,IJB:IJE,JSV) + PEXPL*PSFSVM(IIB:IIE,IJB:IJE,JSV))  &
+                       * PDIRCOSZW(IIB:IIE,IJB:IJE)  
+      !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
     ELSE
-      !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-      ZFLXZ(:,:,IKB) = (PIMPL*PSFSVP(:,:,JSV) + PEXPL*PSFSVM(:,:,JSV))  &
-                       / PDIRCOSZW(:,:)
-     !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+      !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+      ZFLXZ(IIB:IIE,IJB:IJE,IKB) = (PIMPL*PSFSVP(IIB:IIE,IJB:IJE,JSV) + PEXPL*PSFSVM(IIB:IIE,IJB:IJE,JSV))  &
+                       / PDIRCOSZW(IIB:IIE,IJB:IJE)
+     !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
     END IF
     ! extrapolates the flux under the ground so that the vertical average with
     ! the IKB flux gives the ground value
     !
-    !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-    ZFLXZ(:,:,D%NKA) = ZFLXZ(:,:,IKB)
-    !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+    !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+    ZFLXZ(IIB:IIE,IJB:IJE,D%NKA) = ZFLXZ(IIB:IIE,IJB:IJE,IKB)
+    !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
     DO JK=IKTB+1,IKTE-1
-      !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-      PWSV(:,:,JK,JSV)=0.5*(ZFLXZ(:,:,JK)+ZFLXZ(:,:,JK+D%NKL))
-      !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+      !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+      PWSV(IIB:IIE,IJB:IJE,JK,JSV)=0.5*(ZFLXZ(IIB:IIE,IJB:IJE,JK)+ZFLXZ(IIB:IIE,IJB:IJE,JK+D%NKL))
+      !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
     END DO
-    !$mnh_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
-    PWSV(:,:,IKB,JSV)=0.5*(ZFLXZ(:,:,IKB)+ZFLXZ(:,:,IKB+D%NKL))
-    PWSV(:,:,IKE,JSV)=PWSV(:,:,IKE-D%NKL,JSV)
-    !$mnh_end_expand_array(JI=1:D%NIT,JJ=1:D%NJT)
+    !$mnh_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
+    PWSV(IIB:IIE,IJB:IJE,IKB,JSV)=0.5*(ZFLXZ(IIB:IIE,IJB:IJE,IKB)+ZFLXZ(IIB:IIE,IJB:IJE,IKB+D%NKL))
+    PWSV(IIB:IIE,IJB:IJE,IKE,JSV)=PWSV(IIB:IIE,IJB:IJE,IKE-D%NKL,JSV)
+    !$mnh_end_expand_array(JI=IIB:IIE,JJ=IJB:IJE)
  END IF
   !
   IF (OTURB_FLX .AND. TPFILE%LOPENED) THEN
-- 
GitLab