Class | gridset_3d |
In: |
setup/gridset_3d.f90
|
引数に与えられた NAMELIST ファイルから, 格子点情報を取得し, 保管するための変数型モジュール
Subroutine : | |
cfgfile : | character(*), intent(in) |
NAMELIST から情報を得て, 格子点を計算する
This procedure input/output NAMELIST#gridset .
subroutine gridset_init(cfgfile) ! !NAMELIST から情報を得て, 格子点を計算する ! !暗黙の型宣言禁止 implicit none !モジュール読み込み !入力変数 character(*), intent(in) :: cfgfile !内部変数 integer :: i, k, unit integer, parameter :: kind = 8 !精度を表す logical :: DebugOn !----------------------------------------------------------------- ! NAMELIST から情報を取得 !----------------------------------------------------------------- NAMELIST /gridset/ NX, NY, NZ, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax, Xmg, Ymg, Zmg, SpcNum call FileOpen(unit, file=cfgfile, mode='r') read(unit, NML=gridset) close(unit) if ( NX < xmg ) then call MessageNotify( "E", "gridset_init", "NX < Xmg" ) end if if ( NY < ymg ) then call MessageNotify( "E", "gridset_init", "NY < Ymg" ) end if if ( NZ < zmg ) then call MessageNotify( "E", "gridset_init", "NZ < Zmg" ) end if !----------------------------------------------------------------- ! 格子点間隔計算 !----------------------------------------------------------------- call xyz_axis_init(NX, NY, NZ, xmg, ymg, zmg, Xmin, Xmax, Ymin, Ymax, Zmin, Zmax) !----------------------------------------------------------------- ! 物理的に意味のある領域の上限・下限を決める !-----------------------------------------------------------------; RegXMin = 1 RegXMax = NX RegYMin = 1 RegYMax = NY RegZMin = 1 RegZMax = NZ !----------------------------------------------------------------- ! デバッグモードか否かで, ファイル出力に用いる変数の大きさを変える !----------------------------------------------------------------- DebugOn= Debug() if (DebugOn) then FileNX = size(p_X, 1) FileNY = size(q_Y, 1) FileNZ = size(r_Z, 1) FileXMin = DimXMin FileXMax = DimXMax FileYmin = DimYMin FileYMax = DimYMax FileZMin = DimZMin FileZMax = DimZMax else FileNX = NX FileNY = NY FileNZ = NZ FileXMin = RegXMin FileXMax = RegXMax FileYMin = RegYMin FileYMax = RegYMax FileZMin = RegZMin FileZMax = RegZMax end if !----------------------------------------------------------------- ! 確認 !----------------------------------------------------------------- call MessageNotify( "M", "gridset_init", "XMin = %f", d=(/XMin/) ) call MessageNotify( "M", "gridset_init", "XMax = %f", d=(/XMax/) ) call MessageNotify( "M", "gridset_init", "YMin = %f", d=(/YMin/) ) call MessageNotify( "M", "gridset_init", "YMax = %f", d=(/YMax/) ) call MessageNotify( "M", "gridset_init", "ZMin = %f", d=(/ZMin/) ) call MessageNotify( "M", "gridset_init", "ZMax = %f", d=(/ZMax/) ) call MessageNotify( "M", "gridset_init", "NX = %d", i=(/NX/) ) call MessageNotify( "M", "gridset_init", "NY = %d", i=(/NY/) ) call MessageNotify( "M", "gridset_init", "NZ = %d", i=(/NZ/) ) call MessageNotify( "M", "gridset_init", "Xmargin = %d", i=(/Xmg/) ) call MessageNotify( "M", "gridset_init", "Ymargin = %d", i=(/Ymg/) ) call MessageNotify( "M", "gridset_init", "Zmargin = %d", i=(/Zmg/) ) call MessageNotify( "M", "gridset_init", "DimXMin = %d", i=(/DimXMin/) ) call MessageNotify( "M", "gridset_init", "DimXMax = %d", i=(/DimXMax/) ) call MessageNotify( "M", "gridset_init", "DimYMin = %d", i=(/DimYMin/) ) call MessageNotify( "M", "gridset_init", "DimYMax = %d", i=(/DimYMax/) ) call MessageNotify( "M", "gridset_init", "DimZMin = %d", i=(/DimZMin/) ) call MessageNotify( "M", "gridset_init", "DimZMax = %d", i=(/DimZMax/) ) end subroutine gridset_init