crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 14:15:42 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

Merge branch 'develop' into xml

Browse source 
Conflicts:
	src/DEPENDENCIES
	src/Makefile
	src/OBJECTS
	src/utils/build_dependencies.py
This commit is contained in:
Paul Romano 2013-08-15 00:30:04 -05:00
commit 304d33c092
219 changed files with 226809 additions and 1084 deletions
Download patch file Download diff file Expand all files Collapse all files

3
.gitignore vendored
View file

@ -25,3 +25,6 @@ src/xml-fortran/xmlreader
# Modules built from XML templates
src/templates/*.f90
# Test results error file
results_error.dat

5
docs/source/publications.rst
View file

@ -4,6 +4,11 @@
Publications
============
- Paul K. Romano, Andrew R. Siegel, Benoit Forget, and Kord Smith, "Data
decomposition of Monte Carlo particle transport simulations via tally
servers," *J. Comput. Phys.*, **252**, 20--36
(2013). `<http://dx.doi.org/10.1016/j.jcp.2013.06.011>`_
- Paul K. Romano, Bryan R. Herman, Nicholas E. Horelik, Benoit Forget, Kord
Smith, and Andrew R. Siegel, "Progress and Status of the OpenMC Monte Carlo
Code," *Proc. Int. Conf. Mathematics and Computational Methods Applied to

32
src/DEPENDENCIES
View file

@ -91,6 +91,7 @@ cross_section.o: error.o
cross_section.o: fission.o
cross_section.o: global.o
cross_section.o: material_header.o
cross_section.o: particle_header.o
cross_section.o: random_lcg.o
cross_section.o: search.o
@ -104,13 +105,14 @@ eigenvalue.o: math.o
eigenvalue.o: mesh.o
eigenvalue.o: mesh_header.o
eigenvalue.o: output.o
eigenvalue.o: physics.o
eigenvalue.o: particle_header.o
eigenvalue.o: random_lcg.o
eigenvalue.o: search.o
eigenvalue.o: source.o
eigenvalue.o: state_point.o
eigenvalue.o: string.o
eigenvalue.o: tally.o
eigenvalue.o: tracking.o
endf.o: constants.o
endf.o: string.o
@ -138,12 +140,13 @@ fission.o: search.o
fixed_source.o: constants.o
fixed_source.o: global.o
fixed_source.o: output.o
fixed_source.o: physics.o
fixed_source.o: particle_header.o
fixed_source.o: random_lcg.o
fixed_source.o: source.o
fixed_source.o: state_point.o
fixed_source.o: string.o
fixed_source.o: tally.o
fixed_source.o: tracking.o
geometry.o: constants.o
geometry.o: error.o
@ -239,6 +242,7 @@ main.o: particle_restart.o
main.o: plot.o
math.o: constants.o
math.o: random_lcg.o
mesh.o: constants.o
mesh.o: global.o
@ -268,6 +272,7 @@ output_interface.o: mpiio_interface.o
output_interface.o: tally_header.o
particle_header.o: constants.o
particle_header.o: geometry_header.o
particle_restart.o: bank_header.o
particle_restart.o: constants.o
@ -276,13 +281,13 @@ particle_restart.o: global.o
particle_restart.o: hdf5_interface.o
particle_restart.o: output.o
particle_restart.o: particle_header.o
particle_restart.o: physics.o
particle_restart.o: random_lcg.o
particle_restart.o: source.o
particle_restart.o: tracking.o
particle_restart_write.o: bank_header.o
particle_restart_write.o: global.o
particle_restart_write.o: hdf5_interface.o
particle_restart_write.o: particle_header.o
particle_restart_write.o: string.o
physics.o: ace_header.o
@ -291,11 +296,10 @@ physics.o: cross_section.o
physics.o: endf.o
physics.o: error.o
physics.o: fission.o
physics.o: geometry.o
physics.o: geometry_header.o
physics.o: global.o
physics.o: interpolation.o
physics.o: material_header.o
physics.o: math.o
physics.o: mesh.o
physics.o: output.o
physics.o: particle_header.o
@ -303,7 +307,6 @@ physics.o: particle_restart_write.o
physics.o: random_lcg.o
physics.o: search.o
physics.o: string.o
physics.o: tally.o
plot.o: constants.o
plot.o: error.o
@ -314,7 +317,6 @@ plot.o: output.o
plot.o: particle_header.o
plot.o: plot_header.o
plot.o: ppmlib.o
plot.o: source.o
plot.o: string.o
plot_header.o: constants.o
@ -332,9 +334,9 @@ source.o: constants.o
source.o: error.o
source.o: geometry_header.o
source.o: global.o
source.o: math.o
source.o: output.o
source.o: particle_header.o
source.o: physics.o
source.o: random_lcg.o
source.o: string.o
@ -371,7 +373,17 @@ tally_initialize.o: tally_header.o
timer_header.o: constants.o
tracking.o: cross_section.o
tracking.o: error.o
tracking.o: geometry.o
tracking.o: geometry_header.o
tracking.o: global.o
tracking.o: output.o
tracking.o: particle_header.o
tracking.o: physics.o
tracking.o: string.o
tracking.o: tally.o
xml_interface.o: constants.o
xml_interface.o: error.o
xml_interface.o: global.o

8
src/Makefile
View file

@ -1,14 +1,10 @@
program = openmc
prefix = /usr/local
source = $(wildcard *.F90)
objects = $(source:.F90=.o)
xml_lib = -Lxml/lib -lxml
#===============================================================================
# Object Files
#===============================================================================
include OBJECTS
#===============================================================================
# User Options
#===============================================================================

63
src/OBJECTS
View file

@ -1,63 +0,0 @@
objects = \
ace.o \
ace_header.o \
bank_header.o \
cmfd_data.o \
cmfd_execute.o \
cmfd_header.o \
cmfd_input.o \
cmfd_jacobian_operator.o \
cmfd_loss_operator.o \
cmfd_message_passing.o \
cmfd_output.o \
cmfd_power_solver.o \
cmfd_prod_operator.o \
cmfd_snes_solver.o \
constants.o \
cross_section.o \
dict_header.o \
doppler.o \
eigenvalue.o \
endf.o \
endf_header.o \
energy_grid.o \
error.o \
finalize.o \
fission.o \
fixed_source.o \
geometry.o \
geometry_header.o \
global.o \
hdf5_interface.o \
hdf5_summary.o \
initialize.o \
interpolation.o \
input_xml.o \
list_header.o \
main.o \
material_header.o \
math.o \
mesh_header.o \
mesh.o \
mpiio_interface.o \
output.o \
output_interface.o \
particle_header.o \
particle_restart.o \
particle_restart_write.o \
physics.o \
plot.o \
plot_header.o \
ppmlib.o \
random_lcg.o \
search.o \
set_header.o \
source.o \
source_header.o \
state_point.o \
string.o \
tally.o \
tally_header.o \
tally_initialize.o \
timer_header.o \
xml_interface.o

96
src/cross_section.F90
View file

@ -6,10 +6,14 @@ module cross_section
use fission, only: nu_total
use global
use material_header, only: Material
use particle_header, only: Particle
use random_lcg, only: prn
use search, only: binary_search
implicit none
save
integer :: union_grid_index
contains
@ -18,7 +22,9 @@ contains
! particle is currently traveling through.
!===============================================================================
subroutine calculate_xs()
subroutine calculate_xs(p)
type(Particle), intent(inout) :: p
integer :: i ! loop index over nuclides
integer :: i_nuclide ! index into nuclides array
@ -42,7 +48,7 @@ contains
mat => materials(p % material)
! Find energy index on unionized grid
if (grid_method == GRID_UNION) call find_energy_index()
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)
@ -84,7 +90,7 @@ contains
! Calculate microscopic cross section for this nuclide
if (p % E /= micro_xs(i_nuclide) % last_E) then
call calculate_nuclide_xs(i_nuclide, i_sab)
call calculate_nuclide_xs(i_nuclide, i_sab, p % E)
end if
! ========================================================================
@ -125,10 +131,11 @@ contains
! given index in the nuclides array at the energy of the given particle
!===============================================================================
subroutine calculate_nuclide_xs(i_nuclide, i_sab)
subroutine calculate_nuclide_xs(i_nuclide, i_sab, E)
integer, intent(in) :: i_nuclide ! index into nuclides array
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
@ -143,19 +150,19 @@ contains
! 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(p % index_grid)
i_grid = nuc % grid_index(union_grid_index)
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
if (p % E < nuc % energy(1)) then
if (E < nuc % energy(1)) then
i_grid = 1
elseif (p % E > nuc % energy(nuc % n_grid)) then
elseif (E > nuc % energy(nuc % n_grid)) then
i_grid = nuc % n_grid - 1
else
i_grid = binary_search(nuc % energy, nuc % n_grid, p % E)
i_grid = binary_search(nuc % energy, nuc % n_grid, E)
end if
end select
@ -164,7 +171,7 @@ contains
if (nuc % energy(i_grid) == nuc % energy(i_grid+1)) i_grid = i_grid + 1
! calculate interpolation factor
f = (p%E - nuc%energy(i_grid))/(nuc%energy(i_grid+1) - nuc%energy(i_grid))
f = (E - nuc%energy(i_grid))/(nuc%energy(i_grid+1) - nuc%energy(i_grid))
micro_xs(i_nuclide) % index_grid = i_grid
micro_xs(i_nuclide) % interp_factor = f
@ -213,22 +220,22 @@ contains
! need to correct it by subtracting the non-S(a,b) elastic cross section and
! then add back in the calculated S(a,b) elastic+inelastic cross section.
if (i_sab > 0) call calculate_sab_xs(i_nuclide, i_sab)
if (i_sab > 0) call calculate_sab_xs(i_nuclide, i_sab, E)
! if the particle is in the unresolved resonance range and there are
! probability tables, we need to determine cross sections from the table
if (urr_ptables_on .and. nuc % urr_present) then
if (p % E > nuc % urr_data % energy(1) .and. &
p % E < nuc % urr_data % energy(nuc % urr_data % n_energy)) then
call calculate_urr_xs(i_nuclide)
if (E > nuc % urr_data % energy(1) .and. &
E < nuc % urr_data % energy(nuc % urr_data % n_energy)) then
call calculate_urr_xs(i_nuclide, E)
end if
end if
! Set last evaluated energy -- if we're in S(a,b) region, force
! re-calculation of cross-section
if (i_sab == 0) then
micro_xs(i_nuclide) % last_E = p % E
micro_xs(i_nuclide) % last_E = E
else
micro_xs(i_nuclide) % last_E = ZERO
end if
@ -241,10 +248,11 @@ contains
! whatever data were taken from the normal Nuclide table.
!===============================================================================
subroutine calculate_sab_xs(i_nuclide, i_sab)
subroutine calculate_sab_xs(i_nuclide, i_sab, E)
integer, intent(in) :: i_nuclide ! index into nuclides array
integer, intent(in) :: i_sab ! index into sab_tables array
real(8), intent(in) :: E ! energy
integer :: i_grid ! index on S(a,b) energy grid
real(8) :: f ! interp factor on S(a,b) energy grid
@ -259,12 +267,12 @@ contains
sab => sab_tables(i_sab)
! Get index and interpolation factor for inelastic grid
if (p%E < sab % inelastic_e_in(1)) then
if (E < sab % inelastic_e_in(1)) then
i_grid = 1
f = ZERO
else
i_grid = binary_search(sab % inelastic_e_in, sab % n_inelastic_e_in, p%E)
f = (p%E - sab%inelastic_e_in(i_grid)) / &
i_grid = binary_search(sab % inelastic_e_in, sab % n_inelastic_e_in, E)
f = (E - sab%inelastic_e_in(i_grid)) / &
(sab%inelastic_e_in(i_grid+1) - sab%inelastic_e_in(i_grid))
end if
@ -273,32 +281,32 @@ contains
f * sab % inelastic_sigma(i_grid + 1)
! Check for elastic data
if (p % E < sab % threshold_elastic) then
if (E < sab % threshold_elastic) then
! Determine whether elastic scattering is given in the coherent or
! incoherent approximation. For coherent, the cross section is
! represented as P/E whereas for incoherent, it is simply P
if (sab % elastic_mode == SAB_ELASTIC_EXACT) then
if (p % E < sab % elastic_e_in(1)) then
if (E < sab % elastic_e_in(1)) then
! If energy is below that of the lowest Bragg peak, the elastic
! cross section will be zero
elastic = ZERO
else
i_grid = binary_search(sab % elastic_e_in, &
sab % n_elastic_e_in, p % E)
elastic = sab % elastic_P(i_grid) / p % E
sab % n_elastic_e_in, E)
elastic = sab % elastic_P(i_grid) / E
end if
else
! Determine index on elastic energy grid
if (p % E < sab % elastic_e_in(1)) then
if (E < sab % elastic_e_in(1)) then
i_grid = 1
else
i_grid = binary_search(sab % elastic_e_in, &
sab % n_elastic_e_in, p%E)
sab % n_elastic_e_in, E)
end if
! Get interpolation factor for elastic grid
f = (p%E - sab%elastic_e_in(i_grid))/(sab%elastic_e_in(i_grid+1) - &
f = (E - sab%elastic_e_in(i_grid))/(sab%elastic_e_in(i_grid+1) - &
sab%elastic_e_in(i_grid))
! Calculate S(a,b) elastic scattering cross section
@ -325,9 +333,10 @@ contains
! from probability tables
!===============================================================================
subroutine calculate_urr_xs(i_nuclide)
subroutine calculate_urr_xs(i_nuclide, E)
integer, intent(in) :: i_nuclide ! index into nuclides array
real(8), intent(in) :: E ! energy
integer :: i_energy ! index for energy
integer :: i_table ! index for table
@ -350,12 +359,12 @@ contains
! determine energy table
i_energy = 1
do
if (p % E < urr % energy(i_energy + 1)) exit
if (E < urr % energy(i_energy + 1)) exit
i_energy = i_energy + 1
end do
! determine interpolation factor on table
f = (p % E - urr % energy(i_energy)) / &
f = (E - urr % energy(i_energy)) / &
(urr % energy(i_energy + 1) - urr % energy(i_energy))
! sample probability table using the cumulative distribution
@ -377,7 +386,7 @@ contains
f * urr % prob(i_energy + 1, URR_N_GAMMA, i_table)
elseif (urr % interp == LOG_LOG) then
! Get logarithmic interpolation factor
f = log(p % E / urr % energy(i_energy)) / &
f = log(E / urr % energy(i_energy)) / &
log(urr % energy(i_energy + 1) / urr % energy(i_energy))
! Calculate elastic cross section/factor
@ -440,44 +449,31 @@ contains
! Determine nu-fission cross section
if (nuc % fissionable) then
micro_xs(i_nuclide) % nu_fission = nu_total(nuc, p % E) * &
micro_xs(i_nuclide) % nu_fission = nu_total(nuc, E) * &
micro_xs(i_nuclide) % fission
end if
end subroutine calculate_urr_xs
!===============================================================================
! FIND_ENERGY_INDEX determines the index on the union energy grid and the
! interpolation factor for a particle at a certain energy
! FIND_ENERGY_INDEX determines the index on the union energy grid at a certain
! energy
!===============================================================================
subroutine find_energy_index()
subroutine find_energy_index(E)
integer :: i_grid ! index on union energy grid
real(8) :: E ! energy of particle
real(8) :: interp ! interpolation factor
! copy particle's energy
E = p % 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
i_grid = 1
union_grid_index = 1
elseif (E > e_grid(n_grid)) then
i_grid = n_grid - 1
union_grid_index = n_grid - 1
else
i_grid = binary_search(e_grid, n_grid, E)
union_grid_index = binary_search(e_grid, n_grid, E)
end if
! calculate the interpolation factor -- note this will be outside of [0,1)
! for a particle outside the energy range of the union grid
interp = (E - e_grid(i_grid))/(e_grid(i_grid+1) - e_grid(i_grid))
! set particle attributes
p % index_grid = i_grid
p % interp = interp
end subroutine find_energy_index
end module cross_section

12
src/eigenvalue.F90
View file

@ -13,7 +13,7 @@ module eigenvalue
use mesh_header, only: StructuredMesh
use output, only: write_message, header, print_columns, &
print_batch_keff, print_generation
use physics, only: transport
use particle_header, only: Particle
use random_lcg, only: prn, set_particle_seed, prn_skip
use search, only: binary_search
use source, only: get_source_particle
@ -21,6 +21,7 @@ module eigenvalue
use string, only: to_str
use tally, only: synchronize_tallies, setup_active_usertallies, &
reset_result
use tracking, only: transport
private
public :: run_eigenvalue
@ -37,10 +38,9 @@ contains
subroutine run_eigenvalue()
if (master) call header("K EIGENVALUE SIMULATION", level=1)
type(Particle) :: p
! Allocate particle
allocate(p)
if (master) call header("K EIGENVALUE SIMULATION", level=1)
! Display column titles
call print_columns()
@ -74,10 +74,10 @@ contains
PARTICLE_LOOP: do current_work = 1, work
! grab source particle from bank
call get_source_particle(current_work)
call get_source_particle(p, current_work)
! transport particle
call transport()
call transport(p)
end do PARTICLE_LOOP

3
src/error.F90
View file

@ -139,9 +139,6 @@ contains
if (current_batch > 0) then
write(ERROR_UNIT,'(1X,A,I12) ') 'Batch: ', current_batch
write(ERROR_UNIT,'(1X,A,I12) ') 'Generation:', current_gen
write(ERROR_UNIT,'(1X,A,I12)') 'Particle: ', p % id
write(ERROR_UNIT,'(1X,A,3ES12.4)') 'Location: ', p % coord0 % xyz
write(ERROR_UNIT,'(1X,A,3ES12.4)') 'Direction: ', p % coord0 % uvw
write(ERROR_UNIT,*)
end if

34
src/fixed_source.F90
View file

@ -1,15 +1,15 @@
module fixed_source
use constants, only: ZERO
use constants, only: ZERO
use global
use output, only: write_message, header
use physics, only: transport
use random_lcg, only: set_particle_seed
use source, only: initialize_particle, sample_external_source, &
copy_source_attributes
use state_point, only: write_state_point
use string, only: to_str
use tally, only: synchronize_tallies, setup_active_usertallies
use output, only: write_message, header
use particle_header, only: Particle
use random_lcg, only: set_particle_seed
use source, only: sample_external_source, copy_source_attributes
use state_point, only: write_state_point
use string, only: to_str
use tally, only: synchronize_tallies, setup_active_usertallies
use tracking, only: transport
type(Bank), pointer :: source_site => null()
@ -17,12 +17,12 @@ contains
subroutine run_fixedsource()
integer(8) :: i ! index over histories in single cycle
integer(8) :: i ! index over histories in single cycle
type(Particle) :: p
if (master) call header("FIXED SOURCE TRANSPORT SIMULATION", level=1)
! Allocate particle and dummy source site
allocate(p)
allocate(source_site)
! Turn timer and tallies on
@ -61,10 +61,10 @@ contains
call set_particle_seed(p % id)
! grab source particle from bank
call sample_source_particle()
call sample_source_particle(p)
! transport particle
call transport()
call transport(p)
end do PARTICLE_LOOP
@ -120,16 +120,18 @@ contains
! SAMPLE_SOURCE_PARTICLE
!===============================================================================
subroutine sample_source_particle()
subroutine sample_source_particle(p)
type(Particle), intent(inout) :: p
! Set particle
call initialize_particle()
call p % initialize()
! Sample the external source distribution
call sample_external_source(source_site)
! Copy source attributes to the particle
call copy_source_attributes(source_site)
call copy_source_attributes(p, source_site)
end subroutine sample_source_particle

94
src/geometry.F90
View file

@ -5,7 +5,7 @@ module geometry
use geometry_header, only: Cell, Surface, Universe, Lattice
use global
use output, only: write_message
use particle_header, only: LocalCoord, deallocate_coord
use particle_header, only: LocalCoord, deallocate_coord, Particle
use particle_restart_write, only: write_particle_restart
use string, only: to_str
use tally, only: score_surface_current
@ -19,10 +19,11 @@ contains
! particle are inside a cell defined as the intersection of a series of surfaces
!===============================================================================
function simple_cell_contains(c) result(in_cell)
function simple_cell_contains(c, p) result(in_cell)
type(Cell), pointer :: c
logical :: in_cell
type(Cell), pointer :: c
type(Particle), intent(inout) :: p
logical :: in_cell
integer :: i ! index of surfaces in cell
integer :: i_surface ! index in surfaces array (with sign)
@ -47,7 +48,7 @@ contains
! Determine the specified sense of the surface in the cell and the actual
! sense of the particle with respect to the surface
s => surfaces(abs(i_surface))
actual_sense = sense(s)
actual_sense = sense(p, s)
specified_sense = (c % surfaces(i) > 0)
! Compare sense of point to specified sense
@ -68,7 +69,9 @@ contains
! position using simple_cell_contains and the LocalCoord's built up by find_cell
!===============================================================================
subroutine check_cell_overlap()
subroutine check_cell_overlap(p)
type(Particle), intent(inout) :: p
integer :: i ! cell loop index on a level
integer :: n ! number of cells to search on a level
@ -92,7 +95,7 @@ contains
index_cell = univ % cells(i)
c => cells(index_cell)
if (simple_cell_contains(c)) then
if (simple_cell_contains(c, p)) then
! the particle should only be contained in one cell per level
if (index_cell /= coord % cell) then
message = "Overlapping cells detected: " // &
@ -120,10 +123,11 @@ contains
! as it's within the geometry
!===============================================================================
recursive subroutine find_cell(found, search_cells)
recursive subroutine find_cell(p, found, search_cells)
logical, intent(inout) :: found
integer, optional :: search_cells(:)
type(Particle), intent(inout) :: p
logical, intent(inout) :: found
integer, optional :: search_cells(:)
integer :: i ! index over cells
integer :: i_x, i_y, i_z ! indices in lattice
@ -166,7 +170,7 @@ contains
! get pointer to cell
c => cells(index_cell)
if (simple_cell_contains(c)) then
if (simple_cell_contains(c, p)) then
! Set cell on this level
p % coord % cell = index_cell
@ -209,7 +213,7 @@ contains
p % coord % rotated = .true.
end if
call find_cell(found)
call find_cell(p, found)
if (.not. found) exit
elseif (c % type == CELL_LATTICE) then
@ -331,7 +335,7 @@ contains
end if
if (.not. outside_lattice) then
call find_cell(found)
call find_cell(p, found)
if (.not. found) exit
end if
@ -352,9 +356,10 @@ contains
! the geometry, is reflected, or crosses into a new lattice or cell
!===============================================================================
subroutine cross_surface(last_cell)
subroutine cross_surface(p, last_cell)
integer, intent(in) :: last_cell ! last cell particle was in
type(Particle), intent(inout) :: p
integer, intent(in) :: last_cell ! last cell particle was in
real(8) :: x ! x-x0 for sphere
real(8) :: y ! y-y0 for sphere
@ -395,7 +400,7 @@ contains
! physically moving the particle forward slightly
p % coord0 % xyz = p % coord0 % xyz + TINY_BIT * p % coord0 % uvw
call score_surface_current()
call score_surface_current(p)
end if
! Score to global leakage tally
@ -417,7 +422,7 @@ contains
if (.not. associated(p % coord, p % coord0)) then
message = "Cannot reflect particle " // trim(to_str(p % id)) // &
" off surface in a lower universe."
call handle_lost_particle()
call handle_lost_particle(p)
return
end if
@ -427,7 +432,7 @@ contains
if (active_current_tallies % size() > 0) then
p % coord0 % xyz = p % coord0 % xyz - TINY_BIT * p % coord0 % uvw
call score_surface_current()
call score_surface_current(p)
p % coord0 % xyz = p % coord0 % xyz + TINY_BIT * p % coord0 % uvw
end if
@ -568,10 +573,10 @@ contains
if (associated(p % coord0 % next)) then
call deallocate_coord(p % coord0 % next)
call find_cell(found)
call find_cell(p, found)
if (.not. found) then
message = "Couldn't find particle after reflecting from surface."
call handle_lost_particle()
call handle_lost_particle(p)
return
end if
end if
@ -594,14 +599,14 @@ contains
! If coming from negative side of surface, search all the neighboring
! cells on the positive side
call find_cell(found, surf % neighbor_pos)
call find_cell(p, found, surf % neighbor_pos)
if (found) return
elseif (p % surface < 0 .and. allocated(surf % neighbor_neg)) then
! If coming from positive side of surface, search all the neighboring
! cells on the negative side
call find_cell(found, surf % neighbor_neg)
call find_cell(p, found, surf % neighbor_neg)
if (found) return
end if
@ -613,7 +618,7 @@ contains
p % surface = NONE
p % coord => p % coord0
call deallocate_coord(p % coord % next)
call find_cell(found)
call find_cell(p, found)
if (run_mode /= MODE_PLOTTING .and. (.not. found)) then
! If a cell is still not found, there are two possible causes: 1) there is
@ -624,7 +629,7 @@ contains
p % coord => p % coord0
call deallocate_coord(p % coord % next)
p % coord % xyz = p % coord % xyz + TINY_BIT * p % coord % uvw
call find_cell(found)
call find_cell(p, found)
! Couldn't find next cell anywhere! This probably means there is an actual
! undefined region in the geometry.
@ -633,7 +638,7 @@ contains
message = "After particle " // trim(to_str(p % id)) // " crossed surface " &
// trim(to_str(surfaces(i_surface) % id)) // " it could not be &
&located in any cell and it did not leak."
call handle_lost_particle()
call handle_lost_particle(p)
return
end if
end if
@ -644,9 +649,10 @@ contains
! CROSS_LATTICE moves a particle into a new lattice element
!===============================================================================
subroutine cross_lattice(lattice_crossed)
subroutine cross_lattice(p, lattice_crossed)
integer, intent(in) :: lattice_crossed
type(Particle), intent(inout) :: p
integer, intent(in) :: lattice_crossed
integer :: i_x, i_y, i_z ! indices in lattice
integer :: n_x, n_y, n_z ! size of lattice
@ -722,11 +728,11 @@ contains
p % coord => p % coord0
! Search for particle
call find_cell(found)
call find_cell(p, found)
if (.not. found) then
message = "Could not locate particle " // trim(to_str(p % id)) // &
" after crossing a lattice boundary."
call handle_lost_particle()
call handle_lost_particle(p)
return
end if
else
@ -734,7 +740,7 @@ contains
p % coord % universe = lat % universes(i_x, i_y, i_z)
! Find cell in next lattice element
call find_cell(found)
call find_cell(p, found)
if (.not. found) then
! In some circumstances, a particle crossing the corner of a cell may not
! be able to be found in the next universe. In this scenario we cut off
@ -745,11 +751,11 @@ contains
p % coord => p % coord0
! Search for particle
call find_cell(found)
call find_cell(p, found)
if (.not. found) then
message = "Could not locate particle " // trim(to_str(p % id)) // &
" after crossing a lattice boundary."
call handle_lost_particle()
call handle_lost_particle(p)
return
end if
end if
@ -763,11 +769,12 @@ contains
! that has a parent cell, also include the surfaces of the edge of the universe.
!===============================================================================
subroutine distance_to_boundary(dist, surface_crossed, lattice_crossed)
subroutine distance_to_boundary(p, dist, surface_crossed, lattice_crossed)
real(8), intent(out) :: dist
integer, intent(out) :: surface_crossed
integer, intent(out) :: lattice_crossed
type(Particle), intent(inout) :: p
real(8), intent(out) :: dist
integer, intent(out) :: surface_crossed
integer, intent(out) :: lattice_crossed
integer :: i ! index for surface in cell
integer :: index_surf ! index in surfaces array (with sign)
@ -1342,10 +1349,11 @@ contains
! is in.
!===============================================================================
recursive function sense(surf) result(s)
recursive function sense(p, surf) result(s)
type(Surface), pointer :: surf ! surface
logical :: s ! sense of particle
type(Particle), intent(inout) :: p
type(Surface), pointer :: surf ! surface
logical :: s ! sense of particle
real(8) :: x,y,z ! coordinates of particle
real(8) :: func ! surface function evaluated at point
@ -1451,7 +1459,7 @@ contains
! Particle may be coincident with this surface. Artifically move the
! particle forward a tiny bit.
p % coord % xyz = p % coord % xyz + TINY_BIT * p % coord % uvw
s = sense(surf)
s = sense(p, surf)
elseif (func > 0) then
s = .true.
else
@ -1544,11 +1552,13 @@ contains
! HANDLE_LOST_PARTICLE
!===============================================================================
subroutine handle_lost_particle()
subroutine handle_lost_particle(p)
type(Particle), intent(inout) :: p
! Print warning and write lost particle file
call warning()
call write_particle_restart()
call write_particle_restart(p)
! Increment number of lost particles
p % alive = .false.

6
src/global.F90
View file

@ -9,7 +9,6 @@ module global
use geometry_header, only: Cell, Universe, Lattice, Surface
use material_header, only: Material
use mesh_header, only: StructuredMesh
use particle_header, only: Particle
use plot_header, only: ObjectPlot
use set_header, only: SetInt
use source_header, only: ExtSource
@ -27,11 +26,6 @@ module global
implicit none
save
! ============================================================================
! THE PARTICLE
type(Particle), pointer :: p => null()
! ============================================================================
! GEOMETRY-RELATED VARIABLES

51
src/math.F90
View file

@ -1,6 +1,7 @@
module math
use constants, only: PI, ONE, TWO
use constants, only: PI, ONE, TWO
use random_lcg, only: prn
implicit none
@ -162,4 +163,52 @@ contains
end function calc_pn
!===============================================================================
! MAXWELL_SPECTRUM samples an energy from the Maxwell fission distribution based
! on a direct sampling scheme. The probability distribution function for a
! Maxwellian is given as p(x) = 2/(T*sqrt(pi))*sqrt(x/T)*exp(-x/T). This PDF can
! be sampled using rule C64 in the Monte Carlo Sampler LA-9721-MS.
!===============================================================================
function maxwell_spectrum(T) result(E_out)
real(8), intent(in) :: T ! tabulated function of incoming E
real(8) :: E_out ! sampled energy
real(8) :: r1, r2, r3 ! random numbers
real(8) :: c ! cosine of pi/2*r3
r1 = prn()
r2 = prn()
r3 = prn()
! determine cosine of pi/2*r
c = cos(PI/TWO*r3)
! determine outgoing energy
E_out = -T*(log(r1) + log(r2)*c*c)
end function maxwell_spectrum
!===============================================================================
! WATT_SPECTRUM samples the outgoing energy from a Watt energy-dependent fission
! spectrum. Although fitted parameters exist for many nuclides, generally the
! continuous tabular distributions (LAW 4) should be used in lieu of the Watt
! spectrum. This direct sampling scheme is an unpublished scheme based on the
! original Watt spectrum derivation (See F. Brown's MC lectures).
!===============================================================================
function watt_spectrum(a, b) result(E_out)
real(8), intent(in) :: a ! Watt parameter a
real(8), intent(in) :: b ! Watt parameter b
real(8) :: E_out ! energy of emitted neutron
real(8) :: w ! sampled from Maxwellian
w = maxwell_spectrum(a)
E_out = w + a*a*b/4. + (TWO*prn() - ONE)*sqrt(a*a*b*w)
end function watt_spectrum
end module math

8
src/output.F90
View file

@ -11,7 +11,7 @@ module output
use math, only: t_percentile
use mesh_header, only: StructuredMesh
use mesh, only: mesh_indices_to_bin, bin_to_mesh_indices
use particle_header, only: LocalCoord
use particle_header, only: LocalCoord, Particle
use plot_header
use string, only: upper_case, to_str
use tally_header, only: TallyObject
@ -235,7 +235,9 @@ contains
! PRINT_PARTICLE displays the attributes of a particle
!===============================================================================
subroutine print_particle()
subroutine print_particle(p)
type(Particle), intent(in) :: p
integer :: i ! index for coordinate levels
type(Cell), pointer :: c => null()
@ -301,8 +303,6 @@ contains
! Display weight, energy, grid index, and interpolation factor
write(ou,*) ' Weight = ' // to_str(p % wgt)
write(ou,*) ' Energy = ' // to_str(p % E)
write(ou,*) ' Energy grid index = ' // to_str(p % index_grid)
write(ou,*) ' Interpolation factor = ' // to_str(p % interp)
write(ou,*)
end subroutine print_particle

44
src/particle_header.F90
View file

@ -1,6 +1,7 @@
module particle_header
use constants, only: NEUTRON, ONE, NONE, ZERO
use constants, only: NEUTRON, ONE, NONE, ZERO
use geometry_header, only: BASE_UNIVERSE
implicit none
@ -65,10 +66,6 @@ module particle_header
integer :: n_bank ! number of fission sites banked
real(8) :: wgt_bank ! weight of fission sites banked
! Energy grid data
integer :: index_grid ! index on unionized energy grid
real(8) :: interp ! interpolation factor for energy grid
! Indices for various arrays
integer :: surface ! index for surface particle is on
integer :: cell_born ! index for cell particle was born in
@ -78,6 +75,8 @@ module particle_header
! Statistical data
integer :: n_collision ! # of collisions
contains
procedure :: initialize => initialize_particle
end type Particle
contains
@ -102,4 +101,39 @@ contains
end subroutine deallocate_coord
!===============================================================================
! INITIALIZE_PARTICLE sets default attributes for a particle from the source
! bank
!===============================================================================
subroutine initialize_particle(this)
class(Particle) :: this
! Set particle to neutron that's alive
this % type = NEUTRON
this % alive = .true.
! clear attributes
this % surface = NONE
this % cell_born = NONE
this % material = NONE
this % last_material = NONE
this % wgt = ONE
this % last_wgt = ONE
this % absorb_wgt = ZERO
this % n_bank = 0
this % wgt_bank = ZERO
this % n_collision = 0
! remove any original coordinates
call deallocate_coord(this % coord0)
! Set up base level coordinates
allocate(this % coord0)
this % coord0 % universe = BASE_UNIVERSE
this % coord => this % coord0
end subroutine initialize_particle
end module particle_header

25
src/particle_restart.F90
View file

@ -6,11 +6,10 @@ module particle_restart
use constants
use geometry_header, only: BASE_UNIVERSE
use global
use particle_header, only: deallocate_coord
use output, only: write_message, print_particle
use physics, only: transport
use particle_header, only: Particle
use random_lcg, only: set_particle_seed
use source, only: initialize_particle
use tracking, only: transport
#ifdef HDF5
use hdf5_interface
@ -36,7 +35,9 @@ contains
! READ_HDF5_PARTICLE_RESTART
!===============================================================================
subroutine read_hdf5_particle_restart()
subroutine read_hdf5_particle_restart(p)
type(Particle), intent(inout) :: p
! write meessage
message = "Loading particle restart file " // trim(path_particle_restart) &
@ -77,7 +78,9 @@ contains
! READ_BINARY_PARTICLE_RESTART
!===============================================================================
subroutine read_binary_particle_restart()
subroutine read_binary_particle_restart(p)
type(Particle), intent(inout) :: p
integer :: filetype
integer :: revision
@ -121,16 +124,16 @@ contains
subroutine run_particle_restart()
integer(8) :: particle_seed
type(Particle) :: p
! initialize the particle to be tracked
allocate(p)
call initialize_particle()
call p % initialize()
! read in the restart information
#ifdef HDF5
call read_hdf5_particle_restart()
call read_hdf5_particle_restart(p)
#else
call read_binary_particle_restart()
call read_binary_particle_restart(p)
#endif
! set all tallies to 0 for now (just tracking errors)
@ -142,10 +145,10 @@ contains
call set_particle_seed(particle_seed)
! transport neutron
call transport()
call transport(p)
! write output if particle made it
call print_particle()
call print_particle(p)
end subroutine run_particle_restart

21
src/particle_restart_write.F90
View file

@ -2,9 +2,10 @@ module particle_restart_write
use, intrinsic :: ISO_FORTRAN_ENV
use bank_header, only: Bank
use bank_header, only: Bank
use global
use string, only: to_str
use particle_header, only: Particle
use string, only: to_str
#ifdef HDF5
use hdf5_interface
@ -24,16 +25,18 @@ contains
! WRITE_PARTICLE_RESTART
!===============================================================================
subroutine write_particle_restart()
subroutine write_particle_restart(p)
type(Particle), intent(in) :: p
! Dont write another restart file if in particle restart mode
if (run_mode == MODE_PARTICLE) return
! write out binary or HDF5 file
#ifdef HDF5
call write_particle_restart_hdf5()
call write_particle_restart_hdf5(p)
#else
call write_particle_restart_binary()
call write_particle_restart_binary(p)
#endif
end subroutine write_particle_restart
@ -44,7 +47,9 @@ contains
! WRITE_PARTICLE_RESTART_HDF5
!===============================================================================
subroutine write_particle_restart_hdf5()
subroutine write_particle_restart_hdf5(p)
type(Particle), intent(in) :: p
character(MAX_FILE_LEN) :: filename
type(Bank), pointer :: src => null()
@ -88,7 +93,9 @@ contains
! WRITE_PARTICLE_RESTART_BINARY
!===============================================================================
subroutine write_particle_restart_binary()
subroutine write_particle_restart_binary(p)
type(Particle), intent(in) :: p
character(MAX_FILE_LEN) :: filename
type(Bank), pointer :: src => null()

1057
src/physics.F90
View file

File diff suppressed because it is too large Load diff

24
src/plot.F90
View file

@ -6,11 +6,10 @@ module plot
use geometry_header, only: Cell, BASE_UNIVERSE
use global
use output, only: write_message
use particle_header, only: deallocate_coord
use particle_header, only: deallocate_coord, Particle
use plot_header
use ppmlib, only: Image, init_image, allocate_image, &
deallocate_image, set_pixel
use source, only: initialize_particle
use string, only: to_str
implicit none
@ -49,8 +48,9 @@ contains
! current particle's position
!===============================================================================
subroutine position_rgb(pl, rgb, id)
subroutine position_rgb(p, pl, rgb, id)
type(Particle), intent(inout) :: p
type(ObjectPlot), pointer, intent(in) :: pl
integer, intent(out) :: rgb(3)
integer, intent(out) :: id
@ -61,8 +61,8 @@ contains
call deallocate_coord(p % coord0 % next)
p % coord => p % coord0
call find_cell(found_cell)
if (check_overlaps) call check_cell_overlap()
call find_cell(p, found_cell)
if (check_overlaps) call check_cell_overlap(p)
if (.not. found_cell) then
! If no cell, revert to default color
@ -110,6 +110,7 @@ contains
real(8) :: out_pixel
real(8) :: xyz(3)
type(Image) :: img
type(Particle) :: p
! Initialize and allocate space for image
call init_image(img)
@ -142,8 +143,7 @@ contains
end if
! allocate and initialize particle
allocate(p)
call initialize_particle()
call p % initialize()
p % coord % xyz = xyz
p % coord % uvw = (/ 0.5, 0.5, 0.5 /)
p % coord % universe = BASE_UNIVERSE
@ -152,7 +152,7 @@ contains
do x = 1, img % width
! get pixel color
call position_rgb(pl, rgb, id)
call position_rgb(p, pl, rgb, id)
! Create a pixel at (x,y) with color (r,g,b)
call set_pixel(img, x, y, rgb(1), rgb(2), rgb(3))
@ -228,6 +228,7 @@ contains
integer :: id ! id of cell or material
real(8) :: vox(3) ! x, y, and z voxel widths
real(8) :: ll(3) ! lower left starting point for each sweep direction
type(Particle) :: p
! compute voxel widths in each direction
vox = pl % width/dble(pl % pixels)
@ -236,8 +237,7 @@ contains
ll = pl % origin - pl % width / 2.0
! allocate and initialize particle
allocate(p)
call initialize_particle()
call p % initialize()
p % coord0 % xyz = ll
p % coord0 % uvw = (/ 0.5, 0.5, 0.5 /)
p % coord0 % universe = BASE_UNIVERSE
@ -257,7 +257,7 @@ contains
do z = 1, pl % pixels(3)
! get voxel color
call position_rgb(pl, rgb, id)
call position_rgb(p, pl, rgb, id)
! write to plot file
write(UNIT_PLOT) id

51
src/source.F90
View file

@ -5,9 +5,9 @@ module source
use error, only: fatal_error
use geometry_header, only: BASE_UNIVERSE
use global
use math, only: maxwell_spectrum, watt_spectrum
use output, only: write_message
use particle_header, only: deallocate_coord
use physics, only: maxwell_spectrum, watt_spectrum
use particle_header, only: Particle
use random_lcg, only: prn, set_particle_seed
use string, only: to_str
@ -149,19 +149,20 @@ contains
! GET_SOURCE_PARTICLE returns the next source particle
!===============================================================================
subroutine get_source_particle(index_source)
subroutine get_source_particle(p, index_source)
integer(8), intent(in) :: index_source
type(Particle), intent(inout) :: p
integer(8), intent(in) :: index_source
integer(8) :: particle_seed ! unique index for particle
type(Bank), pointer :: src => null()
! set defaults
call initialize_particle()
call p % initialize()
! Copy attributes from source to particle
src => source_bank(index_source)
call copy_source_attributes(src)
call copy_source_attributes(p, src)
! set identifier for particle
p % id = bank_first + index_source - 1
@ -181,9 +182,10 @@ contains
! COPY_SOURCE_ATTRIBUTES
!===============================================================================
subroutine copy_source_attributes(src)
subroutine copy_source_attributes(p, src)
type(Bank), pointer :: src
type(Particle), intent(inout) :: p
type(Bank), pointer :: src
! copy attributes from source bank site
p % wgt = src % wgt
@ -196,37 +198,4 @@ contains
end subroutine copy_source_attributes
!===============================================================================
! INITIALIZE_PARTICLE sets default attributes for a particle from the source
! bank
!===============================================================================
subroutine initialize_particle()
! Set particle to neutron that's alive
p % type = NEUTRON
p % alive = .true.
! clear attributes
p % surface = NONE
p % cell_born = NONE
p % material = NONE
p % last_material = NONE
p % wgt = ONE
p % last_wgt = ONE
p % absorb_wgt = ZERO
p % n_bank = 0
p % wgt_bank = ZERO
p % n_collision = 0
! remove any original coordinates
call deallocate_coord(p % coord0)
! Set up base level coordinates
allocate(p % coord0)
p % coord0 % universe = BASE_UNIVERSE
p % coord => p % coord0
end subroutine initialize_particle
end module source

53
src/tally.F90
View file

@ -10,7 +10,7 @@ module tally
mesh_intersects_2d, mesh_intersects_3d
use mesh_header, only: StructuredMesh
use output, only: header
use particle_header, only: LocalCoord
use particle_header, only: LocalCoord, Particle
use search, only: binary_search
use string, only: to_str
use tally_header, only: TallyResult, TallyMapItem, TallyMapElement
@ -32,7 +32,9 @@ contains
! triggered at every collision, not every event
!===============================================================================
subroutine score_analog_tally()
subroutine score_analog_tally(p)
type(Particle), intent(in) :: p
integer :: i
integer :: i_tally
@ -70,7 +72,7 @@ contains
! =======================================================================
! DETERMINE SCORING BIN COMBINATION
call get_scoring_bins(i_tally, found_bin)
call get_scoring_bins(p, i_tally, found_bin)
if (.not. found_bin) cycle
! =======================================================================
@ -285,7 +287,7 @@ contains
! outgoing energy bins may have been scored to. The following
! logic treats this special case and results to multiple bins
call score_fission_eout(t, score_index)
call score_fission_eout(p, t, score_index)
cycle SCORE_LOOP
else
@ -364,8 +366,9 @@ contains
! neutrons produced with different energies.
!===============================================================================
subroutine score_fission_eout(t, i_score)
subroutine score_fission_eout(p, t, i_score)
type(Particle), intent(in) :: p
type(TallyObject), pointer :: t
integer, intent(in) :: i_score ! index for score
@ -424,9 +427,10 @@ contains
! that require post-collision information.
!===============================================================================
subroutine score_tracklength_tally(distance)
subroutine score_tracklength_tally(p, distance)
real(8), intent(in) :: distance
type(Particle), intent(in) :: p
real(8), intent(in) :: distance
integer :: i
integer :: i_tally
@ -464,14 +468,14 @@ contains
! since multiple bins can be scored to with a single track
if (t % find_filter(FILTER_MESH) > 0) then
call score_tl_on_mesh(i_tally, distance)
call score_tl_on_mesh(p, i_tally, distance)
cycle
end if
! =======================================================================
! DETERMINE SCORING BIN COMBINATION
call get_scoring_bins(i_tally, found_bin)
call get_scoring_bins(p, i_tally, found_bin)
if (.not. found_bin) cycle
! =======================================================================
@ -485,7 +489,7 @@ contains
filter_index = sum((t % matching_bins - 1) * t % stride) + 1
if (t % all_nuclides) then
call score_all_nuclides(i_tally, flux, filter_index)
call score_all_nuclides(p, i_tally, flux, filter_index)
else
NUCLIDE_BIN_LOOP: do k = 1, t % n_nuclide_bins
@ -721,11 +725,12 @@ contains
! the user requests <nuclides>all</nuclides>.
!===============================================================================
subroutine score_all_nuclides(i_tally, flux, filter_index)
subroutine score_all_nuclides(p, i_tally, flux, filter_index)
integer, intent(in) :: i_tally
real(8), intent(in) :: flux
integer, intent(in) :: filter_index
type(Particle), intent(in) :: p
integer, intent(in) :: i_tally
real(8), intent(in) :: flux
integer, intent(in) :: filter_index
integer :: i ! loop index for nuclides in material
integer :: j ! loop index for scoring bin types
@ -946,10 +951,11 @@ contains
! these tallies, it is possible to score to multiple mesh cells for each track.
!===============================================================================
subroutine score_tl_on_mesh(i_tally, d_track)
subroutine score_tl_on_mesh(p, i_tally, d_track)
integer, intent(in) :: i_tally
real(8), intent(in) :: d_track
type(Particle), intent(in) :: p
integer, intent(in) :: i_tally
real(8), intent(in) :: d_track
integer :: i ! loop index for filter/score bins
integer :: j ! loop index for direction
@ -1145,7 +1151,7 @@ contains
if (t % all_nuclides) then
! Score reaction rates for each nuclide in material
call score_all_nuclides(i_tally, flux, filter_index)
call score_all_nuclides(p, i_tally, flux, filter_index)
else
NUCLIDE_BIN_LOOP: do b = 1, t % n_nuclide_bins
@ -1260,10 +1266,11 @@ contains
! for a tally based on the particle's current attributes.
!===============================================================================
subroutine get_scoring_bins(i_tally, found_bin)
subroutine get_scoring_bins(p, i_tally, found_bin)
integer, intent(in) :: i_tally
logical, intent(out) :: found_bin
type(Particle), intent(in) :: p
integer, intent(in) :: i_tally
logical, intent(out) :: found_bin
integer :: i ! loop index for filters
integer :: n ! number of bins for single filter
@ -1370,7 +1377,9 @@ contains
! determining which mesh surfaces were crossed
!===============================================================================
subroutine score_surface_current()
subroutine score_surface_current(p)
type(Particle), intent(in) :: p
integer :: i
integer :: i_tally

191
src/tracking.F90 Normal file
View file

@ -0,0 +1,191 @@
module tracking
use cross_section, only: calculate_xs
use error, only: fatal_error, warning
use geometry, only: find_cell, distance_to_boundary, cross_surface, &
cross_lattice, check_cell_overlap
use geometry_header, only: Universe, BASE_UNIVERSE
use global
use output, only: write_message
use particle_header, only: LocalCoord, Particle
use physics, only: collision
use random_lcg, only: prn
use string, only: to_str
use tally, only: score_analog_tally, score_tracklength_tally, &
score_surface_current
contains
!===============================================================================
! TRANSPORT encompasses the main logic for moving a particle through geometry.
!===============================================================================
subroutine transport(p)
type(Particle), intent(inout) :: p
integer :: surface_crossed ! surface which particle is on
integer :: lattice_crossed ! lattice boundary which particle crossed
integer :: last_cell ! most recent cell particle was in
integer :: n_event ! number of collisions/crossings
real(8) :: d_boundary ! distance to nearest boundary
real(8) :: d_collision ! sampled distance to collision
real(8) :: distance ! distance particle travels
logical :: found_cell ! found cell which particle is in?
type(LocalCoord), pointer :: coord => null()
! Display message if high verbosity or trace is on
if (verbosity >= 9 .or. trace) then
message = "Simulating Particle " // trim(to_str(p % id))
call write_message()
end if
! If the cell hasn't been determined based on the particle's location,
! initiate a search for the current cell
if (p % coord % cell == NONE) then
call find_cell(p, found_cell)
! Particle couldn't be located
if (.not. found_cell) then
message = "Could not locate particle " // trim(to_str(p % id))
call fatal_error()
end if
! set birth cell attribute
p % cell_born = p % coord % cell
end if
! Initialize number of events to zero
n_event = 0
! Add paricle's starting weight to count for normalizing tallies later
total_weight = total_weight + p % wgt
! Force calculation of cross-sections by setting last energy to zero
micro_xs % last_E = ZERO
do while (p % alive)
if (check_overlaps) call check_cell_overlap(p)
! Calculate microscopic and macroscopic cross sections -- note: if the
! material is the same as the last material and the energy of the
! particle hasn't changed, we don't need to lookup cross sections again.
if (p % material /= p % last_material) call calculate_xs(p)
! Find the distance to the nearest boundary
call distance_to_boundary(p, d_boundary, surface_crossed, lattice_crossed)
! Sample a distance to collision
if (material_xs % total == ZERO) then
d_collision = INFINITY
else
d_collision = -log(prn()) / material_xs % total
end if
! Select smaller of the two distances
distance = min(d_boundary, d_collision)
! Advance particle
coord => p % coord0
do while (associated(coord))
coord % xyz = coord % xyz + distance * coord % uvw
coord => coord % next
end do
! Score track-length tallies
if (active_tracklength_tallies % size() > 0) &
call score_tracklength_tally(p, distance)
! Score track-length estimate of k-eff
global_tallies(K_TRACKLENGTH) % value = &
global_tallies(K_TRACKLENGTH) % value + p % wgt * distance * &
material_xs % nu_fission
if (d_collision > d_boundary) then
! ====================================================================
! PARTICLE CROSSES SURFACE
last_cell = p % coord % cell
p % coord % cell = NONE
if (lattice_crossed /= NONE) then
! Particle crosses lattice boundary
p % surface = NONE
call cross_lattice(p, lattice_crossed)
p % event = EVENT_LATTICE
else
! Particle crosses surface
p % surface = surface_crossed
call cross_surface(p, last_cell)
p % event = EVENT_SURFACE
end if
else
! ====================================================================
! PARTICLE HAS COLLISION
! Score collision estimate of keff
global_tallies(K_COLLISION) % value = &
global_tallies(K_COLLISION) % value + p % wgt * &
material_xs % nu_fission / material_xs % total
! 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
! pre-collision direction to figure out what mesh surfaces were crossed
if (active_current_tallies % size() > 0) call score_surface_current(p)
! Clear surface component
p % surface = NONE
call collision(p)
! Score collision estimator tallies -- this is done after a collision
! has occurred rather than before because we need information on the
! outgoing energy for any tallies with an outgoing energy filter
if (active_analog_tallies % size() > 0) call score_analog_tally(p)
! Reset banked weight during collision
p % n_bank = 0
p % wgt_bank = ZERO
! Save coordinates for tallying purposes
p % last_xyz = p % coord0 % xyz
! Set last material to none since cross sections will need to be
! re-evaluated
p % last_material = NONE
! Set all uvws to base level -- right now, after a collision, only the
! base level uvws are changed
coord => p % coord0
do while(associated(coord % next))
if (coord % next % rotated) then
! If next level is rotated, apply rotation matrix
coord % next % uvw = matmul(cells(coord % cell) % &
rotation, coord % uvw)
else
! Otherwise, copy this level's direction
coord % next % uvw = coord % uvw
end if
! Advance coordinate level
coord => coord % next
end do
end if
! If particle has too many events, display warning and kill it
n_event = n_event + 1
if (n_event == MAX_EVENTS) then
message = "Particle " // trim(to_str(p%id)) // " underwent maximum &
&number of events."
call warning()
p % alive = .false.
end if
end do
end subroutine transport
end module tracking

4
src/utils/build_dependencies.py
View file

@ -9,9 +9,9 @@ for src in glob.iglob('*.F90'):
module = src.strip('.F90')
deps = set()
d = re.findall(r'\n\s*use\s+(\w+)',
open(src,'r').read())
open(src, 'r').read())
for name in d:
if name in ['mpi','hdf5','h5lt','fox_dom']:
if name in ['mpi', 'hdf5', 'h5lt', 'fox_dom']:
continue
deps.add(name)
if deps:

39
src/utils/convert_xsdata.py
View file

@ -6,6 +6,7 @@ from xml.dom.minidom import getDOMImplementation
types = {1: "neutron", 2: "dosimetry", 3: "thermal"}
class Xsdata(object):
def __init__(self, filename):
@ -94,7 +95,7 @@ class XsdataTable(object):
continue
string = " ".join(self.alias)
else:
string = "{0}".format(getattr(self,attribute))
string = "{0}".format(getattr(self, attribute))
# Skip metastable and binary if 0
if attribute == "metastable" and self.metastable == 0:
@ -109,7 +110,7 @@ class XsdataTable(object):
# node.appendChild(nodeAttr)
node.setAttribute(attribute, string)
return node
if __name__ == '__main__':
# Read command line arguments
@ -121,25 +122,25 @@ if __name__ == '__main__':
# Read xsdata and create XML document object
xsdataObject = Xsdata(xsdataFile)
doc = xsdataObject.to_xml()
# Reduce number of lines
lines = doc.toprettyxml(indent=' ')
lines = lines.replace('<alias>\n ','<alias>')
lines = lines.replace('\n </alias>','</alias>')
lines = lines.replace('<zaid>\n ','<zaid>')
lines = lines.replace('\n </zaid>','</zaid>')
lines = lines.replace('<type>\n ','<type>')
lines = lines.replace('\n </type>','</type>')
lines = lines.replace('<awr>\n ','<awr>')
lines = lines.replace('\n </awr>','</awr>')
lines = lines.replace('<temperature>\n ','<temperature>')
lines = lines.replace('\n </temperature>','</temperature>')
lines = lines.replace('<path>\n ','<path>')
lines = lines.replace('\n </path>','</path>')
lines = lines.replace('<metastable>\n ','<metastable>')
lines = lines.replace('\n </metastable>','</metastable>')
lines = lines.replace('<binary>\n ','<binary>')
lines = lines.replace('\n </binary>','</binary>')
lines = lines.replace('<alias>\n ', '<alias>')
lines = lines.replace('\n </alias>', '</alias>')
lines = lines.replace('<zaid>\n ', '<zaid>')
lines = lines.replace('\n </zaid>', '</zaid>')
lines = lines.replace('<type>\n ', '<type>')
lines = lines.replace('\n </type>', '</type>')
lines = lines.replace('<awr>\n ', '<awr>')
lines = lines.replace('\n </awr>', '</awr>')
lines = lines.replace('<temperature>\n ', '<temperature>')
lines = lines.replace('\n </temperature>', '</temperature>')
lines = lines.replace('<path>\n ', '<path>')
lines = lines.replace('\n </path>', '</path>')
lines = lines.replace('<metastable>\n ', '<metastable>')
lines = lines.replace('\n </metastable>', '</metastable>')
lines = lines.replace('<binary>\n ', '<binary>')
lines = lines.replace('\n </binary>', '</binary>')
# Write document in pretty XML to specified file
f = open(xmlFile, 'w')

18
src/utils/convert_xsdir.py
View file

@ -16,6 +16,7 @@ elements = [None, "H", "He", "Li", "Be", "B", "C", "N", "O", "F", "Ne", "Na",
"Md", "No", "Lr", "Rf", "Db", "Sg", "Bh", "Hs", "Mt", "Ds", "Rg",
"Cn"]
class Xsdir(object):
def __init__(self, filename):
@ -54,14 +55,14 @@ class Xsdir(object):
# Check for end of second section
if len(words) % 2 != 0 or words[0] == 'directory':
break
for zaid, awr in zip(words[::2], words[1::2]):
self.awr[zaid] = awr
# Read third section
while words[0] != 'directory':
words = self.f.readline().split()
while True:
words = self.f.readline().split()
if not words:
@ -76,7 +77,7 @@ class Xsdir(object):
# Create XsdirTable object and add to line
table = XsdirTable(self.directory)
self.tables.append(table)
# All tables have at least 7 attributes
table.name = words[0]
table.awr = float(words[1])
@ -115,7 +116,6 @@ class Xsdir(object):
self.entries = list(self.entries)[0]
else:
self.recordlength = None
def to_xml(self):
# Create XML document
@ -154,7 +154,7 @@ class Xsdir(object):
# Add a node for each table
for table in self.tables:
if table.name[-1] in ['e', 'p', 'u', 'h', 'g' ,'m', 'd']:
if table.name[-1] in ['e', 'p', 'u', 'h', 'g', 'm', 'd']:
continue
node = table.to_xml_node(doc)
root.appendChild(node)
@ -230,7 +230,7 @@ class XsdirTable(object):
return "{0}-{1}.{2}".format(elements[Z], s, self.xs)
else:
return None
@property
def zaid(self):
if self.name.endswith('c'):
@ -245,10 +245,10 @@ class XsdirTable(object):
def to_xml_node(self, doc):
node = doc.createElement("ace_table")
node.setAttribute("name", self.name)
for attribute in ["alias", "zaid", "type", "metastable", "awr",
for attribute in ["alias", "zaid", "type", "metastable", "awr",
"temperature", "path", "location"]:
if hasattr(self, attribute):
string = str(getattr(self,attribute))
string = str(getattr(self, attribute))
# Skip metastable and binary if 0
if attribute == "metastable" and self.metastable == 0:
@ -274,7 +274,7 @@ if __name__ == '__main__':
# Read xsdata and create XML document object
xsdirObject = Xsdir(xsdirFile)
doc = xsdirObject.to_xml()
# Reduce number of lines
lines = doc.toprettyxml(indent=' ')

1
src/utils/geometry.py
View file

@ -2,6 +2,7 @@
from xml.dom.minidom import parse
class Geometry(object):
def __init__(self, filename):

4
src/utils/memory_usage.py
View file

@ -1,7 +1,7 @@
#!/usr/bin/env python2
# This script reads a cross_sections.out file, adds up the memory usage for each
# nuclide and S(a,b) table, and displays the total memory usage
# This script reads a cross_sections.out file, adds up the memory usage for
# each nuclide and S(a,b) table, and displays the total memory usage
import sys
import os

73
src/utils/particle_restart.py Normal file
View file

@ -0,0 +1,73 @@
#!/usr/bin/env python2
import struct
class Particle(object):
def __init__(self, filename):
if filename.endswith('.h5'):
import h5py
self._f = h5py.File(filename, 'r')
self._hdf5 = True
else:
self._f = open(filename, 'rb')
self._hdf5 = False
# Read all metadata
self._read_data()
def _read_data(self):
# Read filetype
self.filetype = self._get_int(path='filetype')[0]
# Read statepoint revision
self.revision = self._get_int(path='revision')[0]
# Read current batch
self.current_batch = self._get_int(path='current_batch')[0]
# Read run information
self.gen_per_batch = self._get_int(path='gen_per_batch')[0]
self.current_gen = self._get_int(path='current_gen')[0]
self.n_particles = self._get_long(path='n_particles')[0]
# Read particle properties
self.id = self._get_long(path='id')[0]
self.weight = self._get_double(path='weight')[0]
self.energy = self._get_double(path='energy')[0]
self.xyz = self._get_double(3, path='xyz')
self.uvw = self._get_double(3, path='uvw')
def _get_data(self, n, typeCode, size):
return list(struct.unpack('={0}{1}'.format(n, typeCode),
self._f.read(n*size)))
def _get_int(self, n=1, path=None):
if self._hdf5:
return [int(v) for v in self._f[path].value]
else:
return [int(v) for v in self._get_data(n, 'i', 4)]
def _get_long(self, n=1, path=None):
if self._hdf5:
return [long(v) for v in self._f[path].value]
else:
return [long(v) for v in self._get_data(n, 'q', 8)]
def _get_float(self, n=1, path=None):
if self._hdf5:
return [float(v) for v in self._f[path].value]
else:
return [float(v) for v in self._get_data(n, 'f', 4)]
def _get_double(self, n=1, path=None):
if self._hdf5:
return [float(v) for v in self._f[path].value]
else:
return [float(v) for v in self._get_data(n, 'd', 8)]
def _get_string(self, n=1, path=None):
if self._hdf5:
return str(self._f[path].value)
else:
return str(self._get_data(n, 's', 1)[0])

8
src/utils/split_nuclide.py
View file

@ -4,15 +4,13 @@ import sys
if len(sys.argv) != 3:
print("Must supply element and atom/b-cm")
exit
sys.exit(1)
element = sys.argv[1]
ao = float(sys.argv[2])
for line in open('abundances_modified.txt','r'):
for line in open('abundances_modified.txt', 'r'):
words = line.split()
if words[1] == element:
print('<nuclide name="{0}-{1}" ao="{2:10.4e}" />'.format(
element, words[2], float(words[3])*ao))
element, words[2], float(words[3])*ao))

9
src/utils/statepoint_histogram.py
View file

@ -56,14 +56,15 @@ for i, t in enumerate(sp.tallies):
relative_error = t_value*sqrt((s2/n - s*s)/(n-1))/s
uncertainties.append(relative_error)
# Display number of tallies with less than 1% CIs
fraction = (sum([1.0 if re < 0.01 else 0.0 for re in uncertainties])
/n_bins)
# Display number of tallies with less than 1% CIs
fraction = (sum([1.0 if re < 0.01 else 0.0 for re in uncertainties])
/ n_bins)
print("Tally " + str(i+1))
print(" Fraction under 1% = {0}".format(fraction))
print(" Min relative error = {0}".format(min(uncertainties)))
print(" Max relative error = {0}".format(max(uncertainties)))
print(" Non-scoring bins = {0}".format(1.0 - float(len(uncertainties))/n_bins))
print(" Non-scoring bins = {0}".format(
1.0 - float(len(uncertainties))/n_bins))
# Plot histogram
plt.hist(uncertainties, 100)

24
src/utils/statepoint_meshplot.py
View file

@ -60,9 +60,9 @@ for t in sp.tallies:
assert n_bins == nx*ny*nz*ns
# Create lists for tallies
mean = np.zeros((nx,ny))
error = np.zeros((nx,ny))
criteria = np.zeros((nx,ny))
mean = np.zeros((nx, ny))
error = np.zeros((nx, ny))
criteria = np.zeros((nx, ny))
# Determine starting position for data
start = sp._f.tell() + (axial_level-1)*ns*16 + (score - 1)*16
@ -75,19 +75,19 @@ for t in sp.tallies:
# Read sum and sum-squared
s, s2 = sp._get_double(2)
s /= n
mean[x,y] = s
mean[x, y] = s
if s != 0.0:
error[x,y] = t_value*sqrt((s2/n - s*s)/(n-1))/s
criteria[x,y] = 1.0 if error[x,y] < 0.05 else 0.0
error[x, y] = t_value*sqrt((s2/n - s*s)/(n-1))/s
criteria[x, y] = 1.0 if error[x, y] < 0.05 else 0.0
# Make figure
print("Making colorplot...")
plt.imshow(mean, interpolation="nearest")
plt.colorbar()
plt.xlim((0,nx))
plt.ylim((0,ny))
plt.xticks(np.linspace(0,nx,5),
np.linspace(m.lower_left[0],m.upper_right[0],5))
plt.yticks(np.linspace(0,ny,5),
np.linspace(m.lower_left[1],m.upper_right[1],5))
plt.xlim((0, nx))
plt.ylim((0, ny))
plt.xticks(np.linspace(0, nx, 5),
np.linspace(m.lower_left[0], m.upper_right[0], 5))
plt.yticks(np.linspace(0, ny, 5),
np.linspace(m.lower_left[1], m.upper_right[1], 5))
plt.show()

22
tests/test_basic/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_basic/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
0.30408759 0.00088008

10
tests/test_basic/test_basic.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -18,8 +19,15 @@ def test_run():
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_confidence_intervals/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

5
tests/test_confidence_intervals/results_true.dat Normal file
View file

@ -0,0 +1,5 @@
k-combined:
0.29432037 0.01265543
tallies:
65.39238423
539.83246478

11
tests/test_confidence_intervals/test_confidence_intervals.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_created_statepoint():
def test_created_output():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

22
tests/test_density_atombcm/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_density_atombcm/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
1.75925969 0.02728159

13
tests/test_density_atombcm/test_density_atombcm.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -15,8 +16,18 @@ def test_run():
print(proc.communicate()[0])
assert returncode == 0
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

22
tests/test_density_atomcm3/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_density_atomcm3/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
1.12562497 0.01475256

13
tests/test_density_atomcm3/test_density_atomcm3.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -15,8 +16,18 @@ def test_run():
print(proc.communicate()[0])
assert returncode == 0
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

22
tests/test_density_kgm3/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_density_kgm3/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
0.77607961 0.02228188

13
tests/test_density_kgm3/test_density_kgm3.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -15,8 +16,18 @@ def test_run():
print(proc.communicate()[0])
assert returncode == 0
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

22
tests/test_density_sum/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_density_sum/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
0.32600067 0.01356562

10
tests/test_density_sum/test_density_sum.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -18,8 +19,15 @@ def test_run():
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

22
tests/test_eigenvalue_genperbatch/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.7.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_eigenvalue_genperbatch/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
0.28402818 0.02273104

10
tests/test_eigenvalue_genperbatch/test_eigenvalue_genperbatch.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -18,8 +19,15 @@ def test_run():
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.7.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.7.binary']
output = [pwd + '/statepoint.7.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

22
tests/test_eigenvalue_no_inactive/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_eigenvalue_no_inactive/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
0.30896828 0.01118024

10
tests/test_eigenvalue_no_inactive/test_eigenvalue_no_inactive.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -18,8 +19,15 @@ def test_run():
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

22
tests/test_energy_grid/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_energy_grid/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
0.31522753 0.00509959

13
tests/test_energy_grid/test_energy_grid.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -15,8 +16,18 @@ def test_run():
print(proc.communicate()[0])
assert returncode == 0
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

28
tests/test_entropy/results.py Normal file
View file

@ -0,0 +1,28 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'entropy:\n'
for item in sp.entropy:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

13
tests/test_entropy/results_true.dat Normal file
View file

@ -0,0 +1,13 @@
k-combined:
0.30408759 0.00088008
entropy:
6.22453002
7.34335892
7.76336209
7.76885034
7.78607726
7.72581193
7.76121795
7.84101403
7.92560054
7.86038189

10
tests/test_entropy/test_entropy.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -18,8 +19,15 @@ def test_run():
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_filter_cell/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

11
tests/test_filter_cell/results_true.dat Normal file
View file

@ -0,0 +1,11 @@
k-combined:
0.98046852 0.04718832
tallies:
0.00000000
0.00000000
8.78237650
15.55794571
1.84482312
0.70297682
26.52660205
142.33423161

11
tests/test_filter_cell/test_filter_cell.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_statepoint_exists():
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_filter_cellborn/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

11
tests/test_filter_cellborn/results_true.dat Normal file
View file

@ -0,0 +1,11 @@
k-combined:
0.98046852 0.04718832
tallies:
0.00000000
0.00000000
43.50660698
385.10698472
0.00000000
0.00000000
0.00000000
0.00000000

11
tests/test_filter_cellborn/test_filter_cellborn.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_statepoint_exists():
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_filter_energy/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

11
tests/test_filter_energy/results_true.dat Normal file
View file

@ -0,0 +1,11 @@
k-combined:
0.98046852 0.04718832
tallies:
25.67612478
135.06234071
43.69328739
384.31331185
53.66385216
577.08087569
9.12440078
17.45738591

11
tests/test_filter_energy/test_filter_energy.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_statepoint_exists():
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_filter_energyout/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

11
tests/test_filter_energyout/results_true.dat Normal file
View file

@ -0,0 +1,11 @@
k-combined:
0.98046852 0.04718832
tallies:
25.59000000
133.54390000
42.96000000
370.54460000
51.41000000
529.20190000
6.14000000
7.84340000

11
tests/test_filter_energyout/test_filter_energyout.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_statepoint_exists():
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_filter_group_transfer/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

67
tests/test_filter_group_transfer/results_true.dat Normal file
View file

@ -0,0 +1,67 @@
k-combined:
0.98046852 0.04718832
tallies:
22.97000000
108.08230000
0.00000000
0.00000000
0.06000000
0.00140000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
2.62000000
1.37720000
0.00000000
0.00000000
38.33000000
295.26630000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
4.56000000
4.15920000
0.00000000
0.00000000
48.48000000
470.68580000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.01000000
0.00010000
0.00000000
0.00000000
2.93000000
1.78470000
0.00000000
0.00000000
6.14000000
7.84340000
0.00000000
0.00000000

11
tests/test_filter_group_transfer/test_filter_group_transfer.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_statepoint_exists():
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_filter_material/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

11
tests/test_filter_material/results_true.dat Normal file
View file

@ -0,0 +1,11 @@
k-combined:
0.98046852 0.04718832
tallies:
28.88595756
167.70454782
6.68679307
8.94919165
35.19283765
254.25938930
60.96487727
750.99609930

11
tests/test_filter_material/test_filter_material.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_statepoint_exists():
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_filter_mesh_2d/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

581
tests/test_filter_mesh_2d/results_true.dat Normal file
View file

@ -0,0 +1,581 @@
k-combined:
0.98046852 0.04718832
tallies:
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.08588844
0.00737682
1.29395418
0.61345354
3.17605552
2.84377280
2.97221609
2.28098649
2.26795333
1.83191774
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.21020517
0.04418622
0.00000000
0.00000000
0.00000000
0.00000000
0.75237153
0.28460845
0.12767750
0.01563459
1.09635635
0.81124602
1.31009450
0.78424448
4.48862660
4.40176027
2.61476684
2.29263337
2.73483247
2.44147576
0.00222491
0.00000495
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.09474020
0.00897571
0.49862162
0.24862352
0.00000000
0.00000000
0.00000000
0.00000000
0.18606180
0.03461899
0.06427301
0.00413102
1.19911234
0.45648739
1.42142705
0.53678749
3.60239192
2.89778626
0.13297308
0.01768184
1.12151012
0.32281932
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00273377
0.00000747
0.87406271
0.30888052
0.39960895
0.08005836
2.30080722
2.04968651
1.66625586
0.78412054
0.26891589
0.03246764
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.24152882
0.05833617
0.25632992
0.03814006
4.65349777
5.94071891
2.10515795
1.52728004
2.58079707
1.55536111
3.50569432
2.61321114
2.11589300
2.23703372
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.58991041
0.34799429
0.00000000
0.00000000
0.15104665
0.02281509
0.42889841
0.18395384
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.35426490
0.12550362
0.75912683
0.15788880
1.49035861
0.65263822
1.54069090
1.58982707
0.05440226
0.00295961
0.00000000
0.00000000
0.06189472
0.00383096
0.07567228
0.00572629
0.76895320
0.55973718
1.82808885
1.17843830
0.31960577
0.10214785
0.04014117
0.00161131
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00869849
0.00007566
0.19463946
0.03788452
0.19061650
0.03633465
0.35525829
0.12620845
0.00000000
0.00000000
0.00000000
0.00000000
0.14948493
0.02234575
0.26642025
0.04919614
0.18119903
0.01648042
0.33066766
0.04503538
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.10984979
0.01206698
0.30039876
0.05912935
0.04004197
0.00160336
0.17230860
0.02969025
0.57320487
0.17068193
0.00000000
0.00000000
0.00000000
0.00000000
0.04772221
0.00227741
0.01276581
0.00009446
1.37698488
0.71236268
0.11787532
0.01389459
0.00849153
0.00007211
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.02731512
0.00074612
0.67798922
0.21862908
0.56975986
0.07959617
0.90520028
0.37868800
3.15953485
3.14097026
1.22452698
0.53508531
0.03235360
0.00104676
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
1.67247656
0.94402049
0.32669397
0.03903995
0.22350340
0.04995377
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.24456194
0.05981054
3.90325239
4.14221991
3.08364316
2.45517354
2.32748465
1.52115323
4.52477924
5.12896423
2.64295725
1.74944646
0.66721965
0.14668973
0.00000000
0.00000000
0.00000000
0.00000000
1.43623311
0.70640842
0.25035386
0.03446661
1.24322754
0.52709506
1.75784221
1.65092830
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.03373481
0.00113804
1.57605576
0.80561508
2.47246255
1.77773457
2.60036982
1.55878535
3.21120690
2.33678026
4.24458477
4.44119713
0.91438183
0.33422901
0.60473444
0.18557226
0.00000000
0.00000000
1.44569376
0.75022322
0.99247413
0.31036171
0.84162339
0.28807337
0.00270497
0.00000732
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.06002436
0.00360292
1.01351346
0.54569143
0.20762778
0.01641673
0.57527086
0.22871026
0.18057558
0.03011038
1.00598096
0.45972185
0.88036591
0.37710745
0.53086453
0.14090904
0.17632445
0.03109031
0.00000000
0.00000000
0.78089598
0.32633537
2.81720674
2.43564074
0.00942578
0.00008885
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.56039932
0.31404739
0.00450316
0.00002028
0.00000000
0.00000000
0.00000000
0.00000000
0.52140618
0.11042605
1.91904108
1.34050514
0.10411834
0.00745036
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.02751492
0.00075707
0.99249418
0.32868074
0.00856223
0.00007331
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.38273195
0.10389392
0.01182298
0.00013978
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000
0.00000000

11
tests/test_filter_mesh_2d/test_filter_mesh_2d.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_statepoint_exists():
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_filter_mesh_3d/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

9829
tests/test_filter_mesh_3d/results_true.dat Normal file
View file

File diff suppressed because it is too large Load diff

11
tests/test_filter_mesh_3d/test_filter_mesh_3d.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_statepoint_exists():
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

34
tests/test_filter_universe/results.py Normal file
View file

@ -0,0 +1,34 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

11
tests/test_filter_universe/results_true.dat Normal file
View file

@ -0,0 +1,11 @@
k-combined:
0.98046852 0.04718832
tallies:
36.02063446
262.67869132
4.68789690
4.50704437
74.02502234
1111.88560685
8.53811818
14.72177988

11
tests/test_filter_universe/test_filter_universe.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -21,8 +22,16 @@ def test_statepoint_exists():
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

30
tests/test_fixed_source/results.py Normal file
View file

@ -0,0 +1,30 @@
#!/usr/bin/env python
import sys
import numpy as np
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# extract tally results and convert to vector
results = sp.tallies[0].results
shape = results.shape
size = (np.product(shape))
results = np.reshape(results, size)
# set up output string
outstr = ''
# write out tally results
outstr += 'tallies:\n'
for item in results:
outstr += "{0:10.8f}\n".format(item)
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

3
tests/test_fixed_source/results_true.dat Normal file
View file

@ -0,0 +1,3 @@
tallies:
453.85329674
20667.76024307

8
tests/test_fixed_source/tallies.xml Normal file
View file

@ -0,0 +1,8 @@
<?xml version="1.0"?>
<tallies>
<tally id="1">
<scores>flux</scores>
</tally>
</tallies>

16
tests/test_fixed_source/test_fixed_source.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -15,11 +16,22 @@ def test_run():
print(proc.communicate()[0])
assert returncode == 0
def test_created_statepoint():
def test_statepoint_exists():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_output_exists():
assert os.path.exists(pwd + '/tallies.out')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/tallies.out',
pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

20
tests/test_lattice/plots.xml
View file

@ -1,20 +0,0 @@
<?xml version="1.0"?>
<plots>
<!--
===============================================================
Description: Rectangular Array of Water-Moderated MOX Fuel Rods
Case: MIX-COMP-THERM-002, PNL-30
Written By: Paul Romano
Date: 9/6/2012
===============================================================
-->
<plot id="1" color="mat">
<basis>xy</basis>
<origin>25. 25. 4.5</origin>
<width>50. 50.</width>
<pixels>1000 1000</pixels>
</plot>
</plots>

22
tests/test_lattice/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_lattice/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
0.91525327 0.06208235

13
tests/test_lattice/test_lattice.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -15,8 +16,18 @@ def test_run():
print(proc.communicate()[0])
assert returncode == 0
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

22
tests/test_lattice_multiple/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_lattice_multiple/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
0.98046852 0.04718832

13
tests/test_lattice_multiple/test_lattice_multiple.py
View file

@ -2,6 +2,7 @@
import os
from subprocess import Popen, STDOUT, PIPE
import filecmp
pwd = os.path.dirname(__file__)
@ -15,8 +16,18 @@ def test_run():
print(proc.communicate()[0])
assert returncode == 0
def test_created_statepoint():
assert os.path.exists(pwd + '/statepoint.10.binary')
def test_results():
os.system('python results.py')
compare = filecmp.cmp('results_test.dat', 'results_true.dat')
if not compare:
os.rename('results_test.dat', 'results_error.dat')
assert compare
def teardown():
output = [pwd + '/statepoint.10.binary']
output = [pwd + '/statepoint.10.binary', pwd + '/results_test.dat']
for f in output:
if os.path.exists(f):
os.remove(f)

22
tests/test_natural_element/results.py Normal file
View file

@ -0,0 +1,22 @@
#!/usr/bin/env python
import sys
# import statepoint
sys.path.append('../../src/utils')
import statepoint
# read in statepoint file
sp = statepoint.StatePoint('statepoint.10.binary')
sp.read_results()
# set up output string
outstr = ''
# write out k-combined
outstr += 'k-combined:\n'
outstr += "{0:10.8f} {1:10.8f}\n".format(sp.k_combined[0], sp.k_combined[1])
# write results to file
with open('results_test.dat','w') as fh:
fh.write(outstr)

2
tests/test_natural_element/results_true.dat Normal file
View file

@ -0,0 +1,2 @@
k-combined:
0.96912628 0.04173669

Some files were not shown because too many files have changed in this diff Show more