check_ic_common_discrete.f

Go to the documentation of this file.

!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
!                                                                      !
!  Subroutine: CHECK_IC_COMMON_DISCRETE                                !
!  Author:   R.Garg                                   Date: 11-Mar-14  !
!                                                                      !
!  Purpose: check the initial conditions input section common to both  !
!           DEM and MPPIC models                                       !
!     - ensure the first IC is defined over the entire domain with     !
!        ep_g = 1 when more than one IC has solids                     !
!     - ensure the ICs are non-overlapping                             !
!                                                                      !
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
      SUBROUTINE check_ic_common_discrete

! Modules
!---------------------------------------------------------------------//
! Runtime Flag: Generate initial particle configuration.
      USE discretelement, only : gener_part_config
! Simulation dimension (2D/3D)
      USE discretelement, only: dimn
! Number of DEM solids phases.
      USE discretelement, only: des_mmax
! direction wise spans of the domain and grid spacing in each direction
      Use geometry, only: xlength, ylength, zlength
! Flag indicating that the IC region is defined.
      USE ic, only: ic_defined
! IC Region gas volume fraction.
      USE ic, only: ic_ep_g
! IC Region solid volume fraction.
      USE ic, only: ic_ep_s
! IC Region gas volume fraction.
      USE ic, only: ic_theta_m
!
      USE ic, only: ic_x_w, ic_x_e, ic_y_s, ic_y_n, ic_z_b, ic_z_t
! Maximum number of IC regions
      USE param, only: dimension_ic
!
      USE param1, only: undefined, undefined_i, zero, one
!
      use physprop, only: mmax

! Use the error manager for posting error messages.
      use error_manager
      implicit none

! Local variables 
!---------------------------------------------------------------------//
      INTEGER :: ICV, ICV2, M, IDIM
      INTEGER :: COUNT_IC, COUNT_IC_WITH_SOLS
      INTEGER :: FIRST_DEF_IC
      DOUBLE PRECISION :: IC_ORIG(3), IC_END(3), IC2_ORIG(3) , IC2_END(3)
      DOUBLE PRECISION :: IC_MIN, IC_MAX, IC2_MIN, IC2_MAX , TOL_IC_REG
      LOGICAL :: SEP_AXIS, first_ic_ok
!---------------------------------------------------------------------//

      IF (.NOT.gener_part_config) RETURN

! Initialize the error manager.
      CALL init_err_msg("CHECK_IC_COMMON_DISCRETE")

! First check if multiple IC regions are defined for non-zero solids volume
! fraction, then check if the first IC is specified over the whole domain with IC_EP_g = 1

      !total count of defined ICs
      count_ic           = 0
      !total count of defined IC's with solids
      count_ic_with_sols = 0
      first_def_ic = undefined_i
      DO icv = 1, dimension_ic

         IF (ic_defined(icv)) THEN
            count_ic = count_ic + 1
            first_def_ic = min(first_def_ic, icv)

            IF(ic_ep_g(icv).LT.one) count_ic_with_sols &
            = count_ic_with_sols  + 1

         ENDIF ! if(ic_defined(icv))
      end DO

      IF(count_ic_with_sols >= 1 .AND. &
         count_ic > count_ic_with_sols+1) THEN

! If the number of IC's with solids is greater than one, make sure the
! first IC spans the entire domain with voidage of one. This ensures
! that the entire domain has valid ICs defined.
         icv = first_def_ic
         first_ic_ok = .false.
         IF(ic_ep_g(icv).EQ.one &
           .AND.ic_x_w(icv).LE.zero.AND.ic_x_e(icv).GE.xlength         &
           .AND.ic_y_s(icv).LE.zero.AND.ic_y_n(icv).GE.ylength)        &
            first_ic_ok = .true.

         IF (first_ic_ok .AND. ic_z_b(icv) <= zero .AND. &
            ic_z_t(icv) >= zlength) first_ic_ok = .true.

         IF(.NOT.first_ic_ok) THEN
            WRITE(err_msg, 1003)
            CALL flush_err_msg(abort=.true.)
         ENDIF

 1003 FORMAT(' Error 1003: Particle seeding with more than one IC ',   &
         'region requires',/'that IC 1 span the entire domain and ',   &
         'have IC_EP_g(1) = 1.0.',/'Please correct the mfix.dat file.')

      ENDIF

! Check if the ICs are non-overlapping.
      tol_ic_reg  = 1e-04
      icvloop : DO icv = 1, dimension_ic

         IF(.NOT.ic_defined(icv)) cycle icvloop
         IF(ic_ep_g(icv) == 1.d0) cycle icvloop
         ic_orig(1) = ic_x_w(icv)
         ic_orig(2) = ic_y_s(icv)
         ic_orig(3) = ic_z_b(icv)
         ic_end(1)  = ic_x_e(icv)
         ic_end(2)  = ic_y_n(icv)
         ic_end(3)  = ic_z_t(icv)
         icvtwoloop : DO icv2 = icv+1, dimension_ic

            IF(.NOT.ic_defined(icv2)) cycle icvtwoloop
            IF(ic_ep_g(icv2) == 1.0d0) cycle icvtwoloop

            ic2_orig(1) = ic_x_w(icv2)
            ic2_orig(2) = ic_y_s(icv2)
            ic2_orig(3) = ic_z_b(icv2)
            ic2_end(1)  = ic_x_e(icv2)
            ic2_end(2)  = ic_y_n(icv2)
            ic2_end(3)  = ic_z_t(icv2)

            sep_axis  = .false.
            DO idim = 1, dimn

               ic_min = ic_orig(idim)
               ic_max = ic_end(idim)
               ic2_min = ic2_orig(idim)
               ic2_max = ic2_end(idim)

! Check for separating axis. If the separating axis exists, then the IC
! regions can't overlap generally equality implies lack of sep_axis,
! and thus, overlapping. However, doing so will flag all IC's as
! overlapping since IC's have to share common edges. So here the
! equality is considered as existence of a separating axis, and hence,
! no overlap equality is also considered as separating axis which is
               if ((ic_min .ge. ic2_max)  .or. (ic_max .le. ic2_min) ) then
                  sep_axis = .true.
                  exit
               endif
            end DO

! Implies the IC regions could not find a separating axis and are
! therefore overlapping.
            IF(.NOT.sep_axis) THEN
               WRITE(err_msg, 1004) icv, icv2
               CALL flush_err_msg(abort=.true.)
            ENDIF

 1004 FORMAT('Error 1004: Overlapping IC regions with nonzero solids ',&
         'volume',/'fraction detected. This is not supported for ',    &
         'discrete solids.',2/'Overlapping ICs: ',2(2x,i4),2/,         &
         'Please correct the mfix.dat file.')

         end DO icvtwoloop
      end DO icvloop



! Check if IC_theta_M is specified for solids phases wherever IC_EP_g lt 1
      DO icv = 1, dimension_ic

         IF (ic_defined(icv).and.ic_ep_g(icv).LT.one) THEN
            DO m = mmax+1, des_mmax+mmax
               IF(ic_theta_m(icv,m)==undefined) THEN
                  IF(ic_ep_s(icv,m).gt.zero) THEN
                     WRITE(err_msg, 1000) trim(ivar('IC_THETA_M',icv,m))
                     CALL flush_err_msg(abort=.true.)
                  ELSE
                     ic_theta_m(icv,m) = zero
                  ENDIF
               ENDIF
            ENDDO
         ENDIF
      ENDDO

      CALL finl_err_msg

      RETURN

 1000 FORMAT('Error 1000: Required input not specified: ',a,/'Please ',&
         'correct the mfix.dat file.')

      END SUBROUTINE check_ic_common_discrete