From 10de23b56132e35d650ef7aff04d076754d6180b Mon Sep 17 00:00:00 2001 From: Paul Romano Date: Thu, 23 Oct 2014 01:04:13 -0400 Subject: [PATCH] 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 986ea83f65..6144f537d8 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 5d483c5d65..93f5c2db33 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 87ee71b16e..6636538a3c 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 1477dfd3a3..5db2632f4a 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 c9e367200c..fc4af5bc7b 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 a2965755a6..10b178e0a7 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 b16d7e870f..39b9fd8871 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 4d0ef7edd4..835a7d1b06 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