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