Merge remote-tracking branch 'upstream/develop' into hex_lattice

Conflicts:
	src/relaxng/geometry.rnc

docs/source/devguide/statepoint.rst

@@ -4,6 +4,292 @@
 State Point Binary File Specifications
 ======================================
 
+-----------
+Revision 13
+-----------
+
+**integer(4) FILETYPE_STATEPOINT**
+
+    Flags whether this file is a statepoint file or a particle restart file.
+
+**integer(4) REVISION_STATEPOINT**
+
+    Revision of the binary state point file. Any time a change is made in the
+    format of the state-point file, this integer is incremented.
+
+**integer(4) VERSION_MAJOR**
+
+    Major version number for OpenMC
+
+**integer(4) VERSION_MINOR**
+
+    Minor version number for OpenMC
+
+**integer(4) VERSION_RELEASE**
+
+    Release version number for OpenMC
+
+**character(19) time_stamp**
+
+    Date and time the state point was written.
+
+**character(255) path**
+
+    Absolute path to directory containing input files.
+
+**integer(8) seed**
+
+    Pseudo-random number generator seed.
+
+**integer(4) run_mode**
+
+    run mode used. The modes are described in constants.F90.
+
+**integer(8) n_particles**
+
+    Number of particles used per generation.
+
+**integer(4) current_batch**
+
+    The number of batches already simulated.
+
+if (run_mode == MODE_EIGENVALUE)
+
+    **integer(4) n_inactive**
+
+        Number of inactive batches
+
+    **integer(4) gen_per_batch**
+
+        Number of generations per batch for criticality calculations
+
+    *do i = 1, current_batch \* gen_per_batch*
+
+        **real(8) k_generation(i)**
+
+             k-effective for the i-th total generation
+
+    *do i = 1, current_batch \* gen_per_batch*
+
+        **real(8) entropy(i)**
+
+            Shannon entropy for the i-th total generation
+
+    **real(8) k_col_abs**
+
+        Sum of product of collision/absorption estimates of k-effective
+
+    **real(8) k_col_tra**
+
+        Sum of product of collision/track-length estimates of k-effective
+
+    **real(8) k_abs_tra**
+
+        Sum of product of absorption/track-length estimates of k-effective
+
+    **real(8) k_combined(2)**
+
+        Mean and standard deviation of a combined estimate of k-effective
+
+    **integer(4) cmfd_on**
+
+        Flag that cmfd is on
+
+    if (cmfd_on)
+
+        **integer(4) cmfd % indices**
+
+            Indices for cmfd mesh (i,j,k,g)
+
+        **real(8) cmfd % k_cmfd(1:current_batch)**
+
+            CMFD eigenvalues
+
+        **real(8) cmfd % src(1:G,1:I,1:J,1:K)**
+
+            CMFD fission source
+
+        **real(8) cmfd % entropy(1:current_batch)**
+
+            CMFD estimate of Shannon entropy
+
+        **real(8) cmfd % balance(1:current_batch)**
+
+            RMS of the residual neutron balance equation on CMFD mesh
+
+        **real(8) cmfd % dom(1:current_batch)**
+
+            CMFD estimate of dominance ratio
+
+        **real(8) cmfd % scr_cmp(1:current_batch)**
+
+            RMS comparison of difference between OpenMC and CMFD fission source
+
+**integer(4) n_meshes**
+
+    Number of meshes in tallies.xml file
+
+*do i = 1, n_meshes*
+
+    **integer(4) meshes(i) % id**
+
+        Unique ID of mesh.
+
+    **integer(4) meshes(i) % type**
+
+        Type of mesh.
+
+    **integer(4) meshes(i) % n_dimension**
+
+        Number of dimensions for mesh (2 or 3).
+
+    **integer(4) meshes(i) % dimension(:)**
+
+        Number of mesh cells in each dimension.
+
+    **real(8) meshes(i) % lower_left(:)**
+
+        Coordinates of lower-left corner of mesh.
+
+    **real(8) meshes(i) % upper_right(:)**
+
+        Coordinates of upper-right corner of mesh.
+
+    **real(8) meshes(i) % width(:)**
+
+        Width of each mesh cell in each dimension.
+
+**integer(4) n_tallies**
+
+*do i = 1, n_tallies*
+
+    **integer(4) tallies(i) % id**
+
+        Unique ID of tally.
+
+    **integer(4) tallies(i) % n_realizations**
+
+        Number of realizations for the i-th tally.
+
+    **integer(4) size(tallies(i) % scores, 1)**
+
+        Total number of score bins for the i-th tally
+
+    **integer(4) size(tallies(i) % scores, 2)**
+
+        Total number of filter bins for the i-th tally
+
+    **integer(4) tallies(i) % n_filters**
+
+    *do j = 1, tallies(i) % n_filters*
+
+        **integer(4) tallies(i) % filter(j) % type**
+
+            Type of tally filter.
+
+        **integer(4) tallies(i) % filter(j) % n_bins**
+
+            Number of bins for filter.
+
+        **integer(4)/real(8) tallies(i) % filter(j) % bins(:)**
+
+            Value for each filter bin of this type.
+
+    **integer(4) tallies(i) % n_nuclide_bins**
+
+        Number of nuclide bins. If none are specified, this is just one.
+
+    *do j = 1, tallies(i) % n_nuclide_bins*
+
+        **integer(4) tallies(i) % nuclide_bins(j)**
+
+            Values of specified nuclide bins
+
+    **integer(4) tallies(i) % n_score_bins**
+
+        Number of scoring bins.
+
+    *do j = 1, tallies(i) % n_score_bins*
+
+        **integer(4) tallies(i) % score_bins(j)**
+
+            Values of specified scoring bins (e.g. SCORE_FLUX).
+            
+    *do j = 1, tallies(i) % n_score_bins*
+
+        **integer(4) tallies(i) % scatt_order(j)**
+
+            Scattering Order specified scoring bins.
+            
+    **integer(4) tallies(i) % n_score_bins**
+
+        Number of scoring bins without accounting for those added by
+        the scatter-pn command.
+
+**integer(4) source_present**
+
+    Flag indicated if source bank is present in the file
+
+**integer(4) n_realizations**
+
+    Number of realizations for global tallies.
+
+**integer(4) N_GLOBAL_TALLIES**
+
+    Number of global tally scores
+
+*do i = 1, N_GLOBAL_TALLIES*
+
+    **real(8) global_tallies(i) % sum**
+
+        Accumulated sum for the i-th global tally
+
+    **real(8) global_tallies(i) % sum_sq**
+
+        Accumulated sum of squares for the i-th global tally
+
+**integer(4) tallies_on**
+
+    Flag indicated if tallies are present in the file.
+
+if (tallies_on > 0)
+
+    *do i = 1, n_tallies*
+
+        *do k = 1, size(tallies(i) % scores, 2)*
+
+            *do j = 1, size(tallies(i) % scores, 1)*
+
+                **real(8) tallies(i) % scores(j,k) % sum**
+            
+                    Accumulated sum for the j-th score and k-th filter of the
+                    i-th tally
+
+                **real(8) tallies(i) % scores(j,k) % sum_sq**
+
+                    Accumulated sum of squares for the j-th score and k-th
+                    filter of the i-th tally
+
+if (run_mode == MODE_EIGENVALUE and source_present)
+
+    *do i = 1, n_particles*
+
+        **real(8) source_bank(i) % wgt**
+
+            Weight of the i-th source particle
+
+        **real(8) source_bank(i) % xyz(1:3)**
+
+            Coordinates of the i-th source particle.
+
+        **real(8) source_bank(i) % uvw(1:3)**
+
+            Direction of the i-th source particle
+
+        **real(8) source_bank(i) % E**
+
+            Energy of the i-th source particle.
+
 -----------
 Revision 12
 -----------

docs/source/usersguide/input.rst

@@ -757,7 +757,10 @@ Each ``<cell>`` element can have the following attributes or sub-elements:
     is on the negative side of surface 3 and the positive side of surface 5, the
     bounding surfaces would be given as "-3 5".
 
-    *Default*: None
+    .. note:: The surface attribute/element can be omitted to make a cell fill
+              its entire universe.
+
+    *Default*: No surfaces
 
   :rotation:
     If the cell is filled with a universe, this element specifies the angles in
@@ -1245,7 +1248,7 @@ implemented in openMC:
 ``<plot>`` Element
 ------------------
 
-Each plot must contain a combination of the following attributes or
+Each plot is specified by a combination of the following attributes or
 sub-elements:
 
   :id:
@@ -1266,6 +1269,18 @@ sub-elements:
 
     *Default*: ``cell``
 
+  :level:
+    Universe depth to plot at (optional).  This parameter controls how many
+    universe levels deep to pull cell and material ids from when setting plot
+    colors.  If a given location does not have as many levels as specified,
+    colors will be taken from the lowest level at that location. For example, if
+    ``level`` is set to zero colors will be taken from top-level (universe zero)
+    cells only.  However, if ``level`` is set to 1 colors will be taken from
+    cells in universes that fill top-level fill-cells, and from top-level cells
+    that contain materials.
+
+    *Default*: Whatever the deepest universe is in the model
+
   :origin:
     Specifies the (x,y,z) coordinate of the center of the plot.  Should be three
     floats separated by spaces.

src/constants.F90

@@ -11,7 +11,7 @@ module constants
   integer, parameter :: VERSION_RELEASE = 1
 
   ! Revision numbers for binary files
-  integer, parameter :: REVISION_STATEPOINT       = 12
+  integer, parameter :: REVISION_STATEPOINT       = 13
   integer, parameter :: REVISION_PARTICLE_RESTART = 1
 
   ! Binary file types

src/eigenvalue.F90

@@ -23,6 +23,7 @@ module eigenvalue
                           reset_result
   use tracking,     only: transport
 
+  implicit none
   private
   public :: run_eigenvalue
 

src/endf.F90

@@ -3,6 +3,8 @@ module endf
   use constants
   use string, only: to_str
 
+  implicit none
+
 contains
 
 !===============================================================================

src/energy_grid.F90

@@ -5,6 +5,8 @@ module energy_grid
   use list_header,      only: ListReal
   use output,           only: write_message
 
+  implicit none
+
 contains
 
 !===============================================================================

src/fixed_source.F90

@@ -11,6 +11,8 @@ module fixed_source
   use tally,           only: synchronize_tallies, setup_active_usertallies
   use tracking,        only: transport
 
+  implicit none
+
 contains
 
   subroutine run_fixedsource()

src/geometry.F90

@@ -83,7 +83,8 @@ contains
       end if
     end do SURFACE_LOOP
 
-    ! If we've reached here, then the sense matched on every surface
+    ! If we've reached here, then the sense matched on every surface or there
+    ! are no surfaces.
     in_cell = .true.
 
   end function simple_cell_contains

src/initialize.F90

@@ -435,7 +435,6 @@ contains
           i = i + 1
         case ('-t', '-track', '--track')
           write_all_tracks = .true.
-          i = i + 1
         case default
           call fatal_error("Unknown command line option: " // argv(i))
         end select

src/input_xml.F90

@@ -1012,17 +1012,18 @@ contains
         call fatal_error("Cannot specify material and fill simultaneously")
       end if
 
-      ! Check to make sure that surfaces were specified
-      if (.not. check_for_node(node_cell, "surfaces")) then
-        call fatal_error("No surfaces specified for cell " &
-             &// trim(to_str(c % id)))
-      end if
-
       ! Allocate array for surfaces and copy
-      n = get_arraysize_integer(node_cell, "surfaces")
+      if (check_for_node(node_cell, "surfaces")) then
+        n = get_arraysize_integer(node_cell, "surfaces")
+      else
+        n = 0
+      end if
       c % n_surfaces = n
-      allocate(c % surfaces(n))
-      call get_node_array(node_cell, "surfaces", c % surfaces)
+
+      if (n > 0) then
+        allocate(c % surfaces(n))
+        call get_node_array(node_cell, "surfaces", c % surfaces)
+      end if
 
       ! Rotation matrix
       if (check_for_node(node_cell, "rotation")) then
@@ -2916,16 +2917,14 @@ contains
       end select
 
       ! Set output file path
-      filename = "plot"
+      filename = trim(to_str(pl % id)) // "_plot"
       if (check_for_node(node_plot, "filename")) &
         call get_node_value(node_plot, "filename", filename)
       select case (pl % type)
       case (PLOT_TYPE_SLICE)
-        pl % path_plot = trim(path_input) // trim(to_str(pl % id)) // &
-             "_" // trim(filename) // ".ppm"
+        pl % path_plot = trim(path_input) // trim(filename) // ".ppm"
       case (PLOT_TYPE_VOXEL)
-        pl % path_plot = trim(path_input) // trim(to_str(pl % id)) // &
-             "_" // trim(filename) // ".voxel"
+        pl % path_plot = trim(path_input) // trim(filename) // ".voxel"
       end select
 
       ! Copy plot pixel size
@@ -3005,6 +3004,18 @@ contains
         end if
       end if
 
+      ! Copy plot cell universe level
+      if (check_for_node(node_plot, "level")) then
+        call get_node_value(node_plot, "level", pl % level)
+        
+        if (pl % level < 0) then
+          call fatal_error("Bad universe level in plot " &
+               &// trim(to_str(pl % id)))
+        end if
+      else
+        pl % level = PLOT_LEVEL_LOWEST
+      end if
+
       ! Copy plot color type and initialize all colors randomly
       temp_str = "cell"
       if (check_for_node(node_plot, "color")) &

src/output.F90

@@ -1463,6 +1463,12 @@ contains
       ! Plot id
       write(ou,100) "Plot ID:", trim(to_str(pl % id))
 
+      ! Plot filename
+      write(ou,100) "Plot file:", trim(pl % path_plot)
+
+      ! Plot level
+      write(ou,100) "Universe depth:", trim(to_str(pl % level))
+
       ! Plot type
       if (pl % type == PLOT_TYPE_SLICE) then
         write(ou,100) "Plot Type:", "Slice"

src/plot.F90

@@ -7,7 +7,7 @@ module plot
   use global
   use mesh,            only: get_mesh_indices
   use output,          only: write_message
-  use particle_header, only: deallocate_coord, Particle
+  use particle_header, only: deallocate_coord, Particle, LocalCoord
   use plot_header
   use ppmlib,          only: Image, init_image, allocate_image, &
                              deallocate_image, set_pixel
@@ -32,7 +32,7 @@ contains
 
       ! Display output message
       call write_message("Processing plot " // trim(to_str(pl % id)) &
-           &// "...", 5)
+           &// ": " // trim(pl % path_plot) // " ...", 5)
 
       if (pl % type == PLOT_TYPE_SLICE) then
         ! create 2d image
@@ -58,7 +58,9 @@ contains
     integer, intent(out)                  :: id
     
     logical :: found_cell
-    type(Cell), pointer :: c => null()
+    integer :: level
+    type(Cell),       pointer :: c => null()
+    type(LocalCoord), pointer :: coord => null()
     
     call deallocate_coord(p % coord0 % next)
     p % coord => p % coord0
@@ -66,6 +68,16 @@ contains
     call find_cell(p, found_cell)
     if (check_overlaps) call check_cell_overlap(p)
 
+    ! Loop through universes and stop on any specified level
+    level = 0
+    coord => p % coord0
+    do
+      if (level == pl % level) exit
+      if (.not. associated(coord % next)) exit
+      coord => coord % next
+      level = level + 1
+    end do
+    
     if (.not. found_cell) then
       ! If no cell, revert to default color
       rgb = pl % not_found % rgb
@@ -73,19 +85,23 @@ contains
     else
       if (pl % color_by == PLOT_COLOR_MATS) then
         ! Assign color based on material
-        c => cells(p % coord % cell)
+        c => cells(coord % cell)
         if (c % material == MATERIAL_VOID) then
           ! By default, color void cells white
           rgb = 255
           id = -1
+        else if (c % type == CELL_FILL) then
+          ! If we stopped on a middle universe level, treat as if not found
+          rgb = pl % not_found % rgb
+          id = -1
         else
           rgb = pl % colors(c % material) % rgb
           id = materials(c % material) % id
         end if
       else if (pl % color_by == PLOT_COLOR_CELLS) then
         ! Assign color based on cell
-        rgb = pl % colors(p % coord % cell) % rgb
-        id = cells(p % coord % cell) % id
+        rgb = pl % colors(coord % cell) % rgb
+        id = cells(coord % cell) % id
       else
         rgb = 0
         id = -1

src/plot_header.F90

@@ -27,6 +27,7 @@ module plot_header
     integer :: basis                 ! direction of plot slice 
     integer :: pixels(3)             ! pixel width/height of plot slice
     integer :: meshlines_width       ! pixel width of meshlines
+    integer :: level                 ! universe depth to plot the cells of
     type(StructuredMesh), pointer :: meshlines_mesh => null() ! mesh to plot
     type(ObjectColor) :: meshlines_color ! Color for meshlines
     type(ObjectColor) :: not_found   ! color for positions where no cell found
@@ -36,6 +37,9 @@ module plot_header
   ! Plot type
   integer, parameter :: PLOT_TYPE_SLICE = 1
   integer, parameter :: PLOT_TYPE_VOXEL = 2
+  
+  ! Plot level
+  integer, parameter :: PLOT_LEVEL_LOWEST = -1
 
   ! Plot basis plane
   integer, parameter :: PLOT_BASIS_XY = 1

src/relaxng/geometry.rnc

@@ -7,7 +7,7 @@ element geometry {
       (element material { ( xsd:int | "void" ) } | 
         attribute material { ( xsd:int | "void" ) })
     ) &
-    (element surfaces { list { xsd:int+ } } | attribute surfaces { list { xsd:int+ } }) &
+    (element surfaces { list { xsd:int* } } | attribute surfaces { list { xsd:int* } })? &
     (element rotation { list { xsd:double+ } } | attribute rotation { list { xsd:double+ } })? &
     (element translation { list { xsd:double+ } } | attribute translation { list { xsd:double+ } })?
   }*

src/relaxng/plots.rnc

@@ -7,6 +7,7 @@ element plots {
       attribute type { "slice" | "voxel" })? &
     (element color { ( "cell" | "mat" | "material" ) } |
       attribute color { ( "cell" | "mat" | "material" ) })? &
+    (element level { xsd:int } | attribute level { xsd:int })? &
     (element origin { list { xsd:double+ } } |
       attribute origin { list { xsd:double+ } })? &
     (element width { list { xsd:double+ } } |

src/search.F90

@@ -3,6 +3,8 @@ module search
   use constants
   use error,     only: fatal_error
 
+  implicit none
+
   integer, parameter :: MAX_ITERATION = 64
 
   interface binary_search

src/state_point.F90

@@ -83,7 +83,6 @@ contains
       ! Write run information
       call sp % write_data(run_mode, "run_mode")
       call sp % write_data(n_particles, "n_particles")
-      call sp % write_data(n_batches, "n_batches")
 
       ! Write out current batch number
       call sp % write_data(current_batch, "current_batch")
@@ -306,7 +305,7 @@ contains
         ! Set filename
         filename = trim(path_output) // 'source.' // &
             & zero_padded(current_batch, count_digits(n_batches))
- 
+
 #ifdef HDF5
         filename = trim(filename) // '.h5'
 #else
@@ -568,14 +567,15 @@ contains
     ! Read and overwrite run information except number of batches
     call sp % read_data(run_mode, "run_mode")
     call sp % read_data(n_particles, "n_particles")
-    call sp % read_data(int_array(1), "n_batches")
-
-    ! Take maximum of statepoint n_batches and input n_batches
-    n_batches = max(n_batches, int_array(1))
 
     ! Read batch number to restart at
     call sp % read_data(restart_batch, "current_batch")
 
+    if (restart_batch > n_batches) then
+      call fatal_error("The number batches specified in settings.xml is fewer &
+           & than the number of batches in the given statepoint file.")
+    end if
+
     ! Read information specific to eigenvalue run
     if (run_mode == MODE_EIGENVALUE) then
       call sp % read_data(int_array(1), "n_inactive")

src/tracking.F90

@@ -16,6 +16,8 @@ module tracking
   use track_output,    only: initialize_particle_track, write_particle_track, &
                              finalize_particle_track
 
+  implicit none
+
 contains
 
 !===============================================================================

src/utils/statepoint.py

@@ -6,7 +6,7 @@ from collections import OrderedDict
 import numpy as np
 import scipy.stats
 
-REVISION_STATEPOINT = 12
+REVISION_STATEPOINT = 13
 
 filter_types = {1: 'universe', 2: 'material', 3: 'cell', 4: 'cellborn',
                 5: 'surface', 6: 'mesh', 7: 'energyin', 8: 'energyout'}
@@ -182,7 +182,6 @@ class StatePoint(object):
         # Read run information
         self.run_mode = self._get_int(path='run_mode')[0]
         self.n_particles = self._get_long(path='n_particles')[0]
-        self.n_batches = self._get_int(path='n_batches')[0]
 
         # Read current batch
         self.current_batch = self._get_int(path='current_batch')[0]

tests/test_infinite_cell/geometry.xml

@@ -0,0 +1,17 @@
+<?xml version="1.0"?>
+<geometry>
+  <cell id="11" universe="11" material="1"/>
+  <cell id="12" universe="12" material="2" surfaces=""/>
+
+  <lattice id="21" type="rect" dimension="2 2" lower_left="-2.0 -2.0"
+    width="2.0 2.0" outside="2">
+    <universes>
+      11 12
+      12 11
+    </universes>
+  </lattice>
+
+
+  <surface id="101" type="z-cylinder" coeffs="0.0 0.0 5.0" boundary="vacuum"/>
+  <cell id="101" universe="0" fill="21" surfaces="-101"/>
+</geometry>

tests/test_infinite_cell/materials.xml

@@ -0,0 +1,14 @@
+<?xml version="1.0"?>
+<materials>
+
+  <material id="1">
+    <density value="4.5" units="g/cc" />
+    <nuclide name="U-235" xs="71c" ao="1.0" />
+  </material>
+
+  <material id="2">
+    <density value="4.5" units="g/cc" />
+    <nuclide name="U-238" xs="71c" ao="1.0" />
+  </material>
+
+</materials>

tests/test_infinite_cell/results.py

@@ -0,0 +1,25 @@
+#!/usr/bin/env python
+
+import sys
+
+# import statepoint
+sys.path.insert(0, '../../src/utils')
+import statepoint
+
+# read in statepoint file
+if len(sys.argv) > 1:
+    sp = statepoint.StatePoint(sys.argv[1])
+else:
+    sp = statepoint.StatePoint('statepoint.10.binary')
+sp.read_results()
+
+# set up output string
+outstr = ''
+ 
+# write out k-combined
+outstr += 'k-combined:\n'
+outstr += "{0:12.6E} {1:12.6E}\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)

tests/test_infinite_cell/results_true.dat

@@ -0,0 +1,2 @@
+k-combined:
+9.998895E-02 2.846817E-04

tests/test_infinite_cell/settings.xml

@@ -0,0 +1,16 @@
+<?xml version="1.0"?>
+<settings>
+
+  <eigenvalue>
+    <batches>10</batches>
+    <inactive>5</inactive>
+    <particles>1000</particles>
+  </eigenvalue>
+
+  <source>
+    <space type="box">
+      <parameters>-4 -4 -4  4  4  4</parameters>
+    </space>
+  </source>
+
+</settings>

tests/test_infinite_cell/test_infinite_cell.py

@@ -0,0 +1,59 @@
+#!/usr/bin/env python
+
+import os
+from subprocess import Popen, STDOUT, PIPE, call
+import filecmp
+import glob
+from optparse import OptionParser
+
+parser = OptionParser()
+parser.add_option('--mpi_exec', dest='mpi_exec', default='')
+parser.add_option('--mpi_np', dest='mpi_np', default='3')
+parser.add_option('--exe', dest='exe')
+(opts, args) = parser.parse_args()
+cwd = os.getcwd()
+
+def test_run():
+    if opts.mpi_exec != '':
+        proc = Popen([opts.mpi_exec, '-np', opts.mpi_np, opts.exe, cwd],
+               stderr=STDOUT, stdout=PIPE)
+    else:
+        proc = Popen([opts.exe, cwd], stderr=STDOUT, stdout=PIPE)
+    print(proc.communicate()[0])
+    returncode = proc.returncode
+    assert returncode == 0, 'OpenMC did not exit successfully.'
+
+def test_created_statepoint():
+    statepoint = glob.glob(os.path.join(cwd, 'statepoint.10.*'))
+    assert len(statepoint) == 1, 'Either multiple or no statepoint files exist.'
+    assert statepoint[0].endswith('binary') or statepoint[0].endswith('h5'),\
+        'Statepoint file is not a binary or hdf5 file.'
+
+def test_results():
+    statepoint = glob.glob(os.path.join(cwd, 'statepoint.10.*'))
+    call(['python', 'results.py', statepoint[0]])
+    compare = filecmp.cmp('results_test.dat', 'results_true.dat')
+    if not compare:
+      os.rename('results_test.dat', 'results_error.dat')
+    assert compare, 'Results do not agree.'
+
+def teardown():
+    output = glob.glob(os.path.join(cwd, 'statepoint.10.*'))
+    output.append(os.path.join(cwd, 'results_test.dat'))
+    for f in output:
+        if os.path.exists(f):
+            os.remove(f)
+
+if __name__ == '__main__':
+
+    # test for openmc executable
+    if opts.exe is None:
+        raise Exception('Must specify OpenMC executable from command line with --exe.')
+
+    # run tests
+    try:
+        test_run()
+        test_created_statepoint()
+        test_results()
+    finally:
+        teardown()