MFIX  2016-1
conv_rop.f
Go to the documentation of this file.
1 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2 ! C
3 ! Subroutine: CONV_ROP C
4 ! Purpose: Calculate the face value of density used for calculating C
5 ! convection fluxes. Master routine. C
6 ! C
7 ! Author: M. Syamlal Date: 31-MAY-05 C
8 ! C
9 ! C
10 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
11  SUBROUTINE conv_rop()
12 
13 ! Modules
14 !---------------------------------------------------------------------//
15  USE fldvar, only: rop_g, u_g, v_g, w_g
16  USE fldvar, only: rop_s, u_s, v_s, w_s
17  USE mflux, only: rop_ge, rop_gn, rop_gt
18  USE mflux, only: rop_se, rop_sn, rop_st
19  USE physprop, only: mmax
20  USE run, only: discretize
21  IMPLICIT NONE
22 
23 !---------------------------------------------------------------------//
24 ! Local variables
25 !---------------------------------------------------------------------//
26 ! solids phase index
27  INTEGER :: M
28 !---------------------------------------------------------------------//
29 
30  IF (discretize(1) == 0) THEN ! 0 & 1 => first order upwinding
31  CALL conv_rop0 (rop_g, u_g, v_g, w_g, &
32  rop_ge, rop_gn, rop_gt)
33  ELSE
34  CALL conv_rop1 (discretize(1), rop_g, u_g, v_g, w_g, &
35  rop_ge, rop_gn, rop_gt)
36  ENDIF
37 
38  IF (discretize(2) == 0) THEN ! 0 & 1 => first order upwinding
39  DO m = 1, mmax
40  CALL conv_rop0 (rop_s(1,m), u_s(1,m), v_s(1,m), w_s(1,m), &
41  rop_se(1,m), rop_sn(1,m), rop_st(1,m))
42  ENDDO
43  ELSE
44  DO m = 1, mmax
45  CALL conv_rop1 (discretize(2), rop_s(1,m), &
46  u_s(1,m), v_s(1,m), w_s(1,m), &
47  rop_se(1,m), rop_sn(1,m), rop_st(1,m))
48  ENDDO
49  ENDIF
50 
51  RETURN
52  END SUBROUTINE conv_rop
53 
54 
55 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
56 ! C
57 ! Subroutine: CONV_ROP C
58 ! Purpose: Calculate the face value of density used for calculating C
59 ! convection fluxes. FOU routine. C
60 ! C
61 ! Author: M. Syamlal Date: 31-MAY-05 C
62 ! C
63 ! C
64 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
65  SUBROUTINE conv_rop0(ROP, U, V, W, ROP_E, ROP_N, ROP_T)
66 
67 ! Modules
68 !---------------------------------------------------------------------//
69  USE compar, only: ijkstart3, ijkend3
70  USE functions, only: fluid_at
71  USE functions, only: east_of, north_of, top_of
72  USE functions, only: west_of, south_of, bottom_of
73  USE functions, only: im_of, jm_of, km_of
74  USE geometry, only: do_k
75  USE param, only: dimension_3
76  USE param1, only: zero
77  IMPLICIT NONE
78 
79 ! Dummy arguments
80 !---------------------------------------------------------------------//
81 ! macroscopic density (rho_prime)
82  DOUBLE PRECISION, INTENT(IN) :: ROP(dimension_3)
83 ! Velocity components
84  DOUBLE PRECISION, INTENT(IN) :: U(dimension_3)
85  DOUBLE PRECISION, INTENT(IN) :: V(dimension_3)
86  DOUBLE PRECISION, INTENT(IN) :: W(dimension_3)
87 ! Face value of density (for calculating convective fluxes)
88  DOUBLE PRECISION, INTENT(OUT) :: ROP_E(dimension_3)
89  DOUBLE PRECISION, INTENT(OUT) :: ROP_N(dimension_3)
90  DOUBLE PRECISION, INTENT(OUT) :: ROP_T(dimension_3)
91 
92 ! Local variables
93 !---------------------------------------------------------------------//
94 ! Indices
95  INTEGER :: IJK
96  INTEGER :: IJKE, IJKN, IJKT
97  INTEGER :: IJKW, IJKS, IJKB
98  INTEGER :: IMJK, IJMK, IJKM
99 !---------------------------------------------------------------------//
100 
101 
102 !$omp parallel do default(none) &
103 !$omp private(IJK, IJKE, IJKN, IJKT, IJKW, &
104 !$omp IJKS, IJKB, IMJK, IJMK, IJKM) &
105 !$omp shared(ijkstart3, ijkend3, u, v, w, do_k, rop, &
106 !$omp rop_e, rop_n, rop_t)
107  DO ijk = ijkstart3, ijkend3
108 
109  IF (fluid_at(ijk)) THEN
110  ijke = east_of(ijk)
111  ijkn = north_of(ijk)
112  ijkt = top_of(ijk)
113 
114  imjk = im_of(ijk)
115  ijmk = jm_of(ijk)
116 
117 ! East face (i+1/2, j, k)
118  IF (u(ijk) >= zero) THEN
119  rop_e(ijk) = rop(ijk)
120  ELSE
121  rop_e(ijk) = rop(ijke)
122  ENDIF
123 ! West face (i-1/2, j, k)
124  IF (.NOT.fluid_at(imjk)) THEN
125  ijkw = west_of(ijk)
126  IF (u(imjk) >= zero) THEN
127  rop_e(imjk) = rop(ijkw)
128  ELSE
129  rop_e(imjk) = rop(ijk)
130  ENDIF
131  ENDIF
132 
133 
134 ! North face (i, j+1/2, k)
135  IF (v(ijk) >= zero) THEN
136  rop_n(ijk) = rop(ijk)
137  ELSE
138  rop_n(ijk) = rop(ijkn)
139  ENDIF
140 ! South face (i, j-1/2, k)
141  IF (.NOT.fluid_at(ijmk)) THEN
142  ijks = south_of(ijk)
143  IF (v(ijmk) >= zero) THEN
144  rop_n(ijmk) = rop(ijks)
145  ELSE
146  rop_n(ijmk) = rop(ijk)
147  ENDIF
148  ENDIF
149 
150 
151  IF (do_k) THEN
152  ijkm = km_of(ijk)
153 ! Top face (i, j, k+1/2)
154  IF (w(ijk) >= zero) THEN
155  rop_t(ijk) = rop(ijk)
156  ELSE
157  rop_t(ijk) = rop(ijkt)
158  ENDIF
159 ! Bottom face (i, j, k-1/2)
160  IF (.NOT.fluid_at(ijkm)) THEN
161  ijkb = bottom_of(ijk)
162  IF (w(ijkm) >= zero) THEN
163  rop_t(ijkm) = rop(ijkb)
164  ELSE
165  rop_t(ijkm) = rop(ijk)
166  ENDIF
167  ENDIF
168  ENDIF ! end if do_k
169 
170  ENDIF ! end if fluid_at
171  ENDDO ! end do ijk
172 
173  RETURN
174  END SUBROUTINE conv_rop0
175 
176 
177 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
178 ! C
179 ! Subroutine: CONV_ROP1 C
180 ! Purpose: Calculate the face value of density used for calculating C
181 ! convection fluxes. HR routine. Here interpolate the face value of C
182 ! density. C
183 ! C
184 ! Author: M. Syamlal Date: 31-MAY-05 C
185 ! C
186 ! C
187 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
188  SUBROUTINE conv_rop1(DISC, ROP, U, V, W, ROP_E, ROP_N, ROP_T)
189 
190 ! Modules
191 !---------------------------------------------------------------------//
192  USE compar, only: ijkstart3, ijkend3
193  USE functions, only: fluid_at
194  USE functions, only: east_of, north_of, top_of
195  USE functions, only: west_of, south_of, bottom_of
196  USE functions, only: im_of, jm_of, km_of
197  USE geometry, only: do_k
198  USE param, only: dimension_3
199  USE param1, only: one
200  USE xsi, only: calc_xsi
201  IMPLICIT NONE
202 
203 ! Dummy arguments
204 !---------------------------------------------------------------------//
205 ! Discretization scheme
206  INTEGER, INTENT(IN) :: DISC
207 ! macroscopic density (rho_prime)
208  DOUBLE PRECISION, INTENT(IN) :: ROP(dimension_3)
209 ! Velocity components
210  DOUBLE PRECISION, INTENT(IN) :: U(dimension_3)
211  DOUBLE PRECISION, INTENT(IN) :: V(dimension_3)
212  DOUBLE PRECISION, INTENT(IN) :: W(dimension_3)
213 ! Face value of density (for calculating convective fluxes)
214  DOUBLE PRECISION, INTENT(OUT) :: ROP_E(dimension_3)
215  DOUBLE PRECISION, INTENT(OUT) :: ROP_N(dimension_3)
216  DOUBLE PRECISION, INTENT(OUT) :: ROP_T(dimension_3)
217 !
218 ! Local variables
219 !---------------------------------------------------------------------//
220  INTEGER :: IJK, IJKE, IJKN, IJKT
221  INTEGER :: IJKW, IJKS, IJKB, IMJK, IJMK, IJKM
222  Integer :: incr
223 
224  DOUBLE PRECISION, DIMENSION(DIMENSION_3) :: XSI_e, XSI_n, XSI_t
225 
226 !---------------------------------------------------------------------//
227 
228 ! Calculate factors
229  incr=0
230  CALL calc_xsi (disc, rop, u, v, w, xsi_e, xsi_n, xsi_t, incr)
231 
232 !!!$omp parallel do private(IJK, IJKE, IJKN, IJKT, IJKW, IJKS, IJKB, &
233 !!!$omp IMJK, IJMK, IJKM) &
234 !!!$omp& schedule(static)
235  DO ijk = ijkstart3, ijkend3
236 
237  IF (fluid_at(ijk)) THEN
238  ijke = east_of(ijk)
239  ijkn = north_of(ijk)
240  ijkt = top_of(ijk)
241 
242  imjk = im_of(ijk)
243  ijmk = jm_of(ijk)
244 
245 ! East face (i+1/2, j, k)
246  rop_e(ijk) = ((one-xsi_e(ijk))*rop(ijk)+&
247  xsi_e(ijk)*rop(ijke))
248 ! West face (i-1/2, j, k)
249  IF (.NOT.fluid_at(imjk)) THEN
250  ijkw = west_of(ijk)
251  rop_e(imjk) = ((one - xsi_e(imjk))*rop(ijkw)+&
252  xsi_e(imjk)*rop(ijk))
253  ENDIF
254 
255 
256 ! North face (i, j+1/2, k)
257  rop_n(ijk) = ((one-xsi_n(ijk))*rop(ijk)+&
258  xsi_n(ijk)*rop(ijkn))
259 ! South face (i, j-1/2, k)
260  IF (.NOT.fluid_at(ijmk)) THEN
261  ijks = south_of(ijk)
262  rop_n(ijmk) = ((one - xsi_n(ijmk))*rop(ijks)+&
263  xsi_n(ijmk)*rop(ijk))
264  ENDIF
265 
266 
267  IF (do_k) THEN
268  ijkm = km_of(ijk)
269 
270 ! Top face (i, j, k+1/2)
271  rop_t(ijk) = ((one - xsi_t(ijk))*rop(ijk)+&
272  xsi_t(ijk)*rop(ijkt))
273 ! Bottom face (i, j, k-1/2)
274  IF (.NOT.fluid_at(ijkm)) THEN
275  ijkb = bottom_of(ijk)
276  rop_t(ijkm) = ((one - xsi_t(ijkm))*rop(ijkb)+&
277  xsi_t(ijkm)*rop(ijk))
278  ENDIF
279  ENDIF ! end if do_k
280 
281  ENDIF ! end if fluid_at
282  ENDDO ! end do ijk
283 
284  RETURN
285  END SUBROUTINE conv_rop1
mflux::rop_st
double precision, dimension(:,:), allocatable rop_st
Definition: mflux_mod.f:46
fldvar::v_s
double precision, dimension(:,:), allocatable v_s
Definition: fldvar_mod.f:105
compar::ijkend3
integer ijkend3
Definition: compar_mod.f:80
param1::one
double precision, parameter one
Definition: param1_mod.f:29
param::dimension_3
integer dimension_3
Definition: param_mod.f:11
mflux::rop_gn
double precision, dimension(:), allocatable rop_gn
Definition: mflux_mod.f:38
fldvar::w_s
double precision, dimension(:,:), allocatable w_s
Definition: fldvar_mod.f:117
xsi::calc_xsi
subroutine calc_xsi(DISCR, PHI, U, V, W, XSI_E, XSI_N, XSI_T, incr)
Definition: xsi_mod.f:23
xsi
Definition: xsi_mod.f:15
fldvar::u_s
double precision, dimension(:,:), allocatable u_s
Definition: fldvar_mod.f:93
physprop::mmax
integer mmax
Definition: physprop_mod.f:19
conv_rop
subroutine conv_rop()
Definition: conv_rop.f:12
fldvar::v_g
double precision, dimension(:), allocatable v_g
Definition: fldvar_mod.f:99
mflux::rop_ge
double precision, dimension(:), allocatable rop_ge
Definition: mflux_mod.f:36
fldvar::w_g
double precision, dimension(:), allocatable w_g
Definition: fldvar_mod.f:111
run
Definition: run_mod.f:13
param
Definition: param_mod.f:2
conv_rop0
subroutine conv_rop0(ROP, U, V, W, ROP_E, ROP_N, ROP_T)
Definition: conv_rop.f:66
mflux::rop_gt
double precision, dimension(:), allocatable rop_gt
Definition: mflux_mod.f:40
geometry::do_k
logical do_k
Definition: geometry_mod.f:30
mflux::rop_se
double precision, dimension(:,:), allocatable rop_se
Definition: mflux_mod.f:44
compar::ijkstart3
integer ijkstart3
Definition: compar_mod.f:80
fldvar::u_g
double precision, dimension(:), allocatable u_g
Definition: fldvar_mod.f:87
run::discretize
integer, dimension(dim_eqs) discretize
Definition: run_mod.f:67
conv_rop1
subroutine conv_rop1(DISC, ROP, U, V, W, ROP_E, ROP_N, ROP_T)
Definition: conv_rop.f:189
fldvar::rop_s
double precision, dimension(:,:), allocatable rop_s
Definition: fldvar_mod.f:51
mflux::rop_sn
double precision, dimension(:,:), allocatable rop_sn
Definition: mflux_mod.f:42
fldvar::rop_g
double precision, dimension(:), allocatable rop_g
Definition: fldvar_mod.f:38
param1::zero
double precision, parameter zero
Definition: param1_mod.f:27