module shr_frz_mod

  !===============================================================================
  !  This is a module used for the freezing point of salt water
  !===============================================================================

  use shr_kind_mod, only: R8=>SHR_KIND_R8, CS=>SHR_KIND_CS
  use shr_log_mod, only: s_logunit => shr_log_Unit, shr_log_level
  use shr_sys_mod, only: shr_sys_abort

  implicit none

  !----------------------------------------------------------------------------
  ! PUBLIC: Interfaces and global data
  !----------------------------------------------------------------------------
  public :: shr_frz_freezetemp, shr_frz_freezetemp_init

  interface shr_frz_freezetemp
     module procedure shr_frz_freezetemp_0d
     module procedure shr_frz_freezetemp_1d
     module procedure shr_frz_freezetemp_2d
  end interface shr_frz_freezetemp

  integer, public, parameter :: TFREEZE_OPTION_MINUS1P8 = 1
  integer, public, parameter :: TFREEZE_OPTION_LINEAR_SALT = 2
  integer, public, parameter :: TFREEZE_OPTION_MUSHY = 3
  integer, public, parameter :: TFREEZE_OPTION_UNINITIALIZED = -999

  private

  integer :: tfrz_option = TFREEZE_OPTION_UNINITIALIZED

  !===============================================================================
contains
  !===============================================================================

  subroutine shr_frz_freezetemp_init(tfreeze_option, mastertask)

    implicit none

    character(len=*),parameter :: subname = "(shr_frz_freezetemp_init) "
    character(CS),intent(in)  :: tfreeze_option ! option for computing freezing point
    logical, intent(in)  :: mastertask ! for io
    ! minus1p8 is constant -1.8C
    ! linear_salt is linear equation
    ! mushy for CICE mushy-layer nonlinear equation

    !---------------------------------------------------------------
    ! Check tfreeze_option
    !---------------------------------------------------------------
    if (trim(tfreeze_option) == 'minus1p8') then
       if (mastertask .and. shr_log_level>0) write(s_logunit,*) ' tfreeze_option is minus1p8'
       tfrz_option = TFREEZE_OPTION_MINUS1P8
    elseif (trim(tfreeze_option) == 'linear_salt') then
       if (mastertask .and. shr_log_level>0) write(s_logunit,*) ' tfreeze_option is linear_salt'
       tfrz_option = TFREEZE_OPTION_LINEAR_SALT
    elseif (trim(tfreeze_option) == 'mushy') then
       if (mastertask .and. shr_log_level>0) write(s_logunit,*) ' tfreeze_option is mushy'
       tfrz_option = TFREEZE_OPTION_MUSHY
    else
       call shr_sys_abort(subname//' ERROR: not a valid tfreeze_option '//trim(tfreeze_option))
    endif

  end subroutine shr_frz_freezetemp_init

  ! DIMS 0,1,2
  function shr_frz_freezetemp_{DIMS}d(s) result(shr_frz_freezetemp)

    !----------------------------------------------------------------------------
    !
    ! FUNCTION to return the freezing point of salt water in degrees Celsus
    !
    !--------------- Code History -----------------------------------------------
    !
    ! Original Author: David Bailey
    ! Date:            Feb, 2016
    !----------------------------------------------------------------------------

    implicit none

    character(len=*),parameter :: subname = "(shr_frz_freezetemp_{DIMS}d) "

    real   (R8),intent(in) :: s{DIMSTR} ! Salinity in psu
#if ({DIMS}==0)
    real   (R8) :: shr_frz_freezetemp
#elif ({DIMS}==1)
    real   (R8) :: shr_frz_freezetemp(size(s))
#elif ({DIMS}==2)
    real   (R8) :: shr_frz_freezetemp(size(s,1),size(s,2))
#endif

    !----------------------------------------------------------------------------
    shr_frz_freezetemp = -274.0_R8
    if (tfrz_option == TFREEZE_OPTION_MINUS1P8) then
       shr_frz_freezetemp = -1.8_R8
    elseif (tfrz_option == TFREEZE_OPTION_LINEAR_SALT) then
       shr_frz_freezetemp = -0.0544_R8*max(s,0.0_R8)
    elseif (tfrz_option == TFREEZE_OPTION_MUSHY) then
       !     This form is the high temperature part of the liquidus relation (Assur 1958)
       shr_frz_freezetemp = max(s,0.0_R8) &
            / (-18.48_R8 + (0.01848_R8*max(s,0.0_R8)))
    else
       call shr_sys_abort(subname//' ERROR: not intialized correctly with a valid tfreeze_option - &
            &call shr_frz_freezetemp_init first with a valid tfreeze_option')
    endif

#if ({DIMS}==0)
    if (shr_frz_freezetemp < -2.6_R8) &
       call shr_sys_abort(subname//' ERROR: shr_frz_freezetemp is less than -2.6')
#else
    if (minval(shr_frz_freezetemp) < -2.6_R8) &
       call shr_sys_abort(subname//' ERROR: shr_frz_freezetemp is less than -2.6')
#endif

  end function shr_frz_freezetemp_{DIMS}d

  !===============================================================================

end module shr_frz_mod
