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