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

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Subroutine: SET_FLAGS                                               C
4     !  Purpose: This module assigns a flag to a cell to identify its type. C
5     !                                                                      C
6     !  Author: M. Syamlal                                 Date: 29-JAN-92  C
7     !  Reviewer: P. Nicoletti, W. Rogers, S. Venkatesan   Date: 31-JAN-92  C
8     !                                                                      C
9     !  Revision Number: 1                                                  C
10     !  Purpose: In cylindrical geometry, set the b.c at X=0 to free-slip   C
11     !  Author: M. Syamlal                                 Date: 09-APR-92  C
12     !                                                                      C
13     !  Revision Number: 2                                                  C
14     !  Purpose: Initialize FLAG_E, FLAG_N, FLAG_T for IS specifications.   C
15     !           Change definition of Flags.                                C
16     !  Author: M. Syamlal                                 Date: 21-OCT-92  C
17     !  Reviewer: M. Syamlal                               Date: 11-DEC-92  C
18     !                                                                      C
19     !  Revision Number: 3                                                  C
20     !  Purpose: Define FLAG using info from ICBC_FLAG                      C
21     !  Author: M. Syamlal                                 Date: 21-APR-93  C
22     !                                                                      C
23     !  Literature/Document References:                                     C
24     !                                                                      C
25     !  Variables referenced: IMAX2, JMAX2, KMAX2, BC_DEFINED, BC_TYPE,     C
26     !                        BC_K_b, BC_K_t, BC_J_s, BC_J_n, BC_I_w,       C
27     !                        IS_K_b, IS_K_t, IS_J_s, IS_J_n, IS_I_w,       C
28     !                        BC_I_e, NO_I, NO_J, NO_K                      C
29     !  Variables modified: FLAG, FLAG_E, FLAG_N, FLAG_T                    C
30     !                                                                      C
31     !  Local variables: I, J, K, IJK, L, FLAGX                             C
32     !                                                                      C
33     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
34     
35           SUBROUTINE SET_FLAGS
36     
37     !-----------------------------------------------
38     ! Modules
39     !-----------------------------------------------
40           USE bc
41           USE compar
42           USE exit, only: mfix_exit
43           USE fldvar
44           USE function3
45           USE functions
46           USE funits
47           USE geometry
48           USE indices
49           USE is
50           USE parallel
51           USE param
52           USE param1
53           USE physprop
54           USE sendrecv
55           USE sendrecv3
56           use mpi_utility
57           IMPLICIT NONE
58     !-----------------------------------------------
59     ! Local variables
60     !-----------------------------------------------
61     ! Indices
62           INTEGER :: I, J, K, IJK, IJK1
63     ! Local DO loop index for b.c. specification
64           INTEGER :: L
65     ! Temporary storage for FLAG value
66           INTEGER :: FLAGX
67           integer, allocatable :: arr1(:)
68     !-----------------------------------------------
69     
70     !  Cell flag definitions
71     !  FLAG  ICBC_FLAG BC_TYPE        Cell type
72     !  ----- --------- -------        ---------
73     !   1       .        -            Cell containing gas or solids or both
74     !  10       p      P_INFLOW       Specified pressure inflow cell
75     !  11       P      P_OUTFLOW      Specified pressure outflow cell
76     !  20       I      MASS_INFLOW    Specified mass flux inflow cell
77     !  21       O      MASS_OUTFLOW   Specified mass flux outflow cell
78     !  31       o      OUTFLOW        outflow cell
79     ! 100       W      NO_SLIP_WALL   Internal/external wall with no-slip b.c.
80     ! 101       S      FREE_SLIP_WALL Internal/external wall with free-slip
81     ! 102       s      PAR_SLIP_WALL  Internal/external wall with partial-slip b.c.
82     ! 106       c      CYCLIC         Cyclic b.c.
83     ! 107       C      CYCLIC_PD      Cyclic b.c. with pressure drop
84     ! Flag values greater than 100 are considered to be wall cells
85     ! (see function.inc).
86     
87     ! make the wall cells adjacent to flow boundaries free-slip wall to
88     ! avoid unphysical strain rates in fluid cells adjacent to the flow
89     ! boundary
90     ! ---------------------------------------------------------------->>>
91     !!$omp  parallel do private( IJK) &
92     !!$omp  schedule(static)
93           DO i = istart4, iend4
94              DO j = jstart4, jend4
95                 DO k = kstart4, kend4
96     
97                   IJK = funijk(i, j, k)
98                   SELECT CASE (TRIM(ICBC_FLAG(IJK)(1:1)))
99                     CASE ('p', 'P', 'I', 'O', 'o')
100     
101                     ijk1 = bound_funijk(i+1, j, k)
102                     IF(TRIM(ICBC_FLAG(IJK1)(1:1)) == 'W')ICBC_FLAG(IJK1)(1:1)='S'
103     
104                     ijk1 = bound_funijk(i-1, j, k)
105                     IF(TRIM(ICBC_FLAG(IJK1)(1:1)) == 'W')ICBC_FLAG(IJK1)(1:1)='S'
106     
107                     ijk1 = bound_funijk(i, j+1, k)
108                     IF(TRIM(ICBC_FLAG(IJK1)(1:1)) == 'W')ICBC_FLAG(IJK1)(1:1)='S'
109     
110                     ijk1 = bound_funijk(i, j-1, k)
111                     IF(TRIM(ICBC_FLAG(IJK1)(1:1)) == 'W')ICBC_FLAG(IJK1)(1:1)='S'
112     
113                     ijk1 = bound_funijk(i, j, k+1)
114                     IF(TRIM(ICBC_FLAG(IJK1)(1:1)) == 'W')ICBC_FLAG(IJK1)(1:1)='S'
115     
116                     ijk1 = bound_funijk(i, j, k-1)
117                     IF(TRIM(ICBC_FLAG(IJK1)(1:1)) == 'W')ICBC_FLAG(IJK1)(1:1)='S'
118                   END SELECT
119                 ENDDO
120               ENDDO
121           ENDDO
122     ! ----------------------------------------------------------------<<<
123     
124     ! Define the numerical value of the variable flag for all cells based
125     ! on the corresponding character value of icbc_flag.  By this point the
126     ! icbc_flag has been defined in all cells (see check_data_06 and
127     ! check_data07 -> get_flow_bc and get_wall_bc)
128     ! ---------------------------------------------------------------->>>
129     !!$omp  parallel do private( IJK) &
130     !!$omp&  schedule(static)
131           DO IJK = ijkstart3, ijkend3
132              SELECT CASE (TRIM(ICBC_FLAG(IJK)(1:1)))
133              CASE ('.')
134                 FLAG(IJK) = 1
135              CASE ('p')
136                 FLAG(IJK) = 10
137              CASE ('P')
138                 FLAG(IJK) = 11
139              CASE ('I')
140                 FLAG(IJK) = 20
141              CASE ('O')
142                 FLAG(IJK) = 21
143              CASE ('o')
144                 FLAG(IJK) = 31
145              CASE ('W')
146                 FLAG(IJK) = 100
147              CASE ('S')
148                 FLAG(IJK) = 101
149              CASE ('s')
150                 FLAG(IJK) = 102
151              CASE ('c')
152                 FLAG(IJK) = 106
153              CASE ('C')
154                 FLAG(IJK) = 107
155              CASE DEFAULT
156     
157     ! Access to only one thread at a time
158     !!$omp       critical
159                 IF(DMP_LOG)WRITE (UNIT_LOG, 1000) IJK, ICBC_FLAG(IJK)
160                 call mfix_exit(myPE)
161     !!$omp       end critical
162              END SELECT
163     ! ----------------------------------------------------------------<<<
164     
165     ! Initialize cell face flags.  UNDEFINED_I should be a large +ve value.
166              FLAG_E(IJK) = UNDEFINED_I
167              FLAG_N(IJK) = UNDEFINED_I
168              FLAG_T(IJK) = UNDEFINED_I
169           ENDDO
170     
171     
172     
173     ! Setting up flags for the higher-order implementation.
174     ! ---------------------------------------------------------------->>>
175           call send_recv(flag)
176     
177           DO i = istart3, iend3
178              DO j = jstart3, jend3
179                 DO k = kstart3, kend3
180                    Flag3(funijk3(i,j,k)) = Flag(funijk(i,j,k))
181                 ENDDO
182              ENDDO
183           ENDDO
184     
185           DO i = istart4, iend4
186              DO j = jstart4, jend4
187                 DO k = kstart4, kend4
188                    If(i.eq.istart4.and.istart4.ne.istart3) then
189                       Flag3(funijk3(i,j,k)) = Flag3(funijk3(i+1,j,k))
190                    endif
191     
192                    If(j.eq.jstart4.and.kstart4.ne.kstart3) then
193                       Flag3(funijk3(i,j,k)) = Flag3(funijk3(i,j+1,k))
194                    endif
195     
196                    If(k.eq.kstart4.and.kstart4.ne.kstart3) then
197                       Flag3(funijk3(i,j,k)) = Flag3(funijk3(i,j,k+1))
198                    endif
199     
200                    If(i.eq.iend4.and.iend4.ne.iend3) then
201                       Flag3(funijk3(i,j,k)) = Flag3(funijk3(i-1,j,k))
202                    endif
203     
204                    If(j.eq.jend4.and.jend4.ne.jend3) then
205                       Flag3(funijk3(i,j,k)) = Flag3(funijk3(i,j-1,k))
206                    endif
207     
208                    If(k.eq.kend4.and.kend4.ne.kend3) then
209                       Flag3(funijk3(i,j,k)) = Flag3(funijk3(i,j,k-1))
210                    endif
211     
212                 ENDDO
213              ENDDO
214           ENDDO
215     
216           call send_recv3(flag3)
217     ! ----------------------------------------------------------------<<<
218     
219     
220     ! Set flag_e, flag_n and flag_b to indicate any internal surfaces. If
221     ! the flag is greater than or equal to 2000, then there is no internal
222     ! surface.
223     ! ---------------------------------------------------------------->>>
224     
225           DO L = 1, DIMENSION_IS
226     ! Make sure an IS has been specified
227              IF (IS_DEFINED(L)) THEN
228                 IF (IS_TYPE(L)=='IMPERMEABLE' .OR. &
229                     IS_TYPE(L)(3:13)=='IMPERMEABLE') THEN
230                    FLAGX = 0
231                 ELSEIF (IS_TYPE(L)=='SEMIPERMEABLE' .OR. &
232                         IS_TYPE(L)(3:15)=='SEMIPERMEABLE') THEN
233                    FLAGX = 1000 + L
234                 ELSE
235                    IF(DMP_LOG)WRITE (UNIT_LOG, 1100) L
236                    call mfix_exit(myPE)
237                 ENDIF
238     
239                 IF (IS_X_W(L)==IS_X_E(L) .AND. DO_I) THEN
240                    IS_PLANE(L) = 'E'
241                    I = IS_I_W(L)
242                    DO K = IS_K_B(L), IS_K_T(L)
243                       DO J = IS_J_S(L), IS_J_N(L)
244                          IF (.NOT.IS_ON_myPE_plus2layers(I,J,K)) CYCLE
245                             IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
246                          IJK = FUNIJK(I,J,K)
247                          FLAG_E(IJK) = FLAGX
248                       ENDDO
249                    ENDDO
250                 ELSEIF (IS_TYPE(L)(1:1) == 'X') THEN
251                    IS_PLANE(L) = 'E'
252                    DO I = IS_I_W(L), IS_I_E(L)
253                       DO K = IS_K_B(L), IS_K_T(L)
254                          DO J = IS_J_S(L), IS_J_N(L)
255                             IF (.NOT.IS_ON_myPE_plus2layers(I,J,K)) CYCLE
256                             IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
257                             IJK = FUNIJK(I,J,K)
258                             FLAG_E(IJK) = FLAGX
259                          ENDDO
260                       ENDDO
261                    ENDDO
262                 ENDIF
263     
264                 IF (IS_Y_S(L)==IS_Y_N(L) .AND. DO_J) THEN
265                    IS_PLANE(L) = 'N'
266                    J = IS_J_S(L)
267                    DO K = IS_K_B(L), IS_K_T(L)
268                       DO I = IS_I_W(L), IS_I_E(L)
269                          IF (.NOT.IS_ON_myPE_plus2layers(I,J,K)) CYCLE
270                             IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
271                          IJK = FUNIJK(I,J,K)
272                          FLAG_N(IJK) = FLAGX
273                       ENDDO
274                    ENDDO
275                 ELSEIF (IS_TYPE(L)(1:1) == 'Y') THEN
276                    IS_PLANE(L) = 'N'
277                    DO J = IS_J_S(L), IS_J_N(L)
278                       DO K = IS_K_B(L), IS_K_T(L)
279                          DO I = IS_I_W(L), IS_I_E(L)
280                             IF (.NOT.IS_ON_myPE_plus2layers(I,J,K)) CYCLE
281                             IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
282                             IJK = FUNIJK(I,J,K)
283                             FLAG_N(IJK) = FLAGX
284                          ENDDO
285                       ENDDO
286                    ENDDO
287                 ENDIF
288     
289                 IF (IS_Z_B(L)==IS_Z_T(L) .AND. DO_K) THEN
290                    IS_PLANE(L) = 'T'
291                    K = IS_K_B(L)
292                    DO J = IS_J_S(L), IS_J_N(L)
293                       DO I = IS_I_W(L), IS_I_E(L)
294                          IF (.NOT.IS_ON_myPE_plus2layers(I,J,K)) CYCLE
295                             IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
296                          IJK = FUNIJK(I,J,K)
297                          FLAG_T(IJK) = FLAGX
298                       ENDDO
299                    ENDDO
300                 ELSEIF (IS_TYPE(L)(1:1) == 'Z') THEN
301                    IS_PLANE(L) = 'T'
302                    DO K = IS_K_B(L), IS_K_T(L)
303                       DO J = IS_J_S(L), IS_J_N(L)
304                          DO I = IS_I_W(L), IS_I_E(L)
305                             IF (.NOT.IS_ON_myPE_plus2layers(I,J,K)) CYCLE
306                             IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
307                             IJK = FUNIJK(I,J,K)
308                             FLAG_T(IJK) = FLAGX
309                          ENDDO
310                       ENDDO
311                    ENDDO
312                 ENDIF
313     
314              ENDIF
315           call send_recv(flag,2)
316           call send_recv(flag_t,2)
317           call send_recv(flag_n,2)
318           call send_recv(flag_e,2)
319           ENDDO    ! end do loop (l = 1, dimension_is)
320     ! ----------------------------------------------------------------<<<
321     
322           IF (MYPE.EQ.PE_IO) THEN
323              ALLOCATE (ARR1(IJKMAX3))
324           ELSE
325              ALLOCATE (ARR1(1))
326           ENDIF
327     
328           CALL GATHER(FLAG,ARR1,ROOT)
329           CALL SCATTER(FLAG,ARR1,ROOT)
330     
331           DEALLOCATE (ARR1)
332     
333     
334           RETURN
335      1000 FORMAT(/1X,70('*')//' From: SET_FLAGS',/&
336              ' Message: ICBC_FLAG(',I3,') = ',&
337              A3,' is illegal',/1X,70('*')/)
338      1100 FORMAT(/1X,70('*')//' From: SET_FLAGS',/&
339              ' Message: Unknown IS_TYPE(',I3,&
340              ')',/1X,70('*')/)
341     
342           END SUBROUTINE SET_FLAGS
343     
344     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
345     !                                                                      C
346     !  Module name: SET_FLAGS1                                             C
347     !  Purpose: Assign IP flag to the faces of wall cells                  C
348     !                                                                      C
349     !  Notes: This routine may still leave flag_e, flag_n and flag_t       C
350     !         undefined in some boundary cells                             C
351     !                                                                      C
352     !  Author: M. Syamlal                                 Date: 15-MAY-96  C
353     !  Reviewer:                                          Date:            C
354     !                                                                      C
355     !  Literature/Document References:                                     C
356     !                                                                      C
357     !  Variables referenced:                                               C
358     !  Variables modified: FLAG_E, FLAG_N, FLAG_T                          C
359     !  Local variables:                                                    C
360     !                                                                      C
361     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
362     
363           SUBROUTINE SET_FLAGS1
364     
365     !-----------------------------------------------
366     ! Modules
367     !-----------------------------------------------
368           USE param
369           USE param1
370           USE parallel
371           USE fldvar
372           USE geometry
373           USE bc
374           USE is
375           USE indices
376           USE physprop
377           USE funits
378           USE compar
379           USE sendrecv
380           USE mpi_utility
381           USE functions
382           IMPLICIT NONE
383     !-----------------------------------------------
384     ! Local variables
385     !-----------------------------------------------
386     ! Indices
387           INTEGER :: IJK, IMJK, IJMK, IJKM, IPJK, IJPK, IJKP
388           INTEGER :: I, J, K
389     !
390           INTEGER, DIMENSION(:), allocatable :: FLAG_TEMP
391           INTEGER :: flag_size
392     !-----------------------------------------------
393     
394     
395     ! Allocate storage for temporary flag arrays
396           flag_size = ijkmax3
397           if (myPE.eq.root) then
398               flag_size = ijkmax3
399           endif
400           allocate( flag_temp(flag_size) )
401     
402     
403           DO IJK = ijkstart3,ijkend3
404              IMJK = IM_OF(IJK)
405              IJMK = JM_OF(IJK)
406              IJKM = KM_OF(IJK)
407              IPJK = IP_OF(IJK)
408              IJPK = JP_OF(IJK)
409              IJKP = KP_OF(IJK)
410              I = I_OF(IJK)
411              J = J_OF(IJK)
412              K = K_OF(IJK)
413              IF(.NOT.IS_ON_myPE_plus2layers(I,J,K)) CYCLE
414              IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
415     
416     ! If the flag is greater than or equal to 2000, there is no
417     ! internal surface.
418              IF (WALL_AT(IJK)) THEN
419     ! ---------------------------------------------------------------->>>
420     ! the default is equivalent to an impermeable surface and these cells
421     ! will be treated as such in the momentum routines
422                 FLAG_E(IJK) = 0
423                 FLAG_N(IJK) = 0
424                 FLAG_T(IJK) = 0
425                 FLAG_E(IMJK) = 0
426                 FLAG_N(IJMK) = 0
427                 FLAG_T(IJKM) = 0
428     
429                 IF (CYCLIC_AT(IJK)) THEN
430     ! make the upper (E, N, T) boundary permeable
431                    IF (I == IMAX2) THEN
432                       IF ((J/=1.AND.J/=0.) .AND. (J/=JMAX2.AND.J/=JMAX3)) THEN
433                          IF (NO_K) THEN
434                             IF(.NOT.WALL_AT(IMJK)) FLAG_E(IMJK) = 2000
435                          ELSEIF ((K/=1.AND.K/=0) .AND. (K/=KMAX2.AND.K/=KMAX3)) THEN
436                             IF(.NOT.WALL_AT(IMJK)) FLAG_E(IMJK) = 2000
437                          ENDIF
438                       ENDIF
439                    ENDIF
440                    IF (J == JMAX2) THEN
441                       IF ((I/=1.AND.I/=0) .AND. (I/=IMAX2.AND.I/=IMAX3)) THEN
442                          IF (NO_K) THEN
443                             IF(.NOT.WALL_AT(IJMK)) FLAG_N(IJMK) = 2000
444                          ELSE IF ((K/=1.AND.K/=0) .AND. (K/=KMAX2.AND.K/=KMAX3)) THEN
445                             IF(.NOT.WALL_AT(IJMK)) FLAG_N(IJMK) = 2000
446                          ENDIF
447                       ENDIF
448                     ENDIF
449                    IF (K == KMAX2) THEN
450                       IF ((J/=1.AND.J/=0.) .AND. (J/=JMAX2.AND.J/=JMAX3)) THEN
451                          IF ((I/=1.AND.I/=0) .AND. (I/=IMAX2.AND.I/=IMAX3) .AND. &
452                            .NOT.WALL_AT(IJKM)) FLAG_T(IJKM) = 2000
453                       ENDIF
454                    ENDIF
455     
456                 ENDIF   ! end if cyclic_at(ijk)
457     
458     ! ----------------------------------------------------------------<<<
459              ELSEIF (FLUID_AT(IJK)) THEN
460     ! ---------------------------------------------------------------->>>
461     
462                 IF ( .NOT.WALL_AT(IMJK) .AND. FLAG_E(IMJK)==UNDEFINED_I) &
463                    FLAG_E(IMJK) = 2000 + FLAG(IMJK)
464                 IF ( .NOT.WALL_AT(IJMK) .AND. FLAG_N(IJMK)==UNDEFINED_I) &
465                    FLAG_N(IJMK) = 2000 + FLAG(IJMK)
466                 IF ( .NOT.WALL_AT(IJKM) .AND. FLAG_T(IJKM)==UNDEFINED_I) &
467                    FLAG_T(IJKM) = 2000 + FLAG(IJKM)
468                 IF ( .NOT.WALL_AT(IPJK) .AND. FLAG_E(IJK)==UNDEFINED_I) &
469                    FLAG_E(IJK) = 2000 + FLAG(IPJK)
470                 IF ( .NOT.WALL_AT(IJPK) .AND. FLAG_N(IJK)==UNDEFINED_I) &
471                    FLAG_N(IJK) = 2000 + FLAG(IJPK)
472                 IF ( .NOT.WALL_AT(IJKP) .AND. FLAG_T(IJK)==UNDEFINED_I) &
473                    FLAG_T(IJK) = 2000 + FLAG(IJKP)
474     
475              ENDIF   ! end if/else (wall_at(ijk)/fluid_at(ijk))
476     
477           ENDDO    ! end do loop (ijk = ijkstart3,ijkend3)
478     ! ----------------------------------------------------------------<<<
479     
480     ! Fill the ghost layers using gather and scatter
481           call gather( flag_e, flag_temp )
482           call scatter( flag_e, flag_temp )
483           call gather( flag_n, flag_temp )
484           call scatter( flag_n, flag_temp )
485           call gather( flag_t, flag_temp )
486           call scatter( flag_t, flag_temp )
487     
488     ! deallocate storage of temporary flag arrays
489           deallocate( flag_temp )
490           call send_recv(flag_t,2)
491           call send_recv(flag_n,2)
492           call send_recv(flag_e,2)
493     
494           RETURN
495           END SUBROUTINE SET_FLAGS1
496