From 89cfbde9e9ad8dfd534e7bd4dbe87eb52ecd418d Mon Sep 17 00:00:00 2001
From: Philippe WAUTELET <philippe.wautelet@aero.obs-mip.fr>
Date: Mon, 20 Mar 2023 13:32:40 +0100
Subject: [PATCH] Philippe 20/03/2023: ZSOLVER: remove files that are now
 similar to MNH/ files

---
 src/ZSOLVER/dyn_sources.f90 | 607 ------------------------------------
 src/ZSOLVER/update_lm.f90   | 185 -----------
 2 files changed, 792 deletions(-)
 delete mode 100644 src/ZSOLVER/dyn_sources.f90
 delete mode 100644 src/ZSOLVER/update_lm.f90

diff --git a/src/ZSOLVER/dyn_sources.f90 b/src/ZSOLVER/dyn_sources.f90
deleted file mode 100644
index aa8feb2fc..000000000
--- a/src/ZSOLVER/dyn_sources.f90
+++ /dev/null
@@ -1,607 +0,0 @@
-!MNH_LIC Copyright 1994-2021 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.
-!-----------------------------------------------------------------
-!     #######################
-      MODULE MODI_DYN_SOURCES
-!     #######################
-!
-INTERFACE
-!
-      SUBROUTINE DYN_SOURCES( KRR,KRRL, KRRI,                              &
-                              PUT, PVT, PWT, PTHT, PRT,                    &
-                              PCORIOX, PCORIOY, PCORIOZ, PCURVX, PCURVY,   &
-                              PRHODJ, PZZ, PTHVREF, PEXNREF,               &
-                              PRUS, PRVS, PRWS, PRTHS                      )
-!
-INTEGER,                  INTENT(IN)    :: KRR  ! Total number of water var.
-INTEGER,                  INTENT(IN)    :: KRRL ! Number of liquid water var.
-INTEGER,                  INTENT(IN)    :: KRRI ! Number of ice water var.
-!
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PUT, PVT, PWT   ! variables
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT            !     at
-REAL, DIMENSION(:,:,:,:), INTENT(IN)    :: PRT             !  time t
-!
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCORIOX   !    2D horizontal
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCORIOY   !      Coriolis 
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCORIOZ   !       factors
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCURVX    !    2D horizontal
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCURVY    !  curvature factors
-!
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODJ    ! Density * Jacobian
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PZZ       ! Height (z)
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHVREF   ! Virtual Temperature
-                                          ! of the reference state
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PEXNREF   ! Exner function
-                                          ! of the reference state
-!
-REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRUS, PRVS, PRWS ! Sources of Momentum
-REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRTHS            ! Sources of theta
-!
-END SUBROUTINE DYN_SOURCES
-!
-END INTERFACE
-!
-END MODULE MODI_DYN_SOURCES 
-!     ######################################################################
-      SUBROUTINE DYN_SOURCES( KRR,KRRL, KRRI,                              &
-                              PUT, PVT, PWT, PTHT, PRT,                    &
-                              PCORIOX, PCORIOY, PCORIOZ, PCURVX, PCURVY,   &
-                              PRHODJ, PZZ, PTHVREF, PEXNREF,               &
-                              PRUS, PRVS, PRWS, PRTHS                      )
-!     ######################################################################
-!
-!!****  *DYN_SOURCES * - routine to compute the curvature, coriolis and gravity terms
-!!
-!!    PURPOSE
-!!    -------
-!!      The purpose of this routine is to compute the dynamical sources
-!!    (curvature, coriolis and gravity terms) for each component of the
-!!    momentum. The curvature terms arise in case of non-cartesian geometry
-!!    only. For both curvature and coriolis terms, the "thin shell" 
-!!    approximation can be assumed or not. Gravity affects only the vertical
-!!    component of the momentum.
-!!      This routine also adds the Lagrangien derivative of the pressure as a
-!!    source for the theta evolution.
-!!
-!!     
-!!**  METHOD
-!!    ------
-!!      The horizontal curvature and coriolis field arrays are expanded in
-!!    the vertical direction with the SPREAD intrinsics to match with the
-!!    Shuman operators. The source term for theta due to the advection of the
-!!    Exner function is treated in an advective form (or non-conservative form).
-!!    Values of the calorific capacity of the melting and of the potential
-!!    virtual temperature are recovered taking care of the number of water
-!!    variables. The different sources terms are stored for the budget
-!!    computations.
-!!
-!!    EXTERNAL
-!!    --------
-!!      MXM,MYM,MZM : Shuman functions (mean operators)
-!!      MXF,MYF,MZF : Shuman functions (mean operators)
-!!      GZ_M_W      : projection along the vertical direction of the gradient
-!!                    vector. It acts on a field localized in mass point and
-!!                    the result is localized in w point.
-!!      BUDGET      : Stores the different budget components
-!!
-!!    IMPLICIT ARGUMENTS
-!!    ------------------
-!!      Module MODD_CONF     : contains configuration variables for all models
-!!           LTHINSHELL  = .TRUE.  if the THINSHELL approximation is made
-!!           LCARTESIAN  = .TRUE.  if the CARTESIAN approximation is made
-!!      Module MODD_CST      : 
-!!           XG           gravity acceleration
-!!           XRADIUS      Earth radius
-!!           XRD,XRV      Gas constant for dry air and wator vapor
-!!           XCPD,XCPV    Cp for dry air and wator vapor
-!!           XCL,XCI      C (calorific capacity) for liquid and solid water
-!!      Module MODD_DYN      : contains the parameters for the dynamics
-!!           LCORIO      = .FALSE. if the earth rotation is neglected
-!!    
-!!      Module MODD_BUDGET:
-!!         NBUMOD       : model in which budget is calculated
-!!         CBUTYPE      : type of desired budget
-!!                          'CART' for cartesian box configuration
-!!                          'MASK' for budget zone defined by a mask 
-!!                          'NONE'  ' for no budget
-!!         LBU_RTH      : logical for budget of RTH (potential temperature)
-!!                        .TRUE. = budget of RTH        
-!!                        .FALSE. = no budget of RTH
-!!         LBU_RU       : logical for budget of RU (wind component along x)
-!!                        .TRUE. = budget of RU         
-!!                        .FALSE. = no budget of RU 
-!!         LBU_RV       : logical for budget of RV (wind component along y)
-!!                        .TRUE. = budget of RV         
-!!                        .FALSE. = no budget of RV 
-!!         LBU_RW        : logical for budget of RW (wind component along z)
-!!                        .TRUE. = budget of RW         
-!!                        .FALSE. = no budget of RW 
-!!
-!!    REFERENCE
-!!    ---------
-!!      Book2 of documentation ( routine DYN_SOURCE )
-!!
-!!    AUTHOR
-!!    ------
-!!	J.-P. Pinty      * Laboratoire d'Aerologie*
-!!
-!!    MODIFICATIONS
-!!    -------------
-!!      Original    06/06/94 
-!!	Corrections 06/08/94 (J.-P. Lafore) 
-!!	Corrections 17/10/94 (Stein) For LCORIO
-!!	Corrections 22/12/94 (Stein) add the pressure term in the theta evolution
-!!  Corrections 30/12/94 (J.P. Lafore) bug corrections for the pressure term 
-!!  Corrections 16/03/95 (Stein) remove R from the historical variables and
-!!                                   correction of the pressure term   
-!!  Corrections 14/04/95 (Masson Stein) bugs in the curvatureand Coriolis
-!!                                   terms for PRUS,PRVS, PRWS         
-!!  Corrections 01/04/95 (Ph. Hereil J. Nicolau) add the budget computation
-!!  Corrections 16/10/95 (J. Stein)     change the budget calls 
-!!  Corrections 16/02/96 (J.-P. Pinty)  Introduce L1D switches
-!!  Corrections 19/12/96 (J.-P. Pinty)  Update the CALL BUDGET
-!!  Corrections 03/12/02  (P. Jabouille)  add no thinshell condition
-!!  Correction     06/10  (C.Lac) Exclude L1D for Coriolis term 
-!!  Modification   03/11  (C.Lac) Split the gravity term due to buoyancy
-!!   J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1 
-!  P. Wautelet    02/2020: use the new data structures and subroutines for budgets
-!-------------------------------------------------------------------------------
-!
-!*       0.    DECLARATIONS
-!              ------------
-!
-use modd_budget,      only: lbudget_u, lbudget_v, lbudget_w, lbudget_th, &
-                            NBUDGET_U, NBUDGET_V, NBUDGET_W, NBUDGET_TH, &
-                            tbudgets
-USE MODD_CONF
-USE MODD_CST
-USE MODD_DYN
-USE MODD_DYN_n,      ONLY: LOCEAN
-!
-use mode_budget,     only: Budget_store_init, Budget_store_end
-USE MODE_MPPDB
-
-USE MODI_GRADIENT_M
-USE MODI_SHUMAN
-#ifdef MNH_OPENACC
-USE MODI_SHUMAN_DEVICE
-#endif
-#ifdef MNH_OPENACC
-USE MODE_DEVICE
-USE MODE_MNH_ZWORK,      ONLY: MNH_MEM_GET, MNH_MEM_POSITION_PIN, MNH_MEM_RELEASE
-#endif
-
-
-IMPLICIT NONE
-!  
-!*       0.1   Declarations of dummy arguments :
-!
-INTEGER,                  INTENT(IN)    :: KRR  ! Total number of water var.
-INTEGER,                  INTENT(IN)    :: KRRL ! Number of liquid water var.
-INTEGER,                  INTENT(IN)    :: KRRI ! Number of ice water var.
-!
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PUT, PVT, PWT   ! variables
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHT            !     at
-REAL, DIMENSION(:,:,:,:), INTENT(IN)    :: PRT             !  time t
-!
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCORIOX   !    2D horizontal
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCORIOY   !      Coriolis 
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCORIOZ   !       factors
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCURVX    !    2D horizontal
-REAL, DIMENSION(:,:),     INTENT(IN)    :: PCURVY    !  curvature factors
-!
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRHODJ    ! dry Density * Jacobian
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PZZ       ! Height (z)
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PTHVREF   ! Virtual Temperature
-                                          ! of the reference state
-REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PEXNREF   ! Exner function
-                                          ! of the reference state
-!
-REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRUS, PRVS, PRWS ! Sources of Momentum
-REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRTHS            ! Sources of theta
-!
-!
-!*       0.2   Declarations of local variables :
-!
-REAL       ::  ZCPD_OV_RD   ! = CPD / RD
-REAL       ::  ZG_OV_CPD    ! =-XG / XCPD
-INTEGER    ::  JWATER       ! loop index on the different types of water
-INTEGER    ::  JIU, JJU, JKU     
-INTEGER    ::  JI,JJ,JK     
-INTEGER    ::  IKU          ! array size along the k direction 
-REAL       ::  ZD1          ! DELTA1 (switch 0/1) for thinshell approximation
-#ifndef MNH_OPENACC
-REAL, DIMENSION(SIZE(PUT,1),SIZE(PUT,2),SIZE(PUT,3)) ::           &
-                              ZWORK1, ZWORK2, ZWORK3, ZRUT, ZRVT
-#else
-REAL, DIMENSION(:,:,:), pointer , contiguous ::ZWORK1, ZWORK2, ZWORK3, ZRUT, ZRVT
-#endif
-#ifdef MNH_OPENACC
-REAL, DIMENSION(:,:,:), pointer , contiguous :: ZTMP1_DEVICE,ZTMP2_DEVICE,ZTMP3_DEVICE,ZTMP4_DEVICE,ZRVTT,ZRUTT
-REAL, DIMENSION(:,:,:), pointer , contiguous :: ZTMP5_DEVICE,ZTMP6_DEVICE,ZTMP7_DEVICE,ZTMP8_DEVICE
-!!!REAL, DIMENSION(:,:,:), pointer , contiguous :: ZTMP9_DEVICE,ZTMP10_DEVICE,ZTMP11_DEVICE,ZTMP12_DEVICE
-#endif                      
-!
-!-------------------------------------------------------------------------------
-!
-
-JIU = size( PUT, 1 )
-JJU = size( PUT, 2 )
-JKU = size( PUT, 3 )
-
-#ifdef MNH_OPENACC
-CALL MNH_MEM_POSITION_PIN()
-
-CALL MNH_MEM_GET( ztmp1_device, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ztmp2_device, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ztmp3_device, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ztmp4_device, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ztmp5_device, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ztmp6_device, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ztmp7_device, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ztmp8_device, JIU, JJU, JKU )
-!CALL MNH_MEM_GET( ztmp9_device, JIU, JJU, JKU )
-!CALL MNH_MEM_GET( ztmp10_device, JIU, JJU, JKU )
-!CALL MNH_MEM_GET( ztmp11_device, JIU, JJU, JKU )
-!CALL MNH_MEM_GET( ztmp12_device, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ZRVTT, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ZRUTT, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ZWORK1, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ZWORK2, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ZWORK3, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ZRUT, JIU, JJU, JKU )
-CALL MNH_MEM_GET( ZRVT, JIU, JJU, JKU )
-
-#endif
-
-!
-!*       1.     COMPUTES THE TRUE VELOCITY COMPONENTS
-!	        -------------------------------------
-!
-
-!$acc data copyin             ( PUT, PVT, PWT,PTHT, PRT,                     &
-!$acc &                        PCORIOX, PCORIOY, PCORIOZ, PCURVX, PCURVY,   &
-!$acc &                        PZZ, PTHVREF, PEXNREF)               &       
-!$acc &   present               (PRUS,PRVS, PRWS, PRTHS,PRHODJ)        
-!$acc update device(PRUS, PRVS, PRWS,PRTHS,PRHODJ)
-!$acc update device  ( PUT, PVT, PWT,PTHT, PRT,                     &
-!$acc &                        PCORIOX, PCORIOY, PCORIOZ, PCURVX, PCURVY,   &
-!$acc &                        PZZ, PTHVREF, PEXNREF)
-
-
-#ifndef MNH_OPENACC
-ZRUT(:,:,:) = PUT(:,:,:) * MXM(PRHODJ(:,:,:))
-ZRVT(:,:,:) = PVT(:,:,:) * MYM(PRHODJ(:,:,:))
-
-print *,"PUT(1,1,1)",PUT(1,1,1)
-#else
-CALL MXM_DEVICE(PRHODJ(:,:,:), ZTMP1_DEVICE (:,:,:))
-CALL MYM_DEVICE(PRHODJ(:,:,:), ZTMP2_DEVICE (:,:,:))
-!$acc kernels present(ZRUT,ZRVT,PUT,PVT,ZTMP1_DEVICE,ZTMP2_DEVICE)
-ZRUT(:,:,:) = PUT(:,:,:) * ZTMP1_DEVICE(:,:,:)
-ZRVT(:,:,:) = PVT(:,:,:) * ZTMP2_DEVICE(:,:,:)
-!$acc end kernels
-#endif
-!
-IKU = SIZE(PUT,3)
-!
-!-------------------------------------------------------------------------------
-!
-!*       2.     COMPUTES THE CURVATURE TERMS
-!    	        ----------------------------
-!
-! Only when earth rotation is considered but not in 1D and CARTESIAN cases
-!
-IF ((.NOT.L1D).AND.(.NOT.LCARTESIAN) )  THEN
-  if ( lbudget_u ) call Budget_store_init( tbudgets(NBUDGET_U), 'CURV', prus(:, :, :) )
-  if ( lbudget_v ) call Budget_store_init( tbudgets(NBUDGET_V), 'CURV', prvs(:, :, :) )
-
-  IF ( LTHINSHELL ) THEN           !  THINSHELL approximation
-#ifndef  MNH_OPENACC
-    ZWORK1(:,:,:) = SPREAD( PCURVX(:,:),DIM=3,NCOPIES=JKU ) / XRADIUS
-    ZWORK2(:,:,:) = SPREAD( PCURVY(:,:),DIM=3,NCOPIES=JKU ) / XRADIUS
-#else
- !$acc kernels present_cr(ZWORK1,ZWORK2)
-  !$mnh_expand_array(JI=1:JIU  , JJ=1:JJU , JK=1:JKU)   
-    ZWORK1(:,:,:) = PCURVX(:,:) / XRADIUS
-    ZWORK2(:,:,:) = PCURVY(:,:) / XRADIUS
-  !$mnh_end_expand_array()
-!$acc end kernels
-#endif    
-#ifndef MNH_OPENACC
-    PRUS(:,:,:) = PRUS                                   &
-    + ZRUT * MXM(  MYF(PVT) * ZWORK1 )                   &
-    + MXM( MYF(ZRVT*PVT) * ZWORK2 )
-!
-    PRVS(:,:,:) = PRVS                                   &
-    - MYM( MXF(ZRUT*PUT) * ZWORK1 )                      &
-    - ZRVT * MYM( MXF(PUT) * ZWORK2 ) 
-!
-#else
-!$acc kernels
-ZRVTT(:,:,:)=ZRVT(:,:,:)*PVT(:,:,:)
-ZRUTT(:,:,:)=ZRUT(:,:,:)*PUT(:,:,:)
-!$acc end kernels
-CALL MYF_DEVICE(PVT,ZTMP1_DEVICE )   !MYF(PVT)
-CALL MYF_DEVICE(ZRVTT,ZTMP2_DEVICE ) !MYF(ZRVT*PVT)
-CALL MXF_DEVICE(ZRUTT,ZTMP3_DEVICE ) !MXF(ZRUT*PUT)
-CALL MXF_DEVICE(PUT,ZTMP4_DEVICE )   !MXF(PUT)
-
-!$acc kernels present_cr(ZTMP1_DEVICE,ZTMP1_DEVICE,ZTMP1_DEVICE,ZTMP1_DEVICE)
-ZTMP1_DEVICE(:,:,:)=ZTMP1_DEVICE(:,:,:)* ZWORK1(:,:,:)   !MYF(PVT) * ZWORK1
-ZTMP2_DEVICE(:,:,:)=ZTMP2_DEVICE(:,:,:)* ZWORK2(:,:,:)   !MYF(ZRVT*PVT) * ZWORK2
-ZTMP3_DEVICE(:,:,:)=ZTMP3_DEVICE(:,:,:)* ZWORK1(:,:,:)   !MXF(ZRUT*PUT) * ZWORK1
-ZTMP4_DEVICE(:,:,:)=ZTMP4_DEVICE(:,:,:) * ZWORK2(:,:,:)  !MXF(PUT) * ZWORK2
-!$acc end kernels
-
-CALL MXM_DEVICE(ZTMP1_DEVICE,ZTMP5_DEVICE)   !MXM(  MYF(PVT) * ZWORK1 )
-CALL MXM_DEVICE(ZTMP2_DEVICE,ZTMP6_DEVICE)   !MXM( MYF(ZRVT*PVT) * ZWORK2 )
-
-CALL MYM_DEVICE(ZTMP3_DEVICE,ZTMP7_DEVICE)  !MYM( MXF(ZRUT*PUT) * ZWORK1 )
-CALL MYM_DEVICE(ZTMP4_DEVICE,ZTMP8_DEVICE) !MYM( MXF(PUT) * ZWORK2 )
-
-!$acc kernels present_cr(PRUS,PRVS)
-PRUS(:,:,:) = PRUS(:,:,:)                                  &
-+ ZRUT(:,:,:) * ZTMP5_DEVICE(:,:,:)                        &
-    + ZTMP6_DEVICE(:,:,:)
-
-PRVS(:,:,:) = PRVS(:,:,:)                                   &
-    - ZTMP7_DEVICE(:,:,:)                                  &
-    - ZRVT(:,:,:) * ZTMP8_DEVICE(:,:,:)
-
-!$acc end kernels
-#endif
-
-  ELSE                           !  NO THINSHELL approximation
-    if ( lbudget_w ) call Budget_store_init( tbudgets(NBUDGET_W), 'CURV', prws(:, :, :) )
-#ifndef MNH_OPENACC
- ZWORK3(:,:,:) = 1.0 / ( XRADIUS + MZF(PZZ(:,:,:)) )
- ZWORK1(:,:,:) = SPREAD( PCURVX(:,:),DIM=3,NCOPIES=JKU )
- ZWORK2(:,:,:) = SPREAD( PCURVY(:,:),DIM=3,NCOPIES=JKU )
-#else
-CALL MZF_DEVICE(PZZ,ZWORK3) !MZF(PZZ(:,:,:))
-!$acc kernels present_cr(ZWORK3,ZWORK1,ZWORK2)
-    ZWORK3(:,:,:) = 1.0 / ( XRADIUS + ZWORK3(:,:,:) )
-    ZWORK1(:,:,:) = SPREAD( PCURVX(:,:),DIM=3,NCOPIES=JKU )
-    ZWORK2(:,:,:) = SPREAD( PCURVY(:,:),DIM=3,NCOPIES=JKU )
-!$acc end kernels
-#endif
-    IF (MPPDB_INITIALIZED) THEN
-    CALL MPPDB_CHECK3DM("DYN_SOURCES:ZWORK3,ZWORK1,ZWORK2",PRECISION,&
-                      & ZWORK3,ZWORK1,ZWORK2,&
-                      & MXM( MYF(ZRVT*PVT) * ZWORK2 * ZWORK3 ) , &
-                      & MXM( ( MYF(PVT) * ZWORK1 - MZF(PWT) ) * ZWORK3 ) ,&
-                      & MYF(PVT) * ZWORK1 - MZF(PWT) , &
-                      & MYF(PVT) , MZF(PWT) , MXM(PWT) , MYM(PWT) )
-    CALL MPPDB_CHECK3DM("DYN_SOOURCES:SUITE",PRECISION,&
-         &  MXM(ZRVT),MXM(PVT),MXM(PWT),MXM(ZWORK1),MXM(ZWORK2),MXM(ZWORK3),&
-         &  ZRUT,ZRVT,PRUS,PRVS,PRWS )
-    ENDIF
-!
-#ifndef MNH_OPENACC
-    PRUS(:,:,:) = PRUS                                              &
-    + MXM( MYF(ZRVT*PVT) * ZWORK2 * ZWORK3 )                        &
-    + ZRUT * MXM( ( MYF(PVT) * ZWORK1 - MZF(PWT) ) * ZWORK3 )
-!
-    PRVS(:,:,:) = PRVS                                              &
-    - MYM( MXF(ZRUT*PUT) * ZWORK1 * ZWORK3 )                        &
-    - ZRVT * MYM( (MXF(PUT) * ZWORK2 + MZF(PWT) ) * ZWORK3 )
-!
-    PRWS(:,:,:) = PRWS                                              &
-    +MZM( ( MXF(ZRUT*PUT) + MYF(ZRVT*PVT) ) * ZWORK3 )
-#else
-!$acc kernels present_cr(ZRVTT,ZRUTT)
-ZRVTT(:,:,:)=ZRVT(:,:,:)*PVT(:,:,:)
-ZRUTT(:,:,:)=ZRUT(:,:,:)*PUT(:,:,:)
-
-!PRUS(:,:,:) = PRUS                                              &
-!   + MXM( MYF(ZRVT*PVT) * ZWORK2 * ZWORK3 )                        &
-!   + ZRUT * MXM( ( MYF(PVT) * ZWORK1 - MZF(PWT) ) * ZWORK3 )
-!$acc end kernels
-CALL MYF_DEVICE(ZRVTT,ZTMP1_DEVICE )  !MYF(ZRVT*PVT)
-CALL MYF_DEVICE(PVT,ZTMP2_DEVICE )    !MYF(PVT)
-CALL MZF_DEVICE(PWT,ZTMP3_DEVICE )    !MZF(PWT)
-!$acc kernels present_cr(ZTMP1_DEVICE,ZTMP4_DEVICE)
-ZTMP1_DEVICE(:,:,:)=ZTMP1_DEVICE(:,:,:)* ZWORK2(:,:,:) * ZWORK3(:,:,:)  ! MYF(ZRVT*PVT) * ZWORK2 * ZWORK3
-ZTMP4_DEVICE(:,:,:)=(ZTMP2_DEVICE(:,:,:)* ZWORK1(:,:,:) - ZTMP3_DEVICE(:,:,:)) * ZWORK3(:,:,:)  !( MYF(PVT) * ZWORK1 - MZF(PWT) ) * ZWORK3
-!$acc end kernels
-CALL MXM_DEVICE(ZTMP1_DEVICE,ZTMP2_DEVICE) ! MXM( MYF(ZRVT*PVT) * ZWORK2 * ZWORK3 ) 
-CALL MXM_DEVICE(ZTMP4_DEVICE,ZTMP5_DEVICE) ! MXM( ( MYF(PVT) * ZWORK1 - MZF(PWT) ) * ZWORK3 )
-!$acc kernels present_cr(PRUS)
-    PRUS(:,:,:) = PRUS(:,:,:)                                &
-    +ZTMP2_DEVICE(:,:,:)                                     &
-    + ZRUT(:,:,:) * ZTMP5_DEVICE(:,:,:)
-!$acc end kernels
-
-!    PRVS(:,:,:) = PRVS                                              &
-!    - MYM( MXF(ZRUT*PUT) * ZWORK1 * ZWORK3 )                        &
-!    - ZRVT * MYM( (MXF(PUT) * ZWORK2 + MZF(PWT) ) * ZWORK3 )
-CALL MXF_DEVICE(ZRUTT,ZTMP1_DEVICE ) !MXF(ZRUT*PUT)
-CALL MXF_DEVICE(PUT,ZTMP2_DEVICE )   !MXF(PUT)
-!CALL MZF_DEVICE(PWT,ZTMP3_DEVICE )   !MZF(PWT)
-!$acc kernels present_cr(ZTMP1_DEVICE,ZTMP4_DEVICE)
-ZTMP1_DEVICE(:,:,:)=ZTMP1_DEVICE(:,:,:)* ZWORK1(:,:,:) * ZWORK3(:,:,:) !MXF(ZRUT*PUT) * ZWORK1 * ZWORK3
-ZTMP4_DEVICE(:,:,:)=(ZTMP2_DEVICE(:,:,:)* ZWORK2(:,:,:) + ZTMP3_DEVICE(:,:,:)) * ZWORK3(:,:,:) !(MXF(PUT) * ZWORK2 + MZF(PWT) ) * ZWORK3
-!$acc end kernels
-CALL MYM_DEVICE(ZTMP1_DEVICE,ZTMP2_DEVICE) !MYM( MXF(ZRUT*PUT) * ZWORK1 * ZWORK3 )
-CALL MYM_DEVICE(ZTMP4_DEVICE,ZTMP5_DEVICE) !MYM( (MXF(PUT) * ZWORK2 + MZF(PWT) ) * ZWORK3 )
-!$acc kernels present_cr(PRVS)
-    PRVS(:,:,:) = PRVS(:,:,:)                                &
-    - ZTMP2_DEVICE(:,:,:)                                     &
-    - ZRVT(:,:,:) * ZTMP5_DEVICE(:,:,:)
-!$acc end kernels
-
-!    PRWS(:,:,:) = PRWS                                              &
-!    +MZM( ( MXF(ZRUT*PUT) + MYF(ZRVT*PVT) ) * ZWORK3 )
-CALL MYF_DEVICE(ZRVTT,ZTMP1_DEVICE ) !MXF(ZRUT*PUT)
-CALL MXF_DEVICE(ZRUTT,ZTMP2_DEVICE ) !MYF(ZRVT*PVT)
-!$acc kernels present_cr(ZTMP3_DEVICE)
-ZTMP3_DEVICE(:,:,:)=(ZTMP2_DEVICE(:,:,:)+ZTMP1_DEVICE(:,:,:)) * ZWORK3(:,:,:) ! (MXF(ZRUT*PUT) + MYF(ZRVT*PVT) ) * ZWORK3
-!$acc end kernels
-CALL MZM_DEVICE(ZTMP3_DEVICE,ZTMP1_DEVICE) ! MZM( ( MXF(ZRUT*PUT) + MYF(ZRVT*PVT) ) * ZWORK3 )
-!$acc kernels present_cr(PRWS)
-PRWS(:,:,:) = PRWS(:,:,:)+ZTMP1_DEVICE(:,:,:)
-!$acc end kernels
-#endif
-
-
-    if ( lbudget_w ) call Budget_store_end( tbudgets(NBUDGET_W), 'CURV', prws(:, :, :) )
-  END IF
-
-  if ( lbudget_u ) call Budget_store_end( tbudgets(NBUDGET_U), 'CURV', prus(:, :, :) )
-  if ( lbudget_v ) call Budget_store_end( tbudgets(NBUDGET_V), 'CURV', prvs(:, :, :) )
-END IF
-!
-!
-!-------------------------------------------------------------------------------
-!
-!*       3.     COMPUTES THE CORIOLIS TERMS
-!	        ---------------------------
-!
-IF (LCORIO)   THEN 
-  if ( lbudget_u ) call Budget_store_init( tbudgets(NBUDGET_U), 'COR', prus(:, :, :) )
-  if ( lbudget_v ) call Budget_store_init( tbudgets(NBUDGET_V), 'COR', prvs(:, :, :) )
-!$acc kernels present_cr(ZWORK3)
-  ZWORK3(:,:,:) = SPREAD( PCORIOZ(:,:),DIM=3,NCOPIES=JKU ) * PRHODJ(:,:,:)
-!$acc end kernels
-#ifndef MNH_OPENACC
-PRUS(:,:,:) = PRUS + MXM( ZWORK3 * MYF(PVT) )
-PRVS(:,:,:) = PRVS - MYM( ZWORK3 * MXF(PUT) )
-#else
-CALL MYF_DEVICE(PVT,ZTMP1_DEVICE) !MYF(PVT)
-CALL MXF_DEVICE(PUT,ZTMP2_DEVICE) !MXF(PUT)
-!$acc kernels present_cr(ZTMP1_DEVICE,ZTMP2_DEVICE)
-ZTMP1_DEVICE(:,:,:)=ZWORK3(:,:,:) * ZTMP1_DEVICE(:,:,:)
-ZTMP2_DEVICE(:,:,:)=ZWORK3(:,:,:) * ZTMP2_DEVICE(:,:,:)
-!$acc end kernels
-CALL MXM_DEVICE(ZTMP1_DEVICE, ZTMP3_DEVICE) ! MXM( ZWORK3 * MYF(PVT)
-CALL MYM_DEVICE(ZTMP2_DEVICE, ZTMP4_DEVICE) ! MYM( ZWORK3 * MXF(PUT)
-!$acc kernels present_cr(PRUS,PRVS)
-PRUS(:,:,:) = PRUS(:,:,:) + ZTMP3_DEVICE(:,:,:)
-PRVS(:,:,:) = PRVS(:,:,:) - ZTMP4_DEVICE(:,:,:)
-!$acc end kernels
-#endif  
-!
-
-IF ((.NOT.LTHINSHELL) .AND. (.NOT.L1D)) THEN
-if ( lbudget_w ) call Budget_store_init( tbudgets(NBUDGET_W), 'COR', prws(:, :, :) )
-!$acc kernels present_cr(ZWORK1,ZWORK2)
-ZWORK1(:,:,:) = SPREAD( PCORIOX(:,:),DIM=3,NCOPIES=JKU) * PRHODJ(:,:,:) 
-ZWORK2(:,:,:) = SPREAD( PCORIOY(:,:),DIM=3,NCOPIES=JKU) * PRHODJ(:,:,:)
-!$acc end kernels
-
-CALL MPPDB_CHECK3DM("DYN_SOOURCES:CORIOLIS",PRECISION,&
- & ZWORK1,ZWORK2,ZWORK3 )
-!
-#ifndef MNH_OPENACC
-    PRUS(:,:,:) = PRUS - MXM( ZWORK2 * MZF(PWT) )
-!
-    PRVS(:,:,:) = PRVS - MYM( ZWORK1 * MZF(PWT) )
-!
-    PRWS(:,:,:) = PRWS + MZM( ZWORK2 * MXF(PUT) + ZWORK1 * MYF(PVT) )
-#else
-!***PRUS(:,:,:) = PRUS - MXM( ZWORK2 * MZF(PWT) )
-CALL MZF_DEVICE(PWT,ZTMP1_DEVICE) !MZF(PWT)
-!$acc kernels present_cr(ZTMP2_DEVICE)
-ZTMP2_DEVICE(:,:,:)=ZWORK2(:,:,:) * ZTMP1_DEVICE(:,:,:)    !ZWORK2 * MZF(PWT)
-!$acc end kernels
-CALL MXM_DEVICE(ZTMP2_DEVICE,ZTMP3_DEVICE) !MXM( ZWORK2 * MZF(PWT) )
-!$acc kernels present_cr(PRUS)
-    PRUS(:,:,:) = PRUS(:,:,:) - ZTMP3_DEVICE(:,:,:)
-!$acc end kernels
-
-!*** PRVS(:,:,:) = PRVS - MYM( ZWORK1 * MZF(PWT) )
-!$acc kernels present_cr(ZTMP2_DEVICE)
-ZTMP2_DEVICE(:,:,:)=ZWORK1(:,:,:) * ZTMP1_DEVICE(:,:,:)    !ZWORK1 * MZF(PWT)
-!$acc end kernels
-CALL MYM_DEVICE(ZTMP2_DEVICE,ZTMP3_DEVICE) !MYM( ZWORK1 * MZF(PWT) )
-!$acc kernels present_cr(PRVS)
-    PRVS(:,:,:) = PRVS(:,:,:) - ZTMP2_DEVICE(:,:,:)
-!$acc end kernels
-
-!*** PRWS(:,:,:) = PRWS + MZM( ZWORK2 * MXF(PUT) + ZWORK1 * MYF(PVT) )
-CALL MXF_DEVICE(PUT,ZTMP1_DEVICE) !MXF(PUT)
-CALL MYF_DEVICE(PVT,ZTMP2_DEVICE) !MYF(PVT)
-!$acc kernels present_cr(ZTMP3_DEVICE)
-ZTMP3_DEVICE(:,:,:)= ZWORK2(:,:,:) * ZTMP1_DEVICE(:,:,:)+ ZWORK1(:,:,:) * ZTMP2_DEVICE(:,:,:) !ZWORK2 * MXF(PUT) + ZWORK1 * MYF(PVT)
-!$acc end kernels
-CALL MZM_DEVICE(ZTMP3_DEVICE,ZTMP1_DEVICE) !MZM( ZWORK2 * MXF(PUT) + ZWORK1 * MYF(PVT) )
-!$acc kernels present_cr(PRWS)
-    PRWS(:,:,:) = PRWS(:,:,:) + ZTMP1_DEVICE(:,:,:)
-!$acc end kernels
-#endif
-
-    if ( lbudget_w ) call Budget_store_end( tbudgets(NBUDGET_W), 'COR', prws(:, :, :) )
-  END IF
-
-  if ( lbudget_u ) call Budget_store_end( tbudgets(NBUDGET_U), 'COR', prus(:, :, :) )
-  if ( lbudget_v ) call Budget_store_end( tbudgets(NBUDGET_V), 'COR', prvs(:, :, :) )
-END IF
-!!!!$acc end data
-!
-!-------------------------------------------------------------------------------
-!
-!*       4.     COMPUTES THE THETA SOURCE TERM DUE TO THE REFERENCE PRESSURE
-!	        ------------------------------------------------------------
-!
-IF (LCARTESIAN .OR. LTHINSHELL) THEN
-  ZD1=0.
-ELSE
-  ZD1=1.
-ENDIF
-!
-IF( .NOT.L1D ) THEN
- IF (.NOT. LOCEAN) THEN
-!
-  IF (KRR > 0) THEN
-    if ( lbudget_th ) call Budget_store_init( tbudgets(NBUDGET_TH), 'PREF', prths(:, :, :) )
-
-
-    ZCPD_OV_RD = XCPD / XRD
-    ZG_OV_CPD  = -XG / XCPD
-!
-! stores the specific heat capacity (Cph) in ZWORK1
-
-!$acc kernels present_cr( ZWORK1)
-    ZWORK1(:,:,:) = XCPD + XCPV * PRT(:,:,:,1)  ! gas mixing
-!$acc loop seq
-    DO JWATER = 2,1+KRRL             !  loop on the liquid components  
-      ZWORK1(:,:,:) = ZWORK1(:,:,:) + XCL * PRT(:,:,:,JWATER)
-    END DO
-!$acc loop seq
-    DO JWATER = 2+KRRL,1+KRRL+KRRI   !  loop on the solid components   
-      ZWORK1(:,:,:) = ZWORK1(:,:,:) + XCI * PRT(:,:,:,JWATER)
-    END DO
-!$acc end kernels
-! computes the source term
-!
-#ifndef MNH_OPENACC
-    PRTHS(:,:,:) = PRTHS(:,:,:) +  PRHODJ(:,:,:)                             &
-      * ( ( XRD + XRV * PRT(:,:,:,1) ) * ZCPD_OV_RD / ZWORK1(:,:,:)  - 1. )  &
-      * PTHT(:,:,:)/PEXNREF(:,:,:)*MZF(PWT(:,:,:))*(ZG_OV_CPD/PTHVREF(:,:,:) &
-      -ZD1*4./7.*PEXNREF(:,:,:)/( XRADIUS+MZF(PZZ(:,:,:)) ))
-#else
-CALL MZF_DEVICE(PWT,ZTMP1_DEVICE) !MZF(PWT(:,:,:))
-CALL MZF_DEVICE(PZZ(:,:,:),ZTMP2_DEVICE(:,:,:))  !MZF(PZZ(:,:,:))
-!$acc kernels present_cr(PRTHS)
-PRTHS(:,:,:) = PRTHS(:,:,:) +  PRHODJ(:,:,:)                             &
-      * ( ( XRD + XRV * PRT(:,:,:,1) ) * ZCPD_OV_RD / ZWORK1(:,:,:)  - 1. )  &
-      * PTHT(:,:,:)/PEXNREF(:,:,:)*ZTMP1_DEVICE*(ZG_OV_CPD/PTHVREF(:,:,:) &
-      -ZD1*4./7.*PEXNREF(:,:,:)/( XRADIUS+ZTMP2_DEVICE(:,:,:) ))
-!$acc end kernels
-#endif
-    if ( lbudget_th ) call Budget_store_end( tbudgets(NBUDGET_TH), 'PREF', prths(:, :, :) )
-  END IF
- END IF
-!
-END IF
-!
-!-------------------------------------------------------------------------------
-#ifdef MNH_OPENACC
-CALL MNH_MEM_RELEASE()
-#endif
-!$acc end data
-!$acc update self(PRUS, PRVS, PRWS,PRTHS)
-END SUBROUTINE DYN_SOURCES 
diff --git a/src/ZSOLVER/update_lm.f90 b/src/ZSOLVER/update_lm.f90
deleted file mode 100644
index c8908734a..000000000
--- a/src/ZSOLVER/update_lm.f90
+++ /dev/null
@@ -1,185 +0,0 @@
-!MNH_LIC Copyright 2006-2022 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.
-!-----------------------------------------------------------------
-!     ###################
-      MODULE MODI_UPDATE_LM
-!     ###################
-INTERFACE
-!
-SUBROUTINE UPDATE_LM(HLBCX,HLBCY,PLM,PLEPS)
-!
-CHARACTER(LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX ! X boundary type
-CHARACTER(LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY ! Y boundary type
-!
-REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PLM   ! mixing length
-REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PLEPS ! dissipative length
-!
-END SUBROUTINE UPDATE_LM
-!
-END INTERFACE
-!
-END MODULE MODI_UPDATE_LM
-!
-!
-!
-!     #################################################################
-      SUBROUTINE UPDATE_LM(HLBCX,HLBCY,PLM,PLEPS)
-!     #################################################################
-!
-!!****  *UPDATE_LM* - routine to set external points for mixing length
-!!
-!!    PURPOSE
-!!    -------
-!
-!!**  METHOD
-!!    ------
-!!   
-!!    EXTERNAL
-!!    --------   
-!!
-!!    IMPLICIT ARGUMENTS
-!!    ------------------ 
-!!        
-!!    REFERENCE
-!!    ---------
-!!      Book2 of documentation (routine UPDATE_LM)
-!!
-!!    AUTHOR
-!!    ------
-!!	V. Masson        * Meteo France *
-!!
-!!    MODIFICATIONS
-!!    -------------
-!!      Original    april 2006
-!!       V.Masson : Exchange of East and North sides
-!!   J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1 
-!  P. Wautelet 20/05/2019: add name argument to ADDnFIELD_ll + new ADD4DFIELD_ll subroutine
-!-------------------------------------------------------------------------------
-!
-!*       0.    DECLARATIONS
-!         
-USE MODD_ARGSLIST_ll, ONLY : LIST_ll
-USE MODD_CONF
-USE MODD_PARAMETERS
-!
-USE MODE_ll
-use mode_mppdb
-!
-#ifdef MNH_OPENACC
-USE MODI_GET_HALO
-#endif
-!
-IMPLICIT NONE
-!
-!
-!*       0.1   declarations of arguments
-!
-CHARACTER(LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX ! X boundary type
-CHARACTER(LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY ! Y boundary type
-!
-REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PLM   ! mixing length
-REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PLEPS ! dissipative length
-!
-!*       0.2   declarations of local variables
-!
-INTEGER             :: IIB      ! First physical index in x direction
-INTEGER             :: IJB      ! First physical index in y direction
-INTEGER             :: IIE      ! last  physical index in x direction
-INTEGER             :: IJE      ! last  physical index in y direction
-INTEGER             :: JI       ! loop index
-!
-LOGICAL                :: GNORTH, GSOUTH, GWEST, GEAST
-TYPE(LIST_ll), POINTER :: TZLM_ll   ! list of fields to exchange
-INTEGER                :: IINFO_ll       ! return code of parallel routine
-!
-!-------------------------------------------------------------------------------
-
-!$acc data present(PLM,PLEPS)
-
-if ( mppdb_initialized ) then
-  !Check all inout arrays
-  call Mppdb_check( plm,   "Update_lm beg:plm"   )
-  call Mppdb_check( pleps, "Update_lm beg:pleps" )
-end if
-!
-!*       1.    COMPUTE DIMENSIONS OF ARRAYS :
-!              ----------------------------
-CALL GET_INDICE_ll (IIB,IJB,IIE,IJE)
-#ifndef MNH_OPENACC
-NULLIFY(TZLM_ll)
-#endif
-!
-GWEST  = ( HLBCX(1) /= 'CYCL' .AND. LWEST_ll() )
-GEAST  = ( HLBCX(2) /= 'CYCL' .AND. LEAST_ll() )
-GSOUTH = ( HLBCY(1) /= 'CYCL' .AND. LSOUTH_ll() )
-GNORTH = ( HLBCY(2) /= 'CYCL' .AND. LNORTH_ll() )
-
-!-------------------------------------------------------------------------------
-!
-!*       2.  UPDATE HALOs :
-!            -------------
-!
-!
-!!$IF(NHALO == 1) THEN
-#ifndef MNH_OPENACC
-  CALL ADD3DFIELD_ll( TZLM_ll, PLM,   'UPDATE_LM::PLM'   )
-  CALL ADD3DFIELD_ll( TZLM_ll, PLEPS, 'UPDATE_LM::PLEPS' )
-  CALL UPDATE_HALO_ll(TZLM_ll,IINFO_ll)
-  CALL CLEANLIST_ll(TZLM_ll)
-#else
-!
-! /!\ Corner update needed ! -> GET_HALO_DDC
-!
-CALL GET_HALO_DDC( PLM,   HNAME='UPDATE_LM::PLM'   )
-CALL GET_HALO_DDC( PLEPS, HNAME='UPDATE_LM::PLEPS' )
-#endif
-!!$END IF
-!
-!-------------------------------------------------------------------------------
-!
-!*       3.  UPDATE EXTERNAL POINTS OF GLOBAL DOMAIN:
-!            ---------------------------------------
-!
-
-IF ( GWEST ) THEN
-  !$acc kernels async
-  PLM  (IIB-1,:,:) = PLM  (IIB,:,:)
-  PLEPS(IIB-1,:,:) = PLEPS(IIB,:,:)
-  !$acc end kernels
-END IF
-IF ( GEAST ) THEN
-  !$acc kernels async
-  PLM  (IIE+1,:,:) = PLM  (IIE,:,:)
-  PLEPS(IIE+1,:,:) = PLEPS(IIE,:,:)
-  !$acc end kernels
-END IF
-IF ( GSOUTH ) THEN
-  !$acc kernels async 
-  DO JI=1,SIZE(PLM,1)
-    PLM  (JI,IJB-1,:) = PLM  (JI,IJB,:)
-    PLEPS(JI,IJB-1,:) = PLEPS(JI,IJB,:)
-  END DO
-  !$acc end kernels
-END IF
-IF ( GNORTH ) THEN
-  !$acc kernels async
-  DO JI=1,SIZE(PLM,1)
-    PLM  (JI,IJE+1,:) = PLM  (JI,IJE,:)
-    PLEPS(JI,IJE+1,:) = PLEPS(JI,IJE,:)
- END DO
- !$acc end kernels
-END IF
-!$acc wait
-
-if ( mppdb_initialized ) then
-  !Check all inout arrays
-  call Mppdb_check( plm,   "Update_lm end:plm"   )
-  call Mppdb_check( pleps, "Update_lm end:pleps" )
-end if
-
-!$acc end data
-
-!-----------------------------------------------------------------------------
-END SUBROUTINE UPDATE_LM
-- 
GitLab