diffuse_mean_fields.f

Go to the documentation of this file.

!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
!                                                                      !
!  Subroutine: DIFFUSE_MEAN_FIELDS                                     !
!  Author: J.Musser                                   Date: 11-NOV-14  !
!                                                                      !
!  Purpose: Given the field variable PHI (e.g., volume fraction),      !
!  diffuse it across the Eulerian grid such that the Full Width at     !
!  Half Maximum (FWHM) equals the user specified diffusion width.      !
!                                                                      !
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
      SUBROUTINE diffuse_mean_field(PHI, VNAME)

! Global Variables:
!---------------------------------------------------------------------//
! Max bound for array sizes.
      use geometry, only: ijkmax2
! Coefficient matrix and force vector.
      use ambm, only: a_m, b_m
! Method to solve linear system and max iterations
      use leqsol, only: leq_method, leq_it
! Preconditioner, sweep method, convergence tolerance
      use leqsol, only: leq_pc, leq_sweep, leq_tol
! Size of fluid variable arrays.
      use param, only: dimension_3

! Module procedures
!---------------------------------------------------------------------//
! Routines to manage messages to user.
      use error_manager
      use machine, only: wall_time

      IMPLICIT NONE

! Dummy Arguments:
!---------------------------------------------------------------------//
! Variable to diffuse
      DOUBLE PRECISION, INTENT(INOUT) :: PHI(dimension_3)
! Name of variable to diffuse
      CHARACTER(LEN=*), INTENT(IN) :: VNAME

! Local Variables:
!---------------------------------------------------------------------//
! Integer error flag
      INTEGER :: IER
! Linear equation solver method and iterations
      INTEGER :: LEQM, LEQI
! Start, stop and step size of diffusion time march
      DOUBLE PRECISION :: DIF_TIME, DIF_STOP, DIF_DT
! wall time at start
      DOUBLE PRECISION :: WALL_START
! Local flag to print debug messages
      LOGICAL, PARAMETER :: setDBG = .true.
      DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: DIF
!......................................................................!

      IF(setdbg) THEN
         WRITE(err_msg, "(/3x,'Diffusing Variable: ',A)") vname
         CALL flush_err_msg(header=.false., footer=.false.)
         wall_start = wall_time()
      ENDIF

! Populate the diffusion coefficients
      CALL calc_dif_des(dif, setdbg, ier)

      dif_stop = 1.0d0
      dif_time = 0.0d0
      dif_dt = dif_stop/5.0

! Integrate the diffusion equation (time, space)
      DO WHILE(dif_time < dif_stop)
! Initialize the coefficient matrix and force vector
         CALL init_ab_m (a_m, b_m, ijkmax2, 0)
! Calculate the coefficients
         CALL dif_phi_des(0, dif, a_m, b_m)
! Apply zero-flux BC at all walls
         CALL dif_phi_bc_des(phi, 0, a_m, b_m)
! Collect the center coefficient and force vector
         CALL dif_phi_source_des(phi, 0, a_m, b_m, dif_dt)
! Set the local method and iterations.
         CALL adjust_leq(0.0d0, leq_it(10), leq_method(10), leqi, leqm)
! Solve the linear system.
         CALL solve_lin_eq (vname, 10, phi, a_m, b_m, 0, leqi, leqm, &
            leq_sweep(10), leq_tol(10), leq_pc(10), ier)
! Advance time.
         dif_time = dif_time + dif_dt
      ENDDO

! Debugging information
      IF(setdbg) THEN
         WRITE(err_msg, 9001) wall_time() - wall_start
         CALL flush_err_msg(header=.false., footer=.false.)
      ENDIF

 9001 FORMAT(5x,'Wall Time: ',g11.4)

      RETURN
      END SUBROUTINE diffuse_mean_field

!vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
!                                                                      !
!  Subroutine: CALC_DIF_DES                                            !
!  Author: J.Musser                                   Date: 11-NOV-14  !
!                                                                      !
!  Purpose: Calculate the diffusion coefficient for diffusing mean     !
!  fields. Presently the diffusion coefficient is constant.            !
!                                                                      !
!^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
      SUBROUTINE calc_dif_des(DIF, lDBG, IER)

      use param
      use param1

      use compar, only: ijkstart3, ijkend3
      use particle_filter, only: des_diffuse_width
      use functions, only: fluid_at

      use error_manager

      IMPLICIT NONE

      DOUBLE PRECISION, INTENT(INOUT) :: DIF(dimension_3)
      LOGICAL, INTENT(IN) :: lDBG
      INTEGER, INTENT(INOUT) :: IER

! Fluid Cell indices
      INTEGER :: IJK

      DOUBLE PRECISION :: lDIF

      ier = 0

! The diffusion coefficient is set so that over one second, the
! quantity diffuses such that the Full Width at Half Maximum (FWHM)
! equals what the user specified as the "filter wideth."
      ldif = ((0.5*des_diffuse_width)**2) / &
         (2.0*sqrt(2.0*log(2.0)))

! Store the diffusion coefficient in all fluid cells.
      DO ijk = ijkstart3, ijkend3
         dif(ijk) = zero
         IF(fluid_at(ijk)) dif(ijk) = ldif
      ENDDO

! Information included for debugging.
      IF(ldbg) THEN
         WRITE(err_msg, 9100) ival(ldif)
         CALL flush_err_msg(header=.false., footer=.false.)
      ENDIF

 9100 FORMAT(/3x,'Diffusion Coefficient: ',a)

      RETURN
      END SUBROUTINE calc_dif_des