File: /nfs/home/0/users/jenkins/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, PAIR_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( PEA (MAX_PIP, 4) )
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           OLD_PAIR_NUM = 0
145           PAIR_NUM = 0
146           PAIR_MAX = 1024
147     
148           Allocate(  PAIRS (2,PAIR_MAX) )
149           Allocate(  PAIRS_OLD (2,PAIR_MAX) )
150           Allocate(  PV_PAIR (PAIR_MAX) )
151           Allocate(  PV_PAIR_OLD (PAIR_MAX) )
152           Allocate(  PFT_PAIR (3,PAIR_MAX) )
153           Allocate(  PFT_PAIR_OLD (3,PAIR_MAX) )
154           Allocate(  PFN_PAIR (3,PAIR_MAX) )
155           Allocate(  PFN_PAIR_OLD (3,PAIR_MAX) )
156     
157     ! Variable that stores the particle in cell information (ID) on the
158     ! computational fluid grid defined by imax, jmax and kmax in mfix.dat
159           ALLOCATE(PIC(DIMENSION_3))
160           DO IJK=1,DIMENSION_3
161             NULLIFY(pic(ijk)%p)
162           ENDDO
163     
164     ! Particles in a computational fluid cell (for volume fraction)
165           Allocate(  PINC (DIMENSION_3) )
166     
167     ! For each particle track its i,j,k location on computational fluid grid
168     ! defined by imax, jmax and kmax in mfix.dat and phase no.
169           Allocate(  PIJK (MAX_PIP,5) )
170     
171           ALLOCATE(DRAG_AM(DIMENSION_3))
172           ALLOCATE(DRAG_BM(DIMENSION_3, DIMN))
173           ALLOCATE(F_gp(MAX_PIP ))
174           F_gp(1:MAX_PIP)  = ZERO
175     
176     ! Explict drag force acting on a particle.
177           Allocate(DRAG_FC (DIMN,MAX_PIP) )
178     
179     ! force due to gas-pressure gradient
180           ALLOCATE(P_FORCE(DIMN, DIMENSION_3))
181     
182     ! Volume averaged solids volume in a computational fluid cell
183           Allocate(  DES_U_s (DIMENSION_3, DES_MMAX) )
184           Allocate(  DES_V_s (DIMENSION_3, DES_MMAX) )
185           Allocate(  DES_W_s (DIMENSION_3, DES_MMAX) )
186     
187     ! Volume of nodes
188           ALLOCATE(DES_VOL_NODE(DIMENSION_3))
189     
190           ALLOCATE(F_GDS(DIMENSION_3))
191           ALLOCATE(VXF_GDS(DIMENSION_3))
192     
193           SELECT CASE(DES_INTERP_SCHEME_ENUM)
194           CASE(DES_INTERP_DPVM, DES_INTERP_GAUSS)
195              ALLOCATE(FILTER_CELL(FILTER_SIZE, MAX_PIP))
196              ALLOCATE(FILTER_WEIGHT(FILTER_SIZE, MAX_PIP))
197           CASE(DES_INTERP_GARG)
198              ALLOCATE(DES_ROPS_NODE(DIMENSION_3, DES_MMAX))
199              ALLOCATE(DES_VEL_NODE(DIMENSION_3, DIMN, DES_MMAX))
200           END SELECT
201     
202     ! Variables for hybrid model
203           IF (DES_CONTINUUM_HYBRID) THEN
204              ALLOCATE(SDRAG_AM(DIMENSION_3,DIMENSION_M))
205              ALLOCATE(SDRAG_BM(DIMENSION_3, DIMN,DIMENSION_M))
206     
207              ALLOCATE(F_SDS(DIMENSION_3,DIMENSION_M))
208              ALLOCATE(VXF_SDS(DIMENSION_3,DIMENSION_M))
209           ENDIF
210     ! Bulk density in a computational fluid cell / for communication with
211     ! MFIX continuum
212           ALLOCATE( DES_ROP_S(DIMENSION_3, DES_MMAX) )
213           ALLOCATE( DES_ROP_SO(DIMENSION_3, DES_MMAX) )
214     
215     ! MP-PIC related
216           IF(MPPIC) THEN
217              Allocate(PS_FORCE_PIC(DIMENSION_3, DIMN))
218              ALLOCATE(DES_STAT_WT(MAX_PIP))
219              ALLOCATE(DES_VEL_MAX(DIMN))
220              ALLOCATE(PS_GRAD(MAX_PIP, DIMN))
221              ALLOCATE(AVGSOLVEL_P(DIMN, MAX_PIP))
222              ALLOCATE(EPG_P(MAX_PIP))
223     
224              Allocate(PIC_U_S(DIMENSION_3, DES_MMAX))
225              Allocate(PIC_V_S(DIMENSION_3, DES_MMAX))
226              Allocate(PIC_W_S(DIMENSION_3, DES_MMAX))
227     
228              Allocate(PIC_P_s (DIMENSION_3, DES_MMAX) )
229     !         ALLOCATE(MPPIC_VPTAU(MAX_PIP, DIMN))
230              PIC_U_s = zero
231              PIC_V_s = zero
232              PIC_W_s = zero
233              PIC_P_s = zero
234           ENDIF
235     
236     
237     ! Granular temperature in a computational fluid cell
238           Allocate(DES_THETA (DIMENSION_3, DES_MMAX) )
239     
240     ! Averaged velocity obtained by averaging over all the particles
241           ALLOCATE(DES_VEL_AVG(DIMN) )
242     
243     ! Global Granular Energy
244           ALLOCATE(GLOBAL_GRAN_ENERGY(DIMN) )
245           ALLOCATE(GLOBAL_GRAN_TEMP(DIMN) )
246     
247     ! variable for bed height of solids phase M
248           ALLOCATE(BED_HEIGHT(DES_MMAX))
249     
250     ! ---------------------------------------------------------------->>>
251     ! BEGIN COHESION
252           IF(USE_COHESION) THEN
253     ! Matrix location of particle  (should be allocated in case user wishes
254     ! to invoke routines in /cohesion subdirectory
255              Allocate(  PostCohesive (MAX_PIP) )
256           ENDIF
257     ! END COHESION
258     ! ----------------------------------------------------------------<<<
259     
260     ! ---------------------------------------------------------------->>>
261     ! BEGIN Thermodynamic Allocation
262           IF(ENERGY_EQ)THEN
263     ! Particle temperature
264              Allocate( DES_T_s_OLD( MAX_PIP ) )
265              Allocate( DES_T_s_NEW( MAX_PIP ) )
266     ! Specific heat
267              Allocate( DES_C_PS( MAX_PIP ) )
268     ! Species mass fractions comprising a particle. This array may not be
269     ! needed for all thermo problems.
270              Allocate( DES_X_s( MAX_PIP, DIMENSION_N_S))
271     ! Total rate of heat transfer to individual particles.
272              Allocate( Q_Source( MAX_PIP ) )
273     ! Average solids temperature in fluid cell
274              Allocate(avgDES_T_s(DIMENSION_3) )
275     
276              Allocate(DES_ENERGY_SOURCE(DIMENSION_3) )
277     
278     ! Allocate the history variables for Adams-Bashforth integration
279              IF (INTG_ADAMS_BASHFORTH) &
280                 Allocate( Q_Source0( MAX_PIP ) )
281           ENDIF
282     ! End Thermodynamic Allocation
283     ! ----------------------------------------------------------------<<<
284     
285     
286     ! ---------------------------------------------------------------->>>
287     ! BEGIN Species Allocation
288           IF(ANY_SPECIES_EQ)THEN
289     ! Rate of solids phase production for each species
290              Allocate( DES_R_sp( MAX_PIP, DIMENSION_N_s) )
291     ! Rate of solids phase consumption for each species
292              Allocate( DES_R_sc( MAX_PIP, DIMENSION_N_s) )
293     
294     
295              Allocate( DES_R_gp( DIMENSION_3, DIMENSION_N_g ) )
296              Allocate( DES_R_gc( DIMENSION_3, DIMENSION_N_g ) )
297              Allocate( DES_SUM_R_g( DIMENSION_3 ) )
298              Allocate( DES_R_PHASE( DIMENSION_3, DIMENSION_LM+DIMENSION_M-1 ) )
299              Allocate( DES_HOR_g( DIMENSION_3 ) )
300     
301     
302     ! Allocate the history variables for Adams-Bashforth integration
303              IF (INTG_ADAMS_BASHFORTH) THEN
304     ! Rate of change of particle mass
305                 Allocate( dMdt_OLD( MAX_PIP ) )
306     ! Rate of change of particle mass percent species
307                 Allocate( dXdt_OLD( MAX_PIP, DIMENSION_N_s) )
308              ENDIF
309     
310     ! Energy generation from reaction (cal/sec)
311              Allocate( Qint( MAX_PIP ) )
312           ENDIF
313     ! End Species Allocation
314     ! ----------------------------------------------------------------<<<
315     
316           CALL FINL_ERR_MSG
317     
318           RETURN
319           END SUBROUTINE DES_ALLOCATE_ARRAYS
320     
321     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
322     !                                                                      !
323     !  Subroutine: ALLOCATE_DEM_MIO                                        !
324     !                                                                      !
325     !  Purpose:                                                            !
326     !                                                                      !
327     !  Author: J.Musser                                   Date: 17-Aug-09  !
328     !                                                                      !
329     !  Comments:                                                           !
330     !                                                                      !
331     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
332     
333           SUBROUTINE ALLOCATE_DEM_MI
334     
335     !-----------------------------------------------
336     ! Modules
337     !-----------------------------------------------
338           USE des_bc
339           USE discretelement
340           IMPLICIT NONE
341     !-----------------------------------------------
342     
343     ! Particle injection factor
344           Allocate( PI_FACTOR (DEM_BCMI) )
345     ! Particle injection count (injection number)
346           Allocate( PI_COUNT (DEM_BCMI) )
347     ! Particle injection time scale
348           Allocate( DEM_MI_TIME (DEM_BCMI) )
349     ! Array used for polydisperse inlets: stores the particle number
350     ! distribution of an inlet scaled with numfrac_limit
351           Allocate( DEM_BC_POLY_LAYOUT( DEM_BCMI, NUMFRAC_LIMIT ) )
352     ! Data structure for storing BC data.
353           Allocate( DEM_MI(DEM_BCMI) )
354     
355     ! Initializiation
356     ! Integer arrays
357           PI_FACTOR(:) = -1
358           PI_COUNT(:) = -1
359           DEM_BC_POLY_LAYOUT(:,:) = -1
360     ! Double precision arrays
361           DEM_MI_TIME(:) = UNDEFINED
362     
363           allocate( DEM_BCMI_IJKSTART(DEM_BCMI) )
364           allocate( DEM_BCMI_IJKEND(DEM_BCMI) )
365     
366           DEM_BCMI_IJKSTART = -1
367           DEM_BCMI_IJKEND   = -1
368     
369     ! Boundary classification
370     !         Allocate( PARTICLE_PLCMNT (DES_BCMI) )
371     ! Character precision arrays
372     !         PARTICLE_PLCMNT(:) = UNDEFINED_C
373     
374           RETURN
375           END SUBROUTINE ALLOCATE_DEM_MI
376     
377     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
378     !                                                                      !
379     !  Subroutine: ALLOCATE_PIC_MIO                                        !
380     !                                                                      !
381     !  Purpose:                                                            !
382     !                                                                      !
383     !  Author: R. Garg                                    Date: 11-Jun-14  !
384     !                                                                      !
385     !  Comments:                                                           !
386     !                                                                      !
387     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
388     
389           SUBROUTINE ALLOCATE_PIC_MIO
390     
391     !-----------------------------------------------
392     ! Modules
393     !-----------------------------------------------
394           USE pic_bc
395           USE discretelement
396           IMPLICIT NONE
397     !-----------------------------------------------
398     
399     ! Allocate/Initialize for inlets
400           IF(PIC_BCMI /= 0)THEN
401     
402              allocate( PIC_BCMI_IJKSTART(PIC_BCMI) )
403              allocate( PIC_BCMI_IJKEND  (PIC_BCMI) )
404              allocate( PIC_BCMI_NORMDIR (PIC_BCMI,3) )
405     
406              ALLOCATE( PIC_BCMI_OFFSET  (PIC_BCMI,3))
407     
408              ALLOCATE( PIC_BCMI_INCL_CUTCELL(PIC_BCMI) )
409     
410              PIC_BCMI_IJKSTART = -1
411              PIC_BCMI_IJKEND   = -1
412     
413           ENDIF  ! end if PIC_BCMI /= 0
414     
415     
416     
417           IF(PIC_BCMO > 0)THEN
418              allocate( PIC_BCMO_IJKSTART(PIC_BCMO) )
419              allocate( PIC_BCMO_IJKEND(PIC_BCMO) )
420     
421              PIC_BCMO_IJKSTART = -1
422              PIC_BCMO_IJKEND   = -1
423           ENDIF
424     
425     
426           RETURN
427           END SUBROUTINE ALLOCATE_PIC_MIO
428     
429     
430     !``````````````````````````````````````````````````````````````````````!
431     ! Subroutine: ADD_PAIR                                                 !
432     !                                                                      !
433     ! Purpose: Adds a neighbor pair to the pairs array.                    !
434     !                                                                      !
435     !``````````````````````````````````````````````````````````````````````!
436           SUBROUTINE add_pair(ii,jj)
437           USE discretelement
438           USE geometry
439           IMPLICIT NONE
440           INTEGER, INTENT(IN) :: ii,jj
441     
442           pair_num = pair_num +1
443     
444     ! Reallocate to double the size of the arrays if needed.
445           IF(PAIR_NUM > PAIR_MAX) THEN
446              PAIR_MAX = PAIR_MAX*2
447              CALL PAIR_GROW
448           ENDIF
449     
450           pairs(1,pair_num) = ii
451           pairs(2,pair_num) = jj
452     
453           RETURN
454           END SUBROUTINE add_pair
455     
456     !``````````````````````````````````````````````````````````````````````!
457     ! Subroutine: PAIR_GROW                                                !
458     !                                                                      !
459     ! Purpose: Grow pair arrays to pair_max. Note that pair                !
460     ! max should be increased before calling this routine. Also, no        !
461     ! assumption to the previous array size is made as needed for restarts.!
462     !``````````````````````````````````````````````````````````````````````!
463           SUBROUTINE PAIR_GROW
464     
465             USE discretelement
466     
467             IMPLICIT NONE
468     
469             call integer_grow2(pairs,pair_max)
470             call integer_grow2(pairs_old,pair_max)
471     
472             call logical_grow(pv_pair,pair_max)
473             call logical_grow(pv_pair_old,pair_max)
474     
475             call real_grow2(pft_pair,pair_max)
476             call real_grow2(pft_pair_old,pair_max)
477             call real_grow2(pfn_pair,pair_max)
478             call real_grow2(pfn_pair_old,pair_max)
479     
480             RETURN
481     
482           END SUBROUTINE PAIR_GROW
483     
484     !``````````````````````````````````````````````````````````````````````!
485     ! Subroutine: PARTICLE_GROW                                            !
486     !                                                                      !
487     ! Purpose: Grow pair arrays to pair_max. Note that pair                !
488     ! max should be increased before calling this routine. Also, no        !
489     ! assumption to the previous array size is made as needed for restarts.!
490     !``````````````````````````````````````````````````````````````````````!
491           SUBROUTINE PARTICLE_GROW(new_max_pip)
492     
493     !      use discretelement, only: pea, wall_collision_facet_id, wall_collision_pft, MAX_PIP, PIP
494     !      use param1
495     
496             USE des_rxns
497             USE des_thermo
498             USE mfix_pic
499             USE discretelement
500             USE particle_filter
501             USE run
502     
503             IMPLICIT NONE
504     
505             integer, intent(in) :: new_max_pip
506     
507             DO WHILE (MAX_PIP < new_max_pip)
508                MAX_PIP = MAX_PIP*2
509     
510                call real_grow(des_radius,MAX_PIP)
511                call real_grow(RO_Sol,MAX_PIP)
512                call real_grow(PVOL,MAX_PIP)
513                call real_grow(PMASS,MAX_PIP)
514                call real_grow(OMOI,MAX_PIP)
515                call real_grow2(DES_POS_NEW,MAX_PIP)
516                call real_grow2(DES_VEL_NEW,MAX_PIP)
517                call real_grow2(OMEGA_NEW,MAX_PIP)
518                IF(PARTICLE_ORIENTATION) call real_grow2(ORIENTATION,MAX_PIP)
519                call real_grow2(PPOS,MAX_PIP)
520                call logical_grow2_reverse(PEA,MAX_PIP)
521                call integer_grow(iglobal_id,MAX_PIP)
522                call integer_grow2_reverse(pijk,MAX_PIP)
523                call integer_grow(dg_pijk,MAX_PIP)
524                call integer_grow(dg_pijkprv,MAX_PIP)
525                call logical_grow(ighost_updated,MAX_PIP)
526                call real_grow2(FC,MAX_PIP)
527                call real_grow2(TOW,MAX_PIP)
528                call real_grow(F_GP,MAX_PIP)
529                call integer_grow2(WALL_COLLISION_FACET_ID,MAX_PIP)
530                call real_grow3(WALL_COLLISION_PFT,MAX_PIP)
531                call real_grow2(DRAG_FC,MAX_PIP)
532     
533                IF(FILTER_SIZE > 0) THEN
534                   call integer_grow2(FILTER_CELL,MAX_PIP)
535                   call real_grow2(FILTER_WEIGHT,MAX_PIP)
536                ENDIF
537     
538                IF(MPPIC) THEN
539                   call real_grow(DES_STAT_WT,MAX_PIP)
540                   call real_grow2_reverse(PS_GRAD,MAX_PIP)
541                   call real_grow2(AVGSOLVEL_P,MAX_PIP)
542                   call real_grow(EPG_P,MAX_PIP)
543                ENDIF
544     
545                IF(USE_COHESION) THEN
546                   call real_grow(PostCohesive,MAX_PIP)
547                ENDIF
548     
549                IF (DO_OLD) THEN
550                   call real_grow2(DES_POS_OLD,MAX_PIP)
551                   call real_grow2(DES_VEL_OLD,MAX_PIP)
552                   call real_grow2(DES_ACC_OLD,MAX_PIP)
553                   call real_grow2(OMEGA_OLD,MAX_PIP)
554                   call real_grow2(ROT_ACC_OLD,MAX_PIP)
555                ENDIF
556     
557                IF(ENERGY_EQ)THEN
558                   call real_grow(DES_T_s_OLD,MAX_PIP)
559                   call real_grow(DES_T_s_NEW,MAX_PIP)
560                   call real_grow(DES_C_PS,MAX_PIP)
561                   call real_grow2_reverse(DES_X_s,MAX_PIP)
562                   call real_grow(Q_Source,MAX_PIP)
563     
564                   IF (INTG_ADAMS_BASHFORTH) &
565                        call real_grow(Q_Source0,MAX_PIP)
566                ENDIF
567     
568                IF(ANY_SPECIES_EQ)THEN
569                   call real_grow2_reverse( DES_R_sp, MAX_PIP )
570                   call real_grow2_reverse( DES_R_sc, MAX_PIP )
571     
572                   IF (INTG_ADAMS_BASHFORTH) THEN
573                      call real_grow( dMdt_OLD, MAX_PIP )
574                      call real_grow2_reverse( dXdt_OLD, MAX_PIP )
575                   ENDIF
576     
577                   call real_grow( Qint, MAX_PIP )
578                ENDIF
579     
580                IF(DES_USR_VAR_SIZE > 0) &
581                   call real_grow2(DES_USR_VAR,MAX_PIP)
582     
583                CALL DES_INIT_PARTICLE_ARRAYS(MAX_PIP/2+1,MAX_PIP)
584     
585             ENDDO
586     
587           RETURN
588     
589           END SUBROUTINE PARTICLE_GROW
590     
591           SUBROUTINE INTEGER_GROW(integer_array,new_size)
592             IMPLICIT NONE
593     
594             INTEGER, INTENT(IN) :: new_size
595             INTEGER, DIMENSION(:), ALLOCATABLE, INTENT(INOUT) :: integer_array
596             INTEGER, DIMENSION(:), ALLOCATABLE :: integer_tmp
597             INTEGER lSIZE
598     
599             lSIZE = size(integer_array,1)
600             allocate(integer_tmp(new_size))
601             integer_tmp(1:lSIZE) = integer_array(1:lSIZE)
602             call move_alloc(integer_tmp,integer_array)
603     
604           END SUBROUTINE INTEGER_GROW
605     
606           SUBROUTINE INTEGER_GROW2_reverse(integer_array,new_size)
607             IMPLICIT NONE
608     
609             INTEGER, INTENT(IN) :: new_size
610             INTEGER, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: integer_array
611             INTEGER, DIMENSION(:,:), ALLOCATABLE :: integer_tmp
612             INTEGER lSIZE, lSIZE2
613     
614             lSIZE = size(integer_array,1)
615             lSIZE2 = size(integer_array,2)
616             allocate(integer_tmp(new_size,lSIZE2))
617             integer_tmp(1:lSIZE,:) = integer_array(1:lSIZE,:)
618             call move_alloc(integer_tmp,integer_array)
619     
620           END SUBROUTINE INTEGER_GROW2_reverse
621     
622           SUBROUTINE INTEGER_GROW2(integer_array,new_size)
623             IMPLICIT NONE
624     
625             INTEGER, INTENT(IN) :: new_size
626             INTEGER, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: integer_array
627             INTEGER, DIMENSION(:,:), ALLOCATABLE :: integer_tmp
628             INTEGER lSIZE, lSIZE2
629     
630             lSIZE = size(integer_array,1)
631             lSIZE2 = size(integer_array,2)
632             allocate(integer_tmp(lSIZE,new_size))
633             integer_tmp(:,1:lSIZE2) = integer_array(:,1:lSIZE2)
634             call move_alloc(integer_tmp,integer_array)
635     
636           END SUBROUTINE INTEGER_GROW2
637     
638           SUBROUTINE LOGICAL_GROW(logical_array,new_size)
639             IMPLICIT NONE
640     
641             INTEGER, INTENT(IN) :: new_size
642             LOGICAL, DIMENSION(:), ALLOCATABLE, INTENT(INOUT) :: logical_array
643             LOGICAL, DIMENSION(:), ALLOCATABLE :: logical_tmp
644             INTEGER lSIZE2
645     
646             lSIZE2 = size(logical_array,1)
647             allocate(logical_tmp(new_size))
648             logical_tmp(1:lSIZE2) = logical_array(1:lSIZE2)
649             call move_alloc(logical_tmp,logical_array)
650     
651           END SUBROUTINE LOGICAL_GROW
652     
653           SUBROUTINE REAL_GROW(real_array,new_size)
654             IMPLICIT NONE
655     
656             INTEGER, INTENT(IN) :: new_size
657             DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE, INTENT(INOUT) :: real_array
658             DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE :: real_tmp
659             INTEGER lSIZE
660     
661             lSIZE = size(real_array,1)
662             allocate(real_tmp(new_size))
663             real_tmp(1:lSIZE) = real_array(1:lSIZE)
664             call move_alloc(real_tmp,real_array)
665     
666           END SUBROUTINE REAL_GROW
667     
668           SUBROUTINE REAL_GROW2(real_array,new_size)
669             IMPLICIT NONE
670     
671             INTEGER, INTENT(IN) :: new_size
672             DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: real_array
673             DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: real_tmp
674             INTEGER lSIZE, lSIZE2
675     
676             lSIZE = size(real_array,1)
677             lSIZE2 = size(real_array,2)
678             allocate(real_tmp(lSIZE,new_size))
679             real_tmp(:,1:lSIZE2) = real_array(:,1:lSIZE2)
680             call move_alloc(real_tmp,real_array)
681     
682           END SUBROUTINE REAL_GROW2
683     
684           SUBROUTINE REAL_GROW2_reverse(real_array,new_size)
685             IMPLICIT NONE
686     
687             INTEGER, INTENT(IN) :: new_size
688             DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: real_array
689             DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: real_tmp
690             INTEGER lSIZE, lSIZE2
691     
692             lSIZE = size(real_array,1)
693             lSIZE2 = size(real_array,2)
694             allocate(real_tmp(new_size,lSIZE2))
695             real_tmp(1:lSIZE,:) = real_array(1:lSIZE,:)
696             call move_alloc(real_tmp,real_array)
697     
698           END SUBROUTINE REAL_GROW2_REVERSE
699     
700           SUBROUTINE REAL_GROW3(real_array,new_size)
701             IMPLICIT NONE
702     
703             INTEGER, INTENT(IN) :: new_size
704             DOUBLE PRECISION, DIMENSION(:,:,:), ALLOCATABLE, INTENT(INOUT) :: real_array
705             DOUBLE PRECISION, DIMENSION(:,:,:), ALLOCATABLE :: real_tmp
706             INTEGER lSIZE, lSIZE2, lSIZE3
707     
708             lSIZE = size(real_array,1)
709             lSIZE2 = size(real_array,2)
710             lSIZE3 = size(real_array,3)
711             allocate(real_tmp(lSIZE,lSIZE2,new_size))
712             real_tmp(:,:,1:lSIZE3) = real_array(:,:,1:lSIZE3)
713             call move_alloc(real_tmp,real_array)
714     
715           END SUBROUTINE REAL_GROW3
716     
717           SUBROUTINE LOGICAL_GROW2_REVERSE(real_array,new_size)
718             IMPLICIT NONE
719     
720             INTEGER, INTENT(IN) :: new_size
721             LOGICAL, DIMENSION(:,:), ALLOCATABLE, INTENT(INOUT) :: real_array
722             LOGICAL, DIMENSION(:,:), ALLOCATABLE :: real_tmp
723             INTEGER lSIZE, lSIZE2
724     
725             lSIZE = size(real_array,1)
726             lSIZE2 = size(real_array,2)
727             allocate(real_tmp(new_size,lSIZE2))
728             real_tmp(1:lSIZE,:) = real_array(1:lSIZE,:)
729             call move_alloc(real_tmp,real_array)
730     
731           END SUBROUTINE LOGICAL_GROW2_REVERSE
732     
733         END MODULE DES_ALLOCATE
734