Skip to content
Snippets Groups Projects
advecuvw_4th.f90 16.2 KiB
Newer Older
!-----------------------------------------------------------------
!--------------- special set of characters for RCS information
!-----------------------------------------------------------------
! $Source: /home/cvsroot/MNH-VX-Y-Z/src/MNH/advecuvw_4th.f90,v $ $Revision: 1.1.4.1 $
! MASDEV4_7 adiab 2006/05/18 13:07:25
!-----------------------------------------------------------------
!     ###########################
      MODULE MODI_ADVECUVW_4TH
!     ###########################
!
INTERFACE
!
      SUBROUTINE ADVECUVW_4TH ( HLBCX, HLBCY, PRUCT, PRVCT, PRWCT,           &
                                PUT, PVT, PWT, PRUS, PRVS, PRWS, TPHALO2LIST )              
!
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
!
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
!
END SUBROUTINE ADVECUVW_4TH
!
END INTERFACE
!
END MODULE MODI_ADVECUVW_4TH
!
!
!     ######################################################################
      SUBROUTINE ADVECUVW_4TH ( HLBCX, HLBCY, PRUCT, PRVCT, PRWCT,           &
                                PUT, PVT, PWT, PRUS, PRVS, PRWS, TPHALO2LIST ) 

        USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
        USE MODE_MNH_ZWORK  , ONLY : ZT3D, MNH_GET_ZT3D , MNH_REL_ZT3D 
        USE MODE_MNH_ZWORK, ONLY : IIU,IJU,IKU

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
!
!
REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRUCT  ! contravariant
REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRVCT  !  components
REAL, DIMENSION(:,:,:),   INTENT(IN)    :: PRWCT  ! of momentum
!$acc declare present(PRUCT,PRVCT,PRWCT)
!
REAL, DIMENSION(:,:,:),   INTENT(IN) :: PUT, PVT, PWT        ! Variables at t
!$acc declare pcopyin(PUT,PVT,PWT)
!
REAL, DIMENSION(:,:,:),   INTENT(INOUT) :: PRUS, PRVS, PRWS     ! Source terms
!$acc declare pcopy(PRUS, PRVS, PRWS)
!
TYPE(HALO2LIST_ll), POINTER :: TPHALO2LIST ! list for diffusion
!
!*       0.2   Declarations of local variables :
!
INTEGER:: IIB,IJB        ! Begining useful area  in x,y,z directions
INTEGER:: IIE,IJE        ! End useful area in x,y,z directions
!
TYPE(HALO2LIST_ll), POINTER :: TZHALO2LIST

INTEGER   :: IZMEANX, IZMEANY, IZTEMP1,IZTEMP2,IZTEMP3,IZTEMP4 

         CALL  MNH_GET_ZT3D(IZMEANX, IZMEANY,IZTEMP1,IZTEMP2,IZTEMP3,IZTEMP4 )

         CALL  ADVECUVW_4TH_D ( IIU,IJU,IKU,HLBCX, HLBCY, &
                              & PRUCT, PRVCT, PRWCT,        &
                              & PUT, PVT, PWT, PRUS, PRVS, PRWS, TPHALO2LIST, & 
                              & ZT3D(:,:,:,IZMEANX),ZT3D(:,:,:,IZMEANY), & 
                              & ZT3D(:,:,:,IZTEMP1),ZT3D(:,:,:,IZTEMP2), & 
                              & ZT3D(:,:,:,IZTEMP3),ZT3D(:,:,:,IZTEMP4)  & 
                              & )

         CALL  MNH_REL_ZT3D(IZMEANX, IZMEANY, IZTEMP1,IZTEMP2,IZTEMP3,IZTEMP4)

CONTAINS

      SUBROUTINE ADVECUVW_4TH_D ( IIU,IJU,IKU,HLBCX, HLBCY, &
                              & PRUCT, PRVCT, PRWCT,           &
                              & PUT, PVT, PWT, PRUS, PRVS, PRWS, TPHALO2LIST , &
                              & ZMEANX, ZMEANY, ZTEMP1,ZTEMP2,ZTEMP3,ZTEMP4 ) 
         
!     ######################################################################
!
!!****  *ADVECUVW_4TH * - routine to compute the 4th order centered
!!                           advection tendency of momentum (U,V,W)
!!
!!    PURPOSE
!!    -------
!!      The purpose of this routine is to call the ADVEC_4TH_ORDER_ALGO
!!    routine for the horizontal advection and the MZM4 and MZF4 functions for
!!    the vertical advection of momentum. The code is 
!!    parallelized and works for various boundary conditions.
!!
!!**  METHOD
!!    ------
!!      For each wind component the ADVECUVW_4TH routine calls
!!    the ADVEC_4TH_ORDER_ALGO routine which computes the numerical advection
!!    of any 3D field.
!!      The following variables are passed as argument to ADVEC_4TH_ORDER_ALGO :
!!
!!    -- The variable at t
!!    -- The second layer of the halo of the field at t
!!    -- The horizontal advection fluxes
!!    -- The localisation on the model grid :
!!
!!        IGRID = 1 for mass grid point
!!        IGRID = 2 for U grid point
!!        IGRID = 3 for V grid point
!!        IGRID = 4 for W grid point
!!
!!    EXTERNAL
!!    --------
!!      BUDGET      : Stores the different budget components
!!                    (not used in current version)
!!      ADVEC_4TH_ORDER_ALGO : computes the horizontal advection fluxes
!!      MZF4 and MZM4 : computes the vertical advection fluxes
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!    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
!!         NBUPROCCTR   : process counter used for each budget variable
!!         Switches for budgets activations:
!!
!!         LBU_RU       : logical for budget of RU (wind component along x)
!!
!!         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
!!         LBU_RTH      : logical for budget of RTH (potential temperature)
!!                        .TRUE. = budget of RTH
!!                        .FALSE. = no budget of RTH
!!         LBU_RTKE     : logical for budget of RTKE (turbulent kinetic energy)
!!                        .TRUE. = budget of RTKE
!!                        .FALSE. = no budget of RTKE
!!         LBU_RRV      : logical for budget of RRV (water vapor)
!!                        .TRUE. = budget of RRV
!!                        .FALSE. = no budget of RRV
!!         LBU_RRC      : logical for budget of RRC (cloud water)
!!                        .TRUE. = budget of RRC
!!                        .FALSE. = no budget of RRC
!!         LBU_RRR      : logical for budget of RRR (rain water)
!!                        .TRUE. = budget of RRR
!!                        .FALSE. = no budget of RRR
!!         LBU_RRI      : logical for budget of RRI (ice)
!!                        .TRUE. = budget of RRI
!!                        .FALSE. = no budget of RRI
!!         LBU_RRS      : logical for budget of RRS (snow)
!!                        .TRUE. = budget of RRS
!!                        .FALSE. = no budget of RRS
!!         LBU_RRG      : logical for budget of RRG (graupel)
!!                        .TRUE. = budget of RRG
!!                        .FALSE. = no budget of RRG
!!         LBU_RRH      : logical for budget of RRH (hail)
!!                        .TRUE. = budget of RRH
!!                        .FALSE. = no budget of RRH
!!         LBU_RSV      : logical for budget of RSVx (scalar variable)
!!                        .TRUE. = budget of RSVx
!!                        .FALSE. = no budget of RSVx
!!
!!    MODULE MODD_ARGSLIST
!!         HALO2LIST_ll : type for a list of "HALO2_lls"
!!
!!    REFERENCE
!!    ---------
!!      Book2 of documentation ( routine ADVECUVW_4TH )
!!
!!    AUTHOR
!!    ------
!!      J.-P. Pinty      * Laboratoire d'Aerologie*
!!
!!    MODIFICATIONS
!!    -------------
!!      Original   25/10/05
!!
!-------------------------------------------------------------------------------
!
!*       0.    DECLARATIONS
!              ------------
!
USE MODE_ll
!
USE MODD_PARAMETERS
USE MODD_CONF
USE MODD_BUDGET
USE MODD_ARGSLIST_ll, ONLY : HALO2LIST_ll
!
USE MODI_SHUMAN
USE MODI_BUDGET
!
USE MODI_ADVEC_4TH_ORDER_AUX
!
IMPLICIT NONE
!
!*       0.1   Declarations of dummy arguments :
!
!
INTEGER                        , INTENT(IN) :: IIU,IJU,IKU
CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCX ! X direction LBC type
CHARACTER (LEN=4), DIMENSION(2), INTENT(IN) :: HLBCY ! Y direction LBC type
!
!
REAL, DIMENSION(IIU,IJU,IKU),   INTENT(IN)    :: PRUCT  ! contravariant
REAL, DIMENSION(IIU,IJU,IKU),   INTENT(IN)    :: PRVCT  !  components
REAL, DIMENSION(IIU,IJU,IKU),   INTENT(IN)    :: PRWCT  ! of momentum
!$acc declare present(PRUCT,PRVCT,PRWCT)
REAL, DIMENSION(IIU,IJU,IKU),   INTENT(IN) :: PUT, PVT, PWT        ! Variables at t
!$acc declare pcopyin(PUT,PVT,PWT)
REAL, DIMENSION(IIU,IJU,IKU),   INTENT(INOUT) :: PRUS, PRVS, PRWS     ! Source terms
!$acc declare pcopy(PRUS, PRVS, PRWS)
!
TYPE(HALO2LIST_ll), POINTER :: TPHALO2LIST ! list for diffusion
!
!*       0.2   Declarations of local variables :
!
INTEGER:: IIB,IJB        ! Begining useful area  in x,y,z directions
INTEGER:: IIE,IJE        ! End useful area in x,y,z directions
!
TYPE(HALO2LIST_ll), POINTER :: TZHALO2LIST
!
INTEGER :: IGRID ! localisation on the model grid
REAL, DIMENSION(IIU,IJU,IKU) :: ZMEANX, ZMEANY ! fluxes
!acc declare present(ZMEANX, ZMEANY)
!
REAL, DIMENSION(IIU,IJU,IKU) :: ZTEMP1,ZTEMP2,ZTEMP3,ZTEMP4
!$acc declare present(ZTEMP1,ZTEMP2,ZTEMP3,ZTEMP4) 

INTEGER                                              :: II

!
#define dxm(PDXM,PA) PDXM(2:IIU,:,:)   = PA(2:IIU,:,:) - PA(1:IIU-1,:,:)       ; PDXM(1,:,:) = PDXM(IIU-2*JPHEXT+1,:,:) ! DXM(PDXM,PA)
#define mxf(PMXF,PA) PMXF(1:IIU-1,:,:) = 0.5*( PA(2:IIU,:,:)+PA(1:IIU-1,:,:) ) ; PMXF(IIU,:,:) = PMXF(2*JPHEXT,:,:) ! MXF(PMXF,PA)
#define mxm(PMXM,PA) PMXM(2:IIU,:,:) = 0.5*( PA(2:IIU,:,:)+PA(1:IIU-1,:,:) ) ; PMXM(1,:,:) = PMXM(IIU-2*JPHEXT+1,:,:) ! MXM(PMXM,PA)
#define dyf(PDYF,PA) PDYF(:,1:IJU-1,:) = PA(:,2:IJU,:) - PA(:,1:IJU-1,:) ; PDYF(:,IJU,:) = PDYF(:,2*JPHEXT,:) ! DYF(PDYF,PA)
#define dzf(PDZF,PA) PDZF(:,:,1:IKU-1) = PA(:,:,2:IKU) - PA(:,:,1:IKU-1) ; PDZF(:,:,IKU) = -999. ! DZF(PDZF,PA)
#define mzm4(PMZM4,PA) PMZM4(:,:,3:IKU-1) = (7.0*( PA(:,:,3:IKU-1)+PA(:,:,2:IKU-2) ) - (PA(:,:,4:IKU)+PA(:,:,1:IKU-3) ) )/12.0 ; \
  PMZM4(:,:,2) = 0.5*( PA(:,:,2)+PA(:,:,1) ) ; PMZM4(:,:,IKU) = 0.5*( PA(:,:,IKU)+PA(:,:,IKU-1) ) ; PMZM4(:,:,1) = -999.
#define mym(PMYM,PA) PMYM(:,2:IJU,:) = 0.5*( PA(:,2:IJU,:)+PA(:,1:IJU-1,:) ) ; PMYM(:,1,:) = PMYM(:,IJU-2*JPHEXT+1,:) !  MYM(PMYM,PA)
#define dxf(PDXF,PA) PDXF(1:IIU-1,:,:) = PA(2:IIU,:,:) - PA(1:IIU-1,:,:) ; PDXF(IIU,:,:)    = PDXF(2*JPHEXT,:,:) ! DXF(PDXF,PA)
#define myf(PMYF,PA) PMYF(:,1:IJU-1,:) = 0.5*( PA(:,1:IJU-1,:)+PA(:,2:IJU,:) ) ; PMYF(:,IJU,:) = PMYF(:,2*JPHEXT,:) ! MYF(PMYF,PA)
#define dym(PDYM,PA) PDYM(:,2:IJU,:) = PA(:,2:IJU,:) - PA(:,1:IJU-1,:) ; PDYM(:,1,:) = PDYM(:,IJU-2*JPHEXT+1,:) ! DYM(PDYM,PA)
#define mzm(PMZM,PA) PMZM(:,:,2:IKU) = 0.5*( PA(:,:,2:IKU)+PA(:,:,1:IKU-1) ) ; PMZM(:,:,1)    = -999. !  MZM(PMZM,PA)
#define mzf(PMZF,PA) PMZF(:,:,1:IKU-1) = 0.5*( PA(:,:,1:IKU-1)+PA(:,:,2:IKU) ) ; PMZF(:,:,IKU) = -999. ! MZF(PMZF,PA)
#define dzm(PDZM,PA) PDZM(:,:,2:IKU) = PA(:,:,2:IKU) - PA(:,:,1:IKU-1) ; PDZM(:,:,1) = -999. !  DZM(PDZM,PA)
#define mzf4(PMZF4,PA) PMZF4(:,:,2:IKU-2) = (7.0*( PA(:,:,3:IKU-1)+PA(:,:,2:IKU-2) ) - (PA(:,:,4:IKU)+PA(:,:,1:IKU-3) ) )/12.0 ; \
  PMZF4(:,:,1) = 0.5*( PA(:,:,2)+PA(:,:,1) ) ; PMZF4(:,:,IKU-1) = 0.5*( PA(:,:,IKU)+PA(:,:,IKU-1) ) ; PMZF4(:,:,IKU) = -999.
!
!-------------------------------------------------------------------------------
!
!*       1.     COMPUTES THE DOMAIN DIMENSIONS
!               ------------------------------
!$acc data 
! create(ZMEANX, ZMEANY)
!
CALL GET_INDICE_ll(IIB,IJB,IIE,IJE)
!!$CALL GET_DIM_EXT_ll('B',IIU,IJU)     
!!$IKU=SIZE(PUT,3)
!
!-------------------------------------------------------------------------------
!
!*       2.     CALL THE ADVEC_4TH_ORDER_ALGO ROUTINE FOR MOMENTUM
!               --------------------------------------------------
!
IGRID = 2
IF(NHALO == 1) THEN
  TZHALO2LIST => TPHALO2LIST
  CALL ADVEC_4TH_ORDER_ALGO(HLBCX, HLBCY, PUT, IGRID, ZMEANX, ZMEANY, &
                            TZHALO2LIST%HALO2 )
ELSE
  CALL ADVEC_4TH_ORDER_ALGO(HLBCX, HLBCY, PUT, IGRID, ZMEANX, ZMEANY)
ENDIF
!
! pcopy(prus) pcopyin(pruct,ZMEANX) create(ZTEMP1,ZTEMP2,ZTEMP3)
!!$PRUS(:,:,:) = PRUS(:,:,:)                          &
!!$             -DXM( MXF(PRUCT(:,:,:))*ZMEANX(:,:,:) ) 

!$acc kernels   pcopyin(ZMEANX)       
mxf(ZTEMP1,PRUCT)         
ZTEMP2 = ZTEMP1 * ZMEANX
dxm(ZTEMP3,ZTEMP2)
PRUS(:,:,:) = PRUS(:,:,:) -  ZTEMP3
!$acc end kernels   
         
IF (LBUDGET_U)  CALL BUDGET (PRUS,1,'ADVX_BU_RU')

!
!!$PRUS(:,:,:) = PRUS(:,:,:)                          &
!!$             -DYF( MXM(PRVCT(:,:,:))*ZMEANY(:,:,:) )

!$acc kernels   pcopyin(ZMEANY) 
mxm(ZTEMP1,PRVCT)
ZTEMP2 = ZTEMP1 * ZMEANY
dyf(ZTEMP3,ZTEMP2)
PRUS(:,:,:) = PRUS(:,:,:) - ZTEMP3
!$acc end kernels            

IF (LBUDGET_U)  CALL BUDGET (PRUS,1,'ADVY_BU_RU')
!
!!$PRUS(:,:,:) = PRUS(:,:,:)                             &
!!$             -DZF( MXM(PRWCT(:,:,:))*MZM4(PUT(:,:,:)) )

!$acc kernels            
mzm4(ZTEMP1,PUT)
mxm(ZTEMP2,PRWCT)            
ZTEMP3 = ZTEMP1 * ZTEMP2
dzf(ZTEMP4,ZTEMP3)
PRUS(:,:,:) = PRUS(:,:,:) - ZTEMP4
!$acc end kernels            

IF (LBUDGET_U)  CALL BUDGET (PRUS,1,'ADVZ_BU_RU')
!
!
IGRID = 3
IF(NHALO == 1) THEN
  TZHALO2LIST => TZHALO2LIST%NEXT
  CALL ADVEC_4TH_ORDER_ALGO(HLBCX, HLBCY, PVT, IGRID, ZMEANX, ZMEANY, &
                            TZHALO2LIST%HALO2 )
ELSE
  CALL ADVEC_4TH_ORDER_ALGO(HLBCX, HLBCY, PVT, IGRID, ZMEANX, ZMEANY)
ENDIF
!
!!$PRVS(:,:,:) = PRVS(:,:,:)                          &
!!$             -DXF( MYM(PRUCT(:,:,:))*ZMEANX(:,:,:) ) 

!$acc kernels  pcopyin(ZMEANX)          
mym(ZTEMP1,PRUCT)
ZTEMP2 = ZTEMP1 * ZMEANX
dxf(ZTEMP3,ZTEMP2) 
PRVS(:,:,:) = PRVS(:,:,:) -  ZTEMP3   
!$acc end kernels   
                     
IF (LBUDGET_V)  CALL BUDGET (PRVS,2,'ADVX_BU_RV')
!
!!$PRVS(:,:,:) = PRVS(:,:,:)                          &
!!$             -DYM( MYF(PRVCT(:,:,:))*ZMEANY(:,:,:) )

!$acc kernels   pcopyin(ZMEANY)    
myf(ZTEMP1,PRVCT)
ZTEMP2 = ZTEMP1 * ZMEANY
dym(ZTEMP3,ZTEMP2)
PRVS(:,:,:) = PRVS(:,:,:) - ZTEMP3
!$acc end kernels   

IF (LBUDGET_V)  CALL BUDGET (PRVS,2,'ADVY_BU_RV')
!
!!$PRVS(:,:,:) = PRVS(:,:,:)                             &
!!$             -DZF( MYM(PRWCT(:,:,:))*MZM4(PVT(:,:,:)) )

!$acc kernels    
mym(ZTEMP1,PRWCT)
mzm4(ZTEMP2,PVT)
ZTEMP3 = ZTEMP1 * ZTEMP2
dzf(ZTEMP4,ZTEMP3)
PRVS(:,:,:) = PRVS(:,:,:) - ZTEMP4
!$acc end kernels

IF (LBUDGET_V)  CALL BUDGET (PRVS,2,'ADVZ_BU_RV')
!
!
IGRID = 4
!
IF(NHALO == 1) THEN
  TZHALO2LIST => TZHALO2LIST%NEXT
  CALL ADVEC_4TH_ORDER_ALGO(HLBCX, HLBCY, PWT, IGRID, ZMEANX, ZMEANY, &
                            TZHALO2LIST%HALO2 )
ELSE
  CALL ADVEC_4TH_ORDER_ALGO(HLBCX, HLBCY, PWT, IGRID, ZMEANX, ZMEANY)
ENDIF
!
!!$PRWS(:,:,:) = PRWS(:,:,:)                          &
!!$             -DXF( MZM(PRUCT(:,:,:))*ZMEANX(:,:,:) ) 

!$acc kernels    pcopyin(ZMEANX) 
mzm(ZTEMP1,PRUCT)
ZTEMP2 = ZTEMP1 * ZMEANX
dxf(ZTEMP3,ZTEMP2)
PRWS(:,:,:) = PRWS(:,:,:)  - ZTEMP3
!$acc end kernels

IF (LBUDGET_W)  CALL BUDGET (PRWS,3,'ADVX_BU_RW')
!
!!$PRWS(:,:,:) = PRWS(:,:,:)                          &
!!$             -DYF( MZM(PRVCT(:,:,:))*ZMEANY(:,:,:) )

!$acc kernels   pcopyin(ZMEANY)  
mzm(ZTEMP1,PRVCT)
ZTEMP2 = ZTEMP1 * ZMEANY
dyf(ZTEMP3,ZTEMP2)
PRWS(:,:,:) = PRWS(:,:,:) - ZTEMP3
!$acc end kernels

IF (LBUDGET_W)  CALL BUDGET (PRWS,3,'ADVY_BU_RW')
!
!!$PRWS(:,:,:) = PRWS(:,:,:)                             &
!!$             -DZM( MZF(PRWCT(:,:,:))*MZF4(PWT(:,:,:)) )

!$acc kernels
mzf(ZTEMP1,PRWCT)
mzf4(ZTEMP2,PWT)
ZTEMP1 = ZTEMP1 * ZTEMP2 
dzm(ZTEMP4,ZTEMP1)
PRWS(:,:,:) = PRWS(:,:,:) - ZTEMP4
!$acc end kernels
IF (LBUDGET_W)  CALL BUDGET (PRWS,3,'ADVZ_BU_RW')
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE ADVECUVW_4TH_D