From 167f9511ba1a19f3165055f0763046016d2d1c32 Mon Sep 17 00:00:00 2001 From: Nick Horelik Date: Thu, 6 Jun 2013 08:24:13 -0700 Subject: [PATCH 1/7] Added cell overlap checking subroutine --- src/geometry.F90 | 51 ++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 51 insertions(+) diff --git a/src/geometry.F90 b/src/geometry.F90 index 448691ca2d..84edc80930 100644 --- a/src/geometry.F90 +++ b/src/geometry.F90 @@ -62,6 +62,57 @@ contains end function simple_cell_contains + +!=============================================================================== +! CHECK_CELL_OVERLAP checks for overlapping cells at the current particle's +! position using simple_cell_contains and the LocalCoord's built up by +! find_cell. This should be called immediately following calls to +! find_cell. +!=============================================================================== + + subroutine check_cell_overlap() + + integer :: i ! cell loop index on a level + integer :: n ! number of cells to search on a level + integer :: index_cell ! index in cells array + type(Cell), pointer :: c ! pointer to cell + type(Universe), pointer :: univ ! universe to search in + type(LocalCoord), pointer :: coord ! particle coordinate to search on + + coord => p % coord0 + + ! loop through each coordinate level + do while (associated(coord)) + + p % coord => coord + + univ => universes(coord % universe) + n = univ % n_cells + + ! loop through each cell on this level + do i = 1, n + index_cell = univ % cells(i) + c => cells(index_cell) + + if (simple_cell_contains(c)) then + ! that particle should only be contained in one cell per level + if (index_cell /= coord % cell) then + message = "Overlapping cells detected: " // & + trim(to_str(cells(index_cell) % id)) // ", " // & + trim(to_str(cells(coord % cell) % id)) // & + " on universe " // trim(to_str(univ % id)) + call fatal_error() + end if + end if + + end do + + coord => coord % next + + end do + + end subroutine check_cell_overlap + !=============================================================================== ! FIND_CELL determines what cell a source particle is in within a particular ! universe. If the base universe is passed, the particle should be found as long From 6aeb5281b0758bcbd39eda4d1af33a10c3c8ec8b Mon Sep 17 00:00:00 2001 From: Nick Horelik Date: Thu, 6 Jun 2013 08:25:49 -0700 Subject: [PATCH 2/7] Made plots check for cell overlaps by default --- src/plot.F90 | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/plot.F90 b/src/plot.F90 index c38788402a..5158718573 100644 --- a/src/plot.F90 +++ b/src/plot.F90 @@ -2,8 +2,8 @@ module plot use constants use error, only: fatal_error - use geometry, only: find_cell - use geometry_header, only: Universe, BASE_UNIVERSE + use geometry, only: find_cell, check_cell_overlap + use geometry_header, only: Cell, BASE_UNIVERSE use global use output, only: write_message use particle_header, only: deallocate_coord @@ -62,6 +62,7 @@ contains p % coord => p % coord0 call find_cell(found_cell) + call check_cell_overlap() if (.not. found_cell) then ! If no cell, revert to default color From d9165b5a6a23a6383c8bbc1d788bd341b564f25f Mon Sep 17 00:00:00 2001 From: Nick Horelik Date: Thu, 6 Jun 2013 09:08:56 -0700 Subject: [PATCH 3/7] Added command-line flag to enable cell overlap checking during transport --- src/global.F90 | 3 +++ src/initialize.F90 | 3 +++ src/physics.F90 | 5 ++++- 3 files changed, 10 insertions(+), 1 deletion(-) diff --git a/src/global.F90 b/src/global.F90 index e0f0bd5bbf..f2bbff6019 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -274,6 +274,9 @@ module global ! screen and in logs integer :: verbosity = 7 + ! Flag for enabling cell overlap checking during transport + logical :: check_overlaps = .false. + ! Trace for single particle logical :: trace integer :: trace_batch diff --git a/src/initialize.F90 b/src/initialize.F90 index dd28f9bf9c..8cd85a5e90 100644 --- a/src/initialize.F90 +++ b/src/initialize.F90 @@ -297,6 +297,9 @@ contains path_state_point = argv(i) restart_run = .true. + case ('-g', '-geometry-debug', '--geometry-debug') + check_overlaps = .true. + case ('-?', '-help', '--help') call print_usage() stop diff --git a/src/physics.F90 b/src/physics.F90 index 707049b1a8..4dc5e13b29 100644 --- a/src/physics.F90 +++ b/src/physics.F90 @@ -7,7 +7,8 @@ module physics use error, only: fatal_error, warning use fission, only: nu_total, nu_delayed use geometry, only: find_cell, distance_to_boundary, & - cross_surface, cross_lattice + cross_surface, cross_lattice, & + check_cell_overlap use geometry_header, only: Universe, BASE_UNIVERSE use global use interpolation, only: interpolate_tab1 @@ -75,6 +76,8 @@ contains do while (p % alive) + if (check_overlaps) call check_cell_overlap() + ! Calculate microscopic and macroscopic cross sections -- note: if the ! material is the same as the last material and the energy of the ! particle hasn't changed, we don't need to lookup cross sections again. From f96e952c09912f96e3cf828609ae539178f6839f Mon Sep 17 00:00:00 2001 From: Nick Horelik Date: Thu, 6 Jun 2013 11:49:12 -0700 Subject: [PATCH 4/7] Added output summary for cell overlap checking --- src/finalize.F90 | 4 +++- src/geometry.F90 | 5 ++++- src/global.F90 | 6 +++++- src/initialize.F90 | 17 ++++++++++++++--- src/input_xml.F90 | 5 +++++ src/output.F90 | 36 ++++++++++++++++++++++++++++++++++++ src/plot.F90 | 2 +- 7 files changed, 68 insertions(+), 7 deletions(-) diff --git a/src/finalize.F90 b/src/finalize.F90 index d601c2e592..da9d80e2fd 100644 --- a/src/finalize.F90 +++ b/src/finalize.F90 @@ -4,7 +4,8 @@ module finalize use cmfd_output, only: finalize_cmfd # endif use global - use output, only: print_runtime, print_results, write_tallies + use output, only: print_runtime, print_results, & + print_overlap_check, write_tallies use tally, only: tally_statistics #ifdef MPI @@ -50,6 +51,7 @@ contains run_mode /= MODE_PARTICLE)) then call print_runtime() call print_results() + if (check_overlaps) call print_overlap_check() end if ! deallocate arrays diff --git a/src/geometry.F90 b/src/geometry.F90 index 84edc80930..52bd0a43cf 100644 --- a/src/geometry.F90 +++ b/src/geometry.F90 @@ -95,7 +95,7 @@ contains c => cells(index_cell) if (simple_cell_contains(c)) then - ! that particle should only be contained in one cell per level + ! the particle should only be contained in one cell per level if (index_cell /= coord % cell) then message = "Overlapping cells detected: " // & trim(to_str(cells(index_cell) % id)) // ", " // & @@ -103,6 +103,9 @@ contains " on universe " // trim(to_str(univ % id)) call fatal_error() end if + + overlap_check_cnt(index_cell) = overlap_check_cnt(index_cell) + 1 + end if end do diff --git a/src/global.F90 b/src/global.F90 index f2bbff6019..addfaadbed 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -275,7 +275,8 @@ module global integer :: verbosity = 7 ! Flag for enabling cell overlap checking during transport - logical :: check_overlaps = .false. + logical :: check_overlaps = .false. + integer, allocatable :: overlap_check_cnt(:) ! Trace for single particle logical :: trace @@ -399,6 +400,9 @@ contains if (allocated(materials)) deallocate(materials) if (allocated(plots)) deallocate(plots) + ! Deallocate geometry debugging information + if (allocated(overlap_check_cnt)) deallocate(overlap_check_cnt) + ! Deallocate cross section data, listings, and cache if (allocated(nuclides)) then ! First call the clear routines diff --git a/src/initialize.F90 b/src/initialize.F90 index 8cd85a5e90..4870123bb4 100644 --- a/src/initialize.F90 +++ b/src/initialize.F90 @@ -5,14 +5,15 @@ module initialize use constants use dict_header, only: DictIntInt, ElemKeyValueII use energy_grid, only: unionized_grid - use error, only: fatal_error + use error, only: fatal_error, warning use geometry, only: neighbor_lists use geometry_header, only: Cell, Universe, Lattice, BASE_UNIVERSE use global use input_xml, only: read_input_xml, read_cross_sections_xml, & cells_in_univ_dict, read_plots_xml use output, only: title, header, write_summary, print_version, & - print_usage, write_xs_summary, print_plot + print_usage, write_xs_summary, print_plot, & + write_message use random_lcg, only: initialize_prng use source, only: initialize_source use state_point, only: load_state_point @@ -145,9 +146,17 @@ contains end if end if - ! check for particle restart run + ! Check for particle restart run if (particle_restart_run) run_mode = MODE_PARTICLE + ! Warn if overlap checking is on + if (check_overlaps) then + message = "" + call write_message() + message = "Cell overlap checking is ON" + call warning() + end if + ! Stop initialization timer call time_initialize % stop() @@ -272,6 +281,8 @@ contains select case (argv(i)) case ('-p', '-plot', '--plot') run_mode = MODE_PLOTTING + check_overlaps = .true. + case ('-n', '-n_particles', '--n_particles') ! Read number of particles per cycle i = i + 1 diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 75d0518d1f..5c00ed8d0c 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -626,6 +626,11 @@ contains ! Allocate cells array allocate(cells(n_cells)) + if (check_overlaps) then + allocate(overlap_check_cnt(n_cells)) + overlap_check_cnt = 0 + end if + n_universes = 0 do i = 1, n_cells c => cells(i) diff --git a/src/output.F90 b/src/output.F90 index 4ca455624b..fc5d1d6a69 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -1472,6 +1472,42 @@ contains end subroutine print_results +!=============================================================================== +! PRINT_OVERLAP_DEBUG displays information regarding overlap checking results +!=============================================================================== + + subroutine print_overlap_check + + integer :: i, j + integer :: num_sparse = 0 + + ! display header block for geometry debugging section + call header("Cell Overlap Check Summary") + + write(ou,100) 'Cell ID','No. Overlap Checks' + + do i = 1, n_cells + write(ou,101) cells(i) % id, overlap_check_cnt(i) + if (overlap_check_cnt(i) < 10) num_sparse = num_sparse + 1 + end do + write(ou,*) + write(ou,'(1X,A)') 'There were ' // trim(to_str(num_sparse)) // & + ' cells with less than 10 overlap checks' + j = 0 + do i = 1, n_cells + if (overlap_check_cnt(i) < 10) then + j = j + 1 + write(ou,'(1X,A8)', advance='no') trim(to_str(cells(i) % id)) + if (modulo(j,8) == 0) write(ou,*) + end if + end do + write(ou,*) + +100 format (1X,A,T15,A) +101 format (1X,I8,T15,I8) + + end subroutine print_overlap_check + !=============================================================================== ! WRITE_TALLIES creates an output file and writes out the mean values of all ! tallies and their standard deviations diff --git a/src/plot.F90 b/src/plot.F90 index 5158718573..e7b8a235a6 100644 --- a/src/plot.F90 +++ b/src/plot.F90 @@ -62,7 +62,7 @@ contains p % coord => p % coord0 call find_cell(found_cell) - call check_cell_overlap() + if (check_overlaps) call check_cell_overlap() if (.not. found_cell) then ! If no cell, revert to default color From 00286f1e49612ff3e486432fb7b0bf58c3b1b2a8 Mon Sep 17 00:00:00 2001 From: Nick Horelik Date: Thu, 6 Jun 2013 13:04:05 -0700 Subject: [PATCH 5/7] Added support for parallel runs to geometry overlap checking --- src/finalize.F90 | 19 +++++++++++++++++++ src/global.F90 | 4 ++-- src/initialize.F90 | 2 +- src/output.F90 | 2 +- 4 files changed, 23 insertions(+), 4 deletions(-) diff --git a/src/finalize.F90 b/src/finalize.F90 index da9d80e2fd..0236eb4e84 100644 --- a/src/finalize.F90 +++ b/src/finalize.F90 @@ -36,6 +36,7 @@ contains if (output_tallies) then if (master) call write_tallies() end if + if (check_overlaps) call reduce_overlap_count() end if #ifdef PETSC @@ -69,4 +70,22 @@ contains end subroutine finalize_run +!=============================================================================== +! REDUCE_OVERLAP_COUNT accumulates cell overlap check counts to master +!=============================================================================== + + subroutine reduce_overlap_count() + +#ifdef MPI + if (master) then + call MPI_REDUCE(MPI_IN_PLACE, overlap_check_cnt, n_cells, & + MPI_INTEGER8, MPI_SUM, 0, MPI_COMM_WORLD, mpi_err) + else + call MPI_REDUCE(overlap_check_cnt, overlap_check_cnt, n_cells, & + MPI_INTEGER8, MPI_SUM, 0, MPI_COMM_WORLD, mpi_err) + end if +#endif + + end subroutine reduce_overlap_count + end module finalize diff --git a/src/global.F90 b/src/global.F90 index addfaadbed..502968054a 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -275,8 +275,8 @@ module global integer :: verbosity = 7 ! Flag for enabling cell overlap checking during transport - logical :: check_overlaps = .false. - integer, allocatable :: overlap_check_cnt(:) + logical :: check_overlaps = .false. + integer(8), allocatable :: overlap_check_cnt(:) ! Trace for single particle logical :: trace diff --git a/src/initialize.F90 b/src/initialize.F90 index 4870123bb4..c9c572ddda 100644 --- a/src/initialize.F90 +++ b/src/initialize.F90 @@ -150,7 +150,7 @@ contains if (particle_restart_run) run_mode = MODE_PARTICLE ! Warn if overlap checking is on - if (check_overlaps) then + if (master .and. check_overlaps) then message = "" call write_message() message = "Cell overlap checking is ON" diff --git a/src/output.F90 b/src/output.F90 index fc5d1d6a69..727747dff5 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -1504,7 +1504,7 @@ contains write(ou,*) 100 format (1X,A,T15,A) -101 format (1X,I8,T15,I8) +101 format (1X,I8,T15,I12) end subroutine print_overlap_check From d88b2e504016fda7219ba9d5727d3614be17986d Mon Sep 17 00:00:00 2001 From: nhorelik Date: Thu, 6 Jun 2013 15:52:08 -0400 Subject: [PATCH 6/7] Updated documentation with overlap checking --- docs/source/usersguide/processing.rst | 2 +- docs/source/usersguide/troubleshoot.rst | 37 ++++++++++++++++++++++++- man/man1/openmc.1 | 4 +++ src/geometry.F90 | 4 +-- 4 files changed, 42 insertions(+), 5 deletions(-) diff --git a/docs/source/usersguide/processing.rst b/docs/source/usersguide/processing.rst index 874290465d..7105e9f7d1 100644 --- a/docs/source/usersguide/processing.rst +++ b/docs/source/usersguide/processing.rst @@ -288,7 +288,7 @@ SILO file can be created with: import silomesh as sm sm.init_silo("fluxtally.silo") - sm.init_mesh('tally_mesh',*mesh.dimension, *mesh.lower_left, *mesh.upper_right) + sm.init_mesh('tally_mesh', *mesh.dimension, *mesh.lower_left, *mesh.upper_right) sm.init_var('flux_tally_thermal') for x in range(1,nx+1): for y in range(1,ny+1): diff --git a/docs/source/usersguide/troubleshoot.rst b/docs/source/usersguide/troubleshoot.rst index 64b955593e..49cf41c8de 100644 --- a/docs/source/usersguide/troubleshoot.rst +++ b/docs/source/usersguide/troubleshoot.rst @@ -98,8 +98,38 @@ with the :envvar:`CROSS_SECTIONS` environment variable. It is recommended to add a line in your ``.profile`` or ``.bash_profile`` setting the :envvar:`CROSS_SECTIONS` environment variable. +Geometry Debugging +****************** + +Overlapping Cells +^^^^^^^^^^^^^^^^^ + +For fast run times, normal simulations do not check if the geometry is +incorrectly defined to have overlapping cells. This can lead to incorrect +results that may or may not be obvious when there are errors in the geometry +input file. The built-in 2D and 3D plotters will check for cell overlaps at +the center of every pixel or voxel position they process, however this might +not be a sufficient check to ensure correctly defined geometry. For instance, +if an overlap is of small aspect ratio, the plotting resolution might not be +high enough to produce any pixels in the overlapping area. + +To reliably validate a geometry input, it is best to run the problem in +geometry debugging mode with the ``-g``, ``-geometry-debug``, or +``--geometry-debug`` command-line options. This will enable checks for +overlapping cells at every move of esch simulated particle. Depending on the +complexity of the geometry input file, this could add considerable overhead to +the run (these runs can still be done in parallel). As a result, for this run +mode the user will probably want to run fewer particles than a normal +simulation run. In this case it is important to be aware of how much coverage +each area of the geometry is getting. For instance, if certain regions do not +have many particles travelling through them there will not be many locations +where overlaps are checked for in that region. The user should refer to the +output after a geometry debug run to see how many checks were performed in each +cell, and then adjust the number of starting particles or starting source +distributions accordingly to achieve good coverage. + ERROR: After particle __ crossed surface __ it could not be located in any cell and it did not leak. -**************************************************************************************************** +^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ This error can arise either if a problem is specified with no boundary conditions or if there is an error in the geometry itself. First check to ensure @@ -120,5 +150,10 @@ has a collision. For example, if you received this error at cycle 5, generation 5 1 4032 +For large runs it is often advantageous to run only the offending particle by +using particle restart mode with the ``-s``, ``-particle``, or ``--particle`` +command-line options in conjunction with the particle restart files that are +created when particles are lost with this error. + .. _gfortran: http://gcc.gnu.org/wiki/GFortran .. _mailing list: https://groups.google.com/forum/?fromgroups=#!forum/openmc-users diff --git a/man/man1/openmc.1 b/man/man1/openmc.1 index b8a9d9c44c..aefbca9171 100644 --- a/man/man1/openmc.1 +++ b/man/man1/openmc.1 @@ -13,6 +13,10 @@ It is assumed that if no is specified, the XML input files are present in the current directory. .SH OPTIONS .TP +.B "\-g\fR, \fP\-\-geometry-debug" +Run in geometry debugging mode, where cell overlaps are checked for after each +move of a particle +.TP .B "\-p\fR, \fP\-\-plot" Run in plotting mode .TP diff --git a/src/geometry.F90 b/src/geometry.F90 index 52bd0a43cf..cb90997846 100644 --- a/src/geometry.F90 +++ b/src/geometry.F90 @@ -65,9 +65,7 @@ contains !=============================================================================== ! CHECK_CELL_OVERLAP checks for overlapping cells at the current particle's -! position using simple_cell_contains and the LocalCoord's built up by -! find_cell. This should be called immediately following calls to -! find_cell. +! position using simple_cell_contains and the LocalCoord's built up by find_cell !=============================================================================== subroutine check_cell_overlap() From c03da91153b677ce80a21533898b34280c9eaf2b Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sat, 15 Jun 2013 13:39:45 -0400 Subject: [PATCH 7/7] Added -g to print_usage. --- src/output.F90 | 1 + 1 file changed, 1 insertion(+) diff --git a/src/output.F90 b/src/output.F90 index 727747dff5..698338d091 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -165,6 +165,7 @@ contains write(OUTPUT_UNIT,*) 'Usage: openmc [options] [directory]' write(OUTPUT_UNIT,*) write(OUTPUT_UNIT,*) 'Options:' + write(OUTPUT_UNIT,*) ' -g, --geometry-debug Run in geometry debugging mode' write(OUTPUT_UNIT,*) ' -p, --plot Run in plotting mode' write(OUTPUT_UNIT,*) ' -r, --restart Restart a previous run' write(OUTPUT_UNIT,*) ' -s, --particle Run a single particle history'