Class | ReStartFileIO |
In: |
io/restartfileio.f90
|
Copyright (C) GFD Dennou Club, 2004, 2005, 2006. All rights reserved.
* Developer: SUGIYAMA Ko-ichiro (sugiyama@gfd-dennou.org) * Version: $Id: restartfileio.f90,v 1.1.1.1 2006/04/25 03:43:58 deepconv Exp $ * Tag Name: $Name: $ * Change History:
リスタート用の場の情報を netCDF ファイルに出力するためのルーチン
速度を評価する格子位置と, 軸として入力してある座標値は整合的でないことに注意
subroutine ReStartFile_Get ( ReStartTime, pz_VelX_b, xr_VelZ_b, xz_Exner_b, xz_PotTemp_b, xz_Km_b, xz_MixRt_b, pz_VelX_n, xr_VelZ_n, xz_Exner_n, xz_PotTemp_n, xz_Km_n, xz_MixRt_n ) ! !リスタートファイルから情報取得 ! !モジュール読み込み use gt4_history use dc_message ,only: MessageNotify use timeset, only: DelTimeLong use fileset, only: InitFile use gridset, only: DimXMin, DimXMax, DimZMin, DimZMax, SpcNum !暗黙の型宣言禁止 implicit none !変数定義 real(8), intent(out) :: ReStartTime(2) real(8), intent(out) :: pz_VelX_n(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xr_VelZ_n(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xz_Exner_n(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xz_PotTemp_n(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xz_Km_n(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xz_MixRt_n(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) real(8), intent(out) :: pz_VelX_b(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xr_VelZ_b(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xz_Exner_b(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xz_PotTemp_b(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xz_Km_b(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(out) :: xz_MixRt_b(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) real(8) :: DelTime real(8) :: Var3D(DimXMin:DimXMax, DimZMin:DimZMax, 2) ! real(8) :: Var4D(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum, 2) character(30) :: name !変数名 !------------------------------------------------------------- !Get a Value from netCDF File !------------------------------------------------------------- name = "t" call HistoryGet( InitFile, name, ReStartTime ) !時間刻みのチェック DelTime = ReStartTime(2) - ReStartTime(1) if ( DelTime /= real(DelTimeLong, 4) ) then write(*,*) DelTime, DelTimeLong call MessageNotify("Error", "getdistrbvar", & & "DelTime for InitFile is not the same as DelTimeLong") stop end if !------------------------------------------------------------- ! Get a Value from netCDF File !------------------------------------------------------------- name = "VelX" call HistoryGet( InitFile, name, Var3D ) pz_VelX_b = Var3D(:,:,1) pz_VelX_n = Var3D(:,:,2) name = "VelZ" call HistoryGet( InitFile, name, Var3D ) xr_VelZ_b = Var3D(:,:,1) xr_VelZ_n = Var3D(:,:,2) name = "Exner" call HistoryGet( InitFile, name, Var3D ) xz_Exner_b = Var3D(:,:,1) xz_Exner_n = Var3D(:,:,2) name = "PotTemp" call HistoryGet( InitFile, name, Var3D ) xz_PotTemp_b = Var3D(:,:,1) xz_PotTemp_n = Var3D(:,:,2) name = "Km" call HistoryGet( InitFile, name, Var3D ) xz_Km_b = Var3D(:,:,1) xz_Km_n = Var3D(:,:,2) name = "MixRt" ! call HistoryGet( InitFile, name, Var4D ) ! xz_MixRt_b = Var4D(:,:,:,1) ! xz_MixRt_n = Var4D(:,:,:,2) end subroutine
subroutine ReStartFile_Open ( ) ! !リスタートファイルの書き出し ! !モジュール読み込み use gt4_history use gridset, only: s_X, s_Z, SpcNum use fileset, only: exptitle, expsrc, expinst, ReStartFile use basicset, only: xz_ExnerBasicZ, xz_DensBasicZ, & & xz_PotTempBasicZ, xz_VelSoundBasicZ, & & xz_PressBasicZ, xz_TempBasicZ, & & xz_MixRtBasicZ, xz_EffMolWtBasicZ, & & SpcWetID !暗黙の型宣言禁止 implicit none !変数定義 integer :: N, M N = size(s_X, 1) M = size(s_Z, 1) !------------------------------------------------------------- ! ヒストリー作成 !------------------------------------------------------------- call HistoryCreate( & & file = ReStartFile, & & title = exptitle, & & source = expsrc, & & institution = expinst, & & dims=(/'x','z','s','t'/), & & dimsizes=(/N, M, SpcNum, 2/), & !時間方向は配列要素 2 つ & longnames=(/'X-coordinate', & & 'Z-coordinate', & & 'Species Num ', & & 'Time '/), & & units=(/'m','m','1','t'/), origin=0.0, & & interval=0.0 ) !------------------------------------------------------------- ! 変数出力 !------------------------------------------------------------- call HistoryPut('x', s_X ) call HistoryPut('z', s_Z ) call HistoryPut('s', real(SpcWetID, 4) ) !無次元圧力の基本場 call HistoryAddVariable( & & varname='ExnerBasicZ', dims=(/'x','z'/), & & longname='nondimensional pressure', units='1',& & xtype='double' ) !温位の基本場 call HistoryAddVariable( & & varname='PotTempBasicZ', dims=(/'x','z'/), & & longname='potential temperature', & & units='K', xtype='double' ) !密度の基本場 call HistoryAddVariable( & & varname='DensBasicZ', dims=(/'x','z'/), & & longname='density', & & units='Kg/m^3', xtype='double' ) !音波速度の基本場 call HistoryAddVariable( & & varname='VelSoundBasicZ', dims=(/'x','z'/), & & longname='sound velocity', & & units='m/s|2', xtype='double' ) !温度の基本場 call HistoryAddVariable( & & varname='TempBasicZ', dims=(/'x','z'/), & & longname='Temperature of basic state', & & units='K', xtype='double' ) !圧力の基本場 call HistoryAddVariable( & & varname='PressBasicZ', dims=(/'x','z'/), & & longname='Pressure of basic state', & & units='Pa', xtype='double' ) !水蒸気混合比の基本場 call HistoryAddVariable( & & varname='MixRtBasicZ', dims=(/'x','z','s'/), & & longname='Mixing ratio of Condensible volatiles', & & units='kg/kg', xtype='double' ) !分子量効果 call HistoryAddVariable( & & varname='EffMolWtBasicZ', dims=(/'x','z'/), & & longname='Effect of Mole Weight', & & units='1', xtype='double' ) !無次元圧力 call HistoryAddVariable( & & varname='Exner', dims=(/'x','z','t'/), & & longname='nondimensional pressure', & & units='1', & & xtype='double' ) !温位の擾乱 call HistoryAddVariable( & & varname='PotTemp', dims=(/'x','z','t'/), & & longname='virtual potential temperature', & & units='K', & & xtype='double' ) !速度 call HistoryAddVariable( & & varname='VelX', dims=(/'x','z','t'/), & & longname='zonal velocity', & & units='m/s', & & xtype='double' ) !速度 call HistoryAddVariable( & & varname='VelZ', dims=(/'x','z','t'/), & & longname='vertical velocity', & & units='m/s', & & xtype='double' ) !渦粘性係数 call HistoryAddVariable( & & varname='Km', dims=(/'x','z','t'/), & & longname='Km', & & units='1', & & xtype='double' ) !水蒸気混合比 call HistoryAddVariable( & & varname='MixRt', dims=(/'x','z','s','t'/), & & longname='Mixing Ratio of Water Vaper', & & units='1', & & xtype='double' ) !------------------------------------------------------------- ! 基本場のファイル出力 !------------------------------------------------------------- call HistoryPut( 'DensBasicZ', xz_DensBasicZ ) call HistoryPut( 'ExnerBasicZ', xz_ExnerBasicZ ) call HistoryPut( 'PotTempBasicZ', xz_PotTempBasicZ ) call HistoryPut( 'VelSoundBasicZ', xz_VelSoundBasicZ ) call HistoryPut( 'TempBasicZ', xz_TempBasicZ ) call HistoryPut( 'PressBasicZ', xz_PressBasicZ ) call HistoryPut( 'MixRtBasicZ', xz_MixRtBasicZ ) call HistoryPut( 'EffMolWtBasicZ', xz_EffMolWtBasicZ ) end subroutine
subroutine ReStartFile_OutPut ( Time, pz_VelX, xr_VelZ, xz_Exner, xz_PotTemp, xz_Km, xz_MixRt ) ! !リスタートファイルに予報変数を書き出す ! !モジュール読み込み use gt4_history use gridset, only: DimXMin, DimXMax, DimZMin, DimZMax, SpcNum !暗黙の型宣言禁止 implicit none !変数定義 real(8), intent(in) :: Time real(8), intent(in) :: pz_VelX(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(in) :: xr_VelZ(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(in) :: xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(in) :: xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(in) :: xz_Km(DimXMin:DimXMax, DimZMin:DimZMax) real(8), intent(in) :: xz_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum) !------------------------------------------------------------------ ! ファイル出力 !------------------------------------------------------------------ call HistoryPut( 't', Time ) call HistoryPut( 'VelX', pz_VelX ) call HistoryPut( 'VelZ', xr_VelZ ) call HistoryPut( 'Exner', xz_Exner ) call HistoryPut( 'PotTemp', xz_PotTemp ) call HistoryPut( 'Km', xz_Km ) call HistoryPut( 'MixRt', xz_MixRt ) end subroutine