scale_mesh_base1d.F90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
#include "scaleFElib.h"
module scale_mesh_base1d
 
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
 
  use scale_localmesh_1d, only: &
    localmesh1d, localmesh1d_init, localmesh1d_final
  
  use scale_mesh_base, only: meshbase, &
    meshbase_init, meshbase_final, &
    meshbase_setgeometricinfo
  
  use scale_element_base, only: elementbase1d
  use scale_element_line, only: lineelement
 
  !-----------------------------------------------------------------------------
  implicit none
  private
 
  !-----------------------------------------------------------------------------
  !
  !++ Public type & procedure
  ! 
  type, abstract, public, extends(meshbase) :: meshbase1d
    type(localmesh1d), allocatable :: lcmesh_list(:)
    class(elementbase1d), pointer :: refelem1d
    
    integer, public :: neg
    integer, public :: nprc
    real(rp), public :: xmin_gl, xmax_gl
    real(rp), public, allocatable :: fx(:)
  contains
    procedure(meshbase1d_generate), deferred :: generate 
    procedure :: getlocalmesh => meshbase1d_get_localmesh
  end type meshbase1d
  type, abstract, public, extends(meshbase) :: meshbase1d …
 
  interface 
    subroutine meshbase1d_generate(this)
      import meshbase1d
      class(meshbase1d), intent(inout), target :: this
    subroutine meshbase1d_generate(this) …
    end subroutine meshbase1d_generate
  end interface
 
  public :: meshbase1d_init, meshbase1d_final
  public :: meshbase1d_setgeometricinfo, meshbase1d_assigndomid, meshbase1d_setuplocaldom
 
  !-----------------------------------------------------------------------------
  !
  !++ Public parameters & variables
  !
  integer, public :: meshbase1d_dimtype_num   = 2
  integer, public :: meshbase1d_dimtypeid_x   = 1
  integer, public :: meshbase1d_dimtypeid_xt  = 2
  
  !-----------------------------------------------------------------------------
  !
  !++ Private procedure
  !
  
  !-----------------------------------------------------------------------------
  !
  !++ Private parameters & variables
  !
 
contains
!OCL SERIAL
  subroutine meshbase1d_init( this,         &
    NeG,                                    &
    dom_xmin, dom_xmax,                     &
    refElem, NLocalMeshPerPrc,              &
    nprocs, myrank, FX                      )
    
    implicit none
 
    class(meshbase1d), intent(inout) :: this
    integer, intent(in) :: neg
    real(rp), intent(in) :: dom_xmin
    real(rp), intent(in) :: dom_xmax   
    class(elementbase1d), intent(in), target :: refelem
    integer, intent(in) :: nlocalmeshperprc
    integer, intent(in), optional :: nprocs
    integer, intent(in), optional :: myrank
    real(rp), intent(in), optional :: fx(neg+1)
 
    integer :: n
    integer :: k
    real(rp) :: dx
    !-----------------------------------------------------------------------------
    
    this%NeG = neg
 
    this%xmin_gl       = dom_xmin
    this%xmax_gl       = dom_xmax
 
    !- Fx
    allocate( this%FX(neg+1) )
    if ( present(fx) ) then
      this%FX(:) = fx(:)
    else
      this%FX(1    ) = dom_xmin
      this%FX(neg+1) = dom_xmax
      dx = (dom_xmax - dom_xmin) / dble(neg)
      do k=2, neg
        this%FX(k) = this%FX(k-1) + dx
      end do
    end if
 
    this%refElem1D => refelem
    call meshbase_init( this,          &
      meshbase1d_dimtype_num, refelem, &
      nlocalmeshperprc, 2,             &
      nprocs                           )
 
    this%Nprc = this%PRC_NUM
         
    allocate( this%lcmesh_list(this%LOCAL_MESH_NUM) )
    do n=1, this%LOCAL_MESH_NUM
      call localmesh1d_init( this%lcmesh_list(n), n, refelem, myrank )
    end do
 
    call this%SetDimInfo( meshbase1d_dimtypeid_x, "x", "m", "X-coordinate" )
    call this%SetDimInfo( meshbase1d_dimtypeid_xt, "xt", "m", "X-coordinate" )
 
    return
  subroutine meshbase1d_init( this,         & …
  end subroutine meshbase1d_init
 
!OCL SERIAL
  subroutine meshbase1d_final( this )
    implicit none    
    class(meshbase1d), intent(inout) :: this
 
    integer :: n
    !-----------------------------------------------------------------------------
  
    if ( allocated ( this%lcmesh_list ) ) then 
      do n=1, this%LOCAL_MESH_NUM
        call localmesh1d_final( this%lcmesh_list(n), this%isGenerated )
      end do
  
      deallocate( this%lcmesh_list )
    end if
    
    call meshbase_final(this)
 
    return
  subroutine meshbase1d_final( this ) …
  end subroutine meshbase1d_final
  
!OCL SERIAL
  subroutine meshbase1d_get_localmesh( this, id, ptr_lcmesh )
    use scale_localmesh_base, only: localmeshbase
    implicit none
 
    class(meshbase1d), target, intent(in) :: this
    integer, intent(in) :: id
    class(localmeshbase), pointer, intent(out) :: ptr_lcmesh
    !-------------------------------------------------------------
 
    ptr_lcmesh => this%lcmesh_list(id)
    return
  end subroutine meshbase1d_get_localmesh
 
!OCL SERIAL
  subroutine meshbase1d_setgeometricinfo( lcmesh )
    implicit none
    
    type(localmesh1d), intent(inout) :: lcmesh
 
    class(elementbase1d), pointer :: refelem
    integer :: ke
    integer :: f
    integer :: i
 
    integer :: node_ids(lcmesh%refelem%nv)
    real(rp) :: vx(lcmesh%refelem%nv)
    real(rp) :: xr(lcmesh%refelem%np)
    real(dp) :: escale(1,1,lcmesh%refelem%np)
    integer :: fmask(lcmesh%refelem%nfptot)
    integer :: fid(lcmesh%refelem1d%nfp,lcmesh%refelem1d%nfaces)
 
  !-----------------------------------------------------------------------------
 
    call meshbase_setgeometricinfo(lcmesh, 1)
    refelem => lcmesh%refElem1D
    
    fmask(:) = reshape(refelem%Fmask, shape(fmask))
    do f=1, refelem%Nfaces
    do i=1, refelem%Nfp
       fid(i,f) = i + (f-1)*refelem%Nfp
    end do
    end do
  
    do ke=1, lcmesh%Ne
      node_ids(:) = lcmesh%EToV(ke,:)
      vx(:) = lcmesh%pos_ev(node_ids(:),1)
      lcmesh%pos_en(:,ke,1) = vx(1) + 0.5_rp*(refelem%x1(:) + 1.0_rp)*(vx(2) - vx(1))
      
      xr(:) = 0.5_rp*(vx(2) - vx(1)) !matmul(refElem%Dr,mesh%x(:,n))
 
      lcmesh%J(:,ke) = xr
      lcmesh%Escale(:,ke,1,1) =   1.0_rp/lcmesh%J(:,ke)
 
      !* Face
 
      !
      !mesh%fx(:,n) = mesh%x(fmask(:),n)
      !mesh%fy(:,n) = mesh%y(fmask(:),n)
 
      ! Calculate normal vectors
      lcmesh%normal_fn(fid(:,1),ke,1) = - 1.0_rp
      lcmesh%normal_fn(fid(:,2),ke,1) = + 1.0_rp
      lcmesh%sJ(:,ke) = 1.0_rp
      lcmesh%Fscale(:,ke) = lcmesh%sJ(:,ke)/lcmesh%J(fmask(:),ke)
      lcmesh%Gsqrt(:,ke) = 1.0_rp
    end do
 
    return
  subroutine meshbase1d_setgeometricinfo( lcmesh ) …
  end subroutine meshbase1d_setgeometricinfo
 
!OCL SERIAL
  subroutine meshbase1d_assigndomid( this, &
    tileID_table, panelID_table,        &
    pi_table )
  
    use scale_meshutil_1d, only: &       
      meshutil1d_buildglobalmap
    implicit none
    
    class(meshbase1d), intent(inout) :: this    
    integer, intent(out) :: tileid_table(this%local_mesh_num, this%prc_num)
    integer, intent(out) :: panelid_table(this%local_mesh_num*this%prc_num)
    integer, intent(out) :: pi_table(this%local_mesh_num*this%prc_num)
    
    integer :: n
    integer :: p
    integer :: tileid
    !-----------------------------------------------------------------------------
    
    call meshutil1d_buildglobalmap( &
      panelid_table, pi_table,                                                      & ! (out)
      this%tileID_globalMap, this%tileFaceID_globalMap, this%tilePanelID_globalMap, & ! (out)
      this%LOCAL_MESH_NUM_global )                                                    ! (in)
 
    do p=1, this%PRC_NUM
    do n=1, this%LOCAL_MESH_NUM
      tileid = n + (p-1)*this%LOCAL_MESH_NUM
 
      tileid_table(n,p)                   = tileid
      this%tileID_global2localMap(tileid) = n
      this%PRCRank_globalMap(tileid)      = p - 1
    end do
    end do
 
    return
  subroutine meshbase1d_assigndomid( this, & …
  end subroutine meshbase1d_assigndomid
 
!OCL SERIAL
  subroutine meshbase1d_setuplocaldom( lcmesh,  &
    tileID, panelID,                            &
    i, Nprc,                                    &
    dom_xmin, dom_xmax,                         &
    Ne, FX )
 
    use scale_meshutil_1d, only: &
      meshutil1d_genconnectivity,   &
      meshutil1d_genlinedomain,     &
      meshutil1d_buildinteriormap,  &
      meshutil1d_genpatchboundarymap
 
    use scale_localmesh_base, only: bctype_interior
    implicit none
      
    type(localmesh1d), intent(inout) :: lcmesh
    integer, intent(in) :: tileid
    integer, intent(in) :: panelid
    integer, intent(in) :: i
    integer, intent(in) :: nprc
    real(rp) :: dom_xmin, dom_xmax
    integer, intent(in) ::ne
    real(rp), intent(in) :: fx(ne*nprc+1)
    
    class(elementbase1d), pointer :: elem
    real(rp) :: delx
    real(rp) :: fx_lc(ne+1)    
    !-----------------------------------------------------------------------------
 
    elem => lcmesh%refElem1D
    lcmesh%tileID = tileid
    lcmesh%panelID = panelid
    
    !--
 
    lcmesh%Ne  = ne
    lcmesh%Nv  = ne + 1
    lcmesh%NeS = 1
    lcmesh%NeE = lcmesh%Ne
    lcmesh%NeA = lcmesh%Ne + 2
 
    !delx = (dom_xmax - dom_xmin)/dble(Nprc)
    fx_lc(:) = fx((i-1)*ne+1:i*ne)
    lcmesh%xmin = fx_lc(1)
    lcmesh%xmax = fx_lc(ne+1)
 
    allocate(lcmesh%pos_ev(lcmesh%Nv,1))
    allocate( lcmesh%EToV(lcmesh%Ne,2) )
    allocate( lcmesh%EToE(lcmesh%Ne,elem%Nfaces) )
    allocate( lcmesh%EToF(lcmesh%Ne,elem%Nfaces) )
    allocate( lcmesh%BCType(lcmesh%refElem%Nfaces,lcmesh%Ne) )
    allocate( lcmesh%VMapM(elem%NfpTot, lcmesh%Ne) )
    allocate( lcmesh%VMapP(elem%NfpTot, lcmesh%Ne) )
    allocate( lcmesh%MapM(elem%NfpTot, lcmesh%Ne) )
    allocate( lcmesh%MapP(elem%NfpTot, lcmesh%Ne) )
 
    lcmesh%BCType(:,:) = bctype_interior
 
    !----
 
    call meshutil1d_genlinedomain( lcmesh%pos_ev, lcmesh%EToV,   & ! (out)
        lcmesh%Ne, lcmesh%xmin, lcmesh%xmax, fx=fx_lc   )            ! (in)
 
 
    !---
    call meshbase1d_setgeometricinfo( lcmesh )
    
    !---
 
    call meshutil1d_genconnectivity( lcmesh%EToE, lcmesh%EToF, & ! (out)
        & lcmesh%EToV, lcmesh%Ne, elem%Nfaces )                  ! (in)
 
    !---
    call meshutil1d_buildinteriormap( lcmesh%VmapM, lcmesh%VMapP, lcmesh%MapM, lcmesh%MapP,  &
      & lcmesh%pos_en, lcmesh%pos_ev, lcmesh%EToE, lcmesh%EtoF, lcmesh%EtoV,                   &
      & elem%Fmask, lcmesh%Ne, elem%Np, elem%Nfp, elem%Nfaces, lcmesh%Nv )
 
    call meshutil1d_genpatchboundarymap( lcmesh%VMapB, lcmesh%MapB, lcmesh%VMapP, &
      & lcmesh%pos_en, lcmesh%xmin, lcmesh%xmax,                                  &
      & elem%Fmask, lcmesh%Ne, elem%Np, elem%Nfp, elem%Nfaces, lcmesh%Nv)
    
    return
  subroutine meshbase1d_setuplocaldom( lcmesh,  & …
  end subroutine meshbase1d_setuplocaldom
 
end module scale_mesh_base1d