File: /nfs/home/0/users/jenkins/mfix.git/model/mfix.f

1     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
2     !                                                                      !
3     !  Subroutine: MFIX                                                    !
4     !  Author: M. Syamlal                                 Date: 29-JAN-92  !
5     !                                                                      !
6     !  Purpose: The main module in the MFIX program                        !
7     !                                                                      !
8     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
9     !
10     !> \mainpage Multiphase Flow with Interphase eXchanges
11     !!
12     !! MFIX is a general-purpose computer code developed at the National
13     !! Energy Technology Laboratory, NETL, for describing the hydrodynamics,
14     !! heat transfer, and chemical reactions in fluid-solid systems.
15     !!
16     !! It has been used for describing bubbling and circulating fluidized
17     !! beds and spouted beds. MFiX calculations give transient data on the
18     !! three-dimensional distribution of pressure, velocity, temperature,
19     !! and species mass fractions. MFiX code is based on a generally
20     !! accepted set of multiphase flow equations. The code is used as a
21     !! "test-stand" for testing and developing multiphase flow constitutive
22     !!  equations.
23     !!
24     !! \section Notice
25     !! Neither the United States Government nor any agency thereof, nor any
26     !! of their employees, makes any warranty, expressed or implied, or
27     !! assumes any legal liability or responsibility for the accuracy,
28     !! completeness, or usefulness of any information, apparatus, product,
29     !! or process disclosed or represents that its use would not infringe
30     !! privately owned rights.
31     !!
32     !! * MFIX is provided without any user support for applications in the
33     !!   user's immediate organization. It should not be redistributed in
34     !!   whole or in part.
35     !!
36     !! * The use of MFIX is to be acknowledged in any published paper based
37     !!   on computations using this software by citing the MFIX theory
38     !!   manual. Some of the submodels are being developed by researchers
39     !!   outside of NETL. The use of such submodels is to be acknowledged
40     !!   by citing the appropriate papers of the developers of the submodels.
41     !!
42     !! * The authors would appreciate receiving any reports of bugs or other
43     !!   difficulties with the software, enhancements to the software, and
44     !!   accounts of practical applications of this software.
45     !!
46     !! \section Disclaimer
47     !! This report was prepared as an account of work sponsored by an agency
48     !! of the United States Government. Neither the United States Government
49     !! nor any agency thereof, nor any of their employees, makes any
50     !! warranty, express or implied, or assumes any legal liability or
51     !! responsibility for the accuracy, completeness, or usefulness of any
52     !! information, apparatus, product, or process disclosed, or represents
53     !! that its use would not infringe privately owned rights. Reference
54     !! herein to any specific commercial product, process, or service by
55     !! trade name, trademark, manufacturer, or otherwise does not
56     !! necessarily constitute or imply its endorsement, recommendation, or
57     !! favoring by the United States Government or any agency thereof. The
58     !! views and opinions of authors expressed herein do not necessarily
59     !! state or reflect those of the United States Government or any
60     !! agency thereof.
61     
62           PROGRAM MFIX
63     
64     !-----------------------------------------------
65     ! Modules
66     !-----------------------------------------------
67           USE param
68           USE param1
69           USE run
70           USE time_cpu
71           USE funits
72           USE output
73           USE compar
74           USE mpi_utility
75           USE parallel_mpi
76           USE discretelement
77           USE mfix_pic
78           USE cdist
79           USE MFIX_netcdf
80           USE fldvar
81           USE cutcell
82           USE quadric
83           USE dashboard
84           USE qmom_kinetic_equation
85           USE functions
86     
87           USE vtk, only : WRITE_VTK_FILES
88     
89           use error_manager
90     
91           IMPLICIT NONE
92     !-----------------------------------------------
93     ! Local variables
94     !-----------------------------------------------
95     ! Final value of CPU time.
96           DOUBLE PRECISION :: CPU1
97     ! time used for computations.
98           DOUBLE PRECISION :: CPUTIME_USED, WALLTIME_USED
99     ! CPU time unit.
100           CHARACTER(LEN=4) :: TUNIT
101     ! Save TIME in input file for RESTART_2
102           DOUBLE PRECISION :: TIME_SAVE
103     ! Temporary storage for DT
104           DOUBLE PRECISION :: DT_tmp
105     ! loop counter
106           INTEGER :: L
107     ! Error index
108           INTEGER :: IER
109     ! DISTIO variable for specifying the mfix version
110           CHARACTER(LEN=512) :: version
111     ! environment variable
112           CHARACTER(LEN=512) :: omp_num_threads
113           INTEGER :: length
114           INTEGER :: status
115     
116     !$      INTEGER num_threads, threads_specified, omp_id
117     !$      INTEGER omp_get_num_threads
118     !$      INTEGER omp_get_thread_num
119     
120           DOUBLE PRECISION :: WALL_TIME
121     !-----------------------------------------------
122     
123     
124     
125     ! DISTIO
126     ! If you change the value below in this subroutine, you must also
127     ! change it in write_res0.f and the value should also be consistent
128     ! with the check in read_res0
129           version = 'RES = 01.6'
130     
131           bDoing_postmfix = .false.
132     
133     ! Invoke MPI initialization routines and get rank info.
134           CALL PARALLEL_INIT
135           CALL GEN_LOG_BASENAME
136     
137     ! we want only PE_IO to write out common error messages
138           DMP_LOG = (myPE == PE_IO)
139     
140     ! set the version.release of the software
141           ID_VERSION = '2015-1'
142     
143     ! set automatic restart flag to false
144     !      AUTOMATIC_RESTART = .FALSE.
145     !      ITER_RESTART      = 1
146     
147     ! specify the number of processors to be used
148     !$        call get_environment_variable("OMP_NUM_THREADS",omp_num_threads,length,status, .true.)
149     !$      if (status.eq.0 .and. length.ne.0) then
150     !$        read(omp_num_threads,*) threads_specified
151     !$      else
152     !$        WRITE(*,'(A,$)') 'Enter the number of threads to be used for SMP: '
153     !$        READ(*,*) threads_specified
154     !$      endif
155     
156     !$      call omp_set_num_threads(threads_specified)
157     
158     ! Find the number of processors used
159     !$omp  parallel
160     !$      num_threads = omp_get_num_threads()
161     !$      omp_id = omp_get_thread_num()
162     !$      if(omp_id.eq.0) Write(*,*)' Number of threads used for SMP = ',  num_threads
163     !$omp  end parallel
164     
165     
166     ! Set machine dependent constants
167           CALL MACHINE_CONS
168     
169     ! Get the date and time. They give the unique run_id in binary output
170     ! files
171           CALL GET_RUN_ID
172     
173     ! AEOLUS: stop trigger mechanism to terminate MFIX normally before batch
174     ! queue terminates. timestep at the beginning of execution
175           CALL CPU_TIME (CPU00)
176           WALL0 = WALL_TIME()
177     
178     ! Read input data, check data, do computations for IC and BC locations
179     ! and flows, and set geometry parameters such as X, X_E, DToDX, etc.
180           CALL GET_DATA
181     
182     ! Write the initial part of the standard output file
183           CALL WRITE_OUT0
184           IF(.NOT.CARTESIAN_GRID)  CALL WRITE_FLAGS
185     
186     ! Write the initial part of the special output file(s)
187           CALL WRITE_USR0
188     
189     !$    CALL START_LOG
190     !$    IF(DMP_LOG)WRITE (UNIT_LOG, *) ' '
191     !$    IF(DMP_LOG)WRITE (UNIT_LOG, *) ' Number of processors used = ', threads_specified
192     !$    IF(DMP_LOG)WRITE (UNIT_LOG, *) ' '
193     !$    CALL END_LOG
194     
195     !  setup for PC quickwin application
196           CALL PC_QUICKWIN
197     
198     
199           CALL INIT_ERR_MSG('MFIX')
200     
201     
202       101 CONTINUE
203     
204     
205     ! if not netcdf writes asked for ... globally turn off netcdf
206           if (MFIX_usingNETCDF()) then
207              bGlobalNetcdf = .false.
208              do L = 1,20
209                 if (bWrite_netcdf(L)) bGlobalNetcdf = .true.
210              enddo
211           endif
212     
213     
214           IF(AUTOMATIC_RESTART) THEN
215              RUN_TYPE = 'RESTART_1'
216              AUTOMATIC_RESTART = .FALSE.
217              ITER_RESTART = ITER_RESTART + 1
218              CALL CHECK_INITIAL_CONDITIONS
219              CALL CHECK_BOUNDARY_CONDITIONS
220              CALL CHECK_INTERNAL_SURFACES
221              CALL CHECK_POINT_SOURCES
222              CALL CHECK_CHEMICAL_RXNS
223              CALL SET_FLAGS
224              CALL SET_CONSTPROP
225           ENDIF
226     
227           DT_TMP = DT
228           SELECT CASE (TRIM(RUN_TYPE))
229     
230           CASE ('NEW')
231     ! Write the initial part of the restart files
232              CALL WRITE_RES0
233              DO L = 1, N_SPX
234                 CALL WRITE_SPX0 (L, 0)
235              ENDDO
236     
237            CASE ('RESTART_1')
238     ! Read the time-dependant part of the restart file
239              CALL READ_RES1
240              WRITE(ERR_MSG, 1010) TIME, NSTEP
241              CALL FLUSH_ERR_MSG()
242     
243           CASE ('RESTART_2')
244              TIME_SAVE = TIME
245     ! DISTIO
246              if (myPE .ne. PE_IO .and. bDist_IO .and. bStart_with_one_res) then
247                  write (unit_res,rec=1) version
248                  write (unit_res,rec=2) 4
249                  write (unit_res,rec=3) 4
250              endif
251     
252              CALL READ_RES1
253              TIME = TIME_SAVE
254     
255              WRITE(ERR_MSG, 1010) TIME, NSTEP
256              CALL FLUSH_ERR_MSG()
257     
258              CALL WRITE_RES0
259     
260     ! Writing the RES1 and SPX1 can only be done here when re-indexing is turned off
261     ! This will be done after the cell re-indexing is done later in this file.
262     ! This allows restarting independently of the re-indexing setting between
263     ! the previous and current run.
264              IF(.NOT.RE_INDEXING) THEN
265                 CALL WRITE_RES1
266                 DO L = 1, N_SPX
267                    CALL WRITE_SPX0 (L, 0)
268                    CALL WRITE_SPX1 (L, 0)
269                 END DO
270                 call write_netcdf(0,0,time)
271              ENDIF
272     
273           CASE DEFAULT
274              CALL START_LOG
275              IF(DMP_LOG)WRITE (UNIT_LOG, *) &
276                 ' MFIX: Do not know how to process'
277              IF(DMP_LOG)WRITE (UNIT_LOG, *) ' RUN_TYPE in data file'
278              CALL END_LOG
279              call mfix_exit(myPE)
280     
281           END SELECT
282     
283           call MPI_Barrier(MPI_COMM_WORLD,mpierr)
284     
285           IF (DT_TMP /= UNDEFINED) THEN
286              DT = MAX(DT_MIN,MIN(DT_MAX,DT))
287           ELSE
288              DT = DT_TMP
289           ENDIF
290     
291     ! Set arrays for computing indices. A secondary call is made
292     ! after cut cell-preprocessing to update array indices.
293           IF(CARTESIAN_GRID) THEN
294              CALL SET_INCREMENTS
295              CALL SET_INCREMENTS3
296           ENDIF
297     
298     
299     !      IF(.NOT.RE_INDEXING) CALL WRITE_IJK_VALUES
300     
301     
302     ! Set the flags for wall surfaces impermeable and identify flow
303     ! boundaries using FLAG_E, FLAG_N, and FLAG_T
304           CALL SET_FLAGS1
305     
306     !  Update flags for Cartesian_GRID.
307           IF(CARTESIAN_GRID) CALL CHECK_BC_FLAGS
308     
309     ! Calculate cell volumes and face areas
310           IF(.NOT.CARTESIAN_GRID)  THEN
311              CALL SET_GEOMETRY1
312     !       ELSE
313     !         CALL SET_GEOMETRY
314            ENDIF
315     
316     ! Find corner cells and set their face areas to zero
317           IF(.NOT.CARTESIAN_GRID)  THEN
318              CALL GET_CORNER_CELLS()
319           ELSE
320              IF (SET_CORNER_CELLS)  CALL GET_CORNER_CELLS ()
321           ENDIF
322     
323     ! Set constant physical properties
324           CALL SET_CONSTPROP
325     
326     ! Set initial conditions
327           CALL SET_IC
328     
329     ! Set point sources.
330           CALL SET_PS
331     
332     ! Set boundary conditions
333           CALL ZERO_NORM_VEL
334           CALL SET_BC0
335     !      IF(DISCRETE_ELEMENT) CALL MAKE_ARRAYS_DES
336     
337     ! JFD: cartesian grid implementation
338           IF(CARTESIAN_GRID) CALL CG_SET_BC0
339     !      IF(DEM_SOLIDS) CALL SET_BC_DEM
340     
341     ! Set gas mixture molecular weight
342           CALL SET_MW_MIX_G
343     
344     ! Set the pressure field for a fluidized bed
345           IF (RUN_TYPE == 'NEW') CALL SET_FLUIDBED_P
346     
347     ! Initialize densities.
348           IF (RUN_TYPE == 'NEW') CALL SET_RO_G
349           IF (RUN_TYPE == 'NEW') CALL SET_RO_S
350     
351     ! Initialize time dependent boundary conditions
352           CALL SET_BC1
353     
354     ! Check the field variable data and report errors.
355           IF(.NOT.CARTESIAN_GRID)  CALL CHECK_DATA_20
356     
357     
358     !=======================================================================
359     ! JFD: START MODIFICATION FOR RE-INDEXING CELLS
360     !=======================================================================
361           IF(CARTESIAN_GRID.AND.RE_INDEXING) THEN
362     
363              IF(myPE == PE_IO) THEN
364                 WRITE(*,"(72('='))")
365                 WRITE(*,*)' RE-INDEXING CELLS FOR CARTESIAN GRID...'
366              ENDIF
367              CALL RE_INDEX_ARRAYS
368     
369     
370     !      IF(myPE == PE_IO)print*,'Calling REPORT_BEST_IJK_SIZE:'
371     !       CALL REPORT_BEST_IJK_SIZE
372            CALL REPORT_BEST_PROCESSOR_SIZE
373     !      IF(myPE == PE_IO)print*,'Exiting MFIX after REPORT_BEST_IJK_SIZE.'
374     
375     
376              IF(myPE == PE_IO) WRITE(*,"(72('='))")
377     
378     ! In case of a RESTART_2, write the RES1 and SPX1 files here
379     ! This was commented out earlier in this file.
380              IF(RUN_TYPE == 'RESTART_2') THEN
381                 CALL WRITE_RES1
382                 DO L = 1, N_SPX
383                    CALL WRITE_SPX0 (L, 0)
384                    CALL WRITE_SPX1 (L, 0)
385                 END DO
386                 call write_netcdf(0,0,time)
387              ENDIF
388           ENDIF
389     
390     !=======================================================================
391     ! JFD: END MODIFICATION FOR RE-INDEXING CELLS
392     !=======================================================================
393     
394     ! Setup VTK data for regular (no cut cells) grid
395           IF(.NOT.CARTESIAN_GRID.AND.WRITE_VTK_FILES) CALL SETUP_VTK_NO_CUTCELL
396     
397           IF(DISCRETE_ELEMENT) CALL MAKE_ARRAYS_DES
398           IF(QMOMK) CALL QMOMK_MAKE_ARRAYS
399     
400     ! Set the inflow/outflow BCs for DEM solids
401           IF(DEM_SOLIDS) CALL SET_BC_DEM
402     ! Set the inflow/outflow BC for PIC solids
403           IF(PIC_SOLIDS) CALL SET_BC_PIC
404     
405     ! Set the inital properties of each particle.
406           IF(DEM_SOLIDS) CALL SET_IC_DEM
407     
408     ! AEOLUS: debug prints
409           if (DBGPRN_LAYOUT .or. bdist_io) then
410     !     write (*,*) myPE , ' E.4 ... version = ' , version(1:33)
411              call debug_write_layout()
412              call write_parallel_info()
413           endif
414     
415     ! Initializations for CPU time calculations in iterate
416           CPUOS = 0.
417           CALL CPU_TIME (CPU1)
418           CPU_NLOG = CPU1
419           TIME_NLOG = TIME - DT
420     
421     ! Get the initial value of CPU time
422           CALL CPU_TIME (CPU0)
423     
424     ! Find the solution of the equations from TIME to TSTOP at
425     ! intervals of DT
426           CALL TIME_MARCH
427           IF(AUTO_RESTART.AND.AUTOMATIC_RESTART&
428              .AND.ITER_RESTART.LE.10) GOTO 101
429     
430     ! Call user-defined subroutine after time-loop.
431           IF (CALL_USR) CALL USR3
432     
433     ! Get the final value of CPU time.  The difference gives the
434     ! CPU time used for the computations.
435           CALL CPU_TIME (CPU1)
436     
437     ! Compute the CPU time and write it out in the .OUT file.
438           CPUTIME_USED = CPU1 - CPU0 - CPU_IO
439           WALLTIME_USED = WALL_TIME() - WALL0
440           CALL WRITE_OUT3 (CPUTIME_USED, WALLTIME_USED, CPU_IO)
441     
442     ! JFD: cartesian grid implementation
443           IF(WRITE_DASHBOARD) THEN
444              IF(DT>=DT_MIN) THEN
445                 RUN_STATUS = 'Complete.'
446              ELSE
447                 RUN_STATUS = 'DT < DT_MIN.  Recovery not possible!'
448              ENDIF
449              CALL GET_TUNIT(CPUTIME_USED,TUNIT)
450              CALL UPDATE_DASHBOARD(0,CPUTIME_USED,TUNIT)
451           ENDIF
452           IF(CARTESIAN_GRID)  CALL CLOSE_CUT_CELL_FILES
453     
454     ! Finalize and terminate MPI
455           call parallel_fin
456     
457           CALL FINL_ERR_MSG
458     
459           STOP
460     
461      1000 FORMAT(/1X,'MFIX ',A,' Simulation:'/)
462     
463      1010 FORMAT('Message 1010: Read in data from .RES file for TIME = ',&
464              G12.5,/'Time step number (NSTEP) =',I7)
465     
466           END PROGRAM MFIX
467     
468     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv!
469     !                                                                      !
470     !  Subroutien: GEN_LOG_BASENAME                                        !
471     !  Author: Aytekin Gel                                Date: 19-SEP-03  !
472     !                                                                      !
473     !  Purpose: Generate the file base for DMP logs.                       !
474     !                                                                      !
475     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^!
476           SUBROUTINE GEN_LOG_BASENAME
477     
478           use compar, only: myPE
479           use compar, only: fbname
480     
481           implicit none
482     
483     ! Variables for generating file basename with processor id
484           INTEGER :: i1, i10, i100, i1000, i10000
485     
486     ! PAR_I/O Generate file basename for LOG files
487           i10000 = int(myPE/10000)
488           i1000  = int((myPE-i10000*10000)/1000)
489           i100   = int((myPE-i10000*10000-i1000*1000)/100)
490           i10    = int((myPE-i10000*10000-i1000*1000-i100*100)/10)
491           i1     = int((myPE-i10000*10000-i1000*1000-i100*100-i10*10)/1)
492     
493           i10000 = i10000 + 48
494           i1000  = i1000  + 48
495           i100   = i100   + 48
496           i10    = i10    + 48
497           i1     = i1     + 48
498     
499           fbname=char(i10000)//char(i1000)//char(i100)//char(i10)//char(i1)
500     
501           RETURN
502           END SUBROUTINE GEN_LOG_BASENAME
503     
504     
505     
506     
507     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
508     !                                                                      C
509     !  Module name: debug_write()                                          C
510     !  Purpose: Write out full geometry index setup information for the
511     !  case
512     !                                                                      C
513     !  Author: Aytekin Gel                                Date: 19-SEP-03  C
514     !  Reviewer:                                          Date:            C
515     !                                                                      C
516     !                                                                      C
517     !  Literature/Document References:                                     C
518     !                                                                      C
519     !  Variables referenced:                                               C
520     !  Variables modified:                                                 C
521     !                                                                      C
522     !  Local variables:                                                    C
523     !                                                                      C
524     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
525     
526           SUBROUTINE debug_write_layout()
527     
528     !-----------------------------------------------
529     ! Modules
530     !-----------------------------------------------
531           USE param
532           USE param1
533           USE parallel
534           USE matrix
535           USE geometry
536           USE compar
537           USE mpi_utility
538           USE sendrecv
539           USE sendrecv3
540           USE indices
541           USE leqsol
542           USE cdist
543           USE funits
544           USE run
545           USE time_cpu
546           USE functions
547           IMPLICIT NONE
548     !-----------------------------------------------
549     ! Local Variables
550     !-----------------------------------------------
551     ! phase index
552           INTEGER :: M
553     ! indices
554           INTEGER :: i, j, k, ijk, ijk_GL, ijk_PROC, ijk_IO
555     !
556           integer :: indxA, indxA_gl, indxB, indxB_gl, indxC, indxC_gl
557           integer :: indxD, indxD_gl, indxE, indxE_gl, indxF, indxF_gl
558           integer :: indxG, indxG_gl, indxH, indxH_gl
559     !
560           logical :: amgdbg = .TRUE.
561     
562           character(LEN=80) :: fname
563     
564     !DISTIO
565     !      fname = "layout_xxxx.txt"
566     !      write (fname(8:11),'(i4.4)') myPE
567           fname = "layout_xxxxx.txt"
568           write (fname(8:12),'(i5.5)') myPE
569           open (unit=11,file=fname,status='unknown')
570     
571           write (11,*) ' ********************************************'
572           write (11,*) ' ********************************************'
573           write (11,*) ' ********************************************'
574           write (11,*) ' ********************************************'
575           write (11,*) ' '
576           write (11,*) ' '
577           write (11,*) ' myPE =           ' , myPE
578           write (11,*) ' '
579           write (11,*) ' '
580     
581     
582         IF (AMGDBG .OR. bDist_IO) THEN
583           write(11,"('BLK1: Running from istart3,iend3 .AND. jstart3, jend3 .AND. kstart3, kend3')")
584           write(11,"(' (   i ,    j,     k) =>    ijk      ijk_GL     ijk_PROC    ijk_IO')")
585           write(11,"(' ====================      =====     =======    ========    ======')")
586           DO k = kstart3, kend3
587             DO i = istart3,iend3
588               DO j = jstart3, jend3
589                  ijk = FUNIJK(i,j,k)
590                  ijk_GL = FUNIJK_GL(i,j,k)
591                  ijk_PROC = FUNIJK_PROC(i,j,k,myPE)
592                  ijk_IO = FUNIJK_IO(i,j,k)
593                  write(11,"(' (',I4,' , ',I4,' , ',I4,') => ',4(I8,' , '))") &
594                                              i,j,k,ijk,ijk_GL,ijk_PROC,ijk_IO
595               ENDDO
596             ENDDO
597           ENDDO
598     
599           write(11,"(/,/,'BLK2: Print out Bottom, South, West, East, North, Top neighbors')")
600           write(11,"(' (   i ,    j,     k) =>    ijk    ijk_GL    B_of    S_of    W_of    E_of    N_of    T_of')")
601           write(11,"(' ====================      =====   =======  ======  ======  ======  ======  ======  ======')")
602           DO k = kstart3, kend3
603             DO i = istart3,iend3
604               DO j = jstart3, jend3
605                  ijk = FUNIJK(i,j,k)
606                  ijk_GL = FUNIJK_GL(i,j,k)
607                  write(11,"(' (',I4,' , ',I4,' , ',I4,') => ',2(I7,' , '),6(I7,2X))") &
608                                              i,j,k,ijk,ijk_GL,bottom_of(ijk),south_of(ijk),west_of(ijk),&
609                                              east_of(ijk),north_of(ijk),top_of(ijk)
610               ENDDO
611             ENDDO
612           ENDDO
613     
614           write(11,"(/,/,'BLK3: Print out km, jm, im, ip, jp, kp neighbors')")
615           write(11,"(' (   i ,    j,     k) =>    ijk    ijk_GL    km_of   jm_of   im_of   ip_of   jp_of   kp_of')")
616           write(11,"(' ====================      =====   =======  ======  ======  ======  ======  ======  ======')")
617           DO k = kstart3, kend3
618             DO i = istart3,iend3
619               DO j = jstart3, jend3
620                  ijk = FUNIJK(i,j,k)
621                  ijk_GL = FUNIJK_GL(i,j,k)
622                  write(11,"(' (',I4,' , ',I4,' , ',I4,') => ',2(I7,' , '),6(I7,2X))") &
623                                              i,j,k,ijk,ijk_GL,km_of(ijk),jm_of(ijk),im_of(ijk),&
624                                              ip_of(ijk),jp_of(ijk),kp_of(ijk)
625               ENDDO
626             ENDDO
627           ENDDO
628     
629           write(11,"(/,'BLK4a: Active Fluid Cells:FLUID_AT(ijk)=.T.',/,&
630          &           ' (   i ,    j,     k) =>    ijk  [   x ,     ,     z]')")
631           write(11,"(' ====================      =====  ====================')")
632            DO ijk = ijkstart3, ijkend3
633              I = I_OF(IJK)
634              J = J_OF(IJK)
635              K = K_OF(IJK)
636     
637     !         IF (FLOW_AT_E(IJK)) THEN
638              IF (FLUID_AT(IJK)) THEN
639     !          write(11,"(' (',I4,' , ',I4,' , ',I4,') => ',I8)") I,J,K,ijk
640                write(11,"(' (',I4,' , ',I4,' , ',I4,') => ',I8,' [',E12.5,',',E12.5,' ]')") I,J,K,ijk,X(i),Z(k)
641              ENDIF
642           ENDDO
643     
644           write(11,"(/,'BLK4b: Cells that are (.NOT.WALL_AT(IJK)) = .T.',/,&
645          &           ' (   i ,    j,     k) =>    ijk  [   x ,     ,     z]')")
646           write(11,"(' ====================      =====  ====================')")
647            DO ijk = ijkstart3, ijkend3
648              I = I_OF(IJK)
649              J = J_OF(IJK)
650              K = K_OF(IJK)
651     
652              IF (.NOT.WALL_AT(IJK)) THEN
653     !          write(11,"(' (',I4,' , ',I4,' , ',I4,') => ',I8)") I,J,K,ijk
654                write(11,"(' (',I4,' , ',I4,' , ',I4,') => ',I8,' [',E12.5,',',E12.5,' ]')") I,J,K,ijk,X(i),Z(k)
655              ENDIF
656           ENDDO
657     
658           DO k = kstart3, kend3
659             DO i = istart3,iend3
660               DO j = jstart3, jend3
661                  ijk = FUNIJK(i,j,k)
662                  ijk_GL = FUNIJK_GL(i,j,k)
663     
664                  if (i == istart2 .AND. j == jstart2) then
665                      indxA = ijk
666                      indxA_gl = ijk_GL
667                  endif
668                  if (i == istart1 .AND. j == jstart1) then
669                      indxE = ijk
670                      indxE_gl = ijk_GL
671                  endif
672                  if (i == istart2 .AND. j == jend2) then
673                      indxB = ijk
674                      indxB_gl = ijk_GL
675                  endif
676                  if (i == istart1 .AND. j == jend1) then
677                      indxF = ijk
678                      indxF_gl = ijk_GL
679                  endif
680                  if (i == iend1 .AND. j == jstart1) then
681                      indxH = ijk
682                      indxH_gl = ijk_GL
683                  endif
684                  if (i == iend2 .AND. j == jstart2) then
685                      indxD = ijk
686                      indxD_gl = ijk_GL
687                  endif
688                  if (i == iend1 .AND. j == jend1) then
689                      indxG = ijk
690                      indxG_gl = ijk_GL
691                  endif
692                  if (i == iend2 .AND. j == jend2) then
693                      indxC = ijk
694                      indxC_gl = ijk_GL
695                  endif
696               ENDDO
697             ENDDO
698             write(11,"('BLK5:')")
699             write(11,"(57('='))")
700             write(11,"('k= ',I5,/,57('='))") k
701             write(11,"('B= ',I5,' (',I7,')',20X,'C= ',I5,' (',I7,')',/)") indxB, indxB_gl, &
702                             indxC, indxC_gl
703     !        write(UNIT_LOG,"(' \',34X,'/')")
704     !        write(UNIT_LOG,"(2X,'\',32X,'/')")
705             write(11,"(3X,'F= ',I5,' (',I7,')',12X,'G= ',I5,' (',I7,')')") indxF, indxF_gl, &
706                             indxG, indxG_gl
707             write(11,"(4(9X,'|',29X,'|',/))")
708             write(11,"(3X,'E= ',I5,' (',I7,')',12X,'H= ',I5,' (',I7,')',/)") indxE, indxE_gl, &
709                             indxH, indxH_gl
710     !        write(UNIT_LOG,"(2X,'/',32X,'\')")
711     !        write(UNIT_LOG,"('/',34X,'\')")
712             write(11,"('A= ',I5,' (',I7,')',20X,'D= ',I5,' (',I7,')',/,/)") indxA, indxA_gl, &
713                             indxD, indxD_gl
714     
715     !        write(UNIT_LOG,"(' (',I4,' , ',I4,' , ',I4,') => ',2(I7,' , '),6(I7,2X))") &
716     !                                         i,j,k,ijk,ijk_GL,bottom_of(ijk),south_of(ijk),west_of(ijk),&
717     !                                        east_of(ijk),north_of(ijk),top_of(ijk)
718     
719           ENDDO
720     
721     !      write(UNIT_LOG,"(/,' (   i ,    j,     k) =>    ijk (Active Fluid)')")
722     !      write(UNIT_LOG,"(' ====================      =====')")
723     !       DO ijk = ijkstart3, ijkend3
724     !         I = I_OF(IJK)
725     !         J = J_OF(IJK)
726     !         K = K_OF(IJK)
727     
728     !         IF (FLOW_AT_E(IJK)) THEN
729     !         IF (FLUID_AT(IJK)) THEN
730     !           write(UNIT_LOG,"(' (',I4,' , ',I4,' , ',I4,') => ',I8)") I,J,K,ijk
731     !         ENDIF
732     !      END DO
733     
734     
735         endif   ! end if(amgdbg .or. bdist_io)
736     
737           M = 0
738     !      CALL WRITE_AB_M (A_M, B_M, IJKMAX2, M, IER)
739     
740         IF (AMGDBG .OR. bDist_IO) THEN
741           write(11,"(/,/,'BLK6: ========= ORIGINAL MFIX VARIABLES ===========')")
742           write(11,"('PE ',I5,': imin1  = ',I6,3X,'imax1= ',I6,/,'PE ',I5,': jmin1  = ',I6,3X,'jmax1= ',I6)") &
743                  myPE,imin1,imax1,myPE,jmin1,jmax1
744           write(11,"('PE ',I5,': kmin1  = ',I6,3X,'kmax1= ',I6)") myPE,kmin1,kmax1
745           write(11,"('-----')")
746           write(11,"('PE ',I5,': imin2  = ',I6,3X,'imax2= ',I6,/,'PE ',I5,': jmin2  = ',I6,3X,'jmax2= ',I6)") &
747                  myPE,imin2,imax2,myPE,jmin2,jmax2
748           write(11,"('PE ',I5,': kmin2  = ',I6,3X,'kmax2= ',I6)") myPE,kmin2,kmax2
749           write(11,"('----- Below xxx3 set is DMP extension ------------')")
750           write(11,"('PE ',I5,': imin3  = ',I6,3X,'imax3= ',I6,/,'PE ',I5,': jmin3  = ',I6,3X,'jmax3= ',I6)") &
751                  myPE,imin3,imax3,myPE,jmin3,jmax3
752           write(11,"('PE ',I5,': kmin3  = ',I6,3X,'kmax3= ',I6)") myPE,kmin3,kmax3
753           write(11,"('----- End of Below xxx3 set is DMP extension -----')")
754     !      write(11,"('PE ',I5,': ijkmax2= ',I6)") myPE,ijkmax2
755           write(11,"('PE ',I5,': ijmax2 = ',I6)") myPE,ijmax2
756           write(11,"('PE ',I5,': ijkmin1= ',I6,' ijkmax1= ',I12)") myPE,ijkmin1, ijkmax1
757           write(11,"('PE ',I5,':          ',6X,' ijkmax2= ',I12)") myPE,ijkmax2
758           write(11,"('PE ',I5,':          ',6X,' ijkmax3= ',I12)") myPE,ijkmax3
759           write(11,"('PE ',I5,': ijkmin4= ',I6,' ijkmax4= ',I12)") myPE,ijkmin4, ijkmax4
760     
761     
762           write(11,"(/,/,' ========= DMP EXTENSION VARIABLES ===========')")
763     !      write(UNIT_LOG,"('PE ',I5,': ijksize  = ',I6)") myPE,ijksize
764           write(11,"('PE ',I5,': ijksize3 = ',I6,3X,'ijksize3_all = ',I6)") myPE,ijksize3,ijksize3_all(myPE)
765           write(11,"('PE ',I5,': ijksize4 = ',I6,3X,'ijksize4_all = ',I6)") myPE,ijksize4,ijksize4_all(myPE)
766           write(11,"('PE ',I5,': ijkstart3  = ',I6,3X,'ijkend3  = ',I6)") myPE,ijkstart3, ijkend3
767           write(11,"('PE ',I5,': ijkstart3_all = ',I6,3X,'ijkstart4_all = ',I6)") myPE,ijkstart3_all(myPE),ijkstart4_all(myPE)
768           write(11,"('PE ',I5,': istart_all = ',I6,3X,'iend_all = ',I6,/,'PE ',I5,': jstart_all = ',I6,3X,'jend_all = ',I6)") &
769                  myPE,istart_all(myPE),iend_all(myPE),myPE,jstart_all(myPE),jend_all(myPE)
770           write(11,"('PE ',I5,': kstart_all = ',I6,3X,'kend_all = ',I6,/,'----------------------')") &
771                  myPE,kstart_all(myPE),kend_all(myPE)
772     
773           write(11,"('PE ',I5,': istart1_all= ',I6,3X,'iend1_all= ',I6,/,'PE ',I5,': jstart1_all= ',I6,3X,'jend3_all= ',I6)") &
774                  myPE,istart1_all(myPE),iend1_all(myPE),myPE,jstart1_all(myPE),jend1_all(myPE)
775           write(11,"('PE ',I5,': kstart1_all= ',I6,3X,'kend1_all= ',I6,/,'----------------------')") &
776                  myPE,kstart1_all(myPE),kend1_all(myPE)
777     
778           write(11,"('PE ',I5,': istart2_all= ',I6,3X,'iend2_all= ',I6,/,'PE ',I5,': jstart2_all= ',I6,3X,'jend3_all= ',I6)") &
779                  myPE,istart2_all(myPE),iend2_all(myPE),myPE,jstart2_all(myPE),jend2_all(myPE)
780           write(11,"('PE ',I5,': kstart2_all= ',I6,3X,'kend2_all= ',I6,/,'----------------------')") &
781                  myPE,kstart2_all(myPE),kend2_all(myPE)
782     
783           write(11,"('PE ',I5,': istart3_all= ',I6,3X,'iend3_all= ',I6,/,'PE ',I5,': jstart3_all= ',I6,3X,'jend3_all= ',I6)") &
784                  myPE,istart3_all(myPE),iend3_all(myPE),myPE,jstart3_all(myPE),jend3_all(myPE)
785           write(11,"('PE ',I5,': kstart3_all= ',I6,3X,'kend3_all= ',I6,/,'----------------------')") &
786                  myPE,kstart3_all(myPE),kend3_all(myPE)
787     
788           write(11,"('PE ',I5,': istart1= ',I6,3X,'iend1= ',I6,/,'PE ',I5,': jstart1= ',I6,3X,'jend1= ',I6)") &
789                  myPE,istart1,iend1,myPE,jstart1,jend1
790           write(11,"('PE ',I5,': kstart1= ',I6,3X,'kend1= ',I6,/,'----------------------')") &
791                  myPE,kstart1,kend1
792           write(11,"('PE ',I5,': istart2= ',I6,3X,'iend2= ',I6,/,'PE ',I5,': jstart2= ',I6,3X,'jend2= ',I6)") &
793                  myPE,istart2,iend2,myPE,jstart2,jend2
794           write(11,"('PE ',I5,': kstart2= ',I6,3X,'kend2= ',I6,/,'----------------------')") &
795                  myPE,kstart2,kend2
796           write(11,"('PE ',I5,': istart3= ',I6,3X,'iend3= ',I6,/,'PE ',I5,': jstart3= ',I6,3X,'jend3= ',I6)") &
797                  myPE,istart3,iend3,myPE,jstart3,jend3
798           write(11,"('PE ',I5,': kstart3= ',I6,3X,'kend3= ',I6,/,'----------------------')") &
799                  myPE,kstart3,kend3
800     
801         ENDIF   ! end if(amgdbg .or. bdist_io)
802     
803           close(unit=11)
804     
805     
806           RETURN
807           END SUBROUTINE DEBUG_WRITE_LAYOUT
808     
809     
810     
811     !vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvC
812     
813     !^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^C
814     
815           SUBROUTINE write_parallel_info()
816     
817     !-----------------------------------------------
818     !   M o d u l e s
819     !-----------------------------------------------
820           USE param
821           USE param1
822           USE parallel
823           USE matrix
824           USE geometry
825           USE compar
826           USE mpi_utility
827           USE sendrecv
828           USE sendrecv3
829           USE indices
830           USE leqsol
831           USE funits
832           USE run
833           USE time_cpu
834           USE functions
835           IMPLICIT NONE
836     !-----------------------------------------------
837     ! Dummy arguments
838     !-----------------------------------------------
839     ! Local Variables
840     !-----------------------------------------------
841     ! phase index
842           INTEGER :: M
843     ! indices
844           INTEGER :: i, j, k, ijk, ijk_GL, ijk_PROC, ijk_IO
845     !
846           character(LEN=80) :: fname
847     !-----------------------------------------------
848     
849     !DISTIO
850     !      fname = "p_info_xxxx.txt"
851     !      write (fname(8:11),'(i4.4)') myPE
852           fname = "p_info_xxxxx.txt"
853           write (fname(8:12),'(i5.5)') myPE
854           open (unit=11,file=fname,status='unknown')
855     
856           write (11,*) myPe , ' = myPE'
857     
858           write (11,*) myPE , istart3,iend3
859           write (11,*) myPE , jstart3,jend3
860           write (11,*) myPE , kstart3,kend3
861     
862           write(11,"('BLK1: Running from istart3,iend3 .AND. jstart3, jend3 .AND. kstart3, kend3')")
863           write(11,"(' (   i ,    j,     k)       ijk      ijk_GL     ijk_PROC    ijk_IO')")
864           write(11,"(' ====================      =====     =======    ========    ======')")
865           DO k = kstart3, kend3
866             DO i = istart3,iend3
867               DO j = jstart3, jend3
868                  ijk = FUNIJK(i,j,k)
869                  ijk_GL = FUNIJK_GL(i,j,k)
870                  ijk_PROC = FUNIJK_PROC(i,j,k,myPE)
871                  ijk_IO = FUNIJK_IO(i,j,k)
872                  write(11,"('  ',I4,'   ',I4,'   ',I4,'     ',4(I8,'   '))" ) &
873                                              i,j,k,ijk,ijk_GL,ijk_PROC,ijk_IO
874               ENDDO
875             ENDDO
876           ENDDO
877     
878           M = 0
879     !      CALL WRITE_AB_M (A_M, B_M, IJKMAX2, M, IER)
880     
881           close(unit=11)
882     
883           RETURN
884           END SUBROUTINE write_parallel_info
885     
886