File: /nfs/home/0/users/jenkins/mfix.git/model/solve_species_eq.f
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19 SUBROUTINE SOLVE_SPECIES_EQ(IER)
20
21
22
23
24 USE ChiScheme
25 USE bc
26 USE compar
27 USE drag
28 USE energy
29 USE fldvar
30 USE geometry
31 USE indices
32 USE leqsol
33 USE matrix
34 USE mflux
35 USE mpi_utility
36 USE output
37 USE param
38 USE param1
39 USE pgcor
40 USE physprop
41 USE pscor
42 USE residual
43 USE run
44 USE rxns
45 USE sendrecv
46 USE toleranc
47 USE ur_facs
48 Use ambm
49 Use tmp_array, S_p => Array1, S_c => Array2, EPs => Array3, VxGama => Array4
50 use functions
51 use ps
52 use toleranc
53
54 IMPLICIT NONE
55
56
57
58
59 INTEGER, INTENT(INOUT) :: IER
60
61
62
63
64 INTEGER :: M
65
66 INTEGER :: LN
67
68 DOUBLE PRECISION :: apo
69
70 INTEGER :: IJK
71
72 INTEGER :: LEQM, LEQI
73
74 DOUBLE PRECISION :: X_s_temp(DIMENSION_3, DIMENSION_N_s)
75
76
77
78
79
80 INTEGER :: Err_l(0:numPEs-1)
81 INTEGER :: Err_g(0:numPEs-1)
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103 DOUBLE PRECISION , EXTERNAL :: Check_conservation
104
105
106 call lock_ambm
107 call lock_tmp_array
108
109
110
111 = 0
112
113
114
115 IF (SPECIES_EQ(0)) THEN
116 IF(chi_scheme) call set_chi(DISCRETIZE(7), X_g, NMAX(0), &
117 U_g, V_g, W_g)
118
119
120 DO LN = 1, NMAX(0)
121 CALL INIT_AB_M (A_M, B_M, IJKMAX2, 0, IER)
122
123 DO IJK = ijkstart3, ijkend3
124 IF (FLUID_AT(IJK)) THEN
125
126 = ROP_GO(IJK)*VOL(IJK)*ODT
127 S_P(IJK) = APO+&
128 (ZMAX(SUM_R_G(IJK))+ROX_GC(IJK,LN))*VOL(IJK)
129 S_C(IJK) = APO*X_GO(IJK,LN) + &
130 X_G(IJK,LN)*ZMAX((-SUM_R_G(IJK)))*VOL(IJK) +&
131 R_GP(IJK,LN)*VOL(IJK)
132 ELSE
133 S_P(IJK) = ZERO
134 S_C(IJK) = ZERO
135 ENDIF
136 ENDDO
137
138
139 IF(.NOT.ADDED_MASS) THEN
140 CALL CONV_DIF_PHI (X_G(1,LN), DIF_G(1,LN), &
141 DISCRETIZE(7), U_G, V_G, W_G, &
142 Flux_gE, Flux_gN, Flux_gT, 0, A_M, B_M, IER)
143 ELSE
144 CALL CONV_DIF_PHI (X_G(1,LN), DIF_G(1,LN),&
145 DISCRETIZE(7), U_G, V_G, W_G, &
146 Flux_gSE, Flux_gSN, Flux_gST, 0, A_M, B_M, IER)
147 ENDIF
148
149
150 CALL BC_PHI (X_G(1,LN), BC_X_G(1,LN), BC_XW_G(1,LN), &
151 BC_HW_X_G(1,LN), BC_C_X_G(1,LN), 0, A_M, B_M, IER)
152
153
154 CALL SOURCE_PHI (S_P, S_C, EP_G, X_G(1,LN), 0, A_M, B_M)
155
156
157 IF(POINT_SOURCE) CALL POINT_SOURCE_PHI (X_G(1,LN), &
158 PS_X_G(:,LN), PS_MASSFLOW_G, 0, A_M, B_M, IER)
159
160 CALL CALC_RESID_S (X_G(1,LN), A_M, B_M, 0, &
161 NUM_RESID(RESID_X+(LN-1),0), &
162 DEN_RESID(RESID_X+(LN-1),0), RESID(RESID_X+(LN-1),0), &
163 MAX_RESID(RESID_X+(LN-1),0), IJK_RESID(RESID_X+(LN-1),0), &
164 ZERO_X_GS, IER)
165
166 CALL UNDER_RELAX_S (X_G(1,LN), A_M, B_M, 0, UR_FAC(7), IER)
167
168 CALL ADJUST_LEQ (RESID(RESID_X+(LN-1),0), LEQ_IT(7), &
169 LEQ_METHOD(7), LEQI, LEQM, IER)
170
171
172 CALL SOLVE_LIN_EQ ('X_g', 7, X_G(1,LN), A_M, B_M, 0, &
173 LEQI, LEQM, LEQ_SWEEP(7), LEQ_TOL(7), LEQ_PC(7), IER)
174
175
176 IF(ier == -2) Err_l(myPE) = 130
177
178 CALL BOUND_X (X_G(1,LN), IJKMAX2, IER)
179
180 ENDDO
181 IF(chi_scheme) call unset_chi()
182 ENDIF
183
184
185
186
187
188 DO M = 1, SMAX
189 IF (SPECIES_EQ(M)) THEN
190 IF(chi_scheme) THEN
191 DO LN = 1, NMAX(M)
192 DO IJK = ijkstart3, ijkend3
193 X_S_temp(IJK, LN) = X_S(IJK,M,LN)
194 ENDDO
195 ENDDO
196 call set_chi(DISCRETIZE(7), X_S_temp, NMAX(M), &
197 U_S(1,M), V_S(1,M), W_S(1,M))
198 ENDIF
199
200 DO LN = 1, NMAX(M)
201 CALL INIT_AB_M (A_M, B_M, IJKMAX2, M, IER)
202
203
204 DO IJK = ijkstart3, ijkend3
205 IF (FLUID_AT(IJK)) THEN
206 APO = ROP_SO(IJK,M)*VOL(IJK)*ODT
207 S_P(IJK) = APO + &
208 (ZMAX(SUM_R_S(IJK,M))+ROX_SC(IJK,M,LN))*VOL(IJK)
209 S_C(IJK) = APO*X_SO(IJK,M,LN) + &
210 X_S(IJK,M,LN)*ZMAX((-SUM_R_S(IJK,M)))*VOL(IJK) + &
211 R_SP(IJK,M,LN)*VOL(IJK)
212 EPS(IJK) = EP_S(IJK,M)
213 ELSE
214 S_P(IJK) = ZERO
215 S_C(IJK) = ZERO
216 EPS(IJK) = ZERO
217 ENDIF
218 ENDDO
219
220 IF(.NOT.ADDED_MASS .OR. M /= M_AM) THEN
221 CALL CONV_DIF_PHI (X_S(1,M,LN), DIF_S(1,M,LN), &
222 DISCRETIZE(7), U_S(1,M), V_S(1,M), W_S(1,M), &
223 Flux_sE(1,M), Flux_sN(1,M), Flux_sT(1,M), M, A_M, B_M, IER)
224 ELSE
225 CALL CONV_DIF_PHI (X_S(1,M,LN), DIF_S(1,M,LN), &
226 DISCRETIZE(7), U_S(1,M), V_S(1,M), W_S(1,M), &
227 Flux_sSE, Flux_sSN, Flux_sST, M, A_M, B_M, IER)
228 ENDIF
229
230 CALL BC_PHI (X_S(1,M,LN), BC_X_S(1,M,LN), &
231 BC_XW_S(1,M,LN), BC_HW_X_S(1,M,LN), &
232 BC_C_X_S(1,M,LN), M, A_M, B_M, IER)
233
234 CALL SOURCE_PHI (S_P, S_C, EPS, X_S(1,M,LN), M, A_M, B_M)
235
236
237 IF(POINT_SOURCE) CALL POINT_SOURCE_PHI (X_S(1,M,LN), &
238 PS_X_S(:,M,LN), PS_MASSFLOW_S(:,M), M, A_M, B_M, IER)
239
240 CALL CALC_RESID_S (X_S(1,M,LN), A_M, B_M, M, &
241 NUM_RESID(RESID_X+(LN-1),M), &
242 DEN_RESID(RESID_X+(LN-1),M), RESID(RESID_X+(LN-1),&
243 M), MAX_RESID(RESID_X+(LN-1),M), IJK_RESID(RESID_X+(LN-1),M), &
244 ZERO_X_GS, IER)
245
246 CALL UNDER_RELAX_S (X_S(1,M,LN), A_M, B_M, M, UR_FAC(7), IER)
247
248
249
250
251
252
253
254
255
256
257 CALL ADJUST_LEQ (RESID(RESID_X+(LN-1),M), LEQ_IT(7), &
258 LEQ_METHOD(7), LEQI, LEQM, IER)
259
260 CALL SOLVE_LIN_EQ ('X_s', 7, X_S(1,M,LN), A_M, B_M, M,&
261 LEQI, LEQM, LEQ_SWEEP(7), LEQ_TOL(7), LEQ_PC(7), IER)
262
263
264 IF(ier == -2) Err_l(myPE) = 131
265
266 CALL BOUND_X (X_S(1,M,LN), IJKMAX2, IER)
267
268
269 END DO
270
271 if(chi_scheme) call unset_chi()
272 ENDIF
273 END DO
274
275
276
277 call unlock_ambm
278 call unlock_tmp_array
279
280
281
282
283
284
285 CALL global_all_sum(Err_l, Err_g)
286 IER = maxval(Err_g)
287
288 RETURN
289 END SUBROUTINE SOLVE_SPECIES_EQ
290