scale_meshfieldcomm_cubedspheredom3d.F90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
!> module FElib / Data / Communication in 3D cubed-sphere domain
!!
!! @par Description
!!      A module to manage data communication with 3D cubed-sphere domain for element-based methods
!!
!! @author Yuta Kawai, Team SCALE
!<
#include "scaleFElib.h"
module scale_meshfieldcomm_cubedspheredom3d
 
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
  use scale_io
  use scale_prof
 
  use scale_element_base, only: elementbase, elementbase3d
  use scale_mesh_cubedspheredom3d, only: meshcubedspheredom3d
  use scale_meshfield_base, only: meshfield3d
  use scale_meshfieldcomm_base, only: &
    meshfieldcommbase,                               &
    meshfieldcommbase_init, meshfieldcommbase_final, &
    meshfieldcommbase_extract_bounddata,             &
    meshfieldcommbase_extract_bounddata2,            &
    meshfieldcommbase_extract_bounddata_2,           &
    meshfieldcommbase_extract_bounddata_3,           &
    meshfieldcommbase_set_bounddata,                 &
    meshfieldcommbase_set_bounddata_3,               &
    meshfieldcontainer
 
  use scale_localmesh_2d, only: localmesh2d
  use scale_localmesh_3d, only: localmesh3d
 
   
  !-----------------------------------------------------------------------------
  implicit none
  private
 
  !-----------------------------------------------------------------------------
  !
  !++ Public type & procedure
  ! 
 
  !> Derived type to represent covariant vector components
  type :: veccovariantcomp
    type(MeshField3D), pointer :: u1 => null()
    type(MeshField3D), pointer :: u2 => null()
  end type    
 
  !> Base derived type to manage data communication with 3D cubed-sphere domain
  type, public, extends(meshfieldcommbase) :: meshfieldcommcubedspheredom3d
    class(meshcubedspheredom3d), pointer :: mesh3d                          !< Pointer to an object representing 3D cubed-sphere computational mesh
    type(veccovariantcomp), allocatable :: vec_covariant_comp_ptrlist(:)
 
    integer :: halosize_h1d !< Halo size for 1D horizontal direction
    integer :: halosize_v   !< Halo size for vertical direction  
  contains
    procedure, public :: init   => meshfieldcommcubedspheredom3d_init
    procedure, public :: put    => meshfieldcommcubedspheredom3d_put
    procedure, public :: get    => meshfieldcommcubedspheredom3d_get
    procedure, public :: exchange => meshfieldcommcubedspheredom3d_exchange  
    procedure, public :: setcovariantvec => meshfieldcommcubedspheredom3d_set_covariantvec
    procedure, public :: final => meshfieldcommcubedspheredom3d_final
  end type meshfieldcommcubedspheredom3d
 
  !-----------------------------------------------------------------------------
  !
  !++ Public parameters & variables
  !
  
  !-----------------------------------------------------------------------------
  !
  !++ Private type & procedure
  !
  private :: post_exchange_core
  private :: push_localsendbuf
 
  !-----------------------------------------------------------------------------
  !
  !++ Private parameters & variables
  !
  integer, parameter :: comm_face_num = 6  !< Number of faces with data communication
 
contains
!> Initialize an object to manage data communication with 3D cubed-sphere computational mesh
  subroutine meshfieldcommcubedspheredom3d_init( this, &
    sfield_num, hvfield_num, htensorfield_num, mesh3d, &
    haloSize_h1D, haloSize_v )
 
    use scale_meshutil_cubedsphere3d, only: meshutilcubedsphere3d_genpatchboundarymap_wide
    implicit none
    
    class(meshfieldcommcubedspheredom3d), intent(inout), target :: this
    integer, intent(in) :: sfield_num                         !< Number of scalar fields
    integer, intent(in) :: hvfield_num                        !< Number of horizontal vector fields
    integer, intent(in) :: htensorfield_num                   !< Number of horizontal vector fields
    class(meshcubedspheredom3d), intent(in), target :: mesh3d !< Object to manage a 3D cubed-sphere computational mesh
    integer, intent(in), optional :: halosize_h1d             !< Halo size for 1D horizontal direction
    integer, intent(in), optional :: halosize_v               !< Halo size for vertical direction
 
    type(localmesh3d), pointer :: lcmesh
    type(elementbase3d), pointer :: elem
    integer :: n
    integer :: nnode_lcmeshface(comm_face_num,mesh3d%local_mesh_num)
    !-----------------------------------------------------------------------------
    
    this%mesh3d => mesh3d
    lcmesh => mesh3d%lcmesh_list(1)
    elem => lcmesh%refElem3D
 
    !-
    if ( present(halosize_h1d) ) then
      this%haloSize_h1D = halosize_h1d
    else
      this%haloSize_h1D = 1
    end if
    if ( present(halosize_v) ) then
      this%haloSize_v = halosize_v
    else
      this%haloSize_v = 1
    end if
 
    !-
    allocate( this%VMapB_size(this%mesh3d%LOCAL_MESH_NUM) )
 
    this%bufsize_per_field =  2*(lcmesh%NeX + lcmesh%NeY)*lcmesh%NeZ*elem%Nfp_h*this%haloSize_h1D &
                            + 2*lcmesh%NeX*lcmesh%NeY*elem%Nfp_v*this%haloSize_v
 
    do n=1, this%mesh3d%LOCAL_MESH_NUM
      lcmesh => this%mesh3d%lcmesh_list(n)
      nnode_lcmeshface(:,n) = &
          (/ lcmesh%NeX, lcmesh%NeY, lcmesh%NeX, lcmesh%NeY, 0, 0 /) * lcmesh%NeZ * lcmesh%refElem3D%Nfp_h*this%haloSize_h1D &
        + (/ 0, 0, 0, 0, 1, 1 /) * lcmesh%NeX*lcmesh%NeY * lcmesh%refElem3D%Nfp_v*this%haloSize_v
    end do
 
    call meshfieldcommbase_init( this, sfield_num, hvfield_num, htensorfield_num, this%bufsize_per_field, comm_face_num, nnode_lcmeshface, mesh3d )  
  
    !-
    if ( this%haloSize_h1D > 1 .or. this%haloSize_v > 1) then
      this%use_vmap_wide_flag = .true.
      allocate( this%VMapB2(this%bufsize_per_field) )
 
      lcmesh => this%mesh3d%lcmesh_list(1)      
      call meshutilcubedsphere3d_genpatchboundarymap_wide( this%VMapB2, &
        lcmesh%VMapB, this%haloSize_h1D, this%haloSize_v,    &
        lcmesh%NeX, lcmesh%NeY, lcmesh%NeZ,                  &
        elem%Nfp_h, elem%Nfp_v, elem%Nnode_h1D, elem%Nnode_v )
      !$acc enter data copyin(this%VMapB2)
    else
      this%use_vmap_wide_flag = .false.
    end if
    
    do n=1, this%mesh3d%LOCAL_MESH_NUM
      lcmesh => this%mesh3d%lcmesh_list(n)
      if ( this%use_vmap_wide_flag ) then
        this%VMapB_size(n) = size(this%VMapB2)        
      else
        this%VMapB_size(n) = size(lcmesh%VMapB)
      end if
    end do
    !$acc enter data copyin(this%VMapB_size)
 
    !-
    if (hvfield_num > 0) then
      allocate( this%vec_covariant_comp_ptrlist(hvfield_num) )
    end if
 
    return
  end subroutine meshfieldcommcubedspheredom3d_init
 
!> Finalize an object to manage data communication with 3D cubed-sphere computational mesh
  subroutine meshfieldcommcubedspheredom3d_final( this )
    implicit none
    class(meshfieldcommcubedspheredom3d), intent(inout) :: this
    !-----------------------------------------------------------------------------
 
    if ( this%hvfield_num > 0 ) then
      deallocate( this%vec_covariant_comp_ptrlist )
    end if
 
    call meshfieldcommbase_final( this )
    return
  end subroutine meshfieldcommcubedspheredom3d_final
 
  !> Register objects to manage covariant vector components
  subroutine meshfieldcommcubedspheredom3d_set_covariantvec( &
    this, hvfield_ID, u1, u2  )
    implicit none
    class(meshfieldcommcubedspheredom3d), intent(inout) :: this
    integer, intent(in) :: hvfield_id
    type(meshfield3d), intent(in), target :: u1
    type(meshfield3d), intent(in), target :: u2
    !--------------------------------------------------------------
 
    this%vec_covariant_comp_ptrlist(hvfield_id)%u1 => u1
    this%vec_covariant_comp_ptrlist(hvfield_id)%u2 => u2
    return
  end subroutine meshfieldcommcubedspheredom3d_set_covariantvec
 
!> Put field data into temporary buffers
  subroutine meshfieldcommcubedspheredom3d_put(this, field_list, varid_s)
    implicit none
    class(meshfieldcommcubedspheredom3d), intent(inout) :: this
    type(meshfieldcontainer), intent(in) :: field_list(:)  !< Array of objects with 3D mesh field
    integer, intent(in) :: varid_s                         !< Start index with variables when field_list(1) is written to buffers for data communication
  
    integer :: i
    integer :: n
    type(localmesh3d), pointer :: lcmesh
    !-----------------------------------------------------------------------------
    ! call PROF_rapstart( 'comm_put', 2)
 
    if ( this%use_vmap_wide_flag ) then
      call meshfieldcommbase_extract_bounddata_3( field_list, 3, varid_s, this%mesh3d%lcmesh_list, this%VMapB2, this%VMapB_size(1), & ! (in)
        this%send_buf ) ! (out)
    else
      call meshfieldcommbase_extract_bounddata_2( &
        field_list, 3, varid_s, this%mesh3d%lcmesh_list, size(this%mesh3d%lcmesh_list(1)%VMapB), & !(in)
        this%send_buf ) ! (out)
    end if
    ! call PROF_rapend( 'comm_put', 2)
    return
  end subroutine meshfieldcommcubedspheredom3d_put
 
!> Extract field data from temporary buffers
  subroutine meshfieldcommcubedspheredom3d_get(this, field_list, varid_s)
    use scale_meshfieldcomm_base, only: &
      meshfieldcommbase_wait_core
    implicit none
    
    class(meshfieldcommcubedspheredom3d), intent(inout) :: this
    type(meshfieldcontainer), intent(inout) :: field_list(:)  !< Array of objects with 3D mesh field
    integer, intent(in) :: varid_s                            !< Start index with variables when field_list(1) is written to buffers for data communication
 
    integer :: i
    integer :: n
    integer :: ke, ke2d
    integer :: p
    type(localmesh3d), pointer :: lcmesh
    class(elementbase3d), pointer :: elem
 
    integer :: varnum
    integer :: varid_e
    integer :: varid_vec_s
 
    integer, allocatable :: indexh2dto3d(:)
    real(rp), allocatable :: g_ij(:,:,:,:)
 
    integer :: ke_z
    integer, allocatable :: vmapb(:)
    !-----------------------------------------------------------------------------
 
    ! call PROF_rapstart( 'comm_get', 2)
 
    varnum = size(field_list) 
 
    !--
    if ( this%call_wait_flag_sub_get ) then
      ! This workflow should be reconsidered in near future.
      ! The coordinate convresions for horizontal vector/tensor fields in post_exchange_core are applied for the received data (this%recv_buf).
      ! However, MeshFieldCommBase_wait_core does not store the communication data into this%recv_buf, but directly into the field_list.
      call meshfieldcommbase_wait_core( this, this%commdata_list, &
        field_list, 3, varid_s, this%mesh3d%lcmesh_list )
      call post_exchange_core( this )
    else
      do i=1, varnum
      do n=1, this%mesh3d%LOCAL_MESH_NUM
        lcmesh => this%mesh3d%lcmesh_list(n)
        if ( this%use_vmap_wide_flag ) then
          call meshfieldcommbase_set_bounddata_3( this%recv_buf(:,varid_s+i-1,n), lcmesh%refElem, lcmesh, & ! (in)
            this%VMapB2, this%VMapB_size(1),                                                              & ! (in)          
            field_list(i)%field3d%local(n)%val )                                                            !(out)
        else
          call meshfieldcommbase_set_bounddata( this%recv_buf(:,varid_s+i-1,n), lcmesh%refElem, lcmesh, & !(in)
            field_list(i)%field3d%local(n)%val )                                                          !(out)
        end if
      end do
      end do
    end if
 
    varid_e = varid_s + varnum - 1
    if ( varid_e > this%sfield_num ) then
      do i=1, this%hvfield_num
 
        varid_vec_s = this%sfield_num + 2*i - 1 
        if ( varid_vec_s > varid_e ) exit
 
        if (       associated(this%vec_covariant_comp_ptrlist(i)%u1 ) &
            .and.  associated(this%vec_covariant_comp_ptrlist(i)%u2 ) ) then
        
          do n=1, this%mesh3d%LOCAL_MESH_NUM
            lcmesh => this%mesh3d%lcmesh_list(n)
            elem => lcmesh%refElem3D
 
            allocate( g_ij(elem%Np,lcmesh%Ne,2,2), indexh2dto3d(elem%Np) )
            indexh2dto3d(:) = elem%IndexH2Dto3D(:)
 
            if ( this%use_vmap_wide_flag ) then
              allocate(vmapb(this%VMapB_size(n)))
              vmapb(:) = this%VMapB2(:)
            else
              allocate(vmapb(size(lcmesh%VMapB)))
              vmapb(:) = lcmesh%VMapB(:)
            end if
            !$acc enter data create(G_ij) copyin(VMapB, IndexH2Dto3D) async(1)
 
            !$omp parallel do private(ke2D)
            !$acc parallel loop collapse(2) async(1)
            do ke=lcmesh%NeS, lcmesh%NeE
            do p=1, elem%Np
              ke2d = lcmesh%EMap3Dto2D(ke)
              g_ij(p,ke,1,1) = lcmesh%G_ij(indexh2dto3d(p),ke2d,1,1)
              g_ij(p,ke,2,1) = lcmesh%G_ij(indexh2dto3d(p),ke2d,2,1)
              g_ij(p,ke,1,2) = lcmesh%G_ij(indexh2dto3d(p),ke2d,1,2)
              g_ij(p,ke,2,2) = lcmesh%G_ij(indexh2dto3d(p),ke2d,2,2)
            end do
            end do
 
            call set_boundary_data3d_u1u2( &
              this%recv_buf(:,varid_vec_s,n), this%recv_buf(:,varid_vec_s+1,n), & ! (in)
              this%bufsize_per_field, lcmesh%refElem3D, lcmesh, g_ij, vmapb,    & ! (in)
              this%vec_covariant_comp_ptrlist(i)%u1%local(n)%val,               & ! (out)
              this%vec_covariant_comp_ptrlist(i)%u2%local(n)%val                ) ! (out)
            !$acc exit data delete(G_ij, VMapB, IndexH2Dto3D) async(1)
            deallocate( g_ij, vmapb, indexh2dto3d )
          end do
       end if
      end do
    end if
 
    !$acc wait(1)
    ! call PROF_rapend( 'comm_get', 2)
    return
  end subroutine meshfieldcommcubedspheredom3d_get
 
  !> Exchange field data between neighboring MPI processes
  !!
  !! @param do_wait Flag whether MPI_waitall is called and move tmp data of LocalMeshCommData object to a recv buffer
!OCL SERIAL
  subroutine meshfieldcommcubedspheredom3d_exchange( this, do_wait )
    use scale_meshfieldcomm_base, only: &
      meshfieldcommbase_exchange_core,  &
      localmeshcommdata
    use scale_cubedsphere_coord_cnv, only: &
      cubedspherecoordcnv_cs2lonlatvec, &
      cubedspherecoordcnv_lonlat2csvec
    use scale_prc, only: prc_abort
    implicit none
  
    class(meshfieldcommcubedspheredom3d), intent(inout), target :: this
    logical, intent(in), optional :: do_wait
 
    integer :: n, f
    integer :: varid
    integer :: i
 
    real(rp), allocatable :: fpos3d(:,:)
    real(rp), allocatable :: lcfpos3d(:,:)
    real(rp), allocatable :: unity_fac(:)
    real(rp), allocatable :: tmp_svec3d(:,:)
    real(rp), allocatable :: tmp1_htensor3d(:,:,:)
    real(rp), allocatable :: tmp2_htensor3d(:,:,:)
    
    class(elementbase3d), pointer :: elem
    type(localmesh3d), pointer :: lcmesh
    type(localmeshcommdata), pointer :: commdata
 
    integer, pointer :: vmapb(:)
    !-----------------------------------------------------------------------------
 
    ! call PROF_rapstart( 'comm_exchange_1', 2)
 
    if ( present(do_wait) ) then
      if ( .not. do_wait ) then
        log_info("MeshFieldCommCubedSphereDom3D_exchange",*) "do_wait=False is not currently supported. Check!"
        call prc_abort
      end if
    end if
 
    do n=1, this%mesh%LOCAL_MESH_NUM
      lcmesh => this%mesh3d%lcmesh_list(n)
      elem => lcmesh%refElem3D
 
      allocate( fpos3d(this%Nnode_LCMeshAllFace(n),2) )
 
      if ( this%use_vmap_wide_flag ) then
        vmapb => this%VMapB2
      else
        vmapb => lcmesh%VMapB
      end if
 
      !$acc data create(fpos3D) copyin(VmapB)
 
      call meshfieldcommbase_extract_bounddata2( lcmesh%pos_en(:,:,1), vmapb, size(vmapb), elem%Np*lcmesh%Ne, fpos3d(:,1) )
      call meshfieldcommbase_extract_bounddata2( lcmesh%pos_en(:,:,2), vmapb, size(vmapb), elem%Np*lcmesh%Ne, fpos3d(:,2) )
 
      do f=1, this%nfaces_comm
        commdata => this%commdata_list(f,n)
        call push_localsendbuf( commdata%send_buf,                 &  ! (inout)
          this%send_buf(:,:,n), commdata%s_faceID, this%is_f(f,n), &  ! (in)
          commdata%Nnode_LCMeshFace, this%bufsize_per_field,       &  ! (in)
          this%field_num_tot, lcmesh, this%HaloSize_h1D            )  ! (in)
 
        if ( commdata%s_panelID /= lcmesh%panelID &
            .and. ( this%hvfield_num > 0 .or. this%htensorfield_num > 0 ) ) then
          
          allocate( lcfpos3d(commdata%Nnode_LCMeshFace,2), unity_fac(commdata%Nnode_LCMeshFace) )
          unity_fac(:) = 1.0_rp
          !$acc enter data create(lcfpos3D) copyin(unity_fac) async(1)
 
          call push_localsendbuf( lcfpos3d,                            &
            fpos3d, commdata%s_faceID, this%is_f(f,n),                 &
            commdata%Nnode_LCMeshFace, this%Nnode_LCMeshAllFace(n), 2, &
            lcmesh, this%HaloSize_h1D )
        end if
          
        if ( commdata%s_panelID /= lcmesh%panelID ) then
          if ( this%hvfield_num > 0 ) then
            allocate( tmp_svec3d(commdata%Nnode_LCMeshFace,2) )
            !$acc enter data create(tmp_svec3D) async(1)
 
            do varid=this%sfield_num+1, this%sfield_num+2*this%hvfield_num-1, 2
              !$acc parallel loop present(tmp_svec3D, commdata%send_buf) async(1)
              do i=1, size(tmp_svec3d,1)
                tmp_svec3d(i,1) = commdata%send_buf(i,varid  )
                tmp_svec3d(i,2) = commdata%send_buf(i,varid+1)
              end do
  
              call cubedspherecoordcnv_cs2lonlatvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac(:),  & ! (in)
                commdata%Nnode_LCMeshFace,                                   & ! (in)
                tmp_svec3d(:,1), tmp_svec3d(:,2),                            & ! (in)
                commdata%send_buf(:,varid), commdata%send_buf(:,varid+1),    & ! (out) 
                gpu_async_id=1 ) ! (in)
 
            end do
            !$acc exit data delete(tmp_svec3D) async(1)
            deallocate( tmp_svec3d )
          end if
 
          if ( this%htensorfield_num > 0 ) then
            allocate( tmp1_htensor3d(commdata%Nnode_LCMeshFace,2,2) )
            allocate( tmp2_htensor3d(commdata%Nnode_LCMeshFace,2,2) )
            !$acc enter data create(tmp1_htensor3D, tmp2_htensor3D) async(1)
            
            do varid=this%sfield_num+2*this%hvfield_num+1, this%field_num_tot-3, 4
              !$acc parallel loop present(tmp1_htensor3D, commdata%send_buf) async(1)
              do i=1, size(tmp1_htensor3d,1)
                tmp1_htensor3d(i,1,1) = commdata%send_buf(i,varid  )
                tmp1_htensor3d(i,2,1) = commdata%send_buf(i,varid+1)
                tmp1_htensor3d(i,1,2) = commdata%send_buf(i,varid+2)
                tmp1_htensor3d(i,2,2) = commdata%send_buf(i,varid+3)
              end do
  
              call cubedspherecoordcnv_cs2lonlatvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac(:),        & ! (in)
                commdata%Nnode_LCMeshFace,                                         & ! (in)
                tmp1_htensor3d(:,1,1), tmp1_htensor3d(:,2,1),                      & ! (in)
                tmp2_htensor3d(:,1,1), tmp2_htensor3d(:,2,1),                      & ! (out)
                gpu_async_id=1 ) ! (in)
              call cubedspherecoordcnv_cs2lonlatvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac(:),        & ! (in)
                commdata%Nnode_LCMeshFace,                                         & ! (in)
                tmp1_htensor3d(:,1,2), tmp1_htensor3d(:,2,2),                      & ! (in)
                tmp2_htensor3d(:,1,2), tmp2_htensor3d(:,2,2),                      & ! (out)
                gpu_async_id=1 ) ! (in)
              call cubedspherecoordcnv_cs2lonlatvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac(:),        & ! (in)
                commdata%Nnode_LCMeshFace,                                         & ! (in)
                tmp2_htensor3d(:,1,1), tmp2_htensor3d(:,1,2),                      & ! (in)
                commdata%send_buf(:,varid), commdata%send_buf(:,varid+2),          & ! (out)
                gpu_async_id=1 ) ! (in)
              call cubedspherecoordcnv_cs2lonlatvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac(:),        & ! (in)
                commdata%Nnode_LCMeshFace,                                         & ! (in)
                tmp2_htensor3d(:,2,1), tmp2_htensor3d(:,2,2),                      & ! (in)
                commdata%send_buf(:,varid+1), commdata%send_buf(:,varid+3),        & ! (out)
                gpu_async_id=1 ) ! (in)
            end do
            !$acc exit data delete(tmp1_htensor3D, tmp2_htensor3D) async(1)
            deallocate( tmp1_htensor3d, tmp2_htensor3d )
          end if
        
        end if
        
        if ( allocated(lcfpos3d) ) then
          !$acc exit data delete(lcfpos3D, unity_fac) async(1)
          deallocate( lcfpos3d, unity_fac )
        end if
 
      end do
      !$acc end data
      deallocate( fpos3d )
    end do
    !$acc wait(1)
    ! call PROF_rapend( 'comm_exchange_1', 2)
    
    !-----------------------
 
    ! call PROF_rapstart( 'comm_exchange_2', 2)
    call meshfieldcommbase_exchange_core(this, this%commdata_list, do_wait )
    ! call PROF_rapend( 'comm_exchange_2', 2)
 
    ! call PROF_rapstart( 'comm_exchange_3', 2)
    if ( .not. this%call_wait_flag_sub_get ) &
      call post_exchange_core( this )
 
    ! call PROF_rapend( 'comm_exchange_3', 2)
 
    return
  end subroutine meshfieldcommcubedspheredom3d_exchange
 
!----------------------------
 
!OCL SERIAL
  subroutine post_exchange_core( this )
    use scale_meshfieldcomm_base, only: &
      localmeshcommdata
    use scale_cubedsphere_coord_cnv, only: &
      cubedspherecoordcnv_lonlat2csvec
    implicit none
 
    class(meshfieldcommcubedspheredom3d), intent(inout), target :: this
 
    integer :: n, f
    integer :: varid
 
    real(rp), allocatable :: fpos3d(:,:)
    real(rp), allocatable :: lcfpos3d(:,:)
    real(rp), allocatable :: unity_fac(:)
    real(rp), allocatable :: tmp1_htensor3d(:,:,:)
    
    class(elementbase3d), pointer :: elem
    type(localmesh3d), pointer :: lcmesh
    type(localmeshcommdata), pointer :: commdata
 
    integer :: irs, ire
    integer, pointer :: vmapb(:)
    !-----------------------------------------------------------------------------
 
    do n=1, this%mesh%LOCAL_MESH_NUM
      lcmesh => this%mesh3d%lcmesh_list(n)
      elem => lcmesh%refElem3D
 
      allocate( fpos3d(this%Nnode_LCMeshAllFace(n),2) )
      if ( this%use_vmap_wide_flag ) then
        vmapb => this%VMapB2
      else
        vmapb => lcmesh%VMapB
      end if
      !$acc data create(fpos3D) copyin(VMapB)
 
      call meshfieldcommbase_extract_bounddata2( lcmesh%pos_en(:,:,1), vmapb, size(vmapb), elem%Np*lcmesh%Ne, fpos3d(:,1) )
      call meshfieldcommbase_extract_bounddata2( lcmesh%pos_en(:,:,2), vmapb, size(vmapb), elem%Np*lcmesh%Ne, fpos3d(:,2) )
 
      irs = 1
      do f=1, this%nfaces_comm
        commdata => this%commdata_list(f,n)
        ire = irs + commdata%Nnode_LCMeshFace - 1
 
        if ( commdata%s_panelID /= lcmesh%panelID &
            .and. ( this%hvfield_num > 0 .or. this%htensorfield_num > 0 ) ) then
          
          allocate( lcfpos3d(commdata%Nnode_LCMeshFace,2), unity_fac(commdata%Nnode_LCMeshFace) )
          unity_fac(:) = 1.0_rp
          !$acc enter data create(lcfpos3D) copyin(unity_fac) async(1)
 
          call push_localsendbuf( lcfpos3d,                            &
            fpos3d, f, this%is_f(f,n),                                 &
            commdata%Nnode_LCMeshFace, this%Nnode_LCMeshAllFace(n), 2, &
            lcmesh, this%HaloSize_h1D )
        end if
 
        if ( commdata%s_panelID /= lcmesh%panelID ) then
          if ( this%hvfield_num > 0 ) then            
            do varid=this%sfield_num+1, this%sfield_num+2*this%hvfield_num-1, 2
              call cubedspherecoordcnv_lonlat2csvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac,           & ! (in)
                commdata%Nnode_LCMeshFace,                                         & ! (in)
                commdata%recv_buf(:,varid), commdata%recv_buf(:,varid+1),          & ! (in)
                this%recv_buf(irs:ire,varid,n), this%recv_buf(irs:ire,varid+1,n),  & ! (out) 
                gpu_async_id=1 ) ! (in)
            end do
          end if
 
          if ( this%htensorfield_num > 0 ) then
            allocate( tmp1_htensor3d(commdata%Nnode_LCMeshFace,2,2) )
            !$acc enter data create(tmp1_htensor3D) async(1)
            
            do varid=this%sfield_num+2*this%hvfield_num+1, this%field_num_tot-3, 4
              call cubedspherecoordcnv_lonlat2csvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac,            & ! (in)
                commdata%Nnode_LCMeshFace,                                          & ! (in)
                commdata%recv_buf(:,varid), commdata%recv_buf(:,varid+1),           & ! (in)
                tmp1_htensor3d(:,1,1), tmp1_htensor3d(:,2,1),                       & ! (out)
                gpu_async_id=1 ) ! (in)
              call cubedspherecoordcnv_lonlat2csvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac,            & ! (in)
                commdata%Nnode_LCMeshFace,                                          & ! (in)
                commdata%recv_buf(:,varid+2), commdata%recv_buf(:,varid+3),         & ! (in)
                tmp1_htensor3d(:,1,2), tmp1_htensor3d(:,2,2),                       & ! (out)
                gpu_async_id=1 ) ! (in)
              call cubedspherecoordcnv_lonlat2csvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac,            & ! (in)
                commdata%Nnode_LCMeshFace,                                          & ! (in)
                tmp1_htensor3d(:,1,1), tmp1_htensor3d(:,1,2),                       & ! (in)
                this%recv_buf(irs:ire,varid,n), this%recv_buf(irs:ire,varid+2,n),   & ! (out)
                gpu_async_id=1 ) ! (in)
              call cubedspherecoordcnv_lonlat2csvec( &
                lcmesh%panelID, lcfpos3d(:,1), lcfpos3d(:,2), unity_fac,            & ! (in)
                commdata%Nnode_LCMeshFace,                                          & ! (in)
                tmp1_htensor3d(:,2,1), tmp1_htensor3d(:,2,2),                       & ! (in)
                this%recv_buf(irs:ire,varid+1,n), this%recv_buf(irs:ire,varid+3,n), & ! (out)
                gpu_async_id=1 ) ! (in)
            end do
            !$acc exit data delete(tmp1_htensor3D) async(1)
            deallocate( tmp1_htensor3d )
          end if          
        end if
 
        irs = ire + 1
        if ( allocated(lcfpos3d) ) then
          !$acc exit data delete(lcfpos3D, unity_fac) async(1)
          deallocate( lcfpos3d, unity_fac )
        end if
      end do
 
      !$acc end data
      deallocate( fpos3d )
    end do 
    !$acc wait(1)
 
    return
  end subroutine post_exchange_core
 
  !> Push temporary buffer of the local send data to the communication buffer
!OCL SERIAL
  subroutine push_localsendbuf( lc_send_buf, &
    send_buf, s_faceID, is, Nnode_LCMeshFace, bufsize_per_field, var_num, &
    lcmesh, haloSize_h1D )
    implicit none
 
    integer, intent(in) :: var_num
    integer, intent(in) ::  nnode_lcmeshface
    integer, intent(in) :: bufsize_per_field
    real(rp), intent(inout) :: lc_send_buf(nnode_lcmeshface,var_num)
    real(rp), intent(in) :: send_buf(bufsize_per_field,var_num)  
    integer, intent(in) :: s_faceid, is
    type(localmesh3d), pointer :: lcmesh
    integer, intent(in) :: halosize_h1d
 
    integer :: vid, i
    integer :: nnode_v, nnode_h1d, nfp_h
    !-----------------------------------------------------------------------------
 
    if ( s_faceid > 0 ) then
      !$omp parallel do
      !$acc parallel loop collapse(2) present(lc_send_buf, send_buf) async(1)
      do vid=1, var_num
      do i=1, nnode_lcmeshface
        lc_send_buf(i,vid) = send_buf(is+i-1,vid)
      end do
      end do
    else if ( -5 < s_faceid .and. s_faceid < 0) then
 
      nfp_h = lcmesh%refElem3D%Nfp_h
      nnode_h1d = lcmesh%refElem3D%Nnode_h1D
      nnode_v = lcmesh%refElem3D%Nnode_v
 
      call revert_hori( lc_send_buf, &
        send_buf, is, nnode_lcmeshface/(nfp_h * halosize_h1d * lcmesh%NeZ), &
        lcmesh%NeZ, halosize_h1d )      
    end if 
    
    return
  contains
    !> Revert the horizontal ordering of the local send buffer for negative s_faceID
    subroutine revert_hori(revert, ori, is_, Ne_h1D, NeZ, haloSize_h1D_ )
      integer, intent(in) :: is_
      integer, intent(in) :: ne_h1d, nez
      integer, intent(in) :: halosize_h1d_
      real(rp), intent(out) :: revert(nnode_h1d,nnode_v,halosize_h1d_, ne_h1d,nez, var_num)
      real(rp), intent(in)  :: ori(bufsize_per_field, var_num)
      
      integer :: p1, p3, ph, i, k, n
      integer :: i_, p1_
      integer :: id_src
      !-----------------------------------------------------------------------------
      
      !$acc parallel loop collapse(4) present(revert, ori) async(1)
      do n=1, var_num
      do k=1, nez
      do i=1, ne_h1d
      do ph=1, halosize_h1d_ 
        i_ = ne_h1d - i + 1
        do p3=1, nnode_v
        do p1=1, nnode_h1d
          p1_ = nnode_h1d - p1 + 1
          ! source ordering:　(p1_,p3,ph,i_,k)
          id_src = is_ - 1 &
            + p1_ &
            + (p3-1) * nnode_h1d             &
            + (ph-1) * nfp_h                 &
            + (i_-1) * nfp_h * halosize_h1d_ &
            + (k -1) * nfp_h * halosize_h1d_ * ne_h1d
          
          revert(p1,p3,ph,i,k,n) = ori(id_src,n)
        end do
        end do
      end do
      end do
      end do
      end do
      return
    end subroutine revert_hori    
  end subroutine push_localsendbuf
 
  subroutine set_boundary_data3d_u1u2( buf_U, buf_V, &
    bufsize_per_field, elem, mesh, G_ij, VMapB,      &
    u1, u2)
  
    implicit none
    integer, intent(in) :: bufsize_per_field
    type(elementbase3d), intent(in) :: elem
    type(localmesh3d), intent(in) :: mesh
    real(dp), intent(in) :: buf_u(bufsize_per_field)
    real(dp), intent(in) :: buf_v(bufsize_per_field)
    real(dp), intent(in) :: g_ij(elem%np * mesh%ne,2,2)
    integer, intent(in) :: vmapb(bufsize_per_field)
    real(dp), intent(inout) :: u1(elem%np * mesh%nea)
    real(dp), intent(inout) :: u2(elem%np * mesh%nea)
 
    integer :: ii, iis
    integer :: vmap_b
    !------------------------------------------------------------
 
#ifdef _OPENACC
    iis = elem%Np*mesh%NeE
    !$acc parallel loop present(buf_U, buf_V, G_ij, u1, u2, VMapB)
    do ii=1, bufsize_per_field
      vmap_b = vmapb(ii)
      u1(iis+ii) = g_ij(vmap_b,1,1) * buf_u(ii) + g_ij(vmap_b,1,2) * buf_v(ii)
      u2(iis+ii) = g_ij(vmap_b,2,1) * buf_u(ii) + g_ij(vmap_b,2,2) * buf_v(ii)
    end do
#else
    u1(elem%Np*mesh%NeE+1:elem%Np*mesh%NeE+bufsize_per_field) &
      = g_ij(vmapb,1,1) * buf_u(:) + g_ij(vmapb,1,2) * buf_v(:)
    u2(elem%Np*mesh%NeE+1:elem%Np*mesh%NeE+bufsize_per_field) &
      = g_ij(vmapb,2,1) * buf_u(:) + g_ij(vmapb,2,2) * buf_v(:)
#endif
    return
  end subroutine set_boundary_data3d_u1u2
 
end module scale_meshfieldcomm_cubedspheredom3d