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

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Subroutine: set_constprop                                           C
4     !  Purpose: This routine serves two purposes:                          C
5     !    1) initializes various variables everywhere in the domain with    C
6     !       a zero value. the real need for this is unclear. undefined     C
7     !       may be a better approach...                                    C
8     !    2) if defined, sets physical properties to their specified        C
9     !       constant value in the fluid domain. cannot set in flow         C
10     !       boundaries or later checks will also complain (may be          C
11     !       overly strict check)                                           C
12     !                                                                      C
13     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
14           SUBROUTINE SET_CONSTPROP
15     
16     ! Modules
17     !-----------------------------------------------
18           USE param1, only: zero, half, one, undefined
19     
20           USE fldvar, only: ro_g
21           USE fldvar, only: ro_s, d_p
22           use fldvar, only: p_s
23     
24           USE visc_g, only: mu_gt, epmu_gt, lambda_gt, eplambda_gt
25           USE visc_g, only: recalc_visc_g
26           USE visc_s, only: mu_s, epmu_s, lambda_s, eplambda_s, lambda_s_c
27           USE visc_s, only: ep_star_array
28           USE visc_s, only: ep_g_blend_start, ep_g_blend_end
29     
30           USE physprop, only: nmax, mmax, smax
31           USE physprop, only: ro_s0, ro_g0, d_p0
32           USE physprop, only: mu_g0, mu_g, mu_s0
33           USE physprop, only: mw_avg, mw_mix_g
34           USE physprop, only: c_pg0, k_g0, c_pg, k_g, dif_g0, dif_g
35           USE physprop, only: c_ps0, k_s0, c_ps, k_s, dif_s0, dif_s
36           USE physprop, only: cv
37     
38           USE constant, only: ep_s_max_ratio, d_p_ratio, ep_s_max, m_max
39           use constant, only: ep_star, l_scale0
40     
41           USE run, only: energy_eq, k_epsilon, call_dqmom
42           USE run, only: yu_standish, fedors_landel
43           USE run, only: kt_type_enum, ia_2005, gd_1999, gtsh_2012
44           USE run, only: blending_stress, sigm_blend, tanh_blend
45     
46           USE drag, only: f_gs, f_ss
47     
48           use kintheory, only: mu_sm_ip, mu_sl_ip, xi_sm_ip, xi_sl_ip
49           use kintheory, only: fnu_s_ip, ft_sm_ip, ft_sl_ip
50           use kintheory, only: kth_sl_ip, knu_sm_ip, knu_sl_ip, kvel_s_ip
51           use kintheory, only: ed_ss_ip, edvel_sl_ip, edt_s_ip, edvel_sm_ip
52           use kintheory, only: a2_gtsh, xsi_gtsh
53     
54           use mms, only: use_mms
55     
56           USE compar, only: ijkstart3, ijkend3
57           use functions, only: wall_at, fluid_at
58     
59           IMPLICIT NONE
60     !-----------------------------------------------
61     ! Local variables
62     !-----------------------------------------------
63     ! indices
64           INTEGER :: IJK, M, I, J
65           DOUBLE PRECISION :: old_value, DP_TMP(SMAX)
66     !-----------------------------------------------
67     
68     ! First, initialize certain transport coefficients, physical
69     ! properties, and other variable types everywhere in the
70     ! domain to zero. Then, set these to their specified value
71     ! (if defined) in fluid cells. Some are also set in flow cells.
72     ! NOTE: DO NOT simply zero existing field variables.
73           RO_g = ZERO
74           MU_g = ZERO
75           MU_gt = ZERO
76           EPMU_GT = ZERO
77           LAMBDA_GT = ZERO
78           EPLAMBDA_GT = ZERO
79           K_g = ZERO
80           C_PG = ZERO
81           MW_MIX_G = zero
82           DIF_g = ZERO
83           F_GS = ZERO
84     
85           RO_S = ZERO
86           D_P = zero     ! this is done in init_fvars
87           MU_s = ZERO
88           EPMU_s = ZERO
89           LAMBDA_s_c = ZERO
90           LAMBDA_s = ZERO
91           EPLAMBDA_s = ZERO
92     ! not defining p_s in flow cells can become a problem for certain
93     ! cases as gradients in P_s are used (e.g., when a PO is west of a
94     ! MI...?). Simply assigning it a zero value may not be the best
95     ! approach but it is currently used....
96           P_S = ZERO
97           K_s = ZERO
98           C_PS = ZERO
99           DIF_S = ZERO
100           F_SS = ZERO
101     
102     ! Set the flag for recalculating gas viscosity.
103           RECALC_VISC_G = (ENERGY_EQ .OR. L_SCALE0/=ZERO .OR. K_EPSILON)
104     
105     ! Set default value for virtual mass coefficient
106           Cv = HALF
107     
108     ! Variables specific to various kinetic theory models
109           IF (KT_TYPE_ENUM == IA_2005) THEN
110              MU_sM_ip = ZERO
111              MU_sL_ip = ZERO
112              XI_sM_ip = ZERO
113              XI_sL_ip = ZERO
114              Fnu_s_ip = ZERO
115              FT_sM_ip = ZERO
116              FT_sL_ip = ZERO
117              Kth_sL_ip = ZERO
118              Knu_sM_ip = ZERO
119              Knu_sL_ip = ZERO
120              Kvel_s_ip = ZERO
121              ED_ss_ip = ZERO
122              EDvel_sL_ip = ZERO
123           ENDIF
124     
125           IF (KT_TYPE_ENUM == IA_2005 .OR. &
126               KT_TYPE_ENUM == GD_1999 .OR.  &
127               KT_TYPE_ENUM == GTSH_2012) THEN
128              EDT_s_ip = ZERO
129              EDvel_sM_ip = ZERO
130           ENDIF
131     
132           IF(KT_TYPE_ENUM == GTSH_2012) THEN
133              A2_gtsh = ZERO
134              xsi_gtsh = zero
135           ENDIF
136     
137     ! Set specified constant physical properties values
138           DO IJK = ijkstart3, ijkend3
139     
140              IF (.NOT.WALL_AT(IJK)) THEN
141     ! Fluid and inflow/outflow cells: FLAG < 100
142                 IF (RO_G0 /= UNDEFINED) RO_G(IJK) = RO_G0
143                 IF (C_PG0 /= UNDEFINED) C_PG(IJK) = C_PG0
144                 IF (MW_AVG /= UNDEFINED) MW_MIX_G(IJK) = MW_AVG
145              ENDIF
146     
147              IF (FLUID_AT(IJK)) THEN
148     ! Strictly Fluid cells: FLAG = 1
149                 IF (MU_G0 /= UNDEFINED) THEN
150                    MU_G(IJK) = MU_G0
151                    MU_GT(IJK) = MU_G0
152                    EPMU_GT(IJK) = MU_G0
153                    LAMBDA_GT(IJK) = -(2.0d0/3.0d0)*MU_G0
154                    EPLAMBDA_GT(IJK) = -(2.0d0/3.0d0)*MU_G0
155                 ENDIF
156                 IF (K_G0 /= UNDEFINED) K_G(IJK) = K_G0
157                 IF (DIF_G0 /= UNDEFINED) DIF_G(IJK,:NMAX(0)) = DIF_G0
158              ENDIF
159     
160              IF (USE_MMS) THEN
161                 IF (MU_G0 /= UNDEFINED) THEN
162                    MU_G(IJK) = MU_G0
163                    MU_GT(IJK) = MU_G0
164                    EPMU_GT(IJK) = MU_G0
165                    LAMBDA_GT(IJK) = -(2.0d0/3.0d0)*MU_G0
166                    EPLAMBDA_GT(IJK) = -(2.0d0/3.0d0)*MU_G0
167                 ENDIF
168                 IF (K_G0 /= UNDEFINED) K_G(IJK) = K_G0
169                 IF (DIF_G0 /= UNDEFINED) DIF_G(IJK,:NMAX(0)) = DIF_G0
170              ENDIF
171     
172           ENDDO
173     
174     
175           DO M = 1, MMAX
176              DO IJK = ijkstart3, ijkend3
177                 IF(.NOT.WALL_AT(IJK)) THEN
178     ! Fluid and inflow/outflow cells: FLAG < 100
179                    IF (RO_S0(M) /= UNDEFINED) RO_S(IJK,M) = RO_S0(M)
180                    IF (C_PS0(M) /= UNDEFINED) C_PS(IJK,M) = C_PS0(M)
181                    IF (D_P0(M) /= UNDEFINED) D_P(IJK,M) = D_P0(M)
182                 ENDIF
183     
184                 IF (FLUID_AT(IJK)) THEN
185     ! Strictly fluid cells: FLAG = 1
186                    IF (MU_S0(M) /= UNDEFINED) THEN
187                       MU_S(IJK,M) = MU_S0(M)
188                       EPMU_S(IJK,M) = MU_S0(M)
189                       LAMBDA_S(IJK,M) = (-2./3.)*MU_S(IJK,M)
190                       EPLAMBDA_S(IJK,M) = (-2./3.)*MU_S(IJK,M)
191                    ENDIF
192                    IF (K_S0(M) /= UNDEFINED) K_S(IJK,M) = K_S0(M)
193                    IF (DIF_S0 /= UNDEFINED) DIF_S(IJK,M,:NMAX(M)) = DIF_S0
194                 ENDIF
195     
196                 IF (USE_MMS) THEN
197                    IF (MU_S0(M) /= UNDEFINED) THEN
198                       MU_S(IJK,M) = MU_S0(M)
199                       EPMU_S(IJK,M) = MU_S0(M)
200                       LAMBDA_S(IJK,M) = (-2./3.)*MU_S(IJK,M)
201                       EPLAMBDA_S(IJK,M) = (-2./3.)*MU_S(IJK,M)
202                    ENDIF
203                    IF (K_S0(M) /= UNDEFINED) K_S(IJK,M) = K_S0(M)
204                    IF (DIF_S0 /= UNDEFINED) DIF_S(IJK,M,:NMAX(M)) = DIF_S0
205                 ENDIF
206     
207     ! set ep_star_array to user input ep_star in all cells.
208                 EP_star_array(ijk) = ep_star
209     ! initializing blending stress parameters
210                 IF(BLENDING_STRESS.AND.TANH_BLEND) THEN
211                    ep_g_blend_start(ijk) = ep_star_array(ijk) * 0.99d0
212                    ep_g_blend_end(ijk)   = ep_star_array(ijk) * 1.01d0
213                 ELSE IF(BLENDING_STRESS.AND.SIGM_BLEND) THEN
214                    ep_g_blend_start(ijk) = ep_star * 0.97d0
215                    ep_g_blend_end(ijk) = ep_star * 1.01d0
216                 ELSE
217                    ep_g_blend_start(ijk) = ep_star_array(ijk)
218                    ep_g_blend_end(ijk)   = ep_star_array(ijk)
219                 ENDIF
220     
221              ENDDO   ! end loop over ijk
222           ENDDO   ! end loop over MMAX
223     
224     
225     ! Initializing parameters needed if a correlation is used to compute
226     ! ep_star: initializing the indexing system.
227           IF(YU_STANDISH .OR. FEDORS_LANDEL) THEN
228              DO M = 1, SMAX
229                 IF(EP_S_MAX(M) == UNDEFINED) EP_S_MAX(M) = ONE-EP_STAR
230              ENDDO
231     
232              IF (.NOT.CALL_DQMOM) THEN
233     
234     ! refer to Syam's dissertation
235                 IF (SMAX == 2) THEN
236                    ep_s_max_ratio(1,2) = ep_s_max(1)/ &
237                       (ep_s_max(1)+(1.-ep_s_max(1))*ep_s_max(2))
238                 ENDIF
239     
240     ! initialize local variables
241                 DO I = 1, SMAX
242                    DP_TMP(I) = D_P0(I)
243                    M_MAX(I) = I
244                 ENDDO
245     
246     ! Rearrange the indices from coarsest particles to finest to be
247     ! used in CALC_ep_star. Done here because it may need to be done
248     ! for auto_restart
249                 DO I = 1, SMAX
250                    DO J = I, SMAX
251                       IF(DP_TMP(I) < DP_TMP(J)) THEN
252                          old_value = DP_TMP(I)
253                          DP_TMP(I) = DP_TMP(J)
254                          DP_TMP(J) = old_value
255                       ENDIF
256                    ENDDO
257                 ENDDO
258     
259                 DO I = 1, SMAX
260                    DO J = 1, SMAX
261                       IF(DP_TMP(I) == D_P0(J) .AND. D_P0(I) .NE. D_P0(J)) THEN
262                          M_MAX(I) = J
263                       ENDIF
264                    ENDDO
265                 ENDDO
266              ENDIF    ! if .not. call_dqmom
267           ELSE   ! if .not. Yu-standish or Fedors-Landel
268              EP_S_MAX(:) = ZERO
269              EP_S_MAX_RATIO(:,:) = ZERO
270              D_P_RATIO(:,:) = ZERO
271              M_MAX(:) = ZERO
272           ENDIF
273     
274           RETURN
275           END SUBROUTINE SET_CONSTPROP
276