!= 格子設定用モジュール ! ! Authors:: SUGIYAMA Ko-ichiro, ODAKA Masatsugu ! Version:: $Id: gridset.f90,v 1.14 2014/07/08 01:05:33 sugiyama Exp $ ! Tag Name:: $Name: arare5-20160127-2 $ ! Copyright:: Copyright (C) GFD Dennou Club, 2006. All rights reserved. ! License:: See COPYRIGHT[link:../../COPYRIGHT] ! module gridset ! != 格子設定用モジュール ! !引数に与えられた NAMELIST ファイルから, 格子点情報を取得し, !保管するための変数参照型モジュール ! !== Procedures List ! gridset_init :: 初期化ルーチン ! gridset_check :: 設定内容の整合性をチェックするためのルーチン ! !モジュール読み込み use dc_iounit, only: FileOpen use dc_message, only: MessageNotify use mpi_wrapper, only: myrank, nprocs use namelist_util, only: namelist_filename !暗黙の型宣言禁止 implicit none !save 属性 private !公開変数 integer, public, save :: NX = 10 ! x 方向格子点数 integer, public, save :: NY = 10 ! y 方向格子点数 integer, public, save :: NZ = 10 ! z 方向格子点数 integer, public, save :: NCMAX = 1 ! 組成配列要素数 integer, public, save :: Xmg = 2 ! x 方向糊代格子点数 integer, public, save :: Ymg = 2 ! y 方向糊代格子点数 integer, public, save :: Zmg = 3 ! z 方向糊代格子点数 integer, public, save :: imin ! x 方向の配列の下限 integer, public, save :: imax ! x 方向の配列の上限 integer, public, save :: jmin ! y 方向の配列の下限 integer, public, save :: jmax ! y 方向の配列の下限 integer, public, save :: kmin ! z 方向の配列の下限 integer, public, save :: kmax ! z 方向の配列の下限 integer, public, save :: xsub = 1 ! Number of MPI processes (X direction) integer, public, save :: ysub = 1 ! Number of MPI processes (Y direction) logical, public, save :: FlagCalc3D = .true. logical, public, save :: FlagCalcMoist = .true. integer, save :: xdim = 1 integer, save :: ydim = 1 integer, save :: zdim = 1 integer, public, save :: nxny = 100 ! x 方向と y 方向の総格子点数 ! サブルーチンの公開 public gridset_init contains subroutine gridset_init ! != 初期化ルーチン ! ! 設定ファイルから情報を読み込み格子点数を計算する ! !暗黙の型宣言禁止 implicit none !内部変数 integer :: unit !設定ファイル用装置番号 !----------------------------------------------------------------- ! 設定ファイルから情報を読み込み ! NAMELIST /gridset_nml/ xdim, ydim, zdim, NCMAX, Xmg, Ymg, Zmg, xsub, ysub call FileOpen(unit, file=namelist_filename, mode='r') read(unit, NML=gridset_nml) close(unit) !----------------------------------------------------------------- ! NX, NY, NZ を決める ! nx = xdim / xsub imin = 1 - xmg imax = nx + xmg if (ydim == 1) then ! 二次元 ymg = 0 ny = 1 FlagCalc3D = .false. else ! 三次元 ny = ydim / ysub FlagCalc3D = .true. end if jmin = 1 - ymg jmax = ny + ymg nz = zdim kmin = 1 - zmg kmax = nz + zmg nxny = nx * ny ! トレーサー計算のスイッチ ! if ( ncmax == 0 ) then ! FlagCalcTracer = .false. FlagCalcMoist = .false. else ! FlagCalcTracer = .true. FlagCalcMoist = .true. end if !----------------------------------------------------------------- !"myrank == 0" に該当する計算ノードが, 読み込んだ情報を出力 ! if (myrank == 0) then call MessageNotify( "M", "gridset_init", "xsub = %d", i=(/xsub/) ) call MessageNotify( "M", "gridset_init", "ysub = %d", i=(/ysub/) ) call MessageNotify( "M", "gridset_init", "xdim = %d", i=(/xdim/) ) call MessageNotify( "M", "gridset_init", "ydim = %d", i=(/ydim/) ) call MessageNotify( "M", "gridset_init", "zdim = %d", i=(/zdim/) ) call MessageNotify( "M", "gridset_init", "[1 node] NX = %d", i=(/NX/) ) call MessageNotify( "M", "gridset_init", "[1 node] NY = %d", i=(/NY/) ) call MessageNotify( "M", "gridset_init", "[1 node] NZ = %d", i=(/NZ/) ) call MessageNotify( "M", "gridset_init", "[1 node] NCMAX = %d", i=(/NCMAX/) ) call MessageNotify( "M", "gridset_init", "[1 node] xmg = %d", i=(/Xmg/) ) call MessageNotify( "M", "gridset_init", "[1 node] ymg = %d", i=(/Ymg/) ) call MessageNotify( "M", "gridset_init", "[1 node] zmg = %d", i=(/Zmg/) ) call MessageNotify( "M", "gridset_init", "[1 node] imin = %d", i=(/imin/) ) call MessageNotify( "M", "gridset_init", "[1 node] imax = %d", i=(/imax/) ) call MessageNotify( "M", "gridset_init", "[1 node] jmin = %d", i=(/jmin/) ) call MessageNotify( "M", "gridset_init", "[1 node] jmax = %d", i=(/jmax/) ) call MessageNotify( "M", "gridset_init", "[1 node] kmin = %d", i=(/kmin/) ) call MessageNotify( "M", "gridset_init", "[1 node] kmax = %d", i=(/kmax/) ) end if !----------------------------------------------------------------- ! 値のチェック ! call gridset_check end subroutine gridset_init !!!----------------------------------------------------------------- subroutine gridset_check ! != 設定内容の整合性をチェックするためのルーチン ! ! 設定内容に矛盾がないか確認をする. ! 問題が有る場合はエラーメッセージを出して終了する. ! ! nprocs が xsub or ysub で割り切れるかチェック. ! if ( mod( nprocs, xsub ) /= 0) then call MessageNotify( "E", "gridset_init", "mod( nprocs / xsub ) /= 0" ) end if if ( mod( nprocs, ysub ) /= 0) then call MessageNotify( "E", "gridset_init", "mod( nprocs / ysub ) /= 0" ) end if if (xsub * ysub /= nprocs) then call MessageNotify( "E", "gridset_init", "xsub * ysub /= nprocs" ) end if ! if 3D calculation, xsub > 1 and ysub > 1 ! ! if ( FlagCalc3D ) then ! if ( ysub <= 1 ) then ! call MessageNotify( "E", "gridset_init", "ysub has to be > 1" ) ! elseif ( xsub <= 1 ) then ! call MessageNotify( "E", "gridset_init", "xsub has to be > 1" ) ! end if ! end if ! X 方向のマージンの大きさのチェック ! if ( xdim < xmg ) then call MessageNotify( "E", "gridset_init", "xdim < Xmg" ) end if ! xdim が xsub で割り切れるか確認 ! if ( mod(xdim, xsub) /= 0 ) then call MessageNotify( "E", "gridset_init", "mod(xdim, xsub) /= 0" ) end if ! Y 方向のマージンの大きさのチェック ! if ( ydim < ymg ) then call MessageNotify( "E", "gridset_init", "ydim < Ymg" ) end if ! ydim が xsub で割り切れるか確認 ! if ( mod(ydim, ysub) /= 0 ) then call MessageNotify( "E", "gridset_init", "mod(ydim, ysub) /= 0" ) end if ! Z 方向のマージンの大きさのチェック ! if ( zdim < zmg ) then call MessageNotify( "E", "gridset_init", "zdim < Zmg" ) end if end subroutine gridset_check end module gridset