File: N:\mfix\model\des\write_res1_des_mod.f

1           MODULE WRITE_RES1_DES
2     
3           use desmpi
4           use compar, only: myPE
5           use compar, only: PE_IO
6     
7           use cdist, only: bDist_IO
8     
9           use mpi_comm_des, only: DES_GATHER
10           use mpi_comm_des, only: DESMPI_GATHERV
11     
12           IMPLICIT NONE
13     
14           PRIVATE
15     
16           PUBLIC :: INIT_WRITE_RES_DES
17           PUBLIC :: FINL_WRITE_RES_DES
18     
19           PUBLIC :: WRITE_RES_DES
20           PUBLIC :: WRITE_RES_pARRAY
21           PUBLIC :: WRITE_RES_cARRAY
22     
23     ! Write scalar and data WITHOUT MPI data collection.
24           INTERFACE WRITE_RES_DES
25              MODULE PROCEDURE WRITE_RES_DES_0I, WRITE_RES_DES_1I
26              MODULE PROCEDURE WRITE_RES_DES_0D, WRITE_RES_DES_1D
27              MODULE PROCEDURE WRITE_RES_DES_0L, WRITE_RES_DES_1L
28           END INTERFACE
29     
30     ! Write particle array data.
31           INTERFACE WRITE_RES_pARRAY
32              MODULE PROCEDURE WRITE_RES_pARRAY_1B
33              MODULE PROCEDURE WRITE_RES_pARRAY_1I
34              MODULE PROCEDURE WRITE_RES_pARRAY_1D
35              MODULE PROCEDURE WRITE_RES_pARRAY_1L
36           END INTERFACE
37     
38     ! Write neighbor/collision array data.
39           INTERFACE WRITE_RES_cARRAY
40              MODULE PROCEDURE WRITE_RES_cARRAY_1I
41              MODULE PROCEDURE WRITE_RES_cARRAY_1D
42              MODULE PROCEDURE WRITE_RES_cARRAY_1L
43           END INTERFACE
44     
45           INTEGER, PARAMETER :: RDES_UNIT = 901
46     
47     ! Send/Recv parameters for Particle arrays:
48           INTEGER :: pROOTCNT, pPROCCNT
49           INTEGER :: pSEND
50           INTEGER, allocatable :: pGATHER(:)
51           INTEGER, allocatable :: pDISPLS(:)
52     
53     ! Send/Recv parameters for Collision/Neighbor arrays:
54           INTEGER :: cROOTCNT, cPROCCNT
55           INTEGER :: cSEND
56           INTEGER, allocatable :: cGATHER(:)
57           INTEGER, allocatable :: cDISPLS(:)
58     
59           CONTAINS
60     
61     !``````````````````````````````````````````````````````````````````````!
62     ! Subroutine: OPEN_RES_DES                                             !
63     !                                                                      !
64     ! Purpose: Construct the file name and open the DES RES file.          !
65     !``````````````````````````````````````````````````````````````````````!
66           SUBROUTINE OPEN_RES_DES(BASE)
67     
68           use machine, only: OPEN_N1
69     
70           CHARACTER(len=*), INTENT(IN)  :: BASE
71           CHARACTER(len=32) :: lFNAME
72     
73           IF(bDIST_IO) THEN
74              WRITE(lFNAME,'(A,I4.4,A)') BASE//'_DES_',myPE,'.RES'
75              OPEN(CONVERT='BIG_ENDIAN',UNIT=RDES_UNIT, FILE=lFNAME, FORM='UNFORMATTED',         &
76                 STATUS='UNKNOWN', ACCESS='DIRECT', RECL=OPEN_N1)
77     
78           ELSEIF(myPE == PE_IO) THEN
79              WRITE(lFNAME,'(A,A)') BASE//'_DES.RES'
80              OPEN(CONVERT='BIG_ENDIAN',UNIT=RDES_UNIT, FILE=lFNAME, FORM='UNFORMATTED',         &
81                 STATUS='UNKNOWN', ACCESS='DIRECT', RECL=OPEN_N1)
82           ENDIF
83     
84           END SUBROUTINE OPEN_RES_DES
85     
86     !``````````````````````````````````````````````````````````````````````!
87     ! Subroutine: INIT_WRITE_RES_DES                                       !
88     !                                                                      !
89     ! Purpose: Construct the file name and open the DES RES file.          !
90     !``````````````````````````````````````````````````````````````````````!
91           SUBROUTINE INIT_WRITE_RES_DES(BASE, lVERSION, lNEXT_REC)
92     
93           use compar, only: numPEs
94           use mpi_utility, only: GLOBAL_SUM
95     
96           use discretelement, only: PIP, iGHOST_CNT
97           use discretelement, only: NEIGHBORS, NEIGHBOR_INDEX, NEIGH_NUM
98           use functions, only: is_nonexistent
99     
100           CHARACTER(len=*), INTENT(IN)  :: BASE
101           DOUBLE PRECISION, INTENT(IN) :: lVERSION
102           INTEGER, INTENT(OUT) :: lNEXT_REC
103     
104     ! Number of real particles on local rank
105           INTEGER :: lPROC
106     ! Total number of real particles.
107           INTEGER :: lGHOST_CNT
108     ! Local gather counts for send/recv
109           INTEGER :: lGatherCnts(0:NUMPEs-1)
110     ! Loop counters
111           INTEGER :: LC1,part
112     
113           CALL OPEN_RES_DES(BASE)
114     
115           allocate(pGATHER(0:numPEs-1))
116           allocate(pDISPLS(0:numPEs-1))
117     
118           allocate(cGATHER(0:numPEs-1))
119           allocate(cDISPLS(0:numPEs-1))
120     
121           IF(bDIST_IO) THEN
122     
123              pROOTCNT = PIP
124              pPROCCNT = pROOTCNT
125     
126              lGHOST_CNT = iGHOST_CNT
127     
128              cROOTCNT = NEIGH_NUM
129              cPROCCNT = cROOTCNT
130           ELSE
131     
132     ! Setup data for particle array data collection:
133              pROOTCNT = 10
134              pPROCCNT = PIP - iGHOST_CNT
135     
136     ! Rank 0 gets the total number of gloabl particles.
137              CALL GLOBAL_SUM(pPROCCNT, pROOTCNT)
138     
139     ! Serial IO does not store ghost particle data.
140              lGHOST_CNT = 0
141     
142     ! Construct an array for the Root process that states the number of
143     ! (real) particles on each process.
144              pSEND = pPROCCNT
145     
146              lGatherCnts = 0
147              lGatherCnts(myPE) = pPROCCNT
148     
149              CALL GLOBAL_SUM(lGatherCnts, pGATHER)
150     
151     ! Calculate the displacements for each process in the global array.
152              pDISPLS(0) = 0
153              DO lPROC = 1,NUMPES-1
154                 pDISPLS(lPROC) = pDISPLS(lPROC-1) + pGATHER(lPROC-1)
155              ENDDO
156     
157     ! Setup data for neighbor arrays
158              cROOTCNT = 10
159     ! Count the number of real neighbors.
160              cPROCCNT = 0
161              part = 1
162              DO LC1 = 1, NEIGH_NUM
163                 IF (0 .eq. NEIGHBORS(LC1)) EXIT
164                 IF (LC1.eq.NEIGHBOR_INDEX(part)) THEN
165                    part = part + 1
166                 ENDIF
167                 IF(.NOT.IS_NONEXISTENT(part) .AND. .NOT.IS_NONEXISTENT(NEIGHBORS(LC1))) THEN
168                    cPROCCNT = cPROCCNT +1
169                 ENDIF
170     
171              ENDDO
172     
173     ! Rank 0 gets the total number of global particles.
174              CALL GLOBAL_SUM(cPROCCNT, cROOTCNT)
175     
176     ! Construct an array for the Root process that states the number of
177     ! (real) particles on each process.
178              cSEND = cPROCCNT
179     
180              lGatherCnts = 0
181              lGatherCnts(myPE) = cPROCCNT
182     
183              CALL GLOBAL_SUM(lGatherCnts, cGATHER)
184     
185     ! Calculate the displacements for each process in the global array.
186              cDISPLS(0) = 0
187              DO lPROC = 1,NUMPES-1
188                 cDISPLS(lPROC) = cDISPLS(lPROC-1) + cGATHER(lPROC-1)
189              ENDDO
190     
191           ENDIF
192     
193     ! Write out the initial data.
194           lNEXT_REC = 1
195           CALL WRITE_RES_DES(lNEXT_REC, lVERSION)    ! RES file version
196           CALL WRITE_RES_DES(lNEXT_REC, pROOTCNT)    ! Number of Particles
197           CALL WRITE_RES_DES(lNEXT_REC, lGHOST_CNT)  ! Number of Ghosts
198           CALL WRITE_RES_DES(lNEXT_REC, cROOTCNT)    ! Number of neighbors
199     
200           RETURN
201           END SUBROUTINE INIT_WRITE_RES_DES
202     
203     !``````````````````````````````````````````````````````````````````````!
204     ! Subroutine: CLOSE_RES_DES                                            !
205     !                                                                      !
206     ! Purpose: Close the DES RES file.                                     !
207     !``````````````````````````````````````````````````````````````````````!
208           SUBROUTINE FINL_WRITE_RES_DES
209     
210           IF(bDIST_IO .OR. myPE == PE_IO) close(RDES_UNIT)
211     
212           IF(allocated(dPROCBUF)) deallocate(dPROCBUF)
213           IF(allocated(dROOTBUF)) deallocate(dROOTBUF)
214           IF(allocated(iPROCBUF)) deallocate(iPROCBUF)
215           IF(allocated(iROOTBUF)) deallocate(iROOTBUF)
216     
217           if(allocated(pGATHER)) deallocate(pGATHER)
218           if(allocated(pDISPLS)) deallocate(pDISPLS)
219     
220           if(allocated(cGATHER)) deallocate(cGATHER)
221           if(allocated(cDISPLS)) deallocate(cDISPLS)
222     
223           RETURN
224           END SUBROUTINE FINL_WRITE_RES_DES
225     
226     !``````````````````````````````````````````````````````````````````````!
227     ! Subroutine: WRITE_RES_DES_0I                                         !
228     !                                                                      !
229     ! Purpose: Write scalar integers to RES file.                          !
230     !``````````````````````````````````````````````````````````````````````!
231           SUBROUTINE WRITE_RES_DES_0I(lNEXT_REC, INPUT_I)
232     
233           INTEGER, INTENT(INOUT) :: lNEXT_REC
234           INTEGER, INTENT(IN) :: INPUT_I
235     
236           IF(bDIST_IO .OR. myPE == PE_IO) &
237              WRITE(RDES_UNIT, REC=lNEXT_REC) INPUT_I
238     
239           lNEXT_REC = lNEXT_REC + 1
240     
241           RETURN
242           END SUBROUTINE WRITE_RES_DES_0I
243     
244     !``````````````````````````````````````````````````````````````````````!
245     ! Subroutine: WRITE_RES_DES_1I                                         !
246     !                                                                      !
247     ! Purpose: Write an array of integers to RES file. Note that data is   !
248     ! not collected and hsould be on rank 0.                               !
249     !``````````````````````````````````````````````````````````````````````!
250           SUBROUTINE WRITE_RES_DES_1I(lNEXT_REC, INPUT_I)
251     
252           INTEGER, INTENT(INOUT) :: lNEXT_REC
253           INTEGER, INTENT(IN) :: INPUT_I(:)
254     
255           INTEGER :: lSIZE
256     
257           lSIZE = size(INPUT_I)
258     
259           IF(bDIST_IO .OR. myPE == PE_IO) &
260              CALL OUT_BIN_512i(RDES_UNIT, INPUT_I, lSIZE, lNEXT_REC)
261     
262           RETURN
263           END SUBROUTINE WRITE_RES_DES_1I
264     
265     !``````````````````````````````````````````````````````````````````````!
266     ! Subroutine: WRITE_RES_DES_0D                                         !
267     !                                                                      !
268     ! Purpose: Write scalar double percision values to RES file.           !
269     !``````````````````````````````````````````````````````````````````````!
270           SUBROUTINE WRITE_RES_DES_0D(lNEXT_REC, INPUT_D)
271     
272           INTEGER, INTENT(INOUT) :: lNEXT_REC
273           DOUBLE PRECISION, INTENT(IN) :: INPUT_D
274     
275           IF(bDIST_IO .OR. myPE == PE_IO) &
276              WRITE(RDES_UNIT, REC=lNEXT_REC) INPUT_D
277     
278           lNEXT_REC = lNEXT_REC + 1
279     
280           RETURN
281           END SUBROUTINE WRITE_RES_DES_0D
282     
283     !``````````````````````````````````````````````````````````````````````!
284     ! Subroutine: WRITE_RES_DES_1D                                         !
285     !                                                                      !
286     ! Purpose: Write an array of double percision values to RES file. Note !
287     ! that data is not collected and should be on rank 0.                  !
288     !``````````````````````````````````````````````````````````````````````!
289           SUBROUTINE WRITE_RES_DES_1D(lNEXT_REC, INPUT_D)
290     
291           INTEGER, INTENT(INOUT) :: lNEXT_REC
292           DOUBLE PRECISION, INTENT(IN) :: INPUT_D(:)
293     
294           INTEGER :: lSIZE
295     
296           lSIZE = size(INPUT_D)
297     
298           IF(bDIST_IO .OR. myPE == PE_IO) &
299              CALL OUT_BIN_512(RDES_UNIT, INPUT_D, lSIZE, lNEXT_REC)
300     
301           RETURN
302           END SUBROUTINE WRITE_RES_DES_1D
303     
304     !``````````````````````````````````````````````````````````````````````!
305     ! Subroutine: WRITE_RES_DES_0L                                         !
306     !                                                                      !
307     ! Purpose: Write scalar logical values to RES file.                    !
308     !``````````````````````````````````````````````````````````````````````!
309           SUBROUTINE WRITE_RES_DES_0L(lNEXT_REC, INPUT_L)
310     
311           INTEGER, INTENT(INOUT) :: lNEXT_REC
312           LOGICAL, INTENT(IN) :: INPUT_L
313     
314           INTEGER :: INPUT_I
315     
316           INPUT_I = merge(1,0,INPUT_L)
317     
318           IF(bDIST_IO .OR. myPE == PE_IO) &
319              WRITE(RDES_UNIT, REC=lNEXT_REC) INPUT_I
320     
321           lNEXT_REC = lNEXT_REC + 1
322     
323           RETURN
324           END SUBROUTINE WRITE_RES_DES_0L
325     
326     !``````````````````````````````````````````````````````````````````````!
327     ! Subroutine: WRITE_RES_DES_1L                                         !
328     !                                                                      !
329     ! Purpose: Write an array of integers to RES file. Note that data is   !
330     ! not collected and hsould be on rank 0.                               !
331     !``````````````````````````````````````````````````````````````````````!
332           SUBROUTINE WRITE_RES_DES_1L(lNEXT_REC, INPUT_L)
333     
334           INTEGER, INTENT(INOUT) :: lNEXT_REC
335           LOGICAL, INTENT(IN) :: INPUT_L(:)
336     
337           INTEGER, ALLOCATABLE :: INPUT_I(:)
338     
339           INTEGER :: lSIZE, LC1
340     
341           lSIZE = size(INPUT_L)
342           ALLOCATE(INPUT_I(lSIZE))
343     
344           DO LC1=1, lSIZE
345              INPUT_I(LC1) = merge(1,0,INPUT_L(LC1))
346           ENDDO
347     
348           IF(bDIST_IO .OR. myPE == PE_IO) &
349              CALL OUT_BIN_512i(RDES_UNIT, INPUT_I, lSIZE, lNEXT_REC)
350     
351           IF(allocated(INPUT_I)) deallocate(INPUT_I)
352     
353           RETURN
354           END SUBROUTINE WRITE_RES_DES_1L
355     
356     !``````````````````````````````````````````````````````````````````````!
357     ! Subroutine: WRITE_RES_PARRAY_1B                                      !
358     !                                                                      !
359     ! Purpose: Write scalar bytes to RES file.                             !
360     !``````````````````````````````````````````````````````````````````````!
361           SUBROUTINE WRITE_RES_PARRAY_1B(lNEXT_REC, INPUT_B, pLOC2GLB)
362     
363           use desmpi, only: iProcBuf
364           use discretelement, only: MAX_PIP, PIP
365           use discretelement, only: iGLOBAL_ID
366           use functions, only: is_nonexistent
367     
368           INTEGER, INTENT(INOUT) :: lNEXT_REC
369           INTEGER(KIND=1), INTENT(IN) :: INPUT_B(:)
370           LOGICAL, INTENT(IN), OPTIONAL :: pLOC2GLB
371     
372           INTEGER, ALLOCATABLE :: INPUT_I(:)
373           LOGICAL :: lLOC2GLB
374     ! Loop counters
375           INTEGER :: LC1, LC2
376     
377           lLOC2GLB = .FALSE.
378           IF(present(pLOC2GLB)) lLOC2GLB = pLOC2GLB
379     
380           allocate(iPROCBUF(pPROCCNT))
381           allocate(iROOTBUF(pROOTCNT))
382     
383           allocate(input_i(size(input_b)))
384     
385           input_i(:) = input_b(:)
386     
387           iDISPLS = pDISPLS
388           iGath_SendCnt = pSEND
389           iGatherCnts   = pGATHER
390     
391           IF(bDIST_IO) THEN
392              LC1 = 1
393              IF(lLOC2GLB) THEN
394                 DO LC2 = 1, MAX_PIP
395                    IF(LC1 > PIP) EXIT
396                    IF(IS_NONEXISTENT(LC1)) CYCLE
397                    iProcBuf(LC1) = iGLOBAL_ID(INPUT_I(LC2))
398                    LC1 = LC1 + 1
399                 ENDDO
400              ELSE
401                 DO LC2 = 1, MAX_PIP
402                    IF(LC1 > PIP) EXIT
403                    IF(IS_NONEXISTENT(LC1)) CYCLE
404                    iProcBuf(LC1) = INPUT_I(LC2)
405                    LC1 = LC1 + 1
406                 ENDDO
407              ENDIF
408              CALL OUT_BIN_512i(RDES_UNIT, iProcBuf, pROOTCNT, lNEXT_REC)
409     
410           ELSE
411              CALL DES_GATHER(INPUT_I, lLOC2GLB)
412              IF(myPE == PE_IO) &
413                 CALL OUT_BIN_512i(RDES_UNIT,iROOTBUF, pROOTCNT, lNEXT_REC)
414           ENDIF
415     
416           deallocate(iPROCBUF)
417           deallocate(iROOTBUF)
418           deallocate(input_i)
419     
420           RETURN
421           END SUBROUTINE WRITE_RES_PARRAY_1B
422     
423     !``````````````````````````````````````````````````````````````````````!
424     ! Subroutine: WRITE_RES_PARRAY_1I                                      !
425     !                                                                      !
426     ! Purpose: Write scalar integers to RES file.                          !
427     !``````````````````````````````````````````````````````````````````````!
428           SUBROUTINE WRITE_RES_PARRAY_1I(lNEXT_REC, INPUT_I, pLOC2GLB)
429     
430           use desmpi, only: iProcBuf
431           use discretelement, only: MAX_PIP, PIP
432           use discretelement, only: iGLOBAL_ID
433           use functions, only: is_nonexistent
434     
435           INTEGER, INTENT(INOUT) :: lNEXT_REC
436           INTEGER, INTENT(IN) :: INPUT_I(:)
437           LOGICAL, INTENT(IN), OPTIONAL :: pLOC2GLB
438     
439           LOGICAL :: lLOC2GLB
440     ! Loop counters
441           INTEGER :: LC1, LC2
442     
443           lLOC2GLB = .FALSE.
444           IF(present(pLOC2GLB)) lLOC2GLB = pLOC2GLB
445     
446           allocate(iPROCBUF(pPROCCNT))
447           allocate(iROOTBUF(pROOTCNT))
448     
449           iDISPLS = pDISPLS
450           iGath_SendCnt = pSEND
451           iGatherCnts   = pGATHER
452     
453           IF(bDIST_IO) THEN
454              LC1 = 1
455     
456              IF(lLOC2GLB) THEN
457                 DO LC2 = 1, MAX_PIP
458                    IF(LC1 > PIP) EXIT
459                    IF(IS_NONEXISTENT(LC1)) CYCLE
460                    iProcBuf(LC1) = iGLOBAL_ID(INPUT_I(LC2))
461                    LC1 = LC1 + 1
462                 ENDDO
463              ELSE
464                 DO LC2 = 1, MAX_PIP
465                    IF(LC1 > PIP) EXIT
466                    IF(IS_NONEXISTENT(LC1)) CYCLE
467                    iProcBuf(LC1) = INPUT_I(LC2)
468                    LC1 = LC1 + 1
469                 ENDDO
470              ENDIF
471              CALL OUT_BIN_512i(RDES_UNIT, iProcBuf, pROOTCNT, lNEXT_REC)
472     
473           ELSE
474              CALL DES_GATHER(INPUT_I, lLOC2GLB)
475              IF(myPE == PE_IO) &
476                 CALL OUT_BIN_512i(RDES_UNIT,iROOTBUF, pROOTCNT, lNEXT_REC)
477           ENDIF
478     
479           deallocate(iPROCBUF)
480           deallocate(iROOTBUF)
481     
482           RETURN
483           END SUBROUTINE WRITE_RES_PARRAY_1I
484     
485     !``````````````````````````````````````````````````````````````````````!
486     ! Subroutine: WRITE_RES_PARRAY_1D                                      !
487     !                                                                      !
488     ! Purpose: Write scalar integers to RES file.                          !
489     !``````````````````````````````````````````````````````````````````````!
490           SUBROUTINE WRITE_RES_PARRAY_1D(lNEXT_REC, INPUT_D)
491     
492           use discretelement, only: MAX_PIP, PIP
493           use functions, only: is_nonexistent
494     
495           IMPLICIT NONE
496     
497           INTEGER, INTENT(INOUT) :: lNEXT_REC
498           DOUBLE PRECISION, INTENT(IN) :: INPUT_D(:)
499     
500     ! Loop counters
501           INTEGER :: LC1, LC2
502     
503     
504           allocate(dPROCBUF(pPROCCNT))
505           allocate(dROOTBUF(pROOTCNT))
506     
507           iDISPLS = pDISPLS
508           iGath_SendCnt = pSEND
509           iGatherCnts   = pGATHER
510     
511           IF(bDIST_IO) THEN
512              LC1 = 1
513              DO LC2 = 1, MAX_PIP
514                 IF(LC1 > PIP) EXIT
515                 IF(IS_NONEXISTENT(LC1)) CYCLE
516                 dProcBuf(LC1) = INPUT_D(LC2)
517                 LC1 = LC1 + 1
518              ENDDO
519              CALL OUT_BIN_512(RDES_UNIT, dProcBuf, pROOTCNT, lNEXT_REC)
520           ELSE
521              CALL DES_GATHER(INPUT_D)
522              IF(myPE == PE_IO) &
523                 CALL OUT_BIN_512(RDES_UNIT,dRootBuf, pROOTCNT, lNEXT_REC)
524           ENDIF
525     
526           deallocate(dPROCBUF)
527           deallocate(dROOTBUF)
528     
529           RETURN
530           END SUBROUTINE WRITE_RES_PARRAY_1D
531     
532     !``````````````````````````````````````````````````````````````````````!
533     ! Subroutine: WRITE_RES_PARRAY_1D                                      !
534     !                                                                      !
535     ! Purpose: Write scalar integers to RES file.                          !
536     !``````````````````````````````````````````````````````````````````````!
537           SUBROUTINE WRITE_RES_PARRAY_1L(lNEXT_REC, INPUT_L)
538     
539           use desmpi, only: iProcBuf
540           use discretelement, only: MAX_PIP, PIP
541           use functions, only: is_nonexistent
542     
543           INTEGER, INTENT(INOUT) :: lNEXT_REC
544           LOGICAL, INTENT(IN) :: INPUT_L(:)
545     
546     ! Loop counters
547           INTEGER :: LC1, LC2
548     
549           allocate(iPROCBUF(pPROCCNT))
550           allocate(iROOTBUF(pROOTCNT))
551     
552           iDISPLS = pDISPLS
553           iGath_SendCnt = pSEND
554           iGatherCnts   = pGATHER
555     
556           IF(bDIST_IO) THEN
557              LC1 = 1
558              DO LC2 = 1, MAX_PIP
559                 IF(LC1 > PIP) EXIT
560                 IF(IS_NONEXISTENT(LC1)) CYCLE
561                 iProcBuf(LC1) = merge(1,0,INPUT_L(LC2))
562                 LC1 = LC1 + 1
563              ENDDO
564              CALL OUT_BIN_512i(RDES_UNIT, iProcBuf, pROOTCNT, lNEXT_REC)
565           ELSE
566              CALL DES_GATHER(INPUT_L)
567              IF(myPE == PE_IO) &
568                 CALL OUT_BIN_512i(RDES_UNIT,iRootBuf, pROOTCNT, lNEXT_REC)
569           ENDIF
570     
571           deallocate(iPROCBUF)
572           deallocate(iROOTBUF)
573     
574           RETURN
575           END SUBROUTINE WRITE_RES_PARRAY_1L
576     
577     !``````````````````````````````````````````````````````````````````````!
578     ! Subroutine: WRITE_RES_cARRAY_1I                                      !
579     !                                                                      !
580     ! Purpose: Write scalar integers to RES file.                          !
581     !``````````````````````````````````````````````````````````````````````!
582           SUBROUTINE WRITE_RES_cARRAY_1I(lNEXT_REC, INPUT_I, pLOC2GLB)
583     
584           use desmpi, only: iProcBuf
585           use discretelement, only: NEIGHBORS, NEIGHBOR_INDEX, NEIGH_NUM
586           USE functions, only: is_nonexistent
587           use discretelement, only: iGlobal_ID
588     
589           INTEGER, INTENT(INOUT) :: lNEXT_REC
590           INTEGER, INTENT(IN) :: INPUT_I(:)
591           LOGICAL, INTENT(IN), OPTIONAL :: pLOC2GLB
592     
593           LOGICAL :: lLOC2GLB
594     ! Loop counters
595           INTEGER :: LC1, LC2, part
596     
597           lLOC2GLB = .FALSE.
598           IF(present(pLOC2GLB)) lLOC2GLB = pLOC2GLB
599     
600           allocate(iPROCBUF(cPROCCNT))
601           allocate(iROOTBUF(cROOTCNT))
602     
603           iDISPLS = cDISPLS
604           iGath_SendCnt = cSEND
605           iGatherCnts   = cGATHER
606     
607           LC2 = 1
608           part = 1
609     
610           DO LC1 = 1, NEIGH_NUM
611              IF (0 .eq. NEIGHBORS(LC1)) EXIT
612              IF (LC1.eq.NEIGHBOR_INDEX(part)) THEN
613                 part = part + 1
614              ENDIF
615              IF(.NOT.IS_NONEXISTENT(part) .AND. .NOT.IS_NONEXISTENT(NEIGHBORS(LC1))) THEN
616                 IF(lLOC2GLB) THEN
617                    iProcBuf(LC2) = iGLOBAL_ID(INPUT_I(LC1))
618                 ELSE
619                    iProcBuf(LC2) = INPUT_I(LC1)
620                 ENDIF
621                 LC2 = LC2 + 1
622              ENDIF
623           ENDDO
624     
625           IF(bDIST_IO) THEN
626              CALL OUT_BIN_512i(RDES_UNIT, iProcBuf, cPROCCNT, lNEXT_REC)
627     
628           ELSE
629              CALL DESMPI_GATHERV(pTYPE=1)
630              IF(myPE == PE_IO) &
631                 CALL OUT_BIN_512i(RDES_UNIT,iROOTBUF, cROOTCNT, lNEXT_REC)
632           ENDIF
633     
634           deallocate(iPROCBUF)
635           deallocate(iROOTBUF)
636     
637           RETURN
638           END SUBROUTINE WRITE_RES_cARRAY_1I
639     
640     !``````````````````````````````````````````````````````````````````````!
641     ! Subroutine: WRITE_RES_cARRAY_1D                                      !
642     !                                                                      !
643     ! Purpose: Write scalar integers to RES file.                          !
644     !``````````````````````````````````````````````````````````````````````!
645           SUBROUTINE WRITE_RES_cARRAY_1D(lNEXT_REC, INPUT_D)
646     
647           use desmpi, only: dPROCBUF ! Local process buffer
648           use desmpi, only: dROOTBUF ! Root process buffer
649           use discretelement, only: NEIGHBORS, NEIGHBOR_INDEX, NEIGH_NUM
650           USE functions, only: is_nonexistent
651     
652           INTEGER, INTENT(INOUT) :: lNEXT_REC
653           DOUBLE PRECISION, INTENT(IN) :: INPUT_D(:)
654     
655     ! Loop counters
656           INTEGER :: LC1, LC2, part
657     
658           allocate(dPROCBUF(cPROCCNT))
659           allocate(dROOTBUF(cROOTCNT))
660     
661           iDISPLS = cDISPLS
662           iGath_SendCnt = cSEND
663           iGatherCnts   = cGATHER
664     
665           LC2 = 1
666           part = 1
667           DO LC1 = 1, NEIGH_NUM
668              IF (0 .eq. NEIGHBORS(LC1)) EXIT
669              IF (LC1.eq.NEIGHBOR_INDEX(part)) THEN
670                 part = part + 1
671              ENDIF
672              IF(.NOT.IS_NONEXISTENT(part) .AND. .NOT.IS_NONEXISTENT(NEIGHBORS(LC1))) THEN
673                 dProcBuf(LC2) = INPUT_D(LC1)
674                 LC2 = LC2 + 1
675              ENDIF
676           ENDDO
677     
678           IF(bDIST_IO) THEN
679              CALL OUT_BIN_512(RDES_UNIT, dProcBuf, cPROCCNT, lNEXT_REC)
680     
681           ELSE
682              CALL DESMPI_GATHERV(pTYPE=2)
683              IF(myPE == PE_IO) &
684                 CALL OUT_BIN_512(RDES_UNIT, dROOTBUF, cROOTCNT, lNEXT_REC)
685           ENDIF
686     
687           deallocate(dPROCBUF)
688           deallocate(dROOTBUF)
689     
690           RETURN
691           END SUBROUTINE WRITE_RES_cARRAY_1D
692     
693     
694     !``````````````````````````````````````````````````````````````````````!
695     ! Subroutine: WRITE_RES_cARRAY_1L                                      !
696     !                                                                      !
697     ! Purpose: Write scalar integers to RES file.                          !
698     !``````````````````````````````````````````````````````````````````````!
699           SUBROUTINE WRITE_RES_cARRAY_1L(lNEXT_REC, INPUT_L)
700     
701           use desmpi, only: iProcBuf
702           use discretelement, only: NEIGHBORS, NEIGHBOR_INDEX, NEIGH_NUM
703           use functions, only: is_nonexistent
704     
705           INTEGER, INTENT(INOUT) :: lNEXT_REC
706           LOGICAL, INTENT(IN) :: INPUT_L(:)
707     
708     ! Loop counters
709           INTEGER :: LC1, LC2, part
710     
711           allocate(iPROCBUF(cPROCCNT))
712           allocate(iROOTBUF(cROOTCNT))
713     
714           iDISPLS = cDISPLS
715           iGath_SendCnt = cSEND
716           iGatherCnts   = cGATHER
717     
718     ! Pack the local buffer, skipping data for deleted particles.
719           LC2 = 1
720           part = 1
721           DO LC1 = 1, NEIGH_NUM
722              IF (0 .eq. NEIGHBORS(LC1)) EXIT
723              IF (LC1.eq.NEIGHBOR_INDEX(part)) THEN
724                 part = part + 1
725              ENDIF
726              IF(.NOT.IS_NONEXISTENT(part) .AND. .NOT.IS_NONEXISTENT(NEIGHBORS(LC1))) THEN
727                 iProcBuf(LC2) = merge(1,0,INPUT_L(LC1))
728                 LC2 = LC2 + 1
729              ENDIF
730           ENDDO
731     
732           IF(bDIST_IO) THEN
733              CALL OUT_BIN_512i(RDES_UNIT, iProcBuf, cPROCCNT, lNEXT_REC)
734     
735           ELSE
736              CALL DESMPI_GATHERV(pTYPE=1)
737              IF(myPE == PE_IO) &
738                 CALL OUT_BIN_512i(RDES_UNIT,iROOTBUF, cROOTCNT, lNEXT_REC)
739           ENDIF
740     
741           deallocate(iPROCBUF)
742           deallocate(iROOTBUF)
743     
744           RETURN
745           END SUBROUTINE WRITE_RES_cARRAY_1L
746     
747           END MODULE WRITE_RES1_DES
748