Class const_provider
In: shared/const_provider.f90

物理定数のサンプル値提供モジュール

Sample values of physical constants provider

Note that Japanese and English are described in parallel.

物理定数のサンプル値を提供するためのモジュールです. 現在は, 円周率や普遍気体定数などのいわゆる普遍定数と, 地球大気, 木星大気における平均分子量などの値を提供します.

このモジュールはあくまでもメインプログラムに対して いくつかの物理定数を提供するのみであり, モデル全体の物理定数 管理を行うものでないことに注意してください.

This module provides sample values of physical constants. Currently, some universal constants (for example pi and the universal gas constant) are provided. And mean molecular weight etc. in the earth atmosphere and the Jupiter atmosphere are provided.

Note that this module not only manage physical constants of entire models but also only provides some physical constants to main programs.

Procedures List

ConstGet :物理定数のサンプル値を取得します
———— :————
ConstGet :Get sample values of physical constants

Usage

ConstGet を用いて物理定数のサンプル値を取得してください.

Get sample values of physical constants by "ConstGet".

Example

  program const_provider_sample
    use dc_string, only: Printf
    use dc_types, only: DP
    use const_provider, only: ConstGet
    implicit none
    real(DP):: RPlanet, Omega, Grav

    ! print some phycical constants on Earth
    call ConstGet( planet = 'earth', &                 ! (in)
      & RPlanet = RPlanet, Omega = Omega, Grav = Grav) ! (out)
    call Printf( &
      & fmt = 'The radius of Earth is %f (m), ' // &
      &       'the angular velocity of Earth is %f (m/s), ' // &
      &       'the acceleration due to gravity on Earth is %f (m^2/s)', &
      & d=(/RPlanet, Omega, Grav/))

    ! print some phycical constants on Mars
    call ConstGet( planet = 'mars'                     ! (in)
      & RPlanet = RPlanet, Omega = Omega, Grav = Grav) ! (out)
    call Printf( &
      & fmt = 'The radius of Mars is %f (m), ' // &
      &       'the angular velocity of Mars is %f (m/s), ' // &
      &       'the acceleration due to gravity on Mars is %f (m^2/s)', &
      & d=(/RPlanet, Omega, Grav/))

  end program const_provider_sample

Methods

ConstGet   version  

Included Modules

dc_types dc_trace dc_error dc_string dcpam_error

Public Instance methods

Subroutine :
planet :character(*), intent(in), optional
PI :real(DP), intent(out), optional
: $ pi $ . 円周率. Circular constant
GasRUniv :real(DP), intent(out), optional
: $ R^{*} $ [J K-1 mol-1]. 普遍気体定数. Universal gas constant
StB :real(DP), intent(out), optional
: $ sigma_{SB} $ . ステファンボルツマン定数. Stefan-Boltzmann constant
FKarm :real(DP), intent(out), optional
: $ k $ . カルマン定数. Karman constant
RPlanet :real(DP), intent(out), optional
: $ a $ [m]. 惑星半径. Radius of planet
Omega :real(DP), intent(out), optional
: $ Omega $ [s-1]. 回転角速度. Angular velocity
Grav :real(DP), intent(out), optional
: $ g $ [m s-2]. 重力加速度. Gravitational acceleration
MolWtDry :real(DP), intent(out), optional
: $ ?? $ [kg mol-1]. 乾燥大気の平均分子量. Mean molecular weight of dry air
CpDry :real(DP), intent(out), optional
: $ C_p $ [J kg-1 K-1]. 乾燥大気の定圧比熱. Specific heat of air at constant pressure
GasRDry :real(DP), intent(out), optional
: $ R $ [J kg-1 K-1]. 乾燥大気の気体定数. Gas constant of air
MolWtWet :real(DP), intent(out), optional
: $ ?? $ [kg mol-1]. 凝結成分の平均分子量. Mean molecular weight of condensible elements
CpWet :real(DP), intent(out), optional
: $ C_v $ [J kg-1 K-1] . 凝結成分の定圧比熱. Specific heat of condensible elements at constant pressure
GasRWet :real(DP), intent(out), optional
: $ R_v $ [J kg-1 K-1]. 凝結成分の気体定数. Gas constant of condensible elements
LatentHeat :real(DP), intent(out), optional
: $ L $ [J kg-1] . 凝結の潜熱. Latent heat of condensation
EpsV :real(DP), intent(out), optional
: $ epsilon_v $ . 水蒸気分子量比. Molecular weight of water vapor
err :logical, intent(out), optional
: 例外処理用フラグ. デフォルトでは, この手続き内でエラーが 生じた場合, プログラムは強制終了します. 引数 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.

物理定数のサンプル値を返します.

planet に以下の文字列を与えることで, それぞれの惑星 に関する標準的な値が返ります. デフォルトは地球の値です.

  • earth
    • 地球大気 (凝結成分は水のみ)
  • jupiter00
    • 木星大気 (凝結成分は水のみ)

Sample values of physical constants are returned.

If following strings is given to planet, standard values of each planet are returned. By default, values of earth are returned.

  • earth
    • earth atmosphere (condensible element is only wator vapor)
  • jupiter00
    • Jovian atmosphere (condensible element is only wator vapor)

[Source]

  subroutine ConstGet( planet, PI, GasRUniv, StB, FKarm, RPlanet, Omega, Grav, MolWtDry, CpDry, GasRDry, MolWtWet, CpWet, GasRWet, LatentHeat, EpsV, err )
    !
    ! 物理定数のサンプル値を返します. 
    ! 
    ! *planet* に以下の文字列を与えることで, それぞれの惑星
    ! に関する標準的な値が返ります. デフォルトは地球の値です. 
    !
    ! * earth
    !   * 地球大気 (凝結成分は水のみ)
    ! * jupiter00
    !   * 木星大気 (凝結成分は水のみ)
    !
    ! Sample values of physical constants are returned. 
    ! 
    ! If following strings is given to *planet*, 
    ! standard values of each planet are returned. 
    ! By default, values of earth are returned. 
    !
    ! * earth
    !   * earth atmosphere (condensible element is only wator vapor)
    ! * jupiter00
    !   * Jovian atmosphere (condensible element is only wator vapor)
    !
    use dc_trace, only: BeginSub, EndSub
    use dc_types, only: STRING, STDOUT
    use dc_error, only: DC_NOERR
    use dc_string, only: LChar
    use dcpam_error, only: StoreError, DCPAM_ENOPLANET
    implicit none
    character(*), intent(in), optional:: planet
    real(DP), intent(out), optional:: PI
                              ! $ \pi $ .
                              ! 円周率.  Circular constant
    real(DP), intent(out), optional:: GasRUniv
                              ! $ R^{*} $ [J K-1 mol-1].
                              ! 普遍気体定数.  Universal gas constant
    real(DP), intent(out), optional:: StB
                              ! $ \sigma_{SB} $ . 
                              ! ステファンボルツマン定数. 
                              ! Stefan-Boltzmann constant
    real(DP), intent(out), optional:: FKarm
                              ! $ k $ .
                              ! カルマン定数. 
                              ! Karman constant

    real(DP), intent(out), optional:: RPlanet
                              ! $ a $ [m]. 
                              ! 惑星半径. 
                              ! Radius of planet
    real(DP), intent(out), optional:: Omega
                              ! $ \Omega $ [s-1]. 
                              ! 回転角速度. 
                              ! Angular velocity
    real(DP), intent(out), optional:: Grav
                              ! $ g $ [m s-2]. 
                              ! 重力加速度. 
                              ! Gravitational acceleration
    real(DP), intent(out), optional:: MolWtDry
                              ! $ ?? $ [kg mol-1]. 
                              ! 乾燥大気の平均分子量. 
                              ! Mean molecular weight of dry air
    real(DP), intent(out), optional:: CpDry
                              ! $ C_p $ [J kg-1 K-1]. 
                              ! 乾燥大気の定圧比熱. 
                              ! Specific heat of air at constant pressure
    real(DP), intent(out), optional:: GasRDry
                              ! $ R $ [J kg-1 K-1]. 
                              ! 乾燥大気の気体定数. 
                              ! Gas constant of air
    real(DP), intent(out), optional:: MolWtWet
                              ! $ ?? $ [kg mol-1]. 
                              ! 凝結成分の平均分子量. 
                              ! Mean molecular weight of condensible elements
    real(DP), intent(out), optional:: CpWet
                              ! $ C_v $ [J kg-1 K-1] . 
                              ! 凝結成分の定圧比熱. 
                              ! Specific heat of condensible elements at constant pressure
    real(DP), intent(out), optional:: GasRWet
                              ! $ R_v $ [J kg-1 K-1]. 
                              ! 凝結成分の気体定数. 
                              ! Gas constant of condensible elements
    real(DP), intent(out), optional:: LatentHeat
                              ! $ L $ [J kg-1] . 
                              ! 凝結の潜熱. 
                              ! Latent heat of condensation
    real(DP), intent(out), optional:: EpsV
                              ! $ \epsilon_v $ . 
                              ! 水蒸気分子量比. 
                              ! Molecular weight of water vapor
    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 
    character(STRING):: planet_name
    real(DP):: GasRUnivW
                              ! $ R^{*} $ [J K-1 mol-1].
                              ! 普遍気体定数.  Universal gas constant
    real(DP):: MolWtDryW
                              ! $ ?? $ [kg mol-1]. 
                              ! 乾燥大気の平均分子量. 
                              ! Mean molecular weight of dry air
    real(DP):: MolWtWetW
                              ! $ ?? $ [kg mol-1]. 
                              ! 凝結成分の平均分子量. 
                              ! Mean molecular weight of condensible elements


    integer:: stat
    character(STRING):: cause_c
    character(*), parameter:: subname = 'ConstGet'
  continue
    call BeginSub(subname)
    stat = DC_NOERR
    cause_c = ''

    !-----------------------------------------------------------------
    !  初期設定のチェック
    !  Check initialization
    !-----------------------------------------------------------------
    if ( present(planet) ) then
      planet_name = LChar( planet )

      select case ( trim( planet_name ) )
      case ('earth')
      case ('jupiter00')
!!$      case ('mars')
      case ('')
        planet_name = 'earth'
      case default
        stat = DCPAM_ENOPLANET
        cause_c = planet_name
        goto 999
      end select
    else
      planet_name = 'earth'
    end if

    !-----------------------------------------------------------------
    !  普遍定数の提供
    !  Provide universal constants
    !-----------------------------------------------------------------
    GasRUnivW = 8.314_DP

    if ( present(PI)         )  PI        = 3.1415926535897930_DP
    if ( present(GasRUniv)   )  GasRUniv  = GasRUnivW

    if ( present(StB)        )  StB        = 5.67e-8_DP
    if ( present(FKarm)      )  FKarm      = 0.4_DP    

    !-----------------------------------------------------------------
    !  惑星ごとの物理定数の提供
    !  Provide planetary physical constants
    !-----------------------------------------------------------------

    select case ( trim( planet_name ) )
    case ('earth')
      if ( present(RPlanet)    )  RPlanet    = 6.371e6_DP             
      if ( present(Omega)      )  Omega      = 7.29210659088065e-05_DP
                                        ! 2 * \pi / ( 60 * 60 * 23.9345 )
      if ( present(Grav)       )  Grav       = 9.8_DP                 

      if ( present(CpDry)      )  CpDry      = 1004.6_DP
      MolWtDryW = 28.964e-3_DP
      if ( present(MolWtDry)   )  MolWtDry   = MolWtDryW
      if ( present(GasRDry)    )  GasRDry    = GasRUnivW / MolWtDryW

      if ( present(CpWet)      )  CpWet      = 1810.0_DP
      MolWtWetW = 18.01528e-3_DP
      if ( present(MolWtWet)   )  MolWtWet   = MolWtWetW
      if ( present(GasRWet)    )  GasRWet    = GasRUnivW / MolWtWetW
      if ( present(LatentHeat) )  LatentHeat = 2.5e6_DP 

      if ( present(EpsV)       )  EpsV       = MolWtWetW / MolWtDryW

    case ('jupiter00')
      if ( present(RPlanet)    )  RPlanet    = 7.1492e7_DP            
                                        ! 1 気圧の深度での半径. 
                                        ! Radius of depth of 1 atm
      if ( present(Omega)      )  Omega      = 1.75851813802955e-4_DP 
      if ( present(Grav)       )  Grav       = 23.1_DP                
                                        ! For Sugiyama et al. (2008)

      if ( present(CpDry)      )  CpDry      = 11900.9264_DP        
                                        ! For Sugiyama et al. (2008)
      MolWtDryW = 2.3053533e-3
                    ! 分子量の内訳
                    ! * He の分子量と割合 : 4.002602e-3    [kg mol-1], 0.1455
                    ! * H2 の分子量と割合 : 1.00794e-3 * 2 [kg mol-1], 0.8531

      if ( present(MolWtDry)   )  MolWtDry   = MolWtDryW
                                        ! For Sugiyama et al. (2008)
      if ( present(GasRDry)    )  GasRDry    = GasRUnivW / MolWtDryW

      !---------------------------------
      !  凝結成分として水のみ考慮
      !  Condensation component is water only
      if ( present(CpWet)      )  CpWet      = 1810.0_DP
      MolWtWetW = 18.01528e-3_DP
      if ( present(MolWtWet)   )  MolWtWet   = MolWtWetW
      if ( present(GasRWet)    )  GasRWet    = GasRUnivW / MolWtWetW
      if ( present(LatentHeat) )  LatentHeat = 2.5e6_DP 

      if ( present(EpsV)       )  EpsV       = MolWtWetW / MolWtDryW

!!$    case ('mars')
!!$      if ( present(RPlanet)    )  RPlanet    = 3.397e6_DP            
!!$      if ( present(Omega)      )  Omega      = 7.08823595918567e-5_DP
!!$                                        ! 2 * \pi / ( 60 * 60 * 24.6229 )
!!$      if ( present(Grav)       )  Grav       = 3.72_DP                
!!$
!!$      if ( present(CpDry)      )  CpDry      = 
!!$      MolWtDryW = 
!!$      if ( present(MolWtDry)   )  MolWtDry   = MolWtDryW
!!$      if ( present(GasRDry)    )  GasRDry    = GasRUnivW / MolWtDryW
!!$
!!$      if ( present(CpWet)      )  CpWet      = 
!!$      MolWtWetW = 
!!$      if ( present(MolWtWet)   )  MolWtWet   = MolWtWetW
!!$      if ( present(GasRWet)    )  GasRWet    = GasRUnivW / MolWtWetW
!!$      if ( present(LatentHeat) )  LatentHeat = 
!!$
!!$      if ( present(EpsV)       )  EpsV       = MolWtWetW / MolWtDryW

    end select

    !-----------------------------------------------------------------
    !  終了処理, 例外処理
    !  Termination and Exception handling
    !-----------------------------------------------------------------
999 continue
    call StoreError(stat, subname, err, cause_c)
    call EndSub(subname)
  end subroutine ConstGet

Private Instance methods

version
Constant :
version = ’$Name: dcpam4-20080626 $’ // ’$Id: const_provider.f90,v 1.2 2008-04-22 12:41:13 morikawa Exp $’ :character(*), parameter

[Validate]