module AtmBasicState_Builder

  !* gtool ***************
  !*

  ! 種類型パラメータ
  ! Kind type parameters
  !
  use dc_types, only: &
    & DP

  !* IGMBaseLib ***********
  !*


  use Field_Grid1D_Manager, only: &
    & Field_Grid1D

  implicit none
  private

  public :: create_isothermalAtm_BasicState, &
    &       create_constBuoFreqAtm_BasicState

contains
subroutine create_isothermalAtm_BasicState( &
  & z_RhoBasic, z_PressBasic,         &
  & GasConst, Grav, refTemp, refPress &
  & )

  type(Field_Grid1D), intent(out) :: z_RhoBasic
  type(Field_Grid1D), intent(out) :: z_PressBasic
  real(DP), intent(in) :: GasConst
  real(DP), intent(in) :: Grav
  real(DP), intent(in) :: refTemp
  real(DP), intent(in) :: refPress

  real(DP), pointer :: AltZ(:)
  integer :: k, AltZkMax
  real(DP) :: ScaleHeightT0

  Altz => z_RhoBasic%grid1D(:)
  AltZkmax = size(AltZ)
  ScaleHeightT0 = GasConst * refTemp / Grav

  do k=1, AltZkmax
    z_PressBasic%val(k,1) = refPress * exp( - AltZ(k) / ScaleHeightT0 )
    z_RhoBasic%val(k,1) = z_PressBasic%val(k,1) / ( GasConst * refTemp )
  end do

end subroutine create_isothermalAtm_BasicState

subroutine create_constBuoFreqAtm_BasicState( &
  & z_RhoBasic, z_PressBasic,    &
  & N, GasConst, Grav, Cp, Cv,   &
  & refTemp, refPress &
  & )

  type(Field_Grid1D), intent(inout) :: z_RhoBasic
  type(Field_Grid1D), intent(inout) :: z_PressBasic
  real(DP), intent(in) :: N
  real(DP), intent(in) :: GasConst
  real(DP), intent(in) :: Grav
  real(DP), intent(in) :: Cp
  real(DP), intent(in) :: Cv
  real(DP), intent(in) :: refTemp
  real(DP), intent(in) :: refPress

  real(DP), pointer :: AltZ(:)
  integer :: k, AltZkMax

  real(DP) :: SpecHeatRatio, Kapper, N2
  real(DP) :: Coef1, Coef2, tmp

  Altz => z_RhoBasic%grid1D(:)
  AltZkmax = size(AltZ)

  N2 = N**2
  Kapper = GasConst / Cp
  SpecHeatRatio = Cp / Cv
  Coef1 = ( Grav / N )**2 / ( Cp * refTemp )
  Coef2 = refPress**Kapper / ( GasConst * refTemp )

  do k=1, AltZkmax

    tmp = 1.0d0 - Coef1 * ( 1.0d0 - exp( - N2 / Grav * AltZ(k) ) )
    z_PressBasic%val(k,1) = refPress * tmp**( 1.0d0 / kapper )
    z_RhoBasic%val(k,1) = z_PressBasic%val(k,1) / ( GasConst * refTemp * tmp ) * exp(- ( N2 / Grav ) * AltZ(k))

  end do

end subroutine create_constBuoFreqAtm_BasicState

end module AtmBasicState_Builder
