!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
!                                                                      !
!  Subroutine: USR_PROP_Ks                                             !
!  Purpose: User hook for calculating solids phase conductivity.       !
!                                                                      !
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
      SUBROUTINE USR_PROP_Ks(IJK,M)

      use error_manager
      use fldvar, only: ro_s, T_s, X_s, ep_s, ep_g, rop_s
      use param1, only: undefined_i, zero, one, half
      use physprop, only: mw_s, nmax, k_s, k_g
      use run, only: units
      use toleranc, only: dil_ep_s, tmax
      use usr, only: calc_cond_liq
      use usr, only: index_liq
      use usr, only: solvent_absorption

      implicit none

! Dummy arguments
!---------------------------------------------------------------------//
! index
      INTEGER, INTENT(IN) :: IJK
! solids phase index
      INTEGER, INTENT(IN) :: M

! Local Variables:
!---------------------------------------------------------------------//
! error flag
      INTEGER :: IER = undefined_i
      CHARACTER(LEN=40) :: err_prop
! bounded phase temperatures (K)
      DOUBLE PRECISION :: xTl
! liquid phase mass fraction
      DOUBLE PRECISION :: x_lH2O
      DOUBLE PRECISION :: x_lRNH2
! volume fraction of liquid volume fraction of liquid
      DOUBLE PRECISION :: ep_liq
! thermal conductivity
      DOUBLE PRECISION :: K_aqMEA

      INTEGER :: iliquid
      INCLUDE 'species.inc'
!......................................................................!
! if using this quantity then remove definition of ier
!      ier = 1
      IF (.NOT.SOLVENT_ABSORPTION) RETURN

      iliquid=index_liq
! routine will fail if not being called for liquid phase only
      IF (M /= iliquid) ier =1

! volume fraction of liquid
      ep_liq = ep_s(ijk,iliquid)

!      IF (ep_liq > DIL_EP_S) THEN
! bounded liquid phase temperature
         xTl = min(TMAX,T_s(IJK,iliquid))

! mass fraction of MEA
         x_lRNH2 = X_S(IJK,iliquid,lRNH2)

! mass fraction of H2O
         x_lH2O = X_S(IJK,iliquid,lH2O)

END MODULE CALC_K_s_MOD


  real :: K_g, K_s, K_go, ep_g, ep_s, beta, A, B, K

  ! Define the given constants
  ep_g = 0.48 ! Replace with your desired value
  K_go = 0.0257   ! Replace with your desired value
  K_so = 0.29		! Replace with your desired value
  beta = 7.26E-3

  ! Calculate B and A
  B = 1.25 * ((1.0 - ep_g(IJK)) / ep_g(IJK)) ** (10.0 / 9.0)
  A = K_so / K_go

  ! Calculate K
  K = (2.0 / (1.0 - B / A)) * ((A - 1.0) / ((1.0 - B / A) ** 2) * B / A * log(A / B) - (B - 1.0) / (1.0 - B / A) - 0.5 * (B + 1.0))

  ! Calculate ep_s (replace with your desired value)
  ep_s = 0.52

  ! Calculate k_s
  K_s(IJK) = ((beta * A + (1.0 - beta) * K) * k_go) / sqrt(1.0 - ep_s(IJK))

! conductivity of MEA solution
         K_aqMEA = calc_cond_liq(xTl, x_lrnh2, x_lh2o)

! Assign the fluid phase thermal conductivity
         K_S(IJK,M) = K_aqMEA
!      ELSE
!         K_S(IJK,M) = zero
!      ENDIF

      IF (IER /= UNDEFINED_I) THEN
         write(err_prop, '("solids phase ",I2," conductivity")') M
         WRITE(ERR_MSG,9999) trim(err_prop)
         CALL LOG_ERROR()
 9999 FORMAT('ERROR 9999: The user-defined properties routine was ',&
         'invoked for',/,A,' but this generic error',/,'message exi',&
         'sts. Either choose a different model or correct',/,'mfix/,'&
         'model/usr_properties.f')
      ENDIF
      RETURN
      END SUBROUTINE USR_PROP_Ks