File: N:\mfix\model\dbg_mod.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     !  Module: DBG                                                         !
4     !  Author: J.Musser                                   Date: 11-Nov-13  !
5     !                                                                      !
6     !  Purpose: Contains routines used for extracting the Am matrix, Bm    !
7     !  vector in .csv format.  Additionally, routines for extracting array !
8     !  subsets are included.                                               !
9     !                                                                      !
10     !  Comments:                                                           !
11     !                                                                      !
12     !  > initExtract must be invoked prior to invoking any other routine   !
13     !       in this file. The initialization routine needs to set several  !
14     !       values required by other functions.                            !
15     !                                                                      !
16     !                                                                      !
17     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
18           MODULE dbg
19     
20           use exit, only: mfix_exit
21           use mpi_utility
22           use param1, only: UNDEFINED_C
23           use param1, only: UNDEFINED_I
24     
25           IMPLICIT NONE
26     
27           PRIVATE
28     
29     ! Subroutine Access:
30     !---------------------------------------------------------------------//
31           PUBLIC ::  initExtract, matrixExtract, arrayExtract
32     
33     ! Interface
34     !---------------------------------------------------------------------//
35           interface arrayExtract
36              module procedure arrayExtract_int
37              module procedure arrayExtract_dbl
38              module procedure arrayExtract_log
39           end interface
40     
41     ! Module variables (private):
42     !---------------------------------------------------------------------//
43     
44     ! Coeffs for the space filling curve.
45           INTEGER ::  Am_c0, Am_c1, Am_c2
46     
47     ! The extent to extract the Am/Bm matrix/vector.
48           INTEGER :: iUB, iLB
49           INTEGER :: jUB, jLB
50           INTEGER :: kUB, kLB
51     
52     ! Logical flag to inflate A_M matrix to from sparce to full NxN.
53           LOGICAL, parameter :: inflateAM = .FALSE.
54     ! Debug Mode:
55           LOGICAL, parameter :: dbgMode = .FALSE.
56     
57     ! Logical to verify that the initialization routine was called.
58           LOGICAL :: initNotCalled = .TRUE.
59     
60     ! Data storage buffer size.
61           INTEGER :: dbgDIMN
62     
63     ! Buffers for I/J/K indices
64           INTEGER, allocatable :: i_ofBuff(:)
65           INTEGER, allocatable :: j_ofBuff(:)
66           INTEGER, allocatable :: k_ofBuff(:)
67           INTEGER, allocatable :: ijk_Buff(:)
68     
69     ! Buffers for various data types.
70     ! Buffer for generic integer array.
71           INTEGER, allocatable :: outBuff_i(:)
72     ! Buffer for generic double precision array.
73           DOUBLE PRECISION, allocatable :: outBuff_dp(:)
74     ! Buffer for Am matrix.
75           DOUBLE PRECISION, allocatable :: outAm(:)
76     
77     ! Flag stating that the dgb routines are in-use
78           LOGICAL :: dbgLock
79     
80     ! File Unit
81           INTEGER, parameter :: dbgUnit = 9659
82     
83     ! File structure information: (optional)
84     ! Pass count if dumping at multiple locations or during iterations.
85           INTEGER :: fPass
86     ! Flag to append to previous output.
87           LOGICAL :: fApnd
88     ! Flag to write IJK values with output.
89           LOGICAL :: fwIJK, pwIJK
90     
91           contains
92     
93     !``````````````````````````````````````````````````````````````````````!
94     !                                                                      !
95     !                                                                      !
96     !``````````````````````````````````````````````````````````````````````!
97           INTEGER FUNCTION  dbg_funijk(xxi, xxj, xxk)
98           INTEGER, intent(in) :: xxi, xxj, xxk
99     
100           dbg_funijk = xxj + Am_c0 + xxi*Am_c1 + xxk*Am_c2
101     
102           END FUNCTION dbg_funijk
103     
104     !----------------------------------------------------------------------!
105     !                                                                      !
106     !                                                                      !
107     !                                                                      !
108     !----------------------------------------------------------------------!
109           SUBROUTINE initExtract(iLow, iHgh, jLow, jHgh, kLow, kHgh)
110     
111           use compar
112           use funits, only: DMP_LOG
113           use geometry, only: iMin3, iMax3
114           use geometry, only: jMin3, jMax3
115           use geometry, only: kMin3, kMax3, do_K
116     
117           INTEGER, optional, intent(in) :: iLow, iHgh
118           INTEGER, optional, intent(in) :: jLow, jHgh
119           INTEGER, optional, intent(in) :: kLow, kHgh
120     
121           INTEGER :: PROC
122     
123     ! Lock the dbg routines.
124           dbgLock = .TRUE.
125     
126           iLB = merge(iLow, iMin3, present(iLow))
127           iUB = merge(iHgh, iMax3, present(iHgh))
128     
129           jLB = merge(jLow, jMin3, present(jLow))
130           jUB = merge(jHgh, jMax3, present(jHgh))
131     
132           kLB = merge(kLow, kMin3, present(kLow))
133           kUB = merge(kHgh, kMax3, present(kHgh))
134     
135     ! Some basic logical checks:
136           IF(iLB > iUB .OR. iLB<iMin3 .OR. iUB>iMax3) THEN
137              IF(DMP_LOG) WRITE(*,1000)'I',iLB,iMin3,iUB,iMax3
138              CALL MFIX_EXIT(myPE)
139           ELSEIF(jLB > jUB .OR. jLB<jMin3 .OR. jUB>jMax3) THEN
140              IF(DMP_LOG) WRITE(*,1000)'J',jLB,jMin3,jUB,jMax3
141              CALL MFIX_EXIT(myPE)
142           ELSEIF(kLB > kUB .OR. kLB<kMin3 .OR. kUB>kMax3) THEN
143              IF(DMP_LOG) WRITE(*,1000)'K',kLB,kMin3,kUB,kMax3
144              CALL MFIX_EXIT(myPE)
145           ENDIF
146     
147     
148     ! Calculate the coeffs for the space filling curve.
149           Am_c0 = 1 - jLB
150           Am_c1 = (jUB-jLB+1)
151           Am_c2 = (jUB-jLB+1)*(iUB-iLB+1)
152           Am_c0 =  Am_c0  - Am_c1*iLB - Am_c2*kLB
153     
154           dbgDIMN = (1+iUB-iLB) * (1+jUB-jLB) * (1+kUB-kLB)
155     
156     ! Allocate the output vector for Am
157           if(inflateAM) then
158              Allocate( outAm( dbgDIMN) )
159           else
160              if(do_K) then
161                 Allocate( outAm(-3:3) )
162              else
163                 Allocate( outAm(-2:2) )
164              endif
165           endif
166     
167     ! Report some basic data to the screen.
168           IF(DMP_LOG) THEN
169              WRITE(*,"(4/,3x,'Matrix Map:')")
170              write(*,"(/,5X,'Domain Limits >')")
171              write(*,"(7X,'I: ',I4,2x,I4)") iMin3, iMax3
172              write(*,"(7X,'J: ',I4,2x,I4)") jMin3, jMax3
173              if(do_K) write(*,"(7X,'K: ',I4,2x,I4)") kMin3, kMax3
174     
175              write(*,"(/5x,'Local IJK Coeffs >')")
176              write(*,"(7x,'C0: ',I9)") C0
177              write(*,"(7x,'C1: ',I9)") C1
178              write(*,"(7x,'C2: ',I9)") C2
179     
180              write(*,"(/5X,'Extraction Region: >')")
181              write(*,"(7X,'I: ',I4,2x,I4)") iLB, iUB
182              write(*,"(7X,'J: ',I4,2x,I4)") jLB, jUB
183              if(do_K) write(*,"(7X,'K: ',I4,2x,I4)") kLB, kUB
184     
185              write(*,"(/5x,'dbgIJK Coeffs >')")
186              write(*,"(7x,'Am_C0: ',I9)") Am_C0
187              write(*,"(7x,'Am_C1: ',I9)") Am_C1
188              write(*,"(7x,'Am_C2: ',I9)") Am_C2
189     
190              if(inflateAM) then
191                 write(*,"(/5x,'A_M is going to be inflated.')")
192              else
193                 write(*,"(/5x,'A_M is NOT going to be inflated.')")
194              endif
195           ENDIF
196     
197     ! Set the flag indicating that the initialization routine was called.
198           initNotCalled = .FALSE.
199     
200           do proc = 0, numPEs-1
201              if(myPE == proc) then
202                 write(*,"(/3x,'Proc: ',I3)")proc
203                 write(*,"(5x,'I start/end 1:',2(2x,I3))") istart1, iend1
204                 write(*,"(5x,'J start/end 1:',2(2x,I3))") jstart1, jend1
205                 if(do_K)write(*,"(5x,'K start/end 1:',2(2x,I3))") kstart1, kend1
206              endif
207     #ifdef MPI
208              CALL MPI_Barrier(MPI_COMM_WORLD,mpierr)
209     #endif
210           enddo
211     
212     ! Unlock dbg routines.
213           dbgLock = .FALSE.
214     
215           RETURN
216     
217      1000 FORMAT(/1X,70('*')/' From: initExtract',/' Error 1000:',         &
218              ' Invalid parameters.  Axis: ',A1,/8x,'LB =',I4,3x,'Min2 =',  &
219              I12,/8x,'UB =',I4,3x,'Max2 =',I12,/1x,70('*'),2/)
220     
221           END SUBROUTINE initExtract
222     
223     !----------------------------------------------------------------------!
224     !                                                                      !
225     !                                                                      !
226     !                                                                      !
227     !----------------------------------------------------------------------!
228           SUBROUTINE matrixExtract(A_m, B_m, M, VAR, PASS)
229     
230           use compar
231           use fldvar
232     !      use geometry
233     !      use indices
234     !      use run
235     !      use sendrecv
236           use functions
237           use param, only: dimension_3, dimension_m
238           use param1, only: zero
239     
240           IMPLICIT NONE
241     
242     ! Septadiagonal matrix A_m
243           DOUBLE PRECISION, intent(in) :: A_m(DIMENSION_3, -3:3, 0:DIMENSION_M)
244     ! Vector b_m
245           DOUBLE PRECISION, intent(in) :: B_m(DIMENSION_3, 0:DIMENSION_M)
246     
247           INTEGER, intent(in) :: M
248     
249           CHARACTER(len=*), intent(in), optional :: VAR
250     
251           INTEGER, intent(in), optional :: PASS
252     
253     ! Septadiagonal matrix A_m
254           DOUBLE PRECISION :: lA_m(-3:3)
255     ! Vector b_m
256           DOUBLE PRECISION :: lB_m
257     ! Neighbor info for debugging.
258           INTEGER :: NBGHS(-3:3)
259     
260           LOGICAL :: lexist
261           CHARACTER(len=64)  :: AmFName
262           CHARACTER(len=64)  :: BmFName
263     
264           INTEGER, parameter :: AmUnit = 9657
265           INTEGER, parameter :: BmUnit = 9658
266     
267           INTEGER :: IJK, I, J, K, OWNER
268     
269     ! If the initialization routine was not called, flag the error and exit.
270           IF(initNotCalled)THEN
271              IF(DMP_LOG) THEN
272                 WRITE(*,"(3/,' Fatal Error in matrixExtract!')")
273                 WRITE(*,"(' The initialization routine was never called.')")
274                 WRITE(*,"(' USR0 should contain a call to initExtract.')")
275                 WRITE(*,"(' Forcing a hard stop.')")
276              ENDIF
277              CALL MFIX_EXIT(myPE)
278           ENDIF
279     
280           AmFName=''
281           BmFName=''
282           IF(present(VAR) .AND. present(PASS)) THEN
283              write(AmFName,"('Am_',A,'_',I6.6,'.csv')")trim(adjustl(VAR)),PASS
284              write(BmFName,"('Bm_',A,'_',I6.6,'.csv')")trim(adjustl(VAR)),PASS
285           ELSEIF(present(VAR)) THEN
286              write(AmFName,"('Am_',A,'.csv')")trim(adjustl(VAR))
287              write(BmFName,"('Bm_',A,'.csv')")trim(adjustl(VAR))
288           ELSEIF(present(PASS)) THEN
289              write(AmFName,"('Am_',I6.6,'.csv')") PASS
290              write(BmFName,"('Bm_',I6.6,'.csv')") PASS
291           ELSE
292              AmFName='Am.csv'
293              BmFName='Bm.csv'
294           ENDIF
295     
296     ! Rank 0 opens the file.
297           IF(myPE == PE_IO) THEN
298              inquire(file=trim(AmFName),exist=lexist)
299              IF(lexist) THEN
300                 open(AmUnit,file=trim(AmFName), status='replace',convert='big_endian')
301              ELSE
302                 open(AmUnit,file=trim(AmFName), status='new',convert='big_endian')
303              ENDIF
304     
305              inquire(file=trim(BmFName),exist=lexist)
306              IF(lexist) THEN
307                 open(BmUnit,file=trim(BmFName), status='replace',convert='big_endian')
308              ELSE
309                 open(BmUnit,file=trim(BmFName), status='new',convert='big_endian')
310              ENDIF
311           ENDIF
312     
313     ! Everyone waits for the file to open.
314     #ifdef MPI
315           CALL MPI_Barrier(MPI_COMM_WORLD,mpierr)
316     #endif
317     
318           do K = kLB, kUB
319           do I = iLB, iUB
320           do J = jLB, jUB
321     
322     ! All ranks set lAm and lBm to zero.
323              OWNER = 0
324              NBGHS = 0
325              lB_m = ZERO
326              lA_m(-3:3) = ZERO
327     
328     ! Only the rank that owns IJK populates lAm and lBm.
329              if(IS_ON_myPE_owns(I,J,K)) then
330                 IJK = funIJK(I,J,K)
331                 lA_m(-3:3) = A_M(IJK,-3:3,M)
332                 lB_m = B_M(IJK,M)
333     ! Incorporate neighbor information for debugging.
334                 if(dbgMode) then
335                    OWNER = myPE
336                    if(do_K) Nbghs(-3) = Bottom_of(IJK)
337                    NBGHS(-2) = South_of(IJK)
338                    NBGHS(-1) = West_of(IJK)
339                    NBGHS( 0) = IJK
340                    NBGHS( 1) = East_of(IJK)
341                    NBGHS( 2) = North_of(IJK)
342                    if(do_K) Nbghs( 3) = Top_of(IJK)
343                 endif
344              endif
345     
346     ! Global sum updates lAm and lBm for all ranks.
347              CALL global_all_sum(lA_m)
348              CALL global_all_sum(lB_m)
349              if(dbgMode) then
350                 CALL global_all_sum(OWNER)
351                 CALL global_all_sum(NBGHS)
352              endif
353     
354              if(myPE==PE_IO) then
355                 CALL Am_to_Aout(I, J, K, NBGHS, OWNER, lA_m)
356                 CALL write_AoutBm(lB_m)
357              endif
358     
359     ! Everyone waits for the owner to complete this cycle.
360     #ifdef MPI
361              CALL MPI_Barrier(MPI_COMM_WORLD, mpierr)
362     #endif
363     
364           enddo
365           enddo
366           enddo
367     
368     ! Rank 0 closes the file.
369           IF(myPE == PE_IO) then
370              close(AmUnit)
371              close(BmUnit)
372           ENDIF
373     
374           return
375     
376           contains
377     
378     !``````````````````````````````````````````````````````````````````````!
379     !                                                                      !
380     !                                                                      !
381     !``````````````````````````````````````````````````````````````````````!
382           SUBROUTINE Am_to_Aout(I, J, K, lNBGHS, lOWNER, lA_m)
383     
384           INTEGER, intent(in) :: I, J, K
385           INTEGER, intent(in) :: lNBGHS(-3:3), lOWNER
386           DOUBLE PRECISION, intent(in) :: lA_m(-3:3)
387     
388     ! Local process IJK for neighbor cell mapped to global Am output matrix.
389           INTEGER :: nIJK
390     ! Increments from I,J,K for neighbor calculations.
391           INTEGER :: ii, jj, kk
392     
393     ! Bounds for extracted section or walls.
394           INTEGER :: sMin, wMin, bMin
395           INTEGER :: eMax, nMax, tMax
396     
397     ! Initialize:
398           outAm = 0.0d0
399     
400     ! Set up bounds check values.
401           wMin = iLB;   eMax = (1+iUB-iLB)
402           sMin = jLB;   nMax = (1+jUB-jLB)
403           bMin = kLB;   tMax = (1+kUB-kLB)
404     
405     ! Print mapping info.
406           if(dbgMode) write(*,9003) lNBGHS(0), dbg_funijk(i,j,k),           &
407              i, j, k, lOWNER
408     
409     ! Bottom Neighboring Fluid Cell
410     !---------------------------------------------------------------------//
411           if(do_K) then
412              if(k > kLB) then
413                 ii = i;   jj = j;   kk = k-1
414                 nIJK = merge(dbg_funijk(ii,jj,kk), -3, inflateAM)
415                 outAm(nIJK) = lA_m(-3)
416     
417                 if(dbgMode) write(*,9000)'Bottom of ', lNBGHS(-3),         &
418                    nIJK, ii, jj, kk, lA_m(-3)
419              else
420                 if(dbgMode) write(*,9001)'Bottom of ', lA_m(-3)
421              endif
422           endif
423     
424     ! South Neighboring Fluid Cell
425     !---------------------------------------------------------------------//
426           if(j > jLB) then
427     
428              ii = i;   jj = j-1;   kk = k
429              nIJK = merge(dbg_funijk(ii,jj,kk), -2, inflateAM)
430              outAm(nIJK) = lA_m(-2)
431     
432              if(dbgMode) write(*,9000)'South of  ', lNBGHS(-2),            &
433                 nIJK, ii, jj, kk, lA_m(-2)
434           else
435              if(dbgMode) write(*,9001)'South of  ', lA_m(-2)
436           endif
437     
438     ! West Neighboring Fluid Cell
439     !---------------------------------------------------------------------//
440           if(i > iLB) then
441     
442              ii = i-1; jj = j;   kk = k
443              nIJK = merge(dbg_funijk(ii,jj,kk), -1, inflateAM)
444              outAm(nIJK) = lA_m(-1)
445     
446              if(dbgMode) write(*,9000)'West of   ', lNBGHS(-1),            &
447                 nIJK, ii, jj, kk, lA_m(-1)
448           else
449              if(dbgMode) write(*,9001)'West of   ', lA_m(-1)
450           endif
451     
452     ! Center Coefficient
453     !---------------------------------------------------------------------//
454     
455           ii = i;   jj = j;   kk = k
456           nIJK = merge(dbg_funijk(ii,jj,kk), 0, inflateAM)
457           outAm(nIJK) = lA_m(0)
458     
459           if(dbgMode) write(*,9000)'Cntr Coef ', lNBGHS(0),                &
460              nIJK, ii, jj, kk, lA_m(0)
461     
462     ! East Neighboring Fluid Cell
463     !---------------------------------------------------------------------//
464           if(i < iUB) then
465     
466              ii = i+1; jj = j;   kk = k
467              nIJK = merge(dbg_funijk(ii,jj,kk), 1, inflateAM)
468              outAm(nIJK) = lA_m( 1)
469     
470              if(dbgMode) write(*,9000)'East of   ', lNBGHS( 1),            &
471                 nIJK, ii, jj, kk, lA_m( 1)
472           else
473              if(dbgMode) write(*,9001)'East of   ', lA_m( 1)
474           endif
475     
476     ! North Neighboring Fluid Cell
477     !---------------------------------------------------------------------//
478           if(j < jUB) then
479     
480              ii = i; jj = j+1;   kk = k
481              nIJK = merge(dbg_funijk(ii,jj,kk), 2, inflateAM)
482              outAm(nIJK) = lA_m( 2)
483     
484              if(dbgMode) write(*,9000)'North of  ', lNBGHS(2),             &
485                 nIJK, ii, jj, kk, lA_m( 2)
486           else
487              if(dbgMode) write(*,9001)'North of  ', lA_m( 2)
488           endif
489     
490     ! Top Neighboring Fluid Cell
491     !---------------------------------------------------------------------//
492           if(do_K) then
493              if(k < kUB) then
494     
495                 ii = i; jj = j;   kk = k+1
496                 nIJK = merge(dbg_funijk(ii,jj,kk), 3, inflateAM)
497                 outAm(nIJK) = lA_m( 3)
498     
499                 if(dbgMode) write(*,9000)'Top of    ', lNBGHS( 3),         &
500                    nIJK, ii, jj, kk, lA_m( 3)
501              else
502                 if(dbgMode) write(*,9001)'Top of    ', lA_m( 3)
503              endif
504           endif
505     
506           return
507     
508      9000 Format(5x,A,':: ',I4,' --> ',I4,3x,'(',I3,',',I3,',',I3,         &
509              ') = ',F12.4)
510     
511      9001 Format(5x,A,':: ............ OoR ............ = ',F12.4)
512     
513     
514      9003 Format(//3x,'Mapping: ',I4,' --> ',I4,3x,'(',I3,',',I3,',',I3,   &
515              ')',4x,'Rank: ',I5)
516     
517           END SUBROUTINE Am_to_Aout
518     
519     
520     !``````````````````````````````````````````````````````````````````````!
521     !                                                                      !
522     !                                                                      !
523     !``````````````````````````````````````````````````````````````````````!
524           SUBROUTINE write_AoutBm(lB_m)
525     
526           DOUBLE PRECISION, intent(in) :: lB_m
527     
528           INTEGER :: lStart
529           INTEGER :: lEnd
530           INTEGER :: IJK
531     
532           lStart = lbound(outAm,1)
533           lEnd = ubound(outAm,1)
534     
535           if(dbgMode) then
536              do IJK = lStart, lEnd-1
537                 write(AmUnit,"(F12.4,',')",advance='no')outAm(IJK)
538              enddo
539              write(AmUnit,"(F12.4)")outAm(lEnd)
540              write(BmUnit,"(F12.4)")lB_m
541     
542           else
543              do IJK = lStart, lEnd-1
544                 write(AmUnit,"(e14.6,',')",advance='no')outAm(IJK)
545              enddo
546              write(AmUnit,"(e14.6)")outAm(lEnd)
547              write(BmUnit,"(e14.6)")lB_m
548           endif
549     
550           END SUBROUTINE write_AoutBm
551     
552           END SUBROUTINE matrixExtract
553     
554     !----------------------------------------------------------------------!
555     !                                                                      !
556     !                                                                      !
557     !                                                                      !
558     !----------------------------------------------------------------------!
559           SUBROUTINE dbg_write(lMsg, Flush)
560     
561           use compar
562     !      use sendrecv
563     
564           implicit none
565     
566           CHARACTER(len=*), intent(in) :: lMsg
567     
568           LOGICAL, optional, intent(in) :: FLUSH
569           LOGICAL :: lFlush
570     
571           if(dbgMode)then
572              lFlush = merge(FLUSH, .FALSE., present(FLUSH))
573     
574              if(myPE == PE_IO) then
575                 if(lFlush) write(*,"(' ')")
576                 write(*,"(3x,A)") trim(lMsg)
577              endif
578     #ifdef MPI
579                 CALL MPI_Barrier(MPI_COMM_WORLD,mpierr)
580     #endif
581           endif
582     
583           RETURN
584           END SUBROUTINE DBG_WRITE
585     
586     !----------------------------------------------------------------------!
587     !                                                                        !
588     !                                                                      !
589     !                                                                      !
590     !----------------------------------------------------------------------!
591           SUBROUTINE arrayExtract_init(vName)
592     
593           use compar, only: myPE
594           use funits, only: DMP_LOG
595     
596           implicit none
597     
598     ! Variable named used to create file name.
599           CHARACTER(len=*), intent(in) :: vName
600     
601     ! Variables for opening files.
602           LOGICAL :: lexist
603           CHARACTER(len=64)  :: VarFName
604     ! Integer Error Flag.
605           INTEGER :: iErr
606     
607     ! If the initialization routine was not called, flag the error and exit.
608           IF(initNotCalled)THEN
609              IF(DMP_LOG) WRITE(*,1000)
610              CALL MFIX_EXIT(myPE)
611           ENDIF
612     
613     ! Verify that the files are unlocked.
614           IF(dbgLock)THEN
615              IF(DMP_LOG) WRITE(*,1001)
616              CALL MFIX_EXIT(myPE)
617           ELSE
618              dbgLock = .TRUE.
619           ENDIF
620     
621     ! Allocate indice arrays.
622           IF(allocated(i_ofBuff))THEN
623              IF(DMP_LOG) WRITE(*,1002) 'i_ofBuff'
624              CALL MFIX_EXIT(myPE)
625           ELSE
626              allocate( i_ofBuff(dbgDIMN) )
627              i_ofBuff = 0
628           ENDIF
629     
630           IF(allocated(j_ofBuff))THEN
631              IF(DMP_LOG) WRITE(*,1002) 'j_ofBuff'
632              CALL MFIX_EXIT(myPE)
633           ELSE
634              allocate( j_ofBuff(dbgDIMN) )
635              j_ofBuff = 0
636           ENDIF
637     
638           IF(allocated(k_ofBuff))THEN
639              IF(DMP_LOG) WRITE(*,1002) 'k_ofBuff'
640              CALL MFIX_EXIT(myPE)
641           ELSE
642              allocate( k_ofBuff(dbgDIMN) )
643              k_ofBuff = 0
644           ENDIF
645     
646           IF(allocated(ijk_Buff))THEN
647              IF(DMP_LOG) WRITE(*,1002) 'ijk_Buff'
648              CALL MFIX_EXIT(myPE)
649           ELSE
650              allocate( ijk_Buff(dbgDIMN) )
651              ijk_Buff = 0
652           ENDIF
653     
654     ! Construct the file name.
655           VarFName=''
656           IF(fApnd) THEN
657              write(VarFName,"(A,'.csv')")                                  &
658                 trim(adjustl(vName))
659           ELSEIF(fPass /= UNDEFINED_I) THEN
660              write(VarFName,"(A,'_',I6.6,'.csv')")                         &
661                 trim(adjustl(vName)), fPass
662           ELSE
663              write(VarFName,"(A,'.csv')")                                  &
664                 trim(adjustl(vName))
665           ENDIF
666     
667     ! Open the file  -- Rank 0 ONLY.
668           inquire(file=trim(VarFName), exist=lExist)
669           IF(myPE == PE_IO) THEN
670              IF(lExist) THEN
671                 IF(fApnd) THEN
672                    open(dbgUnit,file=trim(VarFName),                       &
673                       status='old', position='append', iostat=iErr,convert='big_endian')
674                 ELSE
675                    open(dbgUnit,file=trim(VarFName),                       &
676                       status='replace', iostat=iErr,convert='big_endian')
677                 ENDIF
678              ELSE
679                 open(dbgUnit,file=trim(VarFName),                          &
680                    status='new', iostat=iErr,convert='big_endian')
681              ENDIF
682           ENDIF
683           CALL BCAST(iErr, PE_IO)
684           IF(iErr /= 0)THEN
685              IF(myPE == PE_IO) write(*,1003) trim(VarFName)
686              CALL MFIX_EXIT(myPE)
687           ENDIF
688     
689     ! Set printIJK Flag. Only print the IJK values if this is the first time
690     ! for an append file.
691           pwIJK = merge(.FALSE., fwIJK, fwIJK.AND.lExist)
692     
693           RETURN
694     
695      1000 FORMAT(//1X,70('*')/' From: dbg_mod -> arrayExtract_init',/,     &
696              ' Error 1000: The initialization routine was never called.',  &
697              ' Include the',/' following in USR0: CALL initExtract.',2/,   &
698              ' These arguments are used to specify a domain subset to',    &
699              ' extract. If',/' not defined, the entire domain (MIN3/MAX3)',&
700              ' is extracted.',2/,' Optional arguments:',/,                 &
701              3x,'iLow - lower I index;  iHgh - Upper I index  (X-axis)',/, &
702              3x,'jLow - lower J index;  jHgh - Upper J index  (Y-axis)',/, &
703              3x,'kLow - lower K index;  kHgh - Upper K index  (Z-axis)',/  &
704              1x,70('*'),2/)
705     
706      1001 FORMAT(//1X,70('*')/' From: dbg_mod -> arrayExtract_init',/,     &
707              ' Error 1001: dgbLock is set. Something must have failed.',/  &
708              1x,70('*'),2/)
709     
710      1002 FORMAT(//1X,70('*')/' From: dbg_mod -> arrayExtract_init',/,     &
711              ' Error 1002: Buffer already allocated: ',A,/1x,70('*'),2/)
712     
713      1003 FORMAT(//1X,70('*')/' From: dbg_mod -> arrayExtract_init',/,     &
714              ' Error 1002: Error opening file: ',A,/1x,70('*'),2/)
715     
716           END SUBROUTINE arrayExtract_init
717     
718     !----------------------------------------------------------------------!
719     !                                                                        !
720     !                                                                      !
721     !                                                                      !
722     !----------------------------------------------------------------------!
723           SUBROUTINE arrayExtract_finl
724     
725           implicit none
726     
727           IF(allocated(i_ofBuff)) deallocate(i_ofBuff)
728           IF(allocated(j_ofBuff)) deallocate(j_ofBuff)
729           IF(allocated(k_ofBuff)) deallocate(k_ofBuff)
730           IF(allocated(ijk_Buff)) deallocate(ijk_Buff)
731     
732           IF(myPE == PE_IO) close(dbgUnit)
733     
734           dbgLock = .FALSE.
735     
736           RETURN
737           END SUBROUTINE arrayExtract_finl
738     
739     !----------------------------------------------------------------------!
740     !                                                                        !
741     !                                                                      !
742     !                                                                      !
743     !----------------------------------------------------------------------!
744           SUBROUTINE arrayExtract_prnt(dType)
745     
746           implicit none
747     
748     ! Data type
749           CHARACTER(len=3), intent(in) :: dType
750     
751     ! Loop counter
752           INTEGER :: IJK
753     
754     ! Rank 0 writes and closes the file.
755           IF(myPE /= PE_IO) RETURN
756     
757           IF(fApnd) THEN
758     ! Print header info.
759              IF(pwIJK) THEN
760                 CALL WRITE_INDEX(i_ofBuff)
761                 CALL WRITE_INDEX(j_ofBuff)
762                 if(do_K) CALL WRITE_INDEX(k_ofBuff)
763                 CALL WRITE_INDEX(ijk_Buff)
764              ENDIF
765     
766     ! Store the pass count.
767              IF(fPass /= UNDEFINED_I)                                      &
768                 WRITE(dbgUnit,"(I14,',')",advance='no') fPass
769     ! Output the data.
770              SELECT CASE(dType)
771                 CASE('INT'); CALL WRITE_INT
772                 CASE('DBL'); CALL WRITE_DBL
773                 CASE('LOG'); CALL WRITE_LOG
774              END SELECT
775     
776           ELSE
777     
778              IF(fPass /= UNDEFINED_I)                                      &
779                 WRITE(dbgUnit,"(2x,'Pass: ',I8)") fPass
780     
781              DO IJK=1, dbgDIMN
782     ! Output the IJK info.
783                 IF(pwIJK)THEN
784                    WRITE(dbgUnit,"(4(I14,','))",advance='no') IJK,         &
785                       ijk_Buff(IJK), i_ofBuff(IJK), j_ofBuff(IJK)
786                    if(do_K)WRITE(dbgUnit,"(I14,',')",advance='no')         &
787                       k_ofBuff(IJK)
788                 ENDIF
789     ! Write the formatted output.
790                 SELECT CASE(dType)
791                    CASE('INT'); CALL WRITE_INT(IJK)
792                    CASE('DBL'); CALL WRITE_DBL(IJK)
793                    CASE('LOG'); CALL WRITE_LOG(IJK)
794                 END SELECT
795              ENDDO
796           ENDIF
797     
798           RETURN
799           END SUBROUTINE arrayExtract_prnt
800     
801     !----------------------------------------------------------------------!
802     !                                                                      !
803     !                                                                      !
804     !                                                                      !
805     !----------------------------------------------------------------------!
806           SUBROUTINE WRITE_INDEX(intArray)
807     
808           implicit none
809     
810           INTEGER, intent(in) :: intArray(dbgDIMN)
811     ! Looped IJK
812           INTEGER :: IJK
813     
814     ! Loop through array entries writing them in one long continuous line.
815           IF(fPass /= UNDEFINED_I) WRITE(dbgUnit,"(14X,',')",advance='no')
816           DO IJK = 1, dbgDIMN-1
817              WRITE(dbgUnit,"(I14,',')",advance='no')intArray(IJK)
818           ENDDO
819           WRITE(dbgUnit,"(I14)")intArray(dbgDIMN)
820     
821           RETURN
822           END SUBROUTINE WRITE_INDEX
823     
824     !----------------------------------------------------------------------!
825     !                                                                      !
826     !                                                                      !
827     !                                                                      !
828     !----------------------------------------------------------------------!
829           SUBROUTINE WRITE_INT(pIJK)
830     
831           implicit none
832     
833     ! Dummy arguments.
834           INTEGER, optional, intent(in) :: pIJK
835     ! Looped IJK
836           INTEGER :: IJK
837     
838           IF(present(pIJK)) THEN
839     ! Write the entry and return.
840              WRITE(dbgUnit,"(I14)") outBuff_i(pIJK)
841     
842           ELSE
843     ! Loop through array entries writing them in one long continuous line.
844              DO IJK = 1, dbgDIMN-1
845                 WRITE(dbgUnit,"(I14,',')",advance='no')outBuff_i(IJK)
846              ENDDO
847              WRITE(dbgUnit,"(I14)")outBuff_i(dbgDIMN)
848           ENDIF
849     
850           RETURN
851           END SUBROUTINE WRITE_INT
852     
853     !----------------------------------------------------------------------!
854     !                                                                      !
855     !                                                                      !
856     !                                                                      !
857     !----------------------------------------------------------------------!
858           SUBROUTINE WRITE_DBL(pIJK)
859     
860           implicit none
861     
862     ! Dummy arguments.
863           INTEGER, optional, intent(in) :: pIJK
864     ! Looped IJK
865           INTEGER :: IJK
866     
867           IF(present(pIJK)) THEN
868     ! Write the entry and return.
869              WRITE(dbgUnit,"(E14.6)") outBuff_dp(pIJK)
870           ELSE
871     ! Loop through array entries writing them in one long continuous line.
872              DO IJK = 1, dbgDIMN-1
873                 WRITE(dbgUnit,"(E14.6,',')",advance='no')outBuff_dp(IJK)
874              ENDDO
875              WRITE(dbgUnit,"(E14.6)")outBuff_dp(dbgDIMN)
876           ENDIF
877     
878           RETURN
879           END SUBROUTINE WRITE_DBL
880     
881     !----------------------------------------------------------------------!
882     !                                                                      !
883     !                                                                      !
884     !                                                                      !
885     !----------------------------------------------------------------------!
886           SUBROUTINE WRITE_LOG(pIJK)
887     
888           implicit none
889     
890     ! Dummy arguments.
891           INTEGER, optional, intent(in) :: pIJK
892     ! Looped IJK
893           INTEGER :: IJK
894           LOGICAL :: INTtoLOG
895     
896           IF(present(pIJK)) THEN
897     ! Write the entry and return.
898              INTtoLOG= (outBuff_i(pIJK) .eq. 1)
899              WRITE(dbgUnit,"(L14)") INTtoLOG
900     
901           ELSE
902     ! Loop through array entries writing them in one long continuous line.
903              DO IJK = 1, dbgDIMN-1
904                 INTtoLOG=(outBuff_i(IJK) .eq. 1)
905                 WRITE(dbgUnit,"(L14,',')",advance='no') INTtoLOG
906              ENDDO
907              INTtoLOG = (outBuff_i(dbgDIMN) .eq. 1)
908              WRITE(dbgUnit,"(L14)") INTtoLOG
909     
910           ENDIF
911     
912           RETURN
913           END SUBROUTINE WRITE_LOG
914     
915     !----------------------------------------------------------------------!
916     !                                                                      !
917     !                                                                      !
918     !                                                                      !
919     !----------------------------------------------------------------------!
920           SUBROUTINE arrayExtract_int(Array, VAR, PASS, APND, withIJK)
921     
922           use compar
923           use fldvar
924           use geometry
925           use indices
926           use mflux
927           use parallel
928           use param
929           use param1
930           use physprop
931           use run
932           use sendrecv
933           USE mpi_utility
934           USE functions
935     
936           implicit none
937     
938     ! Array to be extracted.
939           INTEGER, intent(in) :: Array(DIMENSION_3)
940     ! Variable named used to create file name.
941           CHARACTER(len=*), intent(in) :: VAR
942     ! Pass count if dumping at multiple locations or during iterations.
943           INTEGER, intent(in), optional :: PASS
944     ! Flag to append to previous output.
945           LOGICAL, intent(in), optional :: APND
946     ! Flag to write IJK values with output.
947           LOGICAL, intent(in), optional :: withIJK
948     ! Loop counters:
949           INTEGER :: I, J, K, IJK, dbgIJK
950     ! Debugging message.
951           CHARACTER(len=64) :: MSG
952     
953           MSG='Entered arrayExtract_int'
954           if(dbgMode) CALL DBG_WRITE(trim(MSG),FLUSH=.TRUE.)
955     
956     ! Set local flags based on optional arguments.
957           fPass = merge(PASS, UNDEFINED_I, present(PASS))
958           fApnd = merge(APND, .FALSE., present(APND))
959           fwIJK = merge(withIJK, .FALSE., present(withIJK))
960     
961           MSG='  > Calling arrayExtract_INIT'
962           if(dbgMode) CALL DBG_WRITE(trim(MSG))
963           CALL arrayExtract_init(VAR)
964     
965           MSG='  > Allocating outBuff_i'
966           if(dbgMode) CALL DBG_WRITE(trim(MSG))
967           allocate( outBuff_i(dbgDIMN) ); outBuff_i = 0
968     
969           MSG='  > Extracting array data.'
970           if(dbgMode) CALL DBG_WRITE(trim(MSG))
971     
972           do K = kLB, kUB
973           do I = iLB, iUB
974           do J = jLB, jUB
975              if(IS_ON_myPE_owns(I,J,K)) then
976              IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
977                 IJK = funIJK(I,J,K)
978                 dbgIJK = dbg_funijk(I,J,K)
979     ! Only the rank that owns IJK populates lAm and lBm.
980                 outBuff_i(dbgIJK) = Array(IJK)
981     ! Store indicies.
982                 i_ofBuff(dbgIJK) = I
983                 j_ofBuff(dbgIJK) = J
984                 k_ofBuff(dbgIJK) = K
985                 ijk_Buff(dbgIJK) = IJK
986              endif
987           enddo
988           enddo
989           enddo
990     
991           MSG='  > Collecting array data.'
992           if(dbgMode) CALL DBG_WRITE(trim(MSG))
993     ! Global sum updates lAm and lBm for all ranks.
994           CALL global_all_sum(outBuff_i)
995           if(pwIJK) then
996              CALL global_all_sum(i_ofBuff)
997              CALL global_all_sum(j_ofBuff)
998              CALL global_all_sum(ijk_Buff)
999              if(do_K)CALL global_all_sum(k_ofBuff)
1000           endif
1001     
1002     ! Write the data.
1003           MSG='  > Calling arrayExtract_prnt.'
1004           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1005           CALL arrayExtract_prnt('INT')
1006     
1007     ! Clean up the used memory.
1008           MSG='  > Calling arrayExtract_finl.'
1009           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1010     
1011           if(allocated(outBuff_i)) deallocate(outBuff_i)
1012           CALL arrayExtract_finl
1013     
1014           MSG='Leaving arrayExtract_int.'
1015           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1016     
1017           RETURN
1018           END SUBROUTINE arrayExtract_int
1019     
1020     !----------------------------------------------------------------------!
1021     !                                                                      !
1022     !                                                                      !
1023     !                                                                      !
1024     !----------------------------------------------------------------------!
1025           SUBROUTINE arrayExtract_dbl(Array, VAR, PASS, APND, withIJK)
1026     
1027           use compar
1028           use fldvar
1029           use geometry
1030           use indices
1031           use mflux
1032           use parallel
1033           use param
1034           use param1
1035           use physprop
1036           use run
1037           use sendrecv
1038           USE mpi_utility
1039           USE functions
1040     
1041           implicit none
1042     
1043     ! Array to be extracted.
1044           DOUBLE PRECISION, intent(in) :: Array(DIMENSION_3)
1045     ! Variable named used to create file name.
1046           CHARACTER(len=*), intent(in) :: VAR
1047     ! Pass count if dumping at multiple locations or during iterations.
1048           INTEGER, intent(in), optional :: PASS
1049     ! Flag to append to previous output.
1050           LOGICAL, intent(in), optional :: APND
1051     ! Flag to write IJK values with output.
1052           LOGICAL, intent(in), optional :: withIJK
1053     
1054     ! Loop indices:
1055           INTEGER :: I, J, K, IJK, dbgIJK
1056     ! Debugging message.
1057           CHARACTER(len=64) :: MSG
1058     
1059           MSG='Entered arrayExtract_dbl'
1060           CALL DBG_WRITE(trim(MSG), flush=.TRUE.)
1061     
1062     ! Set local flags based on optional arguments.
1063           fPass = merge(PASS, UNDEFINED_I, present(PASS))
1064           fApnd = merge(APND, .FALSE., present(APND))
1065           fwIJK = merge(withIJK, .FALSE., present(withIJK))
1066     
1067           MSG='  > Calling arrayExtract_INIT'
1068           CALL DBG_WRITE(trim(MSG))
1069           CALL arrayExtract_init(VAR)
1070     
1071     ! Allocate the global storage array.
1072           allocate( outBuff_dp(dbgDIMN) ); outBuff_dp = 0
1073     
1074           MSG='  > Extracting array data.'
1075           CALL DBG_WRITE(trim(MSG))
1076     
1077     ! Extract the data... one cell at a time.
1078           do K = kLB, kUB
1079           do I = iLB, iUB
1080           do J = jLB, jUB
1081              if(IS_ON_myPE_owns(I,J,K)) then
1082              IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
1083                 IJK = funIJK(I,J,K)
1084                 dbgIJK = dbg_funijk(I,J,K)
1085     ! Only the rank that owns IJK populates lAm and lBm.
1086                 outBuff_dp(dbgIJK) = Array(IJK)
1087     ! Store indicies.
1088                 i_ofBuff(dbgIJK) = I
1089                 j_ofBuff(dbgIJK) = J
1090                 k_ofBuff(dbgIJK) = K
1091                 ijk_Buff(dbgIJK) = IJK
1092              endif
1093           enddo
1094           enddo
1095           enddo
1096     
1097     ! Global sum updates lAm and lBm for all ranks.
1098           MSG='  > Collecting array data.'
1099           CALL DBG_WRITE(trim(MSG))
1100     
1101           CALL global_all_sum(outBuff_dp)
1102           if(pwIJK) then
1103              CALL global_all_sum(i_ofBuff)
1104              CALL global_all_sum(j_ofBuff)
1105              CALL global_all_sum(ijk_Buff)
1106              if(do_K)CALL global_all_sum(k_ofBuff)
1107           endif
1108     
1109     ! Write the data.
1110           MSG='  > Calling arrayExtract_prnt.'
1111           CALL DBG_WRITE(trim(MSG))
1112           CALL arrayExtract_prnt('DBL')
1113     
1114     ! Clean up the used memory.
1115           MSG='  > Calling arrayExtract_finl.'
1116           CALL DBG_WRITE(trim(MSG))
1117           if(allocated(outBuff_dp)) deallocate(outBuff_dp)
1118           CALL arrayExtract_finl
1119     
1120           MSG='Leaving arrayExtract_dbl.'
1121           CALL DBG_WRITE(trim(MSG))
1122     
1123           RETURN
1124           END SUBROUTINE arrayExtract_dbl
1125     
1126     !----------------------------------------------------------------------!
1127     !                                                                      !
1128     !                                                                      !
1129     !                                                                      !
1130     !----------------------------------------------------------------------!
1131           SUBROUTINE arrayExtract_log(Array, VAR, PASS, APND, withIJK)
1132     
1133           use compar
1134           use fldvar
1135           use geometry
1136           use indices
1137           use mflux
1138           use parallel
1139           use param
1140           use param1
1141           use physprop
1142           use run
1143           use sendrecv
1144           USE mpi_utility
1145           USE functions
1146     
1147           implicit none
1148     
1149     ! Array to be extracted.
1150           LOGICAL, intent(in) :: Array(DIMENSION_3)
1151     ! Variable named used to create file name.
1152           CHARACTER(len=*), intent(in) :: VAR
1153     ! Pass count if dumping at multiple locations or during iterations.
1154           INTEGER, intent(in), optional :: PASS
1155     ! Flag to append to previous output.
1156           LOGICAL, intent(in), optional :: APND
1157     ! Flag to write IJK values with output.
1158           LOGICAL, intent(in), optional :: withIJK
1159     ! Loop counters:
1160           INTEGER :: I, J, K, IJK, dbgIJK
1161     ! Debugging message.
1162           CHARACTER(len=64) :: MSG
1163     
1164           MSG='Entered arrayExtract_log'
1165           if(dbgMode) CALL DBG_WRITE(trim(MSG),FLUSH=.TRUE.)
1166     
1167     ! Set local flags based on optional arguments.
1168           fPass = merge(PASS, UNDEFINED_I, present(PASS))
1169           fApnd = merge(APND, .FALSE., present(APND))
1170           fwIJK = merge(withIJK, .FALSE., present(withIJK))
1171     
1172           MSG='  > Calling arrayExtract_INIT'
1173           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1174           CALL arrayExtract_init(VAR)
1175     
1176           MSG='  > Allocating outBuff_i'
1177           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1178           allocate( outBuff_i(dbgDIMN) ); outBuff_i = 0
1179     
1180           MSG='  > Extracting array data.'
1181           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1182     
1183           do K = kLB, kUB
1184           do I = iLB, iUB
1185           do J = jLB, jUB
1186              if(IS_ON_myPE_owns(I,J,K)) then
1187              IF (DEAD_CELL_AT(I,J,K)) CYCLE  ! skip dead cells
1188                 IJK = funIJK(I,J,K)
1189                 dbgIJK = dbg_funijk(I,J,K)
1190     ! Only the rank that owns IJK populates lAm and lBm.
1191     ! Convert the logical to interger.
1192                 outBuff_i(dbgIJK) = merge(1,0,Array(IJK))
1193     ! Store indicies.
1194                 i_ofBuff(dbgIJK) = I
1195                 j_ofBuff(dbgIJK) = J
1196                 k_ofBuff(dbgIJK) = K
1197                 ijk_Buff(dbgIJK) = IJK
1198              endif
1199           enddo
1200           enddo
1201           enddo
1202     
1203           MSG='  > Collecting array data.'
1204           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1205     ! Global sum updates lAm and lBm for all ranks.
1206           CALL global_all_sum(outBuff_i)
1207           if(pwIJK) then
1208              CALL global_all_sum(i_ofBuff)
1209              CALL global_all_sum(j_ofBuff)
1210              CALL global_all_sum(ijk_Buff)
1211              if(do_K)CALL global_all_sum(k_ofBuff)
1212           endif
1213     
1214     ! Write the data.
1215           MSG='  > Calling arrayExtract_prnt.'
1216           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1217           CALL arrayExtract_prnt('LOG')
1218     
1219     ! Clean up the used memory.
1220           MSG='  > Calling arrayExtract_finl.'
1221           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1222     
1223           if(allocated(outBuff_i)) deallocate(outBuff_i)
1224           CALL arrayExtract_finl
1225     
1226           MSG='Leaving arrayExtract_log.'
1227           if(dbgMode) CALL DBG_WRITE(trim(MSG))
1228     
1229           RETURN
1230           END SUBROUTINE arrayExtract_log
1231     
1232           END MODULE dbg
1233