File: N:\mfix\model\usr_prop_mod.f

1     ! or put the actual routine calc_usr_source here? move logic variable
2     !  usr_source into run or something...
3           module usr_prop
4           use param, only: dim_m
5     
6     ! If .TRUE. call user-defined physical properties
7           LOGICAL :: USR_ROG
8           LOGICAL :: USR_CPG 
9           LOGICAL :: USR_ROS(DIM_M)
10           LOGICAL :: USR_CPS(DIM_M)
11     
12     ! If .TRUE. call user_defined transport properties
13           LOGICAL :: USR_MUG
14           LOGICAL :: USR_KG
15           LOGICAL :: USR_DIFG
16           LOGICAL :: USR_MUS(DIM_M)
17           LOGICAL :: USR_KS(DIM_M)
18           LOGICAL :: USR_DIFS(DIM_M)
19     
20     ! Interphase momentum exchange coefficients (F_GS, F_SS)
21           LOGICAL :: USR_FGS(DIM_M)
22           LOGICAL :: USR_FSS( (DIM_M*(DIM_M-1)/2)+1 )
23     ! Interphase heat transfer coefficient (GAMA)
24           LOGICAL :: USR_GAMA(DIM_M)
25     
26           ENUM, BIND(C)
27              ENUMERATOR :: Gas_Density, Solids_Density
28              ENUMERATOR :: Gas_SpecificHeat, Solids_SpecificHeat
29              ENUMERATOR :: Gas_Viscosity, Solids_Viscosity
30              ENUMERATOR :: Gas_Conductivity, Solids_Conductivity
31              ENUMERATOR :: Gas_Diffusivity, Solids_Diffusivity 
32              ENUMERATOR :: GasSolids_Drag, SolidsSolids_Drag
33              ENUMERATOR :: GasSolids_HeatTransfer
34              ENUMERATOR :: BLANK
35              ENUMERATOR :: DUMMY
36           END ENUM
37     
38     
39           CONTAINS
40     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
41     !                                                                      C
42     !  Subroutine: CALC_USR_PROP                                           C 
43     !  Purpose: Driver routine to calculate user defined physical,         C
44     !  transport, exchange terms that appear in the governing equations.   C
45     !                                                                      C
46     !                                                                      C
47     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
48           SUBROUTINE CALC_USR_PROP(lprop, lM, lL, lerr)
49     
50     ! Modules
51     !-----------------------------------------------
52           use compar, only: ijkstart3, ijkend3, myPE, numPEs
53           use functions, only: fluid_at, wall_at, fluidorp_flow_at
54           use output, only: REPORT_NEG_DENSITY
55           use output, only: REPORT_NEG_SPECIFICHEAT
56           use param1, only: undefined_i, zero
57           use physprop, only: nmax
58           use physprop, only: k_g, dif_g, c_pg, mu_g
59           use physprop, only: k_s, dif_s, c_ps
60           use fldvar, only: ro_g, rop_g, ep_g
61           use fldvar, only: ro_s, p_s
62           use visc_g, only: mu_gt, lambda_gt
63           use visc_s, only: mu_s, lambda_s
64           use utilities, only: mfix_isnan
65     
66           use error_manager
67     
68           IMPLICIT NONE
69     ! Dummy arguments
70     !-----------------------------------------------
71     ! User properties available
72     !     ROg, ROs (gas and solids density)
73     !     CPg, CPs (gas and solids specific heat)
74     !     Mug, Mus (gas and solids viscosity)
75     !     Kg, Ks (gas and solids conductivity)
76     !     Difg, Difs (gas and solids diffusivity)
77     !     Fgs, Fss (gas-solids and solids-solids drag)
78     !     Gama (gas-solids heat transfere coefficient)
79           INTEGER, INTENT(IN) :: lprop
80     
81     ! Phase index 
82           INTEGER, OPTIONAL, INTENT(IN) :: lM
83     
84     ! Phase index 
85           INTEGER, OPTIONAL, INTENT(IN) :: lL
86     
87     ! Error index
88     ! Arrays for storing errors:
89     ! 100 - Negative gas phase density
90     ! 101 - Negative solids phase density
91     ! 10x - Unclassified
92           INTEGER, OPTIONAL, INTENT(INOUT) :: lerr(0:numPEs-1)
93     
94     ! Local variables
95     !-----------------------------------------------
96     ! indices
97           INTEGER :: IJK
98     ! tmp indices
99           INTEGER :: M, L, N
100     ! Flag to write log header
101           LOGICAL :: wHeader
102     
103     !-----------------------------------------------
104     
105     ! set local values for phase index
106           IF (.NOT.present(lM)) THEN
107              M = UNDEFINED_I
108           ELSE
109              M = lM
110           ENDIF
111     
112     ! set local values for phase index
113           IF (.NOT.present(lL)) THEN
114              L = UNDEFINED_I
115           ELSE
116              L = lL
117           ENDIF
118              
119           N = UNDEFINED_I
120     
121           SELECT CASE(lProp)
122           
123     ! Gas physical properties: density
124           CASE(Gas_Density)
125              wHeader = .TRUE.
126              DO IJK=IJKSTART3,IJKEND3
127                 IF (WALL_AT(IJK)) CYCLE
128     !            CALL USR_PROPERTIES(lprop, IJK, 0, N)
129                 CALL USR_PROP_ROG(IJK)
130                 ROP_G(IJK) = EP_G(IJK)*RO_G(IJK)
131     ! report errors
132                 IF(RO_G(IJK) < ZERO .or. mfix_isnan(RO_G(IJK))) THEN
133                    lErr(myPE) = 100  ! return appropritae error flag
134                    IF(REPORT_NEG_DENSITY) CALL ROgErr_LOG(IJK, wHeader)
135                 ENDIF
136     
137              ENDDO
138     
139     ! Gas physical properties: specific heat
140           CASE(Gas_SpecificHeat)
141              DO IJK=IJKSTART3,IJKEND3
142                 IF (WALL_AT(IJK)) CYCLE
143     !            CALL USR_PROPERTIES(lprop, IJK, 0, N)
144                 CALL USR_PROP_CPG(IJK)
145     ! report errors.
146                 IF(C_PG(IJK) <= ZERO) THEN
147     !               lErr(myPE) = 103
148                    IF(REPORT_NEG_SPECIFICHEAT) CALL CPgErr_LOG(IJK, wHeader)
149                 ENDIF
150              ENDDO
151     
152     ! Gas transport properties: viscosity
153           CASE(Gas_Viscosity)
154              DO IJK=IJKSTART3,IJKEND3
155                 IF (FLUID_AT(IJK)) THEN
156     !               CALL USR_PROPERTIES(lprop, IJK, 0, N)
157                    CALL USR_PROP_MUG(IJK)
158                    MU_GT(IJK) = MU_G(IJK)
159                    LAMBDA_GT(IJK) = -2./3.*MU_GT(IJK)
160                 ELSE ! probably unnecessary
161                    MU_g(IJK) = zero
162                    MU_GT(IJK) = ZERO
163                    LAMBDA_GT(IJK) = ZERO
164                 ENDIF
165              ENDDO
166     
167     ! Gas transport properties: conductivity
168           CASE(Gas_Conductivity)
169              DO IJK=IJKSTART3,IJKEND3
170                 IF (FLUID_AT(IJK)) THEN
171     !               CALL USR_PROPERTIES(lprop, IJK, 0, N)
172                    CALL USR_PROP_KG(IJK)
173                 ELSE ! probably unnecessary
174                    K_G(IJK) = ZERO
175                 ENDIF
176              ENDDO
177     
178     ! Gas transport properties: diffusivity
179           CASE(Gas_Diffusivity)
180              DO N=1,NMAX(0)
181                 DO IJK=IJKSTART3,IJKEND3
182                    IF (FLUID_AT(IJK)) THEN
183     !                  CALL USR_PROPERTIES(lprop, IJK, 0, N)
184                       CALL USR_PROP_DIFG(IJK,N)
185                    ELSE ! probably unnecessary
186                       DIF_G(IJK,N) = ZERO
187                    ENDIF
188                 ENDDO
189              ENDDO
190      
191     ! Solids physical properties: density
192           CASE(Solids_Density)
193              wHeader = .TRUE.
194              DO IJK=IJKSTART3,IJKEND3
195                 IF (WALL_AT(IJK)) CYCLE
196     !            CALL USR_PROPERTIES(lprop, IJK, M, N)
197                 CALL USR_PROP_ROS(IJK,M)
198     ! report errors.
199                 IF(RO_S(IJK,M) <= ZERO) THEN
200                    lErr(myPE) = 101
201                    IF(REPORT_NEG_DENSITY) CALL ROsErr_LOG(IJK, M, wHeader)
202                 ENDIF
203     
204              ENDDO
205     
206     ! Solids physical properties: specific heat
207           CASE(Solids_SpecificHeat)
208              DO IJK=IJKSTART3,IJKEND3
209                 IF (WALL_AT(IJK)) CYCLE
210     !            CALL USR_PROPERTIES(lprop, IJK, M, N)
211                 CALL USR_PROP_CPS(IJK,M)
212     ! report errors.
213                 IF(C_PS(IJK,M) <= ZERO) THEN
214     !               lErr(myPE) = 104
215                    IF(REPORT_NEG_SPECIFICHEAT) CALL CPsErr_LOG(IJK, M, wHeader)
216                 ENDIF
217              ENDDO
218     
219     ! Solids transport properties: viscosity
220           CASE(Solids_Viscosity)
221              DO IJK=IJKSTART3,IJKEND3
222                 IF (FLUID_AT(IJK)) THEN
223     !               CALL USR_PROPERTIES(lprop, IJK, M, N)
224                    CALL USR_PROP_MUS(IJK,M)
225                 ELSE ! probably unnecessary
226                    MU_S(IJK,M) = ZERO
227                    LAMBDA_S(IJK,M) = ZERO
228                 ENDIF
229              ENDDO
230     
231     ! Solids transport properties: conductivity
232           CASE(Solids_Conductivity)
233              DO IJK=IJKSTART3,IJKEND3
234                 IF (FLUID_AT(IJK)) THEN
235     !               CALL USR_PROPERTIES(lprop, IJK, M, N)
236                    CALL USR_PROP_KS(IJK,M)
237                 ELSE ! probably unnecessary
238                    K_S(IJK,M) = ZERO
239                 ENDIF
240              ENDDO
241     
242     ! Solids transport properties: diffusivity
243           CASE(Solids_Diffusivity)
244              DO N=1,NMAX(M)
245                 DO IJK=IJKSTART3,IJKEND3
246                    IF (FLUID_AT(IJK)) THEN
247     !                  CALL USR_PROPERTIES(lprop, IJK, M, N)
248                       CALL USR_PROP_DIFS(IJK,M,N)
249                    ELSE ! probably unnecessary
250                       DIF_S(IJK,M,N) = ZERO
251                    ENDIF
252                 ENDDO
253              ENDDO
254     
255     ! Gas-solids exchange: Interphase heat transfer coefficient (GAMA)
256           CASE(GasSolids_HeatTransfer)
257              DO IJK=IJKSTART3,IJKEND3
258                 IF (FLUID_AT(IJK)) THEN
259     !               CALL USR_PROPERTIES(lprop, IJK, M, N)
260                    CALL USR_PROP_GAMA(IJK,M)
261                 ENDIF
262              ENDDO
263     
264     ! Gas-solids exchange: interphase drag
265           CASE(GasSolids_Drag)
266              DO IJK=IJKSTART3,IJKEND3
267                 IF (FLUIDorP_FLOW_AT(IJK)) THEN
268                    CALL USR_PROPERTIES(lprop, IJK, M, N)
269     ! this hook is not accessed or functioning.
270     ! use drag_type = usr_drag to access this entity
271                 ENDIF
272              ENDDO
273     
274     ! Solids-solids exchange: interphase drag
275           CASE(SolidsSolids_Drag)
276              DO IJK=IJKSTART3,IJKEND3
277                 IF (WALL_AT(IJK)) CYCLE
278     !            CALL USR_PROPERTIES(lprop, IJK, M, L)
279                 CALL USR_PROP_FSS(IJK,L,M)
280              ENDDO
281     
282     ! error out
283           CASE DEFAULT
284     ! should never hit this
285     ! Initialize the error manager.
286              CALL INIT_ERR_MSG("CALC_USR_PROP")
287              WRITE(ERR_MSG, 1001) ival(lprop)
288              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
289      1001 FORMAT('Error 1101: Unknown Property= ', A)
290              CALL FINL_ERR_MSG
291           END SELECT   ! end selection of user property 
292      
293           RETURN
294           END SUBROUTINE CALC_USR_PROP
295     
296           end module usr_prop
297     
298