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