!MNH_LIC Copyright 1994-2014 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.
!-----------------------------------------------------------------
!--------------- special set of characters for RCS information
!-----------------------------------------------------------------
! $Source$ $Revision$
! MASDEV4_7 operators 2006/05/18 13:07:25
!-----------------------------------------------------------------
!     #########################
      MODULE MODI_SHUMAN_DEVICE
!     #########################
!
INTERFACE
!
SUBROUTINE DXF_DEVICE(PA,PDXF)
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDXF   ! result at mass localization
!$acc declare present(PA,PDXM)
END SUBROUTINE DXF_DEVICE
!
SUBROUTINE DXM_DEVICE(PA,PDXM)
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDXM   ! result at flux side
!$acc declare present(PA,PDXM)
END SUBROUTINE DXM_DEVICE
!
SUBROUTINE DYF_DEVICE(PA,PDYF)
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDYF   ! result at mass localization
!$acc declare present(PA,PDYF)
END SUBROUTINE DYF_DEVICE
!
SUBROUTINE DYM_DEVICE(PA,PDYM)
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDYM   ! result at flux side
!$acc declare present(PA,PDYM)
END SUBROUTINE DYM_DEVICE
!
SUBROUTINE DZF_DEVICE(KKA,KKU,KL,PA,PDZF)
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 flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDZF     ! result at mass localization
!$acc declare present(PA,PDZF)
END SUBROUTINE DZF_DEVICE
!
SUBROUTINE DZM_DEVICE(KKA,KKU,KL,PA,PDZM)
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 mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDZM     ! result at flux side
!$acc declare present(PA,PDZM)
END SUBROUTINE DZM_DEVICE
!
SUBROUTINE MXF_DEVICE(PA,PMXF)
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMXF   ! result at mass localization
!$acc declare present(PA,PMXM)
END SUBROUTINE MXF_DEVICE
!
SUBROUTINE MXM_DEVICE(PA,PMXM)
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMXM   ! result at flux localization
!$acc declare present(PA,PMXM)
END SUBROUTINE MXM_DEVICE
!
SUBROUTINE MYF_DEVICE(PA,PMYF)
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMYF   ! result at mass localization
!$acc declare present(PA,PMYF)
END SUBROUTINE MYF_DEVICE
!
SUBROUTINE MYM_DEVICE(PA,PMYM)
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMYM   ! result at flux localization
!$acc declare present(PA,PMYM)
END SUBROUTINE MYM_DEVICE
!
SUBROUTINE MZF_DEVICE(KKA,KKU,KL,PA,PMZF)
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 flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMZF     ! result at mass localization
!$acc declare present(PA,PMZF)
END SUBROUTINE MZF_DEVICE
!
SUBROUTINE MZM_DEVICE(PA,PMZM)
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMZM   ! result at flux localization
!$acc declare present(PA,PMZM)
END SUBROUTINE MZM_DEVICE
!
END INTERFACE
!
END MODULE MODI_SHUMAN_DEVICE
!
!
!     ###############################
      SUBROUTINE MXF_DEVICE(PA,PMXF)
!     ###############################
!
!!****  *MXF* -  Shuman operator : mean operator in x direction for a 
!!                                 variable at a flux side
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a mean 
!     along the x direction (I index) for a field PA localized at a x-flux
!     point (u point). The result is localized at a mass point.
!
!!**  METHOD
!!    ------ 
!!        The result PMXF(i,:,:) is defined by 0.5*(PA(i,:,:)+PA(i+1,:,:))
!!        At i=size(PA,1), PMXF(i,:,:) are replaced by the values of PMXF,
!!    which are the right values in the x-cyclic case
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      Module MODD_PARAMETERS: declaration of parameter variables
!!        JPHEXT: define the number of marginal points out of the 
!!        physical domain along the horizontal directions.
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    04/07/94 
!!      Modification to include the periodic case 13/10/94 J.Stein 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE MODD_PARAMETERS
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument and result
!              ------------------------------------
!
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMXF   ! result at mass localization
!$acc declare present(PA,PMXF)
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
INTEGER :: JI             ! Loop index in x direction
INTEGER :: IIU            ! upper bound in x direction of PA
!
INTEGER :: JJK,IJU,IKU
INTEGER :: JIJK,JIJKOR,JIJKEND
!
!
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF MXF
!              ------------------
!
IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + JPHEXT
JIJKEND = IIU*IJU*IKU
!
!$acc kernels present(PMXF,PA)
!CDIR NODEP
!OCL NOVREC
DO JIJK=JIJKOR , JIJKEND
  PMXF(JIJK-1,1,1) = 0.5*( PA(JIJK-1,1,1)+PA(JIJK,1,1) )
END DO
!
!CDIR NODEP
!OCL NOVREC
DO JJK=1,IJU*IKU
   PMXF(IIU,JJK,1)    = PMXF(2*JPHEXT,JJK,1)
END DO
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE MXF_DEVICE
!
!     ###############################
      SUBROUTINE MXM_DEVICE(PA,PMXM)
!     ###############################
!
!!****  *MXM* -  Shuman operator : mean operator in x direction for a 
!!                                 mass variable 
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a mean 
!     along the x direction (I index) for a field PA localized at a mass
!     point. The result is localized at a x-flux point (u point).
!
!!**  METHOD
!!    ------ 
!!        The result PMXM(i,:,:) is defined by 0.5*(PA(i,:,:)+PA(i-1,:,:))
!!    At i=1, PMXM(1,:,:) are replaced by the values of PMXM,
!!    which are the right values in the x-cyclic case. 
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      Module MODD_PARAMETERS: declaration of parameter variables
!!        JPHEXT: define the number of marginal points out of the 
!!        physical domain along the horizontal directions.
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    04/07/94
!!      Modification to include the periodic case 13/10/94 J.Stein 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE MODD_PARAMETERS
!
IMPLICIT NONE

!
!*       0.1   Declarations of argument and result
!              ------------------------------------
!
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMXM   ! result at flux localization
!$acc declare present(PA,PMXM)

INTEGER :: JI             ! Loop index in x direction
INTEGER :: IIU            ! Size of the array in the x direction
!
INTEGER :: JJK,IJU,IKU
INTEGER :: JIJK,JIJKOR,JIJKEND

!                     
!
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF MXM
!              ------------------
!
IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + JPHEXT
JIJKEND = IIU*IJU*IKU
!
!CDIR NODEP
!OCL NOVREC
!$acc kernels present(PA,PMXM)
DO JIJK=JIJKOR , JIJKEND
   PMXM(JIJK,1,1) = 0.5*( PA(JIJK,1,1)+PA(JIJK-1,1,1) )
END DO
!
!CDIR NODEP
!OCL NOVREC
DO JJK=1,IJU*IKU
   PMXM(1,JJK,1)    = PMXM(IIU-2*JPHEXT+1,JJK,1)
END DO
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE MXM_DEVICE
!
!     ###############################
      SUBROUTINE MYF_DEVICE(PA,PMYF)
!     ###############################
!
!!****  *MYF* -  Shuman operator : mean operator in y direction for a 
!!                                 variable at a flux side
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a mean 
!     along the y direction (J index) for a field PA localized at a y-flux
!     point (v point). The result is localized at a mass point.
!
!!**  METHOD
!!    ------ 
!!        The result PMYF(i,:,:) is defined by 0.5*(PA(:,j,:)+PA(:,j+1,:))
!!        At j=size(PA,2), PMYF(:,j,:) are replaced by the values of PMYF,
!!    which are the right values in the y-cyclic case
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      Module MODD_PARAMETERS: declaration of parameter variables
!!        JPHEXT: define the number of marginal points out of the 
!!        physical domain along the horizontal directions.
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    04/07/94 
!!      Modification to include the periodic case 13/10/94 J.Stein 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE MODD_PARAMETERS
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument and result
!              ------------------------------------
!
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMYF   ! result at mass localization
!$acc declare present(PA,PMYF)
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
INTEGER :: JJ             ! Loop index in y direction
INTEGER :: IJU            ! upper bound in y direction of PA
!
INTEGER :: IIU,IKU
INTEGER :: JIJK,JIJKOR,JIJKEND
!
!
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF MYF
!              ------------------
!
IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + IIU
JIJKEND = IIU*IJU*IKU
!
!$acc kernels present(PA,PMYF)
!CDIR NODEP
!OCL NOVREC
DO JIJK=JIJKOR , JIJKEND
   PMYF(JIJK-IIU,1,1) = 0.5*( PA(JIJK-IIU,1,1)+PA(JIJK,1,1) )
END DO
!
PMYF(:,IJU,:)    = PMYF(:,2*JPHEXT,:)
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE MYF_DEVICE
!
!     ###############################
      SUBROUTINE MYM_DEVICE(PA,PMYM)
!     ###############################
!
!!****  *MYM* -  Shuman operator : mean operator in y direction for a 
!!                                 mass variable 
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a mean 
!     along the y direction (J index) for a field PA localized at a mass
!     point. The result is localized at a y-flux point (v point).
!
!!**  METHOD
!!    ------ 
!!        The result PMYM(:,j,:) is defined by 0.5*(PA(:,j,:)+PA(:,j-1,:))
!!    At j=1, PMYM(:,j,:) are replaced by the values of PMYM,
!!    which are the right values in the y-cyclic case. 
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      Module MODD_PARAMETERS: declaration of parameter variables
!!        JPHEXT: define the number of marginal points out of the 
!!        physical domain along the horizontal directions.
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    04/07/94 
!!      Modification to include the periodic case 13/10/94 J.Stein 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE MODD_PARAMETERS
!
IMPLICIT NONE

REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMYM   ! result at flux localization
!$acc declare present(PA,PMYM)

INTEGER :: JI             ! Loop index in x direction
INTEGER :: IIU            ! Size of the array in the x direction
!
INTEGER :: JJK,IJU,IKU
INTEGER :: JIJK,JIJKOR,JIJKEND

IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + IIU
JIJKEND = IIU*IJU*IKU
!CDIR NODEP
!OCL NOVREC
!$acc kernels present(PA,PMYM)
DO JIJK=JIJKOR , JIJKEND
   PMYM(JIJK,1,1) = 0.5*( PA(JIJK,1,1)+PA(JIJK-IIU,1,1) )
END DO
!
PMYM(:,1,:)    = PMYM(:,IJU-2*JPHEXT+1,:)
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE MYM_DEVICE
!
!     ###############################
      SUBROUTINE MZF_DEVICE(KKA,KKU,KL,PA,PMZF)
!     ###############################
!
!!****  *MZF* -  Shuman operator : mean operator in z direction for a 
!!                                 variable at a flux side
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a mean 
!     along the z direction (K index) for a field PA localized at a z-flux
!     point (w point). The result is localized at a mass point.
!
!!**  METHOD
!!    ------ 
!!        The result PMZF(:,:,k) is defined by 0.5*(PA(:,:,k)+PA(:,:,k+1))
!!        At k=size(PA,3), PMZF(:,:,k) is defined by -999.
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      NONE
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    04/07/94 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument 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 flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMZF     ! result at mass localization 
!$acc declare present(PA,PMZF)
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
INTEGER :: JK             ! Loop index in z direction
INTEGER :: IKU          ! upper bound in z direction of PA 
!     
INTEGER :: IIU,IJU
INTEGER :: JIJ
INTEGER :: JIJK,JIJKOR,JIJKEND
!            
!
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF MZF
!              ------------------
!
IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + IIU*IJU
JIJKEND = IIU*IJU*IKU
!
!$acc kernels present(PA,PMZF)
!CDIR NODEP
!OCL NOVREC
!DO JIJK=JIJKOR , JIJKEND
!   PMZF(JIJK-IIU*IJU,1,1) = 0.5*( PA(JIJK-IIU*IJU,1,1)+PA(JIJK,1,1) )
!END DO
PMZF(:,:,1:IKU-1) = 0.5*( PA(:,:,1:IKU-1)+PA(:,:,2:) )
!
!CDIR NODEP
!OCL NOVREC
!DO JIJ=1,IIU*IJU
!   PMZF(JIJ,1,IKU)    = -999.
!END DO
PMZF(:,:,IKU) = -999.
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE MZF_DEVICE
!
!     ###############################
      SUBROUTINE MZM_DEVICE(PA,PMZM)
!     ###############################
!
!!****  *MZM* -  Shuman operator : mean operator in z direction for a 
!!                                 mass variable 
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a mean
!     along the z direction (K index) for a field PA localized at a mass
!     point. The result is localized at a z-flux point (w point).
!
!!**  METHOD
!!    ------ 
!!        The result PMZM(:,:,k) is defined by 0.5*(PA(:,:,k)+PA(:,:,k-1))
!!        At k=1, PMZM(:,:,1) is defined by -999.
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      NONE
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    04/07/94 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument and result
!              ------------------------------------
!
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PMZM   ! result at flux localization
!$acc declare present(PA,PMZM)
!
!*       0.2   Declarations of local variables
!              -------------------------------
INTEGER :: IKU            ! upper bound in z direction of PA
!
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF MZM
!              ------------------
!
!!$IIU = SIZE(PA,1)
!!$IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!!$!
!!$JIJKOR  = 1 + IIU*IJU
!!$JIJKEND = IIU*IJU*IKU
!!$!
!!$!CDIR NODEP
!!$!OCL NOVREC
!!$DO JIJK=JIJKOR , JIJKEND
!!$   PMZM(JIJK,1,1) = 0.5*( PA(JIJK,1,1)+PA(JIJK-IIU*IJU,1,1) )
!!$END DO
!!$
!!$!
!!$!CDIR NODEP
!!$!OCL NOVREC
!!$DO JIJ=1,IIU*IJU
!!$   PMZM(JIJ,1,1)    = -999.
!!$END DO
!
!$acc kernels present(PA,PMZM)
PMZM(:,:,2:IKU) = 0.5*( PA(:,:,2:IKU)+PA(:,:,1:IKU-1) )
PMZM(:,:,1)    = -999.
!$acc end kernels
!
END SUBROUTINE MZM_DEVICE
!
!     ###############################
      SUBROUTINE DXF_DEVICE(PA,PDXF)
!     ###############################
!
!!****  *DXF* -  Shuman operator : finite difference operator in x direction
!!                                  for a variable at a flux side
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a finite difference 
!     along the x direction (I index) for a field PA localized at a x-flux
!     point (u point). The result is localized at a mass point.
!
!!**  METHOD
!!    ------ 
!!        The result PDXF(i,:,:) is defined by (PA(i+1,:,:)-PA(i,:,:))
!!        At i=size(PA,1), PDXF(i,:,:) are replaced by the values of PDXF,
!!    which are the right values in the x-cyclic case
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      Module MODD_PARAMETERS: declaration of parameter variables
!!        JPHEXT: define the number of marginal points out of the 
!!        physical domain along the horizontal directions.
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    05/07/94 
!!      Modification to include the periodic case 13/10/94 J.Stein 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE MODD_PARAMETERS
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument and result
!              ------------------------------------
!
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDXF   ! result at mass localization 
!$acc declare present(PA,PDXF)
                                              !
!*       0.2   Declarations of local variables
!              -------------------------------
!
INTEGER :: JI             ! Loop index in x direction
INTEGER :: IIU            ! upper bound in x direction of PA 
!             
INTEGER :: JJK,IJU,IKU
INTEGER :: JIJK,JIJKOR,JIJKEND
!             
!
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF DXF
!              ------------------
!
IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + JPHEXT
JIJKEND = IIU*IJU*IKU
!
!$acc kernels present(PA,PDXF)
!CDIR NODEP
!OCL NOVREC
DO JIJK=JIJKOR , JIJKEND
   PDXF(JIJK-1,1,1) = PA(JIJK,1,1) - PA(JIJK-1,1,1) 
END DO
!
!CDIR NODEP
!OCL NOVREC
DO JJK=1,IJU*IKU
   PDXF(IIU,JJK,1) = PDXF(2*JPHEXT,JJK,1) 
END DO
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE DXF_DEVICE
!
!     ###############################
      SUBROUTINE DXM_DEVICE(PA,PDXM)
!     ###############################
!
!!****  *DXM* -  Shuman operator : finite difference operator in x direction
!!                                  for a variable at a mass localization
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a finite difference 
!     along the x direction (I index) for a field PA localized at a mass
!     point. The result is localized at a x-flux point (u point).
!
!!**  METHOD
!!    ------ 
!!        The result PDXM(i,:,:) is defined by (PA(i,:,:)-PA(i-1,:,:))
!!    At i=1, PDXM(1,:,:) are replaced by the values of PDXM,
!!    which are the right values in the x-cyclic case. 
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      Module MODD_PARAMETERS: declaration of parameter variables
!!        JPHEXT: define the number of marginal points out of the 
!!        physical domain along the horizontal directions.
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    05/07/94 
!!      Modification to include the periodic case 13/10/94 J.Stein 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE MODD_PARAMETERS
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument and result
!              ------------------------------------
!
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDXM   ! result at flux side
!$acc declare present(PA,PDXM)
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
INTEGER :: JI             ! Loop index in x direction
INTEGER :: IIU            ! Size of the array in the x direction
!
!
INTEGER :: JJK,IJU,IKU
INTEGER :: JIJK,JIJKOR,JIJKEND
!
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF DXM
!              ------------------
!
IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + 1
JIJKEND = IIU*IJU*IKU
!
!$acc kernels present(PA,PDXM)
!CDIR NODEP
!OCL NOVREC
DO JIJK=JIJKOR , JIJKEND
   PDXM(JIJK,1,1) = PA(JIJK,1,1) - PA(JIJK-1,1,1)
END DO
!
!CDIR NODEP
!OCL NOVREC
DO JJK=1,IJU*IKU
   PDXM(1,JJK,1)    = PDXM(IIU-2*JPHEXT+1,JJK,1)
END DO
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE DXM_DEVICE
!
!     ###############################
      SUBROUTINE DYF_DEVICE(PA,PDYF)
!     ###############################
!
!!****  *DYF* -  Shuman operator : finite difference operator in y direction
!!                                  for a variable at a flux side
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a finite difference 
!     along the y direction (J index) for a field PA localized at a y-flux
!     point (v point). The result is localized at a mass point.
!
!!**  METHOD
!!    ------ 
!!        The result PDYF(:,j,:) is defined by (PA(:,j+1,:)-PA(:,j,:))
!!        At j=size(PA,2), PDYF(:,j,:) are replaced by the values of PDYM,
!!    which are the right values in the y-cyclic case
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      Module MODD_PARAMETERS: declaration of parameter variables
!!        JPHEXT: define the number of marginal points out of the 
!!        physical domain along the horizontal directions.
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    05/07/94 
!!      Modification to include the periodic case 13/10/94 J.Stein 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE MODD_PARAMETERS
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument and result
!              ------------------------------------
!
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDYF   ! result at mass localization 
!$acc declare present(PA,PDYF)
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
INTEGER :: JJ            ! Loop index in y direction
INTEGER :: IJU           ! upper bound in y direction of PA 
!
!          
INTEGER :: IIU,IKU
INTEGER :: JIJK,JIJKOR,JIJKEND
!            
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF DYF
!              ------------------
!
IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + IIU
JIJKEND = IIU*IJU*IKU
!
!$acc kernels present(PA,PDYF)
!CDIR NODEP
!OCL NOVREC
DO JIJK=JIJKOR , JIJKEND
   PDYF(JIJK-IIU,1,1)         = PA(JIJK,1,1)  -  PA(JIJK-IIU,1,1) 
END DO
!
PDYF(:,IJU,:)    = PDYF(:,2*JPHEXT,:)
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE DYF_DEVICE
!
!     ###############################
      SUBROUTINE DYM_DEVICE(PA,PDYM)
!     ###############################
!
!!****  *DYM* -  Shuman operator : finite difference operator in y direction
!!                                  for a variable at a mass localization
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a finite difference 
!     along the y direction (J index) for a field PA localized at a mass
!     point. The result is localized at a y-flux point (v point).
!
!!**  METHOD
!!    ------ 
!!        The result PDYM(:,j,:) is defined by (PA(:,j,:)-PA(:,j-1,:))
!!    At j=1, PDYM(:,1,:) are replaced by the values of PDYM,
!!    which are the right values in the y-cyclic case. 
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      Module MODD_PARAMETERS: declaration of parameter variables
!!        JPHEXT: define the number of marginal points out of the 
!!        physical domain along the horizontal directions.
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    05/07/94 
!!      Modification to include the periodic case 13/10/94 J.Stein 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE MODD_PARAMETERS
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument and result
!              ------------------------------------
!
REAL, DIMENSION(:,:,:), INTENT(IN)  :: PA     ! variable at mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDYM   ! result at flux side
!$acc declare present(PA,PDYM)
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
INTEGER :: JJ             ! Loop index in y direction
INTEGER :: IJU            ! Size of the array in the y direction
!
!    
INTEGER :: IIU,IKU
INTEGER :: JIJK,JIJKOR,JIJKEND
!     
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF DYM
!              ------------------
!
IIU=SIZE(PA,1)
IJU=SIZE(PA,2)
IKU=SIZE(PA,3)
!
JIJKOR  = 1 + IIU
JIJKEND = IIU*IJU*IKU
!
!$acc kernels present(PA,PDYM)
!CDIR NODEP
!OCL NOVREC
DO JIJK=JIJKOR , JIJKEND
   PDYM(JIJK,1,1)           = PA(JIJK,1,1)  -  PA(JIJK-IIU,1,1) 
END DO
!
PDYM(:,1,:)    =  PDYM(:,IJU-2*JPHEXT+1,:)
!$acc end kernels
!
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE DYM_DEVICE
!
!     ###############################
      SUBROUTINE DZF_DEVICE(KKA,KKU,KL,PA,PDZF)
!     ###############################
!
!!****  *DZF* -  Shuman operator : finite difference operator in z direction
!!                                  for a variable at a flux side
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a finite difference 
!     along the z direction (K index) for a field PA localized at a z-flux
!     point (w point). The result is localized at a mass point.
!
!!**  METHOD
!!    ------ 
!!        The result PDZF(:,:,k) is defined by (PA(:,:,k+1)-PA(:,:,k))
!!        At k=size(PA,3), PDZF(:,:,k) is defined by -999.
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      NONE
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    05/07/94 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument 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 flux side
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDZF     ! result at mass localization 
!$acc declare present(PA,PDZF)
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
INTEGER :: JK           ! Loop index in z direction
INTEGER :: IKU          ! upper bound in z direction of PA 
!
!           
INTEGER :: IIU,IJU
INTEGER :: JIJ
INTEGER :: JIJK,JIJKOR,JIJKEND
!         
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF DZF
!              ------------------
!
IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + IIU*IJU
JIJKEND = IIU*IJU*IKU
!
!$acc kernels present(PA,PDZF)
!CDIR NODEP
!OCL NOVREC
DO JIJK=JIJKOR , JIJKEND
   PDZF(JIJK-IIU*IJU,1,1) = PA(JIJK,1,1)-PA(JIJK-IIU*IJU,1,1)
END DO
!
!CDIR NODEP
!OCL NOVREC
DO JIJ=1,IIU*IJU
   PDZF(JIJ,1,IKU) = -999.
END DO
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE DZF_DEVICE
!
!     ###############################
      SUBROUTINE DZM_DEVICE(KKA,KKU,KL,PA,PDZM)
!     ###############################
!
!!****  *DZM* -  Shuman operator : finite difference operator in z direction
!!                                  for a variable at a mass localization
!!
!!    PURPOSE
!!    -------
!       The purpose of this function  is to compute a finite difference 
!     along the z direction (K index) for a field PA localized at a mass
!     point. The result is localized at a z-flux point (w point).
!
!!**  METHOD
!!    ------ 
!!        The result PDZM(:,j,:) is defined by (PA(:,:,k)-PA(:,:,k-1))
!!        At k=1, PDZM(:,:,k) is defined by -999.
!!    
!!
!!    EXTERNAL
!!    --------
!!      NONE
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!      NONE
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation of Meso-NH (SHUMAN operators)
!!      Technical specifications Report of The Meso-NH (chapters 3)  
!!
!!
!!    AUTHOR
!!    ------
!!	V. Ducrocq       * Meteo France *
!!
!!    MODIFICATIONS
!!    -------------
!!      Original    05/07/94 
!!                   optimisation                 20/08/00 J. Escobar
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
IMPLICIT NONE
!
!*       0.1   Declarations of argument 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 mass localization
REAL, DIMENSION(:,:,:), INTENT(OUT) :: PDZM     ! result at flux side
!$acc declare present(PA,PDZM)
!
!*       0.2   Declarations of local variables
!              -------------------------------
!
INTEGER :: JK            ! Loop index in z direction
INTEGER :: IKU           ! upper bound in z direction of PA
!
!         
INTEGER :: IIU,IJU
INTEGER :: JIJ
INTEGER :: JIJK,JIJKOR,JIJKEND
!           
!-------------------------------------------------------------------------------
!
!*       1.    DEFINITION OF DZM
!              ------------------
!
IIU = SIZE(PA,1)
IJU = SIZE(PA,2)
IKU = SIZE(PA,3)
!
JIJKOR  = 1 + IIU*IJU
JIJKEND = IIU*IJU*IKU
!
!$acc kernels present(PA,PDZM)
!CDIR NODEP
!OCL NOVREC
!DO JIJK=JIJKOR , JIJKEND
!   PDZM(JIJK,1,1) = PA(JIJK,1,1)-PA(JIJK-IIU*IJU,1,1)
!END DO
PDZM(:,:,2:) = PA(:,:,2:)-PA(:,:,:IKU-1)
!
!CDIR NODEP
!OCL NOVREC
!DO JIJ=1,IIU*IJU
!   PDZM(JIJ,1,1)    = -999.
!END DO
PDZM(:,:,1) = -999.
!$acc end kernels
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE DZM_DEVICE