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