Class gridset
In: setup/gridset.f90

引数に与えられた NAMELIST ファイルから, 格子点情報を取得し, 保管するための変数型モジュール

Methods

DelX   DelZ   DimXMax   DimXMin   DimZMax   DimZMin   FileNX   FileNZ   FileXMax   FileXMin   FileZMax   FileZMin   MarginX   MarginZ   NX   NZ   RegXMax   RegXMin   RegZMax   RegZMin   SpcNum   Xmax   Xmin   Zmax   Zmin   f_X   f_Z   gridset_init   gridset_init2   gridset_set   p_X   r_Z   s_X   s_Z   x_X   z_Z  

Included Modules

dc_message DebugSet

Public Instance methods

DelX
Variable :
DelX :real(8)
: 格子間隔
DelZ
Variable :
DelZ :real(8)
: 格子間隔
DimXMax
Variable :
DimXMax :integer
: x 方向の配列の上限
DimXMin
Variable :
DimXMin :integer
: x 方向の配列の下限
DimZMax
Variable :
DimZMax :integer
: z 方向の配列の上限
DimZMin
Variable :
DimZMin :integer
: z 方向の配列の下限
FileNX
Variable :
FileNX :integer
: ファイル出力用
FileNZ
Variable :
FileNZ :integer
: ファイル出力用
FileXMax
Variable :
FileXMax :integer
: ファイル出力用
FileXMin
Variable :
FileXMin :integer
: ファイル出力用
FileZMax
Variable :
FileZMax :integer
: ファイル出力用
FileZMin
Variable :
FileZMin :integer
: ファイル出力用
MarginX
Variable :
MarginX :integer
: 境界のグリッド数
MarginZ
Variable :
MarginZ :integer
: 境界のグリッド数
NX
Variable :
NX :integer
: 格子点数
 integer, parameter    :: MarginX = 5   !境界のグリッド数
 integer, parameter    :: MarginZ = 5   !境界のグリッド数
NZ
Variable :
NZ :integer
: 格子点数
 integer, parameter    :: MarginX = 5   !境界のグリッド数
 integer, parameter    :: MarginZ = 5   !境界のグリッド数
RegXMax
Variable :
RegXMax :integer
: x 方向の物理領域の上限
RegXMin
Variable :
RegXMin :integer
: x 方向の物理領域の下限
RegZMax
Variable :
RegZMax :integer
: z 方向の物理領域の上限
RegZMin
Variable :
RegZMin :integer
: z 方向の物理領域の下限
SpcNum
Variable :
SpcNum :integer
: 化学種の数
Xmax
Variable :
Xmax :real(8)
: x 座標の始点・終点
Xmin
Variable :
Xmin :real(8)
: x 座標の始点・終点
Zmax
Variable :
Zmax :real(8)
: z 座標の始点・終点
Zmin
Variable :
Zmin :real(8)
: z 座標の始点・終点
f_X
Variable :
f_X(:) :real(8), allocatable
: X 座標軸(ベクトル格子点)
f_Z
Variable :
f_Z(:) :real(8), allocatable
: Z 座標軸(ベクトル格子点)
Subroutine :
cfgfile :character(*), intent(in)

NAMELIST から情報を得て, 格子点を計算する

This procedure input/output NAMELIST#gridset .

[Source]

  subroutine gridset_init(cfgfile)
    !
    !NAMELIST から情報を得て, 格子点を計算する
    !
    
    !暗黙の型宣言禁止
    implicit none

    !入力変数
    character(*), intent(in) :: cfgfile

    !内部変数
    logical            :: flag

    !-----------------------------------------------------------------
    ! NAMELIST から情報を取得
    !-----------------------------------------------------------------
    NAMELIST /gridset/ NX, NZ, Xmin, Xmax, Zmin, Zmax, SpcNum
    open (10, FILE=cfgfile)
    read(10, NML=gridset)
    close(10)

    MarginX = 5
    MarginZ = 5

    flag = DebugOn
    call gridset_set(flag)

  end subroutine gridset_init
Subroutine :
im :integer, intent(in)
km :integer, intent(in)
X1 :real(8), intent(in)
X2 :real(8), intent(in)
Z1 :real(8), intent(in)
Z2 :real(8), intent(in)
mgn :integer, intent(in)
sNum :integer, intent(in)

格子点を計算する

[Source]

  subroutine gridset_init2( im, km, X1, X2, Z1, Z2, mgn, sNum )
    !
    !格子点を計算する
    !

    !暗黙の型宣言禁止
    implicit none

    !入力変数
    real(8), intent(in) :: X1, X2, Z1, Z2
    integer, intent(in) :: im, km, sNum, mgn
    logical             :: flag

    !内部変数
    integer            :: i, k
    integer, parameter :: kind = 8      !精度を表す

    !-----------------------------------------------------------------
    ! 引数から情報を取得
    !-----------------------------------------------------------------
    NX = im
    NZ = km
    SpcNum = sNum
    XMin = X1
    XMax = X2
    ZMin = Z1
    ZMax = Z2
    MarginX = mgn
    MarginZ = mgn

    flag = .false.
    call gridset_set( flag )
    
  end subroutine gridset_init2
Subroutine :
flag :logical, intent(in)

[Source]

  subroutine gridset_set(flag)
    
    implicit none

    !内部変数
    logical, intent(in) :: flag
    integer, parameter  :: kind = 8      !精度を表す
    integer             :: i,k
    
    !-----------------------------------------------------------------
    ! 格子点間隔計算
    !-----------------------------------------------------------------
    DelX = (Xmax - Xmin) / real(NX, kind)
    DelZ = (Zmax - Zmin) / real(NZ, kind)

    !-----------------------------------------------------------------
    ! 物理的に意味のある領域の上限・下限を決める
    !-----------------------------------------------------------------;
    RegXMin = 0
    RegXMax = NX
    RegZMin = 0
    RegZMax = NZ

!    !-----------------------------------------------------------------
!    ! 物理的に意味のある領域の上限・下限を決める
!    !-----------------------------------------------------------------;
!    RegXMin = 1
!    RegXMax = NX
!    RegZMin = 1
!    RegZMax = NZ
!
!    !-----------------------------------------------------------------
!    ! 壁に対応する領域の配列添え字を保管
!    !-----------------------------------------------------------------;
!    WallXMin = 0
!    WallXMax = NX
!    WallZMin = 0
!    WallZMax = NZ

    !-----------------------------------------------------------------
    ! 配列の上限・下限を決める
    !-----------------------------------------------------------------
    DimXMin = RegXMin - MarginX
    DimXMax = RegXMax + MarginX
    DimZMin = RegZMin - MarginZ
    DimZMax = RegZMax + MarginZ

    !-----------------------------------------------------------------
    ! グリッドの設定
    !-----------------------------------------------------------------
    allocate( f_X(DimXMin:DimXMax), f_Z(DimZMin:DimZMax) )
    allocate( p_X(DimXMin:DimXMax), r_Z(DimZMin:DimZMax) )
    do i = DimXMin, DimXMax
      f_X(i) = Xmin + DelX * ( real(i, kind) - RegXMin )
    end do
    do k = DimZMin, DimZMax
      f_Z(k) = Zmin + DelZ * ( real(k, kind) - RegZMin )
    end do
    
    p_X = f_X
    r_Z = f_Z
    
    !-----------------------------------------------------------------
    ! 半格子ずれたグリッドの設定. 等間隔を前提にしている.
    !-----------------------------------------------------------------
    allocate( s_X(DimXMin: DimXMax), s_Z(DimZMin: DimZMax) )
    allocate( x_X(DimXMin: DimXMax), z_Z(DimZMin: DimZMax) )
    s_X(DimXMin+1:DimXMax) = ( f_X(DimXMin+1:DimXMax) + f_X(DimXMin:DimXMax-1) ) * 5.0d-1
    s_X(DimXMin) = s_X(DimXMin+1) - DelX
    
    s_Z(DimZMin+1:DimZMax) = ( f_Z(DimZMin+1:DimZMax) + f_Z(DimZMin:DimZMax-1) ) * 5.0d-1
    s_Z(DimZMin) = s_Z(DimZMin+1) - DelZ

    x_X = s_X
    z_Z = s_Z
    !-----------------------------------------------------------------
    ! デバッグモードか否かで, ファイル出力に用いる変数の大きさを変える
    !-----------------------------------------------------------------
    if (flag) then
      FileNX   = size(f_X, 1)
      FileNZ   = size(f_Z, 1)
      FileXMin = DimXMin
      FileXMax = DimXMax
      FileZMin = DimZMin
      FileZMax = DimZMax
    else
      FileNX   = NX
      FileNZ   = NZ
      FileXMin = RegXMin + 1
      FileXMax = RegXMax
      FileZMin = RegZMin + 1
      FileZMax = RegZMax
    end if

    !-----------------------------------------------------------------    
    ! 確認
    !-----------------------------------------------------------------
    if (cpurank == 0) then
      call MessageNotify( "M", "gridset_init", "XMin = %f", d=(/XMin/)    )
      call MessageNotify( "M", "gridset_init", "XMax = %f", d=(/XMax/)    )
      call MessageNotify( "M", "gridset_init", "ZMin = %f", d=(/ZMin/)    )
      call MessageNotify( "M", "gridset_init", "ZMax = %f", d=(/ZMax/)    )
      call MessageNotify( "M", "gridset_init", "DelX = %f", d=(/DelX/)    )
      call MessageNotify( "M", "gridset_init", "DelZ = %f", d=(/DelZ/)    )
      call MessageNotify( "M", "gridset_init", "NX = %d",   i=(/NX/)      )
      call MessageNotify( "M", "gridset_init", "NZ = %d",   i=(/NZ/)      )
      call MessageNotify( "M", "gridset_init", "SpcNum = %d", i=(/SpcNum/) )
      call MessageNotify( "M", "gridset_init", "MarginX = %d", i=(/MarginX/) )
      call MessageNotify( "M", "gridset_init", "MarginZ = %d", i=(/MarginZ/) )
      call MessageNotify( "M", "gridset_init", "DimXMin = %d", i=(/DimXMin/) )
      call MessageNotify( "M", "gridset_init", "DimXMax = %d", i=(/DimXMax/) )
      call MessageNotify( "M", "gridset_init", "DimZMin = %d", i=(/DimZMin/) )
      call MessageNotify( "M", "gridset_init", "DimZMax = %d", i=(/DimZMax/) )
      call MessageNotify( "M", "gridset_init", "RegXMin = %d", i=(/RegXMin/) )
      call MessageNotify( "M", "gridset_init", "RegXMax = %d", i=(/RegXMax/) )
      call MessageNotify( "M", "gridset_init", "RegZMin = %d", i=(/RegZMin/) )
      call MessageNotify( "M", "gridset_init", "RegZMax = %d", i=(/RegZMax/) )
      call MessageNotify( "M", "gridset_init", "FileNX = %d",  i=(/FileNX/))
      call MessageNotify( "M", "gridset_init", "FileNZ = %d",  i=(/FileNZ/))
      call MessageNotify( "M", "gridset_init", "FileXMin = %d", i=(/FileXMin/))
      call MessageNotify( "M", "gridset_init", "FileXMax = %d", i=(/FileXMax/))
      call MessageNotify( "M", "gridset_init", "FileZMin = %d", i=(/FileZMin/))
      call MessageNotify( "M", "gridset_init", "FileZMax = %d", i=(/FileZMax/))
    end if

  end subroutine gridset_set
p_X
Variable :
p_X(:) :real(8), allocatable
: X 座標軸(ベクトル格子点)
r_Z
Variable :
r_Z(:) :real(8), allocatable
: Z 座標軸(ベクトル格子点)
s_X
Variable :
s_X(:) :real(8), allocatable
: X 座標軸(スカラー格子点)
s_Z
Variable :
s_Z(:) :real(8), allocatable
: Z 座標軸(スカラー格子点)
x_X
Variable :
x_X(:) :real(8), allocatable
: X 座標軸(スカラー格子点)
z_Z
Variable :
z_Z(:) :real(8), allocatable
: Z 座標軸(スカラー格子点)