File: /nfs/home/0/users/jenkins/mfix.git/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 boundfunijk
33           USE compar
34           USE cutcell
35           USE eos, ONLY: EOSG
36           USE fldvar
37           USE functions
38           USE funits
39           USE geometry
40           USE indices
41           USE mpi_utility
42           USE param
43           USE param1
44           USE physprop
45           USE quadric
46           USE run
47           USE scalars
48           USE scales
49           USE sendrecv
50           USE toleranc
51     
52           IMPLICIT NONE
53     !-----------------------------------------------
54     !   G l o b a l   P a r a m e t e r s
55     !-----------------------------------------------
56     !-----------------------------------------------
57     !   L o c a l   P a r a m e t e r s
58     !-----------------------------------------------
59     !-----------------------------------------------
60     !   L o c a l   V a r i a b l e s
61     !-----------------------------------------------
62     
63     !
64     !                      Local index for boundary condition
65           INTEGER          L
66     !
67     !                      indices
68           INTEGER          IJK, M, N ,IJKW,IJKS,IJKB
69     !
70     !                      Dummy variable for gas pressure
71           DOUBLE PRECISION PJ
72     !
73     !-----------------------------------------------
74     !   E x t e r n a l   F u n c t i o n s
75     !-----------------------------------------------
76           DOUBLE PRECISION , EXTERNAL :: CALC_MW
77     !-----------------------------------------------
78     !----------------------------------------------
79     
80           INTEGER, DIMENSION(8) :: ACCEPTABLE_DEFAULT_WALL=-1
81           LOGICAL :: GLOBAL_CORNER
82     
83     !
84     !  Define global corners as acceptable default walls
85     !  These cells should never be used
86     !
87     
88           IF(.NOT.RE_INDEXING.AND.NumPEs==1) THEN
89     
90           ACCEPTABLE_DEFAULT_WALL(1) = FUNIJK(IMIN3,JMIN3,KMIN3)
91           ACCEPTABLE_DEFAULT_WALL(2) = FUNIJK(IMAX3,JMIN3,KMIN3)
92           ACCEPTABLE_DEFAULT_WALL(3) = FUNIJK(IMIN3,JMAX3,KMIN3)
93           ACCEPTABLE_DEFAULT_WALL(4) = FUNIJK(IMAX3,JMAX3,KMIN3)
94           ACCEPTABLE_DEFAULT_WALL(5) = FUNIJK(IMIN3,JMIN3,KMAX3)
95           ACCEPTABLE_DEFAULT_WALL(6) = FUNIJK(IMAX3,JMIN3,KMAX3)
96           ACCEPTABLE_DEFAULT_WALL(7) = FUNIJK(IMIN3,JMAX3,KMAX3)
97           ACCEPTABLE_DEFAULT_WALL(8) = FUNIJK(IMAX3,JMAX3,KMAX3)
98     
99           ENDIF
100     
101     !      DO N = 1,8
102     !         print*,'acceptable default=',ACCEPTABLE_DEFAULT_WALL(N)
103     !      ENDDO
104     
105     
106           DO IJK = ijkstart3, ijkend3
107     
108              L = BC_ID(IJK)
109     
110              IF(L>0) THEN
111                 IF(BC_TYPE(L)=='CG_PO') THEN
112     
113                    P_STAR(IJK) = ZERO
114                    P_G(IJK) = SCALE(BC_P_G(L))
115        !
116                    IF (BC_EP_G(L) /= UNDEFINED) EP_G(IJK) = BC_EP_G(L)
117                    IF (BC_T_G(L) /= UNDEFINED) then
118                       T_G(IJK) = BC_T_G(L)
119                    ELSE
120                       T_g(IJK) = TMIN
121                    ENDIF
122     
123                    N = 1
124                    IF (NMAX(0) > 0) THEN
125                       WHERE (BC_X_G(L,:NMAX(0)) /= UNDEFINED) X_G(IJK,:&
126                              NMAX(0)) = BC_X_G(L,:NMAX(0))
127                       N = NMAX(0) + 1
128                    ENDIF
129     
130                    IF (NScalar > 0) THEN
131                       WHERE (BC_Scalar(L,:NScalar) /= UNDEFINED)&
132                       Scalar(IJK,:NScalar) = BC_Scalar(L,:NScalar)
133                    ENDIF
134     
135                    IF (K_Epsilon) THEN
136                       IF (BC_K_Turb_G(L) /= UNDEFINED) K_Turb_G(IJK) = BC_K_Turb_G(L)
137                       IF (BC_E_Turb_G(L) /= UNDEFINED) E_Turb_G(IJK) = BC_E_Turb_G(L)
138                    ENDIF
139     
140                    DO M = 1, MMAX
141                       IF (BC_ROP_S(L,M) /= UNDEFINED) ROP_S(IJK,M) = BC_ROP_S(L,M)
142                       IF(BC_T_S(L,M)/=UNDEFINED)T_S(IJK,M)=BC_T_S(L,M)
143                       IF (BC_THETA_M(L,M) /= UNDEFINED) THETA_M(IJK,M)= BC_THETA_M(L,M)
144                       N = 1
145                       IF (NMAX(M) > 0) THEN
146                          WHERE (BC_X_S(L,M,:NMAX(M)) /= UNDEFINED) X_S(&
147                                 IJK,M,:NMAX(M)) = BC_X_S(L,M,:NMAX(M))
148                          N = NMAX(M) + 1
149                       ENDIF
150                    END DO
151     
152                 ELSEIF(BC_TYPE(L)=='CG_MI') THEN
153     
154                    P_STAR(IJK) = ZERO
155        !
156                    EP_G(IJK) = BC_EP_G(L)
157                    P_G(IJK) = SCALE(BC_P_G(L))
158                    T_G(IJK) = BC_T_G(L)
159     
160                    IF (NMAX(0) > 0) THEN
161                       X_G(IJK,:NMAX(0)) = BC_X_G(L,:NMAX(0))
162                    ENDIF
163     
164                    IF (NScalar > 0) THEN
165                       Scalar(IJK,:NScalar) = BC_Scalar(L,:NScalar)
166                    ENDIF
167     
168                    IF (K_Epsilon) THEN
169                       K_Turb_G(IJK) = BC_K_Turb_G(L)
170                       E_Turb_G(IJK) = BC_E_Turb_G(L)
171                    ENDIF
172     
173                    DO M = 1, MMAX
174                       ROP_S(IJK,M) = BC_ROP_S(L,M)
175                       T_S(IJK,M) = BC_T_S(L,M)
176                       THETA_M(IJK,M) = BC_THETA_M(L,M)
177     
178                       IF (NMAX(M) > 0) THEN
179                          X_S(IJK,M,:NMAX(M)) = BC_X_S(L,M,:NMAX(M))
180                       ENDIF
181     
182                    END DO
183     
184                    IF(BC_U_g(L)/=UNDEFINED) THEN
185                       U_G(IJK) =  BC_U_g(L)
186                    ELSE
187                       U_G(IJK) =  BC_VELMAG_g(L)*NORMAL_S(IJK,1)
188                    ENDIF
189     
190     
191                    IF(BC_V_g(L)/=UNDEFINED) THEN
192                       V_G(IJK) =  BC_V_g(L)
193                    ELSE
194                       V_G(IJK) =  BC_VELMAG_g(L)*NORMAL_S(IJK,2)
195                    ENDIF
196     
197                    IF(BC_W_g(L)/=UNDEFINED) THEN
198                       W_G(IJK) =  BC_W_g(L)
199                    ELSE
200                       W_G(IJK) =  BC_VELMAG_g(L)*NORMAL_S(IJK,3)
201                    ENDIF
202     
203                    IJKW = WEST_OF(IJK)
204                    IJKS = SOUTH_OF(IJK)
205                    IJKB = BOTTOM_OF(IJK)
206     
207                    IF(FLUID_AT(IJKW)) THEN
208                       IF(BC_U_g(L)/=UNDEFINED) THEN
209                          U_G(IJKW) =  BC_U_g(L)
210                       ELSE
211                          U_G(IJKW) =  BC_VELMAG_g(L)*NORMAL_S(IJK,1)
212                       ENDIF
213                    ENDIF
214     
215                    IF(FLUID_AT(IJKS)) THEN
216                       IF(BC_V_g(L)/=UNDEFINED) THEN
217                          V_G(IJKS) =  BC_V_g(L)
218                       ELSE
219                          V_G(IJKS) =  BC_VELMAG_g(L)*NORMAL_S(IJK,2)
220                       ENDIF
221                    ENDIF
222     
223                    IF(FLUID_AT(IJKB)) THEN
224                       IF(BC_W_g(L)/=UNDEFINED) THEN
225                          W_G(IJKB) =  BC_W_g(L)
226                       ELSE
227                          W_G(IJKB) =  BC_VELMAG_g(L)*NORMAL_S(IJK,3)
228                       ENDIF
229                    ENDIF
230     
231        !
232                    M = 1
233     
234                    DO M=1,MMAX
235     
236                       IF(BC_U_s(L,M)/=UNDEFINED) THEN
237                          U_S(IJK,M) =  BC_U_S(L,M)
238                       ELSE
239                          U_S(IJK,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,1)
240                       ENDIF
241     
242                       IF(BC_V_S(L,M)/=UNDEFINED) THEN
243                          V_S(IJK,M) =  BC_V_S(L,M)
244                       ELSE
245                          V_S(IJK,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,2)
246                       ENDIF
247     
248                       IF(BC_W_S(L,M)/=UNDEFINED) THEN
249                          W_S(IJK,M) =  BC_W_S(L,M)
250                       ELSE
251                          W_S(IJK,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,3)
252                       ENDIF
253     
254                       IJKW = WEST_OF(IJK)
255                       IJKS = SOUTH_OF(IJK)
256                       IJKB = BOTTOM_OF(IJK)
257     
258                       IF(FLUID_AT(IJKW)) THEN
259                          IF(BC_U_S(L,M)/=UNDEFINED) THEN
260                             U_S(IJKW,M) =  BC_U_S(L,M)
261                          ELSE
262                             U_S(IJKW,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,1)
263                          ENDIF
264                       ENDIF
265     
266                       IF(FLUID_AT(IJKS)) THEN
267                          IF(BC_V_S(L,M)/=UNDEFINED) THEN
268                             V_S(IJKS,M) =  BC_V_S(L,M)
269                          ELSE
270                             V_S(IJKS,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,2)
271                          ENDIF
272                       ENDIF
273     
274                       IF(FLUID_AT(IJKB)) THEN
275                          IF(BC_W_S(L,M)/=UNDEFINED) THEN
276                             W_S(IJKB,M) =  BC_W_S(L,M)
277                          ELSE
278                             W_S(IJKB,M) =  BC_VELMAG_S(L,M)*NORMAL_S(IJK,3)
279                          ENDIF
280                       ENDIF
281     
282                    ENDDO
283     
284        !            IF (MMAX > 0) THEN
285        !               U_S(IJK,:MMAX) = BC_U_S(L,:MMAX)
286        !               V_S(IJK,:MMAX) = BC_V_S(L,:MMAX)
287        !               W_S(IJK,:MMAX) = BC_W_S(L,:MMAX)
288        !
289        !               IF(FLUID_AT(IJKW)) THEN
290        !                  U_S(IJKW,:MMAX) = BC_U_S(L,:MMAX)
291        !               ENDIF
292        !
293        !               IF(FLUID_AT(IJKS)) THEN
294        !                  V_S(IJKS,:MMAX) = BC_V_S(L,:MMAX)
295        !               ENDIF
296        !
297        !               IF(FLUID_AT(IJKB)) THEN
298        !                  W_S(IJKB,:MMAX) = BC_W_S(L,:MMAX)
299        !               ENDIF
300        !               M = MMAX + 1
301        !            ENDIF
302     
303     
304                    IF (MW_AVG == UNDEFINED) THEN
305                       MW_MIX_G(IJK) = CALC_MW(X_G,DIMENSION_3,IJK,NMAX(0),MW_G)
306                    ELSE
307                       MW_MIX_G(IJK) = MW_AVG
308                    ENDIF
309                    IF (RO_G0 == UNDEFINED) RO_G(IJK) = EOSG(MW_MIX_G&
310                       (IJK),P_G(IJK),T_G(IJK))
311                    ROP_G(IJK) = EP_G(IJK)*RO_G(IJK)
312     
313     
314                 ENDIF
315              ENDIF
316     
317     
318     
319     
320              IF(DEFAULT_WALL_AT(IJK)) THEN
321     
322     !            print*,'Default_wall_at IJK=',IJK,I_OF(IJK),J_OF(IJK),K_OF(IJK)
323     
324                 GLOBAL_CORNER = .FALSE.
325                 DO N = 1,8
326                    IF(IJK==ACCEPTABLE_DEFAULT_WALL(N)) GLOBAL_CORNER = .TRUE.
327                 ENDDO
328     
329                 IF(.NOT.GLOBAL_CORNER.AND..NOT.BLOCKED_CELL_AT(IJK)) THEN
330     
331                    ICBC_FLAG(IJK)(2:3) = 'CG'
332     
333                    IF((MyPE == PE_IO).AND.PRINT_WARNINGS) THEN
334                       WRITE(*,*) 'WARNING: DEFAULT WALL DETECTED AT I,J,K = ',I_OF(IJK),J_OF(IJK),K_OF(IJK) ,BLOCKED_CELL_AT(IJK)
335                       WRITE(*,*) '         WHEN USING CARTESIAN GRID CUT-CELL FEATURE.'
336                       WRITE(*,*) '         DEFAULT WALLS ARE NOT ALLOWED WITH CUT-CELLS.'
337                       WRITE(*,*) '         THE DEFAULT WALL WAS REMOVED ALONG THIS CELL.'
338                       WRITE(*,*) ''
339                    ENDIF
340     !               CALL MFIX_EXIT(MYPE)
341     
342                 ENDIF
343     
344              ENDIF
345     
346     
347     
348           ENDDO
349     
350     
351     
352     
353           RETURN
354           END SUBROUTINE CG_SET_BC0
355     
356     !// Comments on the modifications for DMP version implementation
357     !// 001 Include header file and common declarations for parallelization
358     !// 020 New local variables for parallelization: FLAG_G , FLUID_AT_G
359     !// 360 Check if i,j,k resides on current processor
360