!= Module GridSet_3D ! ! Authors:: SUGIYAMA Ko-ichiro, ODAKA Masatsugu ! Version:: $Id: gridset_3d.f90,v 1.2 2007-08-06 15:47:07 odakker Exp $ ! Tag Name:: $Name: $ ! Copyright:: Copyright (C) GFD Dennou Club, 2006. All rights reserved. ! License:: See COPYRIGHT[link:../../COPYRIGHT] ! !== Overview ! !引数に与えられた NAMELIST ファイルから, 格子点情報を取得し, !保管するための変数型モジュール ! !== Error Handling ! !== Known Bugs ! !== Note ! !== Future Plans ! module gridset_3d ! !引数に与えられた NAMELIST ファイルから, 格子点情報を取得し, !保管するための変数型モジュール ! !モジュール読み込み use dc_types, only: DP use dc_iounit, only: FileOpen use dc_trace, only: Debug use dc_message, only: MessageNotify use xyz_base_module, only: DimXmin => imin, DimXMax => imax, & & DimYmin => jmin, DimYMax => jmax, & & DimZmin => kmin, DimZMax => kmax, & & xyz_axis_init, x_X, p_X, y_Y, q_Y, z_Z, r_Z, & & x_dx, p_dx, y_dy, q_dy, z_dz, r_dz, & & xyz_X, xyz_Y, xyz_Z, & & xyz_dX, xyz_dY, xyz_dZ, & & xmargin, ymargin, zmargin !暗黙の型宣言禁止 implicit none !save 属性 save !公開変数 real(DP) :: Xmin, Xmax ! x 座標の始点・終点 real(DP) :: Ymin, Ymax ! x 座標の始点・終点 real(DP) :: Zmin, Zmax ! z 座標の始点・終点 integer :: NX, NY, NZ ! 格子点数 integer :: Xmg, Ymg, Zmg ! 糊代格子点数 integer :: SpcNum ! 化学種の数 integer :: RegXMin ! x 方向の物理領域の下限 integer :: RegXMax ! x 方向の物理領域の上限 integer :: RegYMin ! x 方向の物理領域の下限 integer :: RegYMax ! x 方向の物理領域の上限 integer :: RegZMin ! z 方向の物理領域の下限 integer :: RegZMax ! z 方向の物理領域の上限 integer :: FileNX !ファイル出力用 integer :: FileNY !ファイル出力用 integer :: FileNZ !ファイル出力用 integer :: FileXMin !ファイル出力用 integer :: FileXMax !ファイル出力用 integer :: FileYMin !ファイル出力用 integer :: FileYMax !ファイル出力用 integer :: FileZMin !ファイル出力用 integer :: FileZMax !ファイル出力用 contains 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 end module gridset_3d