File: RELATIVE:/../../../mfix.git/model/chischeme_mod.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Purpose: Determine factors for Chi-scheme of Darwish and Moukalled  C
4     !  to ensure consistency of equation sets (e.g., species mass          C
5     !  fractions add up to 1)                                              C
6     !                                                                      C
7     !  To initiate Chi-Scheme                                              C
8     !     call set_chi( ...)                                               C
9     !  and to terminate Chi-Scheme                                         C
10     !     call unset_chi(ier)                                              C
11     !                                                                      C
12     !                                                                      C
13     !  Author: M. Syamlal                                 Date: 1-AUG-03   C
14     !                                                                      C
15     !                                                                      C
16     !  References:                                                         C
17     !  -Darwish, M. and Moukalled, F., "The Chi-shemes: a new consistent   C
18     !   high-resolution formulation based on the normalized variable       C
19     !   methodology," Comput. Methods Appl. Mech. Engrg., 192, 1711-1730   C
20     !   (2003)                                                             C
21     !                                                                      C
22     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
23     
24           MODULE ChiScheme
25     
26           IMPLICIT NONE
27     
28           DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: Chi_e
29           DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: Chi_n
30           DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: Chi_t
31     
32     
33           LOGICAL :: ChiScheme_allocated = .false.
34           LOGICAL :: Chi_flag = .false.
35     
36     
37           CONTAINS
38     
39     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
40     !                                                                      !
41     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
42     
43           SUBROUTINE Set_Chi(DISCR, PHI, Nmax, U, V, W)
44     
45           USE compar, only: ijkstart3, ijkend3
46           USE param, only: dimension_3
47           USE param1, only: large_number
48           USE sendrecv, only: send_recv
49     
50           USE error_manager, only: err_msg, init_err_msg, finl_err_msg
51           USE error_manager, only: flush_err_msg
52     
53     ! Dummy arguments
54     !---------------------------------------------------------------------//
55     ! discretization method
56           INTEGER, INTENT(IN) :: DISCR
57     ! Second dimension of Phi array
58           INTEGER, INTENT(IN) :: NMax
59     ! convected quantity
60           DOUBLE PRECISION, INTENT(IN) :: PHI(DIMENSION_3, Nmax)
61     ! Velocity components
62           DOUBLE PRECISION, INTENT(IN) :: U(DIMENSION_3)
63           DOUBLE PRECISION, INTENT(IN) :: V(DIMENSION_3)
64           DOUBLE PRECISION, INTENT(IN) :: W(DIMENSION_3)
65     
66     ! Local variables
67     !---------------------------------------------------------------------//
68           DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: Chi_e_temp
69           DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: Chi_n_temp
70           DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: Chi_t_temp
71     ! index
72           INTEGER :: IJK, N
73     !---------------------------------------------------------------------//
74     
75           if(.not.ChiScheme_allocated)then
76              Allocate( Chi_e(DIMENSION_3) , &
77                       Chi_n(DIMENSION_3) , &
78                       Chi_t(DIMENSION_3)     )
79              ChiScheme_allocated = .true.
80           endif
81     
82           if(Chi_flag)then
83     ! Error: Chi-Scheme is already active.  This routine cannot be called
84     ! again before unsetting the flag
85              CALL INIT_ERR_MSG("SET_CHI")
86              WRITE(ERR_MSG, 1102)
87      1102 FORMAT('ERROR 1102: Cannot call Set_Chi again, before Unset_chi')
88              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
89              CALL FINL_ERR_MSG
90           else
91     ! Set Chi_flag to indicate that future calls to calc_Xsi will use
92     ! the Chi-Scheme for discretization
93              Chi_e = large_number
94              Chi_n = large_number
95              Chi_t = large_number
96              Chi_flag = .true.
97           Endif
98     
99           Allocate( Chi_e_temp(DIMENSION_3) , &
100                     Chi_n_temp(DIMENSION_3) , &
101                     Chi_t_temp(DIMENSION_3)  )
102     
103     ! Start Chi calculations
104           DO N = 1, Nmax
105              CALL CALC_CHI(DISCR, PHI(1,N), U, V, W, Chi_e_temp, &
106                            Chi_n_temp, Chi_t_temp)
107     
108     !!!$omp    parallel do private(IJK)
109              DO IJK = ijkstart3, ijkend3
110                 Chi_e(IJK) = MIN(Chi_e(IJK), Chi_e_temp(IJK))
111                 Chi_n(IJK) = MIN(Chi_n(IJK), Chi_n_temp(IJK))
112                 Chi_t(IJK) = MIN(Chi_t(IJK), Chi_t_temp(IJK))
113              ENDDO
114           ENDDO
115     
116           call send_recv(CHI_E,2)
117           call send_recv(CHI_N,2)
118           call send_recv(CHI_T,2)
119     
120           Deallocate( Chi_e_temp , &
121                       Chi_n_temp , &
122                       Chi_t_temp  )
123     
124     
125           RETURN
126           END SUBROUTINE Set_Chi
127     
128     
129     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
130     !                                                                      C
131     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
132           SUBROUTINE Unset_Chi()
133           IMPLICIT NONE
134           Chi_flag = .false.
135           RETURN
136           END SUBROUTINE Unset_Chi
137     
138     
139     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
140     !                                                                      C
141     !  Subroutine: CALC_CHI                                                C
142     !  Purpose: Determine CHI factors for higher order discretization.     C
143     !                                                                      C
144     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
145           SUBROUTINE CALC_CHI(DISCR, PHI, U, V, W, CHI_E, CHI_N, CHI_T)
146     
147     
148     ! Modules
149     !---------------------------------------------------------------------//
150           USE compar, only: ijkstart3, ijkend3
151           USE discretization, only: chi4smart, chi4muscl, phi_c_of
152           USE functions, only: wall_at
153           USE functions, only: east_of, west_of, north_of, south_of
154           USE functions, only: top_of, bottom_of
155           USE geometry, only: do_k
156           USE param, only: dimension_3
157           USE param1, only: zero
158           USE run, only: shear
159           USE sendrecv, only: send_recv
160           USE error_manager, only: err_msg, init_err_msg, finl_err_msg
161           USE error_manager, only: ival, flush_err_msg
162     
163           IMPLICIT NONE
164     
165     ! Dummy arguments
166     !---------------------------------------------------------------------//
167     ! discretization method
168           INTEGER, INTENT(IN) :: DISCR
169     ! convected quantity
170           DOUBLE PRECISION, INTENT(IN) :: PHI(DIMENSION_3)
171     ! Velocity components
172           DOUBLE PRECISION, INTENT(IN) :: U(DIMENSION_3)
173           DOUBLE PRECISION, INTENT(IN) :: V(DIMENSION_3)
174           DOUBLE PRECISION, INTENT(IN) :: W(DIMENSION_3)
175     ! Convection weighting factors
176           DOUBLE PRECISION, INTENT(INOUT) :: CHI_e(DIMENSION_3)
177           DOUBLE PRECISION, INTENT(INOUT) :: CHI_n(DIMENSION_3)
178           DOUBLE PRECISION, INTENT(INOUT) :: CHI_t(DIMENSION_3)
179     
180     ! Local variables
181     !---------------------------------------------------------------------//
182     ! Indices
183           INTEGER :: IJK, IJKC, IJKD, IJKU
184           DOUBLE PRECISION :: PHI_C
185     !---------------------------------------------------------------------//
186     
187     ! calculate CHI_E,CHI_N,CHI_T when periodic shear BCs are used
188           IF (SHEAR) THEN
189     ! this needs implementation...
190     ! note mfix will error before this is hit
191     !         call CXS(incr, DISCR, U, V, W, PHI, CHI_E, CHI_N, CHI_T)
192           ELSE
193     
194     
195           SELECT CASE (DISCR)
196     
197     
198           CASE (:1)                                  !first order upwinding
199     !!!$omp    parallel do private(IJK)
200              DO IJK = ijkstart3, ijkend3
201                 CHI_E(IJK) = ZERO
202                 CHI_N(IJK) = ZERO
203                 IF (DO_K) CHI_T(IJK) = ZERO
204              ENDDO
205     
206     
207     !      CASE (2)                                   !Superbee
208     
209     
210           CASE (3)                                   !SMART
211     !!!$omp    parallel do private(IJK, IJKC,IJKD,IJKU, PHI_C)
212              DO IJK = ijkstart3, ijkend3
213                 IF (.NOT.WALL_AT(IJK)) THEN ! no need to do these calculations for walls
214                    IF (U(IJK) >= ZERO) THEN
215                       IJKC = IJK
216                       IJKD = EAST_OF(IJK)
217                       IJKU = WEST_OF(IJKC)
218                    ELSE
219                       IJKC = EAST_OF(IJK)
220                       IJKD = IJK
221                       IJKU = EAST_OF(IJKC)
222                    ENDIF
223                    PHI_C = PHI_C_OF(PHI(IJKU),PHI(IJKC),PHI(IJKD))
224                    CHI_E(IJK) = CHI4SMART(PHI_C, PHI(IJKU), PHI(IJKC), PHI(IJKD))
225     
226                    IF (V(IJK) >= ZERO) THEN
227                       IJKC = IJK
228                       IJKD = NORTH_OF(IJK)
229                       IJKU = SOUTH_OF(IJKC)
230                    ELSE
231                       IJKC = NORTH_OF(IJK)
232                       IJKD = IJK
233                       IJKU = NORTH_OF(IJKC)
234                    ENDIF
235                    PHI_C = PHI_C_OF(PHI(IJKU),PHI(IJKC),PHI(IJKD))
236                    CHI_N(IJK) = CHI4SMART(PHI_C,PHI(IJKU),PHI(IJKC),PHI(IJKD))
237     
238                    IF (DO_K) THEN
239                       IF (W(IJK) >= ZERO) THEN
240                          IJKC = IJK
241                          IJKD = TOP_OF(IJK)
242                          IJKU = BOTTOM_OF(IJKC)
243                       ELSE
244                          IJKC = TOP_OF(IJK)
245                          IJKD = IJK
246                          IJKU = TOP_OF(IJKC)
247                       ENDIF
248                       PHI_C = PHI_C_OF(PHI(IJKU),PHI(IJKC),PHI(IJKD))
249                       CHI_T(IJK) = CHI4SMART(PHI_C,PHI(IJKU),PHI(IJKC),PHI(IJKD))
250                    ENDIF
251                 ELSE
252                    CHI_E(IJK) = ZERO
253                    CHI_N(IJK) = ZERO
254                    CHI_T(IJK) = ZERO
255                 ENDIF  ! endif (.not. wall_at)
256              ENDDO   ! end do ijk
257     
258     
259     !      CASE (4)                                   !ULTRA-QUICK
260     
261     
262     !      CASE (5)                                   !QUICKEST
263     
264     
265           CASE (6)                                   !MUSCL
266     
267     !!!$omp    parallel do private(IJK, IJKC,IJKD,IJKU, PHI_C )
268              DO IJK = ijkstart3, ijkend3
269                 IF (.NOT.WALL_AT(IJK)) THEN ! no need to do these calculations for walls
270                    IF (U(IJK) >= ZERO) THEN
271                       IJKC = IJK
272                       IJKD = EAST_OF(IJK)
273                       IJKU = WEST_OF(IJKC)
274                    ELSE
275                       IJKC = EAST_OF(IJK)
276                       IJKD = IJK
277                       IJKU = EAST_OF(IJKC)
278                    ENDIF
279                    PHI_C = PHI_C_OF(PHI(IJKU),PHI(IJKC),PHI(IJKD))
280                    CHI_E(IJK) = CHI4MUSCL(PHI_C,PHI(IJKU),PHI(IJKC),PHI(IJKD))
281     
282                    IF (V(IJK) >= ZERO) THEN
283                       IJKC = IJK
284                       IJKD = NORTH_OF(IJK)
285                       IJKU = SOUTH_OF(IJKC)
286                    ELSE
287                       IJKC = NORTH_OF(IJK)
288                       IJKD = IJK
289                       IJKU = NORTH_OF(IJKC)
290                    ENDIF
291                    PHI_C = PHI_C_OF(PHI(IJKU),PHI(IJKC),PHI(IJKD))
292                    CHI_N(IJK) = CHI4MUSCL(PHI_C,PHI(IJKU),PHI(IJKC),PHI(IJKD))
293     
294                    IF (DO_K) THEN
295                       IF (W(IJK) >= ZERO) THEN
296                          IJKC = IJK
297                          IJKD = TOP_OF(IJK)
298                          IJKU = BOTTOM_OF(IJKC)
299                       ELSE
300                          IJKC = TOP_OF(IJK)
301                          IJKD = IJK
302                          IJKU = TOP_OF(IJKC)
303                       ENDIF
304                       PHI_C = PHI_C_OF(PHI(IJKU),PHI(IJKC),PHI(IJKD))
305                       CHI_T(IJK) = CHI4MUSCL(PHI_C,PHI(IJKU),PHI(IJKC),PHI(IJKD))
306                    ENDIF
307                 ELSE
308                    CHI_E(IJK) = ZERO
309                    CHI_N(IJK) = ZERO
310                    CHI_T(IJK) = ZERO
311                 ENDIF   ! endif (.not.wall_at)
312              ENDDO   ! end do ijk
313     
314     
315     !      CASE (7)                                   !Van Leer
316     
317     
318     !      CASE (8)                                   !Minmod
319     
320     
321           CASE DEFAULT                               !Error
322     ! should never hit this
323               CALL INIT_ERR_MSG("CALC_CHI")
324               WRITE(ERR_MSG, 1103) IVAL(DISCR)
325      1103 FORMAT('ERROR 1103: Invalid DISCRETIZE= ', A,' when using '&
326              'chi_scheme',/'The check_data routines should have already ',&
327              'caught this error and ',/,'pevented the simulation from ',&
328              'running. Please notify the MFIX ',/,'MFIX developers.')
329               CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
330               CALL FINL_ERR_MSG
331           END SELECT
332     
333           ENDIF
334     
335           call send_recv(CHI_E,2)
336           call send_recv(CHI_N,2)
337           call send_recv(CHI_T,2)
338     
339           RETURN
340           END SUBROUTINE CALC_CHI
341     
342           END MODULE CHIScheme
343