File: /nfs/home/0/users/jenkins/mfix.git/model/GhdTheory/adjust_eps_ghd.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: ADJUST_EPS_GHD                                         C
4     !  Purpose: Eliminate the solids phases that occupy only very small    C
5     !           fractions of the computational cell volume for GHD theory  C
6     !                                                                      C
7     !  Author: M. Syamlal                                 Date: 28-JAN-92  C
8     !  Reviewer: P. Nicoletti, W. Rogers, S. Venkatesan   Date: 29-JAN-92  C
9     !                                                                      C
10     !  Revision Number: 1                                                  C
11     !  Purpose: Include CONSTANTS.INC                                      C
12     !  Author: M. Syamlal                                 Date: 7-FEB-92   C
13     !  Reviewer: P. Nicoletti                             Date: 11-DEC-92  C
14     !                                                                      C
15     !  Literature/Document References: None                                C
16     !                                                                      C
17     !  Variables referenced: FLAG, RO_s, ROP_s, U_s, V_s, W_s, EP_g        C
18     !                        IMAX1, JMAX1, KMAX1, MMAX, IMIN1, JMIN1, KMIN1C
19     !  Variables modified: ROP_s, U_s, V_s, W_s, EP_g, ROP_g, I, J, K, IJK,C
20     !                      RO_s                                            C
21     !                                                                      C
22     !  Local variables: EPSUM                                              C
23     !                                                                      C
24     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
25     !
26           SUBROUTINE ADJUST_EPS_GHD
27     !...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
28     !...Switches: -xf
29     !
30     !-----------------------------------------------
31     !   M o d u l e s
32     !-----------------------------------------------
33           USE param
34           USE param1
35           USE toleranc
36           USE constant
37           USE fldvar
38           USE geometry
39           USE indices
40           USE physprop
41           USE run
42           USE compar
43           USE sendrecv
44           USE ghdtheory
45           USE fun_avg
46           USE functions
47           IMPLICIT NONE
48     !-----------------------------------------------
49     !   G l o b a l   P a r a m e t e r s
50     !-----------------------------------------------
51     !-----------------------------------------------
52     !   L o c a l   P a r a m e t e r s
53     !-----------------------------------------------
54     !-----------------------------------------------
55     !   L o c a l   V a r i a b l e s
56     !-----------------------------------------------
57     !
58     !                      Indices
59           INTEGER          I, J, K, IJK, IJKE, IJKN, IJKT
60     !
61     !                      Solids phase
62           INTEGER          M
63     !                      Sum of (very small) solids volume fractions that
64     !                      are set to zero.
65           DOUBLE PRECISION epsMixE, epsMixN, epsMixT, epSolid, epSolidE(smax), epSolidN(smax), epSolidT(smax)
66           LOGICAL          DiluteCellE, DiluteCellN, DiluteCellT
67     !
68     ! First set solids cell-center density in very dilute cells to zero
69     !
70           DO IJK = ijkstart3, ijkend3
71                    IF (FLUID_AT(IJK)) THEN
72                       DO M = 1, SMAX
73                          epSolid = ROP_S(IJK,M)/RO_S(IJK,M)
74                          IF (epSolid < ZERO_EP_S) THEN
75     
76     !  Remove solids in very small quantities and set solids velocity to zero
77     !  if there is outflow from the present cell.
78     
79                             ROP_S(IJK,MMAX) = ROP_S(IJK,MMAX) - ROP_S(IJK,M) ! mmax = mixture phase
80                             EP_G(IJK) = EP_G(IJK) + epSolid
81                             ROP_G(IJK) = RO_G(IJK)*EP_G(IJK)
82                             ROP_S(IJK,M) = ZERO
83                          ENDIF
84                       END DO
85                    ENDIF
86           END DO
87     !
88     ! now set velocities at cell faces to zero
89     !
90           DO IJK = ijkstart3, ijkend3
91                    I = I_OF(IJK)
92                    J = J_OF(IJK)
93                    K = K_OF(IJK)
94                    IJKE = EAST_OF(IJK)
95                    IJKN = NORTH_OF(IJK)
96                    IJKT = TOP_OF(IJK)
97                    IF (FLUID_AT(IJK)) THEN
98                       epsMixE = ZERO
99                       epsMixN = ZERO
100                       epsMixT = ZERO
101                       DiluteCellE = .FALSE.
102                       DiluteCellN = .FALSE.
103                       DiluteCellT = .FALSE.
104                       DO M = 1, SMAX
105                          epSolidE(M) = AVG_X(EP_S(IJK,M),EP_S(IJKE,M),I)
106                          epSolidN(M) = AVG_Y(EP_S(IJK,M),EP_S(IJKN,M),J)
107                          epSolidT(M) = AVG_Z(EP_S(IJK,M),EP_S(IJKT,M),K)
108     
109                          epsMixE = epsMixE + epSolidE(M)*RO_S(IJK,M)
110                          epsMixN = epsMixN + epSolidN(M)*RO_S(IJK,M)
111                          epsMixT = epsMixT + epSolidT(M)*RO_S(IJK,M)
112     
113                          IF (epSolidE(M) < ZERO_EP_S) THEN
114                             epsMixE = epsMixE - epSolidE(M)*RO_S(IJK,M)
115                             epSolidE(M) = ZERO
116                             U_S(IJK,M) = ZERO
117                             JoiX(IJK,M) = ZERO
118                             DiluteCellE = .TRUE. ! to make sure at least one species is very dilute before
119                                                  ! re-calculating mixture velocity
120                          ENDIF
121     
122                          IF (epSolidN(M) < ZERO_EP_S) THEN
123                             epsMixN = epsMixN - epSolidN(M)*RO_S(IJK,M)
124                             epSolidN(M) = ZERO
125                             V_S(IJK,M) = ZERO
126                             JoiY(IJK,M) = ZERO
127                             DiluteCellN = .TRUE.
128                          ENDIF
129     
130                          IF (epSolidT(M) < ZERO_EP_S) THEN
131                             epsMixT = epsMixT - epSolidT(M)*RO_S(IJK,M)
132                             epSolidT(M) = ZERO
133                             W_S(IJK,M) = ZERO
134                             JoiZ(IJK,M) = ZERO
135                             DiluteCellT = .TRUE.
136                          ENDIF
137                       END DO
138     !
139     ! compute corrected mixture velcoity and species mass fluxes based on their definition.
140     !
141                          IF (epsMixE > ZERO .AND. DiluteCellE) THEN
142                            U_S(IJK,MMAX) = ZERO
143                            DO M = 1, SMAX
144                              U_S(IJK,MMAX) = U_S(IJK,MMAX) + U_S(IJK,M)*epSolidE(M)*RO_S(IJK,M)
145                            ENDDO
146                            U_S(IJK,MMAX) = U_S(IJK,MMAX)/epsMixE
147                            DO M = 1, SMAX
148                              JoiX(IJK,M) = epSolidE(M)*RO_S(IJK,M)*(U_S(IJK,M)-U_S(IJK,MMAX))
149                            ENDDO
150                          ELSEIF(epsMixE == ZERO) THEN
151                            U_S(IJK,MMAX) = ZERO
152                          ENDIF
153     
154                          IF (epsMixN > ZERO .AND. DiluteCellN) THEN
155                            V_S(IJK,MMAX) = ZERO
156                            DO M = 1, SMAX
157                              V_S(IJK,MMAX) = V_S(IJK,MMAX) + V_S(IJK,M)*epSolidN(M)*RO_S(IJK,M)
158                            ENDDO
159                            V_S(IJK,MMAX) = V_S(IJK,MMAX)/epsMixN
160                            DO M = 1, SMAX
161                              JoiY(IJK,M) = epSolidN(M)*RO_S(IJK,M)*(V_S(IJK,M)-V_S(IJK,MMAX))
162                            ENDDO
163                          ELSEIF(epsMixN == ZERO) THEN
164                            V_S(IJK,MMAX) = ZERO
165                          ENDIF
166     
167                          IF(.NOT.NO_K) THEN  ! for 3D cases only
168                            IF (epsMixT > ZERO .AND. DiluteCellT) THEN
169                              W_S(IJK,MMAX) = ZERO
170                              DO M = 1, SMAX
171                                W_S(IJK,MMAX) = W_S(IJK,MMAX) + W_S(IJK,M)*epSolidT(M)*RO_S(IJK,M)
172                              ENDDO
173                              W_S(IJK,MMAX) = W_S(IJK,MMAX)/epsMixT
174                              DO M = 1, SMAX
175                                JoiZ(IJK,M) = epSolidT(M)*RO_S(IJK,M)*(W_S(IJK,M)-W_S(IJK,MMAX))
176                              ENDDO
177                            ELSEIF(epsMixT == ZERO) THEN
178                              W_S(IJK,MMAX) = ZERO
179                            ENDIF
180                          ENDIF
181                    ENDIF
182           END DO
183     
184           RETURN
185           END SUBROUTINE ADJUST_EPS_GHD
186     
187     !// Comments on the modifications for DMP version implementation
188     !// 001 Include header file and common declarations for parallelization
189     !// 200 1119 Changed the limits for the triple loop, do k=kmin1,kmax1=>kstart1,kend1
190     !// 400 Added sendrecv module and send_recv calls for COMMunication
191