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