From df5cbb325af88ea827f6de860629cddcee6fdb5b Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Wed, 23 Jul 2014 23:48:28 -0400 Subject: [PATCH 01/28] Add new energy grid search based on an equal log-spaced mapping technique as outlined in LA-UR-14-24530. --- src/ace_header.F90 | 2 +- src/constants.F90 | 7 +++--- src/cross_section.F90 | 31 +++++++++++++++++++++----- src/energy_grid.F90 | 51 ++++++++++++++++++++++++++++++++++++++++++- src/global.F90 | 1 + src/initialize.F90 | 4 +++- src/input_xml.F90 | 5 ++--- 7 files changed, 87 insertions(+), 14 deletions(-) diff --git a/src/ace_header.F90 b/src/ace_header.F90 index 27cf64715..11cb109a9 100644 --- a/src/ace_header.F90 +++ b/src/ace_header.F90 @@ -97,7 +97,7 @@ module ace_header ! Energy grid information integer :: n_grid ! # of nuclide grid points - integer, allocatable :: grid_index(:) ! pointers to union grid + integer, allocatable :: grid_index(:) ! union grid pointers / log grid mapping real(8), allocatable :: energy(:) ! energy values corresponding to xs ! Microscopic cross sections diff --git a/src/constants.F90 b/src/constants.F90 index 94edd9d2a..7bee73b1a 100644 --- a/src/constants.F90 +++ b/src/constants.F90 @@ -382,9 +382,10 @@ module constants ! Energy grid methods integer, parameter :: & - GRID_NUCLIDE = 1, & ! non-unionized energy grid - GRID_UNION = 2, & ! union grid with pointers - GRID_LETHARGY = 3 ! lethargy mapping + GRID_NUCLIDE = 1, & ! non-unionized energy grid + GRID_UNION = 2, & ! union grid with pointers + GRID_LOGARITHM = 3 ! logarithmic mapping + integer, parameter :: N_LOG_BINS = 8000 ! Running modes integer, parameter :: & diff --git a/src/cross_section.F90 b/src/cross_section.F90 index 7426519cd..a15067072 100644 --- a/src/cross_section.F90 +++ b/src/cross_section.F90 @@ -141,8 +141,10 @@ contains integer, intent(in) :: i_sab ! index into sab_tables array real(8), intent(in) :: E ! energy - integer :: i_grid ! index on nuclide energy grid - real(8) :: f ! interp factor on nuclide energy grid + integer :: i_grid ! index on nuclide energy grid + integer :: i_low, i_high ! bounding indices from logarithmic mapping + integer :: u ! index into logarithmic mapping array + real(8) :: f ! interp factor on nuclide energy grid type(Nuclide), pointer, save :: nuc => null() !$omp threadprivate(nuc) @@ -157,10 +159,29 @@ contains i_grid = nuc % grid_index(union_grid_index) + case (GRID_LOGARITHM) + ! Determine the energy grid index using a logarithmic mapping to reduce + ! the energy range over which a binary search needs to be performed + + if (E < nuc % energy(1)) then + i_grid = 1 + elseif (E > nuc % energy(nuc % n_grid)) then + i_grid = nuc % n_grid - 1 + else + ! Determine bounding indices based on which equal log-spaced interval + ! the energy is in + u = int(log(E/1.0e-11_8)/log_spacing) + i_low = nuc % grid_index(u) + i_high = nuc % grid_index(u + 1) + 1 + + ! Perform binary search over reduced range + i_grid = binary_search(nuc % energy(i_low:i_high), & + i_high - i_low + 1, E) + i_low - 1 + end if + case (GRID_NUCLIDE) - ! If we're not using the unionized grid, we have to do a binary search on - ! the nuclide energy grid in order to determine which points to - ! interpolate between + ! Perform binary search on the nuclide energy grid in order to determine + ! which points to interpolate between if (E < nuc % energy(1)) then i_grid = 1 diff --git a/src/energy_grid.F90 b/src/energy_grid.F90 index 66b95be0d..934cf2753 100644 --- a/src/energy_grid.F90 +++ b/src/energy_grid.F90 @@ -1,6 +1,6 @@ module energy_grid - use constants, only: MAX_LINE_LEN + use constants, only: MAX_LINE_LEN, N_LOG_BINS use global use list_header, only: ListReal use output, only: write_message @@ -148,4 +148,53 @@ contains end subroutine grid_pointers +!=============================================================================== +! LOGARITHMIC_GRID determines a logarithmic mapping for energies to bounding +! indices on a nuclide energy grid +!=============================================================================== + + subroutine logarithmic_grid() + + integer :: i, j, k ! Loop indices + integer :: M ! Number of equally log-spaced bins + real(8) :: E_max ! Maximum energy in MeV + real(8) :: E_min ! Minimum energy in MeV + real(8), allocatable :: umesh(:) ! Equally log-spaced energy grid + type(Nuclide), pointer :: nuc => null() + + ! Set minimum/maximum energies + E_max = 20.0_8 + E_min = 1.0e-11_8 + + ! Determine equal-logarithmic energy spacing + M = N_LOG_BINS + log_spacing = log(E_max/E_min)/M + + ! Create equally log-spaced energy grid + allocate(umesh(0:M)) + umesh(:) = [(i*log_spacing, i=0, M)] + + do i = 1, n_nuclides_total + ! Allocate logarithmic mapping for nuclide + nuc => nuclides(i) + allocate(nuc % grid_index(0:M)) + + ! Determine corresponding indices in nuclide grid to energies on + ! equal-logarithmic grid + j = 1 + do k = 0, M - 1 + do while (log(nuc%energy(j + 1)/E_min) <= umesh(k)) + j = j + 1 + end do + nuc % grid_index(k) = j + end do + + ! Set the last point explicitly so that we don't have out-of-bounds issues + nuc % grid_index(M) = size(nuc % energy) - 1 + end do + + deallocate(umesh) + + end subroutine logarithmic_grid + end module energy_grid diff --git a/src/global.F90 b/src/global.F90 index 7d728a25d..33674eb2b 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -80,6 +80,7 @@ module global ! Unionized energy grid integer :: grid_method ! how to treat the energy grid integer :: n_grid ! number of points on unionized grid + real(8) :: log_spacing ! spacing on logarithmic grid real(8), allocatable :: e_grid(:) ! energies on unionized grid ! Unreoslved resonance probablity tables diff --git a/src/initialize.F90 b/src/initialize.F90 index c14a8a02a..19586a5e5 100644 --- a/src/initialize.F90 +++ b/src/initialize.F90 @@ -4,7 +4,7 @@ module initialize use bank_header, only: Bank use constants use dict_header, only: DictIntInt, ElemKeyValueII - use energy_grid, only: unionized_grid + use energy_grid, only: unionized_grid, logarithmic_grid use error, only: fatal_error, warning use geometry, only: neighbor_lists use geometry_header, only: Cell, Universe, Lattice, BASE_UNIVERSE @@ -109,6 +109,8 @@ contains call time_unionize % start() call unionized_grid() call time_unionize % stop() + elseif (grid_method == GRID_LOGARITHM) then + call logarithmic_grid() end if ! Allocate and setup tally stride, matching_bins, and tally maps diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 46461b78d..45b6067cb 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -220,9 +220,8 @@ contains grid_method = GRID_NUCLIDE case ('union') grid_method = GRID_UNION - case ('lethargy') - message = "Lethargy mapped energy grid not yet supported." - call fatal_error() + case ('logarithm', 'logarithmic', 'log') + grid_method = GRID_LOGARITHM case default message = "Unknown energy grid method: " // trim(temp_str) call fatal_error() From aeffdb7ce69ae7d28f41c07d0f64a8efcb6f3e66 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Wed, 23 Jul 2014 23:59:12 -0400 Subject: [PATCH 02/28] Change logarithmic energy grid to default. --- src/input_xml.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 45b6067cb..3dda77b3c 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -213,7 +213,7 @@ contains if (check_for_node(doc, "energy_grid")) then call get_node_value(doc, "energy_grid", temp_str) else - temp_str = 'union' + temp_str = 'logarithm' end if select case (trim(temp_str)) case ('nuclide') From 99695b514eedc5ee9745c8a46a733268197488ac Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Wed, 23 Jul 2014 23:59:30 -0400 Subject: [PATCH 03/28] Update documentation for . --- docs/source/usersguide/input.rst | 16 ++++++++++------ src/relaxng/settings.rnc | 3 ++- 2 files changed, 12 insertions(+), 7 deletions(-) diff --git a/docs/source/usersguide/input.rst b/docs/source/usersguide/input.rst index c9aa8cc86..3bf31dcc3 100644 --- a/docs/source/usersguide/input.rst +++ b/docs/source/usersguide/input.rst @@ -134,13 +134,17 @@ should be performed. It has the following attributes/sub-elements: ------------------------- The ```` element determines the treatment of the energy grid during -a simulation. Setting this element to "nuclide" will cause OpenMC to use a -nuclide's energy grid when determining what points to interpolate between for -determining cross sections (i.e. non-unionized energy grid). To use a unionized -energy grid, set this element to "union". Note that the unionized energy grid -treatment is slightly different than that employed in Serpent. +a simulation. The valid options are "nuclide", "union", and "logarithm". Setting +this element to "nuclide" will cause OpenMC to use a nuclide's energy grid when +determining what points to interpolate between for determining cross sections +(i.e. non-unionized energy grid). Setting this element to "union" results in a +unionized energy grid with pointers to nuclide grids. Setting this element to +"logarithm" causes OpenMC to use a logarithmic mapping technique described in +LA-UR-14-24530_. - *Default*: union + *Default*: logarithm + +.. _LA-UR-14-24530: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-ur-14-24530.pdf ```` Element --------------------- diff --git a/src/relaxng/settings.rnc b/src/relaxng/settings.rnc index ea3ecabf5..b083c5692 100644 --- a/src/relaxng/settings.rnc +++ b/src/relaxng/settings.rnc @@ -27,7 +27,8 @@ element settings { (element weight_avg { xsd:double } | attribute weight_avg { xsd:double })? }? & - element energy_grid { ( "nuclide" | "union" | "lethargy" ) }? & + element energy_grid { ( "nuclide" | "union" | "log" | + "logarithm" | "logarithmic" ) }? & element entropy { (element dimension { list { xsd:int+ } } | From 7dee363bc2aed6e1d9d5f9340c04eb0ebb021b2c Mon Sep 17 00:00:00 2001 From: walshjon Date: Sat, 4 Oct 2014 17:05:24 -0700 Subject: [PATCH 04/28] changed function title blocks to reflect function names --- src/string.F90 | 48 ++++++++++++++++++++++++------------------------ 1 file changed, 24 insertions(+), 24 deletions(-) diff --git a/src/string.F90 b/src/string.F90 index d4f8873a0..ea7de1cc9 100644 --- a/src/string.F90 +++ b/src/string.F90 @@ -149,7 +149,7 @@ contains end function concatenate !=============================================================================== -! LOWER_CASE converts a string to all lower case characters +! TO_LOWER converts a string to all lower case characters !=============================================================================== elemental function to_lower(word) result(word_lower) @@ -172,7 +172,7 @@ contains end function to_lower !=============================================================================== -! UPPER_CASE converts a string to all upper case characters +! TO_UPPER converts a string to all upper case characters !=============================================================================== elemental function to_upper(word) result(word_upper) @@ -200,32 +200,32 @@ contains ! integers. !=============================================================================== -function zero_padded(num, n_digits) result(str) - integer, intent(in) :: num - integer, intent(in) :: n_digits - character(11) :: str + function zero_padded(num, n_digits) result(str) + integer, intent(in) :: num + integer, intent(in) :: n_digits + character(11) :: str - character(8) :: zp_form + character(8) :: zp_form - ! Make sure n_digits is reasonable. 10 digits is the maximum needed for the - ! largest integer(4). - if (n_digits > 10) then - message = 'zero_padded called with an unreasonably large n_digits (>10)' - call fatal_error() - end if + ! Make sure n_digits is reasonable. 10 digits is the maximum needed for the + ! largest integer(4). + if (n_digits > 10) then + message = 'zero_padded called with an unreasonably large n_digits (>10)' + call fatal_error() + end if - ! Write a format string of the form '(In.m)' where n is the max width and - ! m is the min width. If a sign is present, then n must be one greater - ! than m. - if (num < 0) then - write(zp_form, '("(I", I0, ".", I0, ")")') n_digits+1, n_digits - else - write(zp_form, '("(I", I0, ".", I0, ")")') n_digits, n_digits - end if + ! Write a format string of the form '(In.m)' where n is the max width and + ! m is the min width. If a sign is present, then n must be one greater + ! than m. + if (num < 0) then + write(zp_form, '("(I", I0, ".", I0, ")")') n_digits+1, n_digits + else + write(zp_form, '("(I", I0, ".", I0, ")")') n_digits, n_digits + end if - ! Format the number. - write(str, zp_form) num -end function zero_padded + ! Format the number. + write(str, zp_form) num + end function zero_padded !=============================================================================== ! IS_NUMBER determines whether a string of characters is all 0-9 characters From 1654659d587fb333a0439b827371d00460ad5a6b Mon Sep 17 00:00:00 2001 From: walshjon Date: Sat, 4 Oct 2014 17:09:13 -0700 Subject: [PATCH 05/28] changed comment to reflect what code is actually doing --- src/cross_section.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/cross_section.F90 b/src/cross_section.F90 index 8e4c19ccc..177956529 100644 --- a/src/cross_section.F90 +++ b/src/cross_section.F90 @@ -198,7 +198,7 @@ contains micro_xs(i_nuclide) % total = (ONE - f) * nuc % total(i_grid) & + f * nuc % total(i_grid+1) - ! Calculate microscopic nuclide total cross section + ! Calculate microscopic nuclide elastic cross section micro_xs(i_nuclide) % elastic = (ONE - f) * nuc % elastic(i_grid) & + f * nuc % elastic(i_grid+1) From c5ead96c164a8e8b7201957549e24654e4a4fb81 Mon Sep 17 00:00:00 2001 From: walshjon Date: Sun, 5 Oct 2014 17:56:39 -0700 Subject: [PATCH 06/28] fixed typo in comment --- src/physics.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/physics.F90 b/src/physics.F90 index 04cc66afa..5ff86deb7 100644 --- a/src/physics.F90 +++ b/src/physics.F90 @@ -234,7 +234,7 @@ contains prob = prob + ((ONE - f)*rxn%sigma(i_grid - rxn%threshold + 1) & + f*(rxn%sigma(i_grid - rxn%threshold + 2))) - ! Create fission bank sites if fission occus + ! Create fission bank sites if fission occurs if (prob > cutoff) exit FISSION_REACTION_LOOP end do FISSION_REACTION_LOOP From 0c82e1f46ce4c55face0ae289baa97ff0be6c95f Mon Sep 17 00:00:00 2001 From: walshjon Date: Sun, 5 Oct 2014 18:41:49 -0700 Subject: [PATCH 07/28] fixed typo in comment --- src/physics.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/physics.F90 b/src/physics.F90 index 5ff86deb7..f339e41e0 100644 --- a/src/physics.F90 +++ b/src/physics.F90 @@ -392,7 +392,7 @@ contains + f*(rxn%sigma(i_grid - rxn%threshold + 2))) end do - ! Perform collision physics for inelastics scattering + ! Perform collision physics for inelastic scattering call inelastic_scatter(nuc, rxn, p % E, p % coord0 % uvw, & p % mu, p % wgt) p % event_MT = rxn % MT From 10de23b56132e35d650ef7aff04d076754d6180b Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 23 Oct 2014 01:04:13 -0400 Subject: [PATCH 08/28] Remove union energy grid treatment. --- src/constants.F90 | 3 +- src/cross_section.F90 | 33 ---------- src/energy_grid.F90 | 144 +----------------------------------------- src/global.F90 | 26 +++----- src/hdf5_summary.F90 | 18 ++---- src/initialize.F90 | 10 +-- src/input_xml.F90 | 2 +- src/output.F90 | 7 -- 8 files changed, 23 insertions(+), 220 deletions(-) diff --git a/src/constants.F90 b/src/constants.F90 index 986ea83f6..6144f537d 100644 --- a/src/constants.F90 +++ b/src/constants.F90 @@ -366,8 +366,7 @@ module constants ! Energy grid methods integer, parameter :: & GRID_NUCLIDE = 1, & ! non-unionized energy grid - GRID_UNION = 2, & ! union grid with pointers - GRID_LOGARITHM = 3 ! logarithmic mapping + GRID_LOGARITHM = 2 ! logarithmic mapping integer, parameter :: N_LOG_BINS = 8000 ! Running modes diff --git a/src/cross_section.F90 b/src/cross_section.F90 index 5d483c5d6..93f5c2db3 100644 --- a/src/cross_section.F90 +++ b/src/cross_section.F90 @@ -14,9 +14,6 @@ module cross_section implicit none save - integer :: union_grid_index -!$omp threadprivate(union_grid_index) - contains !=============================================================================== @@ -50,9 +47,6 @@ contains mat => materials(p % material) - ! Find energy index on unionized grid - if (grid_method == GRID_UNION) call find_energy_index(p % E) - ! Determine if this material has S(a,b) tables check_sab = (mat % n_sab > 0) @@ -154,12 +148,6 @@ contains ! Determine index on nuclide energy grid select case (grid_method) - case (GRID_UNION) - ! If we're using the unionized grid with pointers, finding the index on - ! the nuclide energy grid is as simple as looking up the pointer - - i_grid = nuc % grid_index(union_grid_index) - case (GRID_LOGARITHM) ! Determine the energy grid index using a logarithmic mapping to reduce ! the energy range over which a binary search needs to be performed @@ -517,27 +505,6 @@ contains end subroutine calculate_urr_xs -!=============================================================================== -! FIND_ENERGY_INDEX determines the index on the union energy grid at a certain -! energy -!=============================================================================== - - subroutine find_energy_index(E) - - real(8), intent(in) :: E ! energy of particle - - ! if particle's energy is outside of energy grid range, set to first or last - ! index. Otherwise, do a binary search through the union energy grid. - if (E < e_grid(1)) then - union_grid_index = 1 - elseif (E > e_grid(n_grid)) then - union_grid_index = n_grid - 1 - else - union_grid_index = binary_search(e_grid, n_grid, E) - end if - - end subroutine find_energy_index - !=============================================================================== ! 0K_ELASTIC_XS determines the microscopic 0K elastic cross section ! for a given nuclide at the trial relative energy used in resonance scattering diff --git a/src/energy_grid.F90 b/src/energy_grid.F90 index 87ee71b16..6636538a3 100644 --- a/src/energy_grid.F90 +++ b/src/energy_grid.F90 @@ -1,152 +1,10 @@ module energy_grid - use constants, only: MAX_LINE_LEN, N_LOG_BINS + use constants, only: N_LOG_BINS use global - use list_header, only: ListReal - use output, only: write_message contains -!=============================================================================== -! UNIONIZED_GRID creates a single unionized energy grid combined from each -! nuclide of each material. Right now, the grid for each nuclide is added into a -! linked list one at a time with an effective insertion sort. Could be done with -! a hash for all energy points and then a quicksort at the end (what hash -! function to use?) -!=============================================================================== - - subroutine unionized_grid() - - integer :: i ! index in nuclides array - type(ListReal), pointer :: list => null() - type(Nuclide), pointer :: nuc => null() - - call write_message("Creating unionized energy grid...", 5) - - ! Add grid points for each nuclide in the problem - do i = 1, n_nuclides_total - nuc => nuclides(i) - call add_grid_points(list, nuc % energy) - end do - - ! Set size of unionized energy grid - n_grid = list % size() - - ! create allocated array from linked list - allocate(e_grid(n_grid)) - do i = 1, n_grid - e_grid(i) = list % get_item(i) - end do - - ! delete linked list and dictionary - call list % clear() - deallocate(list) - - ! Set pointers to unionized energy grid for each nuclide - call grid_pointers() - - end subroutine unionized_grid - -!=============================================================================== -! ADD_GRID_POINTS adds energy points from the 'energy' array into a linked list -! of points already stored from previous arrays. -!=============================================================================== - - subroutine add_grid_points(list, energy) - - type(ListReal), pointer :: list - real(8), intent(in) :: energy(:) - - integer :: i ! index in energy array - integer :: n ! size of energy array - integer :: current ! current index - real(8) :: E ! actual energy value - - i = 1 - n = size(energy) - - ! If the original list is empty, we need to allocate the first element and - ! store first energy point - if (.not. associated(list)) then - allocate(list) - do i = 1, n - call list % append(energy(i)) - end do - return - end if - - ! Set current index to beginning of the list - current = 1 - - do while (i <= n) - E = energy(i) - - ! If we've reached the end of the grid energy list, add the remaining - ! energy points to the end - if (current > list % size()) then - ! Finish remaining energies - do while (i <= n) - call list % append(energy(i)) - i = i + 1 - end do - exit - end if - - if (E < list % get_item(current)) then - - ! Insert new energy in this position - call list % insert(current, E) - - ! Advance index in linked list and in new energy grid - i = i + 1 - current = current + 1 - - elseif (E == list % get_item(current)) then - ! Found the exact same energy, no need to store duplicates so just - ! skip and move to next index - i = i + 1 - current = current + 1 - else - current = current + 1 - end if - - end do - - end subroutine add_grid_points - -!=============================================================================== -! GRID_POINTERS creates an array of pointers (ints) for each nuclide to link -! each point on the nuclide energy grid to one on the unionized energy grid -!=============================================================================== - - subroutine grid_pointers() - - integer :: i ! loop index for nuclides - integer :: j ! loop index for nuclide energy grid - integer :: index_e ! index on union energy grid - real(8) :: union_energy ! energy on union grid - real(8) :: energy ! energy on nuclide grid - type(Nuclide), pointer :: nuc => null() - - do i = 1, n_nuclides_total - nuc => nuclides(i) - allocate(nuc % grid_index(n_grid)) - - index_e = 1 - energy = nuc % energy(index_e) - - do j = 1, n_grid - union_energy = e_grid(j) - if (union_energy >= energy .and. index_e < nuc % n_grid) then - index_e = index_e + 1 - energy = nuc % energy(index_e) - end if - nuc % grid_index(j) = index_e - 1 - end do - end do - - end subroutine grid_pointers - !=============================================================================== ! LOGARITHMIC_GRID determines a logarithmic mapping for energies to bounding ! indices on a nuclide energy grid diff --git a/src/global.F90 b/src/global.F90 index 1477dfd3a..5db2632f4 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -62,7 +62,7 @@ module global ! Cross section arrays type(Nuclide), allocatable, target :: nuclides(:) ! Nuclide cross-sections type(SAlphaBeta), allocatable, target :: sab_tables(:) ! S(a,b) tables - type(XsListing), allocatable, target :: xs_listings(:) ! cross_sections.xml listings + type(XsListing), allocatable, target :: xs_listings(:) ! cross_sections.xml listings ! Cross section caches type(NuclideMicroXS), allocatable :: micro_xs(:) ! Cache for each nuclide @@ -79,9 +79,7 @@ module global ! Unionized energy grid integer :: grid_method ! how to treat the energy grid - integer :: n_grid ! number of points on unionized grid real(8) :: log_spacing ! spacing on logarithmic grid - real(8), allocatable :: e_grid(:) ! energies on unionized grid ! Unreoslved resonance probablity tables logical :: urr_ptables_on = .true. @@ -224,7 +222,6 @@ module global type(Timer) :: time_total ! timer for total run type(Timer) :: time_initialize ! timer for initialization type(Timer) :: time_read_xs ! timer for reading cross sections - type(Timer) :: time_unionize ! timer for unionizing energy grid type(Timer) :: time_bank ! timer for fission bank synchronization type(Timer) :: time_bank_sample ! timer for fission bank sampling type(Timer) :: time_bank_sendrecv ! timer for fission bank SEND/RECV @@ -301,7 +298,7 @@ module global logical :: write_initial_source = .false. ! ============================================================================ - ! CMFD VARIABLES + ! CMFD VARIABLES ! Main object type(cmfd_type) :: cmfd @@ -311,11 +308,11 @@ module global ! CMFD communicator integer :: cmfd_comm - + ! Timing objects type(Timer) :: time_cmfd ! timer for whole cmfd calculation type(Timer) :: time_cmfdbuild ! timer for matrix build - type(Timer) :: time_cmfdsolve ! timer for solver + type(Timer) :: time_cmfdsolve ! timer for solver ! Flag for active core map logical :: cmfd_coremap = .false. @@ -391,7 +388,7 @@ module global ! RESONANCE SCATTERING VARIABLES logical :: treat_res_scat = .false. ! is resonance scattering treated? - integer :: n_res_scatterers_total = 0 ! total number of resonant scatterers + integer :: n_res_scatterers_total = 0 ! total number of resonant scatterers type(Nuclide0K), allocatable, target :: nuclides_0K(:) ! 0K nuclides info !$omp threadprivate(micro_xs, material_xs, fission_bank, n_bank, & @@ -400,14 +397,14 @@ module global contains !=============================================================================== -! FREE_MEMORY deallocates and clears all global allocatable arrays in the +! FREE_MEMORY deallocates and clears all global allocatable arrays in the ! program !=============================================================================== subroutine free_memory() - + integer :: i ! Loop Index - + ! Deallocate cells, surfaces, materials if (allocated(cells)) deallocate(cells) if (allocated(universes)) deallocate(universes) @@ -462,9 +459,6 @@ contains if (allocated(matching_bins)) deallocate(matching_bins) if (allocated(tally_maps)) deallocate(tally_maps) - ! Deallocate energy grid - if (allocated(e_grid)) deallocate(e_grid) - ! Deallocate fission and source bank and entropy !$omp parallel if (allocated(fission_bank)) deallocate(fission_bank) @@ -489,7 +483,7 @@ contains ! Deallocate track_identifiers if (allocated(track_identifiers)) deallocate(track_identifiers) - + ! Deallocate dictionaries call cell_dict % clear() call universe_dict % clear() @@ -526,7 +520,7 @@ contains if (allocated(ufs_mesh % width)) deallocate(ufs_mesh % width) deallocate(ufs_mesh) end if - + end subroutine free_memory end module global diff --git a/src/hdf5_summary.F90 b/src/hdf5_summary.F90 index c9e367200..fc4af5bc7 100644 --- a/src/hdf5_summary.F90 +++ b/src/hdf5_summary.F90 @@ -67,7 +67,7 @@ contains end if ! Terminate access to the file. - call su % file_close() + call su % file_close() end subroutine hdf5_write_summary @@ -80,7 +80,7 @@ contains ! Write version information call su % write_data(VERSION_MAJOR, "version_major") call su % write_data(VERSION_MINOR, "version_minor") - call su % write_data(VERSION_RELEASE, "version_release") + call su % write_data(VERSION_RELEASE, "version_release") ! Write current date and time call su % write_data(time_stamp(), "date_and_time") @@ -88,7 +88,7 @@ contains ! Write MPI information call su % write_data(n_procs, "n_procs") call su % write_attribute_string("n_procs", "description", & - "Number of MPI processes") + "Number of MPI processes") end subroutine hdf5_write_header @@ -144,12 +144,12 @@ contains call su % write_data("universe", "fill_type", & group="geometry/cells/cell " // trim(to_str(c % id))) call su % write_data(universes(c % fill) % id, "material", & - group="geometry/cells/cell " // trim(to_str(c % id))) + group="geometry/cells/cell " // trim(to_str(c % id))) case (CELL_LATTICE) call su % write_data("lattice", "fill_type", & group="geometry/cells/cell " // trim(to_str(c % id))) call su % write_data(lattices(c % fill) % id, "lattice", & - group="geometry/cells/cell " // trim(to_str(c % id))) + group="geometry/cells/cell " // trim(to_str(c % id))) end select ! Write list of bounding surfaces @@ -303,7 +303,7 @@ contains else n_z = 1 end if - + ! Write lattice universes allocate(lattice_universes(n_x, n_y, n_z)) do j = 1, n_x @@ -371,7 +371,7 @@ contains group="materials/material " // trim(to_str(m % id))) call su % write_data(m % i_sab_tables, "i_sab_tables", & length=m % n_sab, & - group="materials/material " // trim(to_str(m % id))) + group="materials/material " // trim(to_str(m % id))) end if end do @@ -659,8 +659,6 @@ contains group="timing") call su % write_data(time_read_xs % elapsed, "time_read_xs", & group="timing") - call su % write_data(time_unionize % elapsed, "time_unionize", & - group="timing") call su % write_data(time_transport % elapsed, "time_transport", & group="timing") call su % write_data(time_bank % elapsed, "time_bank", & @@ -685,8 +683,6 @@ contains "Total time elapsed for initialization (s)", group="timing") call su % write_attribute_string("time_read_xs", "description", & "Time reading cross-section libraries (s)", group="timing") - call su % write_attribute_string("time_unionize", "description", & - "Time unionizing energy grid (s)", group="timing") call su % write_attribute_string("time_transport", "description", & "Time in transport only (s)", group="timing") call su % write_attribute_string("time_bank", "description", & diff --git a/src/initialize.F90 b/src/initialize.F90 index a2965755a..10b178e0a 100644 --- a/src/initialize.F90 +++ b/src/initialize.F90 @@ -4,7 +4,7 @@ module initialize use bank_header, only: Bank use constants use dict_header, only: DictIntInt, ElemKeyValueII - use energy_grid, only: unionized_grid, logarithmic_grid + use energy_grid, only: logarithmic_grid use error, only: fatal_error, warning use geometry, only: neighbor_lists use geometry_header, only: Cell, Universe, Lattice, BASE_UNIVERSE @@ -108,12 +108,8 @@ contains ! Create linked lists for multiple instances of the same nuclide call same_nuclide_list() - ! Construct unionized energy grid from cross-sections - if (grid_method == GRID_UNION) then - call time_unionize % start() - call unionized_grid() - call time_unionize % stop() - elseif (grid_method == GRID_LOGARITHM) then + ! Construct logarithmic energy grid for cross-sections + if (grid_method == GRID_LOGARITHM) then call logarithmic_grid() end if diff --git a/src/input_xml.F90 b/src/input_xml.F90 index b16d7e870..39b9fd887 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -213,7 +213,7 @@ contains case ('nuclide') grid_method = GRID_NUCLIDE case ('union') - grid_method = GRID_UNION + call fatal_error("Union energy grid is no longer supported.") case ('logarithm', 'logarithmic', 'log') grid_method = GRID_LOGARITHM case default diff --git a/src/output.F90 b/src/output.F90 index 4d0ef7edd..835a7d1b0 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -1226,12 +1226,6 @@ contains end do end if - ! print summary of unionized energy grid - call header("UNIONIZED ENERGY GRID", unit=UNIT_SUMMARY) - write(UNIT_SUMMARY,*) "Points on energy grid: " // trim(to_str(n_grid)) - write(UNIT_SUMMARY,*) "Extra storage required: " // trim(to_str(& - n_grid*n_nuclides_total*4)) // " bytes" - ! print summary of variance reduction call header("VARIANCE REDUCTION", unit=UNIT_SUMMARY) if (survival_biasing) then @@ -1494,7 +1488,6 @@ contains ! display time elapsed for various sections write(ou,100) "Total time for initialization", time_initialize % elapsed write(ou,100) " Reading cross sections", time_read_xs % elapsed - write(ou,100) " Unionizing energy grid", time_unionize % elapsed write(ou,100) "Total time in simulation", time_inactive % elapsed + & time_active % elapsed write(ou,100) " Time in transport only", time_transport % elapsed From 8dba76302e48c0bddc63b94e36cbe36547311a95 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sat, 25 Oct 2014 13:46:58 -0400 Subject: [PATCH 09/28] Update documentation to reflect no use of union energy grid. --- docs/source/_images/uniongrid.svg | 792 ------------------------- docs/source/methods/cross_sections.rst | 53 +- docs/source/usersguide/input.rst | 10 +- src/relaxng/settings.rnc | 35 +- 4 files changed, 40 insertions(+), 850 deletions(-) delete mode 100644 docs/source/_images/uniongrid.svg diff --git a/docs/source/_images/uniongrid.svg b/docs/source/_images/uniongrid.svg deleted file mode 100644 index 27c3922fd..000000000 --- a/docs/source/_images/uniongrid.svg +++ /dev/null @@ -1,792 +0,0 @@ - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - image/svg+xml - - - - - - - - E1 - - E2 - - E3 - - E4 - - E5 - - E6 - - E7 - - E8 - UnionEnergy Grid - - 0 - NuclidePointers - - 0 - - 1 - - 1 - - 1 - - 2 - - 3 - - 3 - - - - - - - - - - E1 - - E2 - - E3 - NuclideEnergy Grid - - σ1 - - σ2 - - σ3 - NuclideCross Sections - - - - diff --git a/docs/source/methods/cross_sections.rst b/docs/source/methods/cross_sections.rst index a193dde96..ae1c4449e 100644 --- a/docs/source/methods/cross_sections.rst +++ b/docs/source/methods/cross_sections.rst @@ -47,45 +47,30 @@ there would be for burnup calculations. Thus, there is a strong motive to implement a method of reducing the number of energy grid searches in order to speed up the calculation. -Unionized Energy Grid ---------------------- +------------------- +Logarithmic Mapping +------------------- -The most naïve method to reduce the number of energy grid searches is to -construct a new energy grid that consists of the union of the energy points of -each nuclide and use this energy grid for all nuclides. This method is -computationally very efficient as it only requires one energy grid search at -each collision as well as one interpolation between cross section values since -the interpolation factor can be used for all nuclides. However, it requires -redundant storage of cross section values at points which were added to each -nuclide grid. This additional burden on memory storage can become quite -prohibitive. To lessen that burden, the unionized energy grid can be thinned -with cross sections reconstructed on the thinned energy grid. This method is -currently used by default in the Serpent Monte Carlo code. +To speed up energy grid searches, OpenMC uses logarithmic mapping technique +[Brown]_ to limit the range of energies that must be searched for each +nuclide. The entire energy range is divided up into equal-lethargy segments, and +the bounding energies of each segment are mapped to bounding indices on each of +the nuclide energy grids. By default, OpenMC uses 8000 equal-lethargy segments +as recommended by Brown. -Unionized Energy Grid with Nuclide Pointers -------------------------------------------- +------------- +Other Methods +------------- -While having a unionized grid that is used for all nuclides allows for very fast -lookup of cross sections, the burden on memory is in many circumstances -unacceptable. The OpenMC Monte Carlo code utilizes a method that allows for a -single energy grid search to be performed at every collision while avoiding the -redundant storage of cross section values. Instead of using the unionized grid -for every nuclide, the original energy grid of each nuclide is kept and a list -of pointers (of the same length as the unionized energy grid) is constructed for -each nuclide that gives the corresponding grid index on the nuclide grid for a -given grid index on the unionized grid. One must still interpolate on cross -section values for each nuclide since the interpolation factors will generally -be different. The figure below illustrates this method. All values within the -dashed box would need to be stored on a per-nuclide basis, and the union grid -would need to be stored once. This method is also referred to as *double -indexing* and is available as an option in Serpent (see paper by Leppanen_). +A good survey of other energy grid techniques, including unionized energy grids, +can be found in a paper by Leppanen_. -.. figure:: ../_images/uniongrid.* - :width: 600px - :align: center - :figclass: align-center +---------- +References +---------- - Mapping of union energy grid to nuclide energy grid through pointers. +.. [Brown] Forrest B. Brown, "New Hash-based Energy Lookup Algorithm for Monte + Carlo codes," LA-UR-14-24530, Los Alamos National Laboratory (2014). .. _MCNP: http://mcnp.lanl.gov .. _Serpent: http://montecarlo.vtt.fi diff --git a/docs/source/usersguide/input.rst b/docs/source/usersguide/input.rst index 0e9258189..c25a9f146 100644 --- a/docs/source/usersguide/input.rst +++ b/docs/source/usersguide/input.rst @@ -134,13 +134,11 @@ should be performed. It has the following attributes/sub-elements: ------------------------- The ```` element determines the treatment of the energy grid during -a simulation. The valid options are "nuclide", "union", and "logarithm". Setting -this element to "nuclide" will cause OpenMC to use a nuclide's energy grid when +a simulation. The valid options are "nuclide" and "logarithm". Setting this +element to "nuclide" will cause OpenMC to use a nuclide's energy grid when determining what points to interpolate between for determining cross sections -(i.e. non-unionized energy grid). Setting this element to "union" results in a -unionized energy grid with pointers to nuclide grids. Setting this element to -"logarithm" causes OpenMC to use a logarithmic mapping technique described in -LA-UR-14-24530_. +(i.e. non-unionized energy grid). Setting this element to "logarithm" causes +OpenMC to use a logarithmic mapping technique described in LA-UR-14-24530_. *Default*: logarithm diff --git a/src/relaxng/settings.rnc b/src/relaxng/settings.rnc index 458946562..5f793b316 100644 --- a/src/relaxng/settings.rnc +++ b/src/relaxng/settings.rnc @@ -3,19 +3,19 @@ element settings { ( element eigenvalue { - (element batches { xsd:positiveInteger } | + (element batches { xsd:positiveInteger } | attribute batches { xsd:positiveInteger }) & - (element inactive { xsd:nonNegativeInteger } | + (element inactive { xsd:nonNegativeInteger } | attribute inactive { xsd:nonNegativeInteger }) & - (element particles { xsd:positiveInteger } | + (element particles { xsd:positiveInteger } | attribute particles { xsd:positiveInteger }) & - (element generations_per_batch { xsd:positiveInteger } | + (element generations_per_batch { xsd:positiveInteger } | attribute generations_per_batch { xsd:positiveInteger })? } | element fixed_source { - (element batches { xsd:positiveInteger } | + (element batches { xsd:positiveInteger } | attribute batches { xsd:positiveInteger }) & - (element particles { xsd:positiveInteger } | + (element particles { xsd:positiveInteger } | attribute particles { xsd:positiveInteger }) } ) & @@ -27,15 +27,14 @@ element settings { (element weight_avg { xsd:double } | attribute weight_avg { xsd:double })? }? & - element energy_grid { ( "nuclide" | "union" | "log" | - "logarithm" | "logarithmic" ) }? & + element energy_grid { ( "nuclide" | "log" | "logarithm" | "logarithmic" ) }? & element entropy { - (element dimension { list { xsd:int+ } } | + (element dimension { list { xsd:int+ } } | attribute dimension { list { xsd:int+ } })? & - (element lower_left { list { xsd:double+ } } | + (element lower_left { list { xsd:double+ } } | attribute lower_left { list { xsd:double+ } }) & - (element upper_right { list { xsd:double+ } } | + (element upper_right { list { xsd:double+ } } | attribute upper_right { list { xsd:double+ } }) }? & @@ -94,7 +93,7 @@ element settings { ( (element batches { list { xsd:positiveInteger+ } } | attribute batches { list { xsd:positiveInteger+ } }) | - (element interval { xsd:positiveInteger } | + (element interval { xsd:positiveInteger } | attribute interval { xsd:positiveInteger }) ) }? & @@ -103,12 +102,12 @@ element settings { ( (element batches { list { xsd:positiveInteger+ } } | attribute batches { list { xsd:positiveInteger+ } }) | - (element interval { xsd:positiveInteger } | + (element interval { xsd:positiveInteger } | attribute interval { xsd:positiveInteger }) )? & - (element separate { xsd:boolean } | + (element separate { xsd:boolean } | attribute separate { xsd:boolean })? & - (element write { xsd:boolean } | + (element write { xsd:boolean } | attribute write { xsd:boolean })? & (element overwrite_latest { xsd:boolean} | attribute overwrite_latest {xsd:boolean})? @@ -125,11 +124,11 @@ element settings { element verbosity { xsd:positiveInteger }? & element uniform_fs{ - (element dimension { list { xsd:positiveInteger+ } } | + (element dimension { list { xsd:positiveInteger+ } } | attribute dimension { list { xsd:positiveInteger+ } }) & - (element lower_left { list { xsd:double+ } } | + (element lower_left { list { xsd:double+ } } | attribute lower_left { list { xsd:double+ } }) & - (element upper_right { list { xsd:double+ } } | + (element upper_right { list { xsd:double+ } } | attribute upper_right { list { xsd:double+ } }) }? & From 02dc08e274537f3de391a2e885b91b1ce4f80075 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sat, 25 Oct 2014 14:01:28 -0400 Subject: [PATCH 10/28] Fix headings in cross section documentation. --- docs/source/methods/cross_sections.rst | 2 -- 1 file changed, 2 deletions(-) diff --git a/docs/source/methods/cross_sections.rst b/docs/source/methods/cross_sections.rst index ae1c4449e..db2e5156e 100644 --- a/docs/source/methods/cross_sections.rst +++ b/docs/source/methods/cross_sections.rst @@ -47,7 +47,6 @@ there would be for burnup calculations. Thus, there is a strong motive to implement a method of reducing the number of energy grid searches in order to speed up the calculation. -------------------- Logarithmic Mapping ------------------- @@ -58,7 +57,6 @@ the bounding energies of each segment are mapped to bounding indices on each of the nuclide energy grids. By default, OpenMC uses 8000 equal-lethargy segments as recommended by Brown. -------------- Other Methods ------------- From 996a03d5ab90d5a13af442afe7e9fab5421a3aa7 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sat, 25 Oct 2014 16:37:30 -0400 Subject: [PATCH 11/28] Add new evaporation spectrum sampling method outlined in LA-UR-14-27694. Still need documentation. --- src/physics.F90 | 17 ++++++++++------- 1 file changed, 10 insertions(+), 7 deletions(-) diff --git a/src/physics.F90 b/src/physics.F90 index 8fa0575cd..d136d7306 100644 --- a/src/physics.F90 +++ b/src/physics.F90 @@ -797,7 +797,7 @@ contains sampling_scheme = 'cxs' end if - ! otherwise, use free gas model + ! otherwise, use free gas model else if (E >= FREE_GAS_THRESHOLD * kT .and. awr > ONE) then v_target = ZERO @@ -859,7 +859,7 @@ contains m = (nuc % elastic_0K(i_E_up + 1) - xs_up) & & / (nuc % energy_0K(i_E_up + 1) - nuc % energy_0K(i_E_up)) xs_up = xs_up + m * (E_up - nuc % energy_0K(i_E_up)) - + ! get max 0K xs value over range of practical relative energies xs_max = max(xs_low, & & maxval(nuc % elastic_0K(i_E_low + 1 : i_E_up - 1)), xs_up) @@ -972,7 +972,7 @@ contains case default call fatal_error("Not a recognized resonance scattering treatment!") end select - + end subroutine sample_target_velocity !=============================================================================== @@ -1835,14 +1835,15 @@ contains lc = 2 + 2*NR + 2*NE U = edist % data(lc + 1) + y = (E_in - U)/T + v = 1 - exp(-y) + ! sample outgoing energy based on evaporation spectrum probability ! density function n_sample = 0 do - r1 = prn() - r2 = prn() - E_out = -T * log(r1*r2) - if (E_out <= E_in - U) exit + x = -log((1 - v*prn())*(1 - v*prn())) + if (x <= y) exit ! check for large number of rejections n_sample = n_sample + 1 @@ -1852,6 +1853,8 @@ contains end if end do + E_out = x*T + case (11) ! ======================================================================= ! ENERGY-DEPENDENT WATT SPECTRUM From 332970ea2e03d2512e75575c9939b50b8d576f31 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Tue, 28 Oct 2014 21:53:39 -0400 Subject: [PATCH 12/28] Update documentation for evaporation spectrum. --- docs/source/methods/physics.rst | 17 +++++++++++------ 1 file changed, 11 insertions(+), 6 deletions(-) diff --git a/docs/source/methods/physics.rst b/docs/source/methods/physics.rst index 54dc91345..db6cfd89e 100644 --- a/docs/source/methods/physics.rst +++ b/docs/source/methods/physics.rst @@ -682,17 +682,20 @@ nuclear temperature, which is a function of the incoming energy of the neutron. The ACE format contains a list of nuclear temperatures versus incoming energies. The nuclear temperature is interpolated between neighboring incoming energies using a specified interpolation law. Once the temperature :math:`T` is -determined, we then calculate a candidate outgoing energy based on rule C45 in -the `Monte Carlo Sampler`_: +determined, we then calculate a candidate outgoing energy based on the algorithm +given in LA-UR-14-27694_: .. math:: :label: evaporation-E - E' = -T \log (\xi_1 \xi_2) + E' = -T \log ((1 - g\xi_1)(1 - g\xi_2)) -where :math:`\xi_1, \xi_2` are random numbers sampled on the unit -interval. The outgoing energy is only accepted according to a specified -restriction energy as in equation :eq:`maxwell-restriction`. +where :math:`g = 1 - e^{-w}`, :math:`w = (E - U)/T`, :math:`U` is the +restriction energy, and :math:`\xi_1, \xi_2` are random numbers sampled on the +unit interval. The outgoing energy is only accepted according to the restriction +energy as in equation :eq:`maxwell-restriction`. This algorithm has a much +higher rejection efficiency than the standard technique, i.e. rule C45 in the +`Monte Carlo Sampler`_. ACE Law 11 - Energy-Dependent Watt Spectrum +++++++++++++++++++++++++++++++++++++++++++ @@ -1591,6 +1594,8 @@ References .. _Monte Carlo Sampler: https://laws.lanl.gov/vhosts/mcnp.lanl.gov/pdf_files/la-9721_3rdmcsampler.pdf +.. _LA-UR-14-27694: http://permalink.lanl.gov/object/tr?what=info:lanl-repo/lareport/LA-UR-14-27694 + .. _MC21: http://www.osti.gov/bridge/servlets/purl/903083-HT5p1o/903083.pdf .. _Sutton and Brown: http://www.osti.gov/bridge/product.biblio.jsp?osti_id=307911 From 03ffd94e3a3e77ef7e4fcc3760b2d7c6dc4cc344 Mon Sep 17 00:00:00 2001 From: walshjon Date: Tue, 11 Nov 2014 15:01:46 -0800 Subject: [PATCH 13/28] updated copyright year --- src/output.F90 | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/output.F90 b/src/output.F90 index 4d0ef7edd..a66ea46e8 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -157,7 +157,7 @@ contains if (master) then write(UNIT=OUTPUT_UNIT, FMT='(1X,A,1X,I1,".",I1,".",I1)') & "OpenMC version", VERSION_MAJOR, VERSION_MINOR, VERSION_RELEASE - write(UNIT=OUTPUT_UNIT, FMT=*) "Copyright (c) 2011-2013 & + write(UNIT=OUTPUT_UNIT, FMT=*) "Copyright (c) 2011-2014 & &Massachusetts Institute of Technology" write(UNIT=OUTPUT_UNIT, FMT=*) "MIT/X license at & &" @@ -1173,7 +1173,7 @@ contains call header("OpenMC Monte Carlo Code", unit=UNIT_SUMMARY, level=1) write(UNIT=UNIT_SUMMARY, FMT=*) & - "Copyright: 2011-2013 Massachusetts Institute of Technology" + "Copyright: 2011-2014 Massachusetts Institute of Technology" write(UNIT=UNIT_SUMMARY, FMT='(1X,A,7X,2(I1,"."),I1)') & "Version:", VERSION_MAJOR, VERSION_MINOR, VERSION_RELEASE #ifdef GIT_SHA1 From 991cdb9cccf7c07ac4e1d22a169222d42c43171c Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Fri, 14 Nov 2014 16:54:01 -0500 Subject: [PATCH 14/28] Change !$omp critical directives to !$omp atomic where possible. --- src/geometry.F90 | 40 +++++++++++++++++++--------------------- src/physics.F90 | 13 +++++-------- src/tally.F90 | 45 +++++++++++++++------------------------------ src/tracking.F90 | 9 +++------ 4 files changed, 42 insertions(+), 65 deletions(-) diff --git a/src/geometry.F90 b/src/geometry.F90 index 043ded05b..3e124983a 100644 --- a/src/geometry.F90 +++ b/src/geometry.F90 @@ -11,7 +11,7 @@ module geometry use tally, only: score_surface_current implicit none - + contains !=============================================================================== @@ -264,13 +264,13 @@ contains lattice_edge = .true. end if end if - + if (lattice_edge) then - + ! In this case the neutron is leaving the lattice, so we move it ! out, remove all lower coordinate levels and then search from ! universe 0. - + p % coord => p % coord0 call deallocate_coord(p % coord % next) @@ -287,7 +287,7 @@ contains p % last_material = p % material p % material = c % material - ! We'll still make a new coordinate for the particle, as + ! We'll still make a new coordinate for the particle, as ! distance_to_boundary will still need to track through lattice ! widths even though there's nothing in them but this material @@ -406,10 +406,9 @@ contains ! Score to global leakage tally if (tallies_on) then -!$omp critical +!$omp atomic global_tallies(LEAKAGE) % value = & global_tallies(LEAKAGE) % value + p % wgt -!$omp end critical end if ! Display message @@ -644,7 +643,7 @@ contains return end if end if - + end subroutine cross_surface !=============================================================================== @@ -906,7 +905,7 @@ contains if (quad < ZERO) then ! no intersection with cylinder - d = INFINITY + d = INFINITY elseif (on_surface) then ! particle is on the cylinder, thus one distance is @@ -955,7 +954,7 @@ contains if (quad < ZERO) then ! no intersection with cylinder - d = INFINITY + d = INFINITY elseif (on_surface) then ! particle is on the cylinder, thus one distance is @@ -1004,7 +1003,7 @@ contains if (quad < ZERO) then ! no intersection with cylinder - d = INFINITY + d = INFINITY elseif (on_surface) then ! particle is on the cylinder, thus one distance is @@ -1051,7 +1050,7 @@ contains if (quad < ZERO) then ! no intersection with sphere - d = INFINITY + d = INFINITY elseif (on_surface) then ! particle is on the sphere, thus one distance is @@ -1098,7 +1097,7 @@ contains if (quad < ZERO) then ! no intersection with cone - d = INFINITY + d = INFINITY elseif (on_surface) then ! particle is on the cone, thus one distance is positive/negative @@ -1117,7 +1116,7 @@ contains d = (-k - quad)/a b = (-k + quad)/a - ! determine the smallest positive solution + ! determine the smallest positive solution if (d < ZERO) then if (b > ZERO) then d = b @@ -1147,7 +1146,7 @@ contains if (quad < ZERO) then ! no intersection with cone - d = INFINITY + d = INFINITY elseif (on_surface) then ! particle is on the cone, thus one distance is positive/negative @@ -1166,7 +1165,7 @@ contains d = (-k - quad)/a b = (-k + quad)/a - ! determine the smallest positive solution + ! determine the smallest positive solution if (d < ZERO) then if (b > ZERO) then d = b @@ -1196,7 +1195,7 @@ contains if (quad < ZERO) then ! no intersection with cone - d = INFINITY + d = INFINITY elseif (on_surface) then ! particle is on the cone, thus one distance is positive/negative @@ -1215,7 +1214,7 @@ contains d = (-k - quad)/a b = (-k + quad)/a - ! determine the smallest positive solution + ! determine the smallest positive solution if (d < ZERO) then if (b > ZERO) then d = b @@ -1273,7 +1272,7 @@ contains ! logic here checks whether the relative difference is within floating ! point precision. - if (d < dist) then + if (d < dist) then if (abs(d - dist)/dist >= FP_REL_PRECISION) then dist = d if (u > 0) then @@ -1565,9 +1564,8 @@ contains ! Increment number of lost particles p % alive = .false. -!$omp critical +!$omp atomic n_lost_particles = n_lost_particles + 1 -!$omp end critical ! Abort the simulation if the maximum number of lost particles has been ! reached diff --git a/src/physics.F90 b/src/physics.F90 index 8fa0575cd..7988970ef 100644 --- a/src/physics.F90 +++ b/src/physics.F90 @@ -256,22 +256,19 @@ contains p % last_wgt = p % wgt ! Score implicit absorption estimate of keff -!$omp critical +!$omp atomic global_tallies(K_ABSORPTION) % value = & global_tallies(K_ABSORPTION) % value + p % absorb_wgt * & micro_xs(i_nuclide) % nu_fission / micro_xs(i_nuclide) % absorption -!$omp end critical - else ! See if disappearance reaction happens if (micro_xs(i_nuclide) % absorption > & prn() * micro_xs(i_nuclide) % total) then ! Score absorption estimate of keff -!$omp critical +!$omp atomic global_tallies(K_ABSORPTION) % value = & global_tallies(K_ABSORPTION) % value + p % wgt * & micro_xs(i_nuclide) % nu_fission / micro_xs(i_nuclide) % absorption -!$omp end critical p % alive = .false. p % event = EVENT_ABSORB @@ -797,7 +794,7 @@ contains sampling_scheme = 'cxs' end if - ! otherwise, use free gas model + ! otherwise, use free gas model else if (E >= FREE_GAS_THRESHOLD * kT .and. awr > ONE) then v_target = ZERO @@ -859,7 +856,7 @@ contains m = (nuc % elastic_0K(i_E_up + 1) - xs_up) & & / (nuc % energy_0K(i_E_up + 1) - nuc % energy_0K(i_E_up)) xs_up = xs_up + m * (E_up - nuc % energy_0K(i_E_up)) - + ! get max 0K xs value over range of practical relative energies xs_max = max(xs_low, & & maxval(nuc % elastic_0K(i_E_low + 1 : i_E_up - 1)), xs_up) @@ -972,7 +969,7 @@ contains case default call fatal_error("Not a recognized resonance scattering treatment!") end select - + end subroutine sample_target_velocity !=============================================================================== diff --git a/src/tally.F90 b/src/tally.F90 index 87f5a2cab..7c7747ff5 100644 --- a/src/tally.F90 +++ b/src/tally.F90 @@ -281,10 +281,9 @@ contains ! get the score and tally it score = last_wgt * calc_pn(n, mu) -!$omp critical +!$omp atomic t % results(score_index, filter_index) % value = & t % results(score_index, filter_index) % value + score -!$omp end critical end do j = j + t % moment_order(j) cycle SCORE_LOOP @@ -347,10 +346,9 @@ contains ! get the score and tally it score = wgt * calc_pn(n, mu) -!$omp critical +!$omp atomic t % results(score_index, filter_index) % value = & t % results(score_index, filter_index) % value + score -!$omp end critical end do j = j + t % moment_order(j) cycle SCORE_LOOP @@ -542,10 +540,9 @@ contains end select ! Add score to tally -!$omp critical +!$omp atomic t % results(score_index, filter_index) % value = & t % results(score_index, filter_index) % value + score -!$omp end critical end do SCORE_LOOP @@ -617,10 +614,9 @@ contains i_filter = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1 ! Add score to tally -!$omp critical +!$omp atomic t % results(i_score, i_filter) % value = & t % results(i_score, i_filter) % value + score -!$omp end critical end do ! reset outgoing energy bin and score index @@ -1015,10 +1011,9 @@ contains end if ! Add score to tally -!$omp critical +!$omp atomic t % results(score_index, filter_index) % value = & t % results(score_index, filter_index) % value + score -!$omp end critical end do SCORE_LOOP @@ -1214,10 +1209,9 @@ contains end select ! Add score to tally -!$omp critical +!$omp atomic t % results(score_index, filter_index) % value = & t % results(score_index, filter_index) % value + score -!$omp end critical end do SCORE_LOOP @@ -1366,10 +1360,9 @@ contains end select ! Add score to tally -!$omp critical +!$omp atomic t % results(score_index, filter_index) % value = & t % results(score_index, filter_index) % value + score -!$omp end critical end do MATERIAL_SCORE_LOOP @@ -1787,10 +1780,9 @@ contains end if ! Add score to tally -!$omp critical +!$omp atomic t % results(score_index, filter_index) % value = & t % results(score_index, filter_index) % value + score -!$omp end critical end do SCORE_LOOP @@ -2021,10 +2013,9 @@ contains matching_bins(i_filter_mesh) = & mesh_indices_to_bin(m, ijk0 + 1, .true.) filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1 -!$omp critical +!$omp atomic t % results(1, filter_index) % value = & t % results(1, filter_index) % value + p % wgt -!$omp end critical end if end do else @@ -2035,10 +2026,9 @@ contains matching_bins(i_filter_mesh) = & mesh_indices_to_bin(m, ijk0 + 1, .true.) filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1 -!$omp critical +!$omp atomic t % results(1, filter_index) % value = & t % results(1, filter_index) % value + p % wgt -!$omp end critical end if end do end if @@ -2053,10 +2043,9 @@ contains matching_bins(i_filter_mesh) = & mesh_indices_to_bin(m, ijk0 + 1, .true.) filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1 -!$omp critical +!$omp atomic t % results(1, filter_index) % value = & t % results(1, filter_index) % value + p % wgt -!$omp end critical end if end do else @@ -2067,10 +2056,9 @@ contains matching_bins(i_filter_mesh) = & mesh_indices_to_bin(m, ijk0 + 1, .true.) filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1 -!$omp critical +!$omp atomic t % results(1, filter_index) % value = & t % results(1, filter_index) % value + p % wgt -!$omp end critical end if end do end if @@ -2085,10 +2073,9 @@ contains matching_bins(i_filter_mesh) = & mesh_indices_to_bin(m, ijk0 + 1, .true.) filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1 -!$omp critical +!$omp atomic t % results(1, filter_index) % value = & t % results(1, filter_index) % value + p % wgt -!$omp end critical end if end do else @@ -2099,10 +2086,9 @@ contains matching_bins(i_filter_mesh) = & mesh_indices_to_bin(m, ijk0 + 1, .true.) filter_index = sum((matching_bins(1:t%n_filters) - 1) * t % stride) + 1 -!$omp critical +!$omp atomic t % results(1, filter_index) % value = & t % results(1, filter_index) % value + p % wgt -!$omp end critical end if end do end if @@ -2224,10 +2210,9 @@ contains end if ! Add to surface current tally -!$omp critical +!$omp atomic t % results(1, filter_index) % value = & t % results(1, filter_index) % value + p % wgt -!$omp end critical end if ! Calculate new coordinates diff --git a/src/tracking.F90 b/src/tracking.F90 index d108325fa..8eca3871c 100644 --- a/src/tracking.F90 +++ b/src/tracking.F90 @@ -62,9 +62,8 @@ contains n_event = 0 ! Add paricle's starting weight to count for normalizing tallies later -!$omp critical +!$omp atomic total_weight = total_weight + p % wgt -!$omp end critical ! Force calculation of cross-sections by setting last energy to zero micro_xs % last_E = ZERO @@ -112,11 +111,10 @@ contains call score_tracklength_tally(p, distance) ! Score track-length estimate of k-eff -!$omp critical +!$omp atomic global_tallies(K_TRACKLENGTH) % value = & global_tallies(K_TRACKLENGTH) % value + p % wgt * distance * & material_xs % nu_fission -!$omp end critical if (d_collision > d_boundary) then ! ==================================================================== @@ -140,11 +138,10 @@ contains ! PARTICLE HAS COLLISION ! Score collision estimate of keff -!$omp critical +!$omp atomic global_tallies(K_COLLISION) % value = & global_tallies(K_COLLISION) % value + p % wgt * & material_xs % nu_fission / material_xs % total -!$omp end critical ! score surface current tallies -- this has to be done before the collision ! since the direction of the particle will change and we need to use the From 3706b38bf20713cd03ce0c5f1c3406dd8286968d Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Tue, 18 Nov 2014 22:15:35 -0500 Subject: [PATCH 15/28] Make global_tallies allocatable instead of automatic. --- src/global.F90 | 23 ++++++++++++----------- src/tally_initialize.F90 | 3 +++ 2 files changed, 15 insertions(+), 11 deletions(-) diff --git a/src/global.F90 b/src/global.F90 index f4c50b833..874b8d0d1 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -62,7 +62,7 @@ module global ! Cross section arrays type(Nuclide), allocatable, target :: nuclides(:) ! Nuclide cross-sections type(SAlphaBeta), allocatable, target :: sab_tables(:) ! S(a,b) tables - type(XsListing), allocatable, target :: xs_listings(:) ! cross_sections.xml listings + type(XsListing), allocatable, target :: xs_listings(:) ! cross_sections.xml listings ! Cross section caches type(NuclideMicroXS), allocatable :: micro_xs(:) ! Cache for each nuclide @@ -119,7 +119,7 @@ module global ! 2) track-length estimate of k-eff ! 3) leakage fraction - type(TallyResult), target :: global_tallies(N_GLOBAL_TALLIES) + type(TallyResult), allocatable, target :: global_tallies(:) ! Tally map structure type(TallyMap), allocatable :: tally_maps(:) @@ -300,7 +300,7 @@ module global logical :: write_initial_source = .false. ! ============================================================================ - ! CMFD VARIABLES + ! CMFD VARIABLES ! Main object type(cmfd_type) :: cmfd @@ -310,11 +310,11 @@ module global ! CMFD communicator integer :: cmfd_comm - + ! Timing objects type(Timer) :: time_cmfd ! timer for whole cmfd calculation type(Timer) :: time_cmfdbuild ! timer for matrix build - type(Timer) :: time_cmfdsolve ! timer for solver + type(Timer) :: time_cmfdsolve ! timer for solver ! Flag for active core map logical :: cmfd_coremap = .false. @@ -390,7 +390,7 @@ module global ! RESONANCE SCATTERING VARIABLES logical :: treat_res_scat = .false. ! is resonance scattering treated? - integer :: n_res_scatterers_total = 0 ! total number of resonant scatterers + integer :: n_res_scatterers_total = 0 ! total number of resonant scatterers type(Nuclide0K), allocatable, target :: nuclides_0K(:) ! 0K nuclides info !$omp threadprivate(micro_xs, material_xs, fission_bank, n_bank, & @@ -399,14 +399,14 @@ module global contains !=============================================================================== -! FREE_MEMORY deallocates and clears all global allocatable arrays in the +! FREE_MEMORY deallocates and clears all global allocatable arrays in the ! program !=============================================================================== subroutine free_memory() - + integer :: i ! Loop Index - + ! Deallocate cells, surfaces, materials if (allocated(cells)) deallocate(cells) if (allocated(universes)) deallocate(universes) @@ -449,6 +449,7 @@ contains if (allocated(entropy_p)) deallocate(entropy_p) ! Deallocate tally-related arrays + if (allocated(global_tallies)) deallocate(global_tallies) if (allocated(meshes)) deallocate(meshes) if (allocated(tallies)) then ! First call the clear routines @@ -488,7 +489,7 @@ contains ! Deallocate track_identifiers if (allocated(track_identifiers)) deallocate(track_identifiers) - + ! Deallocate dictionaries call cell_dict % clear() call universe_dict % clear() @@ -525,7 +526,7 @@ contains if (allocated(ufs_mesh % width)) deallocate(ufs_mesh % width) deallocate(ufs_mesh) end if - + end subroutine free_memory end module global diff --git a/src/tally_initialize.F90 b/src/tally_initialize.F90 index 23ab5b563..b7dc5ed98 100644 --- a/src/tally_initialize.F90 +++ b/src/tally_initialize.F90 @@ -19,6 +19,9 @@ contains subroutine configure_tallies() + ! Allocate global tallies + allocate(global_tallies(N_GLOBAL_TALLIES)) + call setup_tally_arrays() call setup_tally_maps() From 233cdb3cd00b8a8642d4a740382808fb01500a20 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 11 Dec 2014 11:03:16 -0500 Subject: [PATCH 16/28] Mention CMake in quick install instructions. --- docs/source/quickinstall.rst | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/docs/source/quickinstall.rst b/docs/source/quickinstall.rst index 9dd6be6d1..3fb9036d3 100644 --- a/docs/source/quickinstall.rst +++ b/docs/source/quickinstall.rst @@ -35,9 +35,9 @@ OpenMC from source as is described in :ref:`usersguide_install`. Installing from Source on Linux or Mac OS X ------------------------------------------- -All OpenMC source code is hosted on GitHub_. If you have git_ and the gfortran_ -compiler installed, you can download and install OpenMC be entering the -following commands in a terminal: +All OpenMC source code is hosted on GitHub_. If you have git_, the gfortran_ +compiler, and CMake_ installed, you can download and install OpenMC be entering +the following commands in a terminal: .. code-block:: sh @@ -58,3 +58,4 @@ can be replaced with a local install, e.g. .. _GitHub: https://github.com/mit-crpg/openmc .. _git: http://git-scm.com .. _gfortran: http://gcc.gnu.org/wiki/GFortran +.. _CMake: http://www.cmake.org From 4b9e3f537063c11d28f8a695704fb206b4b22c51 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sun, 14 Dec 2014 21:03:44 -0500 Subject: [PATCH 17/28] Make number of bins for log grid adjustable by user. --- src/constants.F90 | 1 - src/energy_grid.F90 | 3 +-- src/global.F90 | 1 + src/input_xml.F90 | 11 +++++++++++ src/relaxng/settings.rnc | 2 ++ 5 files changed, 15 insertions(+), 3 deletions(-) diff --git a/src/constants.F90 b/src/constants.F90 index 9a9947baa..7e30bd07e 100644 --- a/src/constants.F90 +++ b/src/constants.F90 @@ -367,7 +367,6 @@ module constants integer, parameter :: & GRID_NUCLIDE = 1, & ! non-unionized energy grid GRID_LOGARITHM = 2 ! logarithmic mapping - integer, parameter :: N_LOG_BINS = 8000 ! Running modes integer, parameter :: & diff --git a/src/energy_grid.F90 b/src/energy_grid.F90 index bc111b0e0..98e476002 100644 --- a/src/energy_grid.F90 +++ b/src/energy_grid.F90 @@ -1,6 +1,5 @@ module energy_grid - use constants, only: N_LOG_BINS use global implicit none @@ -26,7 +25,7 @@ contains E_min = 1.0e-11_8 ! Determine equal-logarithmic energy spacing - M = N_LOG_BINS + M = n_log_bins log_spacing = log(E_max/E_min)/M ! Create equally log-spaced energy grid diff --git a/src/global.F90 b/src/global.F90 index 624511912..eb3283ea6 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -79,6 +79,7 @@ module global ! Unionized energy grid integer :: grid_method ! how to treat the energy grid + integer :: n_log_bins ! number of bins for logarithmic grid real(8) :: log_spacing ! spacing on logarithmic grid ! Unreoslved resonance probablity tables diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 6402463d5..281906ff3 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -220,6 +220,17 @@ contains call fatal_error("Unknown energy grid method: " // trim(temp_str)) end select + ! Number of bins for logarithmic grid + if (check_for_node(doc, "log_grid_bins")) then + call get_node_value(doc, "log_grid_bins", n_log_bins) + if (n_log_bins < 1) then + call fatal_error("Number of bins for logarithmic grid must be & + &greater than zero.") + end if + else + n_log_bins = 8000 + end if + ! Verbosity if (check_for_node(doc, "verbosity")) then call get_node_ptr(doc, "verbosity", node_verb) diff --git a/src/relaxng/settings.rnc b/src/relaxng/settings.rnc index 5f793b316..a7b92e173 100644 --- a/src/relaxng/settings.rnc +++ b/src/relaxng/settings.rnc @@ -38,6 +38,8 @@ element settings { attribute upper_right { list { xsd:double+ } }) }? & + element log_grid_bins { xsd:positiveInteger }? & + element natural_elements { xsd:string { maxLength = "20" } }? & element no_reduce { xsd:boolean }? & From 8659a76ea28b9eff0da0e5849a80927c3a0001cb Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Sun, 14 Dec 2014 21:09:09 -0500 Subject: [PATCH 18/28] Add documentation for log_grid_bins. --- docs/source/usersguide/input.rst | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/docs/source/usersguide/input.rst b/docs/source/usersguide/input.rst index 2278e11ff..c8a90c696 100644 --- a/docs/source/usersguide/input.rst +++ b/docs/source/usersguide/input.rst @@ -184,6 +184,16 @@ performed. It has the following attributes/sub-elements: *Default*: None +```` Element +--------------------------- + +The ```` element indicates the number of bins to use for the +logarithmic-mapped energy grid. Using more bins will result in energy grid +searches over a smaller range at the expense of more memory. The default is +based on the recommended value in LA-UR-14-24530_. + + *Default*: 8000 + .. _natural_elements: ```` Element From 9680dbdaa3d01749bae57360f31323cbb9d38024 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Mon, 15 Dec 2014 21:04:23 -0500 Subject: [PATCH 19/28] Move energy grid related variables to energy_grid module. --- src/cross_section.F90 | 1 + src/energy_grid.F90 | 4 ++++ src/global.F90 | 5 ----- src/initialize.F90 | 2 +- src/input_xml.F90 | 3 ++- 5 files changed, 8 insertions(+), 7 deletions(-) diff --git a/src/cross_section.F90 b/src/cross_section.F90 index 93f5c2db3..7c95dfb37 100644 --- a/src/cross_section.F90 +++ b/src/cross_section.F90 @@ -2,6 +2,7 @@ module cross_section use ace_header, only: Nuclide, SAlphaBeta, Reaction, UrrData use constants + use energy_grid, only: grid_method, log_spacing use error, only: fatal_error use fission, only: nu_total use global diff --git a/src/energy_grid.F90 b/src/energy_grid.F90 index 98e476002..234a856c1 100644 --- a/src/energy_grid.F90 +++ b/src/energy_grid.F90 @@ -4,6 +4,10 @@ module energy_grid implicit none + integer :: grid_method ! how to treat the energy grid + integer :: n_log_bins ! number of bins for logarithmic grid + real(8) :: log_spacing ! spacing on logarithmic grid + contains !=============================================================================== diff --git a/src/global.F90 b/src/global.F90 index eb3283ea6..ec1121660 100644 --- a/src/global.F90 +++ b/src/global.F90 @@ -77,11 +77,6 @@ module global type(DictCharInt) :: sab_dict type(DictCharInt) :: xs_listing_dict - ! Unionized energy grid - integer :: grid_method ! how to treat the energy grid - integer :: n_log_bins ! number of bins for logarithmic grid - real(8) :: log_spacing ! spacing on logarithmic grid - ! Unreoslved resonance probablity tables logical :: urr_ptables_on = .true. diff --git a/src/initialize.F90 b/src/initialize.F90 index 2cf5b7eb9..1ab34ca0b 100644 --- a/src/initialize.F90 +++ b/src/initialize.F90 @@ -4,7 +4,7 @@ module initialize use bank_header, only: Bank use constants use dict_header, only: DictIntInt, ElemKeyValueII - use energy_grid, only: logarithmic_grid + use energy_grid, only: logarithmic_grid, grid_method use error, only: fatal_error, warning use geometry, only: neighbor_lists use geometry_header, only: Cell, Universe, Lattice, BASE_UNIVERSE diff --git a/src/input_xml.F90 b/src/input_xml.F90 index 281906ff3..e5c21bdc1 100644 --- a/src/input_xml.F90 +++ b/src/input_xml.F90 @@ -3,6 +3,7 @@ module input_xml use cmfd_input, only: configure_cmfd use constants use dict_header, only: DictIntInt, ElemKeyValueCI + use energy_grid, only: grid_method, n_log_bins use error, only: fatal_error, warning use geometry_header, only: Cell, Surface, Lattice use global @@ -2832,7 +2833,7 @@ contains ! Copy plot cell universe level if (check_for_node(node_plot, "level")) then call get_node_value(node_plot, "level", pl % level) - + if (pl % level < 0) then call fatal_error("Bad universe level in plot " & &// trim(to_str(pl % id))) From 0a46cc209b69f6b7ff28c9ba1b8ec4298e628ab2 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Tue, 16 Dec 2014 06:52:03 -0500 Subject: [PATCH 20/28] Mention ERSN-OpenMC in user's guide section on writing input files. --- docs/source/usersguide/input.rst | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/docs/source/usersguide/input.rst b/docs/source/usersguide/input.rst index 7132a67cd..e9107fb50 100644 --- a/docs/source/usersguide/input.rst +++ b/docs/source/usersguide/input.rst @@ -1572,3 +1572,15 @@ into MATLAB using PETSc-MATLAB utilities. This option can be turned on with "true" and off with "false". *Default*: false + +------------------------------------ +ERSN-OpenMC Graphical User Interface +------------------------------------ + +A third-party Java-based user-friendly graphical user interface for creating XML +input files called ERSN-OpenMC_ is developed and maintained by members of the +Radiation and Nuclear Systems Group at the Faculty of Sciences Tetouan, Morocco. +The GUI also allows one to automatically download prerequisites for installing and +running OpenMC. + +.. _ERSN-OpenMC: https://github.com/EL-Bakkali-Jaafar/ERSN-OpenMC From a157de4624f579f6c5d13f2959ce062c74d72748 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Tue, 16 Dec 2014 19:30:14 -0500 Subject: [PATCH 21/28] Don't allow source energies above 20 MeV for the time being. --- src/source.F90 | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/source.F90 b/src/source.F90 index 07c1bfb41..1f39200b6 100644 --- a/src/source.F90 +++ b/src/source.F90 @@ -211,6 +211,9 @@ contains case (SRC_ENERGY_MONO) ! Monoenergtic source site % E = external_source % params_energy(1) + if (site % E > 20) then + call fatal_error("Source energies above 20 MeV not allowed.") + end if case (SRC_ENERGY_MAXWELL) a = external_source % params_energy(1) From 2c20e54506a55c0e856decf4a33dc7f5157680c0 Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Wed, 17 Dec 2014 09:03:58 -0500 Subject: [PATCH 22/28] Prevent monoenergetic source at exactly 20 MeV --- src/source.F90 | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/source.F90 b/src/source.F90 index 1f39200b6..3a12a3f71 100644 --- a/src/source.F90 +++ b/src/source.F90 @@ -211,7 +211,7 @@ contains case (SRC_ENERGY_MONO) ! Monoenergtic source site % E = external_source % params_energy(1) - if (site % E > 20) then + if (site % E >= 20) then call fatal_error("Source energies above 20 MeV not allowed.") end if From fb75397bcfa454a7fa81706f4cc0fd0fff6b0b8c Mon Sep 17 00:00:00 2001 From: Adam Nelson Date: Sun, 21 Dec 2014 09:43:45 -0500 Subject: [PATCH 23/28] Updated documentation to include estimator option. --- docs/source/usersguide/input.rst | 10 ++++++++++ 1 file changed, 10 insertions(+) diff --git a/docs/source/usersguide/input.rst b/docs/source/usersguide/input.rst index 217535af9..4d96506ec 100644 --- a/docs/source/usersguide/input.rst +++ b/docs/source/usersguide/input.rst @@ -1034,6 +1034,16 @@ The ```` element accepts the following sub-elements: *Default*: total + :estimator: + The estimator element is used to force the us eof either ``analog`` or + ``tracklength`` tally estimation. ''analog'' is generally less efficient + though it can be used with every score type. ''tracklength'' is generally + the most efficient, though its usage is restricted to tallies which do not + score particle information which requires a collision to have occured, such + as a scattering tally which utilizes outgoing energy filters. + + *Default*: ``tracklength`` but will revert to analog if necessary. + :scores: A space-separated list of the desired responses to be accumulated. Accepted options are "flux", "total", "scatter", "absorption", "fission", From 906dfba75446b500675e265231dbffd28dc02a82 Mon Sep 17 00:00:00 2001 From: Adam Nelson Date: Mon, 22 Dec 2014 19:39:48 -0500 Subject: [PATCH 24/28] Fixed typos found by @paulromano --- docs/source/usersguide/input.rst | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/docs/source/usersguide/input.rst b/docs/source/usersguide/input.rst index 4d96506ec..c834bbc92 100644 --- a/docs/source/usersguide/input.rst +++ b/docs/source/usersguide/input.rst @@ -1035,10 +1035,10 @@ The ```` element accepts the following sub-elements: *Default*: total :estimator: - The estimator element is used to force the us eof either ``analog`` or + The estimator element is used to force the use of either ``analog`` or ``tracklength`` tally estimation. ''analog'' is generally less efficient though it can be used with every score type. ''tracklength'' is generally - the most efficient, though its usage is restricted to tallies which do not + the most efficient, though its usage is restricted to tallies that do not score particle information which requires a collision to have occured, such as a scattering tally which utilizes outgoing energy filters. From 3ea5b2fdaae7f1b2e8067e475c2cdc98a13cd873 Mon Sep 17 00:00:00 2001 From: Nick Horelik Date: Wed, 31 Dec 2014 12:48:18 -0500 Subject: [PATCH 25/28] Flipped y pixels for meshline plotting --- src/plot.F90 | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/src/plot.F90 b/src/plot.F90 index ddc4641b2..d9f6b55e2 100644 --- a/src/plot.F90 +++ b/src/plot.F90 @@ -272,11 +272,11 @@ contains outrange(1) = int(frac * real(img % width, 8)) frac = (xyz_ur(outer) - xyz_ll_plot(outer)) / width(outer) outrange(2) = int(frac * real(img % width, 8)) - - frac = (xyz_ll(inner) - xyz_ll_plot(inner)) / width(inner) - inrange(1) = int(frac * real(img % height, 8)) + frac = (xyz_ur(inner) - xyz_ll_plot(inner)) / width(inner) - inrange(2) = int(frac * real(img % height, 8)) + inrange(1) = int((1. - frac) * real(img % height, 8)) + frac = (xyz_ll(inner) - xyz_ll_plot(inner)) / width(inner) + inrange(2) = int((1. - frac) * real(img % height, 8)) ! draw lines do out_ = outrange(1), outrange(2) From c34ab5175ac13253b4589c31ea847b90038a040d Mon Sep 17 00:00:00 2001 From: Sterling Harper Date: Fri, 2 Jan 2015 21:43:02 -0500 Subject: [PATCH 26/28] Fixed doc typo --- docs/source/methods/geometry.rst | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/docs/source/methods/geometry.rst b/docs/source/methods/geometry.rst index 1515ffa89..1ba29b79e 100644 --- a/docs/source/methods/geometry.rst +++ b/docs/source/methods/geometry.rst @@ -421,7 +421,7 @@ satisfy the following equations x^2 + y^2 + z^2 - 10^2 < 0 \\ x - (-3) > 0 \\ - x - 2 < 0 + y - 2 < 0 In order to determine if a point is inside the cell, we would substitute its coordinates into equation :eq:`cell-contains-example`. If the inequalities are From 8a52e5fcdf679cb4fdc1b405d3695ce46376f14e Mon Sep 17 00:00:00 2001 From: John Xia Date: Mon, 5 Jan 2015 18:04:29 -0600 Subject: [PATCH 27/28] Fixing tally.F90 segfault when trying to tally a particle in a void material. --- src/tally.F90 | 6 ++++-- 1 file changed, 4 insertions(+), 2 deletions(-) diff --git a/src/tally.F90 b/src/tally.F90 index 7c7747ff5..c96f1ecd6 100644 --- a/src/tally.F90 +++ b/src/tally.F90 @@ -1837,8 +1837,10 @@ contains p % coord % universe, i_tally) case (FILTER_MATERIAL) - matching_bins(i) = get_next_bin(FILTER_MATERIAL, & - p % material, i_tally) + if (p % material /= MATERIAL_VOID) then + matching_bins(i) = get_next_bin(FILTER_MATERIAL, & + p % material, i_tally) + endif case (FILTER_CELL) ! determine next cell bin From 17466efe2fd44e1a2127f3e905e8039821d61a6e Mon Sep 17 00:00:00 2001 From: walshjon Date: Wed, 7 Jan 2015 23:09:54 -0800 Subject: [PATCH 28/28] updated copyright date to 2015 --- src/output.F90 | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/output.F90 b/src/output.F90 index 26a4a317d..002696671 100644 --- a/src/output.F90 +++ b/src/output.F90 @@ -50,7 +50,7 @@ contains ! Write version information write(UNIT=OUTPUT_UNIT, FMT=*) & - ' Copyright: 2011-2014 Massachusetts Institute of Technology' + ' Copyright: 2011-2015 Massachusetts Institute of Technology' write(UNIT=OUTPUT_UNIT, FMT=*) & ' License: http://mit-crpg.github.io/openmc/license.html' write(UNIT=OUTPUT_UNIT, FMT='(6X,"Version:",8X,I1,".",I1,".",I1)') & @@ -157,7 +157,7 @@ contains if (master) then write(UNIT=OUTPUT_UNIT, FMT='(1X,A,1X,I1,".",I1,".",I1)') & "OpenMC version", VERSION_MAJOR, VERSION_MINOR, VERSION_RELEASE - write(UNIT=OUTPUT_UNIT, FMT=*) "Copyright (c) 2011-2014 & + write(UNIT=OUTPUT_UNIT, FMT=*) "Copyright (c) 2011-2015 & &Massachusetts Institute of Technology" write(UNIT=OUTPUT_UNIT, FMT=*) "MIT/X license at & &" @@ -1173,7 +1173,7 @@ contains call header("OpenMC Monte Carlo Code", unit=UNIT_SUMMARY, level=1) write(UNIT=UNIT_SUMMARY, FMT=*) & - "Copyright: 2011-2014 Massachusetts Institute of Technology" + "Copyright: 2011-2015 Massachusetts Institute of Technology" write(UNIT=UNIT_SUMMARY, FMT='(1X,A,7X,2(I1,"."),I1)') & "Version:", VERSION_MAJOR, VERSION_MINOR, VERSION_RELEASE #ifdef GIT_SHA1