File: N:\mfix\model\cartesian_grid\CG_set_bc0.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: CG_SET_BC0                                             C
4     !  Purpose: This module does the initial setting of boundary           C
5     !           conditions for cut cells only                              C
6     !                                                                      C
7     !  Author: Jeff Dietiker                              Date: 01-Jul-09  C
8     !                                                                      C
9     !  Literature/Document References:                                     C
10     !                                                                      C
11     !  Variables referenced: BC_DEFINED, BC_TYPE, BC_DT_0, TIME, BC_Jet_g0,C
12     !                        BC_K_b, BC_K_t, BC_J_s, BC_J_n, BC_I_w,       C
13     !                        BC_I_e, BC_PLANE, BC_EP_g, BC_P_g, BC_T_g,    C
14     !                        BC_T_s,  BC_U_g, BC_V_g, BC_W_g,              C
15     !                        MMAX, BC_ROP_s, BC_U_s, BC_V_s, BC_W_s        C
16     !  Variables modified: BC_TIME, BC_V_g, I, J, K, IJK, EP_g, P_g, T_g,  C
17     !                      T_s, U_g, V_g, W_g, ROP_s, U_s, V_s, W_s,       C
18     !                      M                                               C
19     !                                                                      C
20     !  Local variables: L, IJK1, IJK2, IJK3                                C
21     !                                                                      C
22     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
23     !
24           SUBROUTINE CG_SET_BC0
25     !...Translated by Pacific-Sierra Research VAST-90 2.06G5  12:17:31  12/09/98
26     !...Switches: -xf
27     !
28     !-----------------------------------------------
29     !   M o d u l e s
30     !-----------------------------------------------
31           USE bc
32           USE compar
33           USE cutcell
34           USE eos, ONLY: EOSG
35           USE fldvar
36           USE functions
37           USE funits
38           USE geometry
39           USE indices
40           USE mpi_utility
41           USE param
42           USE param1
43           USE physprop
44           USE quadric
45           USE run
46           USE scalars
47           USE scales
48           USE sendrecv
49           USE toleranc
50     
51           IMPLICIT NONE
52     !-----------------------------------------------
53     !   G l o b a l   P a r a m e t e r s
54     !-----------------------------------------------
55     !-----------------------------------------------
56     !   L o c a l   P a r a m e t e r s
57     !-----------------------------------------------
58     !-----------------------------------------------
59     !   L o c a l   V a r i a b l e s
60     !-----------------------------------------------
61     
62     !
63     !                      Local index for boundary condition
64           INTEGER          L
65     !
66     !                      indices
67           INTEGER          IJK, M, NN ,IJKW,IJKS,IJKB
68     !
69     !----------------------------------------------
70     
71           INTEGER, DIMENSION(8) :: ACCEPTABLE_DEFAULT_WALL=-1
72           LOGICAL :: GLOBAL_CORNER
73     
74     !
75     !  Define global corners as acceptable default walls
76     !  These cells should never be used
77     !
78     
79           IF(.NOT.RE_INDEXING.AND.NumPEs==1) THEN
80     
81           ACCEPTABLE_DEFAULT_WALL(1) = FUNIJK(IMIN3,JMIN3,KMIN3)
82           ACCEPTABLE_DEFAULT_WALL(2) = FUNIJK(IMAX3,JMIN3,KMIN3)
83           ACCEPTABLE_DEFAULT_WALL(3) = FUNIJK(IMIN3,JMAX3,KMIN3)
84           ACCEPTABLE_DEFAULT_WALL(4) = FUNIJK(IMAX3,JMAX3,KMIN3)
85           ACCEPTABLE_DEFAULT_WALL(5) = FUNIJK(IMIN3,JMIN3,KMAX3)
86           ACCEPTABLE_DEFAULT_WALL(6) = FUNIJK(IMAX3,JMIN3,KMAX3)
87           ACCEPTABLE_DEFAULT_WALL(7) = FUNIJK(IMIN3,JMAX3,KMAX3)
88           ACCEPTABLE_DEFAULT_WALL(8) = FUNIJK(IMAX3,JMAX3,KMAX3)
89     
90           ENDIF
91     
92     !      DO N = 1,8
93     !         print*,'acceptable default=',ACCEPTABLE_DEFAULT_WALL(N)
94     !      ENDDO
95     
96     
97           DO IJK = ijkstart3, ijkend3
98     
99              L = BC_ID(IJK)
100     
101              IF(L>0) THEN
102                 IF(BC_TYPE_ENUM(L)==CG_PO) THEN
103     
104                    P_STAR(IJK) = ZERO
105                    P_G(IJK) = SCALE_PRESSURE(BC_P_G(L))
106        !
107                    IF (BC_EP_G(L) /= UNDEFINED) EP_G(IJK) = BC_EP_G(L)
108                    IF (BC_T_G(L) /= UNDEFINED) then
109                       T_G(IJK) = BC_T_G(L)
110                    ELSE
111                       T_g(IJK) = TMIN
112                    ENDIF
113     
114                    NN = 1
115                    IF (NMAX(0) > 0) THEN
116                       WHERE (BC_X_G(L,:NMAX(0)) /= UNDEFINED) X_G(IJK,:&
117                              NMAX(0)) = BC_X_G(L,:NMAX(0))
118                       NN = NMAX(0) + 1
119                    ENDIF
120     
121                    IF (NScalar > 0) THEN
122                       WHERE (BC_Scalar(L,:NScalar) /= UNDEFINED)&
123                       Scalar(IJK,:NScalar) = BC_Scalar(L,:NScalar)
124                    ENDIF
125     
126                    IF (K_Epsilon) THEN
127                       IF (BC_K_Turb_G(L) /= UNDEFINED) K_Turb_G(IJK) = BC_K_Turb_G(L)
128                       IF (BC_E_Turb_G(L) /= UNDEFINED) E_Turb_G(IJK) = BC_E_Turb_G(L)
129                    ENDIF
130     
131                    DO M = 1, MMAX
132                       IF (BC_ROP_S(L,M) /= UNDEFINED) ROP_S(IJK,M) = BC_ROP_S(L,M)
133                       IF(BC_T_S(L,M)/=UNDEFINED)T_S(IJK,M)=BC_T_S(L,M)
134                       IF (BC_THETA_M(L,M) /= UNDEFINED) THETA_M(IJK,M)= BC_THETA_M(L,M)
135                       NN = 1
136                       IF (NMAX(M) > 0) THEN
137                          WHERE (BC_X_S(L,M,:NMAX(M)) /= UNDEFINED) X_S(&
138                                 IJK,M,:NMAX(M)) = BC_X_S(L,M,:NMAX(M))
139                          NN = NMAX(M) + 1
140                       ENDIF
141                    END DO
142     
143                 ELSEIF(BC_TYPE_ENUM(L)==CG_MI) THEN
144     
145                    P_STAR(IJK) = ZERO
146        !
147                    EP_G(IJK) = BC_EP_G(L)
148                    P_G(IJK) = SCALE_PRESSURE(BC_P_G(L))
149                    T_G(IJK) = BC_T_G(L)
150     
151                    IF (NMAX(0) > 0) THEN
152                       X_G(IJK,:NMAX(0)) = BC_X_G(L,:NMAX(0))
153                    ENDIF
154     
155                    IF (NScalar > 0) THEN
156                       Scalar(IJK,:NScalar) = BC_Scalar(L,:NScalar)
157                    ENDIF
158     
159                    IF (K_Epsilon) THEN
160                       K_Turb_G(IJK) = BC_K_Turb_G(L)
161                       E_Turb_G(IJK) = BC_E_Turb_G(L)
162                    ENDIF
163     
164                    DO M = 1, MMAX
165                       ROP_S(IJK,M) = BC_ROP_S(L,M)
166                       T_S(IJK,M) = BC_T_S(L,M)
167                       THETA_M(IJK,M) = BC_THETA_M(L,M)
168     
169                       IF (NMAX(M) > 0) THEN
170                          X_S(IJK,M,:NMAX(M)) = BC_X_S(L,M,:NMAX(M))
171                       ENDIF
172     
173                    END DO
174     
175                    IF(BC_U_g(L)/=UNDEFINED) THEN
176                       U_G(IJK) =  BC_U_g(L)
177                    ELSE
178                       U_G(IJK) =  BC_VELMAG_g(L)*NORMAL_S(IJK,1)
179                    ENDIF
180     
181     
182                    IF(BC_V_g(L)/=UNDEFINED) THEN
183                       V_G(IJK) =  BC_V_g(L)
184                    ELSE
185                       V_G(IJK) =  BC_VELMAG_g(L)*NORMAL_S(IJK,2)
186                    ENDIF
187     
188                    IF(BC_W_g(L)/=UNDEFINED) THEN
189                       W_G(IJK) =  BC_W_g(L)
190                    ELSE
191                       W_G(IJK) =  BC_VELMAG_g(L)*NORMAL_S(IJK,3)
192                    ENDIF
193     
194                    IJKW = WEST_OF(IJK)
195                    IJKS = SOUTH_OF(IJK)
196                    IJKB = BOTTOM_OF(IJK)
197     
198                    IF(FLUID_AT(IJKW)) THEN
199                       IF(BC_U_g(L)/=UNDEFINED) THEN
200                          U_G(IJKW) =  BC_U_g(L)
201                       ELSE
202                          U_G(IJKW) =  BC_VELMAG_g(L)*NORMAL_S(IJK,1)
203                       ENDIF
204                    ENDIF
205     
206                    IF(FLUID_AT(IJKS)) THEN
207                       IF(BC_V_g(L)/=UNDEFINED) THEN
208                          V_G(IJKS) =  BC_V_g(L)
209                       ELSE
210                          V_G(IJKS) =  BC_VELMAG_g(L)*NORMAL_S(IJK,2)
211                       ENDIF
212                    ENDIF
213     
214                    IF(FLUID_AT(IJKB)) THEN
215                       IF(BC_W_g(L)/=UNDEFINED) THEN
216                          W_G(IJKB) =  BC_W_g(L)
217                       ELSE
218                          W_G(IJKB) =  BC_VELMAG_g(L)*NORMAL_S(IJK,3)
219                       ENDIF
220                    ENDIF
221     
222        !
223                    M = 1
224     
225                    DO M=1,MMAX
226     
227                       IF(BC_U_s(L,M)/=UNDEFINED) THEN
228                          U_S(IJK,M) =  BC_U_S(L,M)
229                       ELSE
230                          U_S(IJK,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,1)
231                       ENDIF
232     
233                       IF(BC_V_S(L,M)/=UNDEFINED) THEN
234                          V_S(IJK,M) =  BC_V_S(L,M)
235                       ELSE
236                          V_S(IJK,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,2)
237                       ENDIF
238     
239                       IF(BC_W_S(L,M)/=UNDEFINED) THEN
240                          W_S(IJK,M) =  BC_W_S(L,M)
241                       ELSE
242                          W_S(IJK,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,3)
243                       ENDIF
244     
245                       IJKW = WEST_OF(IJK)
246                       IJKS = SOUTH_OF(IJK)
247                       IJKB = BOTTOM_OF(IJK)
248     
249                       IF(FLUID_AT(IJKW)) THEN
250                          IF(BC_U_S(L,M)/=UNDEFINED) THEN
251                             U_S(IJKW,M) =  BC_U_S(L,M)
252                          ELSE
253                             U_S(IJKW,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,1)
254                          ENDIF
255                       ENDIF
256     
257                       IF(FLUID_AT(IJKS)) THEN
258                          IF(BC_V_S(L,M)/=UNDEFINED) THEN
259                             V_S(IJKS,M) =  BC_V_S(L,M)
260                          ELSE
261                             V_S(IJKS,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,2)
262                          ENDIF
263                       ENDIF
264     
265                       IF(FLUID_AT(IJKB)) THEN
266                          IF(BC_W_S(L,M)/=UNDEFINED) THEN
267                             W_S(IJKB,M) =  BC_W_S(L,M)
268                          ELSE
269                             W_S(IJKB,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,3)
270                          ENDIF
271                       ENDIF
272     
273                    ENDDO
274     
275        !            IF (MMAX > 0) THEN
276        !               U_S(IJK,:MMAX) = BC_U_S(L,:MMAX)
277        !               V_S(IJK,:MMAX) = BC_V_S(L,:MMAX)
278        !               W_S(IJK,:MMAX) = BC_W_S(L,:MMAX)
279        !
280        !               IF(FLUID_AT(IJKW)) THEN
281        !                  U_S(IJKW,:MMAX) = BC_U_S(L,:MMAX)
282        !               ENDIF
283        !
284        !               IF(FLUID_AT(IJKS)) THEN
285        !                  V_S(IJKS,:MMAX) = BC_V_S(L,:MMAX)
286        !               ENDIF
287        !
288        !               IF(FLUID_AT(IJKB)) THEN
289        !                  W_S(IJKB,:MMAX) = BC_W_S(L,:MMAX)
290        !               ENDIF
291        !               M = MMAX + 1
292        !            ENDIF
293     
294     
295                    IF (MW_AVG == UNDEFINED) THEN
296                       MW_MIX_G(IJK) = CALC_MW(X_G,DIMENSION_3,IJK,NMAX(0),MW_G)
297                    ELSE
298                       MW_MIX_G(IJK) = MW_AVG
299                    ENDIF
300                    IF (RO_G0 == UNDEFINED) RO_G(IJK) = EOSG(MW_MIX_G&
301                       (IJK),P_G(IJK),T_G(IJK))
302                    ROP_G(IJK) = EP_G(IJK)*RO_G(IJK)
303     
304                 ENDIF
305              ENDIF
306     
307              IF(DEFAULT_WALL_AT(IJK)) THEN
308     
309     !            print*,'Default_wall_at IJK=',IJK,I_OF(IJK),J_OF(IJK),K_OF(IJK)
310     
311                 GLOBAL_CORNER = .FALSE.
312                 DO NN = 1,8
313                    IF(IJK==ACCEPTABLE_DEFAULT_WALL(NN)) GLOBAL_CORNER = .TRUE.
314                 ENDDO
315     
316                 IF(.NOT.GLOBAL_CORNER.AND..NOT.BLOCKED_CELL_AT(IJK)) THEN
317     
318                    ICBC_FLAG(IJK)(2:3) = 'CG'
319     
320                    IF((MyPE == PE_IO).AND.PRINT_WARNINGS) THEN
321                       WRITE(*,*) 'WARNING: DEFAULT WALL DETECTED AT I,J,K = ',I_OF(IJK),J_OF(IJK),K_OF(IJK) ,BLOCKED_CELL_AT(IJK)
322                       WRITE(*,*) '         WHEN USING CARTESIAN GRID CUT-CELL FEATURE.'
323                       WRITE(*,*) '         DEFAULT WALLS ARE NOT ALLOWED WITH CUT-CELLS.'
324                       WRITE(*,*) '         THE DEFAULT WALL WAS REMOVED ALONG THIS CELL.'
325                       WRITE(*,*) ''
326                    ENDIF
327     !               CALL MFIX_EXIT(MYPE)
328     
329                 ENDIF
330     
331              ENDIF
332     
333           ENDDO
334     
335           RETURN
336           END SUBROUTINE CG_SET_BC0
337     
338     !// Comments on the modifications for DMP version implementation
339     !// 001 Include header file and common declarations for parallelization
340     !// 020 New local variables for parallelization: FLAG_G , FLUID_AT_G
341     !// 360 Check if i,j,k resides on current processor
342