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