Class BasicSet
In: setup/basicset.f90

    Copyright (C) GFD Dennou Club, 2004, 2005, 2006. All rights reserved.

Module BasicSet

  * Developer: SUGIYAMA Ko-ichiro
  * Version: $Id: basicset.f90,v 1.1.1.1 2006/04/25 03:43:58 deepconv Exp $
  * Tag Name: $Name:  $
  * Change History:

Overview

デフォルトの基本場を設定するための変数参照型モジュール

  * BasicEnvFile_init: 基本場の値を netCDF ファイルから取得
  * BasicEnvCalc_Init: 基本場の情報を Namelist から取得して値を計算

Error Handling

Known Bugs

Note

Future Plans

Methods

Included Modules

dc_message gridset ChemCalc ChemData

Public Instance methods

[Source]

subroutine BasicSetArray_Init (                         
     xz_Press, xz_Exner, xz_Temp, xz_PotTemp, xz_Dens,  
     xz_VelSound, xz_MixRt, xz_EffMolWt )


    !
    !基本場の値を外部から取得する. 
    !
    
    !入力変数
    real(8), intent(in) :: xz_Press(DimXMin:DimXMax, DimZMin:DimZMax)
    real(8), intent(in) :: xz_Exner(DimXMin:DimXMax, DimZMin:DimZMax)
    real(8), intent(in) :: xz_Temp(DimXMin:DimXMax, DimZMin:DimZMax)
    real(8), intent(in) :: xz_PotTemp(DimXMin:DimXMax, DimZMin:DimZMax)
    real(8), intent(in) :: xz_Dens(DimXMin:DimXMax, DimZMin:DimZMax)
    real(8), intent(in) :: xz_VelSound(DimXMin:DimXMax, DimZMin:DimZMax)
    real(8), intent(in) :: xz_MixRt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
    real(8), intent(in) :: xz_EffMolWt(DimXMin:DimXMax, DimZMin:DimZMax)
   
    !---------------------------------------------------------------
    ! *BasicZ な配列の初期化
    !---------------------------------------------------------------
    !配列の割り当て
    allocate( & 
      & xz_DensBasicZ(DimXMin:DimXMax, DimZMin:DimZMax),            &
      & xz_PressBasicZ(DimXMin:DimXMax, DimZMin:DimZMax),           &
      & xz_ExnerBasicZ(DimXMin:DimXMax, DimZMin:DimZMax),           &
      & xz_TempBasicZ(DimXMin:DimXMax, DimZMin:DimZMax),            &
      & xz_PotTempBasicZ(DimXMin:DimXMax, DimZMin:DimZMax),         &
      & xz_VelSoundBasicZ(DimXMin:DimXMax, DimZMin:DimZMax),        &
      & xz_MixRtBasicZ(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum),   &
      & xz_EffMolWtBasicZ(DimXMin:DimXMax, DimZMin:DimZMax)         )
    
    !値の入力
    xz_PressBasicZ    = xz_Press
    xz_ExnerBasicZ    = xz_Exner
    xz_TempBasicZ     = xz_Temp
    xz_PotTempBasicZ  = xz_PotTemp
    xz_DensBasicZ     = xz_Dens
    xz_VelSoundBasicZ = xz_VelSound
    xz_MixRtBasicZ    = xz_MixRt
    xz_EffMolWtBasicZ = xz_EffMolWt

end subroutine

[Source]

subroutine BasicSet_Init (cfgfile)

    !
    !基本場の情報を Namelist から取得して値を計算
    !

    !暗黙の型宣言禁止
    implicit none
    
    !変数定義
    character(*), intent(in) :: cfgfile         !NAMELIST ファイル
!    real(8)                   :: Tropopause      !対流圏圏界面高度
!    real(8)                   :: CpDryMol        !乾燥成分の定圧比熱 [J/K mol]
    integer                   :: SpcDryNum       !乾燥成分の化学種の数
    integer                   :: SpcWetNum       !湿潤成分の化学種の数
    character(15)             :: SpcDrySymbol(5)!乾燥成分の化学種名
!    character(15)             :: SpcWetSymbol(10)!湿潤成分の化学種名
    real(8)                   :: SpcDryMolFr(5) !乾燥成分の化学種の存在度
!    real(8)                   :: SpcWetMolFr(10) !湿潤成分の化学種の存在度
    integer, allocatable      :: SpcDryID(:)     !乾燥成分の化学種のID
!    real(8)                   :: Humidity        !相対湿度
!    character(20)             :: EnvType         !基本場の温度設定, 'Dry' or 'Moist'
    real(8), allocatable      :: PropertyDry(:)  !作業配列
    character(20), allocatable:: Symbol(:)       !作業配列
    integer                   :: s, n1, n2, n3   !作業変数

    !NAMELIST の定義
    NAMELIST /basicset/ &
      & Grav, TempSfc, PressSfc,   Tropopause,  & 
      & SpcDrySymbol, SpcDryMolFr,              &
      & SpcWetSymbol, SpcWetMolFr,              &
      & EnvType, Humidity                       

    SpcDrySymbol = '' 
    SpcDryMolFr  = 0.0d0
    SpcWetSymbol = '' 
    SpcWetMolFr  = 0.0d0
    
    !ファイルオープン. 情報取得. 
    open (10, FILE=cfgfile)
    read(10, NML=basicset)
    close(10)

    !----------------------------------------------------------
    ! 化学種の物性値の初期化
    !----------------------------------------------------------
    !データベースの初期化
    call ChemData_Init
    
    !----------------------------------------------------------
    ! 乾燥成分の物性値の初期化
    !----------------------------------------------------------
    !乾燥成分の個数を数える
    SpcDryNum = count(SpcDrySymbol /= "")
    
    !化学種の ID を取得    
    allocate(SpcDryID(SpcDryNum))    
    SpcDryID = ChemData_SpcID(SpcDryNum, SpcDrySymbol(1:SpcDryNum))
    
    !化学計算モジュールの初期化
    call ChemCalc_Init(DimXMin, DimXMax, DimZMin, DimZMax, SpcDryNum, SpcDryID)
    
    !作業配列の準備
    allocate(PropertyDry(SpcDryNum))

    !定圧比熱
    do s = 1, SpcDryNum
      PropertyDry(s) = CpPerMassRef(SpcDryID(s))
    end do
    CpDry    = dot_product(PropertyDry, SpcDryMolFr(1:SpcDryNum)) 
    
    !定圧比熱(モル当量)
    do s = 1, SpcDryNum    
      PropertyDry(s) = CpPerMolRef(SpcDryID(s))
    end do
    CpDryMol = dot_product(PropertyDry, SpcDryMolFr(1:SpcDryNum)) 
    
    !定積比熱
    do s = 1, SpcDryNum
      PropertyDry(s) = CvPerMassRef(SpcDryID(s))
    end do
    CvDry    = dot_product(PropertyDry, SpcDryMolFr(1:SpcDryNum)) 
    
    !分子量
    do s = 1, SpcDryNum
      PropertyDry(s) = MolWt(SpcDryID(s))
    end do
    MolWtDry = dot_product(PropertyDry, SpcDryMolFr(1:SpcDryNum)) 
    
    !気体定数
    do s = 1, SpcDryNum
      PropertyDry(s) = GasR(SpcDryID(s))
    end do
    GasRDry = dot_product(PropertyDry, SpcDryMolFr(1:SpcDryNum)) 
    
    !----------------------------------------------------------
    ! 湿潤成分の ID を得る
    !----------------------------------------------------------
    !湿潤成分の個数を数える
    SpcWetNum = count(SpcWetSymbol /= "")
    if (SpcWetNum /= SpcNum) then 
      write(*,*) "SpcWetNum /= SpcNum"
      stop
    end if
    
    !配列の割り当て
    allocate(SpcWetID(SpcWetNum), Symbol(SpcWetNum), MolWtWet(SpcWetNum))

    !SpcWetSymbol の文字列から, -Rain, -Cloud を除いたものを Symbol として保管
    do s = 1, SpcWetNum
      n1 = index(SpcWetSymbol(s), '-Cloud' )
      n2 = index(SpcWetSymbol(s), '-Rain' )
      n3 = max(n1, n2)
      if (n3 == 0) then
        Symbol(s) = SpcWetSymbol(s)
      else
        Symbol(s) = SpcWetSymbol(s)(1:n3-1)
      end if
    end do
    
    !化学種の ID を取得
    SpcWetID = ChemData_SpcID(SpcWetNum, Symbol)
    
    !化学計算モジュールの初期化
    call ChemCalc_Init(DimXMin, DimXMax, DimZMin, DimZMax, SpcWetNum, SpcWetID(:))
    
    !分子量を保管
    do s = 1, SpcWetNum
      MolWtWet(s) = MolWt(SpcWetID(s))
    end do
    

    !----------------------------------------------------------
    ! 確認
    !----------------------------------------------------------
    write(*,*) "basicset_init, Grav     ", Grav
    write(*,*) "basicset_init, TempSfc  ", TempSfc
    write(*,*) "basicset_init, PressSfc ", PressSfc

    do s = 1, SpcDryNum
      write(*,*) "basicset_init, SpcDryID ", SpcDryID(s)
      write(*,*) "basicset_init, SpcDrySymbol ", SpcDrySymbol(s)
      write(*,*) "basicset_init, SpcDryMolFr  ", SpcDryMolFr(s)
    end do
    write(*,*) "basicset_init, CpDry    ", CpDry
    write(*,*) "basicset_init, CpDryMol ", CpDryMol
    write(*,*) "basicset_init, CvDry    ", CvDry
    write(*,*) "basicset_init, GasRDry  ", GasRDry
    write(*,*) "basicset_init, MolWtDry ", MolWtDry

    do s = 1, SpcWetNum
      write(*,*) "basicset_init, SpcWetID ", SpcWetID(s)
      write(*,*) "basicset_init, SpcWetSymbol ", SpcWetSymbol(s)
      write(*,*) "basicset_init, SpcWetMolFr  ", SpcWetMolFr(s)
      write(*,*) "basicset_init, MolWtWet ", MolWtWet(s)
    end do

    write(*,*) "basicset_init, Tropopause ", Tropopause
    write(*,*) "basicset_init, EnvType    ", EnvType
    write(*,*) "basicset_init, Humidity   ", Humidity

end subroutine

[Validate]