File: RELATIVE:/../../../mfix.git/model/des/make_arrays_des.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !  Module name: MAKE_ARRAYS_DES                                        !
3     !  Author: Jay Boyalakuntla                           Date: 12-Jun-04  !
4     !                                                                      !
5     !  Purpose: DES - allocating DES arrays
6     !                                                                      !
7     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
8           SUBROUTINE MAKE_ARRAYS_DES
9     
10           USE calc_collision_wall
11           USE compar
12           USE cutcell
13           USE GENERATE_PARTICLES, only: GENERATE_PARTICLE_CONFIG
14           USE des_rxns
15           USE des_stl_functions
16           USE des_thermo
17           USE desgrid
18           USE discretelement
19           USE error_manager
20           USE functions
21           USE funits
22           USE geometry
23           USE mpi_utility
24           USE param1
25           USE run
26           USE stl
27           use mfix_pic, only:  MPPIC
28           use mpi_funs_des, only: DES_PAR_EXCHANGE
29           use constant, only: PI
30           use stl_preproc_des, only: add_facet
31     
32           IMPLICIT NONE
33     !-----------------------------------------------
34     ! Local variables
35     !-----------------------------------------------
36           INTEGER :: I, J, K, L, IJK
37           INTEGER :: I1, I2, J1, J2, K1, K2, II, JJ, KK, IJK2
38           INTEGER :: lcurpar, lpip_all(0:numpes-1), lglobal_id
39           INTEGER :: CELL_ID, I_CELL, J_CELL, K_CELL, COUNT, NF
40           INTEGER :: IMINUS1, IPLUS1, JMINUS1, JPLUS1, KMINUS1, KPLUS1
41     
42     ! MPPIC related quantities
43           CALL INIT_ERR_MSG("MAKE_ARRAYS_DES")
44     
45     ! Check interpolation input.
46           CALL SET_FILTER_DES
47     
48     ! cfassign and des_init_bc called before reading the particle info
49           CALL CFASSIGN
50     
51           VOL_SURR(:) = ZERO
52     
53           ! initialize VOL_SURR array
54           DO K = KSTART2, KEND1
55              DO J = JSTART2, JEND1
56                 DO I = ISTART2, IEND1
57                    IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
58                    IJK = funijk(I,J,K)
59                    I1 = I
60                    I2 = I+1
61                    J1 = J
62                    J2 = J+1
63                    K1 = K
64                    K2 = merge(K, K+1, NO_K)
65     
66     ! looping over stencil points (node values)
67                    DO KK = K1, K2
68                       DO JJ = J1, J2
69                          DO II = I1, I2
70                             IF (DEAD_CELL_AT(II,JJ,KK)) CYCLE  ! skip dead cells
71                             IJK2 = funijk_map_c(II, JJ, KK)
72                             IF(FLUID_AT(IJK2)) VOL_SURR(IJK) = &
73                                VOL_SURR(IJK)+VOL(IJK2)
74                          ENDDO
75                       ENDDO
76                    ENDDO
77                 ENDDO
78              ENDDO
79           ENDDO
80     
81           VERTEX(1,:,WEST_FACEID) = (/zero, zero, zero/)
82           VERTEX(2,:,WEST_FACEID) = (/zero, 2*YLENGTH, zero/)
83           VERTEX(3,:,WEST_FACEID) = (/zero, zero, 2*ZLENGTH/)
84     
85           VERTEX(1,:,EAST_FACEID) = (/XLENGTH, zero, zero/)
86           VERTEX(2,:,EAST_FACEID) = (/XLENGTH, 2*YLENGTH, zero/)
87           VERTEX(3,:,EAST_FACEID) = (/XLENGTH, zero, 2*ZLENGTH/)
88     
89           VERTEX(1,:,SOUTH_FACEID) = (/zero, zero, zero/)
90           VERTEX(2,:,SOUTH_FACEID) = (/2*XLENGTH, zero, zero/)
91           VERTEX(3,:,SOUTH_FACEID) = (/zero, zero, 2*ZLENGTH/)
92     
93           VERTEX(1,:,NORTH_FACEID) = (/zero, YLENGTH, zero/)
94           VERTEX(2,:,NORTH_FACEID) = (/2*XLENGTH, YLENGTH, zero/)
95           VERTEX(3,:,NORTH_FACEID) = (/zero, YLENGTH, 2*ZLENGTH/)
96     
97           VERTEX(1,:,BOTTOM_FACEID) = (/zero, zero, zero/)
98           VERTEX(2,:,BOTTOM_FACEID) = (/2*XLENGTH, zero, zero/)
99           VERTEX(3,:,BOTTOM_FACEID) = (/zero, 2*YLENGTH, zero/)
100     
101           VERTEX(1,:,TOP_FACEID) = (/zero, zero, ZLENGTH/)
102           VERTEX(2,:,TOP_FACEID) = (/2*XLENGTH, zero, ZLENGTH/)
103           VERTEX(3,:,TOP_FACEID) = (/zero, 2*YLENGTH, ZLENGTH/)
104     
105     
106           NORM_FACE(:,WEST_FACEID) = (/one, zero, zero/)
107           NORM_FACE(:,EAST_FACEID) = (/-one, zero, zero/)
108           NORM_FACE(:,SOUTH_FACEID) = (/zero, one, zero/)
109           NORM_FACE(:,NORTH_FACEID) = (/zero, -one, zero/)
110           NORM_FACE(:,BOTTOM_FACEID) = (/zero, zero, one/)
111           NORM_FACE(:,TOP_FACEID) = (/zero, zero, -one/)
112     
113     
114           STL_FACET_TYPE(WEST_FACEID) = FACET_TYPE_NORMAL
115           STL_FACET_TYPE(EAST_FACEID) = FACET_TYPE_NORMAL
116           STL_FACET_TYPE(NORTH_FACEID) = FACET_TYPE_NORMAL
117           STL_FACET_TYPE(SOUTH_FACEID) = FACET_TYPE_NORMAL
118           STL_FACET_TYPE(TOP_FACEID) = FACET_TYPE_NORMAL
119           STL_FACET_TYPE(BOTTOM_FACEID) = FACET_TYPE_NORMAL
120     
121     ! initialize CELLNEIGHBOR_FACET array
122           DO K_CELL = DG_KSTART2, DG_KEND2
123           DO J_CELL = DG_JSTART2, DG_JEND2
124           DO I_CELL = DG_ISTART2, DG_IEND2
125     
126              IPLUS1  =  MIN (I_CELL + 1, DG_IEND2)
127              IMINUS1 =  MAX (I_CELL - 1, DG_ISTART2)
128     
129              JPLUS1  =  MIN (J_CELL + 1, DG_JEND2)
130              JMINUS1 =  MAX (J_CELL - 1, DG_JSTART2)
131     
132              KPLUS1  =  MIN (K_CELL + 1, DG_KEND2)
133              KMINUS1 =  MAX (K_CELL - 1, DG_KSTART2)
134     
135              CELL_ID = DG_FUNIJK(I_CELL,J_CELL,K_CELL)
136     
137              IF(.NOT.DES_PERIODIC_WALLS_X)THEN
138                 IF(I_CELL == DG_IMIN1) CALL ADD_FACET(CELL_ID,WEST_FACEID)
139                 IF(I_CELL == DG_IMAX1) CALL ADD_FACET(CELL_ID,EAST_FACEID)
140              ENDIF
141     
142              IF(.NOT.DES_PERIODIC_WALLS_Y)THEN
143                 IF(J_CELL == DG_JMIN1) CALL ADD_FACET(CELL_ID,SOUTH_FACEID)
144                 IF(J_CELL == DG_JMAX1) CALL ADD_FACET(CELL_ID,NORTH_FACEID)
145              ENDIF
146     
147              IF (DO_K .AND. .NOT.DES_PERIODIC_WALLS_Z)THEN
148                 IF(K_CELL == DG_KMIN1) CALL ADD_FACET(CELL_ID,BOTTOM_FACEID)
149                 IF(K_CELL == DG_KMAX1) CALL ADD_FACET(CELL_ID,TOP_FACEID)
150              ENDIF
151     
152           ENDDO
153           ENDDO
154           ENDDO
155     
156     ! Set the initial particle data.
157           IF(RUN_TYPE == 'NEW') THEN
158              IF(PARTICLES /= 0) THEN
159                 IF(GENER_PART_CONFIG) THEN
160                    CALL GENERATE_PARTICLE_CONFIG
161                 ELSE
162                    CALL READ_PAR_INPUT
163                 ENDIF
164              ENDIF
165     
166     ! Set the global ID for the particles and set the ghost cnt
167              ighost_cnt = 0
168              lpip_all = 0
169              lpip_all(mype) = pip
170              call global_all_sum(lpip_all)
171              lglobal_id = sum(lpip_all(0:mype-1))
172              imax_global_id = 0
173              do lcurpar  = 1,pip
174                 lglobal_id = lglobal_id + 1
175                 iglobal_id(lcurpar) = lglobal_id
176                 imax_global_id = iglobal_id(pip)
177              end do
178              call global_all_max(imax_global_id)
179     
180     ! Initialize old values
181              omega_new(:,:)   = zero
182     
183     ! Particle orientation
184              IF(PARTICLE_ORIENTATION) THEN
185                 ORIENTATION(1,:) = INIT_ORIENTATION(1)
186                 ORIENTATION(2,:) = INIT_ORIENTATION(2)
187                 ORIENTATION(3,:) = INIT_ORIENTATION(3)
188              ENDIF
189     
190     
191              IF (DO_OLD) THEN
192                 omega_old(:,:)   = zero
193                 des_pos_old(:,:) = des_pos_new(:,:)
194                 des_vel_old(:,:) = des_vel_new(:,:)
195              ENDIF
196     
197     ! Read the restart file.
198           ELSEIF(RUN_TYPE == 'RESTART_1' .OR. RUN_TYPE == 'RESTART_2') THEN
199     
200              CALL READ_RES0_DES
201              imax_global_id = maxval(iglobal_id(1:pip))
202              call global_all_max(imax_global_id)
203     
204     ! Initizlie the old values.
205              IF (DO_OLD) THEN
206                 omega_old(:,:)   = omega_new(:,:)
207                 des_pos_old(:,:) = des_pos_new(:,:)
208                 des_vel_old(:,:) = des_vel_new(:,:)
209              ENDIF
210              IF(ENERGY_EQ) DES_T_s_OLD(:) = DES_T_s_NEW(:)
211     
212           ELSE
213     
214              WRITE(ERR_MSG, 1100)
215              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
216      1100 FORMAT('Error 1100: Unsupported RUN_TYPE for DES.')
217     
218           ENDIF
219     
220           IF(RUN_TYPE == 'RESTART_2') VTP_FINDEX=0
221     
222     ! setting additional particle properties now that the particles
223     ! have been identified
224           DO L = 1, MAX_PIP
225     ! Skip 'empty' locations when populating the particle property arrays.
226              IF(IS_NONEXISTENT(L)) CYCLE
227              IF(IS_GHOST(L) .OR. IS_ENTERING_GHOST(L) .OR. IS_EXITING_GHOST(L)) CYCLE
228              PVOL(L) = (4.0D0/3.0D0)*PI*DES_RADIUS(L)**3
229              PMASS(L) = PVOL(L)*RO_SOL(L)
230              OMOI(L) = 2.5D0/(PMASS(L)*DES_RADIUS(L)**2) !ONE OVER MOI
231           ENDDO
232     
233           CALL SET_PHASE_INDEX
234           CALL INIT_PARTICLES_IN_CELL
235     
236     ! do_nsearch should be set before calling particle in cell
237           DO_NSEARCH =.TRUE.
238     ! Bin the particles to the DES grid.
239           CALL DESGRID_PIC(PLOCATE=.TRUE.)
240           CALL DES_PAR_EXCHANGE
241           CALL PARTICLES_IN_CELL
242     
243           IF(DEM_SOLIDS) THEN
244              CALL NEIGHBOUR
245              CALL INIT_SETTLING_DEM
246           ENDIF
247     
248     ! Calculate interpolation weights
249           CALL CALC_INTERP_WEIGHTS
250     ! Calculate mean fields using either interpolation or cell averaging.
251           CALL COMP_MEAN_FIELDS
252     
253           IF(MPPIC) CALL CALC_DTPIC
254     
255           CALL FINL_ERR_MSG
256     
257           RETURN
258           END SUBROUTINE MAKE_ARRAYS_DES
259