File: /nfs/home/0/users/jenkins/mfix.git/model/des/read_res1_des_mod.f

1           MODULE READ_RES1_DES
2     
3           use cdist, only: bDist_IO
4           use compar, only: PE_IO
5           use compar, only: myPE
6           use des_allocate
7           use desmpi
8           use error_manager
9           use mpi_comm_des, only: DESMPI_GATHERV
10           use mpi_comm_des, only: DESMPI_SCATTERV
11     
12           IMPLICIT NONE
13     
14           PRIVATE
15     
16           PUBLIC :: INIT_READ_RES_DES
17           PUBLIC :: FINL_READ_RES_DES
18     
19           PUBLIC :: READ_PAR_POS
20           PUBLIC :: READ_PAR_COL
21     
22           PUBLIC :: READ_RES_DES
23           PUBLIC :: READ_RES_pARRAY
24           PUBLIC :: READ_RES_cARRAY
25     
26           INTERFACE READ_RES_DES
27              MODULE PROCEDURE READ_RES_DES_0I
28              MODULE PROCEDURE READ_RES_DES_1I
29              MODULE PROCEDURE READ_RES_DES_0D
30              MODULE PROCEDURE READ_RES_DES_1D
31              MODULE PROCEDURE READ_RES_DES_0L
32              MODULE PROCEDURE READ_RES_DES_1L
33           END INTERFACE
34     
35           INTERFACE READ_RES_pARRAY
36              MODULE PROCEDURE READ_RES_pARRAY_1I
37              MODULE PROCEDURE READ_RES_pARRAY_1D
38              MODULE PROCEDURE READ_RES_pARRAY_1L
39           END INTERFACE
40     
41           INTERFACE READ_RES_cARRAY
42              MODULE PROCEDURE READ_RES_cARRAY_1I
43              MODULE PROCEDURE READ_RES_cARRAY_1D
44              MODULE PROCEDURE READ_RES_cARRAY_1L
45           END INTERFACE
46     
47     
48           INTEGER, PARAMETER :: RDES_UNIT = 901
49     
50           INTEGER :: pIN_COUNT
51           INTEGER :: cIN_COUNT
52     
53     ! Send/Recv parameters for Particle arrays:
54           INTEGER :: pROOTCNT, pPROCCNT
55           INTEGER :: pRECV
56           INTEGER, allocatable :: pSCATTER(:)
57           INTEGER, allocatable :: pDISPLS(:)
58     
59     ! Variables used for reading restart file
60           INTEGER, ALLOCATABLE :: pRestartMap(:)
61           INTEGER, ALLOCATABLE :: cRestartMap(:)
62     
63     ! Send/Recv parameters for Particle arrays:
64           INTEGER :: cROOTCNT, cPROCCNT
65           INTEGER :: cRECV
66           INTEGER, allocatable :: cSCATTER(:)
67           INTEGER, allocatable :: cDISPLS(:)
68     
69           INTEGER, ALLOCATABLE :: iPAR_COL(:,:)
70     
71           CONTAINS
72     
73     !``````````````````````````````````````````````````````````````````````!
74     ! Subroutine: INIT_READ_RES_DES                                        !
75     !                                                                      !
76     ! Purpose: Construct the file name and open the DES RES file.          !
77     !``````````````````````````````````````````````````````````````````````!
78           SUBROUTINE INIT_READ_RES_DES(BASE, lVERSION, lNEXT_REC)
79     
80           use discretelement, only: MAX_PIP, PIP
81           use discretelement, only: iGHOST_CNT
82           use discretelement, only: PAIR_NUM, PAIR_MAX
83     
84           use compar, only: numPEs
85           use machine, only: OPEN_N1
86           use geometry, only: NO_K
87     
88           use mpi_utility, only: BCAST
89           use mpi_utility, only: GLOBAL_ALL_SUM
90     
91           implicit none
92     
93           CHARACTER(len=*), INTENT(IN)  :: BASE
94           DOUBLE PRECISION, INTENT(OUT) :: lVERSION
95           INTEGER, INTENT(OUT) :: lNEXT_REC
96     
97           CHARACTER(len=32) :: lFNAME
98     
99     ! Integer Error Flag
100           INTEGER :: IER
101     
102     
103           allocate(pSCATTER(0:numPEs-1))
104           allocate(pDISPLS(0:numPEs-1))
105     
106           allocate(cSCATTER(0:numPEs-1))
107           allocate(cDISPLS(0:numPEs-1))
108     
109     
110           IF(bDIST_IO) THEN
111     
112              WRITE(lFNAME,'(A,I4.4,A)') BASE//'_DES_',myPE,'.RES'
113              OPEN(UNIT=RDES_UNIT, FILE=lFNAME, FORM='UNFORMATTED',         &
114                 STATUS='UNKNOWN', ACCESS='DIRECT', RECL=OPEN_N1)
115     
116              READ(RDES_UNIT, REC=1) lVERSION
117              READ(RDES_UNIT, REC=2) pIN_COUNT
118              READ(RDES_UNIT, REC=3) iGHOST_CNT
119              READ(RDES_UNIT, REC=4) cIN_COUNT
120     
121              IF(PIP > MAX_PIP) THEN
122                 write(*,*) "From des_read_restart:"
123                 write(*,*) "Error: The pip is grater than current max_pip"
124                 write(*,*) "pip=" ,pip,"; max_pip =", max_pip
125     
126              ENDIF
127     
128              PIP = pIN_COUNT
129              PAIR_NUM = cIN_COUNT
130     
131              DO WHILE(PAIR_NUM > PAIR_MAX)
132                 PAIR_MAX = 2*PAIR_MAX
133              ENDDO
134              CALL PAIR_GROW
135     
136           ELSE
137     
138              IF(myPE == PE_IO) THEN
139                 WRITE(lFNAME,'(A,A)') BASE//'_DES.RES'
140                 OPEN(UNIT=RDES_UNIT, FILE=lFNAME, FORM='UNFORMATTED',      &
141                    STATUS='UNKNOWN', ACCESS='DIRECT', RECL=OPEN_N1)
142     
143                 READ(RDES_UNIT, REC=1) pIN_COUNT
144     
145                 READ(RDES_UNIT, REC=1) lVERSION
146                 READ(RDES_UNIT, REC=2) pIN_COUNT
147     !           READ(RDES_UNIT, REC=3) -NOTHING-
148                 READ(RDES_UNIT, REC=4) cIN_COUNT
149     
150              ELSE
151                 pIN_COUNT = 10
152              ENDIF
153     
154              IER = 0
155     
156     ! Allocate the particle restart map. This is used in determining were
157     ! particle data is sent. Only process zero needs this array.
158              allocate( pRestartMap(pIN_COUNT), STAT=IER)
159              IF(IER/=0) THEN
160                 WRITE(ERR_MSG, 1200) 'pRestartMap', trim(iVAL(pIN_COUNT))
161                 CALL FLUSH_ERR_MSG
162              ENDIF
163     
164              CALL BCAST(lVERSION, PE_IO)
165     
166     ! Allocate the collision restart map array. All ranks allocatet this
167     ! array so that mapping the collision data can be done in parallel.
168              CALL BCAST(cIN_COUNT, PE_IO)
169              allocate( cRestartMap(cIN_COUNT), STAT=IER)
170              IF(IER/=0) THEN
171                 WRITE(ERR_MSG, 1200) 'cRestartMap', trim(iVAL(cIN_COUNT))
172                 CALL FLUSH_ERR_MSG
173              ENDIF
174     
175      1200 FORMAT('Error 1200: Unable to allocate sufficient memory to ',&
176              'read in DES',/'restart file. size(',A,') = ',A)
177     
178              CALL GLOBAL_ALL_SUM(IER)
179              IF(IER/=0) CALL MFIX_EXIT(myPE)
180     
181           ENDIF
182     
183           lNEXT_REC = 5
184     
185           RETURN
186           END SUBROUTINE INIT_READ_RES_DES
187     
188     
189     !``````````````````````````````````````````````````````````````````````!
190     ! Subroutine: CLOSE_RES_DES                                            !
191     !                                                                      !
192     ! Purpose: Close the DES RES file.                                     !
193     !``````````````````````````````````````````````````````````````````````!
194           SUBROUTINE FINL_READ_RES_DES
195     
196     
197           IF(bDIST_IO .OR. myPE == PE_IO) close(RDES_UNIT)
198     
199           IF(allocated(dPROCBUF)) deallocate(dPROCBUF)
200           IF(allocated(dROOTBUF)) deallocate(dROOTBUF)
201           IF(allocated(iPROCBUF)) deallocate(iPROCBUF)
202           IF(allocated(iROOTBUF)) deallocate(iROOTBUF)
203     
204           IF(allocated(pRestartMap)) deallocate(pRestartMap)
205           IF(allocated(cRestartMap)) deallocate(cRestartMap)
206     
207           IF(allocated(pSCATTER)) deallocate(pSCATTER)
208           IF(allocated(pDISPLS)) deallocate(pDISPLS)
209     
210           IF(allocated(cSCATTER)) deallocate(cSCATTER)
211           IF(allocated(cDISPLS)) deallocate(cDISPLS)
212     
213     
214     
215     
216           RETURN
217           END SUBROUTINE FINL_READ_RES_DES
218     
219     
220     
221     !``````````````````````````````````````````````````````````````````````!
222     ! Subroutine: READ_PAR_POS                                             !
223     !                                                                      !
224     ! Purpose: Generates the mapping used by the scatter routines to send  !
225     ! read data to the correct rank.                                       !
226     !``````````````````````````````````````````````````````````````````````!
227           SUBROUTINE READ_PAR_POS(lNEXT_REC)
228     
229           use discretelement, only: PIP
230           use discretelement, only: DES_POS_NEW
231           use geometry, only: NO_K
232           use compar, only: numPEs
233     
234           use mpi_utility, only: GLOBAL_SUM
235           USE in_binary_512
236     
237           implicit none
238     
239           INTEGER, INTENT(INOUT) :: lNEXT_REC
240     
241           INTEGER :: lDIMN
242           INTEGER :: LC1, lPROC
243           INTEGER :: lScatterCNTS(0:NUMPEs-1)
244     ! The number of particles on each process.
245           INTEGER :: PAR_CNT(0:NUMPEs-1)
246     
247     !-----------------------------------------------
248     
249           CALL INIT_ERR_MSG("READ_PAR_POS")
250     
251           lDIMN = merge(2,3,NO_K)
252     
253     
254     ! All process read positions for distributed IO restarts.
255           IF(bDIST_IO) THEN
256              DO LC1 = 1, lDIMN
257                 CALL READ_RES_DES(lNEXT_REC, DES_POS_NEW(LC1,:))
258              ENDDO
259              RETURN
260           ENDIF
261     
262           allocate( dPAR_POS(lDIMN, pIN_COUNT))
263     
264     ! Only the IO proccess reads positions.
265           IF(myPE == PE_IO) THEN
266              DO LC1=1, merge(2,3,NO_K)
267                 CALL IN_BIN_512(RDES_UNIT, dPAR_POS(LC1,:),                &
268                    pIN_COUNT, lNEXT_REC)
269              ENDDO
270           ENDIF
271     
272     ! Use the particle postions and the domain coverage of each process
273     ! to determine which processor each particle belongs.
274           CALL MAP_pARRAY_TO_PROC(PAR_CNT)
275     
276     ! Send the particle position data to the individual ranks.
277           CALL SCATTER_PAR_POS(PAR_CNT)
278     
279     ! Set up the read/scatter arrary information.
280           pPROCCNT = PIP
281           pROOTCNT = pIN_COUNT
282     
283     ! Set the recv count for this process.
284           pRECV = PIP
285     
286     ! Construct an array for the Root process that states the number of
287     ! (real) particles on each process.
288           lScatterCnts(:) = 0; lScatterCnts(mype) = PIP
289           CALL GLOBAL_SUM(lScatterCnts,pSCATTER)
290     
291     ! Calculate the displacements for each process in the global array.
292           pDispls(0) = 0
293           DO lPROC = 1, NUMPEs-1
294              pDispls(lPROC) = pDispls(lPROC-1) + pSCATTER(lPROC-1)
295           ENDDO
296     
297           IF(allocated(dPAR_POS)) deallocate(dPAR_POS)
298     
299           CALL FINL_ERR_MSG
300     
301           RETURN
302           END SUBROUTINE READ_PAR_POS
303     
304     
305     !``````````````````````````````````````````````````````````````````````!
306     ! Subroutine: MAP_pARRAY_TO_PROC                                       !
307     !                                                                      !
308     ! Purpose: Use the particle positions to determine which processor     !
309     ! they live on and count the number of particles on each process.      !
310     !``````````````````````````````````````````````````````````````````````!
311           SUBROUTINE MAP_pARRAY_TO_PROC(lPAR_CNT)
312     
313           use discretelement, only: PIP, MAX_PIP
314           use discretelement, only: XE, YN, ZT
315           use geometry, only: IMIN1, IMAX1
316           use geometry, only: JMIN1, JMAX1
317           use geometry, only: KMIN1, KMAX1
318           use geometry, only: NO_K, DO_K
319           use compar, only: numPEs
320           use compar, only: ISTART1_ALL, IEND1_ALL
321           use compar, only: JSTART1_ALL, JEND1_ALL
322           use compar, only: KSTART1_ALL, KEND1_ALL
323     
324           use mpi_utility, only: BCAST
325           use mpi_utility, only: GLOBAL_ALL_SUM
326     
327           implicit none
328     
329           INTEGER, INTENT(OUT) :: lPAR_CNT(0:numPEs-1)
330     
331     ! Data dimensionality flag.
332           INTEGER :: lDIMN
333     ! Loop counters.
334           INTEGER :: LC1, lPROC
335     ! Error flag.
336           INTEGER :: IER(0:numPEs-1)
337     ! Local scatter count.
338           INTEGER :: lScatterCNTS(0:NUMPEs-1)
339     ! The X/Y/Z bounds of the physical space "owned" by each process.
340           DOUBLE PRECISION :: lxmin(0:NUMPEs-1), lxmax(0:NUMPEs-1)
341           DOUBLE PRECISION :: lymin(0:NUMPEs-1), lymax(0:NUMPEs-1)
342           DOUBLE PRECISION :: lzmin(0:NUMPEs-1), lzmax(0:NUMPEs-1)
343     !-----------------------------------------------
344     
345           CALL INIT_ERR_MSG("MAP_pARRAY_TO_PROC")
346     
347     ! Initialize the error flag.
348           IER = 0
349     
350           lDIMN = merge(2, 3, NO_K)
351     
352     ! set the domain range for each processor
353           DO lPROC= 0, NUMPEs-1
354              lxmin(lproc) = xe(istart1_all(lproc)-1)
355              lxmax(lproc) = xe(iend1_all(lproc))
356              lymin(lproc) = yn(jstart1_all(lproc)-1)
357              lymax(lproc) = yn(jend1_all(lproc))
358              lzmin(lproc) = zt(kstart1_all(lproc)-1)
359              lzmax(lproc) = zt(kend1_all(lproc))
360     
361     ! modify the range for mass inlet and outlet, as particles injected
362     ! can lie outside the domain and not ghost particles
363              IF(istart1_all(lproc).eq.imin1) &
364                 lxmin(lproc) = xe(istart1_all(lproc)-2)
365              IF(iend1_all(lproc).eq.imax1) &
366                 lxmax(lproc) = xe(iend1_all(lproc)+1)
367              IF(jstart1_all(lproc).eq.jmin1) &
368                 lymin(lproc) = yn(jstart1_all(lproc)-2)
369              IF(jend1_all(lproc).eq.jmax1)  &
370                 lymax(lproc) = yn(jend1_all(lproc)+1)
371              IF(kstart1_all(lproc).eq.kmin1 .AND. DO_K) &
372                 lzmin(lproc) = zt(kstart1_all(lproc)-2)
373              IF(kend1_all(lproc).eq.kmax1 .AND. DO_K) &
374                 lzmax(lproc) = zt(kend1_all(lproc)+1)
375           ENDDO
376     
377     ! build the send buffer in PE_IO proc
378     ! first pass to get the count of particles
379           IER = 0
380           pRestartMap(:) = -1
381           lPAR_CNT(:) = 0
382           IF(myPE == PE_IO) THEN
383              DO LC1 = 1, pIN_COUNT
384                 DO lPROC=0, NUMPEs-1
385                    IF(dPAR_POS(1,LC1) >= lxmin(lproc) .AND. &
386                       dPAR_POS(1,LC1) <  lxmax(lproc) .AND. &
387                       dPAR_POS(2,LC1) >= lymin(lproc) .AND. &
388                       dPAR_POS(2,LC1) <  lymax(lproc)) THEN
389                       IF(NO_K)THEN
390                          lPAR_CNT(lPROC) = lPAR_CNT(lPROC) + 1
391                          pRestartMap(LC1) = lproc
392                          EXIT
393                       ELSE
394                          IF(dPAR_POS(3,LC1) >= lzmin(lproc) .AND. &
395                             dPAR_POS(3,LC1) <  lzmax(lproc)) THEN
396                             lPAR_CNT(lPROC) = lPAR_CNT(lPROC) + 1
397                             pRestartMap(LC1) = lproc
398                             EXIT
399                          ENDIF
400                       ENDIF
401                    ENDIF
402                 ENDDO  ! Loop over processes
403                 IF (pRestartMap(LC1) == -1) then
404                    IER(myPE) = -1
405                    WRITE(ERR_MSG,1000) trim(iVal(LC1))
406                    CALL FLUSH_ERR_MSG(FOOTER=.FALSE.)
407                    IF(NO_K) THEN
408                       WRITE(ERR_MSG,1001) dPAR_POS(1:2,LC1)
409                       CALL FLUSH_ERR_MSG(HEADER=.FALSE.)
410                    ELSE
411                       WRITE(ERR_MSG,1002) dPAR_POS(1:3,LC1)
412                       CALL FLUSH_ERR_MSG(HEADER=.FALSE.)
413                    ENDIF
414                 ENDIF
415              ENDDO  ! Loop over particles
416           ENDIF
417     
418      1000 FORMAT('Error 1000: Unable to locat paritcle inside domain:',/&
419              3x,'Particle Number:',A)
420      1001 FORMAT(3x,'X POS: ',g12.5,/3x,'Y POS: ',g12.5)
421      1002 FORMAT(3x,'X POS: ',g12.5,/3x,'Y POS: ',g12.5,/3x,'Z POS: ',g12.5)
422     
423     ! Send out the error flag and exit if needed.
424           CALL BCAST(IER, PE_IO)
425           IF(IER(PE_IO) /= 0) CALL MFIX_EXIT(myPE)
426     
427     ! PE_IO sends out the number of particles for each process.
428           CALL BCAST(lPAR_CNT(0:NUMPES-1), PE_IO)
429     
430     ! Each process stores the number of particles-on-its-process. The error
431     ! flag is set if that number exceeds the maximum.
432           PIP = lPAR_CNT(myPE)
433           CALL PARTICLE_GROW(PIP)
434     
435     ! Global collection of error flags to abort it the max was exceeded.
436           CALL GLOBAL_ALL_SUM(IER)
437           IF(sum(IER) /= 0) THEN
438              WRITE(ERR_MSG,1100)
439              CALL FLUSH_ERR_MSG(FOOTER=.FALSE.)
440              DO LC1=0, numPEs-1
441                 IF(IER(LC1) /= 0) THEN
442                    WRITE(ERR_MSG,"(3(2x,I10))")LC1,IER(LC1)-1,lPAR_CNT(LC1)
443                    CALL FLUSH_ERR_MSG(HEADER=.FALSE., FOOTER=.FALSE.)
444                 ENDIF
445              ENDDO
446              WRITE(ERR_MSG,"('Aborting.')")
447              CALL FLUSH_ERR_MSG(HEADER=.FALSE.,ABORT=.TRUE.)
448           ENDIF
449     
450      1100 FORMAT('Error 1100: Maximum number of particles exceeded.',2/    &
451              5x,'Process',5x,'Maximum',7x,'Count')
452     
453     
454           CALL FINL_ERR_MSG
455     
456           RETURN
457           END SUBROUTINE MAP_pARRAY_TO_PROC
458     
459     
460     
461     !``````````````````````````````````````````````````````````````````````!
462     ! Subroutine: DES_RESTART_MAP                                          !
463     !                                                                      !
464     ! Purpose: Generates the mapping used by the scatter routines to send  !
465     ! read data to the correct rank.                                       !
466     !``````````````````````````````````````````````````````````````````````!
467           SUBROUTINE SCATTER_PAR_POS(lPAR_CNT)
468     
469           use compar, only: numPEs
470     
471           use discretelement, only: DES_POS_NEW
472           use discretelement, only: PEA
473           use discretelement, only: PIP
474           use geometry, only: NO_K
475     
476           implicit none
477     
478     ! Number of particles on each process.
479           INTEGER, INTENT(INOUT) :: lPAR_CNT(0:numPEs-1)
480     ! Dimensionality flag.
481           INTEGER :: lDIMN
482     ! Loop counters.
483           INTEGER :: LC1, lPROC, lBuf, IER
484     ! Scatter counts.
485           INTEGER :: lScatterCNTS(0:NUMPEs-1)
486     
487     
488           lDIMN = merge(2,3,NO_K)
489     
490     ! Set up the recv count and allocate the local process buffer.
491           iSCR_RECVCNT = PIP*lDIMN
492           allocate (dProcBuf(iscr_recvcnt))
493     
494     ! Allocate the buffer for the root.
495           IF (myPE == PE_IO) THEN
496              allocate (dRootBuf(pIN_COUNT*lDIMN))
497           ELSE
498              allocate (dRootBuf(10))
499           ENDIF
500     
501     ! The IO processor builds drootbuffer and iDISLS
502           IF(myPE == PE_IO) THEN
503     ! Determine the offsets for each process and the amount of data that
504     ! is to be scattered to each.
505              iDISPLS(0) = 0
506              iScatterCnts(0) = lPAR_CNT(0)*lDIMN
507              DO lProc = 1, NUMPES-1
508                 iDispls(lproc) = iDispls(lproc-1) + iScatterCnts(lproc-1)
509                 iScatterCnts(lproc) = lPAR_CNT(lProc)*lDIMN
510              ENDDO
511     ! Copy the position data into the root buffer, mapped to the owner
512     ! process.
513              lPAR_CNT(:) = 0
514              DO LC1 = 1,pIN_COUNT
515                 lPROC = pRestartMap(LC1)
516                 lbuf = iDispls(lProc) + lPAR_CNT(lProc)*lDIMN+1
517                 dRootBuf(lBuf:lBuf+lDIMN-1) = dPAR_POS(1:lDIMN,LC1)
518                 lBuf = lBuf + lDIMN
519                 lPAR_CNT(lProc) = lPAR_CNT(lProc) + 1
520              ENDDO
521           ENDIF
522           CALL DESMPI_SCATTERV(pTYPE=2)
523     
524     ! Unpack the particle data.
525           DO LC1 = 1, PIP
526              lBuf = (LC1-1)*lDIMN+1
527              DES_POS_NEW(1:lDIMN,LC1) = dProcBuf(lBuf:lBuf+lDIMN-1)
528              lBuf = lBuf + lDIMN
529              PEA(LC1,1) = .TRUE.
530           ENDDO
531     
532           IF(allocated(dRootBuf)) deallocate(dRootBuf)
533           IF(allocated(dProcBuf)) deallocate(dProcBuf)
534           IF(allocated(dPAR_POS)) deallocate(dPAR_POS)
535     
536           RETURN
537           END SUBROUTINE SCATTER_PAR_POS
538     
539     
540     !``````````````````````````````````````````````````````````````````````!
541     ! Subroutine: READ_PAR_COL                                             !
542     !                                                                      !
543     ! Purpose: Generates the mapping used by the scatter routines to send  !
544     ! read data to the correct rank.                                       !
545     !``````````````````````````````````````````````````````````````````````!
546           SUBROUTINE READ_PAR_COL(lNEXT_REC)
547     
548           use discretelement, only: PAIRS, PAIR_NUM
549           use compar, only: numPEs
550     
551           use mpi_init_des, only: DES_RESTART_GHOST
552           use mpi_utility, only: BCAST
553           use mpi_utility, only: GLOBAL_SUM
554           use mpi_utility, only: GLOBAL_ALL_SUM
555           use in_binary_512i
556     
557           implicit none
558     
559           INTEGER, INTENT(INOUT) :: lNEXT_REC
560     
561           INTEGER :: LC1, lPROC
562           INTEGER :: lScatterCNTS(0:NUMPEs-1)
563     ! The number of particles on each process.
564           INTEGER :: COL_CNT(0:NUMPEs-1)
565     
566     !-----------------------------------------------
567     
568           CALL INIT_ERR_MSG("READ_PAR_COL")
569     
570     ! All process read positions for distributed IO restarts.
571           IF(bDIST_IO) THEN
572              DO LC1 = 1, 2
573                 CALL READ_RES_DES(lNEXT_REC, PAIRS(LC1,:))
574              ENDDO
575           ENDIF
576     
577           CALL DES_RESTART_GHOST
578     
579           allocate(iPAR_COL(2, cIN_COUNT))
580           iPAR_COL = 0
581     
582     ! Only the IO proccess reads positions.
583           IF(myPE == PE_IO) THEN
584              DO LC1=1, 2
585                 CALL IN_BIN_512i(RDES_UNIT, iPAR_COL(LC1,:),               &
586                    cIN_COUNT, lNEXT_REC)
587              ENDDO
588           ENDIF
589     
590     ! Broadcast collision data to all the other processes.
591            CALL GLOBAL_ALL_SUM(iPAR_COL)
592     
593     ! Determine which process owns the pair datasets. This is done either
594     ! through matching global ids or a search. The actual method depends
595     ! on the ability to allocate a large enough array.
596           CALL MAP_cARRAY_TO_PROC(COL_CNT)
597     
598     ! Send the particle position data to the individual ranks.
599           CALL GLOBAL_TO_LOC_COL
600     
601     ! Set up the read/scatter arrary information.
602           cPROCCNT = PAIR_NUM
603           cROOTCNT = cIN_COUNT
604     
605     ! Set the recv count for this process.
606           cRECV = PAIR_NUM
607     
608     ! Construct an array for the Root process that states the number of
609     ! (real) particles on each process.
610           lScatterCnts(:) = 0; lScatterCnts(mype) = PAIR_NUM
611           CALL GLOBAL_SUM(lScatterCnts,cSCATTER)
612     
613     ! Calculate the displacements for each process in the global array.
614           cDispls(0) = 0
615           DO lPROC = 1, NUMPEs-1
616              cDispls(lPROC) = cDispls(lPROC-1) + cSCATTER(lPROC-1)
617           ENDDO
618     
619           CALL FINL_ERR_MSG
620     
621           RETURN
622           END SUBROUTINE READ_PAR_COL
623     
624     
625     !``````````````````````````````````````````````````````````````````````!
626     ! Subroutine: MAP_cARRAY_TO_PROC                                       !
627     !                                                                      !
628     ! Purpose: Use the particle positions to determine which processor     !
629     ! they live on and count the number of particles on each process.      !
630     !``````````````````````````````````````````````````````````````````````!
631           SUBROUTINE MAP_cARRAY_TO_PROC(lCOL_CNT)
632     
633           use compar, only: numPEs, myPE
634           use discretelement, only: iGLOBAL_ID
635           use discretelement, only: PIP, PEA
636           use discretelement, only: PAIR_MAX, PAIR_NUM
637     
638     !      use mpi_utility, only: BCAST
639           use mpi_utility, only: GLOBAL_ALL_SUM
640           use mpi_utility, only: GLOBAL_ALL_MAX
641     
642           implicit none
643     
644           INTEGER, INTENT(OUT) :: lCOL_CNT(0:numPEs-1)
645     
646     ! Loop counters.
647           INTEGER :: LC1, LC2, lPROC, lBUF
648     ! Error flag.
649           INTEGER :: IER
650     ! Max global id.
651           INTEGER :: MAX_ID, lSTAT
652     
653           INTEGER, ALLOCATABLE :: lGLOBAL_OWNER(:)
654     
655     !-----------------------------------------------
656     
657           CALL INIT_ERR_MSG("MAP_cARRAY_TO_PROC")
658     
659     ! Initialize the error flag.
660           IER = 0
661     
662           MAX_ID = maxval(IGLOBAL_ID(1:PIP))
663           CALL GLOBAL_ALL_MAX(MAX_ID)
664     
665           allocate(lGLOBAL_OWNER(MAX_ID), STAT=lSTAT)
666           CALL GLOBAL_ALL_SUM(lSTAT)
667     
668     ! All ranks successfully allocated the array. This permits a crude
669     ! but much faster collision owner detection.
670           IF(lSTAT == 0) THEN
671     
672              WRITE(ERR_MSG,"('Matching DES pair data by global owner.')")
673              CALL FLUSH_ERR_MSG(HEADER=.FALSE.,FOOTER=.FALSE.)
674     
675              lGLOBAL_OWNER = 0
676              DO LC1=1, PIP
677                 IF(.NOT.PEA(LC1,4)) &
678                    lGLOBAL_OWNER(iGLOBAL_ID(LC1)) = myPE + 1
679              ENDDO
680     
681     ! Loop over the neighbor pair list and match the read global ID to
682     ! one of the global IDs.
683              lCOL_CNT = 0
684              cRestartMap = 0
685              DO LC1=1, cIN_COUNT
686                 IF(lGLOBAL_OWNER(iPAR_COL(1,LC1)) == myPE + 1) THEN
687                    cRestartMap(LC1) = myPE + 1
688                    lCOL_CNT(myPE) = lCOL_CNT(myPE) + 1
689                 ENDIF
690              ENDDO
691     ! One or more ranks could not allocate the memory needed to do the
692     ! quick and dirty match so do a search instead.
693           ELSE
694     
695              WRITE(ERR_MSG,"('Matching DES pair data by search.')")
696              CALL FLUSH_ERR_MSG(HEADER=.FALSE.,FOOTER=.FALSE.)
697     
698     ! Loop over the neighbor pair list and match the read global ID to
699     ! one of the global IDs.
700              lCOL_CNT = 0
701              cRestartMap = 0
702              LC1_LP: DO LC1=1, cIN_COUNT
703                 DO LC2=1, PIP!-iGHOST_CNT
704                    IF(iPAR_COL(1,LC1) == iGLOBAL_ID(LC2)) THEN
705                       cRestartMap(LC1) = myPE + 1
706                       lCOL_CNT(myPE) = lCOL_CNT(myPE) + 1
707                       CYCLE LC1_LP
708                    ENDIF
709                 ENDDO
710              ENDDO LC1_LP
711     
712           ENDIF
713     
714     ! Clean up the large array as it is no longer needed.
715           IF(allocated(lGLOBAL_OWNER)) deallocate(lGLOBAL_OWNER)
716     
717     ! Calculate the number of matched collisions over all processes. Throw
718     ! and error if it doesn't match the number of read collisions.
719           CALL GLOBAL_ALL_SUM(lCOL_CNT)
720           IF(sum(lCOL_CNT) /= cIN_COUNT) THEN
721              WRITE(ERR_MSG,1000) cIN_COUNT, sum(lCOL_CNT)
722              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
723           ENDIF
724     
725     1000 FORMAT('Error 1000: Unable to establish the own of all read ',    &
726              'collision data.',/3x,'Number of Collisions: ',I10,/3x,       &
727              'Matched Collisions:   ',I10)
728     
729     ! Sync the collision restart map arcross all ranks.
730           CALL GLOBAL_ALL_SUM(cRestartMap)
731     
732     ! Error checking and cleanup.
733           DO LC1 = 1, cIN_COUNT
734     ! Verify that each collision is owned by a rank.
735              IF (cRestartMap(LC1) == 0) THEN
736                 IER = -1
737                 WRITE(ERR_MSG,1100) trim(iVal(LC1)), trim(iVal(            &
738                    iPAR_COL(1,LC1))), trim(iVal(iPAR_COL(2,LC1)))
739                 CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
740     
741      1100 FORMAT('Error 1100: Unable to locate process pair owner:',/      &
742              3x,'Pair Number:',A,/3x,'Particles: ',A,' and ',A)
743     
744              ELSEIF(cRestartMap(LC1) > numPEs) THEN
745     
746                 IER = -1
747                 WRITE(ERR_MSG,1101) trim(iVal(LC1)), trim(iVal(            &
748                   iPAR_COL(1,LC1))), trim(iVal(iPAR_COL(2,LC1)))
749                  CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
750     
751      1101 FORMAT('Error 1101: More than one process pair owner:',/         &
752              3x,'Pair Number:',A,/3x,'Particles: ',A,' and ',A)
753     
754     ! Shift the rank ID to the correct value.
755              ELSE
756                 cRestartMap(LC1) = cRestartMap(LC1) - 1
757              ENDIF
758           ENDDO
759     
760     ! Send out the error flag and exit if needed.
761           CALL GLOBAL_ALL_SUM(IER, PE_IO)
762           IF(IER /= 0) CALL MFIX_EXIT(myPE)
763     
764     ! Each process stores the number of particles-on-its-process. The error
765     ! flag is set if that number exceeds the maximum.
766           PAIR_NUM = lCOL_CNT(myPE)
767     
768           DO WHILE(PAIR_NUM > PAIR_MAX)
769              PAIR_MAX = 2*PAIR_MAX
770           ENDDO
771           CALL PAIR_GROW
772     
773           CALL FINL_ERR_MSG
774     
775           RETURN
776           END SUBROUTINE MAP_cARRAY_TO_PROC
777     
778     
779     !``````````````````````````````````````````````````````````````````````!
780     ! Subroutine: GLOBAL_TO_LOC_COL                                        !
781     !                                                                      !
782     ! Purpose: Generates the mapping used by the scatter routines to send  !
783     ! read data to the correct rank.                                       !
784     !``````````````````````````````````````````````````````````````````````!
785           SUBROUTINE GLOBAL_TO_LOC_COL
786     
787           use compar, only: numPEs
788           use discretelement, only: PAIRS, PAIR_NUM
789           use discretelement, only: iGLOBAL_ID
790           use discretelement, only: PIP
791           use discretelement, only: PEA
792     
793           use mpi_utility, only: GLOBAL_ALL_SUM
794           use mpi_utility, only: GLOBAL_ALL_MAX
795     
796           use funits, only: DMP_LOG
797     
798           use error_manager
799     
800           implicit none
801     
802     ! Loop counters.
803           INTEGER :: LC1, LC2, LC3, IER
804           INTEGER :: UNMATCHED
805           INTEGER, ALLOCATABLE :: iLOCAL_ID(:)
806     
807     ! Max global id.
808           INTEGER :: MAX_ID, lSTAT
809     ! Debug flags.
810           LOGICAL :: dFlag
811           LOGICAL, parameter :: setDBG = .FALSE.
812     
813           CALL INIT_ERR_MSG("GLOBAL_TO_LOC_COL")
814     
815     ! Initialize the error flag.
816           IER = 0
817     
818     ! Set the local debug flag.
819           dFlag = (DMP_LOG .AND. setDBG)
820     
821           MAX_ID = maxval(IGLOBAL_ID(1:PIP))
822           CALL GLOBAL_ALL_MAX(MAX_ID)
823     
824           allocate(iLOCAL_ID(MAX_ID), STAT=lSTAT)
825           CALL GLOBAL_ALL_SUM(lSTAT)
826     
827     ! All ranks successfully allocated the array. This permits a crude
828     ! but much faster collision owner detection.
829           IF(lSTAT /= 0) THEN
830              WRITE(ERR_MSG,1000)
831              CALL FLUSH_ERR_MSG(ABORT=.TRUE.)
832           ENDIF
833     
834      1000 FORMAT('Error 1000: Unable to allocate sufficient memory to ',&
835              'generate the',/'map from global to local particle IDs.')
836     
837           iLOCAL_ID = 0
838           DO LC1=1, PIP
839              iLOCAL_ID(iGLOBAL_ID(LC1)) = LC1
840           ENDDO
841     
842     ! Store the particle data.
843           LC3 = 1
844           LC2 = 0
845           UNMATCHED = 0
846           LP1: DO LC1 = 1, cIN_COUNT
847              IF(cRestartMap(LC1) == myPE) THEN
848                 LC2 = LC2 + 1
849                 PAIRS(1,LC2) = iLOCAL_ID(iPAR_COL(1,LC1))
850                 PAIRS(2,LC2) = iLOCAL_ID(iPAR_COL(2,LC1))
851     ! Verify that the local indices are valid. If they do not match it is
852     ! likely because one of the pair was removed via an outlet at the time
853     ! the RES file was written but the ghost data wasn't updated.
854                 IF(PAIRS(1,LC2) == 0 .OR. PAIRS(2,LC2) == 0) THEN
855                    UNMATCHED = UNMATCHED + 1
856                    IF(dFLAG) THEN
857                       WRITE(ERR_MSG,1100) iPAR_COL(1,LC1), PAIRS(1,LC2),   &
858                          iPAR_COL(2,LC1), PAIRS(2,LC2)
859                       CALL FLUSH_ERR_MSG(ABORT=.FALSE.)
860                    ENDIF
861                    DO WHILE(PEA(LC3,1))
862                       LC3 = LC3 + 1
863                    ENDDO
864                    PAIRS(2,LC2) = LC3
865                 ENDIF
866              ENDIF
867           ENDDO LP1
868     
869      1100 FORMAT('Error 1100: Particle pair local indices are invalid.',/  &
870              5x,'Global-ID    Local-ID',/' 1:  ',2(3x,I9),/' 2:  ',2(3x,I9))
871     
872           CALL GLOBAL_ALL_SUM(UNMATCHED)
873           IF(UNMATCHED /= 0) THEN
874              WRITE(ERR_MSG,1101) trim(iVal(UNMATCHED))
875              CALL FLUSH_ERR_MSG
876           ENDIF
877     
878      1101 FORMAT(' Warning: 1101: ',A,' particle pair datasets were ',&
879              'not matched',/' during restart.')
880     
881           IF(allocated(iLOCAL_ID)) deallocate(iLOCAL_ID)
882     
883           CALL FINL_ERR_MSG
884     
885           RETURN
886           END SUBROUTINE GLOBAL_TO_LOC_COL
887     
888     
889     
890     !``````````````````````````````````````````````````````````````````````!
891     ! Subroutine: READ_RES_DES_0I                                          !
892     !                                                                      !
893     ! Purpose: Write scalar integers to RES file.                          !
894     !``````````````````````````````````````````````````````````````````````!
895           SUBROUTINE READ_RES_DES_0I(lNEXT_REC, INPUT_I)
896     
897           use mpi_utility, only: BCAST
898     
899           IMPLICIT NONE
900     
901           INTEGER, INTENT(INOUT) :: lNEXT_REC
902     !      INTEGER, INTENT(OUT) :: INPUT_I
903           INTEGER :: INPUT_I
904     
905           IF(bDIST_IO) THEN
906              READ(RDES_UNIT, REC=lNEXT_REC) INPUT_I
907           ELSE
908              IF(myPE == PE_IO) READ(RDES_UNIT, REC=lNEXT_REC) INPUT_I
909              CALL BCAST(INPUT_I, PE_IO)
910           ENDIF
911     
912           lNEXT_REC = lNEXT_REC + 1
913     
914           RETURN
915           END SUBROUTINE READ_RES_DES_0I
916     
917     
918     !``````````````````````````````````````````````````````````````````````!
919     ! Subroutine: READ_RES_1I                                              !
920     !                                                                      !
921     ! Purpose: Write scalar integers to RES file.                          !
922     !``````````````````````````````````````````````````````````````````````!
923           SUBROUTINE READ_RES_DES_1I(lNEXT_REC, INPUT_I)
924     
925           use mpi_utility, only: BCAST
926           USE in_binary_512i
927     
928           IMPLICIT NONE
929     
930           INTEGER, INTENT(INOUT) :: lNEXT_REC
931           INTEGER, INTENT(OUT) :: INPUT_I(:)
932     
933           INTEGER :: lSIZE
934     
935           lSIZE = size(INPUT_I)
936     
937           IF(bDIST_IO) THEN
938              CALL IN_BIN_512i(RDES_UNIT, INPUT_I, lSIZE, lNEXT_REC)
939           ELSE
940              IF(myPE == PE_IO) &
941                 CALL IN_BIN_512i(RDES_UNIT, INPUT_I, lSIZE, lNEXT_REC)
942              CALL BCAST(INPUT_I, PE_IO)
943           ENDIF
944     
945     
946           RETURN
947           END SUBROUTINE READ_RES_DES_1I
948     
949     
950     !``````````````````````````````````````````````````````````````````````!
951     ! Subroutine: READ_RES_DES_0D                                          !
952     !                                                                      !
953     ! Purpose: Write scalar double percision values to RES file.           !
954     !``````````````````````````````````````````````````````````````````````!
955           SUBROUTINE READ_RES_DES_0D(lNEXT_REC, INPUT_D)
956     
957           use mpi_utility, only: BCAST
958     
959           INTEGER, INTENT(INOUT) :: lNEXT_REC
960           DOUBLE PRECISION, INTENT(OUT) :: INPUT_D
961     
962           IF(bDIST_IO) THEN
963              READ(RDES_UNIT, REC=lNEXT_REC) INPUT_D
964           ELSE
965              IF(myPE == PE_IO) READ(RDES_UNIT, REC=lNEXT_REC) INPUT_D
966              CALL BCAST(INPUT_D, PE_IO)
967           ENDIF
968           lNEXT_REC = lNEXT_REC + 1
969     
970           RETURN
971           END SUBROUTINE READ_RES_DES_0D
972     
973     
974     !``````````````````````````````````````````````````````````````````````!
975     ! Subroutine: READ_RES_DES_1D                                          !
976     !                                                                      !
977     ! Purpose: Write scalar integers to RES file.                          !
978     !``````````````````````````````````````````````````````````````````````!
979           SUBROUTINE READ_RES_DES_1D(lNEXT_REC, INPUT_D)
980     
981           use mpi_utility, only: BCAST
982           USE in_binary_512
983     
984           IMPLICIT NONE
985     
986           INTEGER, INTENT(INOUT) :: lNEXT_REC
987           DOUBLE PRECISION, INTENT(OUT) :: INPUT_D(:)
988     
989           INTEGER :: lSIZE
990     
991           lSIZE = size(INPUT_D)
992     
993           IF(bDIST_IO) THEN
994              CALL IN_BIN_512(RDES_UNIT, INPUT_D, lSIZE, lNEXT_REC)
995           ELSE
996              IF(myPE == PE_IO) &
997                 CALL IN_BIN_512(RDES_UNIT, INPUT_D, lSIZE, lNEXT_REC)
998              CALL BCAST(INPUT_D, PE_IO)
999           ENDIF
1000     
1001     
1002           RETURN
1003           END SUBROUTINE READ_RES_DES_1D
1004     
1005     
1006     !``````````````````````````````````````````````````````````````````````!
1007     ! Subroutine: READ_RES_DES_0L                                          !
1008     !                                                                      !
1009     ! Purpose: Write scalar logical values to RES file.                    !
1010     !``````````````````````````````````````````````````````````````````````!
1011           SUBROUTINE READ_RES_DES_0L(lNEXT_REC, OUTPUT_L)
1012     
1013           use mpi_utility, only: BCAST
1014     
1015           INTEGER, INTENT(INOUT) :: lNEXT_REC
1016           LOGICAL, INTENT(OUT) :: OUTPUT_L
1017     
1018           INTEGER :: OUTPUT_I
1019     
1020           OUTPUT_L = .TRUE.
1021     
1022           IF(bDIST_IO)THEN
1023              READ(RDES_UNIT, REC=lNEXT_REC) OUTPUT_I
1024           ELSE
1025              IF(myPE == PE_IO) READ(RDES_UNIT, REC=lNEXT_REC) OUTPUT_I
1026              CALL BCAST(OUTPUT_I, PE_IO)
1027           ENDIF
1028     
1029           IF(OUTPUT_I == 1) OUTPUT_L = .TRUE.
1030           lNEXT_REC = lNEXT_REC + 1
1031     
1032           RETURN
1033           END SUBROUTINE READ_RES_DES_0L
1034     
1035     
1036     !``````````````````````````````````````````````````````````````````````!
1037     ! Subroutine: READ_RES_DES_1L                                          !
1038     !                                                                      !
1039     ! Purpose: Write scalar integers to RES file.                          !
1040     !``````````````````````````````````````````````````````````````````````!
1041           SUBROUTINE READ_RES_DES_1L(lNEXT_REC, INPUT_L)
1042     
1043           use mpi_utility, only: BCAST
1044           USE in_binary_512i
1045     
1046           IMPLICIT NONE
1047     
1048           INTEGER, INTENT(INOUT) :: lNEXT_REC
1049           LOGICAL, INTENT(OUT) :: INPUT_L(:)
1050     
1051           INTEGER, ALLOCATABLE :: INPUT_I(:)
1052     
1053           INTEGER :: lSIZE, LC1
1054     
1055           lSIZE = size(INPUT_I)
1056           ALLOCATE( INPUT_I(lSIZE))
1057     
1058           IF(bDIST_IO) THEN
1059              CALL IN_BIN_512i(RDES_UNIT, INPUT_I, lSIZE, lNEXT_REC)
1060           ELSE
1061              IF(myPE == PE_IO) &
1062                 CALL IN_BIN_512i(RDES_UNIT, INPUT_I, lSIZE, lNEXT_REC)
1063              CALL BCAST(INPUT_I, PE_IO)
1064           ENDIF
1065     
1066           DO LC1=1, LSIZE
1067              IF(INPUT_I(LC1) == 1) THEN
1068                 INPUT_L(LC1) = .TRUE.
1069              ELSE
1070                 INPUT_L(LC1) = .FALSE.
1071              ENDIF
1072           ENDDO
1073     
1074           IF(allocated(INPUT_I)) deallocate(INPUT_I)
1075     
1076           RETURN
1077           END SUBROUTINE READ_RES_DES_1L
1078     
1079     
1080     !``````````````````````````````````````````````````````````````````````!
1081     ! Subroutine: READ_RES_DES_1I                                          !
1082     !                                                                      !
1083     ! Purpose: Write scalar integers to RES file.                          !
1084     !``````````````````````````````````````````````````````````````````````!
1085           SUBROUTINE READ_RES_pARRAY_1I(lNEXT_REC, OUTPUT_I)
1086     
1087           use discretelement, only: PIP
1088     
1089           use desmpi, only: iRootBuf
1090           use desmpi, only: iProcBuf
1091     
1092           use compar, only: numPEs
1093           USE in_binary_512i
1094     
1095           IMPLICIT NONE
1096     
1097           INTEGER, INTENT(INOUT) :: lNEXT_REC
1098           INTEGER, INTENT(OUT) :: OUTPUT_I(:)
1099     
1100           LOGICAL :: lLOC2GLB
1101     ! Loop counters
1102           INTEGER :: LC1, lc2
1103     
1104           INTEGER :: lPROC
1105     
1106           INTEGER, ALLOCATABLE :: lBUF_I(:)
1107           INTEGER, ALLOCATABLE :: lCOUNT(:)
1108     
1109     
1110           allocate(iPROCBUF(pPROCCNT))
1111           allocate(iROOTBUF(pROOTCNT))
1112     
1113           iDISPLS = pDISPLS
1114           iScr_RecvCNT = pRECV
1115           iScatterCNTS = pSCATTER
1116     
1117           IF(bDIST_IO) THEN
1118              CALL IN_BIN_512i(RDES_UNIT, OUTPUT_I, pIN_COUNT, lNEXT_REC)
1119           ELSE
1120     
1121              IF(myPE == PE_IO) THEN
1122                 allocate(lBUF_I(pIN_COUNT))
1123                 allocate(lCOUNT(0:NUMPEs-1))
1124     
1125                 CALL IN_BIN_512i(RDES_UNIT, lBUF_I, pIN_COUNT, lNEXT_REC)
1126     
1127                 lCOUNT = 0
1128                 DO LC1=1, pIN_COUNT
1129                    lPROC = pRestartMap(LC1)
1130                    lCOUNT(lPROC) = lCOUNT(lPROC) + 1
1131                    iRootBuf(iDispls(lPROC) + lCOUNT(lPROC)) = lBUF_I(LC1)
1132                 ENDDO
1133     
1134                 deallocate(lBUF_I)
1135                 deallocate(lCOUNT)
1136              ENDIF
1137              CALL DESMPI_SCATTERV(ptype=1)
1138              DO LC1=1, PIP
1139                 OUTPUT_I(LC1) = iProcBuf(LC1)
1140              ENDDO
1141     
1142           ENDIF
1143     
1144           deallocate(iPROCBUF)
1145           deallocate(iROOTBUF)
1146     
1147           RETURN
1148           END SUBROUTINE READ_RES_pARRAY_1I
1149     
1150     
1151     
1152     !``````````````````````````````````````````````````````````````````````!
1153     ! Subroutine: READ_RES_pARRAY_1D                                       !
1154     !                                                                      !
1155     ! Purpose: Write scalar integers to RES file.                          !
1156     !``````````````````````````````````````````````````````````````````````!
1157           SUBROUTINE READ_RES_pARRAY_1D(lNEXT_REC, OUTPUT_D)
1158     
1159           use discretelement, only: PIP
1160           use desmpi, only: dRootBuf
1161           use desmpi, only: dProcBuf
1162           use compar, only: numPEs
1163           USE in_binary_512
1164     
1165           IMPLICIT NONE
1166     
1167           INTEGER, INTENT(INOUT) :: lNEXT_REC
1168           DOUBLE PRECISION, INTENT(OUT) :: OUTPUT_D(:)
1169     
1170     ! Loop counters
1171           INTEGER :: LC1
1172     
1173           INTEGER :: lPROC
1174     
1175           DOUBLE PRECISION, ALLOCATABLE :: lBUF_D(:)
1176           INTEGER, ALLOCATABLE :: lCOUNT(:)
1177     
1178     
1179           allocate(dPROCBUF(pPROCCNT))
1180           allocate(dROOTBUF(pROOTCNT))
1181     
1182           iDISPLS = pDISPLS
1183           iScr_RecvCNT = pRECV
1184           iScatterCNTS = pSCATTER
1185     
1186           IF(bDIST_IO) THEN
1187              CALL IN_BIN_512(RDES_UNIT, OUTPUT_D, pIN_COUNT, lNEXT_REC)
1188           ELSE
1189              IF(myPE == PE_IO) THEN
1190                 allocate(lBUF_D(pIN_COUNT))
1191                 allocate(lCOUNT(0:NUMPEs-1))
1192     
1193                 CALL IN_BIN_512(RDES_UNIT, lBUF_D, pIN_COUNT, lNEXT_REC)
1194     
1195                 lCOUNT = 0
1196                 DO LC1=1, pIN_COUNT
1197                    lPROC = pRestartMap(LC1)
1198                    lCOUNT(lPROC) = lCOUNT(lPROC) + 1
1199                    dRootBuf(iDispls(lPROC) + lCOUNT(lPROC)) = lBUF_D(LC1)
1200                 ENDDO
1201     
1202                 deallocate(lBUF_D)
1203                 deallocate(lCOUNT)
1204              ENDIF
1205              CALL DESMPI_SCATTERV(ptype=2)
1206              DO LC1=1, PIP
1207                 OUTPUT_D(LC1) = dProcBuf(LC1)
1208              ENDDO
1209           ENDIF
1210     
1211           deallocate(dPROCBUF)
1212           deallocate(dROOTBUF)
1213     
1214           RETURN
1215           END SUBROUTINE READ_RES_pARRAY_1D
1216     
1217     
1218     !``````````````````````````````````````````````````````````````````````!
1219     ! Subroutine: READ_RES_pARRAY_1L                                       !
1220     !                                                                      !
1221     ! Purpose: Write scalar integers to RES file.                          !
1222     !``````````````````````````````````````````````````````````````````````!
1223           SUBROUTINE READ_RES_pARRAY_1L(lNEXT_REC, OUTPUT_L)
1224     
1225           use discretelement, only: PIP
1226           use desmpi, only: iRootBuf
1227           use desmpi, only: iProcBuf
1228           use compar, only: numPEs
1229           USE in_binary_512i
1230     
1231           IMPLICIT NONE
1232     
1233           INTEGER, INTENT(INOUT) :: lNEXT_REC
1234           LOGICAL, INTENT(OUT) :: OUTPUT_L(:)
1235     
1236           LOGICAL :: lLOC2GLB
1237     ! Loop counters
1238           INTEGER :: LC1
1239     
1240           INTEGER :: lPROC
1241     
1242           INTEGER, ALLOCATABLE :: lBUF_I(:)
1243           INTEGER, ALLOCATABLE :: lCOUNT(:)
1244     
1245           allocate(iPROCBUF(pPROCCNT))
1246           allocate(iROOTBUF(pROOTCNT))
1247     
1248           iDISPLS = pDISPLS
1249           iScr_RecvCNT = pRECV
1250           iScatterCNTS = pSCATTER
1251     
1252           IF(bDIST_IO) THEN
1253              allocate(lBUF_I(pIN_COUNT))
1254              CALL IN_BIN_512i(RDES_UNIT, lBUF_I, pIN_COUNT, lNEXT_REC)
1255              DO LC1=1,pIN_COUNT
1256                 IF(lBUF_I(LC1) == 1) THEN
1257                    OUTPUT_L(LC1) = .TRUE.
1258                 ELSE
1259                    OUTPUT_L(LC1) = .FALSE.
1260                 ENDIF
1261              ENDDO
1262              deallocate(lBUF_I)
1263           ELSE
1264              IF(myPE == PE_IO) THEN
1265                 allocate(lBUF_I(pIN_COUNT))
1266                 allocate(lCOUNT(0:NUMPEs-1))
1267     
1268                 CALL IN_BIN_512i(RDES_UNIT, lBUF_I, pIN_COUNT, lNEXT_REC)
1269     
1270                 lCOUNT = 0
1271                 DO LC1=1, pIN_COUNT
1272                    lPROC = pRestartMap(LC1)
1273                    lCOUNT(lPROC) = lCOUNT(lPROC) + 1
1274                    iRootBuf(iDispls(lPROC) + lCOUNT(lPROC)) = lBUF_I(LC1)
1275                 ENDDO
1276     
1277                 deallocate(lBUF_I)
1278                 deallocate(lCOUNT)
1279              ENDIF
1280              CALL DESMPI_SCATTERV(ptype=1)
1281              DO LC1=1, PIP
1282                 IF(iProcBuf(LC1) == 1) THEN
1283                    OUTPUT_L(LC1) = .TRUE.
1284                 ELSE
1285                    OUTPUT_L(LC1) = .FALSE.
1286                 ENDIF
1287              ENDDO
1288           ENDIF
1289     
1290           deallocate(iPROCBUF)
1291           deallocate(iROOTBUF)
1292     
1293           RETURN
1294           END SUBROUTINE READ_RES_pARRAY_1L
1295     
1296     
1297     !``````````````````````````````````````````````````````````````````````!
1298     ! Subroutine: READ_RES_DES_1I                                          !
1299     !                                                                      !
1300     ! Purpose: Write scalar integers to RES file.                          !
1301     !``````````````````````````````````````````````````````````````````````!
1302           SUBROUTINE READ_RES_cARRAY_1I(lNEXT_REC, OUTPUT_I)
1303     
1304           use desmpi, only: iRootBuf
1305           use desmpi, only: iProcBuf
1306           use compar, only: numPEs
1307           use discretelement, only: PAIR_NUM
1308           USE in_binary_512i
1309     
1310           IMPLICIT NONE
1311     
1312           INTEGER, INTENT(INOUT) :: lNEXT_REC
1313           INTEGER, INTENT(OUT) :: OUTPUT_I(:)
1314     
1315     ! Loop counters
1316           INTEGER :: LC1
1317     
1318           INTEGER :: lPROC
1319     
1320           INTEGER, ALLOCATABLE :: lBUF_I(:)
1321           INTEGER, ALLOCATABLE :: lCOUNT(:)
1322     
1323     
1324           allocate(iPROCBUF(cPROCCNT))
1325           allocate(iROOTBUF(cROOTCNT))
1326     
1327           iDISPLS = cDISPLS
1328           iScr_RecvCNT = cRECV
1329           iScatterCNTS = cSCATTER
1330     
1331           IF(bDIST_IO) THEN
1332              CALL IN_BIN_512i(RDES_UNIT, OUTPUT_I, cIN_COUNT, lNEXT_REC)
1333           ELSE
1334              IF(myPE == PE_IO) THEN
1335                 allocate(lBUF_I(cIN_COUNT))
1336                 allocate(lCOUNT(0:NUMPEs-1))
1337     
1338                 CALL IN_BIN_512i(RDES_UNIT, lBUF_I, cIN_COUNT, lNEXT_REC)
1339     
1340                 lCOUNT = 0
1341                 DO LC1=1, cIN_COUNT
1342                    lPROC = cRestartMap(LC1)
1343                    lCOUNT(lPROC) = lCOUNT(lPROC) + 1
1344                    iRootBuf(iDispls(lPROC) + lCOUNT(lPROC)) = lBUF_I(LC1)
1345                 ENDDO
1346     
1347                 deallocate(lBUF_I)
1348                 deallocate(lCOUNT)
1349              ENDIF
1350              CALL DESMPI_SCATTERV(ptype=1)
1351              DO LC1=1, PAIR_NUM
1352                 OUTPUT_I(LC1) = iProcBuf(LC1)
1353              ENDDO
1354           ENDIF
1355     
1356           deallocate(iPROCBUF)
1357           deallocate(iROOTBUF)
1358     
1359           RETURN
1360           END SUBROUTINE READ_RES_cARRAY_1I
1361     
1362     
1363     !``````````````````````````````````````````````````````````````````````!
1364     ! Subroutine: READ_RES_cARRAY_1D                                       !
1365     !                                                                      !
1366     ! Purpose: Write scalar integers to RES file.                          !
1367     !``````````````````````````````````````````````````````````````````````!
1368           SUBROUTINE READ_RES_cARRAY_1D(lNEXT_REC, OUTPUT_D)
1369     
1370           use compar, only: numPEs
1371           use discretelement, only: PAIR_NUM
1372           use desmpi, only: dRootBuf
1373           use desmpi, only: dProcBuf
1374           USE in_binary_512
1375     
1376           IMPLICIT NONE
1377     
1378           INTEGER, INTENT(INOUT) :: lNEXT_REC
1379           DOUBLE PRECISION, INTENT(OUT) :: OUTPUT_D(:)
1380     
1381     ! Loop counters
1382           INTEGER :: LC1
1383     
1384           INTEGER :: lPROC
1385     
1386           DOUBLE PRECISION, ALLOCATABLE :: lBUF_D(:)
1387           INTEGER, ALLOCATABLE :: lCOUNT(:)
1388     
1389     
1390           allocate(dPROCBUF(cPROCCNT))
1391           allocate(dROOTBUF(cROOTCNT))
1392     
1393           iDISPLS = cDISPLS
1394           iScr_RecvCNT = cRECV
1395           iScatterCNTS = cSCATTER
1396     
1397     
1398           IF(bDIST_IO) THEN
1399              CALL IN_BIN_512(RDES_UNIT, OUTPUT_D, cIN_COUNT, lNEXT_REC)
1400           ELSE
1401              IF(myPE == PE_IO) THEN
1402                 allocate(lBUF_D(cIN_COUNT))
1403                 allocate(lCOUNT(0:NUMPEs-1))
1404     
1405                 CALL IN_BIN_512(RDES_UNIT, lBUF_D, cIN_COUNT, lNEXT_REC)
1406     
1407                 lCOUNT = 0
1408                 DO LC1=1, cIN_COUNT
1409                    lPROC = cRestartMap(LC1)
1410                    lCOUNT(lPROC) = lCOUNT(lPROC) + 1
1411                    dRootBuf(iDispls(lPROC) + lCOUNT(lPROC)) = lBUF_D(LC1)
1412                 ENDDO
1413     
1414                 deallocate(lBUF_D)
1415                 deallocate(lCOUNT)
1416              ENDIF
1417              CALL DESMPI_SCATTERV(ptype=2)
1418              DO LC1=1, PAIR_NUM
1419                 OUTPUT_D(LC1) = dProcBuf(LC1)
1420              ENDDO
1421           ENDIF
1422     
1423           deallocate(dPROCBUF)
1424           deallocate(dROOTBUF)
1425     
1426           RETURN
1427           END SUBROUTINE READ_RES_cARRAY_1D
1428     
1429     
1430     !``````````````````````````````````````````````````````````````````````!
1431     ! Subroutine: READ_RES_pARRAY_1L                                       !
1432     !                                                                      !
1433     ! Purpose: Write scalar integers to RES file.                          !
1434     !``````````````````````````````````````````````````````````````````````!
1435           SUBROUTINE READ_RES_cARRAY_1L(lNEXT_REC, OUTPUT_L)
1436     
1437           use compar, only: numPEs
1438           use discretelement, only: PAIR_NUM
1439           use desmpi, only: iRootBuf
1440           use desmpi, only: iProcBuf
1441           USE in_binary_512i
1442     
1443           IMPLICIT NONE
1444     
1445           INTEGER, INTENT(INOUT) :: lNEXT_REC
1446           LOGICAL, INTENT(OUT) :: OUTPUT_L(:)
1447     
1448           LOGICAL :: lLOC2GLB
1449     ! Loop counters
1450           INTEGER :: LC1
1451     
1452           INTEGER :: lPROC
1453     
1454           INTEGER, ALLOCATABLE :: lBUF_I(:)
1455           INTEGER, ALLOCATABLE :: lCOUNT(:)
1456     
1457           allocate(iPROCBUF(cPROCCNT))
1458           allocate(iROOTBUF(cROOTCNT))
1459     
1460           iDISPLS = cDISPLS
1461           iScr_RecvCNT = cRECV
1462           iScatterCNTS = cSCATTER
1463     
1464           IF(bDIST_IO) THEN
1465              allocate(lBUF_I(cIN_COUNT))
1466              CALL IN_BIN_512i(RDES_UNIT, lBUF_I, cIN_COUNT, lNEXT_REC)
1467              DO LC1=1,cIN_COUNT
1468                 IF(lBUF_I(LC1) == 1) THEN
1469                    OUTPUT_L(LC1) = .TRUE.
1470                 ELSE
1471                    OUTPUT_L(LC1) = .FALSE.
1472                 ENDIF
1473              ENDDO
1474              deallocate(lBUF_I)
1475           ELSE
1476              IF(myPE == PE_IO) THEN
1477                 allocate(lBUF_I(cIN_COUNT))
1478                 allocate(lCOUNT(0:NUMPEs-1))
1479     
1480                 CALL IN_BIN_512i(RDES_UNIT, lBUF_I, cIN_COUNT, lNEXT_REC)
1481     
1482                 lCOUNT = 0
1483                 DO LC1=1, cIN_COUNT
1484                    lPROC = cRestartMap(LC1)
1485                    lCOUNT(lPROC) = lCOUNT(lPROC) + 1
1486                    iRootBuf(iDispls(lPROC) + lCOUNT(lPROC)) = lBUF_I(LC1)
1487                 ENDDO
1488     
1489                 deallocate(lBUF_I)
1490                 deallocate(lCOUNT)
1491              ENDIF
1492              CALL DESMPI_SCATTERV(ptype=1)
1493              DO LC1=1, PAIR_NUM
1494                 IF(iProcBuf(LC1) == 1) THEN
1495                    OUTPUT_L(LC1) = .TRUE.
1496                 ELSE
1497                    OUTPUT_L(LC1) = .FALSE.
1498                 ENDIF
1499              ENDDO
1500           ENDIF
1501     
1502           deallocate(iPROCBUF)
1503           deallocate(iROOTBUF)
1504     
1505           RETURN
1506           END SUBROUTINE READ_RES_cARRAY_1L
1507     
1508     
1509           END MODULE READ_RES1_DES
1510