!= 木星大気を想定した仮想惑星実験
!
!= Virtual planet experiment like Jovian atmosphere
!
! Authors:: Yasuhiro MORIKAWA
! Version:: $Id: phy_jupiter_case00.f90,v 1.6 2008-06-09 16:58:24 morikawa Exp $
! Tag Name:: $Name: dcpam4-20080626 $
! Copyright:: Copyright (C) GFD Dennou Club, 2008. All rights reserved.
! License:: See COPYRIGHT[link:../../../COPYRIGHT]
!
module phy_jupiter_case00
!
!= 木星大気を想定した仮想惑星実験
!
!= Virtual planet experiment like Jovian atmosphere
!
! Note that Japanese and English are described in parallel.
!
! 木星大気を念頭においた, 仮想惑星実験を行います.
! 実際にこのモジュールで計算されている物理過程は以下の通りです.
!
! * 木星大気を模した放射冷却の効果
!
! Sugiyama et al. (2008) の手法を真似て,
! 2.0e5 〜 1.0e4 Pa の範囲に一定の冷却 (-1.0 K/day) を与えます.
!
! * 鉛直拡散 (Mellor and Yamada, 1974, レベル 2)
! * 速度を減衰させるためのスポンジ層
!
! 1.0e4 Pa 以下の層へ, 速度の擾乱成分を潰すための減衰の効果.
! 時定数 1 日で東西平均値に合わせます.
!
! * 湿潤対流調節 (Manabe et al., 1965)
! * 乾燥対流調節
! * 診断型大規模凝結
!
! * 最下層の温度と比湿を固定
!
! 最下層の比湿を 6.11641e-3 に固定します (Sugiyama, 2008 より) .
! 最下層の温度を $ 490 \sigma_1^{R/C_p} $ に固定します (Sugiyama, 2008 より).
! ここで, $ \sigma_1 $ [1] は最下層の中心,
! $ R $ [J kg-1 K-1] は乾燥大気の気体定数,
! $ C_p $ [J kg-1 K-1] は乾燥大気の定圧比熱です.
!
! A virtual planet experiment like Jovian atmosphere is performed.
! Physical processes calculated in this module is as follows.
!
! * Radiation cooling to imitate the Jupiter atmosphere
!
! The method of Sugiyama et al. (2008) is mimicked.
! A constant cooling (-1.0 K/day) is given within the range of
! 2.0e5 -- 1.0e4 Pa.
!
! * Vertical diffusion
! * Sponge layer for damping velocity
!
! Effect of damping perturbation of velocity to layers under 1.0e4 Pa.
! Effect of damping is decided to
! Strength of damping is decided to match to east and west mean
! value by a time constant (1 day).
!
! * Moist convective adjustment
! * Dry convective adjustment
! * Large scale condensation
!
! * Temperature and specific humidity on bottom layer
! are fixed
!
! Specific humidity on bottom layer is fixed to 6.11641e-3 (for Sugiyama, 2008) .
! Temperature on bottom layer is fixed to $ 490 \sigma_1^{R/C_p} $
! (for Sugiyama, 2008),
! where $ \sigma_1 $ [1] is center of botton layer,
! $ R $ [J kg-1 K-1] is specific heat of dry air at constant pressure,
! $ C_p $ [J kg-1 K-1] gas constant of dry air.
!
!== Procedures List
!
! PhyJpCreate :: PHYJP 型変数の初期設定
! PhysicsJupiter :: 木星大気を念頭においた仮想惑星の物理過程
! (放射, 鉛直拡散, スポンジ層) の演算
! PhysicsJupiterAdjust :: 木星大気を念頭においた仮想惑星の物理過程
! (湿潤対流調節, 診断型大規模凝結, 乾燥対流調節,
! 最下層の温度と比湿の固定) の計算
! PhyJpClose :: PHYJP 型変数の終了処理
! PhyJpPutLine :: PHYJP 型変数に格納されている情報の印字
! PhyJpinitialized :: PHYJP 型変数が初期設定されているか否か
! PhyJpSetTime :: 時刻の設定
! ------------ :: ------------
! PhyJpCreate :: Constructor of "PHYJP"
! PhysicsJupiter :: Calculation of physical processes of
! virtual planet experiment like Jovian atmosphere
! (radiation, vertical diffusion, sponge layer)
! PhysicsJupiterAdjust :: Calculate some physical processes
! (moist convective adjustment,
! large scale condensation,
! dry convective adjustment,
! fix of temperature and specific humidity on bottom layer)
! PhyJpClose :: Deconstructor of "PHYJP"
! PhyJpPutLine :: Print information of "PHYJP"
! PhyJpinitialized :: Check initialization of "PHYJP"
! PhyJpSetTime :: Configure time
!
!== Usage
!
! 始めに, PHYJP 型の変数を定義し,
! PhyJpCreate で初期設定を行います.
! PhysicsJupiter は東西風速, 南北風送, 温度, 比湿, 地表面気圧
! を受け取り, 物理過程 (放射など) の効果として
! 東西風速, 南北風送, 温度, 比湿の変化量を返します.
! PhysicsJupiterAdjust に温度, 比湿, 地表面気圧を与える
! ことで, 積雲パラメタリゼーションなどにより調節を行って
! 温度, 比湿, 地表面気圧を返します.
! PHYJP 型の変数の終了処理には
! PhyJpClose を用いてください.
!
! First, initialize "PHYJP" by "PhyJpCreate".
! "PhysicsJupiter" receives zonal and meridional wind,
! and temperature, specific humidity, surface pressure, and
! returns tendency of zonal and meridional wind,
! and temperature, specific humidity with physical processes
! (radiation etc.)
! "PhysicsJupiterAdjust" receives temperature, specific humidity,
! surface pressure, and adjusts them with cumulus parameterization etc.,
! and returns them.
! In order to terminate "PHYJP", use "PhyJpClose".
!
!== References
!
! * Manabe, S., Smagorinsky, J., Strickler, R.F., 1965:
! Simulated climatology of a general circulation model
! with a hydrologic cycle.
! Mon. Weather Rev., 93, 769--798.
!
! * Mellor, G. L, Yamada, T., 1974:
! A hierarchy of turbulence closure models for planetary boundary layers.
! J. Atmos. Sci., 31, 1791--1806.
!
use dc_types, only: DP, TOKEN, STRING
use dc_date_types, only: DC_DIFFTIME
use gt4_history_nmlinfo, only: GTHST_NMLINFO
use phy_ground, only: PHYGRD
use phy_interpolate, only: PHYINTPOL
use phy_cooling_s08, only: PHYCOOL
use phy_implicit, only: PHYIMPL
use phy_verdiff, only: PHYVDIF
use phy_sponge_layer, only: PHYSPOLAY
use phy_cumulus_adjust, only: PHYCUMAD
use phy_lscond, only: PHYLSC
use phy_neg_moist, only: PHYNEGMST
use phy_dryconv_adjust, only: PHYDRYCONV
use phy_adjust_fixed, only: PHYADJFIX
implicit none
private
public:: PHYJP
public:: PhyJpCreate, PhysicsJupiter, PhysicsJupiterAdjust
public:: PhyJpClose, PhyJpPutLine
public:: PhyJpInitialized, PhyJpSetTime
type PHYJP
!
! まず, PhyJpCreate で "PHYJP" 型の
! 変数を初期設定して下さい.
! 初期設定された "PHYJP" 型の変数を再度利用する際には,
! PhyJpClose によって終了処理を行ってください.
!
! Initialize "PHYJP" variable by
! "PhyJpCreate" before usage.
! If you reuse "PHYJP" variable again for another application,
! terminate by "PhyJpClose".
!
logical:: initialized = .false.
! 初期設定フラグ.
! Initialization flag
!-----------------------------------------------------------------
! 格子点数・最大全波数
! Grid points and maximum truncated wavenumber
!-----------------------------------------------------------------
integer:: imax ! 経度格子点数.
! Number of grid points in longitude
integer:: jmax ! 緯度格子点数.
! Number of grid points in latitude
integer:: kmax ! 鉛直層数.
! Number of vertical level
!-----------------------------------------------------------------
! 軸データ
! Axes data
!-----------------------------------------------------------------
real(DP), pointer:: x_Lon (:) =>null()
! 経度. Longitude
real(DP), pointer:: x_Lon_Weight (:) =>null()
! 経度積分用座標重み.
! Weight for integration in longitude
real(DP), pointer:: y_Lat (:) =>null()
! 緯度. Latitude
real(DP), pointer:: y_Lat_Weight (:) =>null()
! 緯度積分用座標重み.
! Weight for integration in latitude
real(DP), pointer:: z_Sigma (:) =>null()
! $ \sigma $ レベル (整数).
! Full $ \sigma $ level
real(DP), pointer:: r_Sigma (:) =>null()
! $ \sigma $ レベル (半整数).
! Half $ \sigma $ level
real(DP), pointer:: z_DelSigma (:) =>null()
! $ \Delta \sigma $ (整数).
! $ \Delta \sigma $ (Full)
!-----------------------------------------------------------------
! 係数
! Coefficients
!-----------------------------------------------------------------
!!$ real(DP):: CoefAlpha ! $ \alpha $ . 係数. Coefficient
!-----------------------------------------------------------------
! 文字データ
! Character data
!-----------------------------------------------------------------
!!$ character(TOKEN):: key00 ! キーワード. Keyword
!!$
!-----------------------------------------------------------------
! 下層のモジュールのオブジェクト
! Objects of low-level modules
!-----------------------------------------------------------------
type(PHYGRD):: phy_grd
! phy_ground モジュールのオブジェクト.
! Object of "phy_ground" module
type(PHYCOOL):: phy_cool
! phy_cooling_s08 モジュールのオブジェクト.
! Object of "phy_cooling_s08" module
type(PHYINTPOL):: phy_intpol
! phy_interpolate モジュールのオブジェクト.
! Object of "phy_interpolate" module
type(PHYVDIF):: phy_vdif
! phy_verdiff モジュールのオブジェクト.
! Object of "phy_verdiff" module
type(PHYIMPL):: phy_impl
! phy_implicit モジュールのオブジェクト.
! Object of "phy_implicit" module
type(PHYSPOLAY):: phy_spo_lay
! phy_sponge_layer モジュールのオブジェクト.
! Object of "phy_sponge_layer" module
type(PHYCUMAD):: phy_cumad
! phy_cumulus_adjust モジュールのオブジェクト.
! Object of "phy_cumulus_adjust" module
type(PHYLSC):: phy_lsc
! phy_lscond モジュールのオブジェクト.
! Object of "phy_lscond" module
type(PHYNEGMST):: phy_neg_mst
! phy_neg_moist モジュールのオブジェクト.
! Object of "phy_neg_moist" module
type(PHYDRYCONV):: phy_dryconv
! phy_dryconv_adjust モジュールのオブジェクト.
! Object of "phy_dryconv_adjust" module
type(PHYADJFIX):: phy_adj_fix
! phy_adjust_fixed モジュールのオブジェクト.
! Object of "phy_adjust_fixed" module
!-----------------------------------------------------------------
! 時刻管理
! Time control
!-----------------------------------------------------------------
type(DC_DIFFTIME):: current_time
! 現在時刻. Current time.
type(DC_DIFFTIME):: delta_time
! $ \Delta t $ . タイムステップ. Time step
!-----------------------------------------------------------------
! ヒストリファイルへのデータ出力設定
! Configure the settings for history data output
!-----------------------------------------------------------------
type(GTHST_NMLINFO):: gthstnml
! NAMELIST#phy_jupiter_case00_history_nml
! から入手される個別のデータ出力情報.
!
! Individual data output information from
! "NAMELIST#phy_jupiter_case00_history_nml".
end type PHYJP
character(*), parameter:: version = &
& '$Name: dcpam4-20080626 $' // &
& '$Id: phy_jupiter_case00.f90,v 1.6 2008-06-09 16:58:24 morikawa Exp $'
!-----------------------------------------------------------------
! 公開手続
! Public procedures
!-----------------------------------------------------------------
interface PhyJpCreate
module procedure PhyJpCreate
end interface
interface PhysicsJupiter
module procedure PhysicsJupiter
end interface
interface PhysicsJupiterAdjust
module procedure PhysicsJupiterAdjust
end interface
!!$ interface PhyJpSample
!!$ module procedure PhyJpSample
!!$ end interface
interface PhyJpClose
module procedure PhyJpClose
end interface
interface PhyJpPutLine
module procedure PhyJpPutLine
end interface
interface PhyJpinitialized
module procedure PhyJpInitialized
end interface
interface PhyJpSetTime
module procedure PhyJpSetTime
end interface
!-----------------------------------------------------------------
! 非公開手続
! Private procedures
!-----------------------------------------------------------------
interface NmlRead
module procedure PhyJpNmlRead
end interface
contains
subroutine PhyJpCreate( phy_jp, &
& imax, jmax, kmax, &
& x_Lon, y_Lat, z_Sigma, r_Sigma, &
& DelTime, &
& GasRDry, Grav, CpDry, &
& LatentHeat, GasRWet, EpsV, ES0, &
& FKarm, &
& x_Lon_Weight, y_Lat_Weight, &
& current_time_value, current_time_unit, &
& history_varlist, &
& history_interval_value, history_interval_unit, &
& history_precision, history_fileprefix, &
& nmlfile, err )
!
! PHYJP 型の変数の初期設定を行います.
! 他のサブルーチンを使用する前に必ずこのサブルーチンによって
! PHYJP 型の変数を初期設定してください.
!
! なお, 与えられた *phy_jp* が既に初期設定されている場合,
! プログラムはエラーを発生させます.
!
! NAMELIST を利用する場合には引数 *nmlfile* に NAMELIST ファイル名
! を与えてください. NAMELIST 変数群の詳細に関しては
! NAMELIST#phy_jupiter_case00_nml を参照してください.
!
! Constructor of "PHYJP".
! Initialize *phy_jp* by this subroutine,
! before other procedures are used,
!
! Note that if *phy_jp* is already initialized
! by this procedure, error is occurred.
!
! In order to use NAMELIST, specify a NAMELIST filename to
! argument *nmlfile*. See "NAMELIST#phy_jupiter_case00_nml"
! for details about a NAMELIST group.
!
use dc_trace, only: BeginSub, EndSub
use dc_string, only: PutLine, Printf, Split, StrInclude, StoA, JoinChar
use dc_types, only: DP, STRING, TOKEN, STDOUT
use dc_present, only: present_and_not_empty, present_and_true, &
& present_select
use dc_message, only: MessageNotify
use dc_error, only: StoreError, DC_NOERR, DC_EALREADYINIT, &
& DC_EARGLACK, DC_ENEGATIVE, DC_ENOFILEREAD, USR_ERRNO, HST_EBADVARNAME
use dc_date, only: DCDiffTimeCreate
use dc_hash, only: HASH, DCHashPut, DCHashRewind, DCHashNext, DCHashDelete
use gt4_history_nmlinfo, only: HstNmlInfoCreate, HstNmlInfoAdd, &
& HstNmlInfoEndDefine, HstNmlInfoPutLine, HstNmlInfoAllNameValid
use gt4_history, only: GT_HISTORY, &
& HistoryAddVariable, HistoryAddAttr
use phy_ground, only: Create, Get
use phy_interpolate, only: Create
use phy_cooling_s08, only: PhyCoolCreate
use phy_implicit, only: Create
use phy_verdiff, only: Create
use phy_sponge_layer, only: PhySpoLayCreate
use phy_cumulus_adjust, only: Create
use phy_lscond, only: Create
use phy_neg_moist, only: Create
use phy_dryconv_adjust, only: Create
use phy_adjust_fixed, only: PhyAdjFixCreate
implicit none
type(PHYJP), intent(inout):: phy_jp
integer, intent(in):: imax ! 経度格子点数.
! Number of grid points in longitude
integer, intent(in):: jmax ! 緯度格子点数.
! Number of grid points in latitude
integer, intent(in):: kmax ! 鉛直層数.
! Number of vertical level
real(DP), intent(in):: x_Lon (0:imax-1)
! 経度. Longitude
real(DP), intent(in):: y_Lat (0:jmax-1)
! 緯度. Latitude
real(DP), intent(in):: z_Sigma (0:kmax-1)
! $ \sigma $ レベル (整数).
! Full $ \sigma $ level
real(DP), intent(in):: r_Sigma (0:kmax)
! $ \sigma $ レベル (半整数).
! Half $ \sigma $ level
real(DP), intent(in):: DelTime ! $ \Delta t $ . タイムステップ. Time step
real(DP), intent(in):: GasRDry
! $ R $ [J kg-1 K-1].
! 乾燥大気の気体定数.
! Gas constant of air
real(DP), intent(in):: Grav
! $ g $ [m s-2].
! 重力加速度.
! Gravitational acceleration
real(DP), intent(in):: CpDry
! $ C_p $ [J kg-1 K-1].
! 乾燥大気の定圧比熱.
! Specific heat of air at constant pressure
real(DP), intent(in):: LatentHeat
! $ L $ [J kg-1] .
! 凝結の潜熱.
! Latent heat of condensation
real(DP), intent(in):: GasRWet
! $ R_v $ [J kg-1 K-1].
! 凝結成分の気体定数.
! Gas constant of condensible elements
real(DP), intent(in):: EpsV
! $ \epsilon_v $ .
! 水蒸気分子量比.
! Molecular weight of water vapor
real(DP), intent(in):: ES0 ! $ e^{*} $ (273K) . 0 ℃での飽和蒸気圧. Saturated vapor pressure at 0 degrees C
real(DP), intent(in):: FKarm
! $ k $ .
! カルマン定数.
! Karman constant
real(DP), intent(in), optional:: x_Lon_Weight (0:imax-1)
! 経度積分用座標重み.
! Weight for integration in longitude
real(DP), intent(in), optional:: y_Lat_Weight (0:jmax-1)
! 緯度積分用座標重み.
! Weight for integration in latitude
real, intent(in), optional:: current_time_value
! 現在時刻の数値. Numerical value of current time
character(*), intent(in), optional:: current_time_unit
! 現在時刻の単位. Unit of current time
character(*), intent(in), optional:: history_varlist
! ヒストリデータの出力変数リスト.
! カンマで区切って並べる.
! (例: "Data1,Data2" ).
!
! List of variables output to history data.
! Delimiter is comma.
! (exp. "Data1,Data2" ).
!
real, intent(in), optional:: history_interval_value
! ヒストリデータの出力間隔の数値.
! Numerical value for interval of history data output
character(*), intent(in), optional:: history_interval_unit
! ヒストリデータの出力間隔の単位.
! Unit for interval of history data output
character(*), intent(in), optional:: history_precision
! ヒストリデータの精度.
! Precision of history data
character(*), intent(in), optional:: history_fileprefix
! ヒストリデータのファイル名の接頭詞.
! Prefix of history data filenames
character(*), intent(in), optional:: nmlfile
! NAMELIST ファイルの名称.
! この引数に空文字以外を与えた場合,
! 指定されたファイルから
! NAMELIST 変数群を読み込みます.
! ファイルを読み込めない場合にはエラーを
! 生じます.
!
! NAMELIST 変数群の詳細に関しては
! NAMELIST#phy_jupiter_case00_nml
! を参照してください.
!
! NAMELIST file name.
! If nonnull character is specified to
! this argument,
! NAMELIST group name is loaded from the
! file.
! If the file can not be read,
! an error occurs.
!
! See "NAMELIST#phy_jupiter_case00_nml"
! for details about a NAMELIST group.
!
logical, intent(out), optional:: err
! 例外処理用フラグ.
! デフォルトでは, この手続き内でエラーが
! 生じた場合, プログラムは強制終了します.
! 引数 *err* が与えられる場合,
! プログラムは強制終了せず, 代わりに
! *err* に .true. が代入されます.
!
! Exception handling flag.
! By default, when error occur in
! this procedure, the program aborts.
! If this *err* argument is given,
! .true. is substituted to *err* and
! the program does not abort.
!-----------------------------------
! 地表面データ
! Data on surface
real(DP):: xy_SurfTemp (0:imax-1, 0:jmax-1)
! 地表面温度.
! Surface temperature
real(DP):: xy_SurfAlbedo (0:imax-1, 0:jmax-1)
! 地表アルベド.
! Surface albedo
real(DP):: xy_SurfHumidCoeff (0:imax-1, 0:jmax-1)
! 地表湿潤度.
! Surface humidity coefficient
real(DP):: xy_SurfRoughLength (0:imax-1, 0:jmax-1)
! 地表粗度長.
! Surface rough length
real(DP):: xy_SurfHeatCapacity (0:imax-1, 0:jmax-1)
! 地表熱容量.
! Surface heat capacity
real(DP):: xy_GroundTempFlux (0:imax-1, 0:jmax-1)
! 地中熱フラックス.
! Ground temperature flux
integer:: xy_SurfCondition (0:imax-1, 0:jmax-1)
! 地表状態.
! Surface condition
!-----------------------------------
! ヒストリファイルへのデータ出力設定
! Configure the settings for history data output
character(STRING):: name = ''
! 変数名. Variable identifier
character(STRING):: longname = ''
! 変数の記述的名称. Descriptive name of variables
character(STRING), allocatable:: dims(:)
! 座標軸の名称. Name of axes
character(STRING):: units = ''
! 単位. Units
type(GT_HISTORY), pointer:: gthist =>null()
! gt4_history モジュール用構造体.
! Derived type for "gt4_history" module
character(TOKEN):: precision
! ヒストリデータの精度.
! Precision of history data
logical:: average
! 出力データの平均化フラグ.
! Flag for average of output data
type(HASH):: registered_varnames
! このモジュールから出力される変数名のリスト.
!
! List of names of variables output
! from this module.
logical:: hash_end
! HASH オブジェクトの終了フラグ.
!
! End flag of "HASH" object
logical:: allvar_invalid
! 無効な変数名のチェックフラグ.
!
! Check flag of invalid variable names
character(STRING):: names_invalid
! 無効な変数名.
!
! Invalid variable names
character(STRING):: nmlfile_msg
! NAMELIST ファイル名に関するメッセージ.
!
! Messages about NAMELIST file name
!-----------------------------------
! 作業変数
! Work variables
integer:: k ! DO ループ用作業変数
! Work variables for DO loop
real(DP), parameter:: PI = 3.1415926535897930_DP
! $ \pi $ . 円周率. Circular constant
integer:: stat
character(STRING):: cause_c
character(*), parameter:: subname = 'PhyJpCreate'
continue
call BeginSub( subname, version )
stat = DC_NOERR
cause_c = ''
!-----------------------------------------------------------------
! 初期設定のチェック
! Check initialization
!-----------------------------------------------------------------
if ( phy_jp % initialized ) then
stat = DC_EALREADYINIT
cause_c = 'PHYJP'
goto 999
end if
!-----------------------------------------------------------------
! 引数の正当性のチェック
! Validate arguments
!-----------------------------------------------------------------
if (imax < 1) then
stat = DC_ENEGATIVE
cause_c = 'imax'
goto 999
end if
if (jmax < 1) then
stat = DC_ENEGATIVE
cause_c = 'jmax'
goto 999
end if
if (kmax < 1) then
stat = DC_ENEGATIVE
cause_c = 'kmax'
goto 999
end if
if (DelTime < 0.0_DP) then
stat = DC_ENEGATIVE
cause_c = 'DelTime'
goto 999
end if
!-----------------------------------------------------------------
! 波数・格子点の設定
! Configure wave number and grid point
!-----------------------------------------------------------------
phy_jp % imax = imax
phy_jp % jmax = jmax
phy_jp % kmax = kmax
!-----------------------------------------------------------------
! 座標軸の設定
! Configure axes
!-----------------------------------------------------------------
allocate( phy_jp % x_Lon (0:imax-1) )
phy_jp % x_Lon = x_Lon
allocate( phy_jp % y_Lat (0:jmax-1) )
phy_jp % y_Lat = y_Lat
allocate( phy_jp % z_Sigma (0:kmax-1) )
phy_jp % z_Sigma = z_Sigma
allocate( phy_jp % r_Sigma (0:kmax) )
phy_jp % r_Sigma = r_Sigma
allocate( phy_jp % x_Lon_Weight (0:imax-1) )
if ( present(x_Lon_Weight) ) then
phy_jp % x_Lon_Weight = x_Lon_Weight
else
phy_jp % x_Lon_Weight = 2.0_DP * PI / imax
end if
allocate( phy_jp % y_Lat_Weight (0:jmax-1) )
if ( present(y_Lat_Weight) ) then
phy_jp % y_Lat_Weight = y_Lat_Weight
else
phy_jp % y_Lat_Weight = 2.0_DP / jmax
end if
allocate( phy_jp % z_DelSigma (0:kmax-1) )
do k = 0, kmax - 1
phy_jp % z_DelSigma(k) = r_Sigma(k) - r_Sigma(k+1)
end do
!-----------------------------------------------------------------
! 時刻管理
! Time control
!-----------------------------------------------------------------
if ( present(current_time_value) .and. present(current_time_unit) ) then
call DCDiffTimeCreate( &
& diff = phy_jp % current_time, & ! (out)
& value = real( current_time_value, DP ), & ! (in)
& unit = current_time_unit ) ! (in)
else
call DCDiffTimeCreate( &
& diff = phy_jp % current_time, & ! (out)
& value = 0.0_DP, & ! (in)
& unit = 'sec' ) ! (in)
end if
call DCDiffTimeCreate( &
& diff = phy_jp % delta_time, & ! (out)
& value = DelTime, & ! (in)
& unit = 'sec' ) ! (in)
!-----------------------------------------------------------------
! "PHYJP" の設定
! Configure the settings for "PHYJP"
!-----------------------------------------------------------------
!!$ phy_jp % CoefAlpha = CoefAlpha
!!$ phy_jp % key00 = ''
!-----------------------------------------------------------------
! 地表条件設定
! Configure ground condition
!-----------------------------------------------------------------
call Create( phy_grd = phy_jp % phy_grd, & ! (inout)
& imax = imax, jmax = jmax, & ! (in)
& SurfTemp = 490.0_DP, & ! (in)
& SurfAlbedo = 0.15_DP, & ! (in)
& SurfHumidCoeff = 1.0_DP, & ! (in)
& SurfRoughLength = 1.0e-4_DP, & ! (in)
& SurfHeatCapacity = 0.0_DP, & ! (in)
& GroundTempFlux = 0.0_DP, & ! (in)
& SurfCondition = 0, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
call Get( phy_grd = phy_jp % phy_grd, & ! (inout)
& xy_SurfTemp = xy_SurfTemp, & ! (out)
& xy_SurfAlbedo = xy_SurfAlbedo, & ! (out)
& xy_SurfHumidCoeff = xy_SurfHumidCoeff, & ! (out)
& xy_SurfRoughLength = xy_SurfRoughLength, & ! (out)
& xy_SurfCondition = xy_SurfCondition, & ! (out)
& xy_SurfHeatCapacity = xy_SurfHeatCapacity, & ! (out)
& xy_GroundTempFlux = xy_GroundTempFlux ) ! (out)
!-----------------------------------------------------------------
! 温度の半整数 $ \sigma $ レベルの補間, 気圧とジオポテンシャルの算出
! Interpolate temperature on half $ \sigma $ level, and
! calculate pressure and geo-potential
!-----------------------------------------------------------------
call Create( phy_intpol = phy_jp % phy_intpol, & ! (inout)
& imax = imax, jmax = jmax, kmax = kmax, & ! (in)
& z_Sigma = z_Sigma, r_Sigma = r_Sigma, & ! (in)
& RAir = GasRDry, Grav = Grav, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 木星を模した大気上層への冷却
! Cooling of atmospheric upper layer to imitate Jupiter
!-----------------------------------------------------------------
call PhyCoolCreate( &
& phy_cool = phy_jp % phy_cool, & ! (inout)
& imax = imax, jmax = jmax, kmax = kmax, & ! (in)
& x_Lon = x_Lon, y_Lat = y_Lat, z_Sigma = z_Sigma, & ! (in)
& DelTime = DelTime, & ! (in)
& day_seconds = 86400.0_DP, & ! (in)
& DTempDdayCooling = -1.0_DP, & ! (in)
& UpperPress = 1.0e4_DP, & ! (in)
& LowerPress = 2.0e5_DP, & ! (in)
& x_Lon_Weight = phy_jp % x_Lon_Weight, & ! (in)
& y_Lat_Weight = phy_jp % y_Lat_Weight, & ! (in)
& z_DelSigma = phy_jp % z_DelSigma, & ! (in)
& current_time_value = current_time_value, & ! (in)
& current_time_unit = current_time_unit, & ! (in)
!!$ & history_varlist = 'DTempDtCooling', & ! (in)
& history_interval_value = history_interval_value, & ! (in)
& history_interval_unit = history_interval_unit, & ! (in)
& history_precision = history_precision, & ! (in)
& history_fileprefix = history_fileprefix, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 陰解配列初期設定
! Initialize implicit matrices
!-----------------------------------------------------------------
call Create( phy_impl = phy_jp % phy_impl, & ! (inout)
& imax = imax, jmax = jmax, kmax = kmax, & ! (in)
& Grav = Grav, Cp = CpDry, EL = LatentHeat, & ! (in)
& DelTime = DelTime, & ! (in)
& xy_SurfHeatCapacity = xy_SurfHeatCapacity, & ! (in)
& xy_SurfCondition = xy_SurfCondition, & ! (in)
& xy_GroundTempFlux = xy_GroundTempFlux, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 鉛直拡散フラックス
! Vertical diffusion flux
!-----------------------------------------------------------------
call Create( phy_vdif = phy_jp % phy_vdif, & ! (inout)
& imax = imax, jmax = jmax, kmax = kmax, & ! (in)
& RAir = GasRDry, Cp = CpDry, Grav = Grav, &! (in)
& EL = LatentHeat, FKarm = FKarm, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 速度を減衰させるためのスポンジ層
! Sponge layer for damping velocity
!-----------------------------------------------------------------
call PhySpoLayCreate( &
& phy_spo_lay = phy_jp % phy_spo_lay, & ! (inout)
& imax = imax, jmax = jmax, kmax = kmax, & ! (in)
& x_Lon = x_Lon, y_Lat = y_Lat, z_Sigma = z_Sigma, & ! (in)
& DelTime = DelTime, & ! (in)
& DampingTime = 86400.0_DP, & ! (in)
& UpperPress = 0.0_DP, & ! (in)
& LowerPress = 1.0e4_DP, & ! (in)
& x_Lon_Weight = phy_jp % x_Lon_Weight, & ! (in)
& y_Lat_Weight = phy_jp % y_Lat_Weight, & ! (in)
& z_DelSigma = phy_jp % z_DelSigma, & ! (in)
& current_time_value = current_time_value, & ! (in)
& current_time_unit = current_time_unit, & ! (in)
!!$ & history_varlist = 'DUDtSpoDamping,DVDtSpoDamping', & ! (in)
& history_interval_value = history_interval_value, & ! (in)
& history_interval_unit = history_interval_unit, & ! (in)
& history_precision = history_precision, & ! (in)
& history_fileprefix = history_fileprefix, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 湿潤過程 (積雲)
! Moist process (cumulus)
!-----------------------------------------------------------------
call Create( phy_cumad = phy_jp % phy_cumad, & ! (inout)
& imax = imax, jmax = jmax, kmax = kmax, & ! (in)
& Grav = Grav, Cp = CpDry, RAir = GasRDry, & ! (in)
& EL = LatentHeat, RVap = GasRWet, & ! (in)
& EpsV = EpsV, ES0 = ES0, & ! (in)
& DelTime = DelTime, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 湿潤過程 (大規模凝結)
! Moist process (large scale condensation)
!-----------------------------------------------------------------
call Create( phy_lsc = phy_jp % phy_lsc, & ! (inout)
& imax = imax, jmax = jmax, kmax = kmax, & ! (in)
& Grav = Grav, Cp = CpDry, & ! (in)
& EL = LatentHeat, RVap = GasRWet, & ! (in)
& EpsV = EpsV, ES0 = ES0, & ! (in)
& DelTime = DelTime, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 負の水蒸気除去
! Remove negative moisture
!-----------------------------------------------------------------
call Create( phy_neg_mst = phy_jp % phy_neg_mst, & ! (inout)
& imax = imax, jmax = jmax, kmax = kmax, & ! (in)
& PI = PI, & ! (in)
& DelTime = DelTime, & ! (in)
& x_Lon_Weight = phy_jp % x_Lon_Weight, & ! (in)
& y_Lat_Weight = phy_jp % y_Lat_Weight, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 乾燥対流調節
! Dry convective adjustment
!-----------------------------------------------------------------
call Create( phy_dryconv = phy_jp % phy_dryconv, & ! (inout)
& imax = imax, jmax = jmax, kmax = kmax, & ! (in)
& RAir = GasRDry, Cp = CpDry, & ! (in)
& DelTime = DelTime, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 最下層の比湿と温度を固定
! Temperature and specific humidity on bottom layer are fixed
!-----------------------------------------------------------------
call PhyAdjFixCreate( &
& phy_adj_fix = phy_jp % phy_adj_fix, & ! (inout)
& imax = imax, jmax = jmax, kmax = 1, & ! (in)
& x_Lon = x_Lon, y_Lat = y_Lat, z_Sigma = z_Sigma, & ! (in)
& DelTime = DelTime, & ! (in)
& TempFixed = 490.0_DP * z_Sigma(0) ** ( GasRDry / CpDry ), & ! (in)
& QVapFixed = 0.00611641_DP, & ! (in)
& x_Lon_Weight = phy_jp % x_Lon_Weight, & ! (in)
& y_Lat_Weight = phy_jp % y_Lat_Weight, & ! (in)
& z_DelSigma = phy_jp % z_DelSigma, & ! (in)
& current_time_value = current_time_value, & ! (in)
& current_time_unit = current_time_unit, & ! (in)
!!$ & history_varlist = 'DTempDtFixed,DQVapDtFixed', & ! (in)
& history_interval_value = history_interval_value, & ! (in)
& history_interval_unit = history_interval_unit, & ! (in)
& history_precision = history_precision, & ! (in)
& history_fileprefix = history_fileprefix, & ! (in)
& nmlfile = nmlfile, & ! (in)
& err = err ) ! (out)
!-----------------------------------------------------------------
! ヒストリファイルへのデータ出力設定
! Configure the settings for history data output
!-----------------------------------------------------------------
!-------------------------
! デフォルト値
! Default values
call HstNmlInfoCreate( gthstnml = phy_jp % gthstnml ) ! (inout)
!-------------------------
! オプショナル引数からの値
! Values from optional arguments
call HstNmlInfoAdd( &
& gthstnml = phy_jp % gthstnml, & ! (inout)
& name = '', & ! (in)
& interval_value = history_interval_value, & ! (in)
& interval_unit = history_interval_unit, & ! (in)
& precision = history_precision, & ! (in)
& average = .false., & ! (in)
& fileprefix = history_fileprefix ) ! (in)
if ( present(history_varlist) ) then
call HstNmlInfoAdd( &
& gthstnml = phy_jp % gthstnml, & ! (inout)
& name = history_varlist ) ! (in)
end if
!-----------------------------------------------------------------
! NAMELIST からの値の読み込み
! Load values from NAMELIST
!-----------------------------------------------------------------
if ( present_and_not_empty(nmlfile) ) then
call MessageNotify( 'M', subname, &
& 'Loading NAMELIST file "%c" ...', &
& c1 = trim(nmlfile) )
call NmlRead ( nmlfile = nmlfile, & ! (in)
!!$ & CoefAlpha = phy_jp % CoefAlpha, & ! (inout)
!!$ & key00_ = phy_jp % key00, & ! (inout)
& gthstnml = phy_jp % gthstnml, & ! (inout)
& err = err ) ! (out)
if ( present_and_true(err) ) then
call MessageNotify( 'W', subname, &
& '"%c" can not be read.', &
& c1 = trim(nmlfile) )
stat = DC_ENOFILEREAD
cause_c = nmlfile
goto 999
end if
end if
call HstNmlInfoEndDefine( &
& gthstnml = phy_jp % gthstnml, & ! (inout)
& err = err ) ! (out)
if ( present_and_true( err ) ) then
stat = USR_ERRNO
goto 999
end if
!-----------------------------------------------------------------
! データ出力の初期設定
! Initialize data output
!-----------------------------------------------------------------
!-------------------------
! xyr_UFlux の出力設定
! Configure the settings for "xyr_UFlux" output
name = 'UFluxJp'
longname = 'surface stress(x)'
units = 'N m-2'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sigm', 'time' )
! 出力ファイルの初期設定.
! * gthist (gt4_history#GT_HISTORY) が設定される.
! Initialize output file.
! * "gthist" (gt4_history#GT_HISTORY) is configured.
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyr_VFlux の出力設定
! Configure the settings for "xyr_VFlux" output
name = 'VFluxJp'
longname = 'surface stress(y)'
units = 'N m-2'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sigm', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyr_TempFlux の出力設定
! Configure the settings for "xyr_TempFlux" output
name = 'TempFluxJp'
longname = 'sensible heat flux'
units = 'W m-2'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sigm', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyr_QVapFlux の出力設定
! Configure the settings for "xyr_QVapFlux" output
name = 'QVapFluxJp'
longname = 'latent heat flux'
units = 'W m-2'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sigm', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyr_RadLFlux の出力設定
! Configure the settings for "xyr_RadLFlux" output
name = 'RadLFluxJp'
longname = 'longwave flux'
units = 'W m-2'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sigm', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyr_RadSFlux の出力設定
! Configure the settings for "xyr_RadSFlux" output
name = 'RadSFlux'
longname = 'shortwave flux'
units = 'W m-2'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sigm', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DTempDtRadL の出力設定
! Configure the settings for "xyz_DTempDtRadL" output
name = 'DTempDtRadL'
longname = 'longwave heating'
units = 'K s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DTempDtRadS の出力設定
! Configure the settings for "xyz_DTempDtRadS" output
name = 'DTempDtRadS'
longname = 'shortwave heating'
units = 'K s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DUDt の出力設定
! Configure the settings for "xyz_DUDt" output
name = 'DUDtJp'
longname = 'diffusive force(x)'
units = 'm s-2'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DVDt の出力設定
! Configure the settings for "xyz_DVDt" output
name = 'DVDtJp'
longname = 'diffusive force(y)'
units = 'm s-2'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DTempDt の出力設定
! Configure the settings for "xyz_DTempDt" output
name = 'DTempDtJp'
longname = 'diffusive heating'
units = 'K s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DQVapDt の出力設定
! Configure the settings for "xyz_DQVapDt" output
name = 'DQVapDtJp'
longname = 'diffusive moistening'
units = 'kg kg-1 s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xy_CumulusRain + xy_LscRain の出力設定
! Configure the settings for "xy_CumulusRain" + "xy_LscRain" output
name = 'Rain'
longname = 'precipitation'
units = 'W m-2'
allocate( dims(3) )
dims = StoA( 'lon', 'lat', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xy_CumulusRain の出力設定
! Configure the settings for "xy_CumulusRain" output
name = 'CumulusRain'
longname = 'precipitation by cumulus scheme'
units = 'W m-2'
allocate( dims(3) )
dims = StoA( 'lon', 'lat', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xy_LscRain の出力設定
! Configure the settings for "xy_LscRain" output
name = 'LscRain'
longname = 'precipitation by large scale condensation'
units = 'W m-2'
allocate( dims(3) )
dims = StoA( 'lon', 'lat', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DNegQVap1Dt + xyz_DNegQVap2Dt の出力設定
! Configure the settings for "xyz_DNegQVap1Dt" + "xyz_DNegQVap2Dt" output
name = 'DNegQVapDt'
longname = 'removed negative moist'
units = 's-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DCumulusTempDt + xyz_DLscTempDt の出力設定
! Configure the settings for "xyz_DCumulusTempDt" + "xyz_DLscTempDt" output
name = 'DCumLscTempDt'
longname = 'cumulus and large-scale condensation heating'
units = 'K s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DCumulusQVapDt + xyz_DLscQVapDt の出力設定
! Configure the settings for "xyz_DCumulusQVapDt" + "xyz_DLscQVapDt" output
name = 'DCumLscQVapDt'
longname = 'cumulus and large-scale condensation moistening'
units = 'kg kg-1 s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DCumulusTempDt の出力設定
! Configure the settings for "xyz_DCumulusTempDt" output
name = 'DCumulusTempDt'
longname = 'cumulus condensation heating'
units = 'K s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DLscTempDt の出力設定
! Configure the settings for "xyz_DLscTempDt" output
name = 'DLscTempDt'
longname = 'large-scale condensation heating'
units = 'K s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DCumulusQVapDt の出力設定
! Configure the settings for "xyz_DCumulusQVapDt" output
name = 'DCumulusQVapDt'
longname = 'cumulus condensation moistening'
units = 'kg kg-1 s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DLscQVapDt の出力設定
! Configure the settings for "xyz_DLscQVapDt" output
name = 'DLscQVapDt'
longname = 'large-scale condensation moistening'
units = 'kg kg-1 s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xyz_DTempDtDryConv の出力設定
! Configure the settings for "xyz_DTempDtDryConv" output
name = 'DTempDtDryConv'
longname = 'temperature tendency by dry convective adjustment'
units = 'K s-1'
allocate( dims(4) )
dims = StoA( 'lon', 'lat', 'sig', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!-------------------------
! xy_DPsDt の出力設定
! Configure the settings for "xy_DPsDt" output
name = 'DPsDt'
longname = 'surface pressure tendency for removal of negative moisture'
units = 'Pa s-1'
allocate( dims(3) )
dims = StoA( 'lon', 'lat', 'time' )
call output_init ! これは内部サブルーチン. This is an internal subroutine
deallocate( dims )
!!$ !-------------------------
!!$ ! x_Data1 の出力設定
!!$ ! Configure the settings for "x_Data1" output
!!$ name = 'Data1'
!!$ longname = 'Sample data (1)'
!!$ units = '1'
!!$ allocate( dims(2) )
!!$ dims = StoA( 'lon', 'time' )
!!$
!!$ ! 出力ファイルの初期設定.
!!$ ! * gthist (gt4_history#GT_HISTORY) が設定される.
!!$ ! Initialize output file.
!!$ ! * "gthist" (gt4_history#GT_HISTORY) is configured.
!!$ call output_init ! これは内部サブルーチン. This is an internal subroutine
!!$
!!$ ! 属性の付加などを行う場合には以下のようにする.
!!$ ! Describe codes as follows in order to add attributes etc.
!!$ if ( associated( gthist ) ) then
!!$ call HistoryAddAttr( &
!!$ & history = gthist, & ! (inout)
!!$ & varname = name, attrname = 'missing_value', & ! (in)
!!$ & value = 9999.0 ) ! (in)
!!$ end if
!!$ deallocate( dims )
!!$
!!$ !-------------------------
!!$ ! y_Data2 の出力設定
!!$ ! Configure the settings for "y_Data2" output
!!$ name = 'Data2'
!!$ longname = 'Sample data (2)'
!!$ units = '1'
!!$ allocate( dims(2) )
!!$ dims = StoA( 'lat', 'time' )
!!$
!!$ ! 出力ファイルの初期設定.
!!$ ! * gthist (gt4_history#GT_HISTORY) が設定される.
!!$ ! Initialize output file.
!!$ ! * "gthist" (gt4_history#GT_HISTORY) is configured.
!!$ call output_init ! これは内部サブルーチン. This is an internal subroutine
!!$ deallocate( dims )
!-----------------------------------------------------------------
! このモジュールから出力される変数名のリストを表示
! Print list of names of variables output from this module
!-----------------------------------------------------------------
call Printf( STDOUT, &
& ' *** MESSAGE *** +---- "%c" module output varnames list -----', &
& c1 = 'phy_jupiter_case00' )
call DCHashRewind( hashv = registered_varnames ) ! (inout)
do
call DCHashNext( hashv = registered_varnames, & ! (inout)
& key = name, value = longname, end = hash_end ) ! (out)
if ( hash_end ) exit
call Printf( STDOUT, &
& ' *** MESSAGE *** | "%c" (%c)', &
& c1 = trim(name), c2 = trim(longname) )
enddo
call DCHashDelete( hashv = registered_varnames ) ! (inout)
call Printf( STDOUT, &
& ' *** MESSAGE *** `----------------------------------------' )
!-----------------------------------------------------------------
! 無効な変数名のチェック
! Check invalid variable names
!-----------------------------------------------------------------
call HstNmlInfoAllNameValid( &
& gthstnml = phy_jp % gthstnml, & ! (inout)
& invalid = allvar_invalid, names = names_invalid, & ! (out)
& err = err ) ! (out)
if ( allvar_invalid ) then
stat = HST_EBADVARNAME
cause_c = names_invalid
phy_jp % initialized = .true.
call PhyJpClose( phy_jp, & ! (inout)
& err = err ) ! (out)
if ( present(nmlfile) ) then
nmlfile_msg = ' or NAMELIST "phy_jupiter_case00_history_nml" in "'// trim(nmlfile) // '"'
else
nmlfile_msg = ''
end if
call MessageNotify( 'W', subname, &
& 'names "%c" from an optional argument "history_varlist"%c are invalid.', &
& c1 = trim(names_invalid), &
& c2 = trim(nmlfile_msg) )
goto 999
end if
!-----------------------------------------------------------------
! 設定値の正当性のチェック
! Validate setting values
!-----------------------------------------------------------------
!!$ if ( phy_jp % CoefAlpha < 0.0_DP ) then
!!$ stat = DC_ENEGATIVE
!!$ cause_c = 'CoefAlpha'
!!$ goto 999
!!$ end if
!-----------------------------------------------------------------
! 終了処理, 例外処理
! Termination and Exception handling
!-----------------------------------------------------------------
phy_jp % initialized = .true.
999 continue
call StoreError( stat, subname, err, cause_c )
call EndSub( subname )
contains
subroutine output_init
!
! 変数 *name* に関して出力ファイルの初期設定を行います.
! 出力ファイル名や出力間隔などの情報は phy_jp % gthstnml
! から取り出されます.
!
! 変数 *name* に関して出力が行われる場合には,
! *gthist* に出力先ファイルの gt4_history#GT_HISTORY
! 型変数を結合させます. そうでない場合は, *gthist* を空状態にします.
!
! また, 出力データの精度を precision に,
! 出力データ平均化の可否を average に設定します.
!
! 標準出力に表示される変数リスト *registered_varnames* に
! *name*, *longname*, *dims*, *units* が登録されます.
!
! An output file is initialized for a variable *name*.
! Information such as the output filename and output intervals
! is taken out of "phy_jp % gthstnml".
!
! When output is done for the variable *name*, *gthist* is
! associated with the "gt4_history#GT_HISTORY" variable of
! the output file. Otherwise, *gthist* is nullified.
!
! Moreover, the accuracy of output data is set to *precision*, and
! right or wrong of averaging the output data is set to *average*.
!
! *name*, *longname*, *dims*, *units* are registered to
! a list of variables *registered_varnames* that is printed to
! standard output.
!
use dc_date, only: DCDiffTimeCreate, EvalSec, EvalByUnit
use gt4_history_nmlinfo, only: GTHST_NMLINFO, &
& HstNmlInfoInitialized, HstNmlInfoInquire, &
& HstNmlInfoOutputValid, HstNmlInfoAssocGtHist, HstNmlInfoPutLine, &
& HstNmlInfoSetValidName
use gt4_history, only: GT_HISTORY, &
& HistoryCreate, HistoryAddVariable, HistoryPut, &
& HistoryAddAttr, HistoryInitialized
!-----------------------------------
! 作業変数
! Work variables
character(STRING):: file
! ヒストリデータのファイル名.
! History data filenames
character(STRING):: dims_str
! 座標軸のリスト.
! List of axes
real:: interval_value
! ヒストリデータの出力間隔の数値.
! Numerical value for interval of history data output
character(TOKEN):: interval_unit
! ヒストリデータの出力間隔の単位.
! Unit for interval of history data output
continue
!-----------------------------------------------------------------
! 標準出力に表示される変数の登録
! Register a variable name for print to standard output
!-----------------------------------------------------------------
if ( allocated(dims) ) then
dims_str = JoinChar( dims, ',' )
else
dims_str = ''
end if
call DCHashPut( hashv = registered_varnames, & ! (inout)
& key = name, & ! (in)
& value = trim( longname ) // ' [' // &
& trim( units ) // '] {' // &
& trim( dims_str ) // '}' ) ! (in)
!-----------------------------------------------------------------
! 変数の初期化
! Initialize variable
!-----------------------------------------------------------------
nullify( gthist )
precision = 'float'
average = .false.
!-----------------------------------------------------------------
! 変数名の有効性を設定
! Set validation of the variable name
!-----------------------------------------------------------------
call HstNmlInfoSetValidName( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name ) ! (in)
!-----------------------------------------------------------------
! 出力が有効かどうかを確認する
! Confirm whether the output is effective
!-----------------------------------------------------------------
if ( .not. HstNmlInfoOutputValid( phy_jp % gthstnml, name ) ) then
return
end if
!-----------------------------------------------------------------
! GT_HISTORY 変数の取得
! Get "GT_HISTORY" variable
!-----------------------------------------------------------------
call HstNmlInfoAssocGtHist( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& history = gthist, & ! (out)
& err = err ) ! (out)
if ( present_and_true( err ) ) return
call HstNmlInfoInquire( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& precision = precision, & ! (out)
& average = average, & ! (out)
& err = err ) ! (out)
if ( present_and_true( err ) ) return
!-----------------------------------------------------------------
! GT_HISTORY 変数の初期設定の確認
! Check initialization of "GT_HISTORY" variable
!-----------------------------------------------------------------
if ( HistoryInitialized( gthist ) ) then
!---------------------------------------------------------------
! HistoryAddVariable による変数作成
! A variable is created by "HistoryAddVariable"
!---------------------------------------------------------------
call HistoryAddVariable( &
& history = gthist, & ! (inout)
& varname = name, dims = dims, & ! (in)
& longname = longname, units = units, & ! (in)
& xtype = precision, average = average ) ! (in)
return
end if
!-----------------------------------------------------------------
! HistoryCreate のための設定値の取得
! Get the settings for "HistoryCreate"
!-----------------------------------------------------------------
call HstNmlInfoInquire( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& file = file, & ! (out)
& interval_unit = interval_unit, & ! (out)
& interval_value = interval_value, & ! (out)
& err = err ) ! (out)
if ( present_and_true( err ) ) return
!-----------------------------------------------------------------
! HistoryCreate によるファイル作成
! Files are created by "HistoryCreate"
!-----------------------------------------------------------------
call HistoryCreate( &
& history = gthist, & ! (out)
& file = file, & ! (in)
& title = 'Virtual planet experiment like Jovian atmosphere', & ! (in)
& source = 'dcpam4 project : ' // trim(version), & ! (in)
& institution = 'GFD Dennou Club', & ! (in)
& dims = StoA( 'lon', 'lat', 'sig', 'sigm', 'time' ), & ! (in)
& dimsizes = (/ phy_jp % imax, phy_jp % jmax, &
& phy_jp % kmax, phy_jp % kmax + 1, 0 /), & ! (in)
& longnames = StoA('longitude', 'latitude', &
& 'sigma at layer midpoints', &
& 'sigma at layer end-points (half level)', &
& 'time'), & ! (in)
& units = StoA( 'degree_east', 'degree_north', &
& '1', '1', interval_unit ), & ! (in)
& origin = real( EvalbyUnit( phy_jp % current_time, &
& interval_unit) ), & ! (in)
& interval = interval_value, & ! (in)
& err = err ) ! (out)
if ( present_and_true( err ) ) then
nullify( gthist )
return
end if
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'lon', attrname = 'standard_name', & ! (in)
& value = 'longitude' ) ! (in)
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'lat', attrname = 'standard_name', & ! (in)
& value = 'latitude' ) ! (in)
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'sig', attrname = 'standard_name', & ! (in)
& value = 'atmosphere_sigma_coordinate' ) ! (in)
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'sigm', attrname = 'standard_name', & ! (in)
& value = 'atmosphere_sigma_coordinate' ) ! (in)
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'time', attrname = 'standard_name', & ! (in)
& value = 'time' ) ! (in)
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'sig', attrname = 'positive', & ! (in)
& value = 'down' ) ! (in)
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'sigm', attrname = 'positive', & ! (in)
& value = 'down' ) ! (in)
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = 'lon', & ! (in)
& array = phy_jp % x_Lon / PI * 180.0_DP ) ! (in)
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = 'lat', & ! (in)
& array = phy_jp % y_Lat / PI * 180.0_DP ) ! (in)
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = 'sig', & ! (in)
& array = phy_jp % z_Sigma ) ! (in)
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = 'sigm', & ! (in)
& array = phy_jp % r_Sigma ) ! (in)
call HistoryAddVariable( &
& history = gthist, & ! (inout)
& varname = 'lon_weight', & ! (in)
& dims = StoA('lon'), & ! (in)
& longname = 'weight for integration in longitude', & ! (in)
& units = 'radian', xtype = 'double' ) ! (in)
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'lon', attrname = 'gt_calc_weight', & ! (in)
& value = 'lon_weight' ) ! (in)
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = 'lon_weight', & ! (in)
& array = phy_jp % x_Lon_Weight ) ! (in)
call HistoryAddVariable( &
& history = gthist, & ! (inout)
& varname = 'lat_weight', & ! (in)
& dims = StoA('lat'), & ! (in)
& longname = 'weight for integration in latitude', & ! (in)
& units = 'radian', xtype = 'double' ) ! (in)
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'lat', attrname = 'gt_calc_weight', & ! (in)
& value = 'lat_weight' ) ! (in)
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = 'lat_weight', & ! (in)
& array = phy_jp % y_Lat_Weight ) ! (in)
call HistoryAddVariable( &
& history = gthist, & ! (inout)
& varname = 'sig_weight', & ! (in)
& dims = StoA('sig'), & ! (in)
& longname = 'weight for integration in sigma', & ! (in)
& units = '1', xtype = 'float' ) ! (in)
call HistoryAddAttr( &
& history = gthist, & ! (inout)
& varname = 'sig', attrname = 'gt_calc_weight', & ! (in)
& value = 'sig_weight' ) ! (in)
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = 'sig_weight', & ! (in)
& array = phy_jp % z_DelSigma ) ! (in)
!-----------------------------------------------------------------
! HistoryAddVariable による変数作成
! A variable is created by "HistoryAddVariable"
!-----------------------------------------------------------------
if ( HistoryInitialized( gthist ) ) then
call HistoryAddVariable( &
& history = gthist, & ! (inout)
& varname = name, dims = dims, & ! (in)
& longname = longname, units = units, & ! (in)
& xtype = precision, average = average ) ! (in)
else
nullify( gthist )
end if
end subroutine output_init
end subroutine PhyJpCreate
subroutine PhysicsJupiter( phy_jp, &
& xyz_U, xyz_V, xyz_Temp, xy_Ps, xyz_QVap, &
& xyz_DUDt, xyz_DVDt, xyz_DTempDt, xyz_DQVapDt, &
& historyput_flag, &
& err )
!
! 放射や鉛直拡散など, いくつかの物理過程
! を計算し, 与えられた風速, 温度, 地表面気圧, 比湿から,
! それぞれの変化量を返します. 詳細は phy_jupiter_case00 を参照してください.
!
! なお, 与えられた *phy_jp* が PhyJpCreate
! によって初期設定されていない場合, プログラムはエラーを発生させます.
!
! Some physical processes (radiation, vertical diffusion, etc)
! are calculated, and tendencies of velocity, temperature,
! surface pressure, specific humidity are returned.
! See "phy_jupiter_case00" for detail.
!
! If *phy_jp* is not initialized
! by "PhyJpCreate" yet, error is occurred.
!
use dc_trace, only: BeginSub, EndSub
use dc_string, only: PutLine, Printf, Split, StrInclude, StoA, JoinChar
use dc_types, only: DP, STRING, TOKEN, STDOUT
use dc_present, only: present_and_true
use dc_date, only: mod, operator(+), operator(==), EvalbyUnit, EvalSec
use dc_error, only: StoreError, DC_NOERR, DC_ENOTINIT, USR_ERRNO
use gt4_history_nmlinfo, only: HstNmlInfoInquire, &
& HstNmlInfoOutputValid, HstNmlInfoAssocGtHist, HstNmlInfoPutLine
use gt4_history, only: GT_HISTORY, HistoryPut, HistoryInitialized
use phy_interpolate, only: InterpolateTemp, InterpolateGeoPot
use phy_cooling_s08, only: Cooling
use phy_implicit, only: GetMatrices, Integrate
use phy_verdiff, only: VerticalDiffusion
use phy_sponge_layer, only: Damping
implicit none
type(PHYJP), intent(inout):: phy_jp
real(DP), intent(in):: xyz_U (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ u $ . 東西風速. Zonal wind
real(DP), intent(in):: xyz_V (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ v $ . 南北風速. Meridional wind
real(DP), intent(in):: xyz_Temp (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ T $ . 温度. Temperature
real(DP), intent(in):: xy_Ps (0:phy_jp%imax-1, 0:phy_jp%jmax-1)
! $ p_s $ . 地表面気圧. Surface pressure
real(DP), intent(in):: xyz_QVap (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ q $ . 比湿. Specific humidity
real(DP), intent(out):: xyz_DUDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ \DP{u}{t} $ .
! 東西風速変化.
! Zonal wind tendency
real(DP), intent(out):: xyz_DVDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ \DP{v}{t} $ .
! 南北風速変化.
! Meridional wind tendency
real(DP), intent(out):: xyz_DTempDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ \DP{T}{t} $ .
! 温度変化.
! Temperature tendency
real(DP), intent(out):: xyz_DQVapDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ \DP{q}{t} $ .
! 比湿変化.
! Temperature tendency
logical, intent(in), optional:: historyput_flag
! データ出力のフラグ.
! SetTime によって時刻を明示的に
! 指定した場合には, この引数に
! .true. または .false. を指定する
! ことでデータ出力のオンオフを
! 明示的に指定する必要があります.
! デフォルトは .false. です.
!
! Data output flag.
! When time is specified by "SetTime",
! explicit specification of data output
! on/off by specifying ".true." or ".false."
! to this argument.
! Default value is ".false.".
!
logical, intent(out), optional:: err
! 例外処理用フラグ.
! デフォルトでは, この手続き内でエラーが
! 生じた場合, プログラムは強制終了します.
! 引数 *err* が与えられる場合,
! プログラムは強制終了せず, 代わりに
! *err* に .true. が代入されます.
!
! Exception handling flag.
! By default, when error occur in
! this procedure, the program aborts.
! If this *err* argument is given,
! .true. is substituted to *err* and
! the program does not abort.
!-----------------------------------
! *phy_jp* から取り出される設定値
! Setting values fetched from *phy_jp*
integer:: kmax ! 鉛直層数.
! Number of vertical level
real(DP):: DelTime ! $ \Delta t $ . タイムステップ. Time step
!-----------------------------------
! 温度 (整数レベル) と地表面気圧から内挿されるデータ
! Data interpolated from temperature (full level) and surface pressure
real(DP):: xyr_Temp (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! $ T $ . 温度 (半整数レベル).
! Temperature (half level)
real(DP):: xyz_Press (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ p $ . 気圧 (整数レベル).
! Air pressure (full level)
real(DP):: xyr_Press (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! $ p $ . 気圧 (半整数レベル).
! Air pressure (half level)
real(DP):: xyz_GeoPot (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ \phi $ . ジオポテンシャル (整数レベル).
! Geo-potential (full level)
real(DP):: xyr_GeoPot (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! $ \phi $ . ジオポテンシャル (半整数レベル).
! Geo-potential (half level)
!-----------------------------------
! 陰解配列
! Implicit matrices
real(DP):: xyza_UVMatrix (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1, -1:1)
! 速度陰解行列.
! Implicit matrix about velocity
real(DP):: xyza_TempMatrix (0:phy_jp%imax-1, 0:phy_jp%jmax-1, -1:phy_jp%kmax-1, -1:1)
! 温度陰解行列.
! Implicit matrix about temperature
real(DP):: xyza_QVapMatrix (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1, -1:1)
! 比湿陰解行列.
! Implicit matrix about specific humidity
!-----------------------------------
! 放射フラックスで計算される量
! Values calculated by radiation flux
real(DP):: xyr_RadLFlux (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! 長波フラックス.
! Longwave flux
real(DP):: xya_SurfRadLMatrix (0:phy_jp%imax-1, 0:phy_jp%jmax-1, -1:1)
! $ T $ 陰解行列: 地表.
! $ T $ implicit matrix: surface
real(DP):: xyra_DelRadLFlux (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax, 0:1)
! 長波地表温度変化.
! Surface temperature tendency with longwave
real(DP):: xyr_RadSFlux (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! 短波 (日射) フラックス.
! Shortwave (insolation) flux
!-----------------------------------
! 鉛直拡散フラックスで計算される量
! Values calculated by vertical diffusion flux
real(DP):: xyr_UFlux (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! 東西風速フラックス.
! Zonal wind flux
real(DP):: xyr_VFlux (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! 南北風速フラックス.
! Meridional wind flux
real(DP):: xyr_TempFlux (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! 温度フラックス.
! Temperature flux
real(DP):: xyr_QVapFlux (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! 比湿フラックス.
! Specific humidity flux
!-----------------------------------
! 地表面フラックスで計算される量
! Values calculated by surface flux
real(DP):: xy_SurfUVMatrix (0:phy_jp%imax-1, 0:phy_jp%jmax-1)
! 速度陰解行列: 地表.
! Implicit matrix about velocity: surface
real(DP):: xyaa_SurfTempMatrix (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:1, -1:1)
! 温度陰解行列: 地表.
! Implicit matrix about temperature: surface
real(DP):: xyaa_SurfQVapMatrix (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:1, -1:1)
! 比湿陰解行列: 地表.
! Implicit matrix about specific humidity: surface
!-----------------------------------
! 温度変換率
! Temperature tendency
real(DP):: xyz_DTempDtRadL (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! 長波加熱率.
! Temperature tendency with longwave
real(DP):: xyz_DTempDtRadS (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! 短波加熱率.
! Temperature tendency with shortwave
real(DP):: xyz_DTempDtWork (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ \DP{T}{t} $ .
! 温度変化. (データ出力用)
! Temperature tendency (for data output)
!-----------------------------------
! 時間変化率の計算 (陰解法)
! Calculate tendency (implicit)
real(DP):: xy_DSurfTempDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1)
! 地表面温度変化率.
! Surface temperature tendency
!-----------------------------------
! 速度を減衰させるためのスポンジ層
! Sponge layer for damping velocity
real(DP):: xyz_DUDtSpoDamping (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ \DP{u}{t} $ .
! 減衰の効果による東西風速変化.
! Zonal wind tendency by damping effect
real(DP):: xyz_DVDtSpoDamping (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ \DP{v}{t} $ .
! 減衰の効果による南北風速変化.
! Meridional wind tendency by damping effect
!-----------------------------------
! ヒストリファイルへのデータ出力設定
! Configure the settings for history data output
character(STRING):: name = ''
! 変数名. Variable identifier
real:: time
! 時刻. Time
type(GT_HISTORY), pointer:: gthist =>null()
! gt4_history モジュール用構造体.
! Derived type for "gt4_history" module
!-----------------------------------
! 作業変数
! Work variables
integer:: k ! DO ループ用作業変数
! Work variables for DO loop
integer:: stat
character(STRING):: cause_c
character(*), parameter:: subname = 'PhysicsJupiter'
continue
call BeginSub( subname )
stat = DC_NOERR
cause_c = ''
!-----------------------------------------------------------------
! 初期設定のチェック
! Check initialization
!-----------------------------------------------------------------
if ( .not. phy_jp % initialized ) then
stat = DC_ENOTINIT
cause_c = 'PHYJP'
goto 999
end if
!-----------------------------------------------------------------
! *phy_jp* に格納されている設定値の取り出し
! Fetch setting values stored in *phy_jp*
!-----------------------------------------------------------------
kmax = phy_jp % kmax
DelTime = EvalSec( phy_jp % delta_time )
!-----------------------------------------------------------------
! 温度の半整数 $ \sigma $ レベルの補間, 気圧とジオポテンシャルの算出 (1)
! Interpolate temperature on half $ \sigma $ level, and
! calculate pressure and geo-potential (1)
!-----------------------------------------------------------------
!-----------------------------------
! 温度の半整数 $ \sigma $ レベルの補間
! Interpolate temperature on half $ \sigma $ level
call InterpolateTemp( phy_intpol = phy_jp % phy_intpol, & ! (inout)
& xyz_Temp = xyz_Temp, & ! (in)
& xyr_Temp = xyr_Temp, & ! (out)
& err = err ) ! (out)
!-----------------------------------
! 気圧とジオポテンシャルの算出
! Calculate pressure and geo-potential
call InterpolateGeoPot( phy_intpol = phy_jp % phy_intpol, & ! (inout)
& xy_Ps = xy_Ps, & ! (in)
& xyz_Temp = xyz_Temp, xyr_Temp = xyr_Temp, & ! (in)
& xyz_Press = xyz_Press, xyr_Press = xyr_Press, & ! (out)
& xyz_GeoPot = xyz_GeoPot, xyr_GeoPot = xyr_GeoPot, & ! (out)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 陰解配列作成
! Create implicit matrices
!-----------------------------------------------------------------
call GetMatrices( phy_impl = phy_jp % phy_impl, & ! (inout)
& xyr_Press = xyr_Press, & ! (in)
& xyr_UFlux = xyr_UFlux, xyr_VFlux = xyr_VFlux, & ! (out)
& xyr_TempFlux = xyr_TempFlux, xyr_QVapFlux = xyr_QVapFlux, & ! (out)
& xyza_UVMatrix = xyza_UVMatrix, & ! (out)
& xyza_TempMatrix = xyza_TempMatrix, & ! (out)
& xyza_QVapMatrix = xyza_QVapMatrix, & ! (out)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 鉛直拡散フラックス
! Vertical diffusion flux
!-----------------------------------------------------------------
! phy_cooling_s08 では, 上空のある範囲のみ冷却を与えるように
! しているので, 短波はゼロで良いが,
! 長波はゼロでは本当はまずいはずなのだが, とりあえず...
! (森川 2008/05/19)
xyr_RadLFlux = 0.0_DP
xya_SurfRadLMatrix = 0.0_DP
xyra_DelRadLFlux = 0.0_DP
xyr_RadSFlux = 0.0_DP
call VerticalDiffusion( phy_vdif = phy_jp % phy_vdif, & ! (inout)
& xyz_U = xyz_U, xyz_V = xyz_V, & ! (in)
& xyz_Temp = xyz_Temp, xyr_Temp = xyr_Temp, & ! (in)
& xyz_QVap = xyz_QVap, & ! (in)
& xyz_Press = xyz_Press, xyr_Press = xyr_Press, & ! (in)
& xyz_GeoPot = xyz_GeoPot, xyr_GeoPot = xyr_GeoPot, & ! (in)
& xyr_UFlux = xyr_UFlux, xyr_VFlux = xyr_VFlux, & ! (inout)
& xyr_TempFlux = xyr_TempFlux, & ! (inout)
& xyr_QVapFlux = xyr_QVapFlux, & ! (inout)
& xyza_UVMatrix = xyza_UVMatrix, & ! (inout)
& xyza_TempMatrix = xyza_TempMatrix, & ! (inout)
& xyza_QVapMatrix = xyza_QVapMatrix, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 地表面フラックスなし
! Surface flux is nothing
!-----------------------------------------------------------------
xy_SurfUVMatrix = 0.0_DP
xyaa_SurfTempMatrix = 0.0_DP
xyaa_SurfQVapMatrix = 0.0_DP
!-----------------------------------------------------------------
! 時間変化率の計算 (陰解法)
! Calculate tendency (implicit)
!-----------------------------------------------------------------
call Integrate( phy_impl = phy_jp % phy_impl, & ! (inout)
& xyr_UFlux = xyr_UFlux, xyr_VFlux = xyr_VFlux, & ! (in)
& xyr_TempFlux = xyr_TempFlux, & ! (in)
& xy_SurfRadSFlux = xyr_RadSFlux(:,:,0), & ! (in)
& xy_SurfRadLFlux = xyr_RadLFlux(:,:,0), & ! (in)
& xyr_QVapFlux = xyr_QVapFlux, & ! (in)
& xyza_UVMatrix = xyza_UVMatrix, & ! (in)
& xyza_TempMatrix = xyza_TempMatrix, & ! (in)
& xyza_QVapMatrix = xyza_QVapMatrix, & ! (in)
& xy_SurfUVMatrix = xy_SurfUVMatrix, & ! (in)
& xyaa_SurfTempMatrix = xyaa_SurfTempMatrix, & ! (in)
& xyaa_SurfQVapMatrix = xyaa_SurfQVapMatrix, & ! (in)
& xya_SurfRadLMatrix = xya_SurfRadLMatrix, & ! (in)
& xyz_DVerdiffUDt = xyz_DUDt, & ! (out)
& xyz_DVerdiffVDt = xyz_DVDt, & ! (out)
& xyz_DVerdiffTempDt = xyz_DTempDt, & ! (out)
& xy_DVerdiffSurfTempDt = xy_DSurfTempDt, & ! (out)
& xyz_DVerdiffQVapDt = xyz_DQVapDt, & ! (out)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 運動量, 温度, 比湿フラックス補正
! Fluxes of momentum, temperature, specifid humidity correction
!-----------------------------------------------------------------
! ※ dcpam3 最新版から導入予定
!!$ call FluxCorrection( phy_impl = phy_jp % phy_impl, & ! (inout)
!!$ & xyz_DVerdiffUDt = xyz_DUDt, & ! (in)
!!$ & xyz_DVerdiffVDt = xyz_DVDt, & ! (in)
!!$ & xyz_DVerdiffTempDt = xyz_DTempDt, & ! (in)
!!$ & xy_DVerdiffSurfTempDt = xy_DSurfTempDt, & ! (in)
!!$ & xyz_DVerdiffQVapDt = xyz_DQVapDt, & ! (in)
!!$ & xyza_UVMatrix = xyza_UVMatrix, & ! (in)
!!$ & xyza_TempMatrix = xyza_TempMatrix, & ! (in)
!!$ & xyza_QVapMatrix = xyza_QVapMatrix, & ! (in)
!!$ & xy_SurfUVMatrix = xy_SurfUVMatrix, & ! (in)
!!$ & xyaa_SurfTempMatrix = xyaa_SurfTempMatrix, & ! (in)
!!$ & xyaa_SurfQVapMatrix = xyaa_SurfQVapMatrix, & ! (in)
!!$ & xyr_UFlux = xyr_UFlux, xyr_VFlux = xyr_VFlux, & ! (inout)
!!$ & xyr_TempFlux = xyr_TempFlux, & ! (inout)
!!$ & xyr_QVapFlux = xyr_QVapFlux ) ! (inout)
!-----------------------------------------------------------------
! 放射フラックス補正
! Radiation flux correction
!-----------------------------------------------------------------
do k = 0, kmax
xyr_RadLFlux(:,:,k) = &
& xyr_RadLFlux(:,:,k) &
& + ( xy_DSurfTempDt * xyra_DelRadLFlux(:,:,k,0) &
& + xyz_DTempDt(:,:,1) * xyra_DelRadLFlux(:,:,k,1) ) &
& * 2.0_DP * DelTime
end do
!-----------------------------------------------------------------
! 木星を模した大気上層への冷却
! Cooling of atmospheric upper layer to imitate Jupiter
!-----------------------------------------------------------------
call Cooling( &
& phy_cool = phy_jp % phy_cool, & ! (inout)
& xy_Ps = xy_Ps, & ! (in)
& xyz_DTempDtCooling = xyz_DTempDtRadL ) ! (out)
xyz_DTempDtRadS = 0.0_DP
!-----------------------------------------------------------------
! 速度を減衰させるためのスポンジ層
! Sponge layer for damping velocity
!-----------------------------------------------------------------
call Damping( &
& phy_spo_lay = phy_jp % phy_spo_lay, & ! (inout)
& xyz_U = xyz_U, xyz_V = xyz_V, xy_Ps = xy_Ps, & ! (in)
& xyz_DUDtSpoDamping = xyz_DUDtSpoDamping, & ! (out)
& xyz_DVDtSpoDamping = xyz_DVDtSpoDamping ) ! (out)
xyz_DUDt = xyz_DUDt + xyz_DUDtSpoDamping
xyz_DVDt = xyz_DVDt + xyz_DVDtSpoDamping
!-----------------------------------------------------------------
! 重力波抵抗
! Gravity wave drag
!-----------------------------------------------------------------
! ※ AGCM5 最新版から導入予定.
!-----------------------------------------------------------------
! 温度変化率へ放射による加熱率を追加
! Radiation heating is added to temperature tendency
!-----------------------------------------------------------------
xyz_DTempDtWork = xyz_DTempDt
xyz_DTempDt = xyz_DTempDt + xyz_DTempDtRadL + xyz_DTempDtRadS
!-----------------------------------------------------------------
! 温度変化分の足し込み
! Add temperature tendency
!-----------------------------------------------------------------
!!$ ! ※ 未導入. (そもそもまともにモジュール化されていない)
!!$ call physics_integrate_surftemp( &
!!$ & xy_SurfTemp, xyr_RadLFlux, xyra_DelRadLFlux, & !(inout)
!!$ & xy_DSurfTempDt, DelTime & !(in)
!!$ & )
!----------------------------------------------------------------
! ヒストリファイルへのデータ出力
! History data output
!----------------------------------------------------------------
!-------------------------
! xyr_UFlux の出力
! Output "xyr_UFlux"
name = 'UFluxJp'
! 出力のチェック.
! * gthist (gt4_history#GT_HISTORY), time (単精度実数型) が設定される.
! Check for output.
! * "gthist" (gt4_history#GT_HISTORY), time (real) are configured.
call output_check ! これは内部サブルーチン. This is an internal subroutine
! 出力データを引数 array に渡す.
! Give output data to argument "array"
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyr_UFlux, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyr_VFlux の出力
! Output "xyr_VFlux"
name = 'VFluxJp'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyr_VFlux, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyr_TempFlux の出力
! Output "xyr_TempFlux"
name = 'TempFluxJp'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyr_TempFlux, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyr_QVapFlux の出力
! Output "xyr_QVapFlux"
name = 'QVapFluxJp'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyr_QVapFlux, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyr_RadLFlux の出力
! Output "xyr_RadLFlux"
name = 'RadLFluxJp'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyr_RadLFlux, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyr_RadSFlux の出力
! Output "xyr_RadSFlux"
name = 'RadSFlux'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyr_RadSFlux, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DTempDtRadL の出力
! Output "xyz_DTempDtRadL"
name = 'DTempDtRadL'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DTempDtRadL, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DTempDtRadS の出力
! Output "xyz_DTempDtRadS"
name = 'DTempDtRadS'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DTempDtRadS, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DUDt の出力
! Output "xyz_DUDt"
name = 'DUDtJp'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DUDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DVDt の出力
! Output "xyz_DVDt"
name = 'DVDtJp'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DVDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DTempDt の出力
! Output "xyz_DTempDt"
name = 'DTempDtJp'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DTempDtWork, & ! (in)
!!$ & varname = name, array = xyz_DTempDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DQVapDt の出力
! Output "xyz_DQVapDt"
name = 'DQVapDtJp'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DQVapDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!!$ !-----------------------------------
!!$ ! x_Data1 の出力
!!$ ! Output "x_Data1"
!!$ name = 'Data1'
!!$
!!$ ! 出力のチェック.
!!$ ! * gthist (gt4_history#GT_HISTORY), time (単精度実数型) が設定される.
!!$ ! Check for output.
!!$ ! * "gthist" (gt4_history#GT_HISTORY), time (real) are configured.
!!$ call output_check ! これは内部サブルーチン. This is an internal subroutine
!!$
!!$ ! 出力データを引数 array に渡す.
!!$ ! Give output data to argument "array"
!!$ if ( associated( gthist ) ) then
!!$ call HistoryPut( &
!!$ & history = gthist, & ! (inout)
!!$ & varname = name, array = x_Data1, & ! (in)
!!$ & time = time, quiet = .false., & ! (in)
!!$ & err = err ) ! (out)
!!$ end if
!!$
!!$ !-----------------------------------
!!$ ! y_Data2 の出力
!!$ ! Output "y_Data2"
!!$ name = 'Data2'
!!$
!!$ ! 出力のチェック.
!!$ ! * gthist (gt4_history#GT_HISTORY), time (単精度実数型) が設定される.
!!$ ! Check for output.
!!$ ! * "gthist" (gt4_history#GT_HISTORY), time (real) are configured.
!!$ call output_check ! これは内部サブルーチン. This is an internal subroutine
!!$
!!$ ! 出力データを引数 array に渡す.
!!$ ! Give output data to argument "array"
!!$ if ( associated( gthist ) ) then
!!$ call HistoryPut( &
!!$ & history = gthist, & ! (inout)
!!$ & varname = name, array = y_Data2, & ! (in)
!!$ & time = time, quiet = .false., & ! (in)
!!$ & err = err ) ! (out)
!!$ end if
!-----------------------------------------------------------------
! 時刻の更新は Adjust にておこなう
! "Adjust" updates time
!-----------------------------------------------------------------
!!$ phy_jp % current_time = &
!!$ & phy_jp % current_time + phy_jp % delta_time
!-----------------------------------------------------------------
! 終了処理, 例外処理
! Termination and Exception handling
!-----------------------------------------------------------------
999 continue
call StoreError( stat, subname, err, cause_c )
call EndSub( subname )
contains
subroutine output_check
!
! 変数 *name* を出力するかどうかをチェックします.
! 出力に関する情報は dcm_sam_cod % gthstnml から取り出されます.
!
! 変数 *name* に関して出力するよう設定されている場合には,
! *gthist* に出力先ファイルの gt4_history#GT_HISTORY
! 型変数を結合させます. そうでない場合は, *gthist* を空状態にします.
!
! また, 現在時刻を *time* に設定します.
!
! Check whether to output variable *name*.
! Information about output is taken out of "dcm_sam_cod % gthstnml".
!
! When output is done for the variable *name*, *gthist* is
! associated with "gt4_history#GT_HISTORY" variable of
! the output file. Otherwise, *gthist* is nullified.
!
! Moreover, current time is set to *time*.
!
character(TOKEN):: interval_unit
! ヒストリデータの出力間隔の単位.
! Unit for interval of history data output
continue
nullify( gthist )
time = 0.0
if ( HstNmlInfoOutputValid( phy_jp % gthstnml, name ) ) then
call HstNmlInfoInquire( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& interval_unit = interval_unit ) ! (out)
time = real( EvalbyUnit( phy_jp % current_time, interval_unit ) )
if ( present_and_true( historyput_flag ) ) time = 0.0
call HstNmlInfoAssocGtHist( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& history = gthist, err = err ) ! (out)
if ( present_and_true( err ) ) then
nullify( gthist )
return
end if
if ( .not. HistoryInitialized( gthist ) ) nullify( gthist )
end if
end subroutine output_check
end subroutine PhysicsJupiter
subroutine PhysicsJupiterAdjust( phy_jp, &
& xyz_Temp, xy_Ps, xyz_QVap, &
& historyput_flag, &
& err )
!
! 積雲パラメタリゼーションや乾燥対流調節など, いくつかの物理過程
! を計算し, 与えられた風速, 温度, 地表面気圧, 比湿を調節して
! 返します. 詳細は phy_jupiter_case00 を参照してください.
!
! なお, 与えられた *phy_jp* が PhyJpCreate
! によって初期設定されていない場合, プログラムはエラーを発生させます.
!
! Some physical processes (cumulus parameterization,
! dry convective adjustment, etc)
! are calculated, and velocity, temperature,
! surface pressure, specific humidity are adjusted and returned.
! See "phy_jupiter_case00" for detail.
!
! If *phy_jp* is not initialized
! by "PhyJpCreate" yet, error is occurred.
!
use dc_trace, only: BeginSub, EndSub
use dc_string, only: PutLine, Printf, Split, StrInclude, StoA, JoinChar
use dc_types, only: DP, STRING, TOKEN, STDOUT
use dc_present, only: present_and_true
use dc_date, only: mod, operator(+), operator(==), EvalbyUnit, EvalSec
use dc_error, only: StoreError, DC_NOERR, DC_ENOTINIT, USR_ERRNO
use gt4_history_nmlinfo, only: HstNmlInfoInquire, &
& HstNmlInfoOutputValid, HstNmlInfoAssocGtHist, HstNmlInfoPutLine
use gt4_history, only: GT_HISTORY, HistoryPut, HistoryInitialized
use phy_interpolate, only: InterpolateTemp, InterpolateGeoPot
use phy_neg_moist, only: RemoveNegQVap
use phy_cumulus_adjust, only: Cumulus
use phy_lscond, only: LScaleCond
use phy_implicit, only: GetMatrices, Integrate
use phy_dryconv_adjust, only: DryConvectAdjust
use phy_adjust_fixed, only: AdjustFixed
implicit none
type(PHYJP), intent(inout):: phy_jp
real(DP), intent(inout):: xyz_Temp (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ T $ . 温度. Temperature
real(DP), intent(inout):: xy_Ps (0:phy_jp%imax-1, 0:phy_jp%jmax-1)
! $ p_s $ . 地表面気圧. Surface pressure
real(DP), intent(inout):: xyz_QVap (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ q $ . 比湿. Specific humidity
logical, intent(in), optional:: historyput_flag
! データ出力のフラグ.
! SetTime によって時刻を明示的に
! 指定した場合には, この引数に
! .true. または .false. を指定する
! ことでデータ出力のオンオフを
! 明示的に指定する必要があります.
! デフォルトは .false. です.
!
! Data output flag.
! When time is specified by "SetTime",
! explicit specification of data output
! on/off by specifying ".true." or ".false."
! to this argument.
! Default value is ".false.".
!
logical, intent(out), optional:: err
! 例外処理用フラグ.
! デフォルトでは, この手続き内でエラーが
! 生じた場合, プログラムは強制終了します.
! 引数 *err* が与えられる場合,
! プログラムは強制終了せず, 代わりに
! *err* に .true. が代入されます.
!
! Exception handling flag.
! By default, when error occur in
! this procedure, the program aborts.
! If this *err* argument is given,
! .true. is substituted to *err* and
! the program does not abort.
!-----------------------------------
! *phy_jp* から取り出される設定値
! Setting values fetched from *phy_jp*
integer:: kmax ! 鉛直層数.
! Number of vertical level
real(DP):: DelTime ! $ \Delta t $ . タイムステップ. Time step
!-----------------------------------
! 負の比湿除去による地表面気圧の変化
! Surface pressure tendency for removal of negative moisture
real(DP):: xy_DPsDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1)
! $ \DP{p_s}{t} $ .
! 負の比湿除去による地表面気圧の変化.
! Surface pressure tendency for removal of negative moisture
!!$ !-----------------------------------
!!$ ! 入力値の保存用変数
!!$ ! Variables for saving input value
!!$ real(DP):: xy_PsB (0:phy_jp%imax-1, 0:phy_jp%jmax-1)
!!$ ! $ p_s $ . 地表面気圧. Surface pressure
!-----------------------------------
! 温度 (整数レベル) と地表面気圧から内挿されるデータ
! Data interpolated from temperature (full level) and surface pressure
real(DP):: xyr_Temp (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! $ T $ . 温度 (半整数レベル).
! Temperature (half level)
real(DP):: xyz_Press (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ p $ . 気圧 (整数レベル).
! Air pressure (full level)
real(DP):: xyr_Press (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! $ p $ . 気圧 (半整数レベル).
! Air pressure (half level)
real(DP):: xyz_GeoPot (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! $ \phi $ . ジオポテンシャル (整数レベル).
! Geo-potential (full level)
real(DP):: xyr_GeoPot (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax)
! $ \phi $ . ジオポテンシャル (半整数レベル).
! Geo-potential (half level)
!-----------------------------------
! 負の水蒸気除去に使用される量
! Values for removal of negative moisture
real(DP):: xyz_DNegQVap1Dt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! 比湿変化率 (1).
! Specific humidity tendency (1)
real(DP):: xyz_DNegQVap2Dt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! 比湿変化率 (2).
! Specific humidity tendency (2)
!-----------------------------------
! 湿潤過程 (積雲) で計算される量
! Values calculated by moist process (cumulus)
real(DP):: xy_CumulusRain (0:phy_jp%imax-1, 0:phy_jp%jmax-1)
! 積雲スキームによる降水量.
! Precipitation by cumulus scheme
real(DP):: xyz_DCumulusTempDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! 積雲スキームによる温度変化率.
! Temperature tendency by cumulus scheme
real(DP):: xyz_DCumulusQVapDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! 積雲スキームによる比湿変化率.
! Specific humidity tendency by cumulus scheme
!-----------------------------------
! 湿潤過程 (大規模凝結) で計算される量
! Values calculated by moist process (large scale condensation)
real(DP):: xyz_DLscTempDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! 大規模凝結による温度変化率.
! Temperature tendency by large scale condensation
real(DP):: xyz_DLscQVapDt (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! 大規模凝結による比湿変化率.
! Specific humidity tendency by large scale condensation
real(DP):: xy_LscRain (0:phy_jp%imax-1, 0:phy_jp%jmax-1)
! 大規模凝結による降水量.
! Precipitation by large scale condensation
!-----------------------------------
! 乾燥対流調節で計算される量
! Values calculated by dry convective adjustment
real(DP):: xyz_DTempDtDryConv (0:phy_jp%imax-1, 0:phy_jp%jmax-1, 0:phy_jp%kmax-1)
! 乾燥対流調節による温度変化率.
! Temperature tendency by dry convective adjustment
!-----------------------------------
! ヒストリファイルへのデータ出力設定
! Configure the settings for history data output
character(STRING):: name = ''
! 変数名. Variable identifier
real:: time
! 時刻. Time
type(GT_HISTORY), pointer:: gthist =>null()
! gt4_history モジュール用構造体.
! Derived type for "gt4_history" module
!-----------------------------------
! 作業変数
! Work variables
integer:: k ! DO ループ用作業変数
! Work variables for DO loop
integer:: stat
character(STRING):: cause_c
character(*), parameter:: subname = 'PhysicsJupiterAdjust'
continue
call BeginSub( subname )
stat = DC_NOERR
cause_c = ''
!-----------------------------------------------------------------
! 初期設定のチェック
! Check initialization
!-----------------------------------------------------------------
if ( .not. phy_jp % initialized ) then
stat = DC_ENOTINIT
cause_c = 'PHYJP'
goto 999
end if
!-----------------------------------------------------------------
! *phy_jp* に格納されている設定値の取り出し
! Fetch setting values stored in *phy_jp*
!-----------------------------------------------------------------
kmax = phy_jp % kmax
DelTime = EvalSec( phy_jp % delta_time )
!!$ !-----------------------------------------------------------------
!!$ ! 入力値の保存
!!$ ! Save input value
!!$ !-----------------------------------------------------------------
!!$ xy_PsB = xy_Ps
!-----------------------------------------------------------------
! 温度の半整数 $ \sigma $ レベルの補間, 気圧とジオポテンシャルの算出 (1)
! Interpolate temperature on half $ \sigma $ level, and
! calculate pressure and geo-potential (1)
!-----------------------------------------------------------------
!-----------------------------------
! 温度の半整数 $ \sigma $ レベルの補間
! Interpolate temperature on half $ \sigma $ level
call InterpolateTemp( phy_intpol = phy_jp % phy_intpol, & ! (inout)
& xyz_Temp = xyz_Temp, & ! (in)
& xyr_Temp = xyr_Temp, & ! (out)
& err = err ) ! (out)
!-----------------------------------
! 気圧とジオポテンシャルの算出
! Calculate pressure and geo-potential
call InterpolateGeoPot( phy_intpol = phy_jp % phy_intpol, & ! (inout)
& xy_Ps = xy_Ps, & ! (in)
& xyz_Temp = xyz_Temp, xyr_Temp = xyr_Temp, & ! (in)
& xyz_Press = xyz_Press, xyr_Press = xyr_Press, & ! (out)
& xyz_GeoPot = xyz_GeoPot, xyr_GeoPot = xyr_GeoPot, & ! (out)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 負の水蒸気除去 (1)
! Remove negative moisture (1)
!-----------------------------------------------------------------
xyz_DNegQVap1Dt = 0.0_DP
call RemoveNegQVap( phy_neg_mst = phy_jp % phy_neg_mst, & ! (inout)
& xyr_Press = xyr_Press, & ! (in)
& xyz_QVap = xyz_QVap, xyz_DNegQVapDt = xyz_DNegQVap1Dt, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 湿潤過程 (積雲)
! Moist process (cumulus)
!-----------------------------------------------------------------
call Cumulus( phy_cumad = phy_jp % phy_cumad, & ! (inout)
& xyz_Temp = xyz_Temp, xyz_QVap = xyz_QVap, & ! (inout)
& xyz_Press = xyz_Press, xyr_Press = xyr_Press, & ! (in)
& xy_Rain = xy_CumulusRain, & ! (out)
& xyz_DTempDt = xyz_DCumulusTempDt, & ! (out)
& xyz_DQVapDt = xyz_DCumulusQVapDt, & ! (out)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 湿潤過程 (大規模凝結)
! Moist process (large scale condensation)
!-----------------------------------------------------------------
call LScaleCond( phy_lsc = phy_jp % phy_lsc, & ! (inout)
& xyz_Temp = xyz_Temp, xyz_QVap = xyz_QVap, & ! (inout)
& xyz_Press = xyz_Press, xyr_Press = xyr_Press, & ! (in)
& xy_Rain = xy_LscRain, & ! (out)
& xyz_DTempDt = xyz_DLscTempDt, & ! (out)
& xyz_DQVapDt = xyz_DLscQVapDt, & ! (out)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 乾燥対流調節
! Dry convective adjustment
!-----------------------------------------------------------------
call DryConvectAdjust( phy_dryconv = phy_jp % phy_dryconv, & ! (inout)
& xyz_Press = xyz_Press, xyr_Press = xyr_Press, & ! (in)
& xyz_Temp = xyz_Temp, & ! (inout)
& xyz_DTempDt = xyz_DTempDtDryConv, & ! (out)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 負の水蒸気除去 (2)
! Remove negative moisture (2)
!-----------------------------------------------------------------
xyz_DNegQVap2Dt = 0.0_DP
call RemoveNegQVap( phy_neg_mst = phy_jp % phy_neg_mst, & ! (inout)
& xyr_Press = xyr_Press, & ! (in)
& xyz_QVap = xyz_QVap, xyz_DNegQVapDt = xyz_DNegQVap2Dt, & ! (inout)
& err = err ) ! (out)
!!$ ! ※ AGCM5 にはなかった, 地表面気圧の再計算
!!$ ! (比湿が減った分だけ気圧を減らす??) をここにいったん退避する.
!!$ ! やまだ由さん作成の dcpam3 に導入された物理過程に入ってた.
!!$
!!$ !----------------------------------- ※ これにいたっては PhysicsJupiter にないといけないが....
!!$ ! 地表面気圧の再計算
!!$ ! Recalculate surface pressure
!!$ do k = 0, kmax - 1
!!$ xy_Ps = xy_Ps &
!!$ & + xyz_DQVapDt(:,:,k) &
!!$ & * ( xyr_Press(:,:,k) - xyr_Press(:,:,k+1) ) &
!!$ & * 2.0_DP * DelTime
!!$ end do
!!$
!!$ !-----------------------------------
!!$ ! 地表面気圧の再計算
!!$ ! Recalculate surface pressure
!!$ do k = 0, kmax - 1
!!$ xy_Ps = xy_Ps &
!!$ & + ( xyz_DLscQVapDt(:,:,k) &
!!$ & + xyz_DCumulusQVapDt(:,:,k) &
!!$ & + xyz_DNegQVap1Dt(:,:,k) ) &
!!$ & * ( xyr_Press(:,:,k) - xyr_Press(:,:,k+1) ) &
!!$ & * 2.0_DP * DelTime
!!$ end do
!!$
!!$ !-----------------------------------
!!$ ! 地表面気圧の再計算
!!$ ! Recalculate surface pressure
!!$ do k = 0, kmax - 1
!!$ xy_Ps = xy_Ps &
!!$ & + xyz_DNegQVap2Dt(:,:,k) &
!!$ & * ( xyr_Press(:,:,k) - xyr_Press(:,:,k+1) ) &
!!$ & * 2.0_DP * DelTime
!!$ end do
xy_DPsDt = 0.0_DP
!-----------------------------------
! 地表面気圧の再計算
! Recalculate surface pressure
do k = 0, kmax - 1
xy_DPsDt = xy_DPsDt &
& + ( xyz_DLscQVapDt(:,:,k) &
& + xyz_DCumulusQVapDt(:,:,k) &
& + xyz_DNegQVap1Dt(:,:,k) ) &
& * ( xyr_Press(:,:,k) - xyr_Press(:,:,k+1) )
end do
do k = 0, kmax - 1
xy_DPsDt = xy_DPsDt &
& + xyz_DNegQVap2Dt(:,:,k) &
& * ( xyr_Press(:,:,k) - xyr_Press(:,:,k+1) )
end do
!-----------------------------------------------------------------
! 最下層の比湿と温度を固定
! Temperature and specific humidity on bottom layer are fixed
!-----------------------------------------------------------------
call AdjustFixed( &
& phy_adj_fix = phy_jp % phy_adj_fix, & ! (inout)
& xyz_Temp = xyz_Temp(:,:,0:0), & ! (inout)
& xyz_QVap = xyz_QVap(:,:,0:0) ) ! (inout)
!----------------------------------------------------------------
! ヒストリファイルへのデータ出力
! History data output
!----------------------------------------------------------------
!-------------------------
! xy_CumulusRain + xy_LscRain の出力
! Output "xy_CumulusRain" + "xy_LscRain"
name = 'Rain'
! 出力のチェック.
! * gthist (gt4_history#GT_HISTORY), time (単精度実数型) が設定される.
! Check for output.
! * "gthist" (gt4_history#GT_HISTORY), time (real) are configured.
call output_check ! これは内部サブルーチン. This is an internal subroutine
! 出力データを引数 array に渡す.
! Give output data to argument "array"
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xy_CumulusRain + xy_LscRain, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xy_CumulusRain の出力
! Output "xy_CumulusRain"
name = 'CumulusRain'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xy_CumulusRain, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xy_LscRain の出力
! Output "xy_LscRain"
name = 'LscRain'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xy_LscRain, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DNegQVap1Dt + xyz_DNegQVap2Dt の出力
! Configure the settings for "xyz_DNegQVap1Dt" + "xyz_DNegQVap2Dt" output
name = 'DNegQVapDt'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DNegQVap1Dt + xyz_DNegQVap2Dt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DCumulusTempDt + xyz_DLscTempDt の出力
! Output "xyz_DCumulusTempDt" + "xyz_DLscTempDt"
name = 'DCumLscTempDt'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DCumulusTempDt + xyz_DLscTempDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DCumulusQVapDt + xyz_DLscQVapDt の出力
! Output "xyz_DCumulusQVapDt" + "xyz_DLscQVapDt"
name = 'DCumLscQVapDt'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DCumulusQVapDt + xyz_DLscQVapDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DCumulusTempDt の出力
! Output "xyz_DCumulusTempDt"
name = 'DCumulusTempDt'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DCumulusTempDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DLscTempDt の出力
! Output "xyz_DLscTempDt"
name = 'DLscTempDt'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DLscTempDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DCumulusQVapDt の出力
! Output "xyz_DCumulusQVapDt"
name = 'DCumulusQVapDt'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DCumulusQVapDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DLscQVapDt の出力
! Output "xyz_DLscQVapDt"
name = 'DLscQVapDt'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DLscQVapDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xyz_DTempDtDryConv の出力
! Output "xyz_DTempDtDryConv"
name = 'DTempDtDryConv'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xyz_DTempDtDryConv, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!-------------------------
! xy_DPsDt の出力
! Output "xy_DPsDt"
name = 'DPsDt'
call output_check ! これは内部サブルーチン. This is an internal subroutine
if ( associated( gthist ) ) then
call HistoryPut( &
& history = gthist, & ! (inout)
& varname = name, array = xy_DPsDt, & ! (in)
& time = time, quiet = .false., & ! (in)
& err = err ) ! (out)
end if
!!$ !-----------------------------------
!!$ ! x_Data1 の出力
!!$ ! Output "x_Data1"
!!$ name = 'Data1'
!!$
!!$ ! 出力のチェック.
!!$ ! * gthist (gt4_history#GT_HISTORY), time (単精度実数型) が設定される.
!!$ ! Check for output.
!!$ ! * "gthist" (gt4_history#GT_HISTORY), time (real) are configured.
!!$ call output_check ! これは内部サブルーチン. This is an internal subroutine
!!$
!!$ ! 出力データを引数 array に渡す.
!!$ ! Give output data to argument "array"
!!$ if ( associated( gthist ) ) then
!!$ call HistoryPut( &
!!$ & history = gthist, & ! (inout)
!!$ & varname = name, array = x_Data1, & ! (in)
!!$ & time = time, quiet = .false., & ! (in)
!!$ & err = err ) ! (out)
!!$ end if
!!$
!!$ !-----------------------------------
!!$ ! y_Data2 の出力
!!$ ! Output "y_Data2"
!!$ name = 'Data2'
!!$
!!$ ! 出力のチェック.
!!$ ! * gthist (gt4_history#GT_HISTORY), time (単精度実数型) が設定される.
!!$ ! Check for output.
!!$ ! * "gthist" (gt4_history#GT_HISTORY), time (real) are configured.
!!$ call output_check ! これは内部サブルーチン. This is an internal subroutine
!!$
!!$ ! 出力データを引数 array に渡す.
!!$ ! Give output data to argument "array"
!!$ if ( associated( gthist ) ) then
!!$ call HistoryPut( &
!!$ & history = gthist, & ! (inout)
!!$ & varname = name, array = y_Data2, & ! (in)
!!$ & time = time, quiet = .false., & ! (in)
!!$ & err = err ) ! (out)
!!$ end if
!-----------------------------------------------------------------
! 時刻の更新
! Update time
!-----------------------------------------------------------------
phy_jp % current_time = &
& phy_jp % current_time + phy_jp % delta_time
!-----------------------------------------------------------------
! 終了処理, 例外処理
! Termination and Exception handling
!-----------------------------------------------------------------
999 continue
call StoreError( stat, subname, err, cause_c )
call EndSub( subname )
contains
subroutine output_check
!
! 変数 *name* を出力するかどうかをチェックします.
! 出力に関する情報は dcm_sam_cod % gthstnml から取り出されます.
!
! 変数 *name* に関して出力するよう設定されている場合には,
! *gthist* に出力先ファイルの gt4_history#GT_HISTORY
! 型変数を結合させます. そうでない場合は, *gthist* を空状態にします.
!
! また, 現在時刻を *time* に設定します.
!
! Check whether to output variable *name*.
! Information about output is taken out of "dcm_sam_cod % gthstnml".
!
! When output is done for the variable *name*, *gthist* is
! associated with "gt4_history#GT_HISTORY" variable of
! the output file. Otherwise, *gthist* is nullified.
!
! Moreover, current time is set to *time*.
!
character(TOKEN):: interval_unit
! ヒストリデータの出力間隔の単位.
! Unit for interval of history data output
continue
nullify( gthist )
time = 0.0
if ( HstNmlInfoOutputValid( phy_jp % gthstnml, name ) ) then
call HstNmlInfoInquire( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& interval_unit = interval_unit ) ! (out)
time = real( EvalbyUnit( phy_jp % current_time, interval_unit ) )
if ( present_and_true( historyput_flag ) ) time = 0.0
call HstNmlInfoAssocGtHist( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& history = gthist, err = err ) ! (out)
if ( present_and_true( err ) ) then
nullify( gthist )
return
end if
if ( .not. HistoryInitialized( gthist ) ) nullify( gthist )
end if
end subroutine output_check
end subroutine PhysicsJupiterAdjust
!!$ subroutine PhyJpSample( phy_jp, err )
!!$ !--
!!$ ! PhyJpSample の要約を記述してください.
!!$ !++
!!$ ! なお, 与えられた *phy_jp* が PhyJpCreate
!!$ ! によって初期設定されていない場合, プログラムはエラーを発生させます.
!!$ !--
!!$ ! Describe brief of "PhyJpSample".
!!$ !++
!!$ ! If *phy_jp* is not initialized
!!$ ! by "PhyJpCreate" yet, error is occurred.
!!$ !
!!$ use dc_trace, only: BeginSub, EndSub
!!$ use dc_string, only: PutLine, Printf, Split, StrInclude, StoA, JoinChar
!!$ use dc_types, only: DP, STRING, TOKEN, STDOUT
!!$ use dc_error, only: StoreError, DC_NOERR, DC_ENOTINIT
!!$ implicit none
!!$ type(PHYJP), intent(inout):: phy_jp
!!$ logical, intent(out), optional:: err
!!$ ! 例外処理用フラグ.
!!$ ! デフォルトでは, この手続き内でエラーが
!!$ ! 生じた場合, プログラムは強制終了します.
!!$ ! 引数 *err* が与えられる場合,
!!$ ! プログラムは強制終了せず, 代わりに
!!$ ! *err* に .true. が代入されます.
!!$ !
!!$ ! Exception handling flag.
!!$ ! By default, when error occur in
!!$ ! this procedure, the program aborts.
!!$ ! If this *err* argument is given,
!!$ ! .true. is substituted to *err* and
!!$ ! the program does not abort.
!!$
!!$ !-----------------------------------
!!$ ! *phy_jp* から取り出される設定値
!!$ ! Setting values fetched from *phy_jp*
!!$!!$ real(DP):: CoefAlpha ! $ \alpha $ . 係数. Coefficient
!!$!!$ real(DP):: DelTime ! $ \Delta t $ . タイムステップ. Time step
!!$
!!$ !-----------------------------------
!!$ ! 作業変数
!!$ ! Work variables
!!$ integer:: stat
!!$ character(STRING):: cause_c
!!$ character(*), parameter:: subname = 'PhyJpSample'
!!$ continue
!!$ call BeginSub( subname )
!!$ stat = DC_NOERR
!!$ cause_c = ''
!!$
!!$ !-----------------------------------------------------------------
!!$ ! 初期設定のチェック
!!$ ! Check initialization
!!$ !-----------------------------------------------------------------
!!$ if ( .not. phy_jp % initialized ) then
!!$ stat = DC_ENOTINIT
!!$ cause_c = 'PHYJP'
!!$ goto 999
!!$ end if
!!$
!!$ !-----------------------------------------------------------------
!!$ ! *phy_jp* に格納されている設定値の取り出し
!!$ ! Fetch setting values stored in *phy_jp*
!!$ !-----------------------------------------------------------------
!!$!!$ CoefAlpha = phy_jp % CoefAlpha
!!$!!$ DelTime = phy_jp % DelTime
!!$
!!$ !-----------------------------------------------------------------
!!$ ! 終了処理, 例外処理
!!$ ! Termination and Exception handling
!!$ !-----------------------------------------------------------------
!!$999 continue
!!$ call StoreError( stat, subname, err, cause_c )
!!$ call EndSub( subname )
!!$
!!$ end subroutine PhyJpSample
subroutine PhyJpClose( phy_jp, err )
!
! PHYJP 型の変数の終了処理を行います.
! なお, 与えられた *phy_jp* が PhyJpCreate
! によって初期設定されていない場合, プログラムはエラーを発生させます.
!
! Deconstructor of "PHYJP".
! Note that if *phy_jp* is not initialized
! by "PhyJpCreate" yet, error is occurred.
!
use dc_trace, only: BeginSub, EndSub
use dc_string, only: PutLine, Printf, Split, StrInclude, StoA, JoinChar
use dc_types, only: DP, STRING, TOKEN, STDOUT
use dc_error, only: StoreError, DC_NOERR, DC_ENOTINIT
use gt4_history_nmlinfo, only: HstNmlInfoClose, HstNmlInfoNames, &
& HstNmlInfoAssocGtHist, HstNmlInfoPutLine
use gt4_history, only: GT_HISTORY, HistoryClose, HistoryInitialized
use phy_ground, only: Close
use phy_interpolate, only: Close
use phy_cooling_s08, only: PhyCoolClose
use phy_implicit, only: Close
use phy_verdiff, only: Close
use phy_sponge_layer, only: PhySpoLayClose
use phy_cumulus_adjust, only: Close
use phy_lscond, only: Close
use phy_neg_moist, only: Close
use phy_dryconv_adjust, only: Close
use phy_adjust_fixed, only: PhyAdjFixClose
implicit none
type(PHYJP), intent(inout):: phy_jp
logical, intent(out), optional:: err
! 例外処理用フラグ.
! デフォルトでは, この手続き内でエラーが
! 生じた場合, プログラムは強制終了します.
! 引数 *err* が与えられる場合,
! プログラムは強制終了せず, 代わりに
! *err* に .true. が代入されます.
!
! Exception handling flag.
! By default, when error occur in
! this procedure, the program aborts.
! If this *err* argument is given,
! .true. is substituted to *err* and
! the program does not abort.
!-----------------------------------
! ヒストリファイルへのデータ出力設定
! Configure the settings for history data output
character(STRING):: name = ''
! 変数名. Variable identifier
character(STRING):: varnames
! 変数名リスト.
! List of variables
character(TOKEN), pointer:: varnames_array(:) =>null()
! 変数名リスト配列.
! List of variables (array)
integer:: i, vnmax
type(GT_HISTORY), pointer:: gthist =>null()
! gt4_history モジュール用構造体.
! Derived type for "gt4_history" module
!-----------------------------------
! 作業変数
! Work variables
integer:: stat
character(STRING):: cause_c
character(*), parameter:: subname = 'PhyJpClose'
continue
call BeginSub( subname )
stat = DC_NOERR
cause_c = ''
!-----------------------------------------------------------------
! 初期設定のチェック
! Check initialization
!-----------------------------------------------------------------
if ( .not. phy_jp % initialized ) then
stat = DC_ENOTINIT
cause_c = 'PHYJP'
goto 999
end if
!-----------------------------------------------------------------
! 軸データの割付解除
! Deallocate axes data
!-----------------------------------------------------------------
deallocate( phy_jp % x_Lon )
deallocate( phy_jp % x_Lon_Weight )
deallocate( phy_jp % y_Lat )
deallocate( phy_jp % y_Lat_Weight )
deallocate( phy_jp % z_Sigma )
deallocate( phy_jp % r_Sigma )
deallocate( phy_jp % z_DelSigma )
!-----------------------------------------------------------------
! 地表条件設定
! Configure ground condition
!-----------------------------------------------------------------
call Close( phy_grd = phy_jp % phy_grd, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 温度の半整数 $ \sigma $ レベルの補間, 気圧とジオポテンシャルの算出
! Interpolate temperature on half $ \sigma $ level, and
! calculate pressure and geo-potential
!-----------------------------------------------------------------
call Close( phy_intpol = phy_jp % phy_intpol, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 木星を模した大気上層への冷却
! Cooling of atmospheric upper layer to imitate Jupiter
!-----------------------------------------------------------------
call PhyCoolClose( &
& phy_cool = phy_jp % phy_cool, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 陰解配列初期設定
! Initialize implicit matrices
!-----------------------------------------------------------------
call Close( phy_impl = phy_jp % phy_impl, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 鉛直拡散フラックス
! Vertical diffusion flux
!-----------------------------------------------------------------
call Close( phy_vdif = phy_jp % phy_vdif, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 速度を減衰させるためのスポンジ層
! Sponge layer for damping velocity
!-----------------------------------------------------------------
call PhySpoLayClose( &
& phy_spo_lay = phy_jp % phy_spo_lay, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 湿潤過程 (積雲)
! Moist process (cumulus)
!-----------------------------------------------------------------
call Close( phy_cumad = phy_jp % phy_cumad, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 湿潤過程 (大規模凝結)
! Moist process (large scale condensation)
!-----------------------------------------------------------------
call Close( phy_lsc = phy_jp % phy_lsc, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 負の水蒸気除去
! Remove negative moisture
!-----------------------------------------------------------------
call Close( phy_neg_mst = phy_jp % phy_neg_mst, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 乾燥対流調節
! Dry convective adjustment
!-----------------------------------------------------------------
call Close( phy_dryconv = phy_jp % phy_dryconv, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 最下層の比湿と温度を固定
! Temperature and specific humidity on bottom layer are fixed
!-----------------------------------------------------------------
call PhyAdjFixClose( &
& phy_adj_fix = phy_jp % phy_adj_fix, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! ヒストリファイルへのデータ出力の終了処理
! Terminate the settings for history data output
!-----------------------------------------------------------------
varnames = HstNmlInfoNames( phy_jp % gthstnml )
call Split( str = varnames, sep = ',', & ! (in)
& carray = varnames_array ) ! (out)
vnmax = size( varnames_array )
do i = 1, vnmax
name = varnames_array(i)
if ( trim( name ) == '' ) exit
nullify( gthist )
call HstNmlInfoAssocGtHist( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& history = gthist, & ! (out)
& err = err ) ! (out)
if ( HistoryInitialized( gthist ) ) then
call HistoryClose( history = gthist, & ! (inout)
& err = err ) ! (out)
end if
end do
!-----------------------------------------------------------------
! ヒストリファイルへのデータ出力設定の割付解除
! Deallocate the settings for history data output
!-----------------------------------------------------------------
call HstNmlInfoClose( &
& phy_jp % gthstnml, & ! (inout)
& err = err ) ! (out)
!-----------------------------------------------------------------
! 終了処理, 例外処理
! Termination and Exception handling
!-----------------------------------------------------------------
phy_jp % initialized = .false.
999 continue
call StoreError( stat, subname, err, cause_c )
call EndSub( subname )
end subroutine PhyJpClose
subroutine PhyJpPutLine( phy_jp, unit, indent, err )
!
! 引数 *phy_jp* に設定されている情報を印字します.
! デフォルトではメッセージは標準出力に出力されます.
! *unit* に装置番号を指定することで, 出力先を変更することが可能です.
!
! Print information of *phy_jp*.
! By default messages are output to standard output.
! Unit number for output can be changed by *unit* argument.
!
use dc_trace, only: BeginSub, EndSub
use dc_string, only: PutLine, Printf, Split, StrInclude, StoA, JoinChar
use dc_types, only: DP, STRING, TOKEN, STDOUT
use dc_date, only: EvalSec
use dc_error, only: StoreError, DC_NOERR, DC_ENOTINIT
use gt4_history_nmlinfo, only: GTHST_NMLINFO, HstNmlInfoPutLine
use phy_ground, only: PutLine
use phy_interpolate, only: PutLine
use phy_cooling_s08, only: PhyCoolPutLine
use phy_implicit, only: PutLine
use phy_verdiff, only: PutLine
use phy_sponge_layer, only: PhySpoLayPutLine
use phy_cumulus_adjust, only: PutLine
use phy_lscond, only: PutLine
use phy_neg_moist, only: PutLine
use phy_dryconv_adjust, only: PutLine
use phy_adjust_fixed, only: PhyAdjFixPutLine
implicit none
type(PHYJP), intent(in):: phy_jp
integer, intent(in), optional:: unit
! 出力先の装置番号.
! デフォルトの出力先は標準出力.
!
! Unit number for output.
! Default value is standard output.
character(*), intent(in), optional:: indent
! 表示されるメッセージの字下げ.
!
! Indent of displayed messages.
logical, intent(out), optional:: err
! 例外処理用フラグ.
! デフォルトでは, この手続き内でエラーが
! 生じた場合, プログラムは強制終了します.
! 引数 *err* が与えられる場合,
! プログラムは強制終了せず, 代わりに
! *err* に .true. が代入されます.
!
! Exception handling flag.
! By default, when error occur in
! this procedure, the program aborts.
! If this *err* argument is given,
! .true. is substituted to *err* and
! the program does not abort.
!-----------------------------------
! 作業変数
! Work variables
integer:: stat
character(STRING):: cause_c
integer:: out_unit
integer:: indent_len
character(STRING):: indent_str
character(*), parameter:: subname = 'PhyJpPutLine'
continue
call BeginSub( subname )
stat = DC_NOERR
cause_c = ''
!-----------------------------------------------------------------
! 出力先装置番号と字下げの設定
! Configure output unit number and indents
!-----------------------------------------------------------------
if ( present(unit) ) then
out_unit = unit
else
out_unit = STDOUT
end if
indent_len = 0
indent_str = ''
if ( present(indent) ) then
if ( len(indent) /= 0 ) then
indent_len = len(indent)
indent_str(1:indent_len) = indent
end if
end if
!-----------------------------------------------------------------
! "PHYJP" の設定の印字
! Print the settings for "PHYJP"
!-----------------------------------------------------------------
if ( phy_jp % initialized ) then
call Printf( unit = out_unit, & ! (in)
& fmt = indent_str(1:indent_len) // &
& '#' ) ! (in)
else
call Printf( unit = out_unit, & ! (in)
& fmt = indent_str(1:indent_len) // & ! (in)
& '#', & ! (in)
& l = (/phy_jp % initialized/) ) ! (in)
end if
!-----------------------------------------------------------------
! 終了処理, 例外処理
! Termination and Exception handling
!-----------------------------------------------------------------
999 continue
call StoreError( stat, subname, err, cause_c )
call EndSub( subname )
end subroutine PhyJpPutLine
logical function PhyJpInitialized( phy_jp ) result(result)
!
! *phy_jp* が初期設定されている場合には .true. が,
! 初期設定されていない場合には .false. が返ります.
!
! If *phy_jp* is initialized, .true. is returned.
! If *phy_jp* is not initialized, .false. is returned.
!
implicit none
type(PHYJP), intent(in):: phy_jp
continue
result = phy_jp % initialized
end function PhyJpInitialized
subroutine PhyJpSetTime( phy_jp, &
& current_time_value, current_time_unit, &
& err )
!
! *phy_jp* に対して時刻の設定を行います.
!
! ヒストリデータを出力している場合には,
! ヒストリデータの出力時刻も設定します.
! 一度でもこのサブルーチンを呼んだ場合には,
! それ以後のヒストリデータ出力前にこのサブルーチンを呼び出し,
! 時刻の設定を行ってください.
! また, データを出力するサブルーチンに対しても
! オプショナル引数 historyput_flag に .true. を与えてください.
!
! なお, 与えられた *phy_jp* が PhyJpCreate
! によって初期設定されていない場合, プログラムはエラーを発生させます.
!
! Set time to *phy_jp*.
!
! When history data are output,
! the output time of history data are specified.
! Once this subroutine is called, the time of history data must be
! specified by this routine before history data output.
! In additional, give ".true." to an optional argument
! "historyput_flag" of a data output subroutine.
!
! If *phy_jp* is not initialized
! by "PhyJpCreate" yet, error is occurred.
!
use dc_trace, only: BeginSub, EndSub
use dc_string, only: PutLine, Printf, Split, StrInclude, StoA, JoinChar
use dc_types, only: DP, STRING, TOKEN, STDOUT
use dc_date, only: DCDiffTimeCreate, EvalbyUnit
use dc_error, only: StoreError, DC_NOERR, DC_ENOTINIT
use gt4_history_nmlinfo, only: HstNmlInfoAdd, HstNmlInfoInquire, &
& HstNmlInfoNames, HstNmlInfoAssocGtHist, &
& HstNmlInfoOutputStepDisable, HstNmlInfoPutLine
use gt4_history, only: GT_HISTORY, HistorySetTime, HistoryInitialized
implicit none
type(PHYJP), intent(inout):: phy_jp
real, intent(in):: current_time_value
! 現在時刻の数値. Numerical value of current time
character(*), intent(in):: current_time_unit
! 現在時刻の単位. Unit of current time
logical, intent(out), optional:: err
! 例外処理用フラグ.
! デフォルトでは, この手続き内でエラーが
! 生じた場合, プログラムは強制終了します.
! 引数 *err* が与えられる場合,
! プログラムは強制終了せず, 代わりに
! *err* に .true. が代入されます.
!
! Exception handling flag.
! By default, when error occur in
! this procedure, the program aborts.
! If this *err* argument is given,
! .true. is substituted to *err* and
! the program does not abort.
!-----------------------------------
! ヒストリファイルへのデータ出力設定
! Configure the settings for history data output
character(STRING):: name = ''
! 変数名. Variable identifier
character(TOKEN):: interval_unit
! ヒストリデータの出力間隔の単位.
! Unit for interval of history data output
character(STRING):: varnames
! 変数名リスト.
! List of variables
character(TOKEN), pointer:: varnames_array(:) =>null()
! 変数名リスト配列.
! List of variables (array)
integer:: i, vnmax
type(GT_HISTORY), pointer:: gthist =>null()
! gt4_history モジュール用構造体.
! Derived type for "gt4_history" module
!-----------------------------------
! 作業変数
! Work variables
integer:: stat
character(STRING):: cause_c
character(*), parameter:: subname = 'PhyJpSetTime'
continue
call BeginSub( subname )
stat = DC_NOERR
cause_c = ''
!-----------------------------------------------------------------
! 初期設定のチェック
! Check initialization
!-----------------------------------------------------------------
if ( .not. phy_jp % initialized ) then
stat = DC_ENOTINIT
cause_c = 'PHYJP'
goto 999
end if
!-----------------------------------------------------------------
! 時刻設定
! Configure time
!-----------------------------------------------------------------
call DCDiffTimeCreate( &
& diff = phy_jp % current_time, & ! (out)
& value = real( current_time_value, DP ), & ! (in)
& unit = current_time_unit ) ! (in)
!-----------------------------------------------------------------
! ヒストリファイルへのデータの時刻設定
! Configure the time of history data
!-----------------------------------------------------------------
varnames = HstNmlInfoNames( phy_jp % gthstnml )
call Split( str = varnames, sep = ',', & ! (in)
& carray = varnames_array ) ! (out)
vnmax = size( varnames_array )
do i = 1, vnmax
name = varnames_array(i)
if ( trim( name ) == '' ) exit
call HstNmlInfoOutputStepDisable( &
& gthstnml = phy_jp % gthstnml, & ! (inout)
& name = name, & ! (in)
& err = err ) ! (out)
nullify( gthist )
call HstNmlInfoAssocGtHist( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& history = gthist, & ! (out)
& err = err ) ! (out)
if ( HistoryInitialized( gthist ) ) then
call HstNmlInfoInquire( &
& gthstnml = phy_jp % gthstnml, & ! (in)
& name = name, & ! (in)
& interval_unit = interval_unit ) ! (out)
call HistorySetTime( &
& history = gthist, & ! (inout)
& time = &
& real( EvalbyUnit( phy_jp % current_time, &
& interval_unit) ) ) ! (in)
end if
end do
!-----------------------------------------------------------------
! 終了処理, 例外処理
! Termination and Exception handling
!-----------------------------------------------------------------
999 continue
call StoreError( stat, subname, err, cause_c )
call EndSub( subname )
end subroutine PhyJpSetTime
subroutine PhyJpNmlRead( nmlfile, &
!!$ & CoefAlpha, &
!!$ & key00_, &
& gthstnml, &
& err )
!
! NAMELIST ファイル *nmlfile* から値を入力するための
! サブルーチンです. PhyJpCreate
! 内で呼び出されることを想定しています.
!
! 値が NAMELIST ファイル内で指定されていない場合には,
! 入力された値がそのまま返ります.
!
! なお, *nmlfile* に空文字が与えられた場合, または
! 与えられた *nmlfile* を読み込むことができない場合,
! プログラムはエラーを発生させます.
!
! This is a subroutine to input values from
! NAMELIST file *nmlfile*. This subroutine is expected to be
! called by "PhyJpCreate".
!
! A value not specified in NAMELIST file is returned
! without change.
!
! If *nmlfile* is empty, or *nmlfile* can not be read,
! error is occurred.
!
use dc_trace, only: BeginSub, EndSub
use dc_string, only: PutLine, Printf, Split, StrInclude, StoA, JoinChar
use dc_types, only: DP, STRING, TOKEN, STDOUT
use dc_iounit, only: FileOpen
use dc_message, only: MessageNotify
use dc_present, only: present_and_true
use dc_date, only: DCDiffTimeCreate
use dc_error, only: StoreError, DC_NOERR, DC_ENOFILEREAD, DC_ENOTINIT
use gt4_history_nmlinfo, only: GTHST_NMLINFO, HstNmlInfoAdd, &
& HstNmlInfoInquire, HstNmlInfoInitialized, HstNmlInfoPutLine
implicit none
character(*), intent(in):: nmlfile
! NAMELIST ファイルの名称.
! NAMELIST file name
!!$ real(DP), intent(inout):: CoefAlpha
!!$ ! $ \alpha $ . 係数. Coefficient
!!$
!!$ character(*), intent(inout):: key00_
!!$ character(TOKEN):: key00
!!$ ! キーワード. Keyword
!!$
type(GTHST_NMLINFO), intent(inout):: gthstnml
! NAMELIST#phy_jupiter_case00_history_nml
! から入手される個別のデータ出力情報.
!
! 初期設定やデフォルト値の設定などを
! 行った後に与えること.
!
! Individual data output information from
! "NAMELIST#phy_jupiter_case00_history_nml".
!
! Before this argument is given to
! this procedure, initialize and
! configure the defaut settings.
logical, intent(out), optional:: err
! 例外処理用フラグ.
! デフォルトでは, この手続き内でエラーが
! 生じた場合, プログラムは強制終了します.
! 引数 *err* が与えられる場合,
! プログラムは強制終了せず, 代わりに
! *err* に .true. が代入されます.
!
! Exception handling flag.
! By default, when error occur in
! this procedure, the program aborts.
! If this *err* argument is given,
! .true. is substituted to *err* and
! the program does not abort.
!!$ namelist /phy_jupiter_case00_nml/ &
!!$ & CoefAlpha, key00
!!$ ! phy_jupiter_case00 モジュール用
!!$ ! NAMELIST 変数群名.
!!$ !
!!$ ! phy_jupiter_case00#PhyJpCreate
!!$ ! を使用する際に, オプショナル引数 *nmlfile*
!!$ ! へ NAMELIST ファイル名を指定することで,
!!$ ! そのファイルからこの NAMELIST 変数群を
!!$ ! 読み込みます.
!!$ !
!!$ ! NAMELIST group name for
!!$ ! "phy_jupiter_case00" module.
!!$ !
!!$ ! If a NAMELIST filename is specified to
!!$ ! an optional argument *nmlfile* when
!!$ ! "phy_jupiter_case00#PhyJpCreate"
!!$ ! is used, this NAMELIST group is
!!$ ! loaded from the file.
character(STRING):: name
! 変数名.
! 空白の場合には, この他の設定値は
! phy_jupiter_case00 モジュールにおいて
! 出力されるデータ全ての
! デフォルト値となります.
!
! "Data1,Data2" のようにカンマで区切って複数
! の変数を指定することも可能です.
!
! Variable identifier.
! If blank is given, other values are
! used as default values of output data
! in "phy_jupiter_case00".
!
! Multiple variables can be specified
! as "Data1,Data2" too. Delimiter is comma.
character(STRING):: file
! 出力ファイル名.
! これはデフォルト値としては使用されません.
! *name* に値が設定されている時のみ有効です.
!
! Output file name.
! This is not used as default value.
! This value is valid only when *name* is
! specified.
real:: interval_value
! ヒストリデータの出力間隔の数値.
! 負の値を与えると, 出力を抑止します.
! Numerical value for interval of history data output
! Negative values suppresses output.
character(TOKEN):: interval_unit
! ヒストリデータの出力間隔の単位.
! Unit for interval of history data output
character(TOKEN):: precision
! ヒストリデータの精度.
! Precision of history data
logical:: average
! 出力データの平均化フラグ.
! Flag for average of output data
character(STRING):: fileprefix
! ヒストリデータのファイル名の接頭詞.
! Prefixes of history data filenames
namelist /phy_jupiter_case00_history_nml/ &
& name, &
& file, &
& interval_value, &
& interval_unit, &
& precision, &
& fileprefix, &
& average
! phy_jupiter_case00 モジュールのヒストリデータ用
! NAMELIST 変数群名.
!
! phy_jupiter_case00#PhyJpCreate
! を使用する際に, オプショナル引数 *nmlfile*
! へ NAMELIST ファイル名を指定することで,
! そのファイルからこの NAMELIST 変数群を
! 読み込みます.
!
! NAMELIST group name for
! history data of "phy_jupiter_case00" module.
!
! If a NAMELIST filename is specified to
! an optional argument *nmlfile* when
! "phy_jupiter_case00#PhyJpCreate"
! is used, this NAMELIST group is
! loaded from the file.
!-----------------------------------
! 作業変数
! Work variables
integer:: stat
character(STRING):: cause_c
integer:: unit_nml ! NAMELIST ファイルオープン用装置番号.
! Unit number for NAMELIST file open
integer:: iostat_nml ! NAMELIST 読み込み時の IOSTAT.
! IOSTAT of NAMELIST read
character(TOKEN):: pos_nml
! NAMELIST 読み込み時のファイル位置.
! File position of NAMELIST read
character(*), parameter:: subname = 'PhyJpNmlRead'
continue
call BeginSub( subname )
stat = DC_NOERR
cause_c = ''
!-----------------------------------------------------------------
! 初期設定のチェック
! Check initialization
!-----------------------------------------------------------------
!-----------------------------------------------------------------
! 文字型引数を NAMELIST 変数群へ代入
! Substitute character arguments to NAMELIST group
!-----------------------------------------------------------------
!!$ key00 = key00_
!----------------------------------------------------------------
! NAMELIST ファイルのオープン
! Open NAMELIST file
!----------------------------------------------------------------
call FileOpen( unit = unit_nml, & ! (out)
& file = nmlfile, mode = 'r', & ! (in)
& err = err ) ! (out)
if ( present_and_true(err) ) then
stat = DC_ENOFILEREAD
cause_c = nmlfile
goto 999
end if
!-----------------------------------------------------------------
! NAMELIST 変数群の取得
! Get NAMELIST group
!-----------------------------------------------------------------
!-------------------------
! 係数などの取得
! Get coefficients etc.
!!$ rewind( unit_nml )
!!$ read( unit = unit_nml, & ! (in)
!!$ & nml = phy_jupiter_case00_nml, & ! (out)
!!$ & iostat = iostat_nml ) ! (out)
!!$ if ( iostat_nml == 0 ) then
!!$ call MessageNotify( 'M', subname, &
!!$ & 'NAMELIST group "%c" is loaded from "%c".', &
!!$ & c1 = 'phy_jupiter_case00_nml', c2 = trim(nmlfile) )
!!$ write(STDOUT, nml = phy_jupiter_case00_nml)
!!$ else
!!$ call MessageNotify( 'W', subname, &
!!$ & 'NAMELIST group "%c" is not found in "%c" (iostat=%d).', &
!!$ & c1 = 'phy_jupiter_case00_nml', c2 = trim(nmlfile), &
!!$ & i = (/iostat_nml/) )
!!$ end if
!!$
!-------------------------
! 出力データの個別情報の取得
! Get individual information of output data
rewind( unit_nml )
iostat_nml = 0
pos_nml = ''
do while ( trim(pos_nml) /= 'APPEND' .and. iostat_nml == 0 )
name = ''
file = ''
call HstNmlInfoInquire( &
& gthstnml = gthstnml, & ! (in)
& interval_value = interval_value, & ! (out)
& interval_unit = interval_unit, & ! (out)
& precision = precision, & ! (out)
& average = average, & ! (out)
& fileprefix = fileprefix ) ! (out)
read( unit = unit_nml, & ! (in)
& nml = phy_jupiter_case00_history_nml, & ! (out)
& iostat = iostat_nml ) ! (out)
inquire( unit = unit_nml, & ! (in)
& position = pos_nml ) ! (out)
if ( iostat_nml == 0 ) then
call MessageNotify( 'M', subname, &
& 'NAMELIST group "%c" is loaded from "%c".', &
& c1='phy_jupiter_case00_history_nml', c2=trim(nmlfile) )
write(STDOUT, nml = phy_jupiter_case00_history_nml)
call HstNmlInfoAdd( &
& gthstnml = gthstnml, & ! (in)
& name = name, & ! (in)
& file = file, & ! (in)
& interval_value = interval_value, & ! (in)
& interval_unit = interval_unit, & ! (in)
& precision = precision, & ! (in)
& average = average, & ! (in)
& fileprefix = fileprefix ) ! (in)
else
call MessageNotify( 'W', subname, &
& 'NAMELIST group "%c" is not found in "%c" any more (iostat=%d).', &
& c1='phy_jupiter_case00_history_nml', c2=trim(nmlfile), &
& i = (/iostat_nml/) )
end if
end do
close( unit_nml )
!-----------------------------------------------------------------
! NAMELIST 変数群を文字型引数へ代入
! Substitute NAMELIST group to character arguments
!-----------------------------------------------------------------
!!$ key00_ = key00
!-----------------------------------------------------------------
! 終了処理, 例外処理
! Termination and Exception handling
!-----------------------------------------------------------------
999 continue
call StoreError( stat, subname, err, cause_c )
call EndSub( subname )
end subroutine PhyJpNmlRead
end module phy_jupiter_case00