usr_prop_mod.f

Go to the documentation of this file.

! or put the actual routine calc_usr_source here? move logic variable
!  usr_source into run or something...
      module usr_prop
      use param, only: dim_m

! If .TRUE. call user-defined physical properties
      LOGICAL :: usr_rog
      LOGICAL :: usr_cpg 
      LOGICAL :: usr_ros(dim_m)
      LOGICAL :: usr_cps(dim_m)

! If .TRUE. call user_defined transport properties
      LOGICAL :: usr_mug
      LOGICAL :: usr_kg
      LOGICAL :: usr_difg
      LOGICAL :: usr_mus(dim_m)
      LOGICAL :: usr_ks(dim_m)
      LOGICAL :: usr_difs(dim_m)

! Interphase momentum exchange coefficients (F_GS, F_SS)
      LOGICAL :: usr_fgs(dim_m)
      LOGICAL :: usr_fss( (dim_m*(dim_m-1)/2)+1 )
! Interphase heat transfer coefficient (GAMA)
      LOGICAL :: usr_gama(dim_m)

      enum, bind(c)
         enumerator :: gas_density, solids_density
         enumerator :: gas_specificheat, solids_specificheat
         enumerator :: gas_viscosity, solids_viscosity
         enumerator :: gas_conductivity, solids_conductivity
         enumerator :: gas_diffusivity, solids_diffusivity 
         enumerator :: gassolids_drag, solidssolids_drag
         enumerator :: gassolids_heattransfer
         enumerator :: blank
         enumerator :: dummy
      END enum


      CONTAINS
!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
!                                                                      C
!  Subroutine: CALC_USR_PROP                                           C 
!  Purpose: Driver routine to calculate user defined physical,         C
!  transport, exchange terms that appear in the governing equations.   C
!                                                                      C
!                                                                      C
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
      SUBROUTINE calc_usr_prop(lprop, lM, lL, lerr)

! Modules
!-----------------------------------------------
      use compar, only: ijkstart3, ijkend3, mype, numpes
      use functions, only: fluid_at, wall_at, fluidorp_flow_at
      use output, only: report_neg_density
      use output, only: report_neg_specificheat
      use param1, only: undefined_i, zero
      use physprop, only: nmax
      use physprop, only: k_g, dif_g, c_pg, mu_g
      use physprop, only: k_s, dif_s, c_ps
      use fldvar, only: ro_g, rop_g, ep_g
      use fldvar, only: ro_s, p_s
      use visc_g, only: mu_gt, lambda_gt
      use visc_s, only: mu_s, lambda_s
      use utilities, only: mfix_isnan

      use error_manager

      IMPLICIT NONE
! Dummy arguments
!-----------------------------------------------
! User properties available
!     ROg, ROs (gas and solids density)
!     CPg, CPs (gas and solids specific heat)
!     Mug, Mus (gas and solids viscosity)
!     Kg, Ks (gas and solids conductivity)
!     Difg, Difs (gas and solids diffusivity)
!     Fgs, Fss (gas-solids and solids-solids drag)
!     Gama (gas-solids heat transfere coefficient)
      INTEGER, INTENT(IN) :: lprop

! Phase index 
      INTEGER, OPTIONAL, INTENT(IN) :: lM

! Phase index 
      INTEGER, OPTIONAL, INTENT(IN) :: lL

! Error index
! Arrays for storing errors:
! 100 - Negative gas phase density
! 101 - Negative solids phase density
! 10x - Unclassified
      INTEGER, OPTIONAL, INTENT(INOUT) :: lerr(0:numpes-1)

! Local variables
!-----------------------------------------------
! indices
      INTEGER :: IJK
! tmp indices
      INTEGER :: M, L, N
! Flag to write log header
      LOGICAL :: wHeader

!-----------------------------------------------

! set local values for phase index
      IF (.NOT.present(lm)) THEN
         m = undefined_i
      ELSE
         m = lm
      ENDIF

! set local values for phase index
      IF (.NOT.present(ll)) THEN
         l = undefined_i
      ELSE
         l = ll
      ENDIF
         
      n = undefined_i

      SELECT CASE(lprop)
      
! Gas physical properties: density
      CASE(gas_density)
         wheader = .true.
         DO ijk=ijkstart3,ijkend3
            IF (wall_at(ijk)) cycle
!            CALL USR_PROPERTIES(lprop, IJK, 0, N)
            CALL usr_prop_rog(ijk)
            rop_g(ijk) = ep_g(ijk)*ro_g(ijk)
! report errors
            IF(ro_g(ijk) < zero .or. mfix_isnan(ro_g(ijk))) THEN
               lerr(mype) = 100  ! return appropritae error flag
               IF(report_neg_density) CALL rogerr_log(ijk, wheader)
            ENDIF

         ENDDO

! Gas physical properties: specific heat
      CASE(gas_specificheat)
         DO ijk=ijkstart3,ijkend3
            IF (wall_at(ijk)) cycle
!            CALL USR_PROPERTIES(lprop, IJK, 0, N)
            CALL usr_prop_cpg(ijk)
! report errors.
            IF(c_pg(ijk) <= zero) THEN
!               lErr(myPE) = 103
               IF(report_neg_specificheat) CALL cpgerr_log(ijk, wheader)
            ENDIF
         ENDDO

! Gas transport properties: viscosity
      CASE(gas_viscosity)
         DO ijk=ijkstart3,ijkend3
            IF (fluid_at(ijk)) THEN
!               CALL USR_PROPERTIES(lprop, IJK, 0, N)
               CALL usr_prop_mug(ijk)
               mu_gt(ijk) = mu_g(ijk)
               lambda_gt(ijk) = -2./3.*mu_gt(ijk)
            ELSE ! probably unnecessary
               mu_g(ijk) = zero
               mu_gt(ijk) = zero
               lambda_gt(ijk) = zero
            ENDIF
         ENDDO

! Gas transport properties: conductivity
      CASE(gas_conductivity)
         DO ijk=ijkstart3,ijkend3
            IF (fluid_at(ijk)) THEN
!               CALL USR_PROPERTIES(lprop, IJK, 0, N)
               CALL usr_prop_kg(ijk)
            ELSE ! probably unnecessary
               k_g(ijk) = zero
            ENDIF
         ENDDO

! Gas transport properties: diffusivity
      CASE(gas_diffusivity)
         DO n=1,nmax(0)
            DO ijk=ijkstart3,ijkend3
               IF (fluid_at(ijk)) THEN
!                  CALL USR_PROPERTIES(lprop, IJK, 0, N)
                  CALL usr_prop_difg(ijk,n)
               ELSE ! probably unnecessary
                  dif_g(ijk,n) = zero
               ENDIF
            ENDDO
         ENDDO
 
! Solids physical properties: density
      CASE(solids_density)
         wheader = .true.
         DO ijk=ijkstart3,ijkend3
            IF (wall_at(ijk)) cycle
!            CALL USR_PROPERTIES(lprop, IJK, M, N)
            CALL usr_prop_ros(ijk,m)
! report errors.
            IF(ro_s(ijk,m) <= zero) THEN
               lerr(mype) = 101
               IF(report_neg_density) CALL roserr_log(ijk, m, wheader)
            ENDIF

         ENDDO

! Solids physical properties: specific heat
      CASE(solids_specificheat)
         DO ijk=ijkstart3,ijkend3
            IF (wall_at(ijk)) cycle
!            CALL USR_PROPERTIES(lprop, IJK, M, N)
            CALL usr_prop_cps(ijk,m)
! report errors.
            IF(c_ps(ijk,m) <= zero) THEN
!               lErr(myPE) = 104
               IF(report_neg_specificheat) CALL cpserr_log(ijk, m, wheader)
            ENDIF
         ENDDO

! Solids transport properties: viscosity
      CASE(solids_viscosity)
         DO ijk=ijkstart3,ijkend3
            IF (fluid_at(ijk)) THEN
!               CALL USR_PROPERTIES(lprop, IJK, M, N)
               CALL usr_prop_mus(ijk,m)
            ELSE ! probably unnecessary
               mu_s(ijk,m) = zero
               lambda_s(ijk,m) = zero
            ENDIF
         ENDDO

! Solids transport properties: conductivity
      CASE(solids_conductivity)
         DO ijk=ijkstart3,ijkend3
            IF (fluid_at(ijk)) THEN
!               CALL USR_PROPERTIES(lprop, IJK, M, N)
               CALL usr_prop_ks(ijk,m)
            ELSE ! probably unnecessary
               k_s(ijk,m) = zero
            ENDIF
         ENDDO

! Solids transport properties: diffusivity
      CASE(solids_diffusivity)
         DO n=1,nmax(m)
            DO ijk=ijkstart3,ijkend3
               IF (fluid_at(ijk)) THEN
!                  CALL USR_PROPERTIES(lprop, IJK, M, N)
                  CALL usr_prop_difs(ijk,m,n)
               ELSE ! probably unnecessary
                  dif_s(ijk,m,n) = zero
               ENDIF
            ENDDO
         ENDDO

! Gas-solids exchange: Interphase heat transfer coefficient (GAMA)
      CASE(gassolids_heattransfer)
         DO ijk=ijkstart3,ijkend3
            IF (fluid_at(ijk)) THEN
!               CALL USR_PROPERTIES(lprop, IJK, M, N)
               CALL usr_prop_gama(ijk,m)
            ENDIF
         ENDDO

! Gas-solids exchange: interphase drag
      CASE(gassolids_drag)
         DO ijk=ijkstart3,ijkend3
            IF (fluidorp_flow_at(ijk)) THEN
               CALL usr_properties(lprop, ijk, m, n)
! this hook is not accessed or functioning.
! use drag_type = usr_drag to access this entity
            ENDIF
         ENDDO

! Solids-solids exchange: interphase drag
      CASE(solidssolids_drag)
         DO ijk=ijkstart3,ijkend3
            IF (wall_at(ijk)) cycle
!            CALL USR_PROPERTIES(lprop, IJK, M, L)
            CALL usr_prop_fss(ijk,l,m)
         ENDDO

! error out
      CASE DEFAULT
! should never hit this
! Initialize the error manager.
         CALL init_err_msg("CALC_USR_PROP")
         WRITE(err_msg, 1001) ival(lprop)
         CALL flush_err_msg(abort=.true.)
 1001 FORMAT('Error 1101: Unknown Property= ', a)
         CALL finl_err_msg
      END SELECT   ! end selection of user property 
 
      RETURN
      END SUBROUTINE calc_usr_prop

      end module usr_prop