!= Module Constants
!
! Authors::   SUGIYAMA Ko-ichiro, ODAKA Masatsugu
! Version::   $Id: constants.f90,v 1.8 2011-06-17 19:04:00 sugiyama Exp $
! Tag Name::  $Name: arare5-20111010 $
! Copyright:: Copyright (C) GFD Dennou Club, 2007. All rights reserved.
! License::   See COPYRIGHT[link:../../COPYRIGHT]
!
!== Overview 
!
!
!== Error Handling
!
!== Known Bugs
!
!== Note
!
!== Future Plans
!

module constants
  !
  !

  !モジュール読み込み
  use dc_types,   only: DP, STRING
  use dc_iounit,  only: FileOpen
  use dc_message, only: MessageNotify   !メッセージ出力
  use ChemData,   only: GasRUniv,             &!普遍気体定数
    &                   ChemData_OneSpcID,    &!化学種の ID
    &                   ChemData_CpPerMolRef, &!標準状態での比熱
    &                   ChemData_MolWt         !分子量
  use mpi_wrapper, only: myrank
  use namelist_util, only: namelist_filename

  !暗黙の型宣言禁止
  implicit none

  private

  !Public Interface
  real(DP), save, public :: Grav = 9.8d0          !重力 [m/s^2]
  real(DP), save, public :: PressBasis = 965.0d0  !温位の基準圧力 [Pa]
  real(DP), save, public :: TempSfc = 300.0d0     !地表面温度 [K]
  real(DP), save, public :: PressSfc = 965.0d0    !地表面圧力 [Pa]
  real(DP), save, public :: CpDry  = 0.0d0        !乾燥成分の定圧比熱 [J/K kg]
  real(DP), save, public :: CpDryMol = 0.0d0      !乾燥成分の定圧比熱 [J/K kg]
  real(DP), save, public :: CvDry = 0.0d0         !乾燥成分の定積比熱 [J/K kg]
  real(DP), save, public :: MolWtDry = 0.0d0      !乾燥成分の分子量   [kg/mol]
  real(DP), save, public :: GasRDry  = 0.0d0      !乾燥成分の気体定数 [J/K kg]
  real(DP), save, public :: DayTime = 86400.0d0   ! 1 日の長さ [s]

  public constants_init

contains

!!!-----------------------------------------------------------------!!!
  subroutine constants_Init()

    !暗黙の型宣言禁止
    implicit none
    
    !変数定義
    integer                  :: SpcDryNum       !乾燥成分の化学種の数
    character(20)            :: SpcDrySymbol(5) !乾燥成分の化学種名
    real(DP)                 :: SpcDryMolFr(5)  !乾燥成分の化学種の存在度
    integer, allocatable     :: SpcDryID(:)     !乾燥成分の化学種のID
    real(DP), allocatable    :: PropertyDry(:)  !作業配列
    integer                  :: s               !作業変数
    integer                  :: unit            !装置番号
    logical                  :: flag = .false.
     
    !NAMELIST の定義
    NAMELIST /constants_nml/ &
      & Grav, PressBasis, TempSfc, PressSfc, & 
      & SpcDrySymbol, SpcDryMolFr, DayTime, CpDry, MolWtDry
 
    SpcDrySymbol = '' 
    SpcDryMolFr  = 0.0d0
    
    !ファイルオープン. 情報取得. 
    call FileOpen(unit, file=namelist_filename, mode='r')
    read(unit, NML=constants_nml)
    close(unit)

    if (CpDry /= 0.0d0 .AND. MolWtDry /= 0.0d0) then 
    ! namelist から CpDry の値が入力された場合

      !定圧モル比熱
      CpDryMol = CpDry *  MolWtDry 
      
      !気体定数
      GasRDry = GasRUniv / MolWtDry
      
      !定積比熱
      CvDry    = CpDry - GasRDry

    elseif (SpcDryMolFr(1) /= 0.0d0) then 
    ! namelist から モル比が入力された場合

      flag = .true.

      !----------------------------------------------------------
      ! 乾燥成分の物性値の初期化
      !
      !乾燥成分の個数を数える
      SpcDryNum = count(SpcDrySymbol /= "")
      
      !化学種の ID を取得    
      allocate(SpcDryID(SpcDryNum))    
      do s = 1, SpcDryNum
        SpcDryID(s) = ChemData_OneSpcID( SpcDrySymbol(s) )
      end do
      
      !作業配列の準備
      allocate(PropertyDry(SpcDryNum))
      
      !分子量
      do s = 1, SpcDryNum
        PropertyDry(s) = ChemData_MolWt(SpcDryID(s))
      end do
      MolWtDry = dot_product(PropertyDry, SpcDryMolFr(1:SpcDryNum)) 
      
      !定圧比熱(モル当量)
      do s = 1, SpcDryNum    
        PropertyDry(s) = ChemData_CpPerMolRef(SpcDryID(s))
      end do
      CpDryMol = dot_product(PropertyDry, SpcDryMolFr(1:SpcDryNum)) 
      
      !定圧比熱
      CpDry    = CpDryMol / MolWtDry
      
      !気体定数
      GasRDry = GasRUniv / MolWtDry
      
      !定積比熱
      CvDry    = CpDry - GasRDry
    end if
    
    !----------------------------------------------------------
    ! 確認
    !----------------------------------------------------------
    if (myrank == 0) then 
      call MessageNotify( "M", &
        & "constants_init", "Grav = %f",     d=(/Grav/) )
      call MessageNotify( "M", &
        & "constants_init", "PressBasis = %f", d=(/PressBasis/))
      call MessageNotify( "M", &
        & "constants_init", "TempSfc = %f",  d=(/TempSfc/) )
      call MessageNotify( "M", &
        & "constants_init", "PressSfc = %f", d=(/PressSfc/) )
      if (flag) then 
        do s = 1, SpcDryNum
          call MessageNotify( "M", &
            &  "constants_init", "SpcDryID = %d",      i=(/SpcDryID(s)/))
          call MessageNotify( "M", &
            &  "constants_init", "SpcDrySymbol = %c", c1=trim(SpcDrySymbol(s)))
          call MessageNotify( "M", &
            &  "constants_init", "SpcDryMolFr = %f",   d=(/SpcDryMolFr(s)/))
        end do
      end if
      
      call MessageNotify( "M", "constants_init", "CpDry    = %f",    d=(/CpDry/) )
      call MessageNotify( "M", "constants_init", "CpDryMol = %f", d=(/CpDryMol/) )
      call MessageNotify( "M", "constants_init", "CvDry    = %f",    d=(/CvDry/) )
      call MessageNotify( "M", "constants_init", "GasRDry  = %f",  d=(/GasRDry/) )
      call MessageNotify( "M", "constants_init", "MolWtDry = %f", d=(/MolWtDry/) )
      call MessageNotify( "M", "constants_init", "DayTime  = %f",  d=(/DayTime/)  )
    end if

  end subroutine Constants_Init

end module Constants