From 44b99f7325077045e78a7a4caa2417647e737c21 Mon Sep 17 00:00:00 2001
From: Philippe WAUTELET <philippe.wautelet@aero.obs-mip.fr>
Date: Thu, 9 Jan 2020 15:20:12 +0100
Subject: [PATCH] Philippe 09/01/2020: remove unused dummy arguments for
 gradient_*.f90 functions

---
 src/MNH/eddyUV_flux_one_wayn.f90       |  6 +--
 src/MNH/eddyUV_fluxn.f90               |  8 ++--
 src/MNH/eddy_flux_one_wayn.f90         |  8 ++--
 src/MNH/eddy_fluxn.f90                 |  6 +--
 src/MNH/elec_fieldn.f90                | 13 +++----
 src/MNH/gradient_m.f90                 | 30 +++++---------
 src/MNH/gradient_u.f90                 | 25 ++++--------
 src/MNH/gradient_uv.f90                | 17 +++-----
 src/MNH/gradient_uw.f90                | 16 ++------
 src/MNH/gradient_v.f90                 | 24 +++---------
 src/MNH/gradient_vw.f90                | 16 ++------
 src/MNH/gradient_w.f90                 | 24 +++---------
 src/MNH/lesn.f90                       | 21 +++++-----
 src/MNH/prandtl.f90                    | 54 +++++++++++++-------------
 src/MNH/turb.f90                       | 10 ++---
 src/MNH/turb_cloud_index.f90           | 14 +++----
 src/MNH/turb_hor_dyn_corr.f90          | 14 +++----
 src/MNH/turb_hor_sv_corr.f90           | 20 +++++-----
 src/MNH/turb_hor_sv_flux.f90           |  8 ++--
 src/MNH/turb_hor_thermo_corr.f90       | 19 +++++----
 src/MNH/turb_hor_thermo_flux.f90       | 24 ++++++------
 src/MNH/turb_hor_uv.f90                | 11 +++---
 src/MNH/turb_hor_uw.f90                |  6 +--
 src/MNH/turb_hor_vw.f90                |  8 ++--
 src/MNH/turb_ver_dyn_flux.f90          | 24 ++++++------
 src/MNH/turb_ver_sv_flux.f90           |  2 +-
 src/MNH/turb_ver_thermo_flux.f90       |  4 +-
 src/MNH/write_lfifm1_for_diag.f90      | 46 +++++++++++-----------
 src/MNH/write_lfifm1_for_diag_supp.f90 | 24 ++++++------
 29 files changed, 211 insertions(+), 291 deletions(-)

diff --git a/src/MNH/eddyUV_flux_one_wayn.f90 b/src/MNH/eddyUV_flux_one_wayn.f90
index 1e8e7c513..61bd98222 100644
--- a/src/MNH/eddyUV_flux_one_wayn.f90
+++ b/src/MNH/eddyUV_flux_one_wayn.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 2011-2019 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2011-2020 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -96,7 +96,6 @@ INTEGER:: IDTRATIO_KMI_1 ! Ratio between the time step of the son and the model
 REAL, DIMENSION(:,:,:), ALLOCATABLE :: ZFLUX2 ! Work array=Dad interpolated flux field
                                               ! on the son grid
 REAL, DIMENSION(:,:,:), ALLOCATABLE :: ZDIV_UV! Work array=DIV of ZFLUX2
-INTEGER :: IKU
 !
 INTEGER :: IDIMX,IDIMY
 INTEGER :: IID, IRESP
@@ -105,7 +104,6 @@ INTEGER :: IID, IRESP
 !
 ! test of temporal synchronisation between the model 1 and the son KMI
 !
-IKU=SIZE(XZHAT)
 IDTRATIO_KMI_1=1
 DO JMI=2,KMI
    IDTRATIO_KMI_1=IDTRATIO_KMI_1*NDTRATIO(JMI)
@@ -138,7 +136,7 @@ IF (ISYNCHRO==1 .OR. IDTRATIO_KMI_1 == 1) THEN
    ZFLUX2(IIE,:,:)    = ZFLUX2(IIE-1,:,:) 
    ZFLUX2(IIE+1,:,:)  = ZFLUX2(IIE,:,:)
 
-   ZDIV_UV(:,:,:) = GX_U_M(1,IKU,1,ZFLUX2,XDXX,XDZZ,XDZX)
+   ZDIV_UV(:,:,:) = GX_U_M(ZFLUX2,XDXX,XDZZ,XDZX)
 
    ! Lateral boundary conditions
    ZDIV_UV(IIB,:,:)  =0.0
diff --git a/src/MNH/eddyUV_fluxn.f90 b/src/MNH/eddyUV_fluxn.f90
index fd153a867..b0110fc06 100644
--- a/src/MNH/eddyUV_fluxn.f90
+++ b/src/MNH/eddyUV_fluxn.f90
@@ -1,6 +1,6 @@
-!MNH_LIC Copyright 1994-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2004-2020 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
-!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
+!MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
 !-----------------------------------------------------------------
 !     #########################
@@ -211,7 +211,7 @@ ZCORIOZ(:,:,:)= SPREAD(XCORIOZ(:,:),3,IKU)
 ZLAT3D(:,:,:) = SPREAD(XLAT(:,:),3,IKU)
 ! 
 ! relative vorticity
-ZVOZ(:,:,:)=GX_V_UV(1,IKU,1,PVM,XDXX,XDZZ,XDZX)
+ZVOZ(:,:,:)=GX_V_UV(PVM,XDXX,XDZZ,XDZX)
 ZVOZ(:,:,2)=ZVOZ(:,:,3)
 ZVOZ(:,:,1)=ZVOZ(:,:,3)
 !
@@ -424,7 +424,7 @@ PVU_FLUX_M(:,:,:) = ZUV_FLUX(:,:,:)
 !             -----------------------------
 !
 ! Take the divergence of the momentum flux 
-ZDIV_UV(:,:,:) = GX_U_M(1,IKU,1,ZUV_FLUX,XDXX,XDZZ,XDZX)
+ZDIV_UV(:,:,:) = GX_U_M(ZUV_FLUX,XDXX,XDZZ,XDZX)
 
 ! Lateral boundary conditions
 ZDIV_UV(IIB,:,:)=0.0
diff --git a/src/MNH/eddy_flux_one_wayn.f90 b/src/MNH/eddy_flux_one_wayn.f90
index 31549eab9..0dfda6a59 100644
--- a/src/MNH/eddy_flux_one_wayn.f90
+++ b/src/MNH/eddy_flux_one_wayn.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 2011-2018 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 2011-2020 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -92,7 +92,6 @@ INTEGER:: IDTRATIO_KMI_1 ! Ratio between the time step of the son and the model
 
 REAL, DIMENSION(:,:,:), ALLOCATABLE :: ZFLUX2 ! Work array=Dad interpolated flux field
                                               ! on the son grid
-INTEGER :: IKU                                              
 !
 INTEGER :: IDIMX,IDIMY
 INTEGER :: IID, IRESP
@@ -101,7 +100,6 @@ INTEGER :: IID, IRESP
 !
 ! test of temporal synchronisation between the model 1 and the son KMI
 !
-IKU=SIZE(XZHAT)
 IDTRATIO_KMI_1=1
 DO JMI=2,KMI
    IDTRATIO_KMI_1=IDTRATIO_KMI_1*NDTRATIO(JMI)
@@ -126,7 +124,7 @@ IF (ISYNCHRO==1 .OR. IDTRATIO_KMI_1 == 1) THEN
                   HLBCX,HLBCY,TFIELDLIST(IID)%TFIELD_X3D(1)%DATA,ZFLUX2)
 
    ! operator GX_U_M used for gradient of v'T' (flux point) placed at a mass point
-   XRTHS_EDDY_FLUX(:,:,:) = - XRHODJ(:,:,:)* GX_U_M(1,IKU,1,ZFLUX2,XDXX,XDZZ,XDZX)
+   XRTHS_EDDY_FLUX(:,:,:) = - XRHODJ(:,:,:)* GX_U_M(ZFLUX2,XDXX,XDZZ,XDZX)
         
    ! w'T' (EDDY_FLUX_MODEL(1)%XWTH_FLUX_M) of model1 interpolation on the son grid put into ZFLUX2
    ZFLUX2 = 0.
@@ -139,7 +137,7 @@ IF (ISYNCHRO==1 .OR. IDTRATIO_KMI_1 == 1) THEN
 
    ! DIV(W'T') put into the source term
    XRTHS_EDDY_FLUX(:,:,:) = XRTHS_EDDY_FLUX(:,:,:) &
-                            - XRHODJ(:,:,:)* GZ_W_M(1,IKU,1,ZFLUX2,XDZZ)
+                            - XRHODJ(:,:,:)* GZ_W_M(ZFLUX2,XDZZ)
 
    DEALLOCATE(ZFLUX2)
 
diff --git a/src/MNH/eddy_fluxn.f90 b/src/MNH/eddy_fluxn.f90
index 9dab9b12c..484bfe0d2 100644
--- a/src/MNH/eddy_fluxn.f90
+++ b/src/MNH/eddy_fluxn.f90
@@ -247,7 +247,7 @@ ENDDO
 ZND(:,:,:) = MXM(ZND(:,:,:))
 !! latitudinal gradient of TH
 ZDTHM_DY(:,:,:) = GX_M_U(1,IKU,1,PTHM,XDXX,XDZZ,XDZX)
-ZDTHM_DZ(:,:,:) = MXM(GZ_M_M(1,IKU,1,PTHM,XDZZ))
+ZDTHM_DZ(:,:,:) = MXM(GZ_M_M(PTHM,XDZZ))
 ! density scale height
 ZH(:,:,:) = PTHM(:,:,:) * XRD * (XG**(-1))
 ZH(:,:,:) = MXM(ZH)
@@ -432,9 +432,9 @@ ENDIF
 !      --------------------      
 ! operator GX_U_M used for gradient of v'T' (flux point) placed at a mass point        
 !
-ZDIV_YTHFLUX(:,:,:) = GX_U_M(1,IKU,1,ZVTH_FLUX,XDXX,XDZZ,XDZX) 
+ZDIV_YTHFLUX(:,:,:) = GX_U_M(ZVTH_FLUX,XDXX,XDZZ,XDZX)
 !
-ZDIV_ZTHFLUX(:,:,:) = GZ_W_M(1,IKU,1,ZWTH_FLUX,XDZZ)
+ZDIV_ZTHFLUX(:,:,:) = GZ_W_M(ZWTH_FLUX,XDZZ)
 
 !
 ! Control test for the sign of the flux 
diff --git a/src/MNH/elec_fieldn.f90 b/src/MNH/elec_fieldn.f90
index ee38a7d92..e6da7d769 100644
--- a/src/MNH/elec_fieldn.f90
+++ b/src/MNH/elec_fieldn.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 1994-2019 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2020 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -99,7 +99,7 @@ INTEGER :: IRESP    ! Return code of FM routines
 INTEGER :: ILENG    ! Length of the data field in LFIFM file
 INTEGER :: IGRID    ! C-grid indicator in LFIFM file
 INTEGER :: ILENCH   ! Length of comment string in LFIFM file
-INTEGER :: IIU,IJU,IKU ! array sizes in I,J,K
+INTEGER :: IIU,IJU  ! array sizes in I,J,K
 INTEGER :: JI,JJ,JK    !loop index on the vertical levels
 INTEGER :: IINFO_ll  
 !
@@ -117,8 +117,7 @@ CALL GET_PHYSICAL_ll(IIB,IJB,IIE,IJE)
 CALL GET_DIM_EXT_ll('B',IIU,IJU)
 !
 IKB = 1 + JPVEXT
-IKU = SIZE(PEFIELDU,3)
-IKE = IKU - JPVEXT
+IKE = SIZE(PEFIELDU,3) - JPVEXT
 !
 ALLOCATE(ZDV_SOURCE(SIZE(PQ_SOURCE,1),SIZE(PQ_SOURCE,2),SIZE(PQ_SOURCE,3)))
 ALLOCATE(ZPHIT(SIZE(PQ_SOURCE,1),SIZE(PQ_SOURCE,2),SIZE(PQ_SOURCE,3)))
@@ -226,9 +225,9 @@ CALL ADD3DFIELD_ll( TZFIELDS_ll, ZPHIT, 'ELEC_FIELD_n::ZPHIT' )
 CALL UPDATE_HALO_ll(TZFIELDS_ll,IINFO_ll)
 CALL CLEANLIST_ll(TZFIELDS_ll)
 !
-PEFIELDU(:,:,:) = PRHODJ(:,:,:) * GX_M_M(1,IKU,1,ZPHIT,XDXX,XDZZ,XDZX)
-PEFIELDV(:,:,:) = PRHODJ(:,:,:) * GY_M_M(1,IKU,1,ZPHIT,XDYY,XDZZ,XDZY)
-PEFIELDW(:,:,:) = PRHODJ(:,:,:) * GZ_M_M(1,IKU,1,ZPHIT,XDZZ) 
+PEFIELDU(:,:,:) = PRHODJ(:,:,:) * GX_M_M(ZPHIT,XDXX,XDZZ,XDZX)
+PEFIELDV(:,:,:) = PRHODJ(:,:,:) * GY_M_M(ZPHIT,XDYY,XDZZ,XDZY)
+PEFIELDW(:,:,:) = PRHODJ(:,:,:) * GZ_M_M(ZPHIT,XDZZ)
 !
 IF (PRESENT(PPHIT)) PPHIT(:,:,:) = - PRHODJ(:,:,:) * ZPHIT(:,:,:) 
 !
diff --git a/src/MNH/gradient_m.f90 b/src/MNH/gradient_m.f90
index 8442b1acf..b5ec025aa 100644
--- a/src/MNH/gradient_m.f90
+++ b/src/MNH/gradient_m.f90
@@ -10,9 +10,8 @@
 INTERFACE
 !
 !
-FUNCTION GX_M_M(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_M_M)
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
+FUNCTION GX_M_M(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_M_M)
+!
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the mass point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -23,9 +22,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGX_M_M ! result mass point
 END FUNCTION GX_M_M
 !
 !
-FUNCTION GY_M_M(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_M_M)
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
+FUNCTION GY_M_M(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_M_M)
+!
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the mass point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -36,10 +34,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGY_M_M ! result mass point
 END FUNCTION GY_M_M
 !
 !
-FUNCTION GZ_M_M(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_M_M)
+FUNCTION GZ_M_M(PA,PDZZ)      RESULT(PGZ_M_M)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the mass point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
@@ -103,7 +99,7 @@ END MODULE MODI_GRADIENT_M
 !
 !
 !     #######################################################
-      FUNCTION GX_M_M(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_M_M)
+      FUNCTION GX_M_M(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_M_M)
 !     #######################################################
 !
 !!****  *GX_M_M* - Cartesian Gradient operator: 
@@ -169,8 +165,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,                 INTENT(IN)  :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,                 INTENT(IN)  :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the mass point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -202,7 +196,7 @@ END FUNCTION GX_M_M
 !
 !
 !     #######################################################
-      FUNCTION GY_M_M(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_M_M)
+      FUNCTION GY_M_M(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_M_M)
 !     #######################################################
 !
 !!****  *GY_M_M* - Cartesian Gradient operator: 
@@ -266,8 +260,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,                 INTENT(IN)  :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,                 INTENT(IN)  :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the mass point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -300,9 +292,9 @@ END FUNCTION GY_M_M
 !
 !
 !
-!     #######################################################
-      FUNCTION GZ_M_M(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_M_M)
-!     #######################################################
+!     #############################################
+      FUNCTION GZ_M_M(PA,PDZZ)      RESULT(PGZ_M_M)
+!     #############################################
 !
 !!****  *GZ_M_M* - Cartesian Gradient operator: 
 !!                          computes the gradient in the cartesian Z
@@ -360,8 +352,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the mass point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
diff --git a/src/MNH/gradient_u.f90 b/src/MNH/gradient_u.f90
index 8b9fa7ff3..3d32ffa80 100644
--- a/src/MNH/gradient_u.f90
+++ b/src/MNH/gradient_u.f90
@@ -10,9 +10,8 @@
 INTERFACE
 !
 !     
-FUNCTION GX_U_M(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_U_M)
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
+FUNCTION GX_U_M(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_U_M)
+!
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the U point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -23,10 +22,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGX_U_M ! result mass point
 END FUNCTION GX_U_M
 !
 !     
-FUNCTION GY_U_UV(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_U_UV)
+FUNCTION GY_U_UV(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_U_UV)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the U point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -37,10 +34,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGY_U_UV ! result UV point
 END FUNCTION GY_U_UV
 !
 !     
-FUNCTION GZ_U_UW(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_U_UW)
+FUNCTION GZ_U_UW(PA,PDZZ)      RESULT(PGZ_U_UW)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the U point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
@@ -56,7 +51,7 @@ END MODULE MODI_GRADIENT_U
 !
 !
 !     #######################################################
-      FUNCTION GX_U_M(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_U_M)
+      FUNCTION GX_U_M(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_U_M)
 !     #######################################################
 !
 !!****  *GX_U_M* - Cartesian Gradient operator: 
@@ -121,8 +116,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,                 INTENT(IN)  :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,                 INTENT(IN)  :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the U point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -154,7 +147,7 @@ END FUNCTION GX_U_M
 !
 ! 
 !     #########################################################
-      FUNCTION GY_U_UV(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_U_UV)
+      FUNCTION GY_U_UV(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_U_UV)
 !     #########################################################
 !
 !!****  *GY_U_UV* - Cartesian Gradient operator: 
@@ -220,8 +213,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,                 INTENT(IN)  :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,                 INTENT(IN)  :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the U point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -252,7 +243,7 @@ END FUNCTION GY_U_UV
 !
 !
 !     #######################################################
-      FUNCTION GZ_U_UW(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_U_UW)
+      FUNCTION GZ_U_UW(PA,PDZZ)      RESULT(PGZ_U_UW)
 !     #######################################################
 !
 !!****  *GZ_U_UW - Cartesian Gradient operator: 
@@ -310,8 +301,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the U point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
diff --git a/src/MNH/gradient_uv.f90 b/src/MNH/gradient_uv.f90
index c350c786f..8e1822ef3 100644
--- a/src/MNH/gradient_uv.f90
+++ b/src/MNH/gradient_uv.f90
@@ -10,9 +10,8 @@
 INTERFACE
 !
 !     
-FUNCTION GX_UV_V(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_UV_V)
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
+FUNCTION GX_UV_V(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_UV_V)
+!
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the UV point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -23,10 +22,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGX_UV_V ! result V point
 END FUNCTION GX_UV_V
 !
 !     
-FUNCTION GY_UV_U(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_UV_U)
+FUNCTION GY_UV_U(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_UV_U)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the UV point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -44,7 +41,7 @@ END MODULE MODI_GRADIENT_UV
 !
 !
 !     #########################################################
-      FUNCTION GX_UV_V(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_UV_V)
+      FUNCTION GX_UV_V(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_UV_V)
 !     #########################################################
 !
 !!****  *GX_UV_V* - Cartesian Gradient operator: 
@@ -107,8 +104,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the UV point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -140,7 +135,7 @@ END FUNCTION GX_UV_V
 !
 ! 
 !     #########################################################
-      FUNCTION GY_UV_U(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_UV_U)
+      FUNCTION GY_UV_U(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_UV_U)
 !     #########################################################
 !
 !!****  *GY_UV_U* - Cartesian Gradient operator: 
@@ -211,8 +206,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the UV point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
diff --git a/src/MNH/gradient_uw.f90 b/src/MNH/gradient_uw.f90
index d9ec4ec0e..946b85c83 100644
--- a/src/MNH/gradient_uw.f90
+++ b/src/MNH/gradient_uw.f90
@@ -10,9 +10,7 @@
 INTERFACE
 !
 !     
-FUNCTION GX_UW_W(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_UW_W)
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
+FUNCTION GX_UW_W(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_UW_W)
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the UW point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -23,10 +21,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGX_UW_W ! result W point
 END FUNCTION GX_UW_W
 !
 !     
-FUNCTION GZ_UW_U(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_UW_U)
+FUNCTION GZ_UW_U(PA,PDZZ)      RESULT(PGZ_UW_U)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the UW point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
@@ -42,7 +38,7 @@ END MODULE MODI_GRADIENT_UW
 !
 !
 !     #########################################################
-      FUNCTION GX_UW_W(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_UW_W)
+      FUNCTION GX_UW_W(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_UW_W)
 !     #########################################################
 !
 !!****  *GX_UW_W* - Cartesian Gradient operator: 
@@ -106,8 +102,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise*
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the UW point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -139,7 +133,7 @@ END FUNCTION GX_UW_W
 !
 ! 
 !     ###############################################
-      FUNCTION GZ_UW_U(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_UW_U)
+      FUNCTION GZ_UW_U(PA,PDZZ)      RESULT(PGZ_UW_U)
 !     ###############################################
 !
 !!****  *GZ_UW_U* - Cartesian Gradient operator: 
@@ -200,8 +194,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the UW point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
diff --git a/src/MNH/gradient_v.f90 b/src/MNH/gradient_v.f90
index 1e609b305..12c1be749 100644
--- a/src/MNH/gradient_v.f90
+++ b/src/MNH/gradient_v.f90
@@ -10,10 +10,8 @@
 INTERFACE
 !
 !           
-FUNCTION GY_V_M(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_V_M)
+FUNCTION GY_V_M(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_V_M)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the V point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -23,10 +21,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGY_V_M ! result mass point
 !
 END FUNCTION GY_V_M
 !           
-FUNCTION GX_V_UV(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_V_UV)
+FUNCTION GX_V_UV(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_V_UV)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the V point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -37,10 +33,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGX_V_UV ! result UV point
 END FUNCTION GX_V_UV
 !
 !           
-FUNCTION GZ_V_VW(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_V_VW)
+FUNCTION GZ_V_VW(PA,PDZZ)      RESULT(PGZ_V_VW)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the V point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
@@ -56,7 +50,7 @@ END MODULE MODI_GRADIENT_V
 !
 !
 !     #######################################################
-      FUNCTION GY_V_M(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_V_M)
+      FUNCTION GY_V_M(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_V_M)
 !     #######################################################
 !
 !!****  *GY_V_M* - Cartesian Gradient operator: 
@@ -120,8 +114,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,                 INTENT(IN)  :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,                 INTENT(IN)  :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the V point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -153,7 +145,7 @@ END FUNCTION GY_V_M
 !
 ! 
 !     #########################################################
-      FUNCTION GX_V_UV(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_V_UV)
+      FUNCTION GX_V_UV(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_V_UV)
 !     #########################################################
 !
 !!****  *GX_V_UV* - Cartesian Gradient operator: 
@@ -218,8 +210,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,                 INTENT(IN)  :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,                 INTENT(IN)  :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the V point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -250,7 +240,7 @@ END FUNCTION GX_V_UV
 !
 !
 !     #######################################################
-      FUNCTION GZ_V_VW(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_V_VW)
+      FUNCTION GZ_V_VW(PA,PDZZ)      RESULT(PGZ_V_VW)
 !     #######################################################
 !
 !!****  *GZ_V_VW - Cartesian Gradient operator: 
@@ -309,8 +299,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the V point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
diff --git a/src/MNH/gradient_vw.f90 b/src/MNH/gradient_vw.f90
index c9016aa32..864be89c8 100644
--- a/src/MNH/gradient_vw.f90
+++ b/src/MNH/gradient_vw.f90
@@ -10,9 +10,7 @@
 INTERFACE
 !
 !     
-FUNCTION GY_VW_W(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_VW_W)
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
+FUNCTION GY_VW_W(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_VW_W)
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the VW point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -23,10 +21,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGY_VW_W ! result W point
 END FUNCTION GY_VW_W
 !
 !     
-FUNCTION GZ_VW_V(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_VW_V)
+FUNCTION GZ_VW_V(PA,PDZZ)      RESULT(PGZ_VW_V)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the VW point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
@@ -42,7 +38,7 @@ END MODULE MODI_GRADIENT_VW
 !
 !
 !     #########################################################
-      FUNCTION GY_VW_W(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_VW_W)
+      FUNCTION GY_VW_W(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_VW_W)
 !     #########################################################
 !
 !!****  *GY_VW_W* - Cartesian Gradient operator: 
@@ -106,8 +102,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the VW point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -139,7 +133,7 @@ END FUNCTION GY_VW_W
 !
 ! 
 !     ###############################################
-      FUNCTION GZ_VW_V(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_VW_V)
+      FUNCTION GZ_VW_V(PA,PDZZ)      RESULT(PGZ_VW_V)
 !     ###############################################
 !
 !!****  *GZ_VW_V* - Cartesian Gradient operator: 
@@ -200,8 +194,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the VW point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
diff --git a/src/MNH/gradient_w.f90 b/src/MNH/gradient_w.f90
index 2547e5a17..1ef8f6916 100644
--- a/src/MNH/gradient_w.f90
+++ b/src/MNH/gradient_w.f90
@@ -10,10 +10,8 @@
 INTERFACE
 !
 !            
-FUNCTION GZ_W_M(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_W_M)
+FUNCTION GZ_W_M(PA,PDZZ)      RESULT(PGZ_W_M)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the W point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
@@ -21,10 +19,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGZ_W_M ! result mass point
 !
 END FUNCTION GZ_W_M
 !            
-FUNCTION GX_W_UW(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_W_UW)
+FUNCTION GX_W_UW(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_W_UW)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the W point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -35,10 +31,8 @@ REAL, DIMENSION(SIZE(PA,1),SIZE(PA,2),SIZE(PA,3)) :: PGX_W_UW ! result UW point
 END FUNCTION GX_W_UW
 !
 !            
-FUNCTION GY_W_VW(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_W_VW)
+FUNCTION GY_W_VW(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_W_VW)
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the W point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dyy
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -56,7 +50,7 @@ END MODULE MODI_GRADIENT_W
 !
 !
 !     #######################################################
-      FUNCTION GZ_W_M(KKA,KKU,KL,PA,PDZZ)      RESULT(PGZ_W_M)
+      FUNCTION GZ_W_M(PA,PDZZ)      RESULT(PGZ_W_M)
 !     #######################################################
 !
 !!****  *GZ_W_M* - Cartesian Gradient operator: 
@@ -109,8 +103,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,              INTENT(IN)     :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,              INTENT(IN)     :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the W point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
 !
@@ -134,7 +126,7 @@ END FUNCTION GZ_W_M
 !
 !
 !     #########################################################
-      FUNCTION GX_W_UW(KKA,KKU,KL,PA,PDXX,PDZZ,PDZX)      RESULT(PGX_W_UW)
+      FUNCTION GX_W_UW(PA,PDXX,PDZZ,PDZX)      RESULT(PGX_W_UW)
 !     #########################################################
 !
 !!****  *GX_W_UW* - Cartesian Gradient operator: 
@@ -189,8 +181,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,                 INTENT(IN)  :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,                 INTENT(IN)  :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the W point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDXX    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
@@ -222,7 +212,7 @@ END FUNCTION GX_W_UW
 !
 !
 !     #########################################################
-      FUNCTION GY_W_VW(KKA,KKU,KL,PA,PDYY,PDZZ,PDZY)      RESULT(PGY_W_VW)
+      FUNCTION GY_W_VW(PA,PDYY,PDZZ,PDZY)      RESULT(PGY_W_VW)
 !     #########################################################
 !
 !!****  *GY_W_VW* - Cartesian Gradient operator: 
@@ -277,8 +267,6 @@ IMPLICIT NONE
 !
 !*       0.1   declarations of arguments and result
 !
-INTEGER,                 INTENT(IN)  :: KKA, KKU ! near ground and uppest atmosphere array indexes
-INTEGER,                 INTENT(IN)  :: KL     ! +1 if grid goes from ground to atmosphere top, -1 otherwise
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PA      ! variable at the W point
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDYY    ! metric coefficient dxx
 REAL, DIMENSION(:,:,:),  INTENT(IN)  :: PDZZ    ! metric coefficient dzz
diff --git a/src/MNH/lesn.f90 b/src/MNH/lesn.f90
index dd9fa4260..aca09170f 100644
--- a/src/MNH/lesn.f90
+++ b/src/MNH/lesn.f90
@@ -199,7 +199,7 @@ REAL, DIMENSION(:),   ALLOCATABLE :: ZWORK   !
 !
 INTEGER :: IRR      ! moist variables counter
 INTEGER :: JSV      ! scalar variables counter
-INTEGER :: IIU, IJU,IKU ! array sizes
+INTEGER :: IIU, IJU ! array sizes
 INTEGER :: IKE,IKB
 INTEGER :: JI, JJ, JK   ! loop counters
 INTEGER :: IIU_ll, IJU_ll    ! total domain I size (fin)
@@ -224,8 +224,7 @@ IIU_ll = IIMAX_ll+JPHEXT
 IJU_ll = IJMAX_ll+JPHEXT
 IIA_ll=JPHEXT+1
 IJA_ll=JPHEXT+1
-IKU=SIZE(XVT,3)
-IKE=IKU-JPVEXT
+IKE=SIZE(XVT,3)-JPVEXT
 IKB=1+JPVEXT
 CALL GET_DIM_EXT_ll('B',IIU,IJU)
 !
@@ -511,25 +510,25 @@ CALL LES_VER_INT( XTHP, ZTP_LES )
 CALL LES_VER_INT( XTR, ZTR_LES )
 CALL LES_VER_INT( XDISS, ZDISS_LES )
 CALL LES_VER_INT( XLEM, ZLM_LES )
-CALL LES_VER_INT( GZ_M_M(1,IKU,1,XPABST,XDZZ), ZDPDZ_LES )
+CALL LES_VER_INT( GZ_M_M(XPABST,XDZZ), ZDPDZ_LES )
 !
 CALL LES_VER_INT( MXF(XUT)  ,ZU_LES  )
 CALL LES_VER_INT( MYF(XVT)  ,ZV_LES  )
 CALL LES_VER_INT( MZF(XWT)  ,ZW_LES  )
 CALL LES_VER_INT( MZF(ZMASSF) ,ZMF_LES)
 CALL LES_VER_INT(     XTHT  ,ZTH_LES )
-CALL LES_VER_INT( MXF(MZF(GZ_U_UW(1,IKU,1,XUT,XDZZ))), ZDUDZ_LES )
-CALL LES_VER_INT( MYF(MZF(GZ_V_VW(1,IKU,1,XVT,XDZZ))), ZDVDZ_LES )
-CALL LES_VER_INT( GZ_W_M(1,IKU,1,XWT,XDZZ), ZDWDZ_LES )
+CALL LES_VER_INT( MXF(MZF(GZ_U_UW(XUT,XDZZ))), ZDUDZ_LES )
+CALL LES_VER_INT( MYF(MZF(GZ_V_VW(XVT,XDZZ))), ZDVDZ_LES )
+CALL LES_VER_INT( GZ_W_M(XWT,XDZZ), ZDWDZ_LES )
 CALL LES_VER_INT( ZEXN, ZEXN_LES)  
 !
-CALL LES_VER_INT( GZ_M_M(1,IKU,1,XTHT,XDZZ), ZDTHDZ_LES )
+CALL LES_VER_INT( GZ_M_M(XTHT,XDZZ), ZDTHDZ_LES )
 !
 CALL LES_VER_INT(ZRHO, ZRHO_LES)
 !
 IF (LUSERV) CALL LES_VER_INT(ZTHV, ZTHV_LES)
 CALL LES_VER_INT(ZTHL, ZTHL_LES)
-CALL LES_VER_INT( GZ_M_M(1,IKU,1,ZTHL,XDZZ), ZDTHLDZ_LES )
+CALL LES_VER_INT( GZ_M_M(ZTHL,XDZZ), ZDTHLDZ_LES )
 !
 CALL LES_VER_INT(     XTKET ,ZTKE_LES)
 IRR = 0
@@ -537,7 +536,7 @@ IF (LUSERV) THEN
   IRR = IRR + 1
   CALL LES_VER_INT(     XRT(:,:,:,IRR)  ,ZRV_LES )
   CALL LES_VER_INT(     ZRT(:,:,:)      ,ZRT_LES )
-  CALL LES_VER_INT( GZ_M_M(1,IKU,1,ZRT,XDZZ), ZDRTDZ_LES )
+  CALL LES_VER_INT( GZ_M_M(ZRT,XDZZ), ZDRTDZ_LES )
   CALL LES_VER_INT(   ZREHU(:,:,:)     ,ZREHU_LES)
 END IF
 IF (LUSERC) THEN
@@ -636,7 +635,7 @@ END IF
 IF (NSV>0) THEN
   DO JSV=1,NSV
     CALL LES_VER_INT(  XSVT(:,:,:,JSV), ZSV_LES(:,:,:,JSV) )
-    CALL LES_VER_INT( GZ_M_M(1,IKU,1,XSVT(:,:,:,JSV),XDZZ), ZDSVDZ_LES(:,:,:,JSV) )
+    CALL LES_VER_INT( GZ_M_M(XSVT(:,:,:,JSV),XDZZ), ZDSVDZ_LES(:,:,:,JSV) )
   END DO
 END IF
 !
diff --git a/src/MNH/prandtl.f90 b/src/MNH/prandtl.f90
index 283f2e02e..e5ea7ee8b 100644
--- a/src/MNH/prandtl.f90
+++ b/src/MNH/prandtl.f90
@@ -376,22 +376,22 @@ ELSE IF (L2D) THEN                      ! 3D case in a 2D model
 !
   IF (KRR /= 0) THEN                 ! moist 3D case
     PRED2TH3(:,:,:)= PREDTH1(:,:,:)**2+(XCTV*PBLL_O_E(:,:,:)*PETHETA(:,:,:) )**2 * &
-      MZM( GX_M_M(KKA,KKU,KKL,PTHLM,PDXX,PDZZ,PDZX)**2 )
+      MZM( GX_M_M(PTHLM,PDXX,PDZZ,PDZX)**2 )
     PRED2TH3(:,:,IKB)=PRED2TH3(:,:,IKB+KKL)
 !
     PRED2R3(:,:,:)= PREDR1(:,:,:)**2 + (XCTV*PBLL_O_E(:,:,:)*PEMOIST(:,:,:))**2 * &
-        MZM( GX_M_M(KKA,KKU,KKL,PRM(:,:,:,1),PDXX,PDZZ,PDZX)**2 )
+        MZM( GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)**2 )
     PRED2R3(:,:,IKB)=PRED2R3(:,:,IKB+KKL)
 !
     PRED2THR3(:,:,:)= PREDR1(:,:,:) * PREDTH1(:,:,:) +  XCTV**2*PBLL_O_E(:,:,:)**2 *   &
                   PEMOIST(:,:,:) * PETHETA(:,:,:) *                         &
-      MZM( GX_M_M(KKA,KKU,KKL,PRM(:,:,:,1),PDXX,PDZZ,PDZX)*     &
-                     GX_M_M(KKA,KKU,KKL,PTHLM,PDXX,PDZZ,PDZX))
+      MZM( GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)*     &
+                     GX_M_M(PTHLM,PDXX,PDZZ,PDZX))
     PRED2THR3(:,:,IKB)=PRED2THR3(:,:,IKB+KKL)
 !
   ELSE                 ! dry 3D case in a 2D model
     PRED2TH3(:,:,:) = PREDTH1(:,:,:)**2 +  XCTV**2*PBLL_O_E(:,:,:)**2 *     &
-      MZM( GX_M_M(KKA,KKU,KKL,PTHLM,PDXX,PDZZ,PDZX)**2 )
+      MZM( GX_M_M(PTHLM,PDXX,PDZZ,PDZX)**2 )
     PRED2TH3(:,:,IKB)=PRED2TH3(:,:,IKB+KKL)
 !
     PRED2R3(:,:,:) = 0.
@@ -404,27 +404,27 @@ ELSE                                 ! 3D case in a 3D model
 !
   IF (KRR /= 0) THEN                 ! moist 3D case
     PRED2TH3(:,:,:)= PREDTH1(:,:,:)**2 +  ( XCTV*PBLL_O_E(:,:,:)*PETHETA(:,:,:) )**2 * &
-      MZM( GX_M_M(KKA,KKU,KKL,PTHLM,PDXX,PDZZ,PDZX)**2 &
-      + GY_M_M(KKA,KKU,KKL,PTHLM,PDYY,PDZZ,PDZY)**2 )
+      MZM( GX_M_M(PTHLM,PDXX,PDZZ,PDZX)**2 &
+      + GY_M_M(PTHLM,PDYY,PDZZ,PDZY)**2 )
     PRED2TH3(:,:,IKB)=PRED2TH3(:,:,IKB+KKL)
 !
     PRED2R3(:,:,:)= PREDR1(:,:,:)**2 + (XCTV*PBLL_O_E(:,:,:)*PEMOIST(:,:,:))**2 *      &
-        MZM( GX_M_M(KKA,KKU,KKL,PRM(:,:,:,1),PDXX,PDZZ,PDZX)**2 + &
-        GY_M_M(KKA,KKU,KKL,PRM(:,:,:,1),PDYY,PDZZ,PDZY)**2 )
+        MZM( GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)**2 + &
+        GY_M_M(PRM(:,:,:,1),PDYY,PDZZ,PDZY)**2 )
     PRED2R3(:,:,IKB)=PRED2R3(:,:,IKB+KKL)
 !
     PRED2THR3(:,:,:)= PREDR1(:,:,:) * PREDTH1(:,:,:) +  XCTV**2*PBLL_O_E(:,:,:)**2 *   &
          PEMOIST(:,:,:) * PETHETA(:,:,:) *                            &
-         MZM( GX_M_M(KKA,KKU,KKL,PRM(:,:,:,1),PDXX,PDZZ,PDZX)*   &
-         GX_M_M(KKA,KKU,KKL,PTHLM,PDXX,PDZZ,PDZX)+                           &
-         GY_M_M(KKA,KKU,KKL,PRM(:,:,:,1),PDYY,PDZZ,PDZY)*                    &
-         GY_M_M(KKA,KKU,KKL,PTHLM,PDYY,PDZZ,PDZY) )
+         MZM( GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)*   &
+         GX_M_M(PTHLM,PDXX,PDZZ,PDZX)+                           &
+         GY_M_M(PRM(:,:,:,1),PDYY,PDZZ,PDZY)*                    &
+         GY_M_M(PTHLM,PDYY,PDZZ,PDZY) )
     PRED2THR3(:,:,IKB)=PRED2THR3(:,:,IKB+KKL)
 !
   ELSE                 ! dry 3D case in a 3D model
     PRED2TH3(:,:,:) = PREDTH1(:,:,:)**2 +  XCTV**2*PBLL_O_E(:,:,:)**2 *                &
-      MZM( GX_M_M(KKA,KKU,KKL,PTHLM,PDXX,PDZZ,PDZX)**2 &
-      + GY_M_M(KKA,KKU,KKL,PTHLM,PDYY,PDZZ,PDZY)**2 )
+      MZM( GX_M_M(PTHLM,PDXX,PDZZ,PDZX)**2 &
+      + GY_M_M(PTHLM,PDYY,PDZZ,PDZY)**2 )
     PRED2TH3(:,:,IKB)=PRED2TH3(:,:,IKB+KKL)
 !
     PRED2R3(:,:,:) = 0.
@@ -461,15 +461,15 @@ ELSE  IF (L2D) THEN ! 3D case in a 2D model
     END IF
     PRED2THS3(:,:,:,JSV) = PREDTH1(:,:,:) * PREDS1(:,:,:,JSV)   +        &
                        ZW1*                                              &
-                       MZM(GX_M_M(KKA,KKU,KKL,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*       &
-                           GX_M_M(KKA,KKU,KKL,PTHLM,PDXX,PDZZ,PDZX)                  &
+                       MZM(GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*       &
+                           GX_M_M(PTHLM,PDXX,PDZZ,PDZX)                  &
                           )
 !
     IF (KRR /= 0) THEN
       PRED2RS3(:,:,:,JSV) = PREDR1(:,:,:) * PREDS1(:,:,:,JSV)   +        &
                        ZW1 * PEMOIST *                                   &
-                       MZM(GX_M_M(KKA,KKU,KKL,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*       &
-                           GX_M_M(KKA,KKU,KKL,PRM(:,:,:,1),PDXX,PDZZ,PDZX)           &
+                       MZM(GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*       &
+                           GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)           &
                           )
     ELSE
       PRED2RS3(:,:,:,JSV) = 0.
@@ -486,19 +486,19 @@ ELSE ! 3D case in a 3D model
     END IF
     PRED2THS3(:,:,:,JSV) = PREDTH1(:,:,:) * PREDS1(:,:,:,JSV)   +        &
                        ZW1*                                              &
-                       MZM(GX_M_M(KKA,KKU,KKL,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*       &
-                           GX_M_M(KKA,KKU,KKL,PTHLM,PDXX,PDZZ,PDZX)                  &
-                          +GY_M_M(KKA,KKU,KKL,PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)*       &
-                           GY_M_M(KKA,KKU,KKL,PTHLM,PDYY,PDZZ,PDZY)                  &
+                       MZM(GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*       &
+                           GX_M_M(PTHLM,PDXX,PDZZ,PDZX)                  &
+                          +GY_M_M(PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)*       &
+                           GY_M_M(PTHLM,PDYY,PDZZ,PDZY)                  &
                           )
 !
     IF (KRR /= 0) THEN
       PRED2RS3(:,:,:,JSV) = PREDR1(:,:,:) * PREDS1(:,:,:,JSV)   +        &
                        ZW1 * PEMOIST *                                   &
-                       MZM(GX_M_M(KKA,KKU,KKL,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*       &
-                           GX_M_M(KKA,KKU,KKL,PRM(:,:,:,1),PDXX,PDZZ,PDZX)           &
-                          +GY_M_M(KKA,KKU,KKL,PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)*       &
-                           GY_M_M(KKA,KKU,KKL,PRM(:,:,:,1),PDYY,PDZZ,PDZY)           &
+                       MZM(GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*       &
+                           GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)           &
+                          +GY_M_M(PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)*       &
+                           GY_M_M(PRM(:,:,:,1),PDYY,PDZZ,PDZY)           &
                           )
     ELSE
       PRED2RS3(:,:,:,JSV) = 0.
diff --git a/src/MNH/turb.f90 b/src/MNH/turb.f90
index 0444bb066..512065c18 100644
--- a/src/MNH/turb.f90
+++ b/src/MNH/turb.f90
@@ -747,8 +747,8 @@ SELECT CASE (HTURBLEN)
 !           ------------------
 
   CASE ('RM17')
-    ZDUDZ = MXF(MZF(GZ_U_UW(1,KKU,1,PUT,PDZZ)))
-    ZDVDZ = MYF(MZF(GZ_V_VW(1,KKU,1,PVT,PDZZ)))
+    ZDUDZ = MXF(MZF(GZ_U_UW(PUT,PDZZ)))
+    ZDVDZ = MYF(MZF(GZ_V_VW(PVT,PDZZ)))
     ZSHEAR = SQRT(ZDUDZ*ZDUDZ + ZDVDZ*ZDVDZ)
     CALL BL89(KKA,KKU,KKL,PZZ,PDZZ,PTHVREF,ZTHLM,KRR,ZRM,PTKET,ZSHEAR,PLEM)
 !
@@ -875,9 +875,9 @@ IF (HTOM=='TM06') THEN
 !
   ZFWTH = -GZ_M_W(KKA,KKU,KKL,ZMWTH,PDZZ)    ! -d(w'2th' )/dz
   !ZFWR  = -GZ_M_W(KKA,KKU,KKL,ZMWR, PDZZ)    ! -d(w'2r'  )/dz
-  ZFTH2 = -GZ_W_M(KKA,KKU,KKL,ZMTH2,PDZZ)    ! -d(w'th'2 )/dz
-  !ZFR2  = -GZ_W_M(KKA,KKU,KKL,ZMR2, PDZZ)    ! -d(w'r'2  )/dz
-  !ZFTHR = -GZ_W_M(KKA,KKU,KKL,ZMTHR,PDZZ)    ! -d(w'th'r')/dz
+  ZFTH2 = -GZ_W_M(ZMTH2,PDZZ)    ! -d(w'th'2 )/dz
+  !ZFR2  = -GZ_W_M(ZMR2, PDZZ)    ! -d(w'r'2  )/dz
+  !ZFTHR = -GZ_W_M(ZMTHR,PDZZ)    ! -d(w'th'r')/dz
 !
   ZFWTH(:,:,IKTE:) = 0.
   ZFWTH(:,:,:IKTB) = 0.
diff --git a/src/MNH/turb_cloud_index.f90 b/src/MNH/turb_cloud_index.f90
index abe8e0f24..811d8ee1d 100644
--- a/src/MNH/turb_cloud_index.f90
+++ b/src/MNH/turb_cloud_index.f90
@@ -1,4 +1,4 @@
-!MNH_LIC Copyright 1994-2019 CNRS, Meteo-France and Universite Paul Sabatier
+!MNH_LIC Copyright 1994-2020 CNRS, Meteo-France and Universite Paul Sabatier
 !MNH_LIC This is part of the Meso-NH software governed by the CeCILL-C licence
 !MNH_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
 !MNH_LIC for details. version 1.
@@ -135,7 +135,6 @@ REAL, DIMENSION(SIZE(PRM,1),SIZE(PRM,2),SIZE(PRM,3),2) :: ZG_RVCI,ZQ_RVCI
 INTEGER             :: JI,JJ,JK     ! loop counters
 INTEGER             :: IIB,IJB,IKB  ! Begin of physical dimensions
 INTEGER             :: IIE,IJE,IKE  ! End   of physical dimensions
-INTEGER             :: IKU          ! array size in k
 INTEGER, DIMENSION(SIZE(PRM,1),SIZE(PRM,2),SIZE(PRM,3)) :: IMASK_CLOUD
                              ! 0 except cloudy points or adjacent points (1)
 TYPE(TFIELDDATA)    :: TZFIELD
@@ -147,8 +146,7 @@ TYPE(TFIELDDATA)    :: TZFIELD
 !
 CALL GET_INDICE_ll (IIB,IJB,IIE,IJE)
 IKB = 1   + JPVEXT
-IKU = SIZE(PRM,3)
-IKE = IKU - JPVEXT
+IKE = SIZE(PRM,3) - JPVEXT
 !
 IMASK_CLOUD(:,:,:) = 0
 PCEI(:,:,:) = 0.
@@ -167,8 +165,8 @@ ZRVCI0(:,:,:) = MAX ( PRRS(:,:,:,1) , 0. ) + MAX ( PRRS(:,:,:,2) , 0. )
 IF (KRRI>=1) ZRVCI0(:,:,:) = ZRVCI0(:,:,:) + MAX ( PRRS(:,:,:,4) , 0. )
 !
 ZRVCI(:,:,:)= PTSTEP *ZRVCI0(:,:,:) /PRHODJ(:,:,:)
-ZG_RVCI(:,:,:,1) = GX_M_M(1,IKU,1,ZRVCI,PDXX,PDZZ,PDZX)
-ZG_RVCI(:,:,:,2) = GY_M_M(1,IKU,1,ZRVCI,PDYY,PDZZ,PDZY)
+ZG_RVCI(:,:,:,1) = GX_M_M(ZRVCI,PDXX,PDZZ,PDZX)
+ZG_RVCI(:,:,:,2) = GY_M_M(ZRVCI,PDYY,PDZZ,PDZY)
 !
 ZGNORM_RVCI(:,:,:) = SQRT( ZG_RVCI(:,:,:,1)*ZG_RVCI(:,:,:,1) +        &
                            ZG_RVCI(:,:,:,2)*ZG_RVCI(:,:,:,2) )
@@ -185,8 +183,8 @@ ZWORK(:,:,:) = ZRVCI0 / PRHODJ(:,:,:) -   &
                ( PRM(:,:,:,1)+ PRM(:,:,:,2) ) / PTSTEP
 IF (KRRI>=1) ZWORK(:,:,:) = ZWORK(:,:,:) - PRM(:,:,:,4) / PTSTEP
 !
-ZQ_RVCI(:,:,:,1) = GX_M_M(1,IKU,1,ZWORK,PDXX,PDZZ,PDZX)
-ZQ_RVCI(:,:,:,2) = GY_M_M(1,IKU,1,ZWORK,PDYY,PDZZ,PDZY)
+ZQ_RVCI(:,:,:,1) = GX_M_M(ZWORK,PDXX,PDZZ,PDZX)
+ZQ_RVCI(:,:,:,2) = GY_M_M(ZWORK,PDYY,PDZZ,PDZY)
 !
 ZQNORM_RVCI(:,:,:) = SQRT( ZQ_RVCI(:,:,:,1)*ZQ_RVCI(:,:,:,1) +         &
                            ZQ_RVCI(:,:,:,2)*ZQ_RVCI(:,:,:,2) )
diff --git a/src/MNH/turb_hor_dyn_corr.f90 b/src/MNH/turb_hor_dyn_corr.f90
index f897784a3..1b81d1ad3 100644
--- a/src/MNH/turb_hor_dyn_corr.f90
+++ b/src/MNH/turb_hor_dyn_corr.f90
@@ -274,9 +274,9 @@ IKU = SIZE(PUM,3)
 !
 ZDIRSINZW(:,:) = SQRT( 1. - PDIRCOSZW(:,:)**2 )
 !
-GX_U_M_PUM = GX_U_M(1,IKU,1,PUM,PDXX,PDZZ,PDZX)
-IF (.NOT. L2D) GY_V_M_PVM = GY_V_M(1,IKU,1,PVM,PDYY,PDZZ,PDZY)
-GZ_W_M_PWM = GZ_W_M(1,IKU,1,PWM,PDZZ)
+GX_U_M_PUM = GX_U_M(PUM,PDXX,PDZZ,PDZX)
+IF (.NOT. L2D) GY_V_M_PVM = GY_V_M(PVM,PDYY,PDZZ,PDZY)
+GZ_W_M_PWM = GZ_W_M(PWM,PDZZ)
 !
 ZMZF_DZZ = MZF(PDZZ)
 !
@@ -582,7 +582,7 @@ PRWS = PRWS(:,:,:) + MZM(PRHODJ(:,:,:))*(ZWP(:,:,:)-PWM(:,:,:))/PTSTEP
 !
 !* recomputes flux using guess of W
 !
-GZ_W_M_ZWP = GZ_W_M(1,IKU,1,ZWP,PDZZ)
+GZ_W_M_ZWP = GZ_W_M(ZWP,PDZZ)
 ZFLX(:,:,IKB+1:)=ZFLX(:,:,IKB+1:) &
   - XCMFS * PK(:,:,IKB+1:) * (4./3.) * (GZ_W_M_ZWP(:,:,IKB+1:) - GZ_W_M_PWM(:,:,IKB+1:))
 !
@@ -605,16 +605,16 @@ IF (LLES_CALL .AND. KSPLT==1) THEN
   CALL SECOND_MNH(ZTIME1)
   CALL LES_MEAN_SUBGRID( ZFLX, X_LES_SUBGRID_W2 ) 
   CALL LES_MEAN_SUBGRID( -ZWORK, X_LES_RES_ddxa_W_SBG_UaW , .TRUE.)
-  CALL LES_MEAN_SUBGRID( GZ_M_M(1,IKU,1,PTHLM,PDZZ)*ZFLX, X_LES_RES_ddxa_Thl_SBG_UaW , .TRUE.)
+  CALL LES_MEAN_SUBGRID( GZ_M_M(PTHLM,PDZZ)*ZFLX, X_LES_RES_ddxa_Thl_SBG_UaW , .TRUE.)
   CALL LES_MEAN_SUBGRID(ZFLX*MZF(GZ_M_W(1,IKU,1,PTHLM,PDZZ)),X_LES_RES_ddz_Thl_SBG_W2)
   IF (KRR>=1) THEN
-    CALL LES_MEAN_SUBGRID( GZ_M_M(1,IKU,1,PRM(:,:,:,1),PDZZ)*ZFLX, &
+    CALL LES_MEAN_SUBGRID( GZ_M_M(PRM(:,:,:,1),PDZZ)*ZFLX, &
                            X_LES_RES_ddxa_Rt_SBG_UaW , .TRUE.)
     CALL LES_MEAN_SUBGRID(ZFLX*MZF(GZ_M_W(1,IKU,1,PRM(:,:,:,1),PDZZ)), &
                            X_LES_RES_ddz_Rt_SBG_W2)
   END IF
   DO JSV=1,NSV
-    CALL LES_MEAN_SUBGRID( GZ_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDZZ)*ZFLX, &
+    CALL LES_MEAN_SUBGRID( GZ_M_M(PSVM(:,:,:,JSV),PDZZ)*ZFLX, &
                            X_LES_RES_ddxa_Sv_SBG_UaW(:,:,:,JSV) , .TRUE.)
     CALL LES_MEAN_SUBGRID(ZFLX*MZF(GZ_M_W(1,IKU,1,PSVM(:,:,:,JSV),PDZZ)), &
                            X_LES_RES_ddz_Sv_SBG_W2(:,:,:,JSV))
diff --git a/src/MNH/turb_hor_sv_corr.f90 b/src/MNH/turb_hor_sv_corr.f90
index bc3b86878..f9e2c7b55 100644
--- a/src/MNH/turb_hor_sv_corr.f90
+++ b/src/MNH/turb_hor_sv_corr.f90
@@ -137,7 +137,6 @@ REAL, DIMENSION(SIZE(PSVM,1),SIZE(PSVM,2),SIZE(PSVM,3))       &
                                      :: ZFLX, ZA
 !
 INTEGER             :: JSV          ! loop counter
-INTEGER             :: IKU
 !
 REAL :: ZTIME1, ZTIME2
 !
@@ -148,7 +147,6 @@ REAL :: ZCQSVD = 2.4  ! constant for humidity - scalar covariance dissipation
 REAL :: ZCSV          !constant for the scalar flux 
 ! ---------------------------------------------------------------------------
 !
-IKU=SIZE(PTKEM,3)
 CALL SECOND_MNH(ZTIME1)
 !
 IF(LBLOWSNOW) THEN
@@ -167,11 +165,11 @@ DO JSV=1,NSV
   IF (LLES_CALL) THEN
     IF (.NOT. L2D) THEN
       ZFLX(:,:,:) =  ZCSV / ZCSVD * PLM(:,:,:) * PLEPS(:,:,:) *   &
-         (  GX_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)**2             &
-          + GY_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)**2 )
+         (  GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)**2             &
+          + GY_M_M(PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)**2 )
     ELSE
       ZFLX(:,:,:) =  ZCSV / ZCSVD * PLM(:,:,:) * PLEPS(:,:,:) *   &
-            GX_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)**2
+            GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)**2
     END IF
     CALL LES_MEAN_SUBGRID( -2.*ZCSVD*SQRT(PTKEM)*ZFLX/PLEPS, &
                            X_LES_SUBGRID_DISS_Sv2(:,:,:,JSV), .TRUE. )
@@ -184,12 +182,12 @@ DO JSV=1,NSV
     ZA(:,:,:)   =  ETHETA(KRR,KRRI,PTHLM,PRM,PLOCPEXNM,PATHETA,PSRCM)
     IF (.NOT. L2D) THEN
       ZFLX(:,:,:)=  PLM(:,:,:) * PLEPS(:,:,:)                                          &
-          *  (  GX_M_M(1,IKU,1,PTHLM,PDXX,PDZZ,PDZX) * GX_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)  &
-              + GY_M_M(1,IKU,1,PTHLM,PDYY,PDZZ,PDZY) * GY_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)  &
+          *  (  GX_M_M(PTHLM,PDXX,PDZZ,PDZX) * GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)  &
+              + GY_M_M(PTHLM,PDYY,PDZZ,PDZY) * GY_M_M(PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)  &
              ) * (XCSHF+ZCSV) / (2.*ZCTSVD)
     ELSE
       ZFLX(:,:,:)=  PLM(:,:,:) * PLEPS(:,:,:)                                          &
-              * GX_M_M(1,IKU,1,PTHLM,PDXX,PDZZ,PDZX) * GX_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)  &
+              * GX_M_M(PTHLM,PDXX,PDZZ,PDZX) * GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)  &
               * (XCSHF+ZCSV) / (2.*ZCTSVD)
     END IF
     CALL LES_MEAN_SUBGRID( ZA*ZFLX, X_LES_SUBGRID_SvThv(:,:,:,JSV) , .TRUE.)
@@ -199,12 +197,12 @@ DO JSV=1,NSV
       ZA(:,:,:)   =  EMOIST(KRR,KRRI,PTHLM,PRM,PLOCPEXNM,PAMOIST,PSRCM)
       IF (.NOT. L2D) THEN
         ZFLX(:,:,:)=  PLM(:,:,:) * PLEPS(:,:,:)                                                 &
-            *  (  GX_M_M(1,IKU,1,PRM(:,:,:,1),PDXX,PDZZ,PDZX) * GX_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)  &
-                + GY_M_M(1,IKU,1,PRM(:,:,:,1),PDYY,PDZZ,PDZY) * GY_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)  &
+            *  (  GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX) * GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)  &
+                + GY_M_M(PRM(:,:,:,1),PDYY,PDZZ,PDZY) * GY_M_M(PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)  &
                ) * (XCHF+ZCSV) / (2.*ZCQSVD)
       ELSE
         ZFLX(:,:,:)=  PLM(:,:,:) * PLEPS(:,:,:)                                                 &
-                * GX_M_M(1,IKU,1,PRM(:,:,:,1),PDXX,PDZZ,PDZX) * GX_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)  &
+                * GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX) * GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)  &
                 * (XCHF+ZCSV) / (2.*ZCQSVD)
       END IF
       CALL LES_MEAN_SUBGRID( ZA*ZFLX, X_LES_SUBGRID_SvThv(:,:,:,JSV) , .TRUE.)
diff --git a/src/MNH/turb_hor_sv_flux.f90 b/src/MNH/turb_hor_sv_flux.f90
index efc04b01e..0bf50d583 100644
--- a/src/MNH/turb_hor_sv_flux.f90
+++ b/src/MNH/turb_hor_sv_flux.f90
@@ -268,9 +268,9 @@ DO JSV=1,ISV
   IF (LLES_CALL .AND. KSPLT==1) THEN
     CALL SECOND_MNH(ZTIME1)
     CALL LES_MEAN_SUBGRID( MXF(ZFLXX), X_LES_SUBGRID_USv(:,:,:,JSV) ) 
-    CALL LES_MEAN_SUBGRID( MZF(MXF(GX_W_UW(1,IKU,1,PWM,PDXX,PDZZ,PDZX)*MZM(ZFLXX))), &
+    CALL LES_MEAN_SUBGRID( MZF(MXF(GX_W_UW(PWM,PDXX,PDZZ,PDZX)*MZM(ZFLXX))), &
                            X_LES_RES_ddxa_W_SBG_UaSv(:,:,:,JSV) , .TRUE. )
-    CALL LES_MEAN_SUBGRID( GX_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*MXF(ZFLXX), &
+    CALL LES_MEAN_SUBGRID( GX_M_M(PSVM(:,:,:,JSV),PDXX,PDZZ,PDZX)*MXF(ZFLXX), &
                            X_LES_RES_ddxa_Sv_SBG_UaSv(:,:,:,JSV), .TRUE. )
     CALL SECOND_MNH(ZTIME2)
     XTIME_LES = XTIME_LES + ZTIME2 - ZTIME1
@@ -323,9 +323,9 @@ DO JSV=1,ISV
   IF (LLES_CALL .AND. KSPLT==1) THEN
     CALL SECOND_MNH(ZTIME1)
     CALL LES_MEAN_SUBGRID( MYF(ZFLXY), X_LES_SUBGRID_VSv(:,:,:,JSV) ) 
-    CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(1,IKU,1,PWM,PDYY,PDZZ,PDZY)*MZM(ZFLXY))), &
+    CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(PWM,PDYY,PDZZ,PDZY)*MZM(ZFLXY))), &
                            X_LES_RES_ddxa_W_SBG_UaSv(:,:,:,JSV) , .TRUE. )
-    CALL LES_MEAN_SUBGRID( GY_M_M(1,IKU,1,PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)*MYF(ZFLXY), &
+    CALL LES_MEAN_SUBGRID( GY_M_M(PSVM(:,:,:,JSV),PDYY,PDZZ,PDZY)*MYF(ZFLXY), &
                            X_LES_RES_ddxa_Sv_SBG_UaSv(:,:,:,JSV) , .TRUE. )
     CALL SECOND_MNH(ZTIME2)
     XTIME_LES = XTIME_LES + ZTIME2 - ZTIME1
diff --git a/src/MNH/turb_hor_thermo_corr.f90 b/src/MNH/turb_hor_thermo_corr.f90
index 8c9462a2c..aa4e22c41 100644
--- a/src/MNH/turb_hor_thermo_corr.f90
+++ b/src/MNH/turb_hor_thermo_corr.f90
@@ -194,7 +194,7 @@ REAL, DIMENSION(SIZE(PTHLM,1),SIZE(PTHLM,2),SIZE(PTHLM,3))       &
                                      :: ZFLX,ZWORK,ZA
     ! work arrays
 !   
-INTEGER             :: IKB,IKE,IKU
+INTEGER             :: IKB,IKE
                                     ! Index values for the Beginning and End
                                     ! mass points of the domain  
 REAL, DIMENSION(SIZE(PDZZ,1),SIZE(PDZZ,2),1+JPVEXT:3+JPVEXT) :: ZCOEFF 
@@ -210,7 +210,6 @@ TYPE(TFIELDDATA) :: TZFIELD
 !
 IKB = 1+JPVEXT               
 IKE = SIZE(PTHLM,3)-JPVEXT   
-IKU = SIZE(PTHLM,3)
 !
 !
 !
@@ -237,10 +236,10 @@ IF ( ( KRRL > 0 .AND. OSUBG_COND) .OR. ( OTURB_FLX .AND. OCLOSE_OUT ) &
   ! Computes the horizontal variance <THl THl>
   IF (.NOT. L2D) THEN
     ZFLX(:,:,:) = XCTV * PLM(:,:,:) * PLEPS(:,:,:) *                           &
-       ( GX_M_M(1,IKU,1,PTHLM,PDXX,PDZZ,PDZX)**2 + GY_M_M(1,IKU,1,PTHLM,PDYY,PDZZ,PDZY)**2 )
+       ( GX_M_M(PTHLM,PDXX,PDZZ,PDZX)**2 + GY_M_M(PTHLM,PDYY,PDZZ,PDZY)**2 )
   ELSE
     ZFLX(:,:,:) = XCTV * PLM(:,:,:) * PLEPS(:,:,:) *                           &
-         GX_M_M(1,IKU,1,PTHLM,PDXX,PDZZ,PDZX)**2
+         GX_M_M(PTHLM,PDXX,PDZZ,PDZX)**2
   END IF
 !
 ! Compute the flux at the first inner U-point with an uncentred vertical  
@@ -306,13 +305,13 @@ IF ( ( KRRL > 0 .AND. OSUBG_COND) .OR. ( OTURB_FLX .AND. OCLOSE_OUT ) &
     IF (.NOT. L2D) THEN
       ZFLX(:,:,:)=                                                               &
             PLM(:,:,:) * PLEPS(:,:,:) *                                          &
-            (GX_M_M(1,IKU,1,PTHLM,PDXX,PDZZ,PDZX) * GX_M_M(1,IKU,1,PRM(:,:,:,1),PDXX,PDZZ,PDZX)  &
-           + GY_M_M(1,IKU,1,PTHLM,PDYY,PDZZ,PDZY) * GY_M_M(1,IKU,1,PRM(:,:,:,1),PDYY,PDZZ,PDZY)  &
+            (GX_M_M(PTHLM,PDXX,PDZZ,PDZX) * GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)  &
+           + GY_M_M(PTHLM,PDYY,PDZZ,PDZY) * GY_M_M(PRM(:,:,:,1),PDYY,PDZZ,PDZY)  &
             ) * (XCHT1+XCHT2)
     ELSE
       ZFLX(:,:,:)=                                                               &
             PLM(:,:,:) * PLEPS(:,:,:) *                                          &
-            (GX_M_M(1,IKU,1,PTHLM,PDXX,PDZZ,PDZX) * GX_M_M(1,IKU,1,PRM(:,:,:,1),PDXX,PDZZ,PDZX)  &
+            (GX_M_M(PTHLM,PDXX,PDZZ,PDZX) * GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)  &
             ) * (XCHT1+XCHT2)
 
     END IF
@@ -393,11 +392,11 @@ IF ( ( KRRL > 0 .AND. OSUBG_COND) .OR. ( OTURB_FLX .AND. OCLOSE_OUT ) &
     ! Computes the horizontal variance <Rnp Rnp>
     IF (.NOT. L2D) THEN
       ZFLX(:,:,:) = XCHV * PLM(:,:,:) * PLEPS(:,:,:) *                      &
-           ( GX_M_M(1,IKU,1,PRM(:,:,:,1),PDXX,PDZZ,PDZX)**2 +                       &
-             GY_M_M(1,IKU,1,PRM(:,:,:,1),PDYY,PDZZ,PDZY)**2 )
+           ( GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)**2 +                       &
+             GY_M_M(PRM(:,:,:,1),PDYY,PDZZ,PDZY)**2 )
     ELSE
       ZFLX(:,:,:) = XCHV * PLM(:,:,:) * PLEPS(:,:,:) *                      &
-           ( GX_M_M(1,IKU,1,PRM(:,:,:,1),PDXX,PDZZ,PDZX)**2  )
+           ( GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)**2  )
     END IF
 !
 ! Compute the flux at the first inner U-point with an uncentred vertical  
diff --git a/src/MNH/turb_hor_thermo_flux.f90 b/src/MNH/turb_hor_thermo_flux.f90
index c2983d6e2..aa88862f5 100644
--- a/src/MNH/turb_hor_thermo_flux.f90
+++ b/src/MNH/turb_hor_thermo_flux.f90
@@ -331,12 +331,12 @@ END IF
 IF (KSPLT==1 .AND. LLES_CALL) THEN
   CALL SECOND_MNH(ZTIME1)
   CALL LES_MEAN_SUBGRID( MXF(ZFLX), X_LES_SUBGRID_UThl ) 
-  CALL LES_MEAN_SUBGRID( MZF(MXF(GX_W_UW(1,IKU,1,PWM,PDXX,PDZZ,PDZX)*MZM(ZFLX))),&
+  CALL LES_MEAN_SUBGRID( MZF(MXF(GX_W_UW(PWM,PDXX,PDZZ,PDZX)*MZM(ZFLX))),&
                          X_LES_RES_ddxa_W_SBG_UaThl , .TRUE. )
-  CALL LES_MEAN_SUBGRID( GX_M_M(1,IKU,1,PTHLM,PDXX,PDZZ,PDZX)*MXF(ZFLX),&
+  CALL LES_MEAN_SUBGRID( GX_M_M(PTHLM,PDXX,PDZZ,PDZX)*MXF(ZFLX),&
                          X_LES_RES_ddxa_Thl_SBG_UaThl , .TRUE. )
   IF (KRR>=1) THEN
-    CALL LES_MEAN_SUBGRID( GX_M_M(1,IKU,1,PRM(:,:,:,1),PDXX,PDZZ,PDZX)*MXF(ZFLX), &
+    CALL LES_MEAN_SUBGRID( GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)*MXF(ZFLX), &
                            X_LES_RES_ddxa_Rt_SBG_UaThl , .TRUE. )
   END IF
   CALL SECOND_MNH(ZTIME2)
@@ -434,11 +434,11 @@ IF (KRR/=0) THEN
   IF (KSPLT==1 .AND. LLES_CALL) THEN
     CALL SECOND_MNH(ZTIME1)
     CALL LES_MEAN_SUBGRID( MXF(ZFLX), X_LES_SUBGRID_URt ) 
-    CALL LES_MEAN_SUBGRID( MZF(MXF(GX_W_UW(1,IKU,1,PWM,PDXX,PDZZ,PDZX)*MZM(ZFLX))),&
+    CALL LES_MEAN_SUBGRID( MZF(MXF(GX_W_UW(PWM,PDXX,PDZZ,PDZX)*MZM(ZFLX))),&
                            X_LES_RES_ddxa_W_SBG_UaRt , .TRUE. )
-    CALL LES_MEAN_SUBGRID( GX_M_M(1,IKU,1,PTHLM,PDXX,PDZZ,PDZX)*MXF(ZFLX),&
+    CALL LES_MEAN_SUBGRID( GX_M_M(PTHLM,PDXX,PDZZ,PDZX)*MXF(ZFLX),&
                            X_LES_RES_ddxa_Thl_SBG_UaRt , .TRUE. )
-    CALL LES_MEAN_SUBGRID( GX_M_M(1,IKU,1,PRM(:,:,:,1),PDXX,PDZZ,PDZX)*MXF(ZFLX),&
+    CALL LES_MEAN_SUBGRID( GX_M_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)*MXF(ZFLX),&
                            X_LES_RES_ddxa_Rt_SBG_UaRt , .TRUE. )
     CALL SECOND_MNH(ZTIME2)
     XTIME_LES = XTIME_LES + ZTIME2 - ZTIME1
@@ -586,12 +586,12 @@ END IF
 IF (KSPLT==1 .AND. LLES_CALL) THEN
   CALL SECOND_MNH(ZTIME1)
   CALL LES_MEAN_SUBGRID( MYF(ZFLX), X_LES_SUBGRID_VThl ) 
-  CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(1,IKU,1,PWM,PDYY,PDZZ,PDZY)*MZM(ZFLX))),&
+  CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(PWM,PDYY,PDZZ,PDZY)*MZM(ZFLX))),&
                          X_LES_RES_ddxa_W_SBG_UaThl , .TRUE. )
-  CALL LES_MEAN_SUBGRID( GY_M_M(1,IKU,1,PTHLM,PDYY,PDZZ,PDZY)*MYF(ZFLX),&
+  CALL LES_MEAN_SUBGRID( GY_M_M(PTHLM,PDYY,PDZZ,PDZY)*MYF(ZFLX),&
                          X_LES_RES_ddxa_Thl_SBG_UaThl , .TRUE. )
   IF (KRR>=1) THEN
-    CALL LES_MEAN_SUBGRID( GY_M_M(1,IKU,1,PRM(:,:,:,1),PDYY,PDZZ,PDZY)*MYF(ZFLX),&
+    CALL LES_MEAN_SUBGRID( GY_M_M(PRM(:,:,:,1),PDYY,PDZZ,PDZY)*MYF(ZFLX),&
                            X_LES_RES_ddxa_Rt_SBG_UaThl , .TRUE. )
   END IF
   CALL SECOND_MNH(ZTIME2)
@@ -698,11 +698,11 @@ IF (KRR/=0) THEN
   IF (KSPLT==1 .AND. LLES_CALL) THEN
     CALL SECOND_MNH(ZTIME1)
     CALL LES_MEAN_SUBGRID( MYF(ZFLX), X_LES_SUBGRID_VRt ) 
-    CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(1,IKU,1,PWM,PDYY,PDZZ,PDZY)*MZM(ZFLX))),&
+    CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(PWM,PDYY,PDZZ,PDZY)*MZM(ZFLX))),&
                            X_LES_RES_ddxa_W_SBG_UaRt , .TRUE. )
-    CALL LES_MEAN_SUBGRID( GY_M_M(1,IKU,1,PTHLM,PDYY,PDZZ,PDZY)*MYF(ZFLX), &
+    CALL LES_MEAN_SUBGRID( GY_M_M(PTHLM,PDYY,PDZZ,PDZY)*MYF(ZFLX), &
                            X_LES_RES_ddxa_Thl_SBG_UaRt , .TRUE. )
-    CALL LES_MEAN_SUBGRID( GY_M_M(1,IKU,1,PRM(:,:,:,1),PDYY,PDZZ,PDZY)*MYF(ZFLX), &
+    CALL LES_MEAN_SUBGRID( GY_M_M(PRM(:,:,:,1),PDYY,PDZZ,PDZY)*MYF(ZFLX), &
                            X_LES_RES_ddxa_Rt_SBG_UaRt , .TRUE. )
     CALL SECOND_MNH(ZTIME2)
     XTIME_LES = XTIME_LES + ZTIME2 - ZTIME1
diff --git a/src/MNH/turb_hor_uv.f90 b/src/MNH/turb_hor_uv.f90
index e73d5bfdf..ae1a38a69 100644
--- a/src/MNH/turb_hor_uv.f90
+++ b/src/MNH/turb_hor_uv.f90
@@ -204,7 +204,7 @@ REAL, DIMENSION(SIZE(PUM,1),SIZE(PUM,2),SIZE(PUM,3))       &
 !   
 REAL, DIMENSION(SIZE(PUM,1),SIZE(PUM,2)) ::ZDIRSINZW 
       ! sinus of the angle between the vertical and the normal to the orography
-INTEGER             :: IKB,IKE,IKU
+INTEGER             :: IKB,IKE
                                     ! Index values for the Beginning and End
                                     ! mass points of the domain  
 !
@@ -220,12 +220,11 @@ TYPE(TFIELDDATA) :: TZFIELD
 !
 IKB = 1+JPVEXT               
 IKE = SIZE(PUM,3)-JPVEXT    
-IKU = SIZE(PUM,3)
 !
 ZDIRSINZW(:,:) = SQRT( 1. - PDIRCOSZW(:,:)**2 )
 !
-GX_V_UV_PVM = GX_V_UV(1,IKU,1,PVM,PDXX,PDZZ,PDZX)
-IF (.NOT. L2D) GY_U_UV_PUM = GY_U_UV(1,IKU,1,PUM,PDYY,PDZZ,PDZY)
+GX_V_UV_PVM = GX_V_UV(PVM,PDXX,PDZZ,PDZX)
+IF (.NOT. L2D) GY_U_UV_PUM = GY_U_UV(PUM,PDYY,PDZZ,PDZY)
 !
 !
 !*      12.   < U'V'>
@@ -350,8 +349,8 @@ END IF
 IF (LLES_CALL .AND. KSPLT==1) THEN
   CALL SECOND_MNH(ZTIME1)
   CALL LES_MEAN_SUBGRID( MXF(MYF(ZFLX)), X_LES_SUBGRID_UV ) 
-  CALL LES_MEAN_SUBGRID( MXF(MYF(GY_U_UV(1,IKU,1,PUM,PDYY,PDZZ,PDZY)*ZFLX)), X_LES_RES_ddxa_U_SBG_UaU , .TRUE.)
-  CALL LES_MEAN_SUBGRID( MXF(MYF(GX_V_UV(1,IKU,1,PVM,PDXX,PDZZ,PDZX)*ZFLX)), X_LES_RES_ddxa_V_SBG_UaV , .TRUE.)
+  CALL LES_MEAN_SUBGRID( MXF(MYF(GY_U_UV(PUM,PDYY,PDZZ,PDZY)*ZFLX)), X_LES_RES_ddxa_U_SBG_UaU , .TRUE.)
+  CALL LES_MEAN_SUBGRID( MXF(MYF(GX_V_UV(PVM,PDXX,PDZZ,PDZX)*ZFLX)), X_LES_RES_ddxa_V_SBG_UaV , .TRUE.)
   CALL SECOND_MNH(ZTIME2)
   XTIME_LES = XTIME_LES + ZTIME2 - ZTIME1
 END IF
diff --git a/src/MNH/turb_hor_uw.f90 b/src/MNH/turb_hor_uw.f90
index 2631fe365..c1ec29645 100644
--- a/src/MNH/turb_hor_uw.f90
+++ b/src/MNH/turb_hor_uw.f90
@@ -202,7 +202,7 @@ IKE = SIZE(PWM,3)-JPVEXT
 IKU = SIZE(PWM,3)
 !
 !
-GX_W_UW_PWM = GX_W_UW(1,IKU,1,PWM,PDXX,PDZZ,PDZX)
+GX_W_UW_PWM = GX_W_UW(PWM,PDXX,PDZZ,PDZX)
 !
 !
 !*      13.   < U'W'>
@@ -256,7 +256,7 @@ IF (KSPLT==1) THEN
   !Contribution to the dynamic production of TKE:
   !
   ZWORK(:,:,:) =-MZF( MXF(                               &
-     ZFLX *( GZ_U_UW(1,IKU,1,PUM,PDZZ) + GX_W_UW_PWM ) ) )
+     ZFLX *( GZ_U_UW(PUM,PDZZ) + GX_W_UW_PWM ) ) )
   !
   !
   ! evaluate the dynamic production at w(IKB+1) in PDP(IKB)
@@ -283,7 +283,7 @@ END IF
 IF (LLES_CALL .AND. KSPLT==1) THEN
   CALL SECOND_MNH(ZTIME1)
   CALL LES_MEAN_SUBGRID( MZF(MXF(ZFLX)), X_LES_SUBGRID_WU , .TRUE. )
-  CALL LES_MEAN_SUBGRID( MZF(MXF(GZ_U_UW(1,IKU,1,PUM,PDZZ)*ZFLX)), X_LES_RES_ddxa_U_SBG_UaU , .TRUE.)
+  CALL LES_MEAN_SUBGRID( MZF(MXF(GZ_U_UW(PUM,PDZZ)*ZFLX)), X_LES_RES_ddxa_U_SBG_UaU , .TRUE.)
   CALL LES_MEAN_SUBGRID( MZF(MXF(GX_W_UW_PWM*ZFLX)), X_LES_RES_ddxa_W_SBG_UaW , .TRUE.)
   CALL LES_MEAN_SUBGRID( MXF(GX_M_U(1,IKU,1,PTHLM,PDXX,PDZZ,PDZX)*MZF(ZFLX)),&
                          X_LES_RES_ddxa_Thl_SBG_UaW , .TRUE.)
diff --git a/src/MNH/turb_hor_vw.f90 b/src/MNH/turb_hor_vw.f90
index f8354f58f..ddcb83c53 100644
--- a/src/MNH/turb_hor_vw.f90
+++ b/src/MNH/turb_hor_vw.f90
@@ -197,7 +197,7 @@ IKE = SIZE(PWM,3)-JPVEXT
 IKU = SIZE(PWM,3)
 !
 !
-IF (.NOT. L2D) GY_W_VW_PWM = GY_W_VW(1,IKU,1,PWM,PDYY,PDZZ,PDZY)
+IF (.NOT. L2D) GY_W_VW_PWM = GY_W_VW(PWM,PDYY,PDZZ,PDZY)
 !
 !
 !*      14.   < V'W'>
@@ -260,7 +260,7 @@ IF (KSPLT==1) THEN
   !Contribution to the dynamic production of TKE:
   !
   IF (.NOT. L2D) THEN
-    ZWORK(:,:,:) =-MZF( MYF( ZFLX *( GZ_V_VW(1,IKU,1,PVM,PDZZ) + GY_W_VW_PWM ) ) )
+    ZWORK(:,:,:) =-MZF( MYF( ZFLX *( GZ_V_VW(PVM,PDZZ) + GY_W_VW_PWM ) ) )
   !
   !
   ! evaluate the dynamic production at w(IKB+1) in PDP(IKB)
@@ -288,9 +288,9 @@ END IF
 IF (LLES_CALL .AND. KSPLT==1) THEN
   CALL SECOND_MNH(ZTIME1)
   CALL LES_MEAN_SUBGRID( MZF(MYF(ZFLX)), X_LES_SUBGRID_WV , .TRUE. )
-  CALL LES_MEAN_SUBGRID( MZF(MYF(GZ_V_VW(1,IKU,1,PVM,PDZZ)*ZFLX)),&
+  CALL LES_MEAN_SUBGRID( MZF(MYF(GZ_V_VW(PVM,PDZZ)*ZFLX)),&
                          X_LES_RES_ddxa_V_SBG_UaV , .TRUE.)
-  CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(1,IKU,1,PWM,PDYY,PDZZ,PDZY)*ZFLX)),&
+  CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(PWM,PDYY,PDZZ,PDZY)*ZFLX)),&
                          X_LES_RES_ddxa_W_SBG_UaW , .TRUE.)
   CALL LES_MEAN_SUBGRID( MXF(GY_M_V(1,IKU,1,PTHLM,PDYY,PDZZ,PDZY)*MZF(ZFLX)),&
                          X_LES_RES_ddxa_Thl_SBG_UaW , .TRUE.)
diff --git a/src/MNH/turb_ver_dyn_flux.f90 b/src/MNH/turb_ver_dyn_flux.f90
index 58b0edcf1..d7ef5a5df 100644
--- a/src/MNH/turb_ver_dyn_flux.f90
+++ b/src/MNH/turb_ver_dyn_flux.f90
@@ -531,7 +531,7 @@ PWU(:,:,:) = ZFLXZ(:,:,:)
 ! Contribution to the dynamic production of TKE
 ! compute the dynamic production at the mass point
 !
-PDP(:,:,:) = - MZF( MXF ( ZFLXZ * GZ_U_UW(KKA,KKU,KKL,PUM,PDZZ) )  )
+PDP(:,:,:) = - MZF( MXF ( ZFLXZ * GZ_U_UW(PUM,PDZZ) )  )
 !
 ! evaluate the dynamic production at w(IKB+KKL) in PDP(IKB)
 PDP(:,:,IKB:IKB) = - MXF (                                                      &
@@ -544,7 +544,7 @@ PDP(:,:,IKB:IKB) = - MXF (
 IF (LLES_CALL) THEN
   CALL SECOND_MNH(ZTIME1)
   CALL LES_MEAN_SUBGRID( MZF(MXF(ZFLXZ)), X_LES_SUBGRID_WU )
-  CALL LES_MEAN_SUBGRID( MZF(MXF(GZ_U_UW(KKA,KKU,KKL,PUM,PDZZ) &
+  CALL LES_MEAN_SUBGRID( MZF(MXF(GZ_U_UW(PUM,PDZZ) &
                           & *ZFLXZ)), X_LES_RES_ddxa_U_SBG_UaU )
   CALL LES_MEAN_SUBGRID( XCMFS * ZKEFF, X_LES_SUBGRID_Km )
   CALL SECOND_MNH(ZTIME2)
@@ -571,7 +571,7 @@ IF(HTURBDIM=='3DIM') THEN
   !
   ! Complete the Dynamical production with the W wind component 
   !
-  ZA(:,:,:)=-MZF( MXF ( ZFLXZ * GX_W_UW(KKA,KKU,KKL, PWM,PDXX,PDZZ,PDZX) )  )
+  ZA(:,:,:)=-MZF( MXF ( ZFLXZ * GX_W_UW( PWM,PDXX,PDZZ,PDZX) )  )
   !
   !
   ! evaluate the dynamic production at w(IKB+KKL) in PDP(IKB)
@@ -593,16 +593,16 @@ IF(HTURBDIM=='3DIM') THEN
   ! 
   IF (LLES_CALL) THEN
     CALL SECOND_MNH(ZTIME1)
-    CALL LES_MEAN_SUBGRID( MZF(MXF(GX_W_UW(KKA,KKU,KKL,PWM,PDXX,&
+    CALL LES_MEAN_SUBGRID( MZF(MXF(GX_W_UW(PWM,PDXX,&
       PDZZ,PDZX)*ZFLXZ)), X_LES_RES_ddxa_W_SBG_UaW )
     CALL LES_MEAN_SUBGRID( MXF(GX_M_U(KKA,KKU,KKL,PTHLM,PDXX,PDZZ,PDZX)&
       * MZF(ZFLXZ)), X_LES_RES_ddxa_Thl_SBG_UaW )
     IF (KRR>=1) THEN
-      CALL LES_MEAN_SUBGRID(MXF(GX_U_M(KKA,KKU,KKL,PRM(:,:,:,1),PDXX,PDZZ,PDZX)&
+      CALL LES_MEAN_SUBGRID(MXF(GX_U_M(PRM(:,:,:,1),PDXX,PDZZ,PDZX)&
       *MZF(ZFLXZ)),X_LES_RES_ddxa_Rt_SBG_UaW )
     END IF
     DO JSV=1,NSV
-      CALL LES_MEAN_SUBGRID( MXF(GX_U_M(KKA,KKU,KKL,PSVM(:,:,:,JSV),PDXX,PDZZ,&
+      CALL LES_MEAN_SUBGRID( MXF(GX_U_M(PSVM(:,:,:,JSV),PDXX,PDZZ,&
       PDZX)*MZF(ZFLXZ)),X_LES_RES_ddxa_Sv_SBG_UaW(:,:,:,JSV) )
     END DO
     CALL SECOND_MNH(ZTIME2)
@@ -705,7 +705,7 @@ PWV(:,:,:) = ZFLXZ(:,:,:)
 !  Contribution to the dynamic production of TKE
 ! compute the dynamic production contribution at the mass point
 !
-ZA(:,:,:) = - MZF( MYF ( ZFLXZ * GZ_V_VW(KKA,KKU,KKL,PVM,PDZZ) ) )
+ZA(:,:,:) = - MZF( MYF ( ZFLXZ * GZ_V_VW(PVM,PDZZ) ) )
 !
 ! evaluate the dynamic production at w(IKB+KKL) in PDP(IKB)
 ZA(:,:,IKB:IKB)  =                                                 &
@@ -721,7 +721,7 @@ PDP(:,:,:)=PDP(:,:,:)+ZA(:,:,:)
 IF (LLES_CALL) THEN
   CALL SECOND_MNH(ZTIME1)
   CALL LES_MEAN_SUBGRID( MZF(MYF(ZFLXZ)), X_LES_SUBGRID_WV )
-  CALL LES_MEAN_SUBGRID( MZF(MYF(GZ_V_VW(KKA,KKU,KKL,PVM,PDZZ)*&
+  CALL LES_MEAN_SUBGRID( MZF(MYF(GZ_V_VW(PVM,PDZZ)*&
                     & ZFLXZ)), X_LES_RES_ddxa_V_SBG_UaV )
   CALL SECOND_MNH(ZTIME2)
   XTIME_LES = XTIME_LES + ZTIME2 - ZTIME1
@@ -748,7 +748,7 @@ IF(HTURBDIM=='3DIM') THEN
   ! 
   ! Complete the Dynamical production with the W wind component 
   IF (.NOT. L2D) THEN
-    ZA(:,:,:) = - MZF( MYF ( ZFLXZ * GY_W_VW(KKA,KKU,KKL, PWM,PDYY,PDZZ,PDZY) )  )
+    ZA(:,:,:) = - MZF( MYF ( ZFLXZ * GY_W_VW( PWM,PDYY,PDZZ,PDZY) )  )
   !
   ! evaluate the dynamic production at w(IKB+KKL) in PDP(IKB)
     ZA(:,:,IKB:IKB) = - MYF (                                              &
@@ -771,12 +771,12 @@ IF(HTURBDIM=='3DIM') THEN
   !
   IF (LLES_CALL) THEN
     CALL SECOND_MNH(ZTIME1)
-    CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(KKA,KKU,KKL,PWM,PDYY,&
+    CALL LES_MEAN_SUBGRID( MZF(MYF(GY_W_VW(PWM,PDYY,&
     PDZZ,PDZY)*ZFLXZ)), X_LES_RES_ddxa_W_SBG_UaW , .TRUE. )
     CALL LES_MEAN_SUBGRID( MYF(GY_M_V(KKA,KKU,KKL,PTHLM,PDYY,PDZZ,PDZY)&
     *MZF(ZFLXZ)), X_LES_RES_ddxa_Thl_SBG_UaW , .TRUE. )
     IF (KRR>=1) THEN
-      CALL LES_MEAN_SUBGRID( MYF(GY_V_M(KKA,KKU,KKL,PRM(:,:,:,1),PDYY,PDZZ,&
+      CALL LES_MEAN_SUBGRID( MYF(GY_V_M(PRM(:,:,:,1),PDYY,PDZZ,&
       PDZY)*MZF(ZFLXZ)),X_LES_RES_ddxa_Rt_SBG_UaW , .TRUE. )
     END IF
     CALL SECOND_MNH(ZTIME2)
@@ -793,7 +793,7 @@ END IF
 !
 IF ( OTURB_FLX .AND. OCLOSE_OUT .AND. HTURBDIM == '1DIM') THEN
   ZFLXZ(:,:,:)= (2./3.) * PTKEM(:,:,:)                     &
-     -XCMFS*PLM(:,:,:)*SQRT(PTKEM(:,:,:))*GZ_W_M(KKA,KKU,KKL,PWM,PDZZ)
+     -XCMFS*PLM(:,:,:)*SQRT(PTKEM(:,:,:))*GZ_W_M(PWM,PDZZ)
   ! to be tested &
   !   +XCMFB*(4./3.)*PLM(:,:,:)/SQRT(PTKEM(:,:,:))*PTP(:,:,:) 
   ! stores the W variance
diff --git a/src/MNH/turb_ver_sv_flux.f90 b/src/MNH/turb_ver_sv_flux.f90
index 8c51e8e42..6603b8ad6 100644
--- a/src/MNH/turb_ver_sv_flux.f90
+++ b/src/MNH/turb_ver_sv_flux.f90
@@ -477,7 +477,7 @@ DO JSV=1,ISV
   IF (LLES_CALL) THEN
     CALL SECOND_MNH(ZTIME1)
     CALL LES_MEAN_SUBGRID( MZF(ZFLXZ), X_LES_SUBGRID_WSv(:,:,:,JSV) )
-    CALL LES_MEAN_SUBGRID( GZ_W_M(KKA,KKU,KKL,PWM,PDZZ)*MZF(ZFLXZ), &
+    CALL LES_MEAN_SUBGRID( GZ_W_M(PWM,PDZZ)*MZF(ZFLXZ), &
                            X_LES_RES_ddxa_W_SBG_UaSv(:,:,:,JSV) )
     CALL LES_MEAN_SUBGRID( MZF(GZ_M_W(KKA,KKU,KKL,PSVM(:,:,:,JSV),PDZZ)*ZFLXZ), &
                            X_LES_RES_ddxa_Sv_SBG_UaSv(:,:,:,JSV) )
diff --git a/src/MNH/turb_ver_thermo_flux.f90 b/src/MNH/turb_ver_thermo_flux.f90
index cbec5013a..9ee8de8bf 100644
--- a/src/MNH/turb_ver_thermo_flux.f90
+++ b/src/MNH/turb_ver_thermo_flux.f90
@@ -656,7 +656,7 @@ IF (LLES_CALL) THEN
   CALL SECOND_MNH(ZTIME1)
   CALL LES_MEAN_SUBGRID( MZF(ZFLXZ), X_LES_SUBGRID_WThl )
   CALL LES_MEAN_SUBGRID( MZF(PWM*ZFLXZ), X_LES_RES_W_SBG_WThl )
-  CALL LES_MEAN_SUBGRID( GZ_W_M(KKA,KKU,KKL,PWM,PDZZ)*MZF(ZFLXZ),&
+  CALL LES_MEAN_SUBGRID( GZ_W_M(PWM,PDZZ)*MZF(ZFLXZ),&
       & X_LES_RES_ddxa_W_SBG_UaThl )
   CALL LES_MEAN_SUBGRID( MZF(PDTH_DZ*ZFLXZ), X_LES_RES_ddxa_Thl_SBG_UaThl )
   CALL LES_MEAN_SUBGRID( -XCTP*PSQRT_TKE/PLM*MZF(ZFLXZ), X_LES_SUBGRID_ThlPz )
@@ -831,7 +831,7 @@ IF (KRR /= 0) THEN
     CALL SECOND_MNH(ZTIME1)
     CALL LES_MEAN_SUBGRID( MZF(ZFLXZ), X_LES_SUBGRID_WRt )
     CALL LES_MEAN_SUBGRID( MZF(PWM*ZFLXZ), X_LES_RES_W_SBG_WRt )
-    CALL LES_MEAN_SUBGRID( GZ_W_M(KKA,KKU,KKL,PWM,PDZZ)*MZF(ZFLXZ),&
+    CALL LES_MEAN_SUBGRID( GZ_W_M(PWM,PDZZ)*MZF(ZFLXZ),&
     & X_LES_RES_ddxa_W_SBG_UaRt )
     CALL LES_MEAN_SUBGRID( MZF(PDTH_DZ*ZFLXZ), X_LES_RES_ddxa_Thl_SBG_UaRt )
     CALL LES_MEAN_SUBGRID( MZF(PDR_DZ*ZFLXZ), X_LES_RES_ddxa_Rt_SBG_UaRt )
diff --git a/src/MNH/write_lfifm1_for_diag.f90 b/src/MNH/write_lfifm1_for_diag.f90
index 8f16be807..90eb9db9f 100644
--- a/src/MNH/write_lfifm1_for_diag.f90
+++ b/src/MNH/write_lfifm1_for_diag.f90
@@ -481,16 +481,16 @@ IF (INDEX(CISO,'TK') /= 0) THEN
 END IF
 !
 ZCORIOZ(:,:,:)=SPREAD( XCORIOZ(:,:),DIM=3,NCOPIES=IKU )
-ZVOX(:,:,:)=GY_W_VW(1,IKU,1,XWT,XDYY,XDZZ,XDZY)-GZ_V_VW(1,IKU,1,XVT,XDZZ)
+ZVOX(:,:,:)=GY_W_VW(XWT,XDYY,XDZZ,XDZY)-GZ_V_VW(XVT,XDZZ)
 ZVOX(:,:,2)=ZVOX(:,:,3)
-ZVOY(:,:,:)=GZ_U_UW(1,IKU,1,XUT,XDZZ)-GX_W_UW(1,IKU,1,XWT,XDXX,XDZZ,XDZX)
+ZVOY(:,:,:)=GZ_U_UW(XUT,XDZZ)-GX_W_UW(XWT,XDXX,XDZZ,XDZX)
 ZVOY(:,:,2)=ZVOY(:,:,3)
-ZVOZ(:,:,:)=GX_V_UV(1,IKU,1,XVT,XDXX,XDZZ,XDZX)-GY_U_UV(1,IKU,1,XUT,XDYY,XDZZ,XDZY)
+ZVOZ(:,:,:)=GX_V_UV(XVT,XDXX,XDZZ,XDZX)-GY_U_UV(XUT,XDYY,XDZZ,XDZY)
 ZVOZ(:,:,2)=ZVOZ(:,:,3)
 ZVOZ(:,:,1)=ZVOZ(:,:,3)
-ZWORK31(:,:,:)=GX_M_M(1,IKU,1,XTHT,XDXX,XDZZ,XDZX)
-ZWORK32(:,:,:)=GY_M_M(1,IKU,1,XTHT,XDYY,XDZZ,XDZY)
-ZWORK33(:,:,:)=GZ_M_M(1,IKU,1,XTHT,XDZZ)
+ZWORK31(:,:,:)=GX_M_M(XTHT,XDXX,XDZZ,XDZX)
+ZWORK32(:,:,:)=GY_M_M(XTHT,XDYY,XDZZ,XDZY)
+ZWORK33(:,:,:)=GZ_M_M(XTHT,XDZZ)
 ZPOVO(:,:,:)= ZWORK31(:,:,:)*MZF(MYF(ZVOX(:,:,:)))     &
              + ZWORK32(:,:,:)*MZF(MXF(ZVOY(:,:,:)))     &
              + ZWORK33(:,:,:)*(MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:))
@@ -743,13 +743,13 @@ IF (LHU_FLX) THEN
   ZWORK35(:,:,:) = XRHODREF(:,:,:) * XRT(:,:,:,1)
   ZWORK31(:,:,:) = MXM(ZWORK35(:,:,:)) * XUT(:,:,:)
   ZWORK32(:,:,:) = MYM(ZWORK35(:,:,:)) * XVT(:,:,:)
-  ZWORK35(:,:,:) = GX_U_M(1,IKU,1,ZWORK31,XDXX,XDZZ,XDZX) + GY_V_M(1,IKU,1,ZWORK32,XDYY,XDZZ,XDZY)
+  ZWORK35(:,:,:) = GX_U_M(ZWORK31,XDXX,XDZZ,XDZX) + GY_V_M(ZWORK32,XDYY,XDZZ,XDZY)
   IF  (CCLOUD(1:3) == 'ICE' .OR. CCLOUD == 'LIMA') THEN
     ZWORK36(:,:,:) = ZWORK35(:,:,:) + XRHODREF(:,:,:) * (XRT(:,:,:,2) + &
     XRT(:,:,:,3) + XRT(:,:,:,4) + XRT(:,:,:,5) + XRT(:,:,:,6))
     ZWORK33(:,:,:) = MXM(ZWORK36(:,:,:)) * XUT(:,:,:)
     ZWORK34(:,:,:) = MYM(ZWORK36(:,:,:)) * XVT(:,:,:)
-    ZWORK36(:,:,:) = GX_U_M(1,IKU,1,ZWORK33,XDXX,XDZZ,XDZX) + GY_V_M(1,IKU,1,ZWORK34,XDYY,XDZZ,XDZY)
+    ZWORK36(:,:,:) = GX_U_M(ZWORK33,XDXX,XDZZ,XDZX) + GY_V_M(ZWORK34,XDYY,XDZZ,XDZY)
   ENDIF
   !
   ! Integration sur 3000 m
@@ -2700,9 +2700,9 @@ END IF
 !
 ! Virtual Potential Vorticity in PV units
 IF (LMOIST_V .AND. (NRR>0) ) THEN
-  ZWORK31(:,:,:)=GX_M_M(1,IKU,1,ZTHETAV,XDXX,XDZZ,XDZX)
-  ZWORK32(:,:,:)=GY_M_M(1,IKU,1,ZTHETAV,XDYY,XDZZ,XDZY)
-  ZWORK33(:,:,:)=GZ_M_M(1,IKU,1,ZTHETAV,XDZZ)
+  ZWORK31(:,:,:)=GX_M_M(ZTHETAV,XDXX,XDZZ,XDZX)
+  ZWORK32(:,:,:)=GY_M_M(ZTHETAV,XDYY,XDZZ,XDZY)
+  ZWORK33(:,:,:)=GZ_M_M(ZTHETAV,XDZZ)
   ZWORK34(:,:,:)= ZWORK31(:,:,:)*MZF(MYF(ZVOX(:,:,:)))     &
                + ZWORK32(:,:,:)*MZF(MXF(ZVOY(:,:,:)))     &
                + ZWORK33(:,:,:)*(MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:))
@@ -2746,9 +2746,9 @@ END IF
 ! Equivalent Potential Vorticity in PV units
 IF (LMOIST_E .AND. (NRR>0) ) THEN
 !
-  ZWORK31(:,:,:)=GX_M_M(1,IKU,1,ZTHETAE,XDXX,XDZZ,XDZX)
-  ZWORK32(:,:,:)=GY_M_M(1,IKU,1,ZTHETAE,XDYY,XDZZ,XDZY)
-  ZWORK33(:,:,:)=GZ_M_M(1,IKU,1,ZTHETAE,XDZZ)
+  ZWORK31(:,:,:)=GX_M_M(ZTHETAE,XDXX,XDZZ,XDZX)
+  ZWORK32(:,:,:)=GY_M_M(ZTHETAE,XDYY,XDZZ,XDZY)
+  ZWORK33(:,:,:)=GZ_M_M(ZTHETAE,XDZZ)
   ZWORK34(:,:,:)= ZWORK31(:,:,:)*MZF(MYF(ZVOX(:,:,:)))     &
                 + ZWORK32(:,:,:)*MZF(MXF(ZVOY(:,:,:)))     &
                 + ZWORK33(:,:,:)*(MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:))
@@ -2793,9 +2793,9 @@ END IF
 !
 ! Equivalent Saturated Potential Vorticity in PV units
 IF (LMOIST_ES .AND. (NRR>0) ) THEN
-  ZWORK31(:,:,:)=GX_M_M(1,IKU,1,ZTHETAES,XDXX,XDZZ,XDZX)
-  ZWORK32(:,:,:)=GY_M_M(1,IKU,1,ZTHETAES,XDYY,XDZZ,XDZY)
-  ZWORK33(:,:,:)=GZ_M_M(1,IKU,1,ZTHETAES,XDZZ)
+  ZWORK31(:,:,:)=GX_M_M(ZTHETAES,XDXX,XDZZ,XDZX)
+  ZWORK32(:,:,:)=GY_M_M(ZTHETAES,XDYY,XDZZ,XDZY)
+  ZWORK33(:,:,:)=GZ_M_M(ZTHETAES,XDZZ)
   ZWORK34(:,:,:)= ZWORK31(:,:,:)*MZF(MYF(ZVOX(:,:,:)))     &
                 + ZWORK32(:,:,:)*MZF(MXF(ZVOY(:,:,:)))     &
                 + ZWORK33(:,:,:)*(MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:))
@@ -2820,7 +2820,7 @@ ENDIF
 !
 IF (LDIV) THEN
 !
-  ZWORK31=GX_U_M(1,IKU,1,XUT,XDXX,XDZZ,XDZX) + GY_V_M(1,IKU,1,XVT,XDYY,XDZZ,XDZY)
+  ZWORK31=GX_U_M(XUT,XDXX,XDZZ,XDZX) + GY_V_M(XVT,XDYY,XDZZ,XDZY)
   TZFIELD%CMNHNAME   = 'HDIV'
   TZFIELD%CSTDNAME   = ''
   TZFIELD%CLONGNAME  = 'HDIV'
@@ -2846,7 +2846,7 @@ IF (LDIV) THEN
     TZFIELD%LTIMEDEP   = .TRUE.
     ZWORK31=MXM(XRHODREF*XRT(:,:,:,1))*XUT
     ZWORK32=MYM(XRHODREF*XRT(:,:,:,1))*XVT
-    ZWORK33=GX_U_M(1,IKU,1,ZWORK31,XDXX,XDZZ,XDZX) + GY_V_M(1,IKU,1,ZWORK32,XDYY,XDZZ,XDZY)
+    ZWORK33=GX_U_M(ZWORK31,XDXX,XDZZ,XDZX) + GY_V_M(ZWORK32,XDYY,XDZZ,XDZY)
     CALL IO_Field_write(TPFILE,TZFIELD,ZWORK33)
   END IF
 !
@@ -2919,13 +2919,13 @@ IF (LGEO .OR. LAGEO) THEN
     ZPHI(1,IJU,:)=2*ZPHI(1,IJU-1,:)-ZPHI(1,IJU-2,:)
     ZPHI(IIU,1,:)=2*ZPHI(IIU,2,:)-ZPHI(IIU,3,:)
     ZPHI(IIU,IJU,:)=2*ZPHI(IIU,IJU-1,:)-ZPHI(IIU,IJU-2,:)
-    ZWORK31(:,:,:)=-MXM(GY_M_M(1,IKU,1,ZPHI,XDYY,XDZZ,XDZY)*XCPD*XTHVREF/ZCORIOZ)
+    ZWORK31(:,:,:)=-MXM(GY_M_M(ZPHI,XDYY,XDZZ,XDZY)*XCPD*XTHVREF/ZCORIOZ)
     !
     ZPHI(1,1,:)=2*ZPHI(2,1,:)-ZPHI(3,1,:)
     ZPHI(IIU,1,:)=2*ZPHI(IIU-1,1,:)-ZPHI(IIU-2,1,:)
     ZPHI(1,IJU,:)=2*ZPHI(2,IJU,:)-ZPHI(3,IJU,:)
     ZPHI(IIU,IJU,:)=2*ZPHI(IIU-1,IJU,:)-ZPHI(IIU-2,IJU,:)
-    ZWORK32(:,:,:)=MYM(GX_M_M(1,IKU,1,ZPHI,XDXX,XDZZ,XDZX)*XCPD*XTHVREF/ZCORIOZ)
+    ZWORK32(:,:,:)=MYM(GX_M_M(ZPHI,XDXX,XDZZ,XDZX)*XCPD*XTHVREF/ZCORIOZ)
   !
   ELSE IF(CEQNSYS=='LHE') THEN
     ZPHI(:,:,:)= ((XPABST(:,:,:)/XP00)**(XRD/XCPD)-XEXNREF(:,:,:))   &
@@ -2935,13 +2935,13 @@ IF (LGEO .OR. LAGEO) THEN
     ZPHI(1,IJU,:)=2*ZPHI(1,IJU-1,:)-ZPHI(1,IJU-2,:)
     ZPHI(IIU,1,:)=2*ZPHI(IIU,2,:)-ZPHI(IIU,3,:)
     ZPHI(IIU,IJU,:)=2*ZPHI(IIU,IJU-1,:)-ZPHI(IIU,IJU-2,:)
-    ZWORK31(:,:,:)=-MXM(GY_M_M(1,IKU,1,ZPHI,XDYY,XDZZ,XDZY)/ZCORIOZ)
+    ZWORK31(:,:,:)=-MXM(GY_M_M(ZPHI,XDYY,XDZZ,XDZY)/ZCORIOZ)
     !
     ZPHI(1,1,:)=2*ZPHI(2,1,:)-ZPHI(3,1,:)
     ZPHI(IIU,1,:)=2*ZPHI(IIU-1,1,:)-ZPHI(IIU-2,1,:)
     ZPHI(1,IJU,:)=2*ZPHI(2,IJU,:)-ZPHI(3,IJU,:)
     ZPHI(IIU,IJU,:)=2*ZPHI(IIU-1,IJU,:)-ZPHI(IIU-2,IJU,:)
-    ZWORK32(:,:,:)=MYM(GX_M_M(1,IKU,1,ZPHI,XDXX,XDZZ,XDZX)/ZCORIOZ)
+    ZWORK32(:,:,:)=MYM(GX_M_M(ZPHI,XDXX,XDZZ,XDZX)/ZCORIOZ)
   END IF
   DEALLOCATE(ZPHI)
 !
diff --git a/src/MNH/write_lfifm1_for_diag_supp.f90 b/src/MNH/write_lfifm1_for_diag_supp.f90
index a55567a21..03881cc17 100644
--- a/src/MNH/write_lfifm1_for_diag_supp.f90
+++ b/src/MNH/write_lfifm1_for_diag_supp.f90
@@ -1193,16 +1193,16 @@ ALLOCATE(ZWORK34(IIU,IJU,IKU))
 ! Potential Vorticity
 ! *********************
   ZCORIOZ(:,:,:)=SPREAD( XCORIOZ(:,:),DIM=3,NCOPIES=IKU )
-  ZVOX(:,:,:)=GY_W_VW(1,IKU,1,XWT,XDYY,XDZZ,XDZY)-GZ_V_VW(1,IKU,1,XVT,XDZZ)
+  ZVOX(:,:,:)=GY_W_VW(XWT,XDYY,XDZZ,XDZY)-GZ_V_VW(XVT,XDZZ)
   ZVOX(:,:,2)=ZVOX(:,:,3)
-  ZVOY(:,:,:)=GZ_U_UW(1,IKU,1,XUT,XDZZ)-GX_W_UW(1,IKU,1,XWT,XDXX,XDZZ,XDZX)
+  ZVOY(:,:,:)=GZ_U_UW(XUT,XDZZ)-GX_W_UW(XWT,XDXX,XDZZ,XDZX)
   ZVOY(:,:,2)=ZVOY(:,:,3)
-  ZVOZ(:,:,:)=GX_V_UV(1,IKU,1,XVT,XDXX,XDZZ,XDZX)-GY_U_UV(1,IKU,1,XUT,XDYY,XDZZ,XDZY)
+  ZVOZ(:,:,:)=GX_V_UV(XVT,XDXX,XDZZ,XDZX)-GY_U_UV(XUT,XDYY,XDZZ,XDZY)
   ZVOZ(:,:,2)=ZVOZ(:,:,3)
   ZVOZ(:,:,1)=ZVOZ(:,:,3)
-  ZWORK31(:,:,:)=GX_M_M(1,IKU,1,XTHT,XDXX,XDZZ,XDZX)
-  ZWORK32(:,:,:)=GY_M_M(1,IKU,1,XTHT,XDYY,XDZZ,XDZY)
-  ZWORK33(:,:,:)=GZ_M_M(1,IKU,1,XTHT,XDZZ)
+  ZWORK31(:,:,:)=GX_M_M(XTHT,XDXX,XDZZ,XDZX)
+  ZWORK32(:,:,:)=GY_M_M(XTHT,XDYY,XDZZ,XDZY)
+  ZWORK33(:,:,:)=GZ_M_M(XTHT,XDZZ)
   ZPOVO(:,:,:)= ZWORK31(:,:,:)*MZF(MYF(ZVOX(:,:,:)))     &
   + ZWORK32(:,:,:)*MZF(MXF(ZVOY(:,:,:)))     &
    + ZWORK33(:,:,:)*(MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:))
@@ -1342,16 +1342,16 @@ IF (LISOAL .AND.XISOAL(1)/=0.) THEN
 ! Potential Vorticity
 ! *********************
   ZCORIOZ(:,:,:)=SPREAD( XCORIOZ(:,:),DIM=3,NCOPIES=IKU )
-  ZVOX(:,:,:)=GY_W_VW(1,IKU,1,XWT,XDYY,XDZZ,XDZY)-GZ_V_VW(1,IKU,1,XVT,XDZZ)
+  ZVOX(:,:,:)=GY_W_VW(XWT,XDYY,XDZZ,XDZY)-GZ_V_VW(XVT,XDZZ)
   ZVOX(:,:,2)=ZVOX(:,:,3)
-  ZVOY(:,:,:)=GZ_U_UW(1,IKU,1,XUT,XDZZ)-GX_W_UW(1,IKU,1,XWT,XDXX,XDZZ,XDZX)
+  ZVOY(:,:,:)=GZ_U_UW(XUT,XDZZ)-GX_W_UW(XWT,XDXX,XDZZ,XDZX)
   ZVOY(:,:,2)=ZVOY(:,:,3)
-  ZVOZ(:,:,:)=GX_V_UV(1,IKU,1,XVT,XDXX,XDZZ,XDZX)-GY_U_UV(1,IKU,1,XUT,XDYY,XDZZ,XDZY)
+  ZVOZ(:,:,:)=GX_V_UV(XVT,XDXX,XDZZ,XDZX)-GY_U_UV(XUT,XDYY,XDZZ,XDZY)
   ZVOZ(:,:,2)=ZVOZ(:,:,3)
   ZVOZ(:,:,1)=ZVOZ(:,:,3)
-  ZWORK31(:,:,:)=GX_M_M(1,IKU,1,XTHT,XDXX,XDZZ,XDZX)
-  ZWORK32(:,:,:)=GY_M_M(1,IKU,1,XTHT,XDYY,XDZZ,XDZY)
-  ZWORK33(:,:,:)=GZ_M_M(1,IKU,1,XTHT,XDZZ)
+  ZWORK31(:,:,:)=GX_M_M(XTHT,XDXX,XDZZ,XDZX)
+  ZWORK32(:,:,:)=GY_M_M(XTHT,XDYY,XDZZ,XDZY)
+  ZWORK33(:,:,:)=GZ_M_M(XTHT,XDZZ)
   ZPOVO(:,:,:)= ZWORK31(:,:,:)*MZF(MYF(ZVOX(:,:,:)))     &
   + ZWORK32(:,:,:)*MZF(MXF(ZVOY(:,:,:)))     &
    + ZWORK33(:,:,:)*(MYF(MXF(ZVOZ(:,:,:))) + ZCORIOZ(:,:,:))
-- 
GitLab