MFIX  2016-1
calc_p_star.f
Go to the documentation of this file.
1 ! TODO:
2 ! p_star calculation should be based on the sum of volume fractions of
3 ! close-packed solids.
4 
5 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
6 ! C
7 ! Subroutine: CALC_P_star C
8 ! Purpose: Calculate P_star in cells where solids continuity is C
9 ! solved C
10 ! C
11 ! Author: M. Syamlal Date: 21-AUG-96 C
12 ! Reviewer: Date: C
13 ! C
14 ! C
15 ! Literature/Document References: C
16 ! C
17 ! Variables referenced: C
18 ! Variables modified: P_STAR, C
19 ! if yu_standish or fedors_landel: ep_star_array, C
20 ! ep_g_blend_start, C
21 ! ep_g_blend_end C
22 ! Local variables: C
23 ! C
24 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
25 
26  SUBROUTINE calc_p_star(EP_G, P_STAR)
27 
28 !-----------------------------------------------
29 ! Modules
30 !-----------------------------------------------
31  USE param
32  USE param1
33  USE parallel
34  USE geometry
35  USE indices
36  USE physprop
37  USE constant
38  USE pgcor
39  USE pscor
40  USE ur_facs
41  USE residual
42  USE compar
43  USE run
44  USE visc_s
45  USE solids_pressure
46  USE functions
47  IMPLICIT NONE
48 !-----------------------------------------------
49 ! Dummy arguments
50 !-----------------------------------------------
51 ! Gas volume fraction
52  DOUBLE PRECISION, INTENT(IN) :: EP_g(dimension_3)
53 ! Solids pressure
54  DOUBLE PRECISION, INTENT(INOUT) :: P_star(dimension_3)
55 !-----------------------------------------------
56 ! Local variables
57 !-----------------------------------------------
58 !! HPF$ align P_star(:) with TT(:)
59 !! HPF$ align EP_g(:) with TT(:)
60 
61 ! Indices
62  INTEGER :: IJK
63 ! Blend factor
64  DOUBLE PRECISION :: blend
65 !-----------------------------------------------
66 ! External functions
67 !-----------------------------------------------
68  DOUBLE PRECISION :: CALC_EP_STAR
69 !-----------------------------------------------
70 
71 !!$omp parallel do private(ijk)
72 !! HPF$ independent
73 
74  DO ijk = ijkstart3, ijkend3
75  IF (fluid_at(ijk)) THEN
76 
77 
78  IF (yu_standish .OR. fedors_landel) THEN
79 ! if Yu_Standish or Fedors_Landel correlations are used, then
80 ! ep_star_array is modified. this is the only time ep_star_array is
81 ! modified (see set_constprop). (sof Nov-16-2005)
82  ep_star_array(ijk) = calc_ep_star(ijk)
83 
84 ! now the values of ep_g_blend_start and ep_g_blend_end need to be
85 ! reassigned based on the new values of ep_star_array
86  IF(blending_stress.AND.tanh_blend) THEN
87  ep_g_blend_start(ijk) = ep_star_array(ijk) * 0.99d0
88  ep_g_blend_end(ijk) = ep_star_array(ijk) * 1.01d0
89  ELSEIF(blending_stress.AND.sigm_blend) THEN
90  ep_g_blend_start(ijk) = ep_star_array(ijk) * 0.97d0
91  ep_g_blend_end(ijk) = ep_star_array(ijk) * 1.01d0
92  ELSE
93  ep_g_blend_start(ijk) = ep_star_array(ijk)
94  ep_g_blend_end(ijk) = ep_star_array(ijk)
95  ENDIF
96  ENDIF
97 
98  IF (ep_g(ijk) < ep_g_blend_end(ijk)) THEN
99  p_star(ijk) = neg_h(ep_g(ijk),ep_g_blend_end(ijk))
100  IF(blending_stress) THEN
101  blend = blend_function(ijk)
102  p_star(ijk) = (1.0d0-blend) * p_star(ijk)
103  ENDIF
104  ELSE
105  p_star(ijk) = zero
106  ENDIF
107  ENDIF
108  ENDDO
109 
110  RETURN
111  END SUBROUTINE calc_p_star
112 
113 
114 !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
115 ! C
116 ! FUNCTION: CALC_ep_star C
117 ! Purpose: calculate the local value of maximum packing C
118 ! C
119 ! Author: D. Gera and M. Syamlal Date: 31-DEC-02 C
120 ! Reviewer: Date: C
121 ! Modified: S. Benyahia Date: 02-May-05 C
122 ! C
123 ! Literature/Document References: C
124 ! A.B. Yu and N. Standish. Powder Tech, 52 (1987) 233-241 C
125 ! R.F. Fedors and R.F. Landel. Powder Tech, 23 (1979) 225-231 C
126 ! C
127 ! Variables referenced: C
128 ! Variables modified: C
129 ! Local variables: C
130 ! C
131 !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
132 
133  DOUBLE PRECISION FUNCTION calc_ep_star(IJK)
134 
135 !-----------------------------------------------
136 ! Modules
137 !-----------------------------------------------
138  USE param
139  USE param1
140  USE fldvar
141  USE geometry
142  USE indices
143  USE physprop
144  USE constant
145  USE toleranc
146  USE compar
147  USE run
148  USE fun_avg
149  USE functions
150  IMPLICIT NONE
151 !-----------------------------------------------
152 ! Dummy arguments
153 !-----------------------------------------------
154 ! IJK index
155  INTEGER, INTENT(IN) :: IJK
156 !-----------------------------------------------
157 ! Local variables
158 !-----------------------------------------------
159 ! Indices
160  INTEGER :: I, J
161 
162 ! DOUBLE PRECISION :: xbar
163 
164 ! start sof modifications (02-May-05)
165 ! maximum packing for the mixture
166  DOUBLE PRECISION :: P_IT(mmax)
167 ! true maximum packing for the mixture
168  DOUBLE PRECISION :: EPs_max_local
169 ! maximum packing fraction for a binary mixture
170  DOUBLE PRECISION :: P_IJ(mmax, mmax)
171 ! particle diameter ratio
172  DOUBLE PRECISION :: R_IJ(mmax, mmax)
173 ! fractional solids volume corresponding to P_IJ
174  DOUBLE PRECISION :: X_IJ(mmax, mmax)
175 ! fractional solids volume in a mixture
176 ! this is Xj in eq. 22 of Yu-Standish
177  DOUBLE PRECISION :: COMP_X_I(mmax), SUM_LOCAL
178 ! local aliases for particle diameter, solids volume fraction and the
179 ! maximum solids volume fraction which are used to rearrange solids
180 ! phases from coarsest to finest
181  DOUBLE PRECISION :: DP_TMP(mmax), EPs_TMP(mmax), &
182  EPs_max_TMP(mmax), old_value
183 !-----------------------------------------------
184 
185 
186  IF (call_dqmom) THEN
187 ! sort particles to start from coarsest to finest particles
188 ! assigning values to local aliases
189  DO i = 1, smax
190  dp_tmp(i) = d_p(ijk,i)
191  eps_tmp(i) = ep_s(ijk,i)
192  eps_max_tmp(i) = ep_s_max(i)
193  ENDDO
194 ! sorting particles from coarse to fine
195  DO i = 1, smax
196  DO j = i , smax
197 ! check if phase J is larger than phase I
198  IF(dp_tmp(i) < dp_tmp(j)) THEN
199 ! temporarily store phase i diameter
200  old_value = dp_tmp(i)
201 ! overwrite phase i diameter with smaller phase j diameter
202  dp_tmp(i) = dp_tmp(j)
203 ! overwrite phase j diameter with stoired phase i diameter
204  dp_tmp(j) = old_value
205 
206  old_value = eps_tmp(i)
207  eps_tmp(i) = eps_tmp(j)
208  eps_tmp(j) = old_value
209 
210  old_value = eps_max_tmp(i)
211  eps_max_tmp(i) = eps_max_tmp(j)
212  eps_max_tmp(j) = old_value
213  ENDIF
214  ENDDO
215  ENDDO
216 
217  ELSE ! not dqmom
218 
219 ! assigning values to local aliases
220  DO i = 1, smax
221  dp_tmp(i) = d_p(ijk,m_max(i))
222  eps_tmp(i) = ep_s(ijk,m_max(i))
223  eps_max_tmp(i) = ep_s_max(m_max(i))
224  ENDDO
225  ENDIF ! end if/else (call_dqmom)
226 
227 
228 ! this is the way the algorithm was written by Yu and Standish (sof).
229 ! compute equations 25 in Yu-Standish
230 ! (this is also needed by Fedors_Landel)
231  DO i = 1, smax
232  sum_local = zero
233  DO j = 1, smax
234  IF(i .GE. j) THEN
235  r_ij(i,j) = dp_tmp(i)/dp_tmp(j)
236  ELSE
237  r_ij(i,j) = dp_tmp(j)/dp_tmp(i)
238  ENDIF
239  sum_local = sum_local + eps_tmp(j)
240  ENDDO ! end do (j=1,smax)
241 
242  IF(sum_local > dil_ep_s) THEN
243 ! fractional solids volume see eq. 20
244  comp_x_i(i) = eps_tmp(i)/sum_local
245  ELSE
246 ! return first phase ep_s_max in case very dilute
247  calc_ep_star = one - eps_max_tmp(1)
248  RETURN
249  ENDIF
250  ENDDO ! end do (i=1,smax)
251 
252 ! Begin YU_STANDISH section
253 ! ---------------------------------------------------------------->>>
254  IF(yu_standish) THEN
255 ! compute equation 23-24 in Yu-Standish
256  DO i = 1, smax
257  DO j = 1, smax
258  IF(r_ij(i,j) .LE. 0.741d0) THEN
259  IF(j .LT. i) THEN
260  x_ij(i,j) = (one - r_ij(i,j)*r_ij(i,j))/&
261  (2.0d0 - eps_max_tmp(i))
262  ELSE
263  x_ij(i,j) = one - (one - r_ij(i,j)*r_ij(i,j))/&
264  (2.0d0 - eps_max_tmp(i))
265  ENDIF
266  p_ij(i, j) = eps_max_tmp(i) + eps_max_tmp(i)*&
267  (one-eps_max_tmp(i)) * (one - 2.35d0*r_ij(i,j) + &
268  1.35d0*r_ij(i,j)*r_ij(i,j) )
269  ELSE
270  p_ij(i, j) = eps_max_tmp(i)
271  ENDIF
272  ENDDO ! end do (j=1,smax)
273  ENDDO ! end do (i=1,smax)
274 
275 ! Compute equation 22
276  eps_max_local = one
277  DO i = 1, smax
278  sum_local = zero
279 
280  IF(i .GE. 2) THEN
281  DO j = 1, (i-1)
282  IF(p_ij(i,j) == eps_max_tmp(i)) THEN
283  sum_local = sum_local
284  ELSE
285  sum_local = sum_local + (one - eps_max_tmp(i)/&
286  p_ij(i,j))*comp_x_i(j)/x_ij(i,j)
287  ENDIF
288  ENDDO
289  ENDIF
290 
291  IF((i+1) .LE. smax) THEN
292  DO j = (i+1), smax
293  IF( p_ij(i, j) == eps_max_tmp(i) ) THEN
294  sum_local = sum_local
295  ELSE
296  sum_local = sum_local + (one - eps_max_tmp(i)/&
297  p_ij(i, j))*comp_x_i(j)/x_ij(i, j)
298  ENDIF
299  ENDDO
300  ENDIF
301 
302  IF (sum_local .NE. zero) THEN
303  p_it(i) = eps_max_tmp(i)/(one - sum_local)
304  ELSE
305 ! do nothing if particles have same diameter
306  p_it(i) = one
307  ENDIF
308 
309  eps_max_local = min(p_it(i), eps_max_local)
310  ENDDO ! end do (i=1,smax)
311 
312 ! for the case of all phases having same diameter
313 
314  IF (eps_max_local == one) eps_max_local = eps_max_tmp(1)
315  calc_ep_star = one - eps_max_local
316 ! end YU_STANDISH section
317 ! ----------------------------------------------------------------<<<
318 
319 ! Part implemented by Dinesh for binary mixture, uncomment to use (Sof)
320 
321 ! if ((EP_s(IJK,1)+EP_s(IJK,2)) .NE. ZERO) THEN
322 ! xbar = EP_s(IJK,1)/(EP_s(IJK,1)+EP_s(IJK,2))
323 
324 ! if (xbar .LE. ep_s_max_ratio(1,2)) THEN
325 ! CALC_EP_star =MAX(0.36d0, (ONE-(((ep_s_max(1)-ep_s_max(2))+&
326 ! (ONE-d_p_ratio(1,2))*(ONE-ep_s_max(1))*ep_s_max(2))*(ep_s_max(1)+&
327 ! (ONE-ep_s_max(1)) *ep_s_max(2))*xbar/ep_s_max(1)+ep_s_max(2))))
328 ! else
329 ! CALC_EP_star =MAX(0.36d0, (ONE-((ONE -d_p_ratio(1,2))*(ep_s_max(1)&
330 ! +(ONE-ep_s_max(1))*ep_s_max(2))*(ONE -xbar) +ep_s_max(1))))
331 ! end if
332 ! else
333 ! CALC_EP_star = ONE - MIN(ep_s_max(1), ep_s_max(2)) !corrected by sof
334 ! end if
335 
336 ! Use the code (below) instead of the above commented code because the
337 ! phases were not rearranged and I didn't want to modify it (sof)
338 ! If you don't understand what's going on, contact me: sof@fluent.com
339 
340 ! In the case of binary mixture (Fedors-Landel empirical correlation)
341 ! ---------------------------------------------------------------->>>
342  ELSEIF(fedors_landel) THEN
343 
344  IF(comp_x_i(1) .LE. (eps_max_tmp(1)/(eps_max_tmp(1)+ &
345  (one - eps_max_tmp(1))*eps_max_tmp(2))) ) THEN
346 
347  calc_ep_star = (eps_max_tmp(1) - eps_max_tmp(2) + &
348  (1 - sqrt(r_ij(2,1))) * (one - eps_max_tmp(1)) * &
349  eps_max_tmp(2) )*&
350  (eps_max_tmp(1) + (one - eps_max_tmp(1)) * &
351  eps_max_tmp(2)) * comp_x_i(1)/eps_max_tmp(1) + &
352  eps_max_tmp(2)
353  ELSE
354  calc_ep_star = (one-sqrt(r_ij(2,1))) * (eps_max_tmp(1)+&
355  (one-eps_max_tmp(1)) * eps_max_tmp(2)) * &
356  (one - comp_x_i(1)) + eps_max_tmp(1)
357  ENDIF
358 ! this is gas volume fraction at packing
359  calc_ep_star = one - calc_ep_star
360  ENDIF ! for Yu_Standish and Fedors_Landel correlations
361 ! end FEDORS_LANDEL correlation
362 ! ----------------------------------------------------------------<<<
363 
364  RETURN
365  END FUNCTION calc_ep_star
calc_ep_star
double precision function calc_ep_star(IJK)
Definition: calc_p_star.f:134
visc_s
Definition: visc_s_mod.f:1
compar::ijkend3
integer ijkend3
Definition: compar_mod.f:80
physprop::blend_function
double precision function blend_function(IJK)
Definition: physprop_mod.f:244
pgcor
Definition: pgcor_mod.f:1
param1::one
double precision, parameter one
Definition: param1_mod.f:29
param::dimension_3
integer dimension_3
Definition: param_mod.f:11
run::sigm_blend
logical sigm_blend
Definition: run_mod.f:163
constant::ep_s_max
double precision, dimension(dim_m) ep_s_max
Definition: constant_mod.f:38
visc_s::ep_star_array
double precision, dimension(:), allocatable ep_star_array
Definition: visc_s_mod.f:54
constant::m_max
integer, dimension(dim_m) m_max
Definition: constant_mod.f:41
fun_avg
Definition: fun_avg_mod.f:1
solids_pressure
Definition: solids_pressure_mod.f:25
fldvar::d_p
double precision, dimension(:,:), allocatable d_p
Definition: fldvar_mod.f:57
visc_s::ep_g_blend_end
double precision, dimension(:), allocatable ep_g_blend_end
Definition: visc_s_mod.f:60
physprop::mmax
integer mmax
Definition: physprop_mod.f:19
solids_pressure::neg_h
double precision function neg_h(XXX, YYY)
Definition: solids_pressure_mod.f:45
ur_facs
Definition: ur_facs_mod.f:1
visc_s::ep_g_blend_start
double precision, dimension(:), allocatable ep_g_blend_start
Definition: visc_s_mod.f:57
pscor
Definition: pscor_mod.f:1
run::call_dqmom
logical call_dqmom
Definition: run_mod.f:127
run
Definition: run_mod.f:13
run::yu_standish
logical yu_standish
Definition: run_mod.f:180
param
Definition: param_mod.f:2
toleranc::dil_ep_s
double precision, parameter dil_ep_s
Definition: toleranc_mod.f:24
run::blending_stress
logical blending_stress
Definition: run_mod.f:161
run::fedors_landel
logical fedors_landel
Definition: run_mod.f:183
compar::ijkstart3
integer ijkstart3
Definition: compar_mod.f:80
fldvar::ep_s
double precision function ep_s(IJK, xxM)
Definition: fldvar_mod.f:178
calc_p_star
subroutine calc_p_star(EP_G, P_STAR)
Definition: calc_p_star.f:27
physprop::smax
integer smax
Definition: physprop_mod.f:22
run::tanh_blend
logical tanh_blend
Definition: run_mod.f:162
param1::zero
double precision, parameter zero
Definition: param1_mod.f:27