!= Subroutine BasicEnvInit
!
! Authors::   SUGIYAMA Koichiro, ODAKA Masatsugu
! Version::   $Id: basicenv_init.f90,v 1.3 2009-08-24 03:50:12 sugiyama Exp $
! Tag Name::  $Name:  $
! Copyright:: Copyright (C) GFD Dennou Club, 2006. All rights reserved.
! License::   See COPYRIGHT[link:../../COPYRIGHT]
!
!== Overview
!
!デフォルトの基本場を設定するための変数参照型モジュール
!   * BasicEnvFile_init: 基本場の値を netCDF ファイルから取得
!   * BasicEnvCalc_Init: 基本場の情報を Namelist から取得して値を計算
!
!== Error Handling
!
!== Known Bugs
!
!== Note
!
!== Future Plans
!
subroutine BasicEnv_Init( myrank, nprocs, cfgfile)
  !
  !デフォルトの基本場を設定するためのサブルーチン.
  !基本場を計算し, BasicSet モジュールの値を初期化する.
  !
  !コンパイルの順序の問題から, 基本場の値(hogeBasicZ な変数)を
  !計算する部分をBasicSet モジュールから切り離している.
  !ECCM 始め, BasicSet 自体に依存するが hogeBasicZ は use しない
  !外部サブルーチンを利用するためである.
  !

  !モジュール読み込み
  use dc_message, only: MessageNotify

  use gridset,  only: DimXMin,       &!配列の X 方向の下限
    &                 DimXMax,       &!配列の X 方向の上限
    &                 DimZMin,       &!配列の Z 方向の下限
    &                 DimZMax,       &!配列の Z 方向の上限
    &                 SpcNum          !凝縮成分の数
  use basicset, only: BasicSetArray_Init, &!
    &                 PressBasis,    &!温位の基準圧力
    &                 GasRDry,       &!乾燥成分の定圧比熱
    &                 CpDry,         &!乾燥成分の定圧比熱
    &                 CvDry,         &!乾燥成分の定積比熱
    &                 MolWtDry,      &!乾燥成分の分子量
    &                 SpcWetMolFr,   &!凝縮成分の初期モル比
    &                 MolWtWet,      &!凝縮成分の分子量
    &                 GasRUniv        !気体定数
  use Boundary,only:  BoundaryXCyc_xz,  &!
    &                 BoundaryZSym_xz,  &!
    &                 BoundaryXCyc_xza, &!
    &                 BoundaryZSym_xza   !
  use ECCM,     only: ECCM_Dry,         &!
    &                 ECCM_Wet

  !暗黙の型宣言禁止
  implicit none

  !変数の定義
  integer, intent(in)      :: myrank, nprocs
  character(*), intent(in) :: cfgfile
  real(8)                  :: xz_DensBasicZ(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)                  :: xz_PressBasicZ(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)                  :: xz_ExnerBasicZ(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)                  :: xz_TempBasicZ(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)                  :: xz_PotTempBasicZ(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)                  :: xz_VelSoundBasicZ(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)                  :: xza_MixRtBasicZ(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
  real(8)                  :: xz_EffMolWtBasicZ(DimXMin:DimXMax, DimZMin:DimZMax)
  real(8)                  :: z_TempBasicZ(DimZMin:DimZMax)
  real(8)                  :: z_PressBasicZ(DimZMin:DimZMax)
  real(8)                  :: MolFrIni(SpcNum)
  real(8)                  :: xza_MolFr(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
  real(8)                  :: za_MolFr(DimZMin:DimZMax, SpcNum)
  real(8)                  :: xza_MixRtDivMolWt(DimXMin:DimXMax, DimZMin:DimZMax, SpcNum)
  real(8)                  :: Humidity              !相対湿度
  character(20)            :: EnvType               !基本場の温度設定, 'Dry' or 'Moist'
  real(8)                  :: TempStrat             !成層圏の温度 [k]
  real(8)                  :: Dhight                !重み関数のパラメータ [m]
  real(8)                  :: z_llhum(DimZMin:DimZMax)
  integer                  :: i, s

  !---------------------------------------------------------------
  ! NAMELIST の定義
  !---------------------------------------------------------------
  NAMELIST /basicenv/ EnvType
  NAMELIST /basicenv_dry/ Humidity, TempStrat, Dhight
  NAMELIST /basicenv_wet/ Humidity, TempStrat, Dhight

  !---------------------------------------------------------------
  ! 配列の初期化
  !---------------------------------------------------------------
  xz_PressBasicZ    = 0.0d0
  xz_ExnerBasicZ    = 0.0d0
  xz_TempBasicZ     = 0.0d0
  xz_PotTempBasicZ  = 0.0d0
  xz_VelSoundBasicZ = 0.0d0
  xza_MixRtBasicZ   = 0.0d0
  xz_EffMolWtBasicZ = 0.0d0
  z_TempBasicZ      = 0.0d0
  z_PressBasicZ     = 0.0d0
  za_MolFr          = 0.0d0

  !---------------------------------------------------------------
  ! EnvType を元に, 温度, 圧力, 組成を決める
  !---------------------------------------------------------------
  MolFrIni = SpcWetMolFr(1:SpcNum)
  open (10, FILE=cfgfile)
  read(10, NML=basicenv)
  close(10)

  select case(EnvType)
  case("Dry") ! 乾燥断熱的な初期場
     ! NAMELIST ファイルの読み込み
     open (10, FILE=cfgfile)
     read(10, NML=basicenv_dry)
     close(10)

     call ECCM_Dry( MolFrIni, Humidity, z_TempBasicZ, z_PressBasicZ, za_MolFr )
     call basicset_upper(Humidity, TempStrat, Dhight, z_TempBasicZ, z_PressBasicZ)

  case("Wet") ! 湿潤断熱的な初期場
     ! NAMELIST ファイルの読み込み
     open (10, FILE=cfgfile)
     read(10, NML=basicenv_wet)
     close(10)

     call ECCM_Wet( MolFrIni, Humidity, z_TempBasicZ, z_PressBasicZ, za_MolFr )
     call basicset_upper(Humidity, TempStrat, Dhight, z_TempBasicZ, z_PressBasicZ)

  case("N1994") ! Yamasaki(1983)の温度と相対湿度の観測値を使用する場合
     call ECCM_N1994( z_TempBasicZ, z_PressBasicZ, za_MolFr )

  case("Takemi2007") ! Takemi(2007)の基本場を使用する場合
     call ECCM_Takemi2007( z_TempBasicZ, z_PressBasicZ, za_MolFr, z_llhum )
  end select

  ! 2 次元配列に格納
  do i = DimXMin, DimXMax
    xz_TempBasicZ(i,:)  = z_TempBasicZ
    xz_PressBasicZ(i,:) = z_PressBasicZ
  end do

  !境界条件
  call BoundaryXCyc_xz( xz_TempBasicZ )
  call BoundaryZSym_xz( xz_TempBasicZ )
  call BoundaryXCyc_xz( xz_PressBasicZ )
  call BoundaryZSym_xz( xz_PressBasicZ )

  !---------------------------------------------------------------
  ! 混合比
  !---------------------------------------------------------------
  !水平方向には一様
  do i = DimXMin, DimXMax
    xza_MolFr(i,:,:) = za_MolFr
  end do

  !気相のモル比を混合比に変換
  do s = 1, SpcNum
    xza_MixRtBasicZ(:,:,s) = xza_MolFr(:,:,s) * MolWtWet(s) / MolWtDry
  end do

  !境界条件
  call BoundaryXCyc_xza( xza_MixRtBasicZ )
  call BoundaryZSym_xza( xza_MixRtBasicZ )

  !---------------------------------------------------------------
  ! 分子量の効果
  !---------------------------------------------------------------
  do s = 1, SpcNum
    xza_MixRtDivMolWt(:,:,s) = xza_MixRtBasicZ(:,:,s) / MolWtWet(s)
  end do

  xz_EffMolWtBasicZ = &
    & (1.0d0 + sum(xza_MixRtBasicZ,3) ) &
    & / ( MolWtDry * ((1.0d0 / MolWtDry) + sum(xza_MixRtDivMolWt,3)) )

  !境界条件
  call BoundaryXCyc_xz( xz_EffMolWtBasicZ )
  call BoundaryZSym_xz( xz_EffMolWtBasicZ )

  !---------------------------------------------------------------
  ! 温位
  !---------------------------------------------------------------
  xz_PotTempBasicZ = &
    & xz_TempBasicZ * (PressBasis / xz_PressBasicZ) ** (GasRDry / CpDry)

  !境界条件
  call BoundaryXCyc_xz( xz_PotTempBasicZ )
  call BoundaryZSym_xz( xz_PotTempBasicZ )

  !---------------------------------------------------------------
  ! エクスナー関数
  !---------------------------------------------------------------
  xz_ExnerBasicZ = xz_TempBasicZ / xz_PotTempBasicZ

  !境界条件
  call BoundaryXCyc_xz( xz_ExnerBasicZ )
  call BoundaryZSym_xz( xz_ExnerBasicZ )

  !---------------------------------------------------------------
  ! 密度
  !---------------------------------------------------------------
  xz_DensBasicZ = &
    & PressBasis * (xz_ExnerBasicZ ** (CvDry / GasRDry)) &
    &  / (GasRDry * xz_PotTempBasicZ / xz_EffMolWtBasicZ)

  !境界条件
  call BoundaryXCyc_xz( xz_DensBasicZ )
  call BoundaryZSym_xz( xz_DensBasicZ )

  !---------------------------------------------------------------
  ! 音速
  !---------------------------------------------------------------
  xz_VelSoundBasicZ = &
    & sqrt ( &
    &   CpDry * GasRDry * xz_ExnerBasicZ * xz_PotTempBasicZ &
    &   / (CvDry * xz_EffMolWtBasicZ) &
    & )

  !境界条件
  call BoundaryXCyc_xz( xz_VelSoundBasicZ )
  call BoundaryZSym_xz( xz_VelSoundBasicZ )

  !----------------------------------------------------------
  ! BasicSet モジュールに値を設定
  !----------------------------------------------------------
  call BasicSetArray_Init(                                  &
    & xz_PressBasicZ,    xz_ExnerBasicZ, xz_TempBasicZ,     &
    & xz_PotTempBasicZ,  xz_DensBasicZ,  xz_VelSoundBasicZ, &
    & xza_MixRtBasicZ, xz_EffMolWtBasicZ )

end subroutine BasicEnv_Init