File: /nfs/home/0/users/jenkins/mfix.git/model/des/diffuse_mean_fields.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     !  Subroutine: DIFFUSE_MEAN_FIELDS                                     !
4     !  Author: J.Musser                                   Date: 11-NOV-14  !
5     !                                                                      !
6     !  Purpose: Given the field variable PHI (e.g., volume fraction),      !
7     !  diffuse it across the Eulerian grid such that the Full Width at     !
8     !  Half Maximum (FWHM) equals the user specified diffusion width.      !
9     !                                                                      !
10     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
11           SUBROUTINE DIFFUSE_MEAN_FIELD(PHI, VNAME)
12     
13     ! Global Variables:
14     !---------------------------------------------------------------------//
15     ! Max bound for array sizes.
16           use geometry, only: IJKMAX2
17     ! Coefficient matrix and force vector.
18           use ambm, only: A_M, B_M
19     ! Method to solve linear system and max iterations
20           use leqsol, only: LEQ_METHOD, LEQ_IT
21     ! Preconditioner, sweep method, convergence tolerance
22           use leqsol, only: LEQ_PC, LEQ_SWEEP, LEQ_TOL
23     ! Cell-center gas velocities.
24           use tmp_array, only: DIF => ARRAY1
25     ! Size of fluid variable arrays.
26           use param, only: DIMENSION_3
27     
28     ! Module procedures
29     !---------------------------------------------------------------------//
30     ! Lock/Unlock the temp arrays to prevent double usage.
31           use tmp_array, only: LOCK_TMP_ARRAY
32           use tmp_array, only: UNLOCK_TMP_ARRAY
33     ! Routines to mange messages to user.
34           use error_manager
35     
36           IMPLICIT NONE
37     
38     ! Dummy Arguments:
39     !---------------------------------------------------------------------//
40     ! Variable to diffuse
41           DOUBLE PRECISION, INTENT(INOUT) :: PHI(DIMENSION_3)
42     ! Name of variable to diffuse
43           CHARACTER(LEN=*), INTENT(IN) :: VNAME
44     
45     ! Local Variables:
46     !---------------------------------------------------------------------//
47     ! Integer error flag
48           INTEGER :: IER
49     ! Linear equation solver method and iterations
50           INTEGER :: LEQM, LEQI
51     ! Start, stop and step size of diffusion time march
52           DOUBLE PRECISION :: DIF_TIME, DIF_STOP, DIF_DT
53     ! External function for wall time and wall time at start
54           DOUBLE PRECISION :: WALL_TIME, WALL_START
55     ! Local flag to print debug messages
56           LOGICAL, PARAMETER :: setDBG = .TRUE.
57     !......................................................................!
58     
59           IF(setDBG) THEN
60              WRITE(ERR_MSG, "(/3x,'Diffusing Variable: ',A)") VNAME
61              CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
62              WALL_START = WALL_TIME()
63           ENDIF
64     
65     ! Lock the temp arrays.
66           CALL LOCK_TMP_ARRAY
67     
68     ! Populate the diffusion coefficients
69           CALL CALC_DIF_DES(DIF, setDBG, IER)
70     
71     
72           DIF_STOP = 1.0d0
73           DIF_TIME = 0.0d0
74           DIF_DT = DIF_STOP/5.0
75     
76     ! Integrate the diffusion equation (time, space)
77           DO WHILE(DIF_TIME < DIF_STOP)
78     ! Initialize the coefficient matrix and force vector
79              CALL INIT_AB_M (A_M, B_M, IJKMAX2, 0, IER)
80     ! Calculate the coefficients
81              CALL DIF_PHI_DES(0, DIF, A_M, B_M, IER)
82     ! Apply zero-flux BC at all walls
83              CALL DIF_PHI_BC_DES(PHI, 0, A_M, B_M, IER)
84     ! Collect the center coefficient and force vector
85              CALL DIF_PHI_SOURCE_DES(PHI, 0, A_M, B_M, DIF_DT, IER)
86     ! Set the local method and iterations.
87              CALL ADJUST_LEQ(0.0d0, LEQ_IT(10), LEQ_METHOD(10), LEQI, LEQM)
88     ! Solve the linear system.
89              CALL SOLVE_LIN_EQ (VNAME, 10, PHI, A_M, B_M, 0, LEQI, LEQM, &
90                 LEQ_SWEEP(10), LEQ_TOL(10), LEQ_PC(10), IER)
91     ! Advance time.
92              DIF_TIME = DIF_TIME + DIF_DT
93           ENDDO
94     
95     ! Debugging information
96           IF(setDBG) THEN
97              WRITE(ERR_MSG, 9001) WALL_TIME() - WALL_START
98              CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
99           ENDIF
100     
101      9001 FORMAT(5x,'Wall Time: ',g11.4)
102     
103     ! Unlock the temp arrays.
104           CALL UNLOCK_TMP_ARRAY
105     
106     
107           RETURN
108           END SUBROUTINE DIFFUSE_MEAN_FIELD
109     
110     
111     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
112     !                                                                      !
113     !  Subroutine: CALC_DIF_DES                                            !
114     !  Author: J.Musser                                   Date: 11-NOV-14  !
115     !                                                                      !
116     !  Purpose: Calculate the diffusion coefficient for diffusing mean     !
117     !  fields. Presently the diffusion coefficient is constant.            !
118     !                                                                      !
119     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
120           SUBROUTINE CALC_DIF_DES(DIF, lDBG, IER)
121     
122           use param
123           use param1
124     
125           use compar, only: IJKStart3, IJKEnd3
126           use particle_filter, only: DES_DIFFUSE_WIDTH
127           use functions, only: FLUID_AT
128     
129           use error_manager
130     
131           IMPLICIT NONE
132     
133           DOUBLE PRECISION, INTENT(INOUT) :: DIF(DIMENSION_3)
134           LOGICAL, INTENT(IN) :: lDBG
135           INTEGER, INTENT(INOUT) :: IER
136     
137     ! Fluid Cell indices
138           INTEGER :: IJK
139     
140           DOUBLE PRECISION :: lDIF
141     
142           IER = 0
143     
144     ! The diffusion coefficient is set so that over one second, the
145     ! quantity diffuses such that the Full Width at Half Maximum (FWHM)
146     ! equals what the user specified as the "filter wideth."
147           lDIF = ((0.5*DES_DIFFUSE_WIDTH)**2) / &
148              (2.0*sqrt(2.0*log(2.0)))
149     
150     ! Store the diffusion coefficient in all fluid cells.
151           DO IJK = IJKStart3, IJKEnd3
152              DIF(IJK) = ZERO
153              IF(FLUID_AT(IJK)) DIF(IJK) = lDIF
154           ENDDO
155     
156     ! Information included for debugging.
157           IF(lDBG) THEN
158              WRITE(ERR_MSG, 9100) iVal(lDIF)
159              CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
160           ENDIF
161     
162      9100 FORMAT(/3x,'Diffusion Coefficient: ',A)
163     
164           RETURN
165           END SUBROUTINE CALC_DIF_DES
166