!> \~japanese
!! \brief 正二十面体格子における格子点およびそのコントロールボリュームの頂点の座標データを管理する構造型と
!!それらのデータを操作する手続きを提供するモジュール.
!!
!! \par IcGrid2D_FVM_Manager モジュールについて
!!
!! \~english
!! \brief This module provides a derived type, which manages the coordinate data for grid points and their control volume in icosahedral gird,
!!and some subroutines to manupulate their data.
!!
!! \par About IcGrid2D_FVM_Manager module
!!
!! \~
!! \since IGMBaseLib 1.0
!!
!! <br><br>
!! \b Copyright (C) GFD Dennou Club, 2011-2012. All rights reserved. <br>
!! \b License ?? <br>
!! \author Yuta Kawai
!!
module IcGrid2D_FVM_Manager

  ! モジュール引用 ; Use statements 
  !

  ! 種類型パラメータ
  ! Kind type parameter
  !
  use dc_types, only: DP ! 倍精度実数型 ; Double recision
    
  ! 数学定数
  ! Mathematic
  !
  use igmcore_math, only: PI     ! 円周率. PI.

  ! 線形代数
  ! Linear algebra
  !
  use igmcore_linear_algebra, only: vec_length, &
    &                       rotateZ

  ! 座標変換
  ! Coordinate conversion
  !
  use igmcore_coordinate_conversion, only: geo_to_orth_pos

  ! 正二十面格子における広域矩形領域のサブドメインである局所矩形領域(各 MPI プロセスと対応)の情報管理
  ! Manager of information about a local rectangle region (corresponding to each MPI process)
  !
  use LocalProcRect_Manager, only: &
    & LocalProcRectInfo, &
    & LocalProcRectInfo_Init, LocalProcRectInfo_Final

  ! 宣言文 ; Declaration statemenets
  !

  implicit none
  private
  
  ! 公開変数
  ! Public variables 
  !

  !> \~japanese 全矩形領域の数.
  !! \~english Number of all rectangle regions.
  integer, parameter, public :: RC_REGIONS_NUM = 10     

  !> \~japanese 正二十面体の頂点数.
  !! \~english The number of vertices of the icosahedron.
  integer, parameter, public :: ICOSAHEDRON_VERTEX_NUM = 12     

  !> \~japanese 格子点に付随するコントロールボリュームの頂点数(特異点は例外的に 5 ).
  !! \~english
  integer, parameter, public :: CELL_POINTS_NUM = 6     

  !> \~japanese 指定された格子点が両極と一致していないことを示す定数.
  !! \~english The constant indicating that the specified grid point coincides with both poles.
  integer, parameter, public :: NOT_POLE_FLAG  =  0

  !> \~japanese 指定された格子点が北極と一致していることを示す定数.
  !! \~english The constant indicating that the specified grid point concides with the north pole.
  integer, parameter, public :: NORTH_POLE_FLAG = 1

  !> \~japanese 指定された格子点が南極と一致していることを示す定数.
  !! \~english The constant indicating that the specified grid point concides with the south pole.
  integer, parameter, public :: SOUTH_POLE_FLAG = 2

  ! 格子タイプを表す定数
  !

  !>
  !!
  integer, parameter, public :: GTYPE_STDGRID = 1

  !>
  !!
  integer, parameter, public :: GTYPE_STDGCGRID = 2

  !>
  !!
  integer, parameter, public :: GTYPE_SPRGCGRID = 3

  !> \~japanse 
  !! \brief 正二十面格子のある格子点において定義される 1 個のコントロールボリュームの情報を管理する構造型.
  !! \~english
  !! \brief 正二十面格子のある格子点において定義される 1 個のコントロールボリュームの情報を管理する構造型.
  !!
  type, public :: IcGrid2D_FVM_CV
    
    ! \~japanese コントロールボリュームの各頂点の位置ベクトルのデータを保持する配列. 
    ! \~english  Array to store the position vectors of 6 (or 5) vertecies of a control volume.
    real(DP) :: cv_pts(3, CELL_POINTS_NUM)

    ! \~japanese コントロールボリュームの面積. 
    ! \~english The area of a control volume.
    real(DP) :: cv_area

  end type IcGrid2D_FVM_CV

  !
  !> \~japanese
  !! \brief 正二十面体格子における格子点およびそのコントロールボリュームの頂点の座標データを管理する構造型.
  !! \~english
  !! \brief Derived type which manages the coordinate datas for grid points and their control volume in icosahedral gird.
  !!
  type, public :: IcGrid2D_FVM

    !> \~japanese
    !! 正二十面格子の格子点(A-grid)の座標データを保持する配列.
    !! 配列の各次元は, それぞれ,  
    !! 3 次元位置ベクトルの各成分, 水平方向の添え字, 水平方向の添え字, 広域な矩形領域の ID を意味する.
    !!
    !! \~english
    !! 正二十面格子の格子点(A-grid)の座標データを保持する配列.
    !! 配列の各次元は, それぞれ,  
    !! 3 次元位置ベクトルの各成分, 水平方向の添え字, 水平方向の添え字, 矩形領域の ID を意味する.
    !!
    real(DP), pointer :: rcs_AGrid(:,:,:,:) => null()

    !> \~japanese
    !! 正二十面格子の各格子点のコントーロルボリュームに関する情報を管理する構造体 IcGrid2D_FVM_CV の変数の配列. 
    !! 配列の各次元は, それぞれ,  
    !! 水平方向の添え字, 水平方向の添え字, 広域な矩形領域の ID を意味する.
    !!
    !! \~english
    !! 正二十面格子の各格子点のコントーロルボリュームに関する情報を管理する構造体 IcGrid2D_FVM_CV の変数の配列. 
    !! 配列の各次元は, それぞれ,  
    !! 水平方向の添え字, 水平方向の添え字, 矩形領域の ID を意味する.
    !!
    type(IcGrid2D_FVM_CV), pointer :: rcs_CVInfo(:,:,:) => null()
    

    ! 非公開メンバ変数
    ! Private member variables
    !

    !> \~japanese 水平分割レベル.
    !! \~english Horizontal division level.
    integer :: h_glevel
    ! f2003: integer, private :: h_glevel

    !> \~japanese 正二十面体を内包する球の半径.
    !! \~english The radius of a sphere which an icosahedron is embedded into.
    real(DP) :: icRadius
    ! f2003: real(DP), private :: icRadius


    !> \~japanese ノリシロ部分を除いた水平一方向の配列サイズ.
    !! \~english ノリシロ部分を除いた水平一方向の配列サイズ.
    integer :: grid_EffSize
    ! f2003: integer, private :: grid_EffSize

    !> \~japanese ノリシロの大きさ(デフォルトは 1 ).
    !! \~english ノリシロの大きさ(デフォルトは 1 ).
    integer :: grid_margine = 1
    ! f2003: integer, private :: grid_margine = 1

    !> \~japanese 格子点座標を保持する配列インデックス(水平一方向)の下限.
    !! \~english 格子点座標を保持する配列インデックス(水平一方向)の下限.
    integer :: idMin
    ! f2003: integer, private :: idMin

    !> \~japanese 格子点座標を保持する配列インデックス(水平一方向)の下限.
    !! \~english 格子点座標を保持する配列インデックス(水平一方向)の下限.
    integer :: idMax
    ! f2003: integer, private :: idMax

    !> \~japanese 格子点座標を保持する配列ののりしろ部分除いたインデックス(水平一方向)の下限(デフォルトは 1 ).
    !> \~english 格子点座標を保持する配列ののりしろ部分除いたインデックス(水平一方向)の下限(デフォルトは 1 ).
    integer :: idEfMin = 1
    ! f2003: integer, private :: idEfMin = 1

    !> \~japanese 格子点座標を保持する配列ののりしろ部分除いたインデックス(水平一方向)の上限.
    !! \~english 格子点座標を保持する配列ののりしろ部分除いたインデックス(水平一方向)の上限.
    integer :: idEfMax
    ! f2003: integer, private :: idEfMax
    
    !> \~japanese
     !! 正二十面格子データを管理する当 IcGrid2D_FVM オブジェクトと関係づけられる局所矩形領域の情報を管理するための
     !! 構造型 LocalProcRectInfo の変数.
    !! \~english 
    !! A variable of derived type LocalProcRectInfo which manages the information about local rectangle region
    !! associated with this IcGrid2D_FVM object managing icosahedral grid data.
    type(LocalProcRectInfo) :: lcproc_RCInfo

  end type IcGrid2D_FVM
  

  ! 公開関数
  ! Public procedures
  !
  public :: IcGrid2D_FVM_Init, IcGrid2D_FVM_Final

  public :: malloc_GP_array, malloc_CV_array,                       &    
    &       generate_CV5_GPindex, generate_CV6_GPindex,             & 
    &       paste_margin_width, check_pole, calc_icosahedron_vertex

  public :: set_rcregion_bounds 
  public :: get_EffSize_Min, get_EffSize_Max, get_IdMin, get_IdMax, &
    &       get_glevel, get_IcRadius, get_IcGridMeanInterval,       &
    &       get_LocalProcRectInfo


contains

!
!> \~japanese
!! \brief 構造型 IcGrid2D_FVM の変数の初期化を行う.
!!
!! 正二十面格子の属性(ex 水平分割レベル)を設定する.
!! この初期化では, 正二十面格子やそのコントロールボリュームの情報を保持する配列のメモリは確保されない.
!! これらは, それぞれ本モジュールが提供する手続き malloc_GP_array, malloc_CV_array を用いて行うことができる.
!!
!! (optional 属性の)　4 番目の引数 mpi_flag が省略された場合は, MPI 並列は無効になる(default).
!!
!! @param[in,out] self      構造型 IcGrid2D_FVM の変数.
!! @param[in]     glevel    水平分割レベル.
!! @param[in]     ic_radius 正二十面体を内包する球の半径.
!! @param[in]     mpi_flag  MPI による並列化を有効にするかどうかを指定するフラグ.
!!
!! \~english
!! \brief Initializes the variable of derived type FVM_IcGrid2D.
!!
!! Some icosahedral grid attributes(for example, horizontal division level) are set.
!! In this initialization, memory for array to store data of icosahedral grid and the control volumes is NOT allocated.
!! If you want to allocate memory for the array, please use the malloc_GP_array and malloc_CV_array which are procedures in this module.
!!
!! If forth (optional) argument is ommitted, the parallelization with MPI will be disabled(default).
!!
!! @param[in,out] self      A variable of derived type IcGrid2D_FVM.
!! @param[in]     glevel    Horizontal division level.
!! @param[in]     ic_radius Radius of sphere embedded in a icosahedron.
!! @param[in]     mpi_flag  A flag which specifies whether to enable or disable the parallelization with MPI.
!!
subroutine IcGrid2D_FVM_Init( &
  & self,              &  ! (inout)
  & glevel, ic_radius, &  ! (in)
  & mpi_flag           &
  )

  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(inout) :: self
  integer, intent(in) :: glevel
  real(DP), intent(in) :: ic_radius
  logical, optional :: mpi_flag
  
  ! 作業変数
  ! Work variables
  !
  logical :: mpi_enable_flag

  ! 実行文 ; Executable statements
  !
  self%h_glevel = glevel
  self%icRadius = ic_radius

  ! 指定された水平格子の分割レベルに対応する要素数を計算する. 
  self%grid_EffSize = 2**self%h_glevel + 1 
  
  ! IcGrid2D_FVM クラスの配列変数のインデックスののりしろを含めた下限・上限,
  ! 物理的に意味のある下限・上限を求める. 
  self%idEfMax = self%idEfMin + self%grid_EffSize - 1
  self%idMin = self%idEfMin - self%grid_margine
  self%idMax = self%idEfMax + self%grid_margine
  
  ! MPI による並列化が有効かフラグを判定する. 
  if ( present(mpi_flag) ) then
    mpi_enable_flag = mpi_flag
  else
    mpi_enable_flag = .false.
  end if

  ! 局所矩形領域のハロー領域に関する情報を保持する構造型 LocalProcRectInfo の変数を初期化する. 
  call LocalProcRectInfo_Init(self%lcproc_RCInfo, &
    & self%idEfMin, self%idEfMax, self%idMin, self%idMax, &
    & mpi_enable_flag )

end subroutine IcGrid2D_FVM_Init

!
!> \~japanese
!! \~brief 正二十面格子の座標を保持する配列のメモリを確保する.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!!
!! \~english
!! \brief Allocates memory for array to stroe the vertices of icosahedral grid.
!!
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!!
subroutine malloc_GP_array( &
  & self                    &  ! (inout)
  & )

  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(inout) :: self

  ! 作業変数
  ! Work variables
  !
  
  ! 実行文 ; Executable statements
  !
  
  allocate(self%rcs_AGrid(3, self%idMin:self%idMax, self%idMin:self%idMax, RC_REGIONS_NUM))

end subroutine malloc_GP_array

!
!> \~japanese
!! \brief 各格子点に付随するコントロールボリュームの情報を保持する構造型 IcGrid2D_FVM_CV の配列のメモリを確保する.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!!
!! \~english
!! \brief Allocates memory for array of derived type IcGrid2D_FVM_CV to store some information about control volumes 
!!        surrounding each grid point.
!!
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!!
!!
subroutine malloc_CV_array( &
  & self                    &  ! (inout)
  & )
  
  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(inout) :: self
  
  ! 実行文 ; Executable statements
  !

  allocate( self%rcs_CVInfo(self%idEfMin:self%idEfMax, self%idEfMin:self%idEfMax, RC_REGIONS_NUM) )
  
end subroutine malloc_CV_array
 
!
!> \~japanese
!! \brief 構造型 IcGrid2D_FVM の変数を最終化する.
!!
!! 構造型 IcGrid2D_FVM において動的に確保された全ての配列のメモリを解放する.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!!
!! \~english
!! \brief Finalizes a variable of derived type IcGrid2D_FVM.
!!
!! This subroutine releases memory for all dynamically allocated arrays in derived type IcGrid2D_FVM.
!!
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!!
subroutine IcGrid2D_FVM_Final( &
  & self             &  ! (inout)
  & )

  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(inout) :: self

  ! 実行文 ; Executable statements
  !

  ! 本モジュールにおいて, 動的に確保した配列を解放する.

  if ( associated(self%rcs_AGrid) ) then
    deallocate(self%rcs_AGrid)
  end if

  if ( associated(self%rcs_CVInfo) ) then
     deallocate(self%rcs_CVInfo)
  end if

  call LocalProcRectInfo_Final(self%lcproc_RCInfo)

end subroutine IcGrid2D_FVM_Final

function get_LocalProcRectInfo(self) result(lcprocRCInfo)

  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(in), target :: self
  type(LocalProcRectInfo), pointer :: lcprocRCInfo

  !
  !
  !
  lcprocRCInfo => self%lcproc_RCInfo

end function get_LocalProcRectInfo
!
!> \~japanese
!! \brief 格子点座標を保持する配列において, 物理的意味がある配列の最小インデックスを取得する.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!! @return    gridEffMin      格子点座標を保持する配列ののりしろ部分除いたインデックス(水平一方向)の下限.
!!
!! \~english
!! \brief 格子点座標を保持する配列において, 物理的意味がある配列の最小インデックスを取得する.
!!
!! @param[in,out] self  A variable of derived type IcGrid2D_FVM.
!! @return    gridEffMin      格子点座標を保持する配列ののりしろ部分除いたインデックス(水平一方向)の下限. 
!!
function get_EffSize_Min( &
  & self                  &  ! (in)
  & ) result(gridEffMin)
 
  type(IcGrid2D_FVM), intent(in) :: self
  integer gridEffMin
  
  gridEffMin = self%idEfMin

end function get_EffSize_Min

!
!> \~japanese
!! \brief 格子点座標を保持する配列において, 物理的意味がある配列の最大インデックスを取得する.
!!
!! @param[in,out] self　 構造型 IcGrid2D_FVM の変数.
!! @return    gridEffMin 格子点座標を保持する配列ののりしろ部分除いたインデックス(水平一方向)の上限.
!!
!! \~english
!! \brief 格子点座標を保持する配列において, 物理的意味がある配列の最大インデックスを取得する.
!!
!! @param[in,out] self　 A variable of derived type IcGrid2D_FVM.
!! @return    gridEffMin 格子点座標を保持する配列ののりしろ部分除いたインデックス(水平一方向)の上限.
!!
function get_EffSize_Max( &
  & self                  &  ! (in)
  & ) result(gridEffMax)

  type(IcGrid2D_FVM), intent(in) :: self
  integer gridEffMax
  
  gridEffMax = self%idEfMax

end function get_EffSize_Max

!
!> \~japanese
!! \brief 格子点座標を保持する配列(構造体 IcGrid2D_FVM の成分 rcs_AGrid)において, 袖領域も含めた配列の最大インデックスを取得する.
!!
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!! @return rcs_AGrid の袖領域も含めた配列の最大インデックス.
!!
!! \~english
!! \brief 格子点座標を保持する配列(構造体 IcGrid2D_FVM の成分 rcs_AGrid)において, 袖領域も含めた配列の最大インデックスを取得する.
!!
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!! @return rcs_AGrid の袖領域も含めた配列の最大インデックス.
!!
function get_IdMax( &
  & self            &  ! (in)
  & ) result(id_Max)

  type(IcGrid2D_FVM), intent(in) :: self
  integer id_Max
  
  id_Max = self%idMax

end function get_IdMax

!
!> \~japanese
!! \brief 格子点座標を保持する配列(構造体 IcGrid2D_FVM の成分 rcs_AGrid)において, 袖領域も含めた配列の最小インデックスを取得する.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!! @return rcs_AGrid の袖領域を含めた配列の最小インデックス.
!!
!! \~english
!! \brief 格子点座標を保持する配列(構造体 IcGrid2D_FVM の成分 rcs_AGrid)において, のりしろも含めた配列の最小インデックスを取得する.
!!
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!! @return rcs_AGrid の袖領域を含めた配列の最小インデックス.
!!
function get_IdMin( &
  & self            &  ! (in)
  & ) result(id_Min)

  type(IcGrid2D_FVM), intent(in) :: self
  integer id_Min
  
  id_Min = self%idMin

end function get_IdMin

!
!> \~japanese
!! \brief 正二十面体を内包する球の半径を取得する.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!! @return 正二十面体を内包する球の半径.
!!
!! \~english
!! \brief Gets the radius of a sphere which an icosahedron is embedded into.
!! 
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!! @return The radius of a sphere which an icosahedron is embedded into.
!!
function get_IcRadius( &
  & self               &  ! (in)
  & ) result(ic_radius)
 
  type(IcGrid2D_FVM), intent(in) :: self
  real(DP) ic_radius
  
  ic_radius = self%icRadius
end function get_IcRadius

!
!> \~japanese
!! \brief 水平格子分割レベル(glevel)を取得する.
!!
!! @param[in] self 構造型 IcGrid2D_FVM の変数.
!! @return 水平格子分割レベル. 
!!
!! \~english
!! \brief Get the horizontal division level (glevel).
!!
!! @param[in] self A variable of derived type IcGrid2D_FVM.
!! @return Horizontal division level.
!!
function get_glevel( &
  & self             &  ! (in)
  & ) result(glevel)
 
  type(IcGrid2D_FVM), intent(in) :: self  ! IcGrid2D_FVM クラスのインスタンスの参照
  integer glevel
  
  glevel = self%h_glevel
end function get_glevel

!
!> \~japanese
!! \brief 矩形領域の端の 4 個の頂点座標を設定する.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!! @param[in]     rcID   頂点座標を設定する矩形領域の ID.
!! @param[in]     top    上側頂点に対応する格子点座標.
!! @param[in]     left   左側頂点に対応する格子点座標.
!! @param[in]     right  右側頂点に対応する格子点座標.
!! @param[in]     bottom 下側頂点に対応する格子点座標.
!!
!! \~english
!! \brief 矩形領域の端の 4 個の頂点座標を設定する.
!!
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!! @param[in]     rcID   頂点座標を設定する矩形領域の ID. 
!! @param[in]     top    上側頂点に対応する格子点座標. 
!! @param[in]     left   左側頂点に対応する格子点座標. 
!! @param[in]     right  右側頂点に対応する格子点座標. 
!! @param[in]     bottom 下側頂点に対応する格子点座標. 
!!
subroutine set_rcregion_bounds(    &
  & self,                          &  ! (inout)
  & rcID, top, left, bottom, right &  ! (in)
  & )

  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(inout) :: self
  integer, intent(in) :: rcID
  real(DP), intent(in) :: top(3)
  real(DP), intent(in) :: left(3)
  real(DP), intent(in) :: bottom(3)
  real(DP), intent(in) :: right(3)
  
  self%rcs_AGrid(:, self%idEfMin, self%idEfMin, rcID) = top
  self%rcs_AGrid(:, self%idEfMax, self%idEfMin, rcID) = left
  self%rcs_AGrid(:, self%idEfMax, self%idEfMax, rcID) = bottom
  self%rcs_AGrid(:, self%idEfMin, self%idEfMax, rcID) = right

end subroutine set_rcregion_bounds

!
!> \~japanese
!! \brief のりしろ部分をとなりの矩形領域の内側の格子座標データで埋める.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の参照.
!!
!! \~english
!! \brief のりしろ部分をとなりの矩形領域の内側の格子座標データで埋める.
!!
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!!
subroutine paste_margin_width( &
  & self                       &  ! (inout) 
  & )

  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(inout) :: self

  ! 作業変数
  ! Work variables
  !
  integer :: NHS_glRCIdList(7) = (/ 5, 1, 2, 3, 4, 5, 1 /)
  integer :: SHS_glRCIDList(7) = (/ 10, 6, 7, 8, 9, 10, 6 /)
  integer :: EMin, EMax, marginw
  integer rcID, k

  ! 実行文 ; Executable statements
  !

  EMin = get_EffSize_Min(self)
  EMax = get_EffSize_Max(self)
  marginw = EMin - get_IdMin(self)

  do rcID=1, RC_REGIONS_NUM/2
    do k=1, marginw
      ! 北半球側の矩形領域ののりしろに値をコピーする.
      self%rcs_AGrid(:, EMin:EMax, EMin-k, rcID) = self%rcs_AGrid(:, EMin+k, EMin:EMax, NHS_glRCIdList(rcID))
      self%rcs_AGrid(:, EMax+k, EMin:EMax, rcID) = self%rcs_AGrid(:, EMin+k, EMin:EMax, SHS_glRCIDList(rcID))
      self%rcs_AGrid(:, EMin:EMax, EMax+k, rcID) = self%rcs_AGrid(:, EMin:EMax, EMin+k, SHS_glRCIDList(rcID+1))
      self%rcs_AGrid(:, EMin-k, EMin:EMax, rcID) = self%rcs_AGrid(:, EMin:EMax, EMin+k, NHS_glRCIdList(rcID+2))

      ! 南半球側の矩形領域ののりしろに値をコピーする.
      self%rcs_AGrid(:, EMin:EMax, EMin-k, rcID+5) = self%rcs_AGrid(:, EMin:EMax, EMax-k, NHS_glRCIdList(rcID+1))
      self%rcs_AGrid(:, EMax+k, EMin:EMax, rcID+5) = self%rcs_AGrid(:, EMin:EMax, EMax-k, SHS_glRCIDList(rcID))
      self%rcs_AGrid(:, EMin:EMax, EMax+k, rcID+5) = self%rcs_AGrid(:, EMax-k, EMin:EMax, SHS_glRCIDList(rcID+2))
      self%rcs_AGrid(:, EMin-k, EMin:EMax, rcID+5) = self%rcs_AGrid(:, EMax-k, EMin:EMax, NHS_glRCIdList(rcID+2))
    end do
  end do

end subroutine paste_margin_width

!
!> \~japanese
!! \brief (特異格子点(rcID,GP_i,GP_j)の)コントールボリュームの各頂点(全 5 個)の計算に必要な格子点の id の組を取得する.
!!
!! ただし, 北半球および南半球の極にあたる格子点のコントロールボリュームに対しては, 適用できない.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!! @param[in]     GP_i 対象となる格子点(rcID, GP_i, GP_j) の GP_i.
!! @param[in]     GP_i 対象となる格子点(rcID, GP_i, GP_j) の GP_j.
!! @param[in]     rcID 対象となる格子点(rcID, GP_i, GP_j) の rcID.
!! @return CV_GPId コントールボリュームの各頂点の計算に必要な格子点の id の組.
!!
!! \~english
!! \brief (特異格子点(rcID,GP_i,GP_j)の)コントールボリュームの各頂点(全 5 個)の計算に必要な格子点の id の組を取得する.
!!
!! ただし, 北半球および南半球の極にあたる格子点のコントロールボリュームに対しては, 適用できない. 
!!
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!! @param[in]    i    i of target grid point(rcID, i, j).
!! @param[in]    j    j of target grid point(rcID, i, j).
!! @param[in]    rcID The rectangle region index rcID of target grid point(rcID, i, j).
!! @return CV_GPId The set of grid-point ID nessesary for calculating the vertices of the control volume.
!!
function generate_CV5_GPindex( &
  & self,                      &  ! (inout)
  & GP_i, GP_j, rcID           &  ! (in)
  & ) result(CV_GPId)
  
  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(inout) :: self
  integer, intent(in) :: GP_i
  integer, intent(in) :: GP_j
  integer, intent(in) :: rcID
  integer CV_GPId(5,3,2)  

  ! 作業変数
  ! Work variables
  !
  integer :: glEMin, glEMax
  integer :: i,j
  integer :: CV_GPIdShape(2) = (/ 3, 2 /)

  ! 実行文 ; Executable statements
  !
  
  glEMin = get_EffSize_Min(self)
  glEMax = get_EffSize_Max(self)

  i = GP_i; j = GP_j
  if( rcID <= 5 )then
    ! 北半球に対して
    if( i==glEMax .and. j == glEMax ) then
      ! 下端特異点
      CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape ) 
      CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape ) 
      CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape ) 
      CV_GPId(4,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape ) 
      CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape ) 
    else if( j == glEMin ) then
      ! 左端特異点
      CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape )
      CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape )
      CV_GPId(3,:,:) = reshape( (/ i,i-1,i-1, j,j,j-1 /), CV_GPIdShape )
      CV_GPId(4,:,:) = reshape( (/ i,i-1,i, j,j-1,j-1 /), CV_GPIdShape )
      CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape )
      return
    else if( j == glEMax ) then
      ! 右端特異点
      CV_GPId(1,:,:) = reshape( (/ i,i-1,i-1, j,j,j-1 /), CV_GPIdShape )
      CV_GPId(2,:,:) = reshape( (/ i,i-1,i, j,j-1,j-1 /), CV_GPIdShape )
      CV_GPId(3,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape )
      CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape )
      CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape )
      return
    end if
  else
    ! 南半球に対して
    if ( i==glEMin .and. j==glEMin ) then
      ! 上端特異点
      CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape )
      CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape )
      CV_GPId(3,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape )
      CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape )
      CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape ) 
      return
    else if( j == glEMin ) then
      ! 左端特異点
      CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape )
      CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape )
      CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape )
      CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j,j+1 /), CV_GPIdShape )
      CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j+1,j+1 /), CV_GPIdShape )
      return
    else if( j == glEMax ) then
      ! 右端特異点
      CV_GPId(1,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape )
      CV_GPId(2,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape )
      CV_GPId(3,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape )
      CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j,j+1 /), CV_GPIdShape )
      CV_GPId(5,:,:) = reshape( (/ i,i+1,i, j,j+1,j+1 /), CV_GPIdShape )
      return
    end if
  end if

end function generate_CV5_GPindex

!
!> \~japanese
!! \brief (特異格子点(rcID,GP_i,GP_j)の)コントールボリュームの各頂点(全 6 個)の計算に必要な格子点の id の組を取得する.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!! @param[in]     GP_i 対象となる格子点(rcID, GP_i, GP_j) の GP_i.
!! @param[in]     GP_i 対象となる格子点(rcID, GP_i, GP_j) の GP_j.
!! @param[in]     rcID 対象となる格子点(rcID, GP_i, GP_j) の矩形領域のインデックス rcID.
!! @return CV_GPId コントールボリュームの各頂点の計算に必要な格子点の id の組.
!!
!! \~english
!! \brief (特異格子点(rcID,GP_i,GP_j)の)コントールボリュームの各頂点(全 6 個)の計算に必要な格子点の id の組を取得する.
!!
!! @param[in,out] self A variable of derived type IcGrid2D_FVM.
!! @param[in]    i    i of target grid point(rcID, i, j).
!! @param[in]    j    j of target grid point(rcID, i, j).
!! @param[in]    rcID The rectangle region index rcID of target grid point(rcID, i, j).
!! @return CV_GPId The set of grid-point ID nessesary for calculating the vertices of the control volume.
!!
function generate_CV6_GPindex( &
  & self,                      &  ! (inout)
  & GP_i, GP_j, rcID           &  ! (in)
  & ) result(CV_GPId)
 
  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(inout) :: self
  integer, intent(in) :: GP_i
  integer, intent(in) :: GP_j
  integer, intent(in) :: rcID
  integer :: CV_GPId(6,3,2)
  
  ! 作業変数
  ! Work variables
  !
  integer :: i,j, glEMin, glEMax 
  integer :: CV_GPIdShape(2) = (/ 3, 2 /)

  ! 実行文 ; Executable statements
  ! 

  glEMin = get_EffSize_Min(self)
  glEMax = get_EffSize_Max(self)
  
  i = GP_i; j = GP_j

  ! 先に, インデックス i,j が矩形領域の境界に当たる場合を処理する. 
  if( rcID <= 5 )then
    ! 北半球の矩形領域に対して
    if(  i == glEMin ) then
      ! 右上端
      CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j /), CV_GPIdShape )
      CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j,j-1 /), CV_GPIdShape )
      CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j-1,j-1 /), CV_GPIdShape )
      CV_GPId(4,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape )
      CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape )
      CV_GPId(6,:,:) = reshape( (/ i,i+1,i, j,j,j+1  /), CV_GPIdShape ) 
      return
    else if( j == glEMin ) then
      ! 左上端
      CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape )
      CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape )
      CV_GPId(3,:,:) = reshape( (/ i,i-1,i-1, j,j,j-1 /), CV_GPIdShape )
      CV_GPId(4,:,:) = reshape( (/ i,i-1,i, j,j-1,j-1 /), CV_GPIdShape )
      CV_GPId(5,:,:) = reshape( (/ i,i,i+1, j,j-1,j /), CV_GPIdShape )
      CV_GPId(6,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape )           
      return
    end if

  else
    ! 南半球の矩形領域に対して
    if( i == glEMax ) then
      ! 左下端
      CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape )
      CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape )
      CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape )
      CV_GPId(4,:,:) = reshape( (/ i,i,i+1, j,j-1,j /), CV_GPIdShape )
      CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j,j+1 /), CV_GPIdShape )
      CV_GPId(6,:,:) = reshape( (/ i,i+1,i, j,j+1,j+1 /), CV_GPIdShape ) 
      return
    else if( j == glEMax ) then
      ! 右下端
      CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j /), CV_GPIdShape )
      CV_GPId(2,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape )
      CV_GPId(3,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape )
      CV_GPId(4,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape )
      CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j,j+1 /), CV_GPIdShape )
!      CV_GPId(6,:,:) = reshape( (/ i,i,i+1, j,j+1,j+1 /), CV_GPIdShape )          
      CV_GPId(6,:,:) = reshape( (/ i,i+1,i, j,j+1,j+1 /), CV_GPIdShape )          
      return
    end if
  end if
  
  ! インデックス i,j が矩形領域の内部(境界以外)に対応する場合
  CV_GPId(1,:,:) = reshape( (/ i,i,i-1, j,j+1,j+1 /), CV_GPIdShape )
  CV_GPId(2,:,:) = reshape( (/ i,i-1,i-1, j,j+1,j /), CV_GPIdShape )
  CV_GPId(3,:,:) = reshape( (/ i,i-1,i, j,j,j-1 /), CV_GPIdShape )
  CV_GPId(4,:,:) = reshape( (/ i,i,i+1, j,j-1,j-1 /), CV_GPIdShape )
  CV_GPId(5,:,:) = reshape( (/ i,i+1,i+1, j,j-1,j /), CV_GPIdShape )
  CV_GPId(6,:,:) = reshape( (/ i,i+1,i, j,j,j+1 /), CV_GPIdShape ) 

end function generate_CV6_GPindex

!
!> \~japanese
!! \brief 格子点(rcID, i, j) が極に対応するかを調べる.
!!
!! 北極であった場合 NORTH_POLE_FLAG, 南極であった場合 SOUTH_POLE_FLAG, いずれでもない場合 NOT_POLE_FLAG
!! を返す.
!!
!! @param[in,out] self 構造型 IcGrid2D_FVM の変数.
!! @param[in]    rcID 対象となる格子点(rcID, i, j) の矩形領域のインデックス rcID.
!! @param[in]    i    対象となる格子点(rcID, i, j) の i.
!! @param[in]    j    対象となる格子点(rcID, i, j) の j.
!! @return       pole_flag 指定された格子点が極と一致しているか否かを表す定数.
!!
!! \~english
!! \brief Check whether the grid point(rcID, i, j) concides with the pole.
!!
!! 北極であった場合 NORTH_POLE_FLAG, 南極であった場合 SOUTH_POLE_FLAG, いずれでもない場合 NOT_POLE_FLAG
!! を返す
!!
!! @param[in,out] self  A variable of derived type IcGrid2D_FVM.
!! @param[in]     rcID  The rectangle region index rcID of target grid point(rcID, i, j).
!! @param[in]     i     i of target grid point(rcID, i, j).
!! @param[in]     j     j of target grid point(rcID, i, j).
!! @return        pole_flag The constant indicating whether the specified grid point concides with the pole.
!!
function check_pole(self, rcID, i, j) result (pole_flag)

  ! 宣言文 ; Declaration statements
  !
  type(IcGrid2D_FVM), intent(in) :: self
  integer, intent(in) :: rcID
  integer, intent(in) :: i
  integer, intent(in) :: j
  integer pole_flag

  ! 実行文 ; Executable statements
  !

  if( rcID <= 5 .and. i==self%idEfMin .and. j== self%idEfMin ) then
    pole_flag = NORTH_POLE_FLAG
  else if( rcID > 5 .and. i==self%idEfMax .and. j== self%idEfMax ) then
    pole_flag = SOUTH_POLE_FLAG
  else
    pole_flag = NOT_POLE_FLAG
  end if

end function check_pole

!
!> \~japanese
!! \brief 単位球に内包される正二十面体の頂点 12 個の座標を計算する.
!!
!! @param[in,out] orth_icvertex 12 個の頂点座標を保存する配列.
!!
!! \~english
!! \brief Calculates the coordinates of twelve vertices of an icosahedron
!! embedded into a unit sphere.
!!
!! @param[in,out] orth_icvertex The array in which the coordinates of twelve vertices will be stored.
!!
subroutine calc_icosahedron_vertex( &
  & orth_icvertex                   &  ! (inout)
  & )
  
  ! 宣言文 ; Declaration statements
  !

  real(DP), intent(inout) :: orth_icvertex(ICOSAHEDRON_VERTEX_NUM, 3)   

  ! 作業変数
  ! Work variable
  !
  real(DP) :: unit_radius = 1.0d0
  real(DP) :: del_lambda = 2.0d0 * PI / 5.0d0
  real(DP) :: pos2_theta
  real(DP) :: pos2_lambda 
  integer i

  ! 実行文 ; Executable statements
  ! 
  
  !　頂点1, 12 を z 軸上に配置する. 
  orth_icvertex(1, :) = (/ 0.0d0, 0.0d0, unit_radius /)
  orth_icvertex(12, :) = (/ 0.0d0, 0.0d0, - unit_radius /)
   
  ! 頂点 2 の緯度・経度を設定する. 
  pos2_lambda = 0.0d0 
  pos2_theta = 2.0d0 * dasin( 1.0d0 / (2.0d0 * dsin(0.2d0 * PI)) ) - PI/2.0d0
  
  ! 地理座標系における頂点 2, 7 の座標をデカルト座標系における座標へ変換する. 
  ! なお, 頂点 7 は頂点 2 と xy 平面を挟んで反対側に存在する. 
  orth_icvertex(2, :) = geo_to_orth_pos( (/ pos2_lambda, pos2_theta, unit_radius /) )
  orth_icvertex(7, :) = geo_to_orth_pos( (/ pos2_lambda + del_lambda/2.0d0, - pos2_theta, unit_radius /) )

  ! z 軸を軸にして, 座標が既知である頂点を del_lambda づつ回転させて, 残りの頂点座標を計算する.  
  do i=3, 6
    orth_icvertex(i, :) = rotateZ( orth_icvertex(i-1, :), del_lambda )
    orth_icvertex(i+5, :) = rotateZ( orth_icvertex( (i+5)-1, : ), del_lambda )
  end do

end subroutine calc_icosahedron_vertex

!
!> \~japanese 
!! \brief 正二十面格子の間隔の特徴的な長さを取得する. 
!!
!! @param[in] self
!! @return
!!
!! \~english 
!! \brief Get a characteristic length of icosahedral grid interval.
!!
!! The length is expressed using the grid division level m as
!! \f[
!!     \lambda = \dfrac{ 2 \Pi a}{10 \times 2^{m-1}},
!! \f]
!! which \f$a\f$ is the radius of sphere.
!! The denominator on right hand side in the above equation are the number of
!! grid points on the equator.
!!
!! @param[in] self
!! @return
!!
function get_IcGridMeanInterval( &
  & self                         &  ! (in)
  & ) result(meanDist)

  ! 宣言文 ; Declaration statements
  !

  type(IcGrid2D_FVM), intent(in) :: self
  real(DP) meanDist

  ! 宣言文 ; Declaration statements
  !

  meanDist = ( 2.0d0 * PI * get_IcRadius(self) ) &
    &          / ( 10.0d0 * 2.0d0**( get_glevel(self) - 1 ) )

end function get_IcGridMeanInterval

end module IcGrid2D_FVM_Manager
