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

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
2     !                                                                      C
3     !  Module name: DES_ALLOCATE                                           C
4     !                                                                      C
5     !  Purpose: subroutines to allocate all DEM arrays                     C
6     !                                                                      C
7     !  Author: Rahul Garg                               Date: 1-Dec-2013   C
8     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
9     
10     MODULE DES_ALLOCATE
11     
12       PUBLIC:: DES_ALLOCATE_ARRAYS, ADD_PAIR, PARTICLE_GROW, ALLOCATE_DEM_MI
13     
14     CONTAINS
15     
16     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
17     !                                                                      C
18     !  Subroutine: DES_ALLOCATE_ARRAYS                                     C
19     !  Purpose: Original allocate arrays subroutines for DES               C
20     !                                                                      C
21     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
22           SUBROUTINE DES_ALLOCATE_ARRAYS
23     
24     !-----------------------------------------------
25     ! Modules
26     !-----------------------------------------------
27           USE param
28           USE param1
29           USE constant
30           USE discretelement
31           Use indices
32           Use geometry
33           Use compar
34           Use physprop
35           Use des_bc
36           Use pic_bc
37           use funits
38           USE mfix_pic
39           Use des_thermo
40           Use des_rxns
41           USE cutcell
42           USE functions
43     
44           use run, only: ENERGY_EQ
45           use run, only: ANY_SPECIES_EQ
46     
47           use particle_filter, only: DES_INTERP_SCHEME_ENUM
48           use particle_filter, only: DES_INTERP_GARG
49           use particle_filter, only: DES_INTERP_DPVM
50           use particle_filter, only: DES_INTERP_GAUSS
51           use particle_filter, only: FILTER_SIZE
52           use particle_filter, only: FILTER_CELL
53           use particle_filter, only: FILTER_WEIGHT
54     
55     ! Use the error manager for posting error messages.
56     !---------------------------------------------------------------------//
57           use error_manager
58     
59           IMPLICIT NONE
60     !-----------------------------------------------
61     ! Local variables
62     !-----------------------------------------------
63     ! indices
64           INTEGER :: IJK
65     !-----------------------------------------------
66     
67           CALL INIT_ERR_MSG("DES_ALLOCATE_ARRAYS")
68     
69     ! For parallel processing the array size required should be either
70     ! specified by the user or could be determined from total particles
71     ! with some factor.
72           MAX_PIP = merge(0, PARTICLES/numPEs, PARTICLES==UNDEFINED_I)
73           MAX_PIP = MAX(MAX_PIP,4)
74     
75           WRITE(ERR_MSG,1000) trim(iVal(MAX_PIP))
76           CALL FLUSH_ERR_MSG(HEADER = .FALSE., FOOTER = .FALSE.)
77     
78      1000 FORMAT('Initial DES Particle array size: ',A)
79     
80     ! DES Allocatable arrays
81     !-----------------------------------------------
82     ! Dynamic particle info including another index for parallel
83     ! processing for ghost
84           ALLOCATE( PARTICLE_STATE (MAX_PIP) )
85           ALLOCATE (iglobal_id(max_pip))
86     
87     ! R.Garg: Allocate necessary arrays for PIC mass inlet/outlet BCs
88           IF(PIC_BCMI /= 0 .OR. PIC_BCMO /=0) CALL ALLOCATE_PIC_MIO
89     
90     ! Particle attributes
91     ! Radius, density, mass, moment of inertia
92           Allocate(  DES_RADIUS (MAX_PIP) )
93           Allocate(  RO_Sol (MAX_PIP) )
94           Allocate(  PVOL (MAX_PIP) )
95           Allocate(  PMASS (MAX_PIP) )
96           Allocate(  OMOI (MAX_PIP) )
97     
98     ! Old and new particle positions, velocities (translational and
99     ! rotational)
100           Allocate(  DES_POS_NEW (DIMN,MAX_PIP) )
101           Allocate(  DES_VEL_NEW (DIMN,MAX_PIP) )
102           Allocate(  OMEGA_NEW (DIMN,MAX_PIP) )
103     
104           IF(PARTICLE_ORIENTATION) Allocate(  ORIENTATION (DIMN,MAX_PIP) )
105     
106           IF (DO_OLD) THEN
107              Allocate(  DES_POS_OLD (DIMN,MAX_PIP) )
108              Allocate(  DES_VEL_OLD (DIMN,MAX_PIP) )
109              Allocate(  DES_ACC_OLD (DIMN,MAX_PIP) )
110              Allocate(  OMEGA_OLD (DIMN,MAX_PIP) )
111              Allocate(  ROT_ACC_OLD (DIMN,MAX_PIP))
112           ENDIF
113     
114     ! Allocating user defined array
115           IF(DES_USR_VAR_SIZE > 0) &
116              Allocate( DES_USR_VAR(DES_USR_VAR_SIZE,MAX_PIP) )
117     
118     ! Particle positions at the last call neighbor search algorithm call
119           Allocate(  PPOS (DIMN,MAX_PIP) )
120     
121     ! Total, normal and tangetial forces
122           Allocate(  FC (DIMN,MAX_PIP) )
123     
124     ! Torque
125           Allocate(  TOW (DIMN,MAX_PIP) )
126     
127     
128     ! allocate variable for des grid binning
129           allocate(dg_pijk(max_pip)); dg_pijk=0
130           allocate(dg_pijkprv(max_pip)); dg_pijkprv=0
131     
132     ! allocate variables related to ghost particles
133           allocate(ighost_updated(max_pip))
134     
135     
136     
137           Allocate(  wall_collision_facet_id (COLLISION_ARRAY_MAX, MAX_PIP) )
138           wall_collision_facet_id(:,:) = -1
139           Allocate(  wall_collision_PFT (DIMN, COLLISION_ARRAY_MAX, MAX_PIP) )
140     
141     ! Temporary variables to store wall position, velocity and normal vector
142           Allocate(  WALL_NORMAL  (NWALLS,DIMN) )
143     
144           NEIGH_MAX = MAX_PIP
145     
146           Allocate(  NEIGHBOR_INDEX (MAX_PIP) )
147           Allocate(  NEIGHBOR_INDEX_OLD (MAX_PIP) )
148           Allocate(  NEIGHBORS (NEIGH_MAX) )
149           NEIGHBORS(:) = 0
150           Allocate(  NEIGHBORS_OLD (NEIGH_MAX) )
151           Allocate(  PFT_NEIGHBOR (3,NEIGH_MAX) )
152           Allocate(  PFT_NEIGHBOR_OLD (3,NEIGH_MAX) )
153     
154     ! Variable that stores the particle in cell information (ID) on the
155     ! computational fluid grid defined by imax, jmax and kmax in mfix.dat
156           ALLOCATE(PIC(DIMENSION_3))
157           DO IJK=1,DIMENSION_3
158             NULLIFY(pic(ijk)%p)
159           ENDDO
160     
161     ! Particles in a computational fluid cell (for volume fraction)
162           Allocate(  PINC (DIMENSION_3) )
163     
164     ! For each particle track its i,j,k location on computational fluid grid
165     ! defined by imax, jmax and kmax in mfix.dat and phase no.
166           Allocate(  PIJK (MAX_PIP,5) )
167     
168           ALLOCATE(DRAG_AM(DIMENSION_3))
169           ALLOCATE(DRAG_BM(DIMENSION_3, DIMN))
170           ALLOCATE(F_gp(MAX_PIP ))
171           F_gp(1:MAX_PIP)  = ZERO
172     
173     ! Explicit drag force acting on a particle.
174           Allocate(DRAG_FC (DIMN,MAX_PIP) )
175     
176     ! force due to gas-pressure gradient
177           ALLOCATE(P_FORCE(DIMN, DIMENSION_3))
178     
179     ! Volume averaged solids volume in a computational fluid cell
180           Allocate(  DES_U_s (DIMENSION_3, DES_MMAX) )
181           Allocate(  DES_V_s (DIMENSION_3, DES_MMAX) )
182           Allocate(  DES_W_s (DIMENSION_3, DES_MMAX) )
183     
184     ! Volume of nodes
185           ALLOCATE(DES_VOL_NODE(DIMENSION_3))
186     
187           ALLOCATE(F_GDS(DIMENSION_3))
188           ALLOCATE(VXF_GDS(DIMENSION_3))
189     
190           SELECT CASE(DES_INTERP_SCHEME_ENUM)
191           CASE(DES_INTERP_DPVM, DES_INTERP_GAUSS)
192              ALLOCATE(FILTER_CELL(FILTER_SIZE, MAX_PIP))
193              ALLOCATE(FILTER_WEIGHT(FILTER_SIZE, MAX_PIP))
194           CASE(DES_INTERP_GARG)
195              ALLOCATE(DES_ROPS_NODE(DIMENSION_3, DES_MMAX))
196              ALLOCATE(DES_VEL_NODE(DIMENSION_3, DIMN, DES_MMAX))
197           END SELECT
198     
199     ! Variables for hybrid model
200           IF (DES_CONTINUUM_HYBRID) THEN
201              ALLOCATE(SDRAG_AM(DIMENSION_3,DIMENSION_M))
202              ALLOCATE(SDRAG_BM(DIMENSION_3, DIMN,DIMENSION_M))
203     
204              ALLOCATE(F_SDS(DIMENSION_3,DIMENSION_M))
205              ALLOCATE(VXF_SDS(DIMENSION_3,DIMENSION_M))
206           ENDIF
207     ! Bulk density in a computational fluid cell / for communication with
208     ! MFIX continuum
209           ALLOCATE( DES_ROP_S(DIMENSION_3, DES_MMAX) )
210           ALLOCATE( DES_ROP_SO(DIMENSION_3, DES_MMAX) )
211     
212     ! MP-PIC related
213           IF(MPPIC) THEN
214              Allocate(PS_FORCE_PIC(DIMENSION_3, DIMN))
215              ALLOCATE(DES_STAT_WT(MAX_PIP))
216              ALLOCATE(DES_VEL_MAX(DIMN))
217              ALLOCATE(PS_GRAD(MAX_PIP, DIMN))
218              ALLOCATE(AVGSOLVEL_P(DIMN, MAX_PIP))
219              ALLOCATE(EPG_P(MAX_PIP))
220     
221              Allocate(PIC_U_S(DIMENSION_3, DES_MMAX))
222              Allocate(PIC_V_S(DIMENSION_3, DES_MMAX))
223              Allocate(PIC_W_S(DIMENSION_3, DES_MMAX))
224     
225              Allocate(PIC_P_s (DIMENSION_3, DES_MMAX) )
226     !         ALLOCATE(MPPIC_VPTAU(MAX_PIP, DIMN))
227              PIC_U_s = zero
228              PIC_V_s = zero
229              PIC_W_s = zero
230              PIC_P_s = zero
231           ENDIF
232     
233     
234     ! Granular temperature in a computational fluid cell
235           Allocate(DES_THETA (DIMENSION_3, DES_MMAX) )
236     
237     ! Averaged velocity obtained by averaging over all the particles
238           ALLOCATE(DES_VEL_AVG(DIMN) )
239     
240     ! Global Granular Energy
241           ALLOCATE(GLOBAL_GRAN_ENERGY(DIMN) )
242           ALLOCATE(GLOBAL_GRAN_TEMP(DIMN) )
243     
244     ! variable for bed height of solids phase M
245           ALLOCATE(BED_HEIGHT(DES_MMAX))
246     
247     ! ---------------------------------------------------------------->>>
248     ! BEGIN COHESION
249           IF(USE_COHESION) THEN
250     ! Matrix location of particle  (should be allocated in case user wishes
251     ! to invoke routines in /cohesion subdirectory
252              Allocate(  PostCohesive (MAX_PIP) )
253           ENDIF
254     ! END COHESION
255     ! ----------------------------------------------------------------<<<
256     
257     ! ---------------------------------------------------------------->>>
258     ! BEGIN Thermodynamic Allocation
259           IF(ENERGY_EQ)THEN
260     ! Particle temperature
261              Allocate( DES_T_s_OLD( MAX_PIP ) )
262              Allocate( DES_T_s_NEW( MAX_PIP ) )
263     ! Specific heat
264              Allocate( DES_C_PS( MAX_PIP ) )
265     ! Species mass fractions comprising a particle. This array may not be
266     ! needed for all thermo problems.
267              Allocate( DES_X_s( MAX_PIP, DIMENSION_N_S))
268     ! Total rate of heat transfer to individual particles.
269              Allocate( Q_Source( MAX_PIP ) )
270     ! Average solids temperature in fluid cell
271              Allocate(avgDES_T_s(DIMENSION_3) )
272     
273              Allocate(DES_ENERGY_SOURCE(DIMENSION_3) )
274     
275     ! Allocate the history variables for Adams-Bashforth integration
276              IF (INTG_ADAMS_BASHFORTH) &
277                 Allocate( Q_Source0( MAX_PIP ) )
278           ENDIF
279     ! End Thermodynamic Allocation
280     ! ----------------------------------------------------------------<<<
281     
282     
283     ! ---------------------------------------------------------------->>>
284     ! BEGIN Species Allocation
285           IF(ANY_SPECIES_EQ)THEN
286     ! Rate of solids phase production for each species
287              Allocate( DES_R_sp( MAX_PIP, DIMENSION_N_s) )
288     ! Rate of solids phase consumption for each species
289              Allocate( DES_R_sc( MAX_PIP, DIMENSION_N_s) )
290     
291     
292              Allocate( DES_R_gp( DIMENSION_3, DIMENSION_N_g ) )
293              Allocate( DES_R_gc( DIMENSION_3, DIMENSION_N_g ) )
294              Allocate( DES_SUM_R_g( DIMENSION_3 ) )
295              Allocate( DES_R_PHASE( DIMENSION_3, DIMENSION_LM+DIMENSION_M-1 ) )
296              Allocate( DES_HOR_g( DIMENSION_3 ) )
297     
298     
299     ! Allocate the history variables for Adams-Bashforth integration
300              IF (INTG_ADAMS_BASHFORTH) THEN
301     ! Rate of change of particle mass
302                 Allocate( dMdt_OLD( MAX_PIP ) )
303     ! Rate of change of particle mass percent species
304                 Allocate( dXdt_OLD( MAX_PIP, DIMENSION_N_s) )
305              ENDIF
306     
307     ! Energy generation from reaction (cal/sec)
308              Allocate( Qint( MAX_PIP ) )
309           ENDIF
310     ! End Species Allocation
311     ! ----------------------------------------------------------------<<<
312     
313           CALL FINL_ERR_MSG
314     
315           RETURN
316           END SUBROUTINE DES_ALLOCATE_ARRAYS
317     
318     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
319     !                                                                      !
320     !  Subroutine: ALLOCATE_DEM_MIO                                        !
321     !                                                                      !
322     !  Purpose:                                                            !
323     !                                                                      !
324     !  Author: J.Musser                                   Date: 17-Aug-09  !
325     !                                                                      !
326     !  Comments:                                                           !
327     !                                                                      !
328     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
329     
330           SUBROUTINE ALLOCATE_DEM_MI
331     
332     !-----------------------------------------------
333     ! Modules
334     !-----------------------------------------------
335           USE des_bc
336           USE discretelement
337           IMPLICIT NONE
338     !-----------------------------------------------
339     
340     ! Particle injection factor
341           Allocate( PI_FACTOR (DEM_BCMI) )
342     ! Particle injection count (injection number)
343           Allocate( PI_COUNT (DEM_BCMI) )
344     ! Particle injection time scale
345           Allocate( DEM_MI_TIME (DEM_BCMI) )
346     ! Array used for polydisperse inlets: stores the particle number
347     ! distribution of an inlet scaled with numfrac_limit
348           Allocate( DEM_BC_POLY_LAYOUT( DEM_BCMI, NUMFRAC_LIMIT ) )
349     ! Data structure for storing BC data.
350           Allocate( DEM_MI(DEM_BCMI) )
351     
352     ! Initializiation
353     ! Integer arrays
354           PI_FACTOR(:) = -1
355           PI_COUNT(:) = -1
356           DEM_BC_POLY_LAYOUT(:,:) = -1
357     ! Double precision arrays
358           DEM_MI_TIME(:) = UNDEFINED
359     
360           allocate( DEM_BCMI_IJKSTART(DEM_BCMI) )
361           allocate( DEM_BCMI_IJKEND(DEM_BCMI) )
362     
363           DEM_BCMI_IJKSTART = -1
364           DEM_BCMI_IJKEND   = -1
365     
366     ! Boundary classification
367     !         Allocate( PARTICLE_PLCMNT (DES_BCMI) )
368     ! Character precision arrays
369     !         PARTICLE_PLCMNT(:) = UNDEFINED_C
370     
371           RETURN
372           END SUBROUTINE ALLOCATE_DEM_MI
373     
374     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
375     !                                                                      !
376     !  Subroutine: ALLOCATE_PIC_MIO                                        !
377     !                                                                      !
378     !  Purpose:                                                            !
379     !                                                                      !
380     !  Author: R. Garg                                    Date: 11-Jun-14  !
381     !                                                                      !
382     !  Comments:                                                           !
383     !                                                                      !
384     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
385     
386           SUBROUTINE ALLOCATE_PIC_MIO
387     
388     !-----------------------------------------------
389     ! Modules
390     !-----------------------------------------------
391           USE pic_bc
392           USE discretelement
393           IMPLICIT NONE
394     !-----------------------------------------------
395     
396     ! Allocate/Initialize for inlets
397           IF(PIC_BCMI /= 0)THEN
398     
399              allocate( PIC_BCMI_IJKSTART(PIC_BCMI) )
400              allocate( PIC_BCMI_IJKEND  (PIC_BCMI) )
401              allocate( PIC_BCMI_NORMDIR (PIC_BCMI,3) )
402     
403              ALLOCATE( PIC_BCMI_OFFSET  (PIC_BCMI,3))
404     
405              ALLOCATE( PIC_BCMI_INCL_CUTCELL(PIC_BCMI) )
406     
407              PIC_BCMI_IJKSTART = -1
408              PIC_BCMI_IJKEND   = -1
409     
410           ENDIF  ! end if PIC_BCMI /= 0
411     
412     
413     
414           IF(PIC_BCMO > 0)THEN
415              allocate( PIC_BCMO_IJKSTART(PIC_BCMO) )
416              allocate( PIC_BCMO_IJKEND(PIC_BCMO) )
417     
418              PIC_BCMO_IJKSTART = -1
419              PIC_BCMO_IJKEND   = -1
420           ENDIF
421     
422     
423           RETURN
424           END SUBROUTINE ALLOCATE_PIC_MIO
425     
426     !``````````````````````````````````````````````````````````````````````!
427     ! Subroutine: ADD_PAIR                                                 !
428     !                                                                      !
429     ! Purpose: Adds a neighbor pair to the pairs array.                    !
430     !                                                                      !
431     !``````````````````````````````````````````````````````````````````````!
432           DOUBLE PRECISION FUNCTION add_pair(ii,jj)
433           USE discretelement
434           IMPLICIT NONE
435           INTEGER, INTENT(IN) :: ii,jj
436     
437           if (NEIGHBOR_INDEX(ii) > NEIGH_MAX) then
438              NEIGH_MAX = 2*NEIGH_MAX
439              CALL NEIGHBOR_GROW(NEIGH_MAX)
440           endif
441     
442           NEIGHBORS(NEIGHBOR_INDEX(ii)) = jj
443           NEIGHBOR_INDEX(ii) = NEIGHBOR_INDEX(ii) + 1
444           add_pair = NEIGHBOR_INDEX(ii)
445     
446           RETURN
447           END FUNCTION add_pair
448     
449     !``````````````````````````````````````````````````````````````````````!
450     ! Subroutine: NEIGHBOR_GROW                                            !
451     !                                                                      !
452     ! Purpose: Grow neighbors arrays to new_neigh_max. Note that neighbor      !
453     ! max should be increased before calling this routine. Also, no        !
454     ! assumption to the previous array size is made as needed for restarts.!
455     !``````````````````````````````````````````````````````````````````````!
456           SUBROUTINE NEIGHBOR_GROW(new_neigh_max)
457             USE discretelement
458             USE geometry
459             IMPLICIT NONE
460     
461             integer, intent(in) :: new_neigh_max
462     
463             INTEGER :: lSIZE1
464             INTEGER, DIMENSION(:), ALLOCATABLE :: neigh_tmp
465             DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: pf_tmp
466     
467             lSIZE1 = size(neighbors,1)
468     
469             allocate(neigh_tmp(new_neigh_max))
470             neigh_tmp(1:lSIZE1) = neighbors(1:lSIZE1)
471             neigh_tmp(lSIZE1+1:) = 0
472             call move_alloc(neigh_tmp,neighbors)
473     
474             allocate(neigh_tmp(new_neigh_max))
475             neigh_tmp(1:lSIZE1) = neighbors_old(1:lSIZE1)
476             neigh_tmp(lSIZE1+1:) = 0
477             call move_alloc(neigh_tmp,neighbors_old)
478     
479             allocate(pf_tmp(3,new_neigh_max))
480             pf_tmp(:,1:lSIZE1) = pft_neighbor(:,1:lSIZE1)
481             pf_tmp(:,lSIZE1+1:) = 0
482             call move_alloc(pf_tmp,pft_neighbor)
483     
484             allocate(pf_tmp(3,new_neigh_max))
485             pf_tmp(:,1:lSIZE1) = pft_neighbor_old(:,1:lSIZE1)
486             pf_tmp(:,lSIZE1+1:) = 0
487             call move_alloc(pf_tmp,pft_neighbor_old)
488     
489     
490           END SUBROUTINE NEIGHBOR_GROW
491     
492     !``````````````````````````````````````````````````````````````````````!
493     ! Subroutine: PARTICLE_GROW                                            !
494     !                                                                      !
495     ! Purpose: Grow particle arrays to new_max_pip. Note that pair         !
496     ! max should be increased before calling this routine. Also, no        !
497     ! assumption to the previous array size is made as needed for restarts.!
498     !``````````````````````````````````````````````````````````````````````!
499           SUBROUTINE PARTICLE_GROW(new_max_pip)
500     
501             USE des_rxns
502             USE des_thermo
503             USE mfix_pic
504             USE discretelement
505             USE particle_filter
506             USE run
507     
508             IMPLICIT NONE
509     
510             integer, intent(in) :: new_max_pip
511     
512             DO WHILE (MAX_PIP < new_max_pip)
513                MAX_PIP = MAX_PIP*2
514     
515                call real_grow(des_radius,MAX_PIP)
516                call real_grow(RO_Sol,MAX_PIP)
517                call real_grow(PVOL,MAX_PIP)
518                call real_grow(PMASS,MAX_PIP)
519                call real_grow(OMOI,MAX_PIP)
520                call real_grow2(DES_POS_NEW,MAX_PIP)
521                call real_grow2(DES_VEL_NEW,MAX_PIP)
522                call real_grow2(OMEGA_NEW,MAX_PIP)
523                call real_grow2(PPOS,MAX_PIP)
524                call byte_grow(PARTICLE_STATE,MAX_PIP)
525                call integer_grow(iglobal_id,MAX_PIP)
526                call integer_grow2_reverse(pijk,MAX_PIP)
527                call integer_grow(dg_pijk,MAX_PIP)
528                call integer_grow(dg_pijkprv,MAX_PIP)
529                call logical_grow(ighost_updated,MAX_PIP)
530                call real_grow2(FC,MAX_PIP)
531                call real_grow2(TOW,MAX_PIP)
532                call real_grow(F_GP,MAX_PIP)
533                call integer_grow2(WALL_COLLISION_FACET_ID,MAX_PIP)
534                call real_grow3(WALL_COLLISION_PFT,MAX_PIP)
535                call real_grow2(DRAG_FC,MAX_PIP)
536     
537                call integer_grow(NEIGHBOR_INDEX,MAX_PIP)
538                call integer_grow(NEIGHBOR_INDEX_OLD,MAX_PIP)
539     
540                IF(PARTICLE_ORIENTATION) call real_grow2(ORIENTATION,MAX_PIP)
541     
542                IF(FILTER_SIZE > 0) THEN
543                   call integer_grow2(FILTER_CELL,MAX_PIP)
544                   call real_grow2(FILTER_WEIGHT,MAX_PIP)
545                ENDIF
546     
547                IF(MPPIC) THEN
548                   call real_grow(DES_STAT_WT,MAX_PIP)
549                   call real_grow2_reverse(PS_GRAD,MAX_PIP)
550                   call real_grow2(AVGSOLVEL_P,MAX_PIP)
551                   call real_grow(EPG_P,MAX_PIP)
552                ENDIF
553     
554                IF(USE_COHESION) THEN
555                   call real_grow(PostCohesive,MAX_PIP)
556                ENDIF
557     
558                IF (DO_OLD) THEN
559                   call real_grow2(DES_POS_OLD,MAX_PIP)
560                   call real_grow2(DES_VEL_OLD,MAX_PIP)
561                   call real_grow2(DES_ACC_OLD,MAX_PIP)
562                   call real_grow2(OMEGA_OLD,MAX_PIP)
563                   call real_grow2(ROT_ACC_OLD,MAX_PIP)
564                ENDIF
565     
566                IF(ENERGY_EQ)THEN
567                   call real_grow(DES_T_s_OLD,MAX_PIP)
568                   call real_grow(DES_T_s_NEW,MAX_PIP)
569                   call real_grow(DES_C_PS,MAX_PIP)
570                   call real_grow2_reverse(DES_X_s,MAX_PIP)
571                   call real_grow(Q_Source,MAX_PIP)
572     
573                   IF (INTG_ADAMS_BASHFORTH) &
574                        call real_grow(Q_Source0,MAX_PIP)
575                ENDIF
576     
577                IF(ANY_SPECIES_EQ)THEN
578                   call real_grow2_reverse( DES_R_sp, MAX_PIP )
579                   call real_grow2_reverse( DES_R_sc, MAX_PIP )
580     
581                   IF (INTG_ADAMS_BASHFORTH) THEN
582                      call real_grow( dMdt_OLD, MAX_PIP )
583                      call real_grow2_reverse( dXdt_OLD, MAX_PIP )
584                   ENDIF
585     
586                   call real_grow( Qint, MAX_PIP )
587                ENDIF
588     
589                IF(DES_USR_VAR_SIZE > 0) &
590                   call real_grow2(DES_USR_VAR,MAX_PIP)
591     
592                CALL DES_INIT_PARTICLE_ARRAYS(MAX_PIP/2+1,MAX_PIP)
593     
594             ENDDO
595     
596           RETURN
597     
598           END SUBROUTINE PARTICLE_GROW
599     
600           SUBROUTINE BYTE_GROW(byte_array,new_size)
601             IMPLICIT NONE
602     
603             INTEGER, INTENT(IN) :: new_size
604             INTEGER(KIND=1), DIMENSION(:), ALLOCATABLE, INTENT(INOUT) :: byte_array
605             INTEGER(KIND=1), DIMENSION(:), ALLOCATABLE :: byte_tmp
606             INTEGER lSIZE
607     
608             lSIZE = size(byte_array,1)
609             allocate(byte_tmp(new_size))
610             byte_tmp(1:lSIZE) = byte_array(1:lSIZE)
611             call move_alloc(byte_tmp,byte_array)
612     
613           END SUBROUTINE BYTE_GROW
614     
615           SUBROUTINE INTEGER_GROW(integer_array,new_size)
616             IMPLICIT NONE
617     
618             INTEGER, INTENT(IN) :: new_size
619             INTEGER, DIMENSION(:), ALLOCATABLE, INTENT(INOUT) :: integer_array
620             INTEGER, DIMENSION(:), ALLOCATABLE :: integer_tmp
621             INTEGER lSIZE
622     
623             lSIZE = size(integer_array,1)
624             allocate(integer_tmp(new_size))
625             integer_tmp(1:lSIZE) = integer_array(1:lSIZE)
626             call move_alloc(integer_tmp,integer_array)
627     
628           END SUBROUTINE INTEGER_GROW
629     
630           SUBROUTINE INTEGER_GROW2_reverse(integer_array,new_size)
631             IMPLICIT NONE
632     
633             INTEGER, INTENT(IN) :: new_size
634             INTEGER, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: integer_array
635             INTEGER, DIMENSION(:,:), ALLOCATABLE :: integer_tmp
636             INTEGER lSIZE, lSIZE2
637     
638             lSIZE = size(integer_array,1)
639             lSIZE2 = size(integer_array,2)
640             allocate(integer_tmp(new_size,lSIZE2))
641             integer_tmp(1:lSIZE,:) = integer_array(1:lSIZE,:)
642             call move_alloc(integer_tmp,integer_array)
643     
644           END SUBROUTINE INTEGER_GROW2_reverse
645     
646           SUBROUTINE INTEGER_GROW2(integer_array,new_size)
647             IMPLICIT NONE
648     
649             INTEGER, INTENT(IN) :: new_size
650             INTEGER, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: integer_array
651             INTEGER, DIMENSION(:,:), ALLOCATABLE :: integer_tmp
652             INTEGER lSIZE, lSIZE2
653     
654             lSIZE = size(integer_array,1)
655             lSIZE2 = size(integer_array,2)
656             allocate(integer_tmp(lSIZE,new_size))
657             integer_tmp(:,1:lSIZE2) = integer_array(:,1:lSIZE2)
658             call move_alloc(integer_tmp,integer_array)
659     
660           END SUBROUTINE INTEGER_GROW2
661     
662           SUBROUTINE LOGICAL_GROW(logical_array,new_size)
663             IMPLICIT NONE
664     
665             INTEGER, INTENT(IN) :: new_size
666             LOGICAL, DIMENSION(:), ALLOCATABLE, INTENT(INOUT) :: logical_array
667             LOGICAL, DIMENSION(:), ALLOCATABLE :: logical_tmp
668             INTEGER lSIZE
669     
670             lSIZE = size(logical_array,1)
671             allocate(logical_tmp(new_size))
672             logical_tmp(1:lSIZE) = logical_array(1:lSIZE)
673             call move_alloc(logical_tmp,logical_array)
674     
675           END SUBROUTINE LOGICAL_GROW
676     
677           SUBROUTINE LOGICAL_GROW2(logical_array,new_size)
678             IMPLICIT NONE
679     
680             INTEGER, INTENT(IN) :: new_size
681             LOGICAL, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: logical_array
682             LOGICAL, DIMENSION(:,:), ALLOCATABLE :: logical_tmp
683             INTEGER lSIZE, lSIZE2
684     
685             lSIZE = size(logical_array,1)
686             lSIZE2 = size(logical_array,2)
687             allocate(logical_tmp(lSIZE,new_size))
688             logical_tmp(:,1:lSIZE2) = logical_array(:,1:lSIZE2)
689             call move_alloc(logical_tmp,logical_array)
690     
691           END SUBROUTINE LOGICAL_GROW2
692     
693           SUBROUTINE REAL_GROW(real_array,new_size)
694             IMPLICIT NONE
695     
696             INTEGER, INTENT(IN) :: new_size
697             DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE, INTENT(INOUT) :: real_array
698             DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: real_tmp
699             INTEGER lSIZE
700     
701             lSIZE = size(real_array,1)
702             allocate(real_tmp(new_size))
703             real_tmp(1:lSIZE) = real_array(1:lSIZE)
704             call move_alloc(real_tmp,real_array)
705     
706           END SUBROUTINE REAL_GROW
707     
708           SUBROUTINE REAL_GROW2(real_array,new_size)
709             IMPLICIT NONE
710     
711             INTEGER, INTENT(IN) :: new_size
712             DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: real_array
713             DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: real_tmp
714             INTEGER lSIZE, lSIZE2
715     
716             lSIZE = size(real_array,1)
717             lSIZE2 = size(real_array,2)
718             allocate(real_tmp(lSIZE,new_size))
719             real_tmp(:,1:lSIZE2) = real_array(:,1:lSIZE2)
720             call move_alloc(real_tmp,real_array)
721     
722           END SUBROUTINE REAL_GROW2
723     
724           SUBROUTINE REAL_GROW2_reverse(real_array,new_size)
725             IMPLICIT NONE
726     
727             INTEGER, INTENT(IN) :: new_size
728             DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: real_array
729             DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: real_tmp
730             INTEGER lSIZE, lSIZE2
731     
732             lSIZE = size(real_array,1)
733             lSIZE2 = size(real_array,2)
734             allocate(real_tmp(new_size,lSIZE2))
735             real_tmp(1:lSIZE,:) = real_array(1:lSIZE,:)
736             call move_alloc(real_tmp,real_array)
737     
738           END SUBROUTINE REAL_GROW2_REVERSE
739     
740           SUBROUTINE REAL_GROW3(real_array,new_size)
741             IMPLICIT NONE
742     
743             INTEGER, INTENT(IN) :: new_size
744             DOUBLE PRECISION, DIMENSION(:,:,:), ALLOCATABLE, INTENT(INOUT) :: real_array
745             DOUBLE PRECISION, DIMENSION(:,:,:), ALLOCATABLE :: real_tmp
746             INTEGER lSIZE, lSIZE2, lSIZE3
747     
748             lSIZE = size(real_array,1)
749             lSIZE2 = size(real_array,2)
750             lSIZE3 = size(real_array,3)
751             allocate(real_tmp(lSIZE,lSIZE2,new_size))
752             real_tmp(:,:,1:lSIZE3) = real_array(:,:,1:lSIZE3)
753             call move_alloc(real_tmp,real_array)
754     
755           END SUBROUTINE REAL_GROW3
756     
757           SUBROUTINE LOGICAL_GROW2_REVERSE(real_array,new_size)
758             IMPLICIT NONE
759     
760             INTEGER, INTENT(IN) :: new_size
761             LOGICAL, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: real_array
762             LOGICAL, DIMENSION(:,:), ALLOCATABLE :: real_tmp
763             INTEGER lSIZE, lSIZE2
764     
765             lSIZE = size(real_array,1)
766             lSIZE2 = size(real_array,2)
767             allocate(real_tmp(new_size,lSIZE2))
768             real_tmp(1:lSIZE,:) = real_array(1:lSIZE,:)
769             call move_alloc(real_tmp,real_array)
770     
771           END SUBROUTINE LOGICAL_GROW2_REVERSE
772     
773         END MODULE DES_ALLOCATE
774