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

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: CALC_MFLUX(IER)                                        C
4     !  Purpose: Calculate the convection fluxes. Master routine.           C
5     !                                                                      C
6     !  Author: M. Syamlal                                 Date: 31-MAY-05  C
7     !  Reviewer:                                          Date:            C
8     !                                                                      C
9     !                                                                      C
10     !  Literature/Document References:                                     C
11     !                                                                      C
12     !  Variables referenced:                                               C
13     !  Variables modified:                                                 C
14     !                                                                      C
15     !  Local variables:                                                    C
16     !                                                                      C
17     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
18     !
19           SUBROUTINE CALC_MFLUX(IER)
20     !
21     !-----------------------------------------------
22     !   M o d u l e s
23     !-----------------------------------------------
24           USE param
25           USE param1
26           USE fldvar
27           USE mflux
28           USE physprop
29           USE run
30           IMPLICIT NONE
31     !-----------------------------------------------
32     !   G l o b a l   P a r a m e t e r s
33     !-----------------------------------------------
34     !-----------------------------------------------
35     !   D u m m y   A r g u m e n t s
36     !-----------------------------------------------
37     !
38     !
39     !
40     !                      Error index
41           INTEGER          IER
42     !
43     !                      solids phase index
44           INTEGER          M
45     !
46     !
47     !
48     !
49           IF(.NOT.Added_Mass) THEN
50              CALL CALC_MFLUX0 (U_g, V_g, W_g, ROP_gE, ROP_gN, ROP_gT, &
51                                Flux_gE, Flux_gN, Flux_gT)
52              DO M = 1, MMAX
53                CALL CALC_MFLUX0 (U_s(1, M), V_s(1, M), W_s(1, M), &
54                                  ROP_sE(1, M), ROP_sN(1, M), ROP_sT(1, M), &
55                                  Flux_sE(1, M), Flux_sN(1, M), Flux_sT(1, M))
56              ENDDO
57           ELSE
58     
59              CALL CALC_MFLUX1 (U_g, V_g, W_g, ROP_gE, ROP_gN, ROP_gT, &
60                                ROP_sE(1,M_AM), ROP_sN(1,M_AM), ROP_sT(1,M_AM), &
61                                Flux_gE, Flux_gN, Flux_gT, M_AM, IER)
62     
63              CALL CALC_MFLUX0 (U_g, V_g, W_g, ROP_gE, ROP_gN, ROP_gT, &
64                                Flux_gSE, Flux_gSN, Flux_gST)
65              DO M = 1, MMAX ! New fluxes are defined only for M = M_am where virtual mass force is added.
66                IF(M==M_AM) THEN
67                   CALL CALC_MFLUX1 (U_s(1, M), V_s(1, M), W_s(1, M), &
68                                  ROP_sE(1, M), ROP_sN(1, M), ROP_sT(1, M), &
69                                  ROP_gE, ROP_gN, ROP_gT, &
70                                  Flux_sE(1, M), Flux_sN(1, M), Flux_sT(1, M), M_AM, IER)
71                   CALL CALC_MFLUX0 (U_s(1, M), V_s(1, M), W_s(1, M), &
72                                  ROP_sE(1, M), ROP_sN(1, M), ROP_sT(1, M), &
73                                  Flux_sSE, Flux_sSN, Flux_sST)
74                ELSE
75                   CALL CALC_MFLUX0 (U_s(1, M), V_s(1, M), W_s(1, M), &
76                                  ROP_sE(1, M), ROP_sN(1, M), ROP_sT(1, M), &
77                                  Flux_sE(1, M), Flux_sN(1, M), Flux_sT(1, M))
78                ENDIF
79              ENDDO
80           ENDIF
81     
82           RETURN
83           END SUBROUTINE CALC_MFLUX
84     !
85     !
86     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
87     !                                                                      C
88     !                                                                      C
89     !  Module name: CALC_MFLUX0(U, V, W, ROP_E, ROP_N, ROP_T, Flux_E,      C
90     !                                         Flux_N, Flux_T, IER)         C
91     !  Purpose: Calculate the convection fluxes.                           C
92     !                                                                      C
93     !  Author: M. Syamlal                                 Date: 31-MAY-05  C
94     !  Reviewer:                                          Date:            C
95     !                                                                      C
96     !                                                                      C
97     !  Literature/Document References:                                     C
98     !                                                                      C
99     !  Variables referenced:                                               C
100     !  Variables modified:                                                 C
101     !                                                                      C
102     !  Local variables:                                                    C
103     !                                                                      C
104     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
105     !
106           SUBROUTINE CALC_MFLUX0(U, V, W, ROP_E, ROP_N, ROP_T, Flux_E, Flux_N, Flux_T)
107     
108     !-----------------------------------------------
109     !   M o d u l e s
110     !-----------------------------------------------
111           USE param
112           USE param1
113           USE parallel
114           USE geometry
115           USE indices
116           USE compar
117           USE functions
118           IMPLICIT NONE
119     !-----------------------------------------------
120     !   G l o b a l   P a r a m e t e r s
121     !-----------------------------------------------
122     !-----------------------------------------------
123     !   D u m m y   A r g u m e n t s
124     !-----------------------------------------------
125     !
126     !
127     !
128     !                      Velocity components
129           DOUBLE PRECISION U(DIMENSION_3), V(DIMENSION_3), W(DIMENSION_3)
130     !
131     !                      Face value of density (for calculating convective fluxes)
132           DOUBLE PRECISION ROP_E(DIMENSION_3), ROP_N(DIMENSION_3), ROP_T(DIMENSION_3)
133     !
134     !                      Convective mass fluxes
135           DOUBLE PRECISION Flux_E(DIMENSION_3), Flux_N(DIMENSION_3), Flux_T(DIMENSION_3)
136     !
137     !                      Indices
138           INTEGER          IJK, IMJK, IJMK, IJKM
139     !-----------------------------------------------
140     
141     !
142     !  Interpolate the face value of density for calculating the convection fluxes
143     !$omp  parallel do default(none) private( IJK, IMJK, IJMK, IJKM) &
144     !$omp              shared(ijkstart3, ijkend3, flux_e, flux_n, flux_t, rop_e, rop_n, rop_t, axz, axy, ayz, u, v, w, do_k)
145           DO IJK = ijkstart3, ijkend3
146     !
147              IF (FLUID_AT(IJK)) THEN
148     !
149     !         East face (i+1/2, j, k)
150                 Flux_E(IJK) = ROP_E(IJK)*AYZ(IJK)*U(IJK)
151     !
152     !         North face (i, j+1/2, k)
153                 Flux_N(IJK) = ROP_N(IJK)*AXZ(IJK)*V(IJK)
154     !
155     !         Top face (i, j, k+1/2)
156                 IF (DO_K) THEN
157                    Flux_T(IJK) = ROP_T(IJK)*AXY(IJK)*W(IJK)
158                 ENDIF
159     !
160     !         West face (i-1/2, j, k)
161                 IMJK = IM_OF(IJK)
162                 IF (.NOT.FLUID_AT(IMJK)) THEN
163                    Flux_E(IMJK) = ROP_E(IMJK)*AYZ(IMJK)*U(IMJK)
164                 ENDIF
165     !
166     !         South face (i, j-1/2, k)
167                 IJMK = JM_OF(IJK)
168                 IF (.NOT.FLUID_AT(IJMK)) THEN
169                   Flux_N(IJMK) = ROP_N(IJMK)*AXZ(IJMK)*V(IJMK)
170                 ENDIF
171     !
172     !         Bottom face (i, j, k-1/2)
173                 IF (DO_K) THEN
174                    IJKM = KM_OF(IJK)
175                    IF (.NOT.FLUID_AT(IJKM)) THEN
176                      Flux_T(IJKM) = ROP_T(IJKM)*AXY(IJKM)*W(IJKM)
177                    ENDIF
178                 ENDIF
179              ENDIF
180           END DO
181     
182           RETURN
183           END SUBROUTINE CALC_MFLUX0
184     !
185     !
186     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
187     !                                                                      C
188     !                                                                      C
189     !  Module name: CALC_MFLUX1(U, V, W, ROP_E, ROP_N, ROP_T, Flux_E,      C
190     !                                         Flux_N, Flux_T, IER)         C
191     !  Purpose: Calculate the convection fluxes.                           C
192     !                                                                      C
193     !  Author: M. Syamlal                                 Date: 31-MAY-05  C
194     !  Reviewer:                                          Date:            C
195     !                                                                      C
196     !                                                                      C
197     !  Literature/Document References:                                     C
198     !                                                                      C
199     !  Variables referenced:                                               C
200     !  Variables modified:                                                 C
201     !                                                                      C
202     !  Local variables:                                                    C
203     !                                                                      C
204     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
205     !
206           SUBROUTINE CALC_MFLUX1(U, V, W, ROP_E, ROP_N, ROP_T, ROPa_E, ROPa_N, &
207                                           ROPa_T, Flux_E, Flux_N, Flux_T, M_AM, IER)
208     !...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
209     !...Switches: -xf
210     !
211     !  Include param.inc file to specify parameter values
212     !
213     !-----------------------------------------------
214     !   M o d u l e s
215     !-----------------------------------------------
216           USE param
217           USE param1
218           USE parallel
219           USE physprop
220           USE geometry
221           USE indices
222           USE compar
223           use fldvar, only: RO_S
224           USE functions
225           IMPLICIT NONE
226     !-----------------------------------------------
227     !   G l o b a l   P a r a m e t e r s
228     !-----------------------------------------------
229     !-----------------------------------------------
230     !   D u m m y   A r g u m e n t s
231     !-----------------------------------------------
232     !
233     !
234     !
235     !                      Velocity components
236           DOUBLE PRECISION U(DIMENSION_3), V(DIMENSION_3), W(DIMENSION_3)
237     !
238     !                      Face value of density (for calculating convective fluxes)
239           DOUBLE PRECISION ROP_E(DIMENSION_3), ROP_N(DIMENSION_3), ROP_T(DIMENSION_3), &
240                            ROPa_E(DIMENSION_3), ROPa_N(DIMENSION_3), ROPa_T(DIMENSION_3)
241     !
242     !                      Convective mass fluxes
243           DOUBLE PRECISION Flux_E(DIMENSION_3), Flux_N(DIMENSION_3), Flux_T(DIMENSION_3)
244     !
245     !                      Error index
246           INTEGER          IER
247     !
248     !                      Indices
249           INTEGER          IJK, IMJK, IJMK, IJKM, M_AM
250     !-----------------------------------------------
251     
252     !
253     !  Interpolate the face value of density for calculating the convection fluxes
254     !!!$omp  parallel do private( IJK, IMJK, IJMK, IJKM) &
255     !!!$omp&  schedule(static)
256           DO IJK = ijkstart3, ijkend3
257     !
258              IF (FLUID_AT(IJK)) THEN
259     !
260     !         East face (i+1/2, j, k)
261     !QX RO_S changed
262                 Flux_E(IJK) = ROP_E(IJK)*(ONE + Cv*ROPa_E(IJK)/RO_S(IJK,M_AM))*AYZ(IJK)*U(IJK)
263     !
264     !         North face (i, j+1/2, k)
265                 Flux_N(IJK) = ROP_N(IJK)*(ONE + Cv*ROPa_N(IJK)/RO_S(IJK,M_AM))*AXZ(IJK)*V(IJK)
266     !
267     !         Top face (i, j, k+1/2)
268                 IF (DO_K) THEN
269                    Flux_T(IJK) = ROP_T(IJK)*(ONE + Cv*ROPa_T(IJK)/RO_S(IJK,M_AM))*AXY(IJK)*W(IJK)
270                 ENDIF
271     !
272     !         West face (i-1/2, j, k)
273                 IMJK = IM_OF(IJK)
274                 IF (.NOT.FLUID_AT(IMJK)) THEN
275                    Flux_E(IMJK) = ROP_E(IMJK)*(ONE + Cv*ROPa_E(IMJK)/RO_S(IJK,M_AM))*AYZ(IMJK)*U(IMJK)
276                 ENDIF
277     !
278     !         South face (i, j-1/2, k)
279                 IJMK = JM_OF(IJK)
280                 IF (.NOT.FLUID_AT(IJMK)) THEN
281                   Flux_N(IJMK) = ROP_N(IJMK)*(ONE + Cv*ROPa_N(IJMK)/RO_S(IJK,M_AM))*AXZ(IJMK)*V(IJMK)
282                 ENDIF
283     !
284     !         Bottom face (i, j, k-1/2)
285                 IF (DO_K) THEN
286                    IJKM = KM_OF(IJK)
287                    IF (.NOT.FLUID_AT(IJKM)) THEN
288                      Flux_T(IJKM) = ROP_T(IJKM)*(ONE + Cv*ROPa_T(IJKM)/RO_S(IJK,M_AM))*AXY(IJKM)*W(IJKM)
289                    ENDIF
290                 ENDIF
291              ENDIF
292           END DO
293     
294           RETURN
295           END SUBROUTINE CALC_MFLUX1
296