diff --git a/src/ZSOLVER/advecuvw_rk.f90 b/src/ZSOLVER/advecuvw_rk.f90
new file mode 100644
index 0000000000000000000000000000000000000000..4a71554b59ef8b3e814b61ca58240aff0b902409
--- /dev/null
+++ b/src/ZSOLVER/advecuvw_rk.f90
@@ -0,0 +1,562 @@
+!MNH_LIC Copyright 1994-2019 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_ADVECUVW_RK
+! #####################
+!
+INTERFACE
+ SUBROUTINE ADVECUVW_RK (HUVW_ADV_SCHEME, &
+ HTEMP_SCHEME, KWENO_ORDER, &
+ HLBCX, HLBCY, PTSTEP, &
+ PU, PV, PW, &
+ PUT, PVT, PWT, &
+ PMXM_RHODJ, PMYM_RHODJ, PMZM_RHODJ, &
+ PRUCT, PRVCT, PRWCT, &
+ PRUS_ADV, PRVS_ADV, PRWS_ADV, &
+ PRUS_OTHER, PRVS_OTHER, PRWS_OTHER &
+#ifndef MNH_OPENACC
+ )
+#else
+ ,ZUT, ZVT, ZWT, ZRUS, ZRVS, ZRWS)
+#endif
+!
+CHARACTER(LEN=6), INTENT(IN) :: HUVW_ADV_SCHEME! to the selected
+CHARACTER(LEN=4), INTENT(IN) :: HTEMP_SCHEME ! Temporal scheme
+!
+INTEGER, INTENT(IN) :: KWENO_ORDER ! Order of the WENO
+ ! scheme (3 or 5)
+!
+CHARACTER(LEN=4),DIMENSION(2),INTENT(IN):: HLBCX, HLBCY ! X- and Y-direc LBC
+!
+REAL, INTENT(IN) :: PTSTEP
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PU , PV , PW
+ ! Variables to advect
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT, PVT , PWT
+ ! Variables at t
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PMXM_RHODJ
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PMYM_RHODJ
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PMZM_RHODJ
+ ! metric coefficients
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRUCT , PRVCT, PRWCT
+ ! Contravariant wind components
+REAL, DIMENSION(:,:,:), INTENT(OUT) :: PRUS_ADV , PRVS_ADV, PRWS_ADV
+ ! Tendency due to advection
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRUS_OTHER , PRVS_OTHER, PRWS_OTHER
+! ! tendencies from other processes
+#ifdef MNH_OPENACC
+! Work arrays
+REAL, DIMENSION(:,:,:) :: ZUT, ZVT, ZWT
+REAL, DIMENSION(:,:,:,:) :: ZRUS,ZRVS,ZRWS
+#endif
+!
+!
+END SUBROUTINE ADVECUVW_RK
+!
+END INTERFACE
+!
+END MODULE MODI_ADVECUVW_RK
+! ##########################################################################
+ SUBROUTINE ADVECUVW_RK (HUVW_ADV_SCHEME, &
+ HTEMP_SCHEME, KWENO_ORDER, &
+ HLBCX, HLBCY, PTSTEP, &
+ PU, PV, PW, &
+ PUT, PVT, PWT, &
+ PMXM_RHODJ, PMYM_RHODJ, PMZM_RHODJ, &
+ PRUCT, PRVCT, PRWCT, &
+ PRUS_ADV, PRVS_ADV, PRWS_ADV, &
+ PRUS_OTHER, PRVS_OTHER, PRWS_OTHER &
+#ifndef MNH_OPENACC
+ )
+#else
+ ,ZUT, ZVT, ZWT, ZRUS, ZRVS, ZRWS)
+#endif
+! ##########################################################################
+!
+!!**** *ADVECUVW_RK * - routine to call the specialized advection routines for wind
+!!
+!! PURPOSE
+!! -------
+!!
+!!** METHOD
+!! ------
+!!
+!! EXTERNAL
+!! --------
+!!
+!! IMPLICIT ARGUMENTS
+!! ------------------
+!! NONE
+!!
+!! REFERENCE
+!! ---------
+!! Book1 and book2 ( routine ADVECTION )
+!!
+!! AUTHOR
+!! ------
+!! J.-P. Pinty * Laboratoire d'Aerologie*
+!! J.-P. Lafore * Meteo France *
+!!
+!! MODIFICATIONS
+!! -------------
+!! Original 06/07/94
+!! 01/04/95 (Ph. Hereil J. Nicolau) add the model number
+!! 23/10/95 (J. Vila and JP Lafore) advection schemes scalar
+!! 16/01/97 (JP Pinty) change presentation
+!! 30/04/98 (J. Stein P Jabouille) extrapolation for the cyclic
+!! case and parallelisation
+!! 24/06/99 (P Jabouille) case of NHALO>1
+!! 25/10/05 (JP Pinty) 4th order scheme
+!! 24/04/06 (C.Lac) Split scalar and passive
+!! tracer routines
+!! 08/06 (T.Maric) PPM scheme
+!! 04/2011 (V. Masson & C. Lac) splits the routine and adds
+!! time splitting
+!! J.Escobar 21/03/2013: for HALOK comment all NHALO=1 test
+!! J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1
+!! F.Auguste and C.Lac : 08/16 : CEN4TH with RKC4
+!! C.Lac 10/16 : Correction on RK loop
+! P. Wautelet 10/04/2019: replace ABORT and STOP calls by Print_msg
+! P. Wautelet 20/05/2019: add name argument to ADDnFIELD_ll + new ADD4DFIELD_ll subroutine
+!!
+!-------------------------------------------------------------------------------
+!
+!* 0. DECLARATIONS
+! ------------
+!
+USE MODD_ARGSLIST_ll, ONLY: LIST_ll, HALO2LIST_ll
+USE MODD_CONF, ONLY: NHALO
+USE MODD_PARAMETERS, ONLY: JPVEXT
+!
+USE MODE_ll
+USE MODE_MPPDB
+use mode_msg
+!
+USE MODI_ADV_BOUNDARIES
+USE MODI_ADVECUVW_4TH
+USE MODI_ADVECUVW_WENO_K
+USE MODI_GET_HALO
+USE MODI_SHUMAN
+!
+!
+#ifdef MNH_OPENACC
+USE MODE_DEVICE
+USE MODE_MNH_ZWORK, ONLY : ZT3D, MNH_GET_ZT3D , MNH_REL_ZT3D
+#endif
+!
+!-------------------------------------------------------------------------------
+!
+IMPLICIT NONE
+!
+!* 0.1 Declarations of dummy arguments :
+!
+CHARACTER(LEN=6), INTENT(IN) :: HUVW_ADV_SCHEME! to the selected
+CHARACTER(LEN=4), INTENT(IN) :: HTEMP_SCHEME ! Temporal scheme
+!
+INTEGER, INTENT(IN) :: KWENO_ORDER ! Order of the WENO
+ ! scheme (3 or 5)
+!
+CHARACTER(LEN=4),DIMENSION(2),INTENT(IN):: HLBCX, HLBCY ! X- and Y-direc LBC
+!
+REAL, INTENT(IN) :: PTSTEP
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PU , PV , PW
+ ! Variables to advect
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT, PVT , PWT
+ ! Variables at t
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PMXM_RHODJ
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PMYM_RHODJ
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PMZM_RHODJ
+ ! metric coefficients
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRUCT , PRVCT, PRWCT
+ ! Contravariant wind components
+REAL, DIMENSION(:,:,:), INTENT(OUT) :: PRUS_ADV , PRVS_ADV, PRWS_ADV
+ ! Tendency due to advection
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRUS_OTHER , PRVS_OTHER, PRWS_OTHER
+! ! tendencies from other processes
+#ifdef MNH_OPENACC
+REAL, DIMENSION(:,:,:) :: ZUT, ZVT, ZWT
+REAL, DIMENSION(:,:,:,:) :: ZRUS,ZRVS,ZRWS
+#endif
+!
+!
+!
+!* 0.2 declarations of local variables
+!
+!
+!
+character(len=3) :: ynum
+INTEGER :: IKE ! indice K End in z direction
+!
+#ifndef MNH_OPENACC
+REAL, DIMENSION(:,:,:), ALLOCATABLE :: ZUT, ZVT, ZWT
+! Intermediate Guesses inside the RK loop
+!
+REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: ZRUS,ZRVS,ZRWS
+#endif
+! Momentum tendencies due to advection
+REAL, DIMENSION(:,:), ALLOCATABLE :: ZBUT ! Butcher array coefficients
+ ! at the RK sub time step
+REAL, DIMENSION(:), ALLOCATABLE :: ZBUTS! Butcher array coefficients
+ ! at the end of the RK loop
+!JUAN
+TYPE(LIST_ll), POINTER :: TZFIELDMT_ll ! list of fields to exchange
+TYPE(HALO2LIST_ll), POINTER :: TZHALO2MT_ll ! momentum variables
+INTEGER :: INBVAR
+INTEGER :: IIU, IJU, IKU ! array sizes
+!JUAN
+
+#ifdef MNH_OPENACC
+INTEGER :: IZMEAN, IZWORK
+#endif
+! Momentum tendencies due to advection
+INTEGER :: ISPL ! Number of RK splitting loops
+INTEGER :: JI, JS ! Loop index
+!
+INTEGER :: IINFO_ll ! return code of parallel routine
+TYPE(LIST_ll), POINTER :: TZFIELD_ll ! list of fields to exchange
+TYPE(LIST_ll), POINTER :: TZFIELDS_ll ! list of fields to exchange
+TYPE(LIST_ll), POINTER :: TZFIELDS0_ll ! list of fields to exchange
+TYPE(LIST_ll), POINTER :: TZFIELDS4_ll ! list of fields to exchange
+!
+!-------------------------------------------------------------------------------
+!$acc data present( PU, PV, PW, PUT, PVT, PWT, PMXM_RHODJ, PMYM_RHODJ, PMZM_RHODJ, &
+!$acc & PRUCT, PRVCT, PRWCT, PRUS_ADV, PRVS_ADV, PRWS_ADV, &
+!$acc & PRUS_OTHER, PRVS_OTHER, PRWS_OTHER, ZUT, ZVT, ZWT, ZRUS, ZRVS, ZRWS )
+
+IF (MPPDB_INITIALIZED) THEN
+ !Check all IN arrays
+ CALL MPPDB_CHECK(PU,"ADVECUVW_RK beg:PU")
+ CALL MPPDB_CHECK(PV,"ADVECUVW_RK beg:PV")
+ CALL MPPDB_CHECK(PW,"ADVECUVW_RK beg:PW")
+ CALL MPPDB_CHECK(PUT,"ADVECUVW_RK beg:PUT")
+ CALL MPPDB_CHECK(PVT,"ADVECUVW_RK beg:PVT")
+ CALL MPPDB_CHECK(PWT,"ADVECUVW_RK beg:PWT")
+ CALL MPPDB_CHECK(PMXM_RHODJ,"ADVECUVW_RK beg:PMXM_RHODJ")
+ CALL MPPDB_CHECK(PMYM_RHODJ,"ADVECUVW_RK beg:PMYM_RHODJ")
+ CALL MPPDB_CHECK(PMZM_RHODJ,"ADVECUVW_RK beg:PMZM_RHODJ")
+ CALL MPPDB_CHECK(PRUCT,"ADVECUVW_RK beg:PRUCT")
+ CALL MPPDB_CHECK(PRVCT,"ADVECUVW_RK beg:PRVCT")
+ CALL MPPDB_CHECK(PRWCT,"ADVECUVW_RK beg:PRWCT")
+ CALL MPPDB_CHECK(PRUS_OTHER,"ADVECUVW_RK beg:PRUS_OTHER")
+ CALL MPPDB_CHECK(PRVS_OTHER,"ADVECUVW_RK beg:PRVS_OTHER")
+ CALL MPPDB_CHECK(PRWS_OTHER,"ADVECUVW_RK beg:PRWS_OTHER")
+END IF
+!
+!* 0. INITIALIZATION
+! --------------
+!
+#ifndef MNH_OPENACC
+allocate(ZUT(SIZE(PUT,1),SIZE(PUT,2),SIZE(PUT,3)))
+allocate(ZVT(SIZE(PUT,1),SIZE(PUT,2),SIZE(PUT,3)))
+allocate(ZWT(SIZE(PUT,1),SIZE(PUT,2),SIZE(PUT,3)))
+#endif
+
+#ifdef MNH_OPENACC
+CALL INIT_ON_HOST_AND_DEVICE(ZUT,4e99,'ADVECUVW_RK::ZUT')
+CALL INIT_ON_HOST_AND_DEVICE(ZVT,5e99,'ADVECUVW_RK::ZVT')
+CALL INIT_ON_HOST_AND_DEVICE(ZWT,6e99,'ADVECUVW_RK::ZWT')
+!
+CALL MNH_GET_ZT3D(IZMEAN,IZWORK)
+#endif
+!
+IKE = SIZE(PWT,3) - JPVEXT
+IIU=SIZE(PUT,1)
+IJU=SIZE(PUT,2)
+IKU=SIZE(PUT,3)
+!
+SELECT CASE (HTEMP_SCHEME)
+ CASE('RK11')
+ ISPL = 1
+ CASE('RK21')
+ ISPL = 2
+ CASE('NP32')
+ ISPL = 3
+ CASE('SP32')
+ ISPL = 3
+ CASE('RK33')
+ ISPL = 3
+ CASE('RKC4')
+ ISPL = 4
+ CASE('RK4B')
+ ISPL = 4
+ CASE('RK53')
+ ISPL = 5
+ CASE('RK62')
+ ISPL = 6
+ CASE('RK65')
+ ISPL = 6
+ CASE DEFAULT
+ call Print_msg(NVERB_FATAL,'GEN','ADVECUVW_RK','unknown HTEMP_SCHEME')
+END SELECT
+!
+!
+ALLOCATE(ZBUT(ISPL-1,ISPL-1))
+ALLOCATE(ZBUTS(ISPL))
+
+!$acc data create(ZBUT,ZBUTS)
+
+SELECT CASE (HTEMP_SCHEME)
+ CASE('RK11')
+ ZBUTS = (/ 1. /)
+ CASE('RK21')
+ ZBUTS = (/ 0. , 1. /)
+ ZBUT(1,1) = 3./4.
+ CASE('RK33')
+ ZBUTS = (/ 1./6. , 1./6. , 2./3. /)
+ ZBUT(1,1) = 1.
+ ZBUT(1,2) = 0.
+ ZBUT(2,1) = 1./4.
+ ZBUT(2,2) = 1./4.
+ CASE('NP32')
+ ZBUTS = (/ 1./2. , 0., 1./2. /)
+ ZBUT(1,1) = 1./3.
+ ZBUT(1,2) = 0.
+ ZBUT(2,1) = 0.
+ ZBUT(2,2) = 1.
+ CASE('SP32')
+ ZBUTS = (/ 1./3. , 1./3. , 1./3. /)
+ ZBUT(1,1) = 1./2.
+ ZBUT(1,2) = 0.
+ ZBUT(2,1) = 1./2.
+ ZBUT(2,2) = 1./2.
+ CASE('RKC4')
+ ZBUTS = (/ 1./6. , 1./3. , 1./3. , 1./6./)
+ ZBUT = 0.
+ ZBUT(1,1) = 1./2.
+ ZBUT(2,2) = 1./2.
+ ZBUT(3,3) = 1.
+ CASE('RK4B')
+ ZBUTS = (/ 1./8. , 3./8. , 3./8. , 1./8./)
+ ZBUT = 0.
+ ZBUT(1,1) = 1./3.
+ ZBUT(2,1) = -1./3.
+ ZBUT(2,2) = 1.
+ ZBUT(3,1) = 1.
+ ZBUT(3,2) = -1.
+ ZBUT(3,3) = 1.
+ CASE('RK53')
+ ZBUTS = (/ 1./4. , 0. , 0. , 0. , 3./4. /)
+ ZBUT = 0.
+ ZBUT(1,1) = 1./7.
+ ZBUT(2,2) = 3./16.
+ ZBUT(3,3) = 1./3.
+ ZBUT(4,4) = 2./3.
+ CASE('RK62')
+ ZBUTS = (/ 1./6. , 1./6. , 1./6. , 1./6. , 1./6. , 1./6. /)
+ ZBUT = 0.
+ ZBUT(1,1) = 1./5.
+ ZBUT(2,1) = 1./5.
+ ZBUT(2,2) = 1./5.
+ ZBUT(3,1) = 1./5.
+ ZBUT(3,2) = 1./5.
+ ZBUT(3,3) = 1./5.
+ ZBUT(4,1) = 1./5.
+ ZBUT(4,2) = 1./5.
+ ZBUT(4,3) = 1./5.
+ ZBUT(4,4) = 1./5.
+ ZBUT(5,1) = 1./5.
+ ZBUT(5,2) = 1./5.
+ ZBUT(5,3) = 1./5.
+ ZBUT(5,4) = 1./5.
+ ZBUT(5,5) = 1./5.
+CASE('RK65')
+ ZBUTS= (/ 7./90. , 0. , 16./45. , 2./15. , 16./45. , 7./90. /)
+ ZBUT= 0.
+ ZBUT(1,1) = 1./4.
+ ZBUT(2,1) = 1./8.
+ ZBUT(2,2) = 1./8.
+ ZBUT(3,1) = 0
+ ZBUT(3,2) = -1./2.
+ ZBUT(3,3) = 1
+ ZBUT(4,1) = 3./16.
+ ZBUT(4,2) = 0
+ ZBUT(4,3) = 0
+ ZBUT(4,4) = 9./16.
+ ZBUT(5,1) = -3./7.
+ ZBUT(5,2) = 2./7.
+ ZBUT(5,3) = 12./7.
+ ZBUT(5,4) = -12./7.
+ ZBUT(5,5) = 8./7.
+END SELECT
+!$acc update device(ZBUTS,ZBUT)
+!
+#ifndef MNH_OPENACC
+ALLOCATE(ZRUS(SIZE(PUT,1),SIZE(PUT,2),SIZE(PWT,3),ISPL))
+ALLOCATE(ZRVS(SIZE(PUT,1),SIZE(PUT,2),SIZE(PWT,3),ISPL))
+ALLOCATE(ZRWS(SIZE(PUT,1),SIZE(PUT,2),SIZE(PWT,3),ISPL))
+#endif
+!
+!$acc kernels present(PRUS_ADV,PRVS_ADV,PRWS_ADV) present(ZUT,ZVT,ZWT) present(PU,PV,PW)
+PRUS_ADV = 0.
+PRVS_ADV = 0.
+PRWS_ADV = 0.
+!
+!-------------------------------------------------------------------------------
+!
+!* 2. Wind guess before RK loop
+! -------------------------
+!
+ZUT = PU
+ZVT = PV
+ZWT = PW
+!$acc end kernels
+!
+#ifndef MNH_OPENACC
+CALL ADV_BOUNDARIES (HLBCX, HLBCY, ZUT, PUT, 'U' )
+CALL ADV_BOUNDARIES (HLBCX, HLBCY, ZVT, PVT, 'V' )
+CALL ADV_BOUNDARIES (HLBCX, HLBCY, ZWT, PWT, 'W' )
+#else
+CALL ADV_BOUNDARIES_DEVICE (HLBCX, HLBCY, ZUT, PUT, 'U' )
+CALL ADV_BOUNDARIES_DEVICE (HLBCX, HLBCY, ZVT, PVT, 'V' )
+CALL ADV_BOUNDARIES_DEVICE (HLBCX, HLBCY, ZWT, PWT, 'W' )
+#endif
+!
+NULLIFY(TZFIELDMT_ll)
+CALL ADD3DFIELD_ll( TZFIELDMT_ll, ZUT, 'ADVECUVW_RK::ZUT' )
+CALL ADD3DFIELD_ll( TZFIELDMT_ll, ZVT, 'ADVECUVW_RK::ZVT' )
+CALL ADD3DFIELD_ll( TZFIELDMT_ll, ZWT, 'ADVECUVW_RK::ZWT' )
+INBVAR = 3
+CALL INIT_HALO2_ll(TZHALO2MT_ll,INBVAR,SIZE(PUT,1),SIZE(PUT,2),SIZE(PWT,3))
+!
+!$acc kernels
+ZRUS(:, :, :, : ) = 0.
+ZRVS(:, :, :, : ) = 0.
+ZRWS(:, :, :, : ) = 0.
+!$acc end kernels
+
+!Necessary to work around a PGI bug (19.10)
+!because following update ZRUS(:,:,:,JS) are done on the WHOLE array
+!$acc update self(ZRUS,ZRVS,ZRWS)
+
+!-------------------------------------------------------------------------------
+!
+!* 3. BEGINNING of Runge-Kutta loop
+! -----------------------------
+!
+ DO JS = 1, ISPL
+!
+#ifndef MNH_OPENACC
+ CALL ADV_BOUNDARIES (HLBCX, HLBCY, ZUT, PUT, 'U' )
+ CALL ADV_BOUNDARIES (HLBCX, HLBCY, ZVT, PVT, 'V' )
+ CALL ADV_BOUNDARIES (HLBCX, HLBCY, ZWT, PWT, 'W' )
+#else
+ CALL ADV_BOUNDARIES_DEVICE (HLBCX, HLBCY, ZUT, PUT, 'U' )
+ CALL ADV_BOUNDARIES_DEVICE (HLBCX, HLBCY, ZVT, PVT, 'V' )
+ CALL ADV_BOUNDARIES_DEVICE (HLBCX, HLBCY, ZWT, PWT, 'W' )
+#endif
+ !
+#ifndef MNH_OPENACC
+ CALL UPDATE_HALO_ll(TZFIELDMT_ll,IINFO_ll)
+ CALL UPDATE_HALO2_ll(TZFIELDMT_ll, TZHALO2MT_ll, IINFO_ll)
+#else
+! acc update self(ZUT,ZVT,ZWT)
+ CALL GET_HALO_D(ZUT,HNAME='ZUT')
+ CALL GET_HALO_D(ZVT,HNAME='ZVT')
+ CALL GET_HALO_D(ZWT,HNAME='ZWT')
+ CALL UPDATE_HALO2_ll(TZFIELDMT_ll, TZHALO2MT_ll, IINFO_ll)
+! acc update device(ZUT,ZVT,ZWT)
+#endif
+
+
+
+!
+!* 4. Advection with WENO
+! -------------------
+!
+
+ IF (HUVW_ADV_SCHEME=='WENO_K') THEN
+ CALL ADVECUVW_WENO_K (HLBCX, HLBCY, KWENO_ORDER, ZUT, ZVT, ZWT, &
+ PRUCT, PRVCT, PRWCT, &
+ ZRUS(:,:,:,JS), ZRVS(:,:,:,JS), ZRWS(:,:,:,JS), &
+#ifndef MNH_OPENACC
+ TZHALO2MT_ll )
+#else
+ TZHALO2MT_ll, ZT3D(:,:,:,IZMEAN), ZT3D(:,:,:,IZWORK) )
+#endif
+ ELSE IF ((HUVW_ADV_SCHEME=='CEN4TH') .AND. (HTEMP_SCHEME=='RKC4')) THEN
+ CALL ADVECUVW_4TH (HLBCX, HLBCY, PRUCT, PRVCT, PRWCT, &
+ ZUT, ZVT, ZWT, &
+ ZRUS(:,:,:,JS), ZRVS(:,:,:,JS), ZRWS(:,:,:,JS), &
+ TZHALO2MT_ll )
+ ENDIF
+!
+
+!
+ write ( ynum, '( I3 )' ) js
+#ifndef MNH_OPENACC
+ NULLIFY(TZFIELDS4_ll)
+ CALL ADD3DFIELD_ll( TZFIELDS4_ll, ZRUS(:,:,:,JS), 'ADVECUVW_RK::ZRUS(:,:,:,'//trim( adjustl( ynum ) )//')' )
+ CALL ADD3DFIELD_ll( TZFIELDS4_ll, ZRVS(:,:,:,JS), 'ADVECUVW_RK::ZRVS(:,:,:,'//trim( adjustl( ynum ) )//')' )
+ CALL ADD3DFIELD_ll( TZFIELDS4_ll, ZRWS(:,:,:,JS), 'ADVECUVW_RK::ZRWS(:,:,:,'//trim( adjustl( ynum ) )//')' )
+ CALL UPDATE_HALO_ll(TZFIELDS4_ll,IINFO_ll)
+ CALL CLEANLIST_ll(TZFIELDS4_ll)
+#else
+ CALL GET_HALO_D(ZRUS(:,:,:,JS),HNAME='ADVECUVW_RK::ZRUS(:,:,:,'//trim( adjustl( ynum ) )//')' )
+ CALL GET_HALO_D(ZRVS(:,:,:,JS),HNAME='ADVECUVW_RK::ZRVS(:,:,:,'//trim( adjustl( ynum ) )//')' )
+ CALL GET_HALO_D(ZRWS(:,:,:,JS),HNAME='ADVECUVW_RK::ZRWS(:,:,:,'//trim( adjustl( ynum ) )//')' )
+! acc update device(ZRUS(:,:,:,JS),ZRVS(:,:,:,JS),ZRWS(:,:,:,JS))
+#endif
+!
+!
+! Guesses at the end of the RK loop
+!
+!$acc kernels present(PRUS_ADV,PRVS_ADV,PRWS_ADV,ZBUTS) present(ZRUS,ZRVS,ZRWS)
+ PRUS_ADV(:,:,:) = PRUS_ADV(:,:,:) + ZBUTS(JS) * ZRUS(:,:,:,JS)
+ PRVS_ADV(:,:,:) = PRVS_ADV(:,:,:) + ZBUTS(JS) * ZRVS(:,:,:,JS)
+ PRWS_ADV(:,:,:) = PRWS_ADV(:,:,:) + ZBUTS(JS) * ZRWS(:,:,:,JS)
+!$acc end kernels
+!
+ IF ( JS < ISPL ) THEN
+!$acc kernels present(ZUT,ZVT,ZWT) present(ZBUT) present(PU,PV,PW) &
+!$acc & present(ZRUS,ZRVS,ZRWS) present(PRUS_OTHER,PRVS_OTHER,PRWS_OTHER) &
+!$acc & present(PMXM_RHODJ,PMYM_RHODJ,PMZM_RHODJ)
+!
+ ZUT = PU
+ ZVT = PV
+ ZWT = PW
+!
+ DO JI = 1, JS
+!
+! Intermediate guesses inside the RK loop
+!
+ ZUT(:,:,:) = ZUT(:,:,:) + ZBUT(JS,JI) * PTSTEP * &
+ ( ZRUS(:,:,:,JI) + PRUS_OTHER(:,:,:) ) / PMXM_RHODJ
+ ZVT(:,:,:) = ZVT(:,:,:) + ZBUT(JS,JI) * PTSTEP * &
+ ( ZRVS(:,:,:,JI) + PRVS_OTHER(:,:,:) ) / PMYM_RHODJ
+ ZWT(:,:,:) = ZWT(:,:,:) + ZBUT(JS,JI) * PTSTEP * &
+ ( ZRWS(:,:,:,JI) + PRWS_OTHER(:,:,:) ) / PMZM_RHODJ
+!
+ END DO
+!$acc end kernels
+!$acc update self(ZUT,ZVT,ZWT)
+ END IF
+!
+! End of the RK loop
+ END DO
+!
+!
+#ifdef MNH_OPENACC
+CALL MNH_REL_ZT3D(IZMEAN,IZWORK)
+#endif
+!
+CALL CLEANLIST_ll(TZFIELDMT_ll)
+CALL DEL_HALO2_ll(TZHALO2MT_ll)
+!$acc update self(PRUS_ADV,PRVS_ADV,PRWS_ADV)
+!-------------------------------------------------------------------------------
+!
+IF (MPPDB_INITIALIZED) THEN
+ !Check all OUT arrays
+ CALL MPPDB_CHECK(PRUS_ADV,"ADVECUVW_RK end:PRUS_ADV")
+ CALL MPPDB_CHECK(PRVS_ADV,"ADVECUVW_RK end:PRVS_ADV")
+ CALL MPPDB_CHECK(PRWS_ADV,"ADVECUVW_RK end:PRWS_ADV")
+END IF
+
+!$acc end data
+
+!$acc end data
+
+END SUBROUTINE ADVECUVW_RK
diff --git a/src/ZSOLVER/advecuvw_weno_k.f90 b/src/ZSOLVER/advecuvw_weno_k.f90
new file mode 100644
index 0000000000000000000000000000000000000000..96c6a85e5e6ea000a1a2508d7b2bc815c3fca85a
--- /dev/null
+++ b/src/ZSOLVER/advecuvw_weno_k.f90
@@ -0,0 +1,662 @@
+!MNH_LIC Copyright 2013-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.
+!-----------------------------------------------------------------
+! ###########################
+ MODULE MODI_ADVECUVW_WENO_K
+! ###########################
+!
+INTERFACE
+!
+ SUBROUTINE ADVECUVW_WENO_K(HLBCX, HLBCY, KWENO_ORDER, PUT, PVT, PWT, &
+ PRUCT, PRVCT, PRWCT, PRUS, PRVS, PRWS, TPHALO2LIST &
+#ifndef MNH_OPENACC
+ )
+#else
+ , ZMEAN, ZWORK)
+#endif
+!
+USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
+!
+CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX ! X direction LBC type
+CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY ! Y direction LBC type
+INTEGER, INTENT(IN) :: KWENO_ORDER ! Order of the WENO
+ ! scheme (3 or 5)
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRUCT ! contravariant
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRVCT ! components
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRWCT ! of momentum
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT, PVT, PWT ! U,V,W at t
+!
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRUS, PRVS, PRWS ! Source terms
+!
+TYPE(HALO2LIST_ll), POINTER :: TPHALO2LIST ! list for diffusion
+!
+#ifdef MNH_OPENACC
+! Work arrays
+REAL, DIMENSION(:,:,:) :: ZMEAN, ZWORK
+#endif
+!
+END SUBROUTINE ADVECUVW_WENO_K
+!
+END INTERFACE
+!
+END MODULE MODI_ADVECUVW_WENO_K
+!
+! ##########################################################################
+ SUBROUTINE ADVECUVW_WENO_K(HLBCX, HLBCY, KWENO_ORDER, PUT, PVT, PWT, &
+ PRUCT, PRVCT, PRWCT, PRUS, PRVS, PRWS, TPHALO2LIST &
+#ifndef MNH_OPENACC
+ )
+#else
+ , ZMEAN, ZWORK)
+#endif
+! ##########################################################################
+!
+!! AUTHOR
+!! ------
+!!
+!!
+!! MODIFICATIONS
+!! -------------
+!! J.Escobar 21/03/2013: for HALOK comment all NHALO=1 tests
+!! T.Lunet 02/10/2014: add get_halo for WENO 5
+!! suppress comment of NHALO=1 tests
+!! J.Escobar : 15/09/2015 : WENO5 & JPHEXT <> 1
+!!
+!-------------------------------------------------------------------------------
+!
+!* 0. DECLARATIONS
+! ------------
+!
+USE MODE_ll
+!
+USE MODD_PARAMETERS
+USE MODD_CONF
+USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
+!
+#ifndef MNH_OPENACC
+USE MODI_SHUMAN
+#else
+USE MODI_SHUMAN_DEVICE
+#endif
+USE MODI_ADVEC_WENO_K_1_AUX
+USE MODI_ADVEC_WENO_K_2_AUX
+USE MODI_ADVEC_WENO_K_3_AUX
+!
+USE MODD_CONF, ONLY : NHALO
+USE MODE_MPPDB
+USE MODI_GET_HALO
+!
+#ifdef MNH_OPENACC
+USE MODE_DEVICE
+USE MODE_MNH_ZWORK, ONLY : ZT3D, MNH_GET_ZT3D , MNH_REL_ZT3D
+#endif
+!
+IMPLICIT NONE
+!
+!* 0.1 Declarations of dummy arguments :
+!
+CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX ! X direction LBC type
+CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY ! Y direction LBC type
+INTEGER, INTENT(IN) :: KWENO_ORDER ! Order of the WENO
+ ! scheme (3 or 5)
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRUCT ! contravariant
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRVCT ! components
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PRWCT ! of momentum
+!
+REAL, DIMENSION(:,:,:), INTENT(IN) :: PUT, PVT, PWT ! Variables at t
+!
+REAL, DIMENSION(:,:,:), INTENT(INOUT) :: PRUS, PRVS, PRWS ! Source terms
+!
+TYPE(HALO2LIST_ll), POINTER :: TPHALO2LIST ! list for diffusion
+!
+!* 0.2 Declarations of local variables :
+!
+TYPE(HALO2LIST_ll), POINTER :: TZHALO2_UT,TZHALO2_VT,TZHALO2_WT
+
+TYPE(LIST_ll), POINTER :: TZHALO2_ZMEAN
+INTEGER :: IINFO_ll ! return code of parallel routine
+!
+#ifndef MNH_OPENACC
+REAL, DIMENSION(:,:,:), allocatable :: ZMEAN, ZWORK
+#else
+REAL, DIMENSION(:,:,:) :: ZMEAN, ZWORK
+#endif
+!
+INTEGER :: IKU
+#ifdef MNH_OPENACC
+INTEGER :: IZFPOS1, IZFPOS2, IZFPOS3
+INTEGER :: IZFNEG1, IZFNEG2, IZFNEG3
+INTEGER :: IZBPOS1, IZBPOS2, IZBPOS3
+INTEGER :: IZBNEG1, IZBNEG2, IZBNEG3
+INTEGER :: IZOMP1, IZOMP2, IZOMP3
+INTEGER :: IZOMN1, IZOMN2, IZOMN3
+#endif
+!
+!$acc data present( PRUCT, PRVCT, PRWCT, PUT, PVT, PWT, PRUS, PRVS, PRWS, ZMEAN, ZWORK )
+!
+IF (MPPDB_INITIALIZED) THEN
+ !Check all IN arrays
+ CALL MPPDB_CHECK(PRUCT,"ADVECUVW_WENO_K beg:PRUCT")
+ CALL MPPDB_CHECK(PRVCT,"ADVECUVW_WENO_K beg:PRVCT")
+ CALL MPPDB_CHECK(PRWCT,"ADVECUVW_WENO_K beg:PRWCT")
+ CALL MPPDB_CHECK(PUT,"ADVECUVW_WENO_K beg:PUT")
+ CALL MPPDB_CHECK(PVT,"ADVECUVW_WENO_K beg:PVT")
+ CALL MPPDB_CHECK(PWT,"ADVECUVW_WENO_K beg:PWT")
+ !Check all INOUT arrays
+ CALL MPPDB_CHECK(PRUS,"ADVECUVW_WENO_K beg:PRUS")
+ CALL MPPDB_CHECK(PRVS,"ADVECUVW_WENO_K beg:PRVS")
+ CALL MPPDB_CHECK(PRWS,"ADVECUVW_WENO_K beg:PRWS")
+END IF
+
+#ifndef MNH_OPENACC
+allocate(ZMEAN(SIZE(PUT,1), SIZE(PUT,2), SIZE(PUT,3)))
+allocate(ZWORK(SIZE(PUT,1), SIZE(PUT,2), SIZE(PUT,3)))
+#endif
+
+#ifdef MNH_OPENACC
+CALL INIT_ON_HOST_AND_DEVICE(ZMEAN,1e90,'ADVECUVW_WENO_K::ZMEAN')
+CALL INIT_ON_HOST_AND_DEVICE(ZWORK,2e90,'ADVECUVW_WENO_K::ZWORK')
+#endif
+!
+!------------------------- ADVECTION OF MOMENTUM ------------------------------
+!
+!
+TZHALO2_UT => TPHALO2LIST ! 1rst add3dfield in model_n
+TZHALO2_VT => TPHALO2LIST%NEXT ! 2nd add3dfield in model_n
+TZHALO2_WT => TPHALO2LIST%NEXT%NEXT ! 3rst add3dfield in model_n
+!
+IKU=SIZE(PUT,3)
+! -------------------------------------------------------
+!
+SELECT CASE(KWENO_ORDER)
+!
+CASE(1) ! WENO 1
+#ifndef MNH_OPENACC
+!
+! U component
+!
+ PRUS = PRUS - DXM(UP_UX(PUT,MXF(PRUCT)))
+!
+ PRUS = PRUS - DYF(UP_MY(PUT,MXM(PRVCT)))
+!
+ PRUS = PRUS - DZF(UP_MZ(PUT,MXM(PRWCT)))
+!
+! V component
+!
+ PRVS = PRVS - DXF(UP_MX(PVT,MYM(PRUCT)))
+!
+ PRVS = PRVS - DYM(UP_VY(PVT,MYF(PRVCT)))
+!
+ PRVS = PRVS - DZF(UP_MZ(PVT,MYM(PRWCT)))
+!
+! W component
+!
+ PRWS = PRWS - DXF(UP_MX(PWT,MZM(PRUCT)))
+!
+ PRWS = PRWS - DYF(UP_MY(PWT,MZM(PRVCT)))
+!
+ PRWS = PRWS - DZM(UP_WZ(PWT,MZF(PRWCT)))
+#else
+!
+! U component
+!
+ !PRUS = PRUS - DXM(UP_UX(PUT,MXF(PRUCT)))
+ CALL MXF_DEVICE(PRUCT,ZWORK)
+ CALL UP_UX_DEVICE(PUT,ZWORK,ZMEAN)
+ CALL DXM_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRUS = PRUS - ZWORK
+!$acc end kernels
+!
+ !PRUS = PRUS - DYF(UP_MY(PUT,MXM(PRVCT)))
+ CALL MXM_DEVICE(PRVCT,ZWORK)
+ CALL UP_MY_DEVICE(PUT,ZWORK,ZMEAN)
+ CALL DYF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRUS = PRUS - ZWORK
+!$acc end kernels
+!
+ !PRUS = PRUS - DZF(1,IKU,1,UP_MZ(PUT,MXM(PRWCT)))
+ CALL MXM_DEVICE(PRWCT,ZWORK)
+ CALL UP_MZ_DEVICE(PUT,ZWORK,ZMEAN)
+ CALL DZF_DEVICE(1,IKU,1,ZMEAN,ZWORK)
+!$acc kernels
+ PRUS = PRUS - ZWORK
+!$acc end kernels
+!
+! V component
+!
+ !PRVS = PRVS - DXF(UP_MX(PVT,MYM(PRUCT)))
+ CALL MYM_DEVICE(PRUCT,ZWORK)
+ CALL UP_MX_DEVICE(PVT,ZWORK,ZMEAN)
+ CALL DXF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRVS = PRVS - ZWORK
+!$acc end kernels
+!
+ !PRVS = PRVS - DYM(UP_VY(PVT,MYF(PRVCT)))
+ CALL MYF_DEVICE(PRVCT,ZWORK)
+ CALL UP_VY_DEVICE(PVT,ZWORK,ZMEAN)
+ CALL DYM_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRVS = PRVS - ZWORK
+!$acc end kernels
+!
+ !PRVS = PRVS - DZF(1,IKU,1,UP_MZ(PVT,MYM(PRWCT)))
+ CALL MYM_DEVICE(PRWCT,ZWORK)
+ CALL UP_MZ_DEVICE(PVT,ZWORK,ZMEAN)
+ CALL DZF_DEVICE(1,IKU,1,ZMEAN,ZWORK)
+!$acc kernels
+ PRVS = PRVS - ZWORK
+!$acc end kernels
+!
+! W component
+!
+ !PRWS = PRWS - DXF(UP_MX(PWT,MZM(1,IKU,1,PRUCT)))
+ CALL MZM_DEVICE(PRUCT,ZWORK)
+ CALL UP_MX_DEVICE(PWT,ZWORK,ZMEAN)
+ CALL DXF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRWS = PRWS - ZWORK
+!$acc end kernels
+!
+ !PRWS = PRWS - DYF(UP_MY(PWT,MZM(1,IKU,1,PRVCT)))
+ CALL MZM_DEVICE(PRVCT,ZWORK)
+ CALL UP_MY_DEVICE(PWT,ZWORK,ZMEAN)
+ CALL DYF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRWS = PRWS - ZWORK
+!$acc end kernels
+!
+ !PRWS = PRWS - DZM(1,IKU,1,UP_WZ(PWT,MZF(1,IKU,1,PRWCT)))
+ CALL MZF_DEVICE(1,IKU,1,PRWCT,ZWORK)
+ CALL UP_WZ_DEVICE(PWT,ZWORK,ZMEAN)
+ CALL DZM_DEVICE(1,IKU,1,ZMEAN,ZWORK)
+!$acc kernels
+ PRWS = PRWS - ZWORK
+!$acc end kernels
+#endif
+!
+!
+CASE(3) ! WENO 3
+#ifndef MNH_OPENACC
+!
+! U component
+!
+ ZWORK = MXF(PRUCT)
+ CALL ADVEC_WENO_K_2_UX(HLBCX, PUT, ZWORK, ZMEAN, TZHALO2_UT%HALO2)
+ PRUS = PRUS - DXM(ZMEAN)
+
+!
+ IF (.NOT.L2D) THEN
+ ZWORK = MXM(PRVCT)
+ CALL ADVEC_WENO_K_2_MY(HLBCY, PUT, ZWORK, ZMEAN, TZHALO2_UT%HALO2)
+ PRUS = PRUS - DYF(ZMEAN)
+ END IF
+!
+ PRUS = PRUS - DZF(WENO_K_2_MZ(PUT, MXM(PRWCT)))
+!
+! V component
+!
+ IF (.NOT.L2D) THEN
+ ZWORK = MYM(PRUCT)
+ CALL ADVEC_WENO_K_2_MX(HLBCX, PVT, ZWORK, ZMEAN, TZHALO2_VT%HALO2)
+ PRVS = PRVS - DXF(ZMEAN)
+!
+ ZWORK = MYF(PRVCT)
+ CALL ADVEC_WENO_K_2_VY(HLBCY, PVT, ZWORK, ZMEAN, TZHALO2_VT%HALO2)
+ PRVS = PRVS - DYM(ZMEAN)
+!
+ PRVS = PRVS - DZF(WENO_K_2_MZ(PVT, MYM(PRWCT)))
+ END IF
+!
+! W component
+!
+ ZWORK = MZM(PRUCT)
+ CALL ADVEC_WENO_K_2_MX(HLBCX, PWT, ZWORK, ZMEAN, TZHALO2_WT%HALO2)
+ PRWS = PRWS - DXF(ZMEAN)
+!
+ IF (.NOT.L2D) THEN
+ ZWORK = MZM(PRVCT)
+ CALL ADVEC_WENO_K_2_MY(HLBCY, PWT, ZWORK, ZMEAN, TZHALO2_WT%HALO2)
+ PRWS = PRWS - DYF(ZMEAN)
+ END IF
+!
+ PRWS = PRWS - DZM(WENO_K_2_WZ(PWT,MZF(PRWCT)))
+#else
+ CALL MNH_GET_ZT3D(IZFPOS1,IZFPOS2,IZFNEG1,IZFNEG2,IZBPOS1,IZBPOS2,IZBNEG1,IZBNEG2,IZOMP1,IZOMP2,IZOMN1,IZOMN2)
+!
+! U component
+!
+ CALL MXF_DEVICE(PRUCT,ZWORK)
+ CALL ADVEC_WENO_K_2_UX(HLBCX, PUT, ZWORK, ZMEAN, TZHALO2_UT%HALO2%WEST, TZHALO2_UT%HALO2%EAST, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2) )
+ CALL DXM_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRUS = PRUS - ZWORK
+!$acc end kernels
+!
+ IF (.NOT.L2D) THEN
+ CALL MXM_DEVICE(PRVCT,ZWORK)
+ CALL ADVEC_WENO_K_2_MY(HLBCY, PUT, ZWORK, ZMEAN, TZHALO2_UT%HALO2%NORTH, TZHALO2_UT%HALO2%SOUTH, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2) )
+ CALL DYF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRUS = PRUS - ZWORK
+!$acc end kernels
+ END IF
+!
+! PRUS = PRUS - DZF(1,IKU,1,WENO_K_2_MZ(PUT, MXM(PRWCT)))
+ CALL MXM_DEVICE(PRWCT,ZWORK)
+ CALL WENO_K_2_MZ(PUT, ZWORK, ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2) )
+ CALL DZF_DEVICE(1,IKU,1,ZMEAN,ZWORK)
+!$acc kernels
+ PRUS = PRUS - ZWORK
+!$acc end kernels
+!
+! V component
+!
+ IF (.NOT.L2D) THEN
+ CALL MYM_DEVICE(PRUCT,ZWORK)
+ CALL ADVEC_WENO_K_2_MX(HLBCX, PVT, ZWORK, ZMEAN, TZHALO2_VT%HALO2%WEST, TZHALO2_VT%HALO2%EAST, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2) )
+ CALL DXF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRVS = PRVS - ZWORK
+!$acc end kernels
+!
+ CALL MYF_DEVICE(PRVCT,ZWORK)
+ CALL ADVEC_WENO_K_2_VY(HLBCY, PVT, ZWORK, ZMEAN, TZHALO2_VT%HALO2%NORTH, TZHALO2_VT%HALO2%SOUTH, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2) )
+ CALL DYM_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRVS = PRVS - ZWORK
+!$acc end kernels
+!
+! PRVS = PRVS - DZF(1,IKU,1,WENO_K_2_MZ(PVT, MYM(PRWCT)))
+ CALL MYM_DEVICE(PRWCT,ZWORK)
+ CALL WENO_K_2_MZ(PVT, ZWORK, ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2) )
+ CALL DZF_DEVICE(1,IKU,1,ZMEAN,ZWORK)
+!$acc kernels
+ PRVS = PRVS - ZWORK
+!$acc end kernels
+ END IF
+!
+! W component
+!
+! ZWORK = MZM(1,IKU,1,PRUCT)
+ CALL MZM_DEVICE(PRUCT,ZWORK)
+ CALL ADVEC_WENO_K_2_MX(HLBCX, PWT, ZWORK, ZMEAN, TZHALO2_WT%HALO2%WEST, TZHALO2_WT%HALO2%EAST, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2) )
+ CALL DXF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRWS = PRWS - ZWORK
+!$acc end kernels
+!
+ IF (.NOT.L2D) THEN
+! ZWORK = MZM(1,IKU,1,PRVCT)
+ CALL MZM_DEVICE(PRVCT,ZWORK)
+ CALL ADVEC_WENO_K_2_MY(HLBCY, PWT, ZWORK, ZMEAN, TZHALO2_WT%HALO2%NORTH, TZHALO2_WT%HALO2%SOUTH, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2) )
+ CALL DYF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRWS = PRWS - ZWORK
+!$acc end kernels
+ END IF
+!
+! PRWS = PRWS - DZM(1,IKU,1,WENO_K_2_WZ(PWT,MZF(1,IKU,1,PRWCT)))
+ CALL MZF_DEVICE(1,IKU,1,PRWCT,ZWORK)
+ CALL WENO_K_2_WZ(PWT, ZWORK, ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2) )
+ CALL DZM_DEVICE(1,IKU,1,ZMEAN,ZWORK)
+!$acc kernels
+ PRWS = PRWS - ZWORK
+!$acc end kernels
+!
+ CALL MNH_REL_ZT3D(IZFPOS1,IZFPOS2,IZFNEG1,IZFNEG2,IZBPOS1,IZBPOS2,IZBNEG1,IZBNEG2,IZOMP1,IZOMP2,IZOMN1,IZOMN2)
+#endif
+!
+!
+CASE(5) ! WENO 5
+#ifndef MNH_OPENACC
+!
+! U component
+!
+ ZWORK = MXF(PRUCT)
+ CALL ADVEC_WENO_K_3_UX(HLBCX, PUT, ZWORK, ZMEAN)
+ CALL GET_HALO(ZMEAN)! Update HALO
+ PRUS = PRUS - DXM(ZMEAN)
+!
+ IF (.NOT.L2D) THEN! 3D Case
+ ZWORK = MXM(PRVCT)
+ CALL ADVEC_WENO_K_3_MY(HLBCY, PUT, ZWORK, ZMEAN)
+ CALL GET_HALO(ZMEAN)! Update HALO
+ PRUS = PRUS - DYF(ZMEAN)
+ END IF
+!
+ ZMEAN = WENO_K_3_MZ(PUT, MXM(PRWCT))
+ CALL GET_HALO(ZMEAN)! Update HALO - maybe not necessary (T.Lunet)
+ PRUS = PRUS - DZF(ZMEAN)
+!
+! V component, only called in 3D case
+!
+ IF (.NOT.L2D) THEN
+!
+ ZWORK = MYM(PRUCT)
+ CALL ADVEC_WENO_K_3_MX(HLBCX, PVT, ZWORK, ZMEAN)
+ CALL GET_HALO(ZMEAN)! Update HALO
+ PRVS = PRVS - DXF(ZMEAN)
+!
+ ZWORK = MYF(PRVCT)
+ CALL ADVEC_WENO_K_3_VY(HLBCY, PVT, ZWORK, ZMEAN)
+ CALL GET_HALO(ZMEAN)! Update HALO
+ PRVS = PRVS - DYM(ZMEAN)
+!
+ ZMEAN = WENO_K_3_MZ(PVT, MYM(PRWCT))
+ CALL GET_HALO(ZMEAN)! Update HALO - maybe not necessary (T.Lunet)
+ PRVS = PRVS - DZF(ZMEAN)
+!
+ END IF
+!
+! W component
+!
+ ZWORK = MZM(PRUCT)
+ CALL ADVEC_WENO_K_3_MX(HLBCX, PWT, ZWORK, ZMEAN)
+ CALL GET_HALO(ZMEAN)! Update HALO
+ PRWS = PRWS - DXF(ZMEAN)
+!
+ IF (.NOT.L2D) THEN! 3D Case
+ ZWORK = MZM(PRVCT)
+ CALL ADVEC_WENO_K_3_MY(HLBCY, PWT, ZWORK, ZMEAN)
+ CALL GET_HALO(ZMEAN)! Update HALO
+ PRWS = PRWS - DYF(ZMEAN)
+ END IF
+!
+ ZMEAN = WENO_K_3_WZ(PWT,MZF(PRWCT))
+ CALL GET_HALO(ZMEAN)! Update HALO - maybe not necessary (T.Lunet)
+ PRWS = PRWS - DZM(ZMEAN)
+#else
+ CALL MNH_GET_ZT3D(IZFPOS1,IZFPOS2,IZFPOS3,IZFNEG1,IZFNEG2,IZFNEG3,IZBPOS1, &
+ IZBPOS2,IZBPOS3,IZBNEG1,IZBNEG2,IZBNEG3,IZOMP1,IZOMP2,IZOMP3,IZOMN1,IZOMN2,IZOMN3)
+!
+! U component
+!
+ CALL MXF_DEVICE(PRUCT,ZWORK)
+ CALL ADVEC_WENO_K_3_UX(HLBCX, PUT, ZWORK, ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFPOS3), &
+ ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), ZT3D(:,:,:,IZFNEG3), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBPOS3), &
+ ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), ZT3D(:,:,:,IZBNEG3), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMP3), &
+ ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2), ZT3D(:,:,:,IZOMN3) )
+ CALL GET_HALO_D(ZMEAN)! Update HALO
+ CALL DXM_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRUS = PRUS - ZWORK
+!$acc end kernels
+!
+ IF (.NOT.L2D) THEN
+ CALL MXM_DEVICE(PRVCT,ZWORK)
+ CALL ADVEC_WENO_K_3_MY(HLBCY, PUT, ZWORK, ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFPOS3), &
+ ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), ZT3D(:,:,:,IZFNEG3), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBPOS3), &
+ ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), ZT3D(:,:,:,IZBNEG3), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMP3), &
+ ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2), ZT3D(:,:,:,IZOMN3) )
+ CALL GET_HALO_D(ZMEAN)! Update HALO
+ CALL DYF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRUS = PRUS - ZWORK
+!$acc end kernels
+ END IF
+!
+ CALL MXM_DEVICE(PRWCT,ZWORK)
+ CALL WENO_K_3_MZ(PUT,ZWORK,ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFPOS3), &
+ ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), ZT3D(:,:,:,IZFNEG3), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBPOS3), &
+ ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), ZT3D(:,:,:,IZBNEG3), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMP3), &
+ ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2), ZT3D(:,:,:,IZOMN3) )
+ CALL GET_HALO_D(ZMEAN)! Update HALO - maybe not necessary (T.Lunet)
+ CALL DZF_DEVICE(1,IKU,1,ZMEAN,ZWORK)
+!$acc kernels
+ PRUS = PRUS - ZWORK
+!$acc end kernels
+!
+! V component
+!
+ IF (.NOT.L2D) THEN
+ CALL MYM_DEVICE(PRUCT,ZWORK)
+ CALL ADVEC_WENO_K_3_MX(HLBCX, PVT, ZWORK, ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFPOS3), &
+ ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), ZT3D(:,:,:,IZFNEG3), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBPOS3), &
+ ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), ZT3D(:,:,:,IZBNEG3), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMP3), &
+ ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2), ZT3D(:,:,:,IZOMN3) )
+ CALL GET_HALO_D(ZMEAN)! Update HALO
+ CALL DXF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRVS = PRVS - ZWORK
+!$acc end kernels
+!
+ CALL MYF_DEVICE(PRVCT,ZWORK)
+ CALL ADVEC_WENO_K_3_VY(HLBCY, PVT, ZWORK, ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFPOS3), &
+ ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), ZT3D(:,:,:,IZFNEG3), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBPOS3), &
+ ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), ZT3D(:,:,:,IZBNEG3), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMP3), &
+ ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2), ZT3D(:,:,:,IZOMN3) )
+ CALL GET_HALO_D(ZMEAN)! Update HALO
+ CALL DYM_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRVS = PRVS - ZWORK
+!$acc end kernels
+!
+ CALL MYM_DEVICE(PRWCT,ZWORK)
+ CALL WENO_K_3_MZ(PVT,ZWORK,ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFPOS3), &
+ ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), ZT3D(:,:,:,IZFNEG3), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBPOS3), &
+ ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), ZT3D(:,:,:,IZBNEG3), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMP3), &
+ ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2), ZT3D(:,:,:,IZOMN3) )
+ CALL GET_HALO_D(ZMEAN)! Update HALO - maybe not necessary (T.Lunet)
+ CALL DZF_DEVICE(1,IKU,1,ZMEAN,ZWORK)
+!$acc kernels
+ PRVS = PRVS - ZWORK
+!$acc end kernels
+ END IF
+!
+! W component
+!
+ CALL MZM_DEVICE(PRUCT,ZWORK)
+ CALL ADVEC_WENO_K_3_MX(HLBCX, PWT, ZWORK, ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFPOS3), &
+ ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), ZT3D(:,:,:,IZFNEG3), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBPOS3), &
+ ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), ZT3D(:,:,:,IZBNEG3), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMP3), &
+ ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2), ZT3D(:,:,:,IZOMN3) )
+ CALL GET_HALO_D(ZMEAN)! Update HALO
+ CALL DXF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRWS = PRWS - ZWORK
+!$acc end kernels
+!
+ IF (.NOT.L2D) THEN
+ CALL MZM_DEVICE(PRVCT,ZWORK)
+ CALL ADVEC_WENO_K_3_MY(HLBCY, PWT, ZWORK, ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFPOS3), &
+ ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), ZT3D(:,:,:,IZFNEG3), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBPOS3), &
+ ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), ZT3D(:,:,:,IZBNEG3), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMP3), &
+ ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2), ZT3D(:,:,:,IZOMN3) )
+ CALL GET_HALO_D(ZMEAN)! Update HALO
+ CALL DYF_DEVICE(ZMEAN,ZWORK)
+!$acc kernels
+ PRWS = PRWS - ZWORK
+!$acc end kernels
+ END IF
+!
+ CALL MZF_DEVICE(1,IKU,1,PRWCT,ZWORK)
+ CALL WENO_K_3_WZ(PWT,ZWORK,ZMEAN, &
+ ZT3D(:,:,:,IZFPOS1), ZT3D(:,:,:,IZFPOS2), ZT3D(:,:,:,IZFPOS3), &
+ ZT3D(:,:,:,IZFNEG1), ZT3D(:,:,:,IZFNEG2), ZT3D(:,:,:,IZFNEG3), &
+ ZT3D(:,:,:,IZBPOS1), ZT3D(:,:,:,IZBPOS2), ZT3D(:,:,:,IZBPOS3), &
+ ZT3D(:,:,:,IZBNEG1), ZT3D(:,:,:,IZBNEG2), ZT3D(:,:,:,IZBNEG3), &
+ ZT3D(:,:,:,IZOMP1), ZT3D(:,:,:,IZOMP2), ZT3D(:,:,:,IZOMP3), &
+ ZT3D(:,:,:,IZOMN1), ZT3D(:,:,:,IZOMN2), ZT3D(:,:,:,IZOMN3) )
+ CALL GET_HALO_D(ZMEAN)! Update HALO - maybe not necessary (T.Lunet)
+ CALL DZM_DEVICE(1,IKU,1,ZMEAN,ZWORK)
+!$acc kernels
+ PRWS = PRWS - ZWORK
+!$acc end kernels
+!
+ CALL MNH_REL_ZT3D(IZFPOS1,IZFPOS2,IZFPOS3,IZFNEG1,IZFNEG2,IZFNEG3,IZBPOS1, &
+ IZBPOS2,IZBPOS3,IZBNEG1,IZBNEG2,IZBNEG3,IZOMP1,IZOMP2,IZOMP3,IZOMN1,IZOMN2,IZOMN3)
+#endif
+!
+!
+END SELECT
+! ---------------------------------
+!$acc update self(PRUS,PRVS,PRWS)
+!
+IF (MPPDB_INITIALIZED) THEN
+ CALL MPPDB_CHECK(PRUS,"ADVECUVW_WENO_K end:PRUS")
+ CALL MPPDB_CHECK(PRVS,"ADVECUVW_WENO_K end:PRVS")
+ CALL MPPDB_CHECK(PRWS,"ADVECUVW_WENO_K end:PRWS")
+END IF
+
+!$acc end data
+
+END SUBROUTINE ADVECUVW_WENO_K