diff --git a/ChangeLog b/ChangeLog
index fde6957f29..80e5465e1a 100644
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,27 @@
+2011-02-03 Paul Romano <romano7@mit.edu>
+
+	* ace.f90: Added initial routines for reading continuous-energy
+	ACE cross section tables. No secondary energy/angle distributions
+	or nu data is read.
+	* cross_section.f90: Added read_xsdata routine which reads
+	information from an xsdata file from SERPENT.
+	* data_structures.f90: Moved derived types to the types module and
+	change dict_create so that one doesn't need to add a (key,value)
+	pair on initialization.
+	* energy_grid.f90: First routines for energy grid unionization
+	based on insertion sorts.
+	* fileio.f90: Added routines for reading integer and real data,
+	particularly from cross-section libraries. Made an interface so
+	that a single routine can be called regardless of type. Ability to
+	read material cards and added subroutine for normalization of
+	atom/weight percents.
+	* Makefile: Added 'unittest' target
+	* output.f90: Fixed bug in get_today subroutine whereby the day
+	was displayed wrong.
+	* search.f90: Added binary search which will eventually be used
+	for energy grid searches
+	* unittest.f90: Added this program for unit testing.
+
 2011-01-26 Paul Romano <romano7@mit.edu>
 
 	* physics.f90: Added collision routine, still segfaults though.
diff --git a/src/Makefile b/src/Makefile
index b0c122fd82..bb1eb4218f 100644
--- a/src/Makefile
+++ b/src/Makefile
@@ -1,18 +1,24 @@
 program = openmc
 
-src = 	data_structures.f90 \
+main = main.f90
+test = unittest.f90
+modules = 	ace.f90 \
+		cross_section.f90 \
+		data_structures.f90 \
+		energy_grid.f90 \
 		fileio.f90 \
 		geometry.f90 \
 		global.f90 \
-		main.f90 \
 		mcnp_random.f90 \
 		output.f90 \
 		physics.f90 \
+		search.f90 \
 		source.f90 \
 		string.f90 \
-		types.f90
+		types.f90 \
 
-objects = $(src:.f90=.o)
+main_objects = $(modules:.f90=.o) $(main:.f90=.o)
+test_objects = $(modules:.f90=.o) $(test:.f90=.o)
 
 #--------------------------------------------------------------------
 # Compiler Options
@@ -22,12 +28,14 @@ F90FLAGS = -g
 #--------------------------------------------------------------------
 # Targets
 all:	$(program)
-$(program): $(objects)
-	$(F90) $(F90FLAGS) $(objects) -o $@
+$(program): $(main_objects)
+	$(F90) $(F90FLAGS) $(main_objects) -o $@
 clean:
 	@rm -f *.o *.mod $(program)
 neat:
 	@rm -f *.o *.mod
+unittest: $(test_objects)
+	$(F90) $(F90FLAGS) $(test_objects) -o $@
 
 #--------------------------------------------------------------------
 # Rules & misc
@@ -39,14 +47,20 @@ neat:
 
 #--------------------------------------------------------------------
 # Dependencies
+ace.o: global.o output.o string.o fileio.o string.o
+cross_section.o: global.o string.o data_structures.o output.o
 data_structures.o: global.o
-fileio.o: types.o global.o string.o output.o
+energy_grid.o: global.o
+fileio.o: types.o global.o string.o output.o data_structures.o
 geometry.o: types.o global.o output.o string.o
 global.o: types.o
 main.o: global.o fileio.o output.o geometry.o mcnp_random.o \
-		source.o physics.o
+		source.o physics.o cross_section.o data_structures.o \
+		ace.o energy_grid.o
 output.o: global.o
 physics.o: types.o global.o mcnp_random.o geometry.o output.o
+search.o:
 source.o: global.o mcnp_random.o
 string.o: global.o output.o
 types.o:
+unittest.o: global.o energy_grid.o
diff --git a/src/ace.f90 b/src/ace.f90
new file mode 100644
index 0000000000..0227ce2216
--- /dev/null
+++ b/src/ace.f90
@@ -0,0 +1,328 @@
+module ace
+
+  use global
+  use output, only: error, message
+  use string, only: lower_case
+  use fileio, only: read_line, read_data, skip_lines
+  use string, only: split_string, str_to_real
+  use data_structures, only: dict_create, dict_add_key, dict_has_key, &
+       &                     dict_get_key
+
+  integer :: NXS(16)
+  integer :: JXS(32)
+  real(8), allocatable :: XSS(:)
+  integer :: XSS_index
+
+  private :: NXS
+  private :: JXS
+  private :: XSS
+
+contains
+
+!=====================================================================
+! READ_XS reads all the cross sections for the problem and stores them
+! in xs_continuous and xs_thermal arrays
+!=====================================================================
+
+  subroutine read_xs()
+
+    type(Material),      pointer :: mat => null()
+    type(xsData),        pointer :: iso => null()
+    type(AceContinuous), pointer :: ace_cont => null()
+    type(AceThermal),    pointer :: ace_thermal => null()
+    integer :: i, j
+    integer :: index
+    character(10) :: table
+    character(250) :: msg
+    integer :: n
+    integer :: n_continuous
+    integer :: n_thermal
+    integer :: index_continuous
+    integer :: index_thermal
+
+    ! determine how many continuous-energy tables and how many S(a,b)
+    ! thermal scattering tables there are
+    n_continuous = 0
+    n_thermal = 0
+    do i = 1, n_materials
+       mat => materials(i)
+       do j = 1, mat%n_isotopes
+          index = mat%isotopes(j)
+          table = xsdatas(index)%id
+          n = len_trim(table)
+          call lower_case(table)
+          select case (table(n:n))
+          case ('c')
+             n_continuous = n_continuous + 1
+          case ('t')
+             n_thermal = n_thermal + 1
+          case default
+             msg = "Unknown cross section table type: " // table
+             call error(msg)
+          end select
+       end do
+    end do
+
+    ! allocate arrays for ACE table storage
+    allocate(xs_continuous(n_continuous))
+    allocate(xs_thermal(n_thermal))
+
+    ! loop over all materials
+    call dict_create(ace_dict)
+    index_continuous = 0
+    index_thermal = 0
+    do i = 1, n_materials
+       mat => materials(i)
+       do j = 1, mat%n_isotopes
+          index = mat%isotopes(j)
+          table = xsdatas(index)%id
+          n = len_trim(table)
+          call lower_case(table)
+          select case (table(n:n))
+          case ('c')
+             if (.not. dict_has_key(ace_dict, table)) then
+                index_continuous = index_continuous + 1
+                call dict_add_key(ace_dict, table, index_continuous)
+                call read_ACE_continuous(index_continuous, index)
+                mat%table(j) = index_continuous
+             end if
+          case ('t')
+             index_thermal = index_thermal + 1
+             n_thermal = n_thermal + 1
+          case default
+             msg = "Unknown cross section table type: " // table
+             call error(msg)
+          end select
+       end do
+    end do
+       
+  end subroutine read_xs
+
+!=====================================================================
+! READ_ACE_CONTINUOUS reads in a single ACE continuous-energy cross
+! section table
+!=====================================================================
+
+  subroutine read_ACE_continuous(index_table, index)
+
+    integer, intent(in) :: index_table
+    integer, intent(in) :: index
+
+    type(AceContinuous), pointer :: table => null()
+    integer :: i
+    integer :: in = 7
+    integer :: ioError
+    integer :: words_per_line
+    integer :: lines
+    integer :: n
+    logical :: file_exists
+    logical :: found_xs
+    character(7) :: readable
+    character(250) :: msg, line
+    character(32) :: words(max_words)
+    character(100) :: filename
+    character(10) :: tablename
+    real(8) :: kT
+
+    filename = xsdatas(index)%path
+    tablename = xsdatas(index)%id
+
+    table => xs_continuous(index_table)
+
+    ! Check if input file exists and is readable
+    inquire( FILE=filename, EXIST=file_exists, READ=readable )
+    if ( .not. file_exists ) then
+       msg = "ACE library '" // trim(filename) // "' does not exist!"
+       call error( msg )
+    elseif ( readable(1:3) /= 'YES' ) then
+       msg = "ACE library '" // trim(filename) // "' is not readable! &
+            &Change file permissions with chmod command."
+       call error( msg )
+    end if
+
+    ! display message
+    msg = "Loading ACE cross section table: " // tablename
+    call message(msg, 6)
+
+    ! open file
+    open(file=filename, unit=in, status='old', action='read')
+
+    found_xs = .false.
+    do while (.not. found_xs)
+       call read_line(in, line, ioError)
+       if (ioError < 0) then
+          msg = "Could not find ACE table " // tablename // "."
+          call error(msg)
+       end if
+       call split_string(line, words, n)
+       if (trim(words(1)) == trim(tablename)) then
+          found_xs = .true.
+          table%name = words(1)
+          table%awr = str_to_real(words(2))
+          kT = str_to_real(words(3))
+          table%temp = kT / K_BOLTZMANN
+       end if
+       
+       ! Skip 5 lines
+       call skip_lines(in, 5, ioError)
+
+       ! Read NXS data
+       lines = 2
+       words_per_line = 8
+       call read_data(in, NXS, 16, lines, words_per_line)
+
+       ! Read JXS data
+       lines = 4
+       call read_data(in, JXS, 32, lines, words_per_line)
+
+       ! Calculate how many data points and lines in the XSS array
+       n = NXS(1)
+       lines = (n + 3)/4
+
+       if (found_xs) then
+          ! allocate storage for XSS array
+          allocate(XSS(n))
+          
+          ! Read XSS
+          words_per_line = 4
+          call read_data(in, XSS, n, lines, words_per_line)
+       else
+          call skip_lines(in, lines, ioError)
+       end if
+
+    end do
+
+    call parse_ESZ(table)
+    call parse_MTR(table)
+
+    ! Free memory from XSS array
+    if(allocated(XSS)) deallocate(XSS)
+    if(associated(table)) nullify(table)
+
+    close(unit=in)
+
+  end subroutine read_ACE_continuous
+
+!=====================================================================
+! PARSE_ESZ - reads through the ESZ block. This block contains the
+! energy grid, total xs, absorption xs, elastic scattering xs, and
+! heating numbers.
+!=====================================================================
+
+  subroutine parse_ESZ(table)
+
+    type(AceContinuous), pointer, intent(inout) :: table
+
+    integer :: NE
+
+    ! determine number of energy points
+    NE = NXS(3)
+
+    ! allocate storage for arrays
+    allocate(table%energy(NE))
+    allocate(table%sigma_t(NE))
+    allocate(table%sigma_a(NE))
+    allocate(table%sigma_el(NE))
+    allocate(table%heating(NE))
+
+    ! read data from XSS
+    XSS_index = 1
+    table%energy = get_real(NE)
+    table%sigma_t = get_real(NE)
+    table%sigma_a = get_real(NE)
+    table%sigma_el = get_real(NE)
+    table%heating = get_real(NE)
+    
+  end subroutine parse_ESZ
+
+!=====================================================================
+! PARSE_MTR - Get the list of reaction MTs for this cross-section
+! table. The MT values are somewhat arbitrary. Also read in Q-values,
+! neutron multiplicities, and cross-sections.
+! =====================================================================
+
+  subroutine parse_MTR(table)
+
+    type(AceContinuous), pointer :: table
+
+    type(AceReaction), pointer :: rxn
+    integer :: LMT   ! index of MT list in XSS
+    integer :: NMT   ! Number of reactions
+    integer :: JXS4  ! index of Q values in XSS
+    integer :: JXS5  ! index of neutron multiplicities in XSS
+    integer :: JXS7  ! index of reactions cross-sections in XSS
+    integer :: LXS   ! 
+    integer :: LOCA  ! 
+    integer :: NE    ! number of energies for reaction
+    
+    LMT = JXS(3)
+    JXS4 = JXS(4)
+    JXS5 = JXS(5)
+    LXS  = JXS(6)
+    JXS7 = JXS(7)
+    NMT = NXS(4)
+
+    ! allocate array of reactions
+    allocate(table%reactions(NMT))
+
+    do i = 1, NMT
+       rxn => table%reactions(i)
+
+       ! read MT number, Q-value, and neutrons produced
+       rxn%MT      = XSS(LMT+i-1)
+       rxn%Q_value = XSS(JXS4+i-1)
+       rxn%TY      = XSS(JXS5+i-1)
+
+       ! read cross section values
+       LOCA = XSS(LXS+i-1)
+       rxn%energy_index = XSS(JXS7 + LOCA - 1)
+       NE = XSS(JXS7 + LOCA)
+       allocate(rxn%sigma(NE))
+       XSS_index = JXS7 + LOCA + 1
+       rxn%sigma = get_real(NE)
+    end do
+
+  end subroutine parse_MTR
+
+!=====================================================================
+! GET_INT returns an array of integers read from the current position
+! in the XSS array
+!=====================================================================
+
+    function get_int(n_values) result(array)
+
+      integer, intent(in) :: n_values
+      integer :: array(n_values)
+
+      array = int(XSS(XSS_index:XSS_index + n_values - 1))
+      XSS_index = XSS_index + n_values
+      
+    end function get_int
+
+!=====================================================================
+! GET_REAL returns an array of real(8)s read from the current position
+! in the XSS array
+!=====================================================================
+
+    function get_real(n_values) result(array)
+
+      integer, intent(in) :: n_values
+      real(8) :: array(n_values)
+
+      array = XSS(XSS_index:XSS_index + n_values - 1)
+      XSS_index = XSS_index + n_values
+
+    end function get_real
+
+!=====================================================================
+! READ_ACE_THERMAL reads in a single ACE S(a,b) thermal scattering
+! cross section table
+!=====================================================================
+
+  subroutine read_ACE_thermal( filename )
+
+    character(*), intent(in) :: filename
+
+  end subroutine read_ACE_thermal
+
+end module ace
diff --git a/src/cross_section.f90 b/src/cross_section.f90
new file mode 100644
index 0000000000..79d629d938
--- /dev/null
+++ b/src/cross_section.f90
@@ -0,0 +1,113 @@
+module cross_section
+
+  use string, only: split_string, str_to_real
+  use global, only: str_to_int, max_words, xsdata_dict, xsdatas
+  use data_structures, only: Dictionary, dict_create, dict_add_key, dict_get_key
+  use output, only: error
+  use types, only: xsData
+
+  implicit none
+
+contains
+
+!=====================================================================
+! READ_XSDATA reads the data in a SERPENT xsdata file and builds a
+! dictionary to find cross-section information later on.
+!=====================================================================
+
+  subroutine read_xsdata( path )
+
+    character(*), intent(in) :: path
+
+    type(xsData), pointer :: iso => null()
+    character(250) :: line, msg
+    character(100) :: words(max_words)
+    character(100) :: filename
+    integer :: i, n
+    integer :: in = 7
+    logical :: file_exists
+    character(7) :: readable
+    integer :: count
+    integer :: index
+    integer :: ioError
+
+    ! Construct filename
+    n = len(path)
+    if (path(n:n) == '/') then
+       filename = trim(path) // 'xsdata'
+    else
+       filename = trim(path) // '/xsdata'
+    end if
+
+    ! Check if xsdata exists and is readable
+    inquire( FILE=filename, EXIST=file_exists, READ=readable )
+    if ( .not. file_exists ) then
+       msg = "Cross section summary '" // trim(filename) // "' does not exist!"
+       call error( msg )
+    elseif ( readable(1:3) /= 'YES' ) then
+       msg = "Cross section summary '" // trim(filename) // "' is not readable!" &
+            & // "Change file permissions with chmod command."
+       call error( msg )
+    end if
+
+    ! open file
+    open(file=filename, unit=in, status='old', action='read')
+
+    ! determine how many lines
+    count = 0
+    do
+       read(unit=in, fmt='(A250)', iostat=ioError) line
+       if ( ioError < 0 ) then
+          ! reached end of file
+          exit
+       elseif ( ioError > 0 ) then
+          msg = "Unknown error while reading file: " // filename
+          close(unit=in)
+          call error( msg )
+       end if
+       count = count + 1
+    end do
+    allocate(xsdatas(count))
+
+    ! create dictionary
+    call dict_create(xsdata_dict)
+
+    ! read actual lines
+    index = 0
+    rewind(in)
+    do
+       read(unit=in, fmt='(A250)', iostat=ioError) line
+       if ( ioError < 0 ) exit
+       index = index + 1
+       call split_string(line, words, n)
+       if ( n == 0 ) cycle ! skip blank line
+
+       ! Check to make sure there are enough arguments
+       if ( n < 9 ) then
+          msg = "Not enough arguments on xsdata line: " // line
+          close(unit=in)
+          call error( msg )
+       end if
+
+       iso => xsdatas(index)
+
+       ! store data
+       iso%alias = words(1)
+       iso%id = words(2)
+       iso%type = str_to_int(words(3))
+       iso%zaid = str_to_int(words(4))
+       iso%isomeric = str_to_int(words(5))
+       iso%awr = str_to_real(words(6))
+       iso%temp = str_to_real(words(7))
+       iso%binary = str_to_int(words(8))
+       iso%path = words(9)
+
+       ! create dictionary entry
+       call dict_add_key(xsdata_dict, iso%alias, index)
+    end do
+
+    close(unit=in)
+
+  end subroutine read_xsdata
+
+end module cross_section
diff --git a/src/data_structures.f90 b/src/data_structures.f90
index faca8b751d..1882f570b0 100644
--- a/src/data_structures.f90
+++ b/src/data_structures.f90
@@ -14,34 +14,10 @@ module data_structures
 ! to ListData, dict_create, dict_add_key, and dict_get_key).
 !=====================================================================
 
-  use global, only: DICT_KEY_LENGTH
+  use types, only: ListData, LinkedList, HashList, Dictionary
 
   implicit none
 
-  ! Data stored in a linked list. In this case, we store the
-  ! (key,value) pair for a dictionary. Note that we need to store the
-  ! key in addition to the value for collision resolution.
-  type ListData
-     character(len=DICT_KEY_LENGTH) :: key
-     integer                        :: value
-  end type ListData
-
-  ! A simple linked list
-  type LinkedList
-     type(LinkedList), pointer :: next
-     type(ListData)            :: data
-  end type LinkedList
-
-  type HashList
-     type(LinkedList), pointer :: list
-  end type HashList
-
-  ! A dictionary of (key,value) pairs
-  type Dictionary
-     private
-     type(HashList), pointer, dimension(:) :: table
-  end type Dictionary
-
   ! Hide objects that do not need to be publicly accessible
   private :: ListData
   private :: HashList
@@ -281,15 +257,15 @@ contains
 !     destroy up an old list.
 !=====================================================================
 
-  subroutine dict_create( dict, key, value )
+  subroutine dict_create( dict ) !, key, value )
 
     type(Dictionary), pointer   :: dict
-    character(len=*), intent(in) :: key
-    integer,  intent(in) :: value
+!!$    character(len=*), intent(in) :: key
+!!$    integer,  intent(in) :: value
 
-    type(ListData)              :: data
+!!$    type(ListData)              :: data
     integer                      :: i
-    integer                      :: hash
+!!$    integer                      :: hash
 
     allocate( dict )
     allocate( dict%table(hash_size) )
@@ -298,11 +274,11 @@ contains
        dict%table(i)%list => null()
     enddo
 
-    data%key   = key
-    data%value = value
-
-    hash = dict_hashkey( trim(key ) )
-    call list_create( dict%table(hash)%list, data )
+!!$    data%key   = key
+!!$    data%value = value
+!!$
+!!$    hash = dict_hashkey( trim(key ) )
+!!$    call list_create( dict%table(hash)%list, data )
 
   end subroutine dict_create
 
diff --git a/src/energy_grid.f90 b/src/energy_grid.f90
new file mode 100644
index 0000000000..6357e4bdf8
--- /dev/null
+++ b/src/energy_grid.f90
@@ -0,0 +1,167 @@
+module energy_grid
+
+  use global
+  use output, only: message
+
+contains
+
+!=====================================================================
+! UNIONIZED_GRID creates a single unionized energy grid combined from
+! each nuclide of each material. Right now, the grid for each nuclide
+! is added into a linked list one at a time with an effective
+! insertion sort. Could be done with a hash for all energy points and
+! then a quicksort at the end (what hash function to use?)
+!=====================================================================
+
+  subroutine unionized_grid()
+
+    type(LinkedListGrid), pointer :: list => null()
+    type(LinkedListGrid), pointer :: current => null()
+
+    type(Material),       pointer :: mat => null()
+    type(AceContinuous),  pointer :: table => null()
+    type(AceReaction),    pointer :: rxn => null()
+    integer :: i, j
+    integer :: n
+    character(100) :: msg
+
+    msg = "Creating unionized energy grid"
+    call message(msg, 5)
+
+    ! loop over all materials
+    do i = 1, n_materials
+       mat => materials(i)
+       
+       ! loop over all isotopes
+       do j = 1, mat%n_isotopes
+          table => xs_continuous(mat%table(j))
+
+          ! loop over energy points
+          n = size(table%energy)
+          call add_grid_points(list, table%energy, n)
+       end do
+    end do
+
+  end subroutine unionized_grid
+
+!=====================================================================
+! ADD_GRID_POINTS adds energy points from the 'energy' array into a
+! linked list of points already stored from previous arrays. 
+!=====================================================================
+
+  subroutine add_grid_points(list, energy, n_energy)
+
+    type(LinkedListGrid), pointer :: list
+    real(8), intent(in) :: energy(n_energy)
+    integer, intent(in) :: n_energy
+
+    type(LinkedListGrid), pointer :: current => null()
+    type(LinkedListGrid), pointer :: previous => null()
+    type(LinkedListGrid), pointer :: head => null()
+    type(LinkedListGrid), pointer :: tmp => null()
+    integer :: index
+    real(8) :: E
+
+    index = 1
+
+    ! if the original list is empty, we need to allocate the first
+    ! element and store first energy point
+    if (grid_count(list) == 0) then
+       allocate(list)
+       current => list
+       do index = 1, n_energy
+          current%energy = energy(index)
+          if (index == n_energy) then
+             current%next => null()
+             return
+          end if
+          allocate(current%next)
+          current => current%next
+       end do
+    end if
+
+    current => list
+    head => list
+
+    E = energy(index)
+    do while (index <= n_energy)
+
+       ! If we've reached the end of the grid energy list, add the
+       ! remaining energy points to the end
+       if (.not. associated(current)) then
+          ! finish remaining energies
+          do while (index <= n_energy)
+             allocate(previous%next)
+             current => previous%next
+             current%energy = energy(index)
+             previous => current
+             index = index + 1
+          end do
+          current%next => null()
+          exit
+       end if
+       
+       if (E < current%energy) then
+          ! create new element and insert it in energy grid list
+          allocate(tmp)
+          tmp%energy = E
+          tmp%next => current
+          if (associated(previous)) then
+             previous%next => tmp
+             previous => tmp
+          else
+             previous => tmp
+             head => previous
+          end if
+          nullify(tmp)
+
+          ! advance index
+          index = index + 1
+          E = energy(index)
+
+       elseif (E == current%energy) then
+          ! found the exact same energy, no need to store duplicates
+          ! so just skip and move to next index
+          index = index + 1
+          E = energy(index)
+       else
+          previous => current
+          current => current%next
+       end if
+       
+    end do
+
+    ! It's possible that an element was inserted at the front of the
+    ! list, so we need to move the list pointer back to the start of
+    ! the list
+    list => head
+
+  end subroutine add_grid_points
+
+!=====================================================================
+! GRID_COUNT gives the number of energy points in list. This should
+! eventually be replaced with a generic list function.
+!=====================================================================
+
+  function grid_count(list) result(count)
+
+    type(LinkedListGrid), pointer :: list
+    integer :: count
+
+    type(LinkedListGrid), pointer :: current
+
+    ! determine size of list
+    if (associated(list)) then
+       count = 1
+       current => list
+       do while (associated(current%next))
+          current => current%next
+          count = count + 1
+       enddo
+    else
+       count = 0
+    end if
+
+  end function grid_count
+
+end module energy_grid
diff --git a/src/fileio.f90 b/src/fileio.f90
index d343c8989f..8ae87126b2 100644
--- a/src/fileio.f90
+++ b/src/fileio.f90
@@ -4,9 +4,21 @@ module fileio
   use types,  only: Cell, Surface, ExtSource
   use string, only: split_string_wl, lower_case, str_to_real, split_string
   use output, only: message, warning, error
+  use data_structures, only: Dictionary, dict_create, dict_add_key, dict_get_key
 
   implicit none
 
+!=====================================================================
+! READ_DATA interface allows data to be read with one function
+! regardless of whether it is integer or real data. E.g. NXS and JXS
+! can be read with the integer version and XSS can be read with the
+! real version
+!=====================================================================
+
+  interface read_data
+     module procedure read_data_int, read_data_real
+  end interface
+
 contains
 
 !=====================================================================
@@ -22,21 +34,21 @@ contains
 
     argc = COMMAND_ARGUMENT_COUNT()
     if ( argc > 0 ) then
-       call GET_COMMAND_ARGUMENT(1,inputfile)
+       call GET_COMMAND_ARGUMENT(1, path_input)
     else
        msg = "No input file specified!"
-       call error( msg )
+       call error(msg)
     end if
 
     ! Check if input file exists and is readable
-    inquire( FILE=inputfile, EXIST=file_exists, READ=readable )
-    if ( .not. file_exists ) then
-       msg = "Input file '" // trim(inputfile) // "' does not exist!"
-       call error( msg )
-    elseif ( readable(1:3) /= 'YES' ) then
-       msg = "Input file '" // trim(inputfile) // "' is not readable! &
+    inquire(FILE=path_input, EXIST=file_exists, READ=readable)
+    if (.not. file_exists) then
+       msg = "Input file '" // trim(path_input) // "' does not exist!"
+       call error(msg)
+    elseif (readable(1:3) /= 'YES') then
+       msg = "Input file '" // trim(path_input) // "' is not readable! &
             &Change file permissions with chmod command."
-       call error( msg )
+       call error(msg)
     end if
 
   end subroutine read_command_line
@@ -65,40 +77,49 @@ contains
 
     do
        read(unit=in, fmt='(A250)', iostat=readError) line
-       if ( readError /= 0 ) exit
+       if (readError /= 0) exit
+       ! TODO: don't really need to split entire string, this could be
+       ! replaced by just finding first word
        call split_string(line, words, n)
-       if ( n == 0 ) cycle
+       if (n == 0) cycle
 
-       select case ( trim(words(1)) )
-       case ( 'cell' )
+       select case (trim(words(1)))
+       case ('cell')
           n_cells = n_cells + 1
-       case ( 'surface' )
+       case ('surface')
           n_surfaces = n_surfaces + 1
-       case ( 'material' ) 
+       case ('material') 
           n_materials = n_materials + 1
        end select
     end do
 
     ! Check to make sure there are cells, surface, and materials
     ! defined for the problem
-    if ( n_cells == 0 ) then
+    if (n_cells == 0) then
        msg = "No cells specified!"
-       call error( msg )
-    elseif ( n_surfaces == 0 ) then
+       close(unit=in)
+       call error(msg)
+    elseif (n_surfaces == 0) then
        msg = "No surfaces specified!"
-       call error( msg )
-    elseif ( n_materials == 0 ) then
+       close(unit=in)
+       call error(msg)
+    elseif (n_materials == 0) then
        msg = "No materials specified!"
-       call error( msg )
+       close(unit=in)
+       call error(msg)
     end if
 
     close(unit=in)
 
     ! Allocate arrays for cells, surfaces, etc.
-    allocate( cells(n_cells)         )
-    allocate( surfaces(n_surfaces)   )
-    allocate( materials(n_materials) )
-    
+    allocate(cells(n_cells))
+    allocate(surfaces(n_surfaces))
+    allocate(materials(n_materials))
+
+    ! Set up dictionaries
+    call dict_create(cell_dict)
+    call dict_create(surface_dict)
+
   end subroutine read_count
 
 !=====================================================================
@@ -134,26 +155,31 @@ contains
        call split_string_wl(line, words, n)
 
        ! skip comments
-       if ( n==0 .or. words(1)(1:1) == '#' ) cycle
+       if (n==0 .or. words(1)(1:1) == '#') cycle
 
-       select case ( trim(words(1)) )
+       select case (trim(words(1)))
        case ('cell')
           ! Read data from cell entry
           index_cell = index_cell + 1
-          call read_cell( index_cell, words, n )
+          call read_cell(index_cell, words, n)
           
        case ('surface')
           ! Read data
           index_surface = index_surface + 1
-          call read_surface( index_surface, words, n )
+          call read_surface(index_surface, words, n)
 
-       case ( 'material' )
+       case ('material')
+          index_material = index_material + 1
+          call read_material(index_material, words, n)
 
-       case ( 'source' )
-          call read_source( words, n )
+       case ('source')
+          call read_source(words, n)
 
-       case ( 'criticality' )
-          call read_criticality( words, n )
+       case ('criticality')
+          call read_criticality(words, n)
+
+       case ('xs_data')
+          path_xsdata = words(2)
 
        case default
           ! do nothing
@@ -163,13 +189,15 @@ contains
 
     close(unit=in)
 
+    ! call adjust_cell_index
+
   end subroutine read_input
 
 !=====================================================================
 ! READ_CELL parses the data on a cell card.
 !=====================================================================
 
-  subroutine read_cell( index, words, n_words )
+  subroutine read_cell(index, words, n_words)
 
     integer,      intent(in) :: index
     character(*), intent(in) :: words(max_words)
@@ -185,25 +213,26 @@ contains
 
     ! Read cell identifier
     read(words(2), fmt='(I8)', iostat=readError) this_cell%uid
-    if ( readError > 0 ) then
+    if (readError > 0) then
        msg = "Invalid cell name: " // words(2)
-       call error( msg )
+       call error(msg)
     end if
+    call dict_add_key(cell_dict, words(2), index)
 
     ! Read cell material
     read(words(3), fmt='(I8)', iostat=readError) this_cell%material
 
     ! Read list of surfaces
-    allocate( this_cell%boundary_list(n_words-3) )
+    allocate(this_cell%boundary_list(n_words-3))
     do i = 1, n_words-3
        word = words(i+3)
-       if ( word(1:1) == '(' ) then
+       if (word(1:1) == '(') then
           this_cell%boundary_list(i) = OP_LEFT_PAREN
-       elseif ( word(1:1) == ')' ) then
+       elseif (word(1:1) == ')') then
           this_cell%boundary_list(i) = OP_RIGHT_PAREN
-       elseif ( word(1:1) == ':' ) then
+       elseif (word(1:1) == ':') then
           this_cell%boundary_list(i) = OP_UNION
-       elseif ( word(1:1) == '#' ) then
+       elseif (word(1:1) == '#') then
           this_cell%boundary_list(i) = OP_DIFFERENCE
        else
           read(word, fmt='(I8)', iostat=readError) this_cell%boundary_list(i)
@@ -216,7 +245,7 @@ contains
 ! READ_SURFACE parses the data on a surface card.
 !=====================================================================
 
-  subroutine read_surface( index, words, n_words )
+  subroutine read_surface(index, words, n_words)
 
     integer, intent(in) :: index
     character(*), intent(in) :: words(max_words)
@@ -233,67 +262,67 @@ contains
 
     ! Read surface identifier
     read(words(2), fmt='(I8)', iostat=readError) this_surface%uid
-    if ( readError > 0 ) then
+    if (readError > 0) then
        msg = "Invalid surface name: " // words(2)
-       call error( msg )
+       call error(msg)
     end if
 
     ! Read surface type
     word = words(3)
     call lower_case(word)
-    select case ( trim(word) )
-       case ( 'px' )
+    select case (trim(word))
+       case ('px')
           this_surface%type = SURF_PX
           coeffs_reqd  = 1
-       case ( 'py' )
+       case ('py')
           this_surface%type = SURF_PY
           coeffs_reqd  = 1
-       case ( 'pz' )
+       case ('pz')
           this_surface%type = SURF_PZ
           coeffs_reqd  = 1
-       case ( 'plane' )
+       case ('plane')
           this_surface%type = SURF_PLANE
           coeffs_reqd  = 4
-       case ( 'cylx' )
+       case ('cylx')
           this_surface%type = SURF_CYL_X
           coeffs_reqd  = 3
-       case ( 'cyly' )
+       case ('cyly')
           this_surface%type = SURF_CYL_Y
           coeffs_reqd  = 3
-       case ( 'cylz' )
+       case ('cylz')
           this_surface%type = SURF_CYL_Z
           coeffs_reqd  = 3
-       case ( 'sph' )
+       case ('sph')
           this_surface%type = SURF_SPHERE
           coeffs_reqd  = 4
-       case ( 'boxx' )
+       case ('boxx')
           this_surface%type = SURF_BOX_X
           coeffs_reqd  = 4
-       case ( 'boxy' )
+       case ('boxy')
           this_surface%type = SURF_BOX_Y
           coeffs_reqd  = 4
-       case ( 'boxz' )
+       case ('boxz')
           this_surface%type = SURF_BOX_Z
           coeffs_reqd  = 4
-       case ( 'box' ) 
+       case ('box') 
           this_surface%type = SURF_BOX
           coeffs_reqd  = 6
-       case ( 'gq' )
+       case ('gq')
           this_surface%type = SURF_GQ
           coeffs_reqd  = 10
        case default
           msg = "Invalid surface type: " // words(3)
-          call error( msg )
+          call error(msg)
     end select
 
     ! Make sure there are enough coefficients for surface type
-    if ( n_words-3 < coeffs_reqd ) then
+    if (n_words-3 < coeffs_reqd) then
        msg = "Not enough coefficients for surface: " // words(2)
-       call error( msg )
+       call error(msg)
     end if
     
     ! Read list of surfaces
-    allocate( this_surface%coeffs(n_words-3) )
+    allocate(this_surface%coeffs(n_words-3))
     do i = 1, n_words-3
        this_surface%coeffs(i) = str_to_real(words(i+3))
     end do
@@ -304,7 +333,7 @@ contains
 ! READ_SOURCE parses the data on a source card.
 !=====================================================================
 
-  subroutine read_source( words, n_words )
+  subroutine read_source(words, n_words)
 
     character(*), intent(in) :: words(max_words)
     integer, intent(in) :: n_words
@@ -318,29 +347,135 @@ contains
     ! Read source type
     word = words(2)
     call lower_case(word)
-    select case ( trim(word) )
-    case ( 'box' )
+    select case (trim(word))
+    case ('box')
        external_source%type = SRC_BOX
        values_reqd = 6
     case default
        msg = "Invalid source type: " // words(2)
-       call error( msg )
+       call error(msg)
     end select
 
     ! Make sure there are enough values for this source type
-    if ( n_words-2 < values_reqd ) then
+    if (n_words-2 < values_reqd) then
        msg = "Not enough values for source of type: " // words(2)
-       call error( msg )
+       call error(msg)
     end if
     
     ! Read list of surfaces
-    allocate( external_source%values(n_words-2) )
+    allocate(external_source%values(n_words-2))
     do i = 1, n_words-2
        external_source%values(i) = str_to_real(words(i+2))
     end do
     
   end subroutine read_source
 
+!=====================================================================
+! READ_MATERIAL parses a material card, creating the appropriate
+! <Isotope>s and <Material>s.
+!=====================================================================
+
+  subroutine read_material(index, words, n_words)
+
+    integer, intent(in) :: index
+    character(*), intent(in) :: words(n_words)
+    integer, intent(in) :: n_words
+
+    character(100) :: msg
+    integer :: i
+    integer :: n_isotopes
+    type(Material), pointer :: mat
+
+    if (mod(n_words,2) /= 0) then
+       msg = "Invalid number of arguments for material: " // words(2)
+       call error(msg)
+    end if
+
+    n_isotopes = (n_words-2)/2
+    mat => materials(index)
+    mat%n_isotopes = n_isotopes
+
+    ! allocate isotope and density list
+    allocate(mat%names(n_isotopes))
+    allocate(mat%isotopes(n_isotopes))
+    allocate(mat%table(n_isotopes))
+    allocate(mat%atom_percent(n_isotopes))
+
+    ! read names and percentage
+    do i = 1, n_isotopes
+       mat%names(i) = words(2*i+1)
+       mat%atom_percent(i) = str_to_real(words(2*i+2))
+    end do
+
+  end subroutine read_material
+
+!=====================================================================
+! NORMALIZE_AO normalizes the atom or weight percentages for each
+! material
+!=====================================================================
+
+  subroutine normalize_ao()
+
+    type(xsData), pointer :: iso
+    type(Material), pointer :: mat
+    integer :: index
+    integer :: i, j
+    character(10) :: key
+    real(8) :: sum_percent
+    real(8) :: awr, w
+    character(100) :: msg
+    
+    ! first find the index in the xsdata array for each isotope in
+    ! each material
+    do i = 1, n_materials
+       mat => materials(i)
+
+       ! Check to make sure either all atom percents or all weight
+       ! percents are given
+       if (.not. (all(mat%atom_percent > 0) .or. & 
+            & all(mat%atom_percent < 0))) then
+          msg = "Cannot mix atom and weight percents in material " // &
+               & int_to_str(mat%uid)
+          call error(msg)
+       end if
+
+       ! If data given in a/o, normalize it by sum
+       if (mat%atom_percent(1) > 0) then
+          sum_percent = sum(mat%atom_percent)
+          mat%atom_percent = mat%atom_percent / sum_percent
+          do j = 1, mat%n_isotopes
+             key = mat%names(j)
+             index = dict_get_key(xsdata_dict, key)
+             mat%isotopes(j) = index
+          end do
+          cycle
+       end if
+
+       ! Otherwise, need to find atomic weight ratios for each isotope
+       ! in order to normalize to a/o. To convert weight percent to
+       ! atom percent, the formula is: x = (w/awr) / sum(w/awr)
+       sum_percent = 0
+       do j = 1, mat%n_isotopes
+          key = mat%names(j)
+          index = dict_get_key(xsdata_dict, key)
+          mat%isotopes(j) = index
+          awr = xsdatas(index)%awr
+          w = mat%atom_percent(j)
+
+          sum_percent = sum_percent + w/awr
+       end do
+       
+       ! change each isotope from w/o to a/o
+       do j = 1, mat%n_isotopes
+          index = mat%isotopes(j)
+          awr = xsdatas(index)%awr
+          w = mat%atom_percent(j)
+          mat%atom_percent(j) = (w/awr)/sum_percent
+       end do
+    end do
+
+  end subroutine normalize_ao
+
 !=====================================================================
 ! READ_XS_LIBRARY parses the data on a xs_library card. This card
 ! specifies what cross section library should be used by default
@@ -373,29 +508,121 @@ contains
 
     ! Read number of cycles
     read(words(2), fmt='(I8)', iostat=readError) n_cycles
-    if ( readError > 0 ) then
+    if (readError > 0) then
        msg = "Invalid number of cycles: " // words(2)
-       call error( msg )
+       call error(msg)
     end if
 
     ! Read number of inactive cycles
     read(words(3), fmt='(I8)', iostat=readError) n_inactive
-    if ( readError > 0 ) then
+    if (readError > 0) then
        msg = "Invalid number of inactive cycles: " // words(2)
-       call error( msg )
+       call error(msg)
     end if
 
     ! Read number of particles
     read(words(4), fmt='(I8)', iostat=readError) n_particles
     if ( readError > 0 ) then
        msg = "Invalid number of particles: " // words(2)
-       call error( msg )
+       call error(msg)
     end if
 
   end subroutine read_criticality
 
+!=====================================================================
+! READ_LINE reads a line from a file open on a unit
+!=====================================================================
+
+  subroutine read_line(unit, line, ioError)
+
+    integer,             intent(in)  :: unit
+    character(max_line), intent(out) :: line
+    integer,             intent(out) :: ioError
+
+    read(UNIT=unit, FMT='(A250)', IOSTAT=ioError) line
+
+  end subroutine read_line
+
+!=====================================================================
+! SKIP_LINES skips 'n_lines' lines from a file open on a unit
+!=====================================================================
+
+  subroutine skip_lines(unit, n_lines, ioError)
+
+    integer, intent(in) :: unit
+    integer, intent(in) :: n_lines
+    integer, intent(out) :: ioError
+
+    integer :: i
+    character(max_line) :: tmp
+
+    do i = 1, n_lines
+       read(UNIT=unit, fmt='(A250)', IOSTAT=ioError) tmp
+    end do
+
+  end subroutine skip_lines
+
+!=====================================================================
+! READ_DATA_INT reads integer data into an arrya from a file open
+!=====================================================================
+
+  subroutine read_data_int(in, array, n, lines, words_per_line)
+
+    integer, intent(in) :: in
+    integer, intent(out) :: array(n)
+    integer, intent(in) :: n
+    integer, intent(in) :: lines
+    integer, intent(in) :: words_per_line
+
+    integer :: i, j
+    character(250) :: line
+    character(32) :: words(max_words)
+    integer :: ioError
+    integer :: n_words
+    integer :: index
+    integer :: n_last
+
+    index = 1
+    do i = 1, lines
+       if (i == lines) then
+          read(in,*) array(index:n)
+       else
+          read(in,*) array(index:index+words_per_line-1)
+          index = index + words_per_line
+       end if
+    end do
+
+  end subroutine read_data_int
+
+!=====================================================================
+! READ_DATA_REAL reads real(8) data into an arrya from a file open
+!=====================================================================
+
+  subroutine read_data_real(in, array, n, lines, words_per_line)
+
+    integer, intent(in)  :: in
+    real(8), intent(out) :: array(n)
+    integer, intent(in)  :: n
+    integer, intent(in)  :: lines
+    integer, intent(in)  :: words_per_line
+
+    integer :: i, j
+    character(250) :: line
+    character(32) :: words(max_words)
+    integer :: ioError
+    integer :: n_words
+    integer :: index
+
+    index = 1
+    do i = 1, lines
+       if (i == lines) then
+          read(in,*) array(index:n)
+       else
+          read(in,*) array(index:index+words_per_line-1)
+          index = index + words_per_line
+       end if
+    end do
+
+  end subroutine read_data_real
+
 end module fileio
-
-    
-
-    
diff --git a/src/global.f90 b/src/global.f90
index 03c865bdc5..5b6c3ae431 100644
--- a/src/global.f90
+++ b/src/global.f90
@@ -5,13 +5,28 @@ module global
   implicit none
 
   ! Main arrays for cells, surfaces, materials
-  type(Cell),      allocatable, target :: cells(:)
-  type(Surface),   allocatable, target :: surfaces(:)
-  type(Material),  allocatable, target :: materials(:)
+  type(Cell),     allocatable, target :: cells(:)
+  type(Surface),  allocatable, target :: surfaces(:)
+  type(Material), allocatable, target :: materials(:)
+  type(Isotope),  allocatable, target :: isotopes(:)
+  type(xsData),   allocatable, target :: xsdatas(:)
   integer :: n_cells     ! # of cells
   integer :: n_surfaces  ! # of surfaces
   integer :: n_materials ! # of materials
 
+  type(Dictionary), pointer :: cell_dict
+  type(Dictionary), pointer :: surface_dict
+  type(Dictionary), pointer :: xsdata_dict
+  type(Dictionary), pointer :: isotope_dict
+  type(Dictionary), pointer :: ace_dict
+
+  ! Cross section arrays
+  type(AceContinuous), pointer :: xs_continuous(:)
+  type(AceThermal),    allocatable, target :: xs_thermal(:)
+
+  ! unionized energy grid
+  real(8), allocatable :: energy_grid(:)
+
   ! Histories/cycles/etc for both external source and criticality
   integer :: n_particles ! # of particles (per cycle for criticality)
   integer :: n_cycles    ! # of cycles
@@ -24,10 +39,19 @@ module global
   type(Neutron), allocatable, target :: source_bank(:)
   type(Bank),    allocatable, target :: fission_bank(:)
 
+  ! Paths to input file, cross section data, etc
+  character(100) :: & 
+       & path_input,    &
+       & path_xsdata
+
   ! Physical constants
-  real(8), parameter ::           &
-       & PI       = 2.*acos(0.0), &  ! pi
-       & INFINITY = huge(0.0_8)      ! positive infinity
+  real(8), parameter ::            &
+       & PI           = 2.*acos(0.0),   & ! pi
+       & MASS_NEUTRON = 1.0086649156,   & ! mass of a neutron
+       & MASS_PROTON  = 1.00727646677,  & ! mass of a proton
+       & AMU          = 1.66053873e-27, & ! 1 amu in kg
+       & K_BOLTZMANN  = 8.617342e-5,    & ! Boltzmann constant in eV/K
+       & INFINITY  = huge(0.0_8)          ! positive infinity
 
   ! Boundary conditions
   integer, parameter ::    &
@@ -81,18 +105,14 @@ module global
        & PHOTON_   = 2, &
        & ELECTRON_ = 3
 
-  character(32) :: inputfile
-
   integer :: verbosity = 5
   integer, parameter :: max_words = 100
+  integer, parameter :: max_line  = 250
 
   ! Versioning numbers
   integer, parameter :: VERSION_MAJOR = 0
   integer, parameter :: VERSION_MINOR = 1
-  integer, parameter :: VERSION_RELEASE = 1
-
-  ! Key length for dictionary
-  integer, parameter :: DICT_KEY_LENGTH = 20
+  integer, parameter :: VERSION_RELEASE = 2
 
 contains
 
@@ -118,23 +138,16 @@ contains
   subroutine free_memory()
 
     ! Deallocate cells, surfaces, materials
-    if (allocated(cells)) then
-       deallocate(cells)
-    end if
-    if (allocated(surfaces)) then
-       deallocate(surfaces)
-    end if
-    if (allocated(materials)) then
-       deallocate(materials)
-    end if
+    if (allocated(cells)) deallocate(cells)
+    if (allocated(surfaces)) deallocate(surfaces)
+    if (allocated(materials)) deallocate(materials)
+
+    ! Deallocate xsdata list
+    if (allocated(xsdatas)) deallocate(xsdatas)
 
     ! Deallocate fission and source bank
-    if (allocated(fission_bank)) then
-       deallocate(fission_bank)
-    end if
-    if (allocated(source_bank)) then
-       deallocate(source_bank)
-    end if
+    if (allocated(fission_bank)) deallocate(fission_bank)
+    if (allocated(source_bank)) deallocate(source_bank)
 
     ! End program
     stop
@@ -146,14 +159,36 @@ contains
 ! to integers less than 10 billion.
 !=====================================================================
 
-  function int_to_str( num )
+  function int_to_str(num) result(str)
 
     integer, intent(in) :: num
-    character(10) :: int_to_str
+    character(10) :: str
 
-    write ( int_to_str, '(I10)' ) num
-    int_to_str = adjustl(int_to_str)
+    write ( str, '(I10)' ) num
+    str = adjustl(str)
 
   end function int_to_str
 
+!=====================================================================
+! INT_TO_STR converts a string to an integer. 
+!=====================================================================
+
+  function str_to_int(str) result(num)
+
+    character(*), intent(in) :: str
+    integer :: num
+    
+    character(5) :: fmt
+    integer :: w
+
+    w = len_trim(str)
+    
+    ! Create format specifier for reading string
+    write(fmt, '("(I",I2,")")') w
+
+    ! read string into integer
+    read(str,fmt) num
+
+  end function str_to_int
+
 end module global
diff --git a/src/input_sample b/src/input_sample
index fc2720427c..0988914339 100644
--- a/src/input_sample
+++ b/src/input_sample
@@ -6,9 +6,10 @@ cell 2 1    2 -1
 surface 1 sph 0 0 0 5
 surface 2 sph 0 0.0 0.0 3
 
-material 1 U235 1.00
+material 1   1001.03c 2.0   8016.03c 1.0 92235.12c 3.0
 
 source box -3 -3 -3  3 3 3
 
 xs_library endfb7
+xs_data /opt/serpent/xsdata/endfb7
 criticality 1 1 10
diff --git a/src/main.f90 b/src/main.f90
index d301fbc65c..7370abdbc8 100644
--- a/src/main.f90
+++ b/src/main.f90
@@ -1,7 +1,8 @@
 program main
 
   use global
-  use fileio,      only: read_input, read_command_line, read_count
+  use fileio,      only: read_input, read_command_line, read_count, &
+       &                 normalize_ao
   use output,      only: title, echo_input, message, warning, error
   use geometry,    only: sense, cell_contains
   use mcnp_random, only: RN_init_problem, rang, RN_init_particle
@@ -9,6 +10,9 @@ program main
   use physics,     only: transport
   use data_structures, only: Dictionary, dict_create, dict_add_key, & 
        &                     dict_get_key
+  use cross_section, only: read_xsdata
+  use ace, only: read_xs
+  use energy_grid, only: unionized_grid
 
   implicit none
 
@@ -17,6 +21,7 @@ program main
 
   ! Print the OpenMC title and version/date/time information
   call title()
+  verbosity = 10
 
   ! Initialize random number generator
   call RN_init_problem( 3, 0_8, 0_8, 0_8, 0 )
@@ -30,15 +35,21 @@ program main
   ! Read input file -- make a first pass through the file to count
   ! cells, surfaces, etc in order to allocate arrays, then do a second
   ! pass to actually read values
-  call read_count(inputfile)
-  call read_input(inputfile)
+  call read_count(path_input)
+  call read_input(path_input)
+  call read_xsdata(path_xsdata)
+  call normalize_ao()
+  call read_xs()
+  call unionized_grid()
+
+  stop
+  
   call echo_input()
 
   ! create source particles
   call init_source()
 
   ! start problem
-  verbosity = 10
   surfaces(1)%bc = BC_VACUUM
   call run_problem()
 
diff --git a/src/output.f90 b/src/output.f90
index d08b7fd104..e864ca432b 100644
--- a/src/output.f90
+++ b/src/output.f90
@@ -203,13 +203,13 @@ module output
       
       if ( val(2) < 10 ) then
          if ( val(3) < 10 ) then
-            today_date = date(6:6) // "/" // date(7:7) // "/" // date(1:4)
+            today_date = date(6:6) // "/" // date(8:8) // "/" // date(1:4)
          else
             today_date = date(6:6) // "/" // date(7:8) // "/" // date(1:4)
          end if
       else
          if ( val(3) < 10 ) then
-            today_date = date(5:6) // "/" // date(7:7) // "/" // date(1:4)
+            today_date = date(5:6) // "/" // date(8:8) // "/" // date(1:4)
          else
             today_date = date(5:6) // "/" // date(7:8) // "/" // date(1:4)
          end if
diff --git a/src/search.f90 b/src/search.f90
new file mode 100644
index 0000000000..6763d593c3
--- /dev/null
+++ b/src/search.f90
@@ -0,0 +1,54 @@
+module search
+
+contains
+
+!=====================================================================
+! BINARY_SEARCH performs a binary search of an array to find where a
+! specific value lies in the array. This is used extensively for
+! energy grid searching
+!=====================================================================
+
+  function binary_search( array, n, val ) result( index )
+
+    real(8), intent(in) :: array(n)
+    integer, intent(in) :: n
+    real(8), intent(in) :: val
+    integer :: index
+
+    integer :: L
+    integer :: R
+    real(8) :: testval
+
+    L = 1
+    R = n
+
+    if (val < array(L) .or. val > array(R)) then
+       ! error
+    end if
+    
+    do while (R - L > 1)
+       
+       ! Check boundaries
+       if (val > array(L) .and. val < array(L+1)) then
+          index = L
+          return
+       elseif (val > array(R+1) .and. val < array(R)) then
+          index = R
+          return
+       end if
+
+       ! Find values at midpoint
+       index = L + (R - L)/2
+       testval = array(index)
+       if (val > testval) then
+          L = index + 1
+       elseif (val < testval) then
+          R = index - 1
+       end if
+    end do
+
+    index = L
+
+  end function binary_search
+
+end module search
diff --git a/src/string.f90 b/src/string.f90
index d0f3911238..8dc2f0dc51 100644
--- a/src/string.f90
+++ b/src/string.f90
@@ -1,7 +1,7 @@
 module string
 
   use global, only: max_words
-  use output, only: error
+  use output, only: error, warning
 
   implicit none
 
@@ -27,6 +27,7 @@ contains
     integer       :: i       ! current index
     integer       :: i_start ! starting index of word
     integer       :: i_end   ! ending index of word
+    character(250) :: msg
 
     i_start = 0
     i_end = 0
@@ -43,6 +44,11 @@ contains
           if (i == len(string))   i_end = i
           if (i_end > 0) then
              n = n + 1
+             if ( i_end - i_start + 1 > len(words(n)) ) then
+                msg = "The word '" // string(i_start:i_end) // "' is longer than " &
+                     & // "the space allocated for it."
+                call warning( msg )
+             end if
              words(n) = string(i_start:i_end)
              ! reset indices
              i_start = 0
diff --git a/src/types.f90 b/src/types.f90
index a66c55156f..f3f9777af2 100644
--- a/src/types.f90
+++ b/src/types.f90
@@ -2,7 +2,12 @@ module types
 
   implicit none
 
-  type :: Surface
+!=====================================================================
+! SURFACE type defines a first- or second-order surface that can be
+! used to construct closed volumes (cells)
+!=====================================================================
+
+  type Surface
      integer :: uid
      integer :: type
      real(8), allocatable :: coeffs(:)
@@ -11,14 +16,23 @@ module types
      integer :: bc
   end type Surface
 
-  type :: Cell
+!=====================================================================
+! CELL defines a closed volume by its bounding surfaces
+!=====================================================================
+
+  type Cell
      integer :: uid
      integer :: n_items
      integer, allocatable :: boundary_list(:)
      integer :: material
   end type Cell
 
-  type :: Neutron
+!=====================================================================
+! NEUTRON describes the state of a neutron being transported through
+! the geometry
+!=====================================================================
+
+  type Neutron
     integer  :: uid     ! Unique ID
     real(8)  :: xyz(3)  ! location
     real(8)  :: uvw(3)  ! directional cosines
@@ -28,27 +42,188 @@ module types
     logical  :: alive   ! is particle alive?
   end type Neutron
 
-  type :: Bank
+!=====================================================================
+! BANK is used for storing fission sites in criticality
+! calculations. Since all the state information of a neutron is not
+! needed, this type allows sites to be stored with less memory
+!=====================================================================
+
+  type Bank
      integer :: uid    ! Unique ID
      real(8) :: xyz(3) ! Location of bank particle
   end type Bank
 
-  type :: Isotope
+!=====================================================================
+! ISOTOPE describes an isotope, e.g. U-235, within a material. Note
+! that two separate variables must be used for the same isotope in two
+! different materials since they will generally have different
+! densities
+!=====================================================================
+
+  type Isotope
      integer :: uid     ! unique identifier
      integer :: zaid    ! ZAID, e.g. 92235
      character(3) :: xs ! cross section identifier, e.g. 70c
      real(8) :: density ! density in atom/b-cm
   end type Isotope
 
-  type :: Material
+!=====================================================================
+! MATERIAL describes a material by its constituent isotopes
+!=====================================================================
+
+  type Material
      integer :: uid
      integer :: n_isotopes
-     type(Isotope), allocatable :: isotopes(:)
+     character(10), allocatable :: names(:)
+     integer, allocatable :: isotopes(:)
+     integer, allocatable :: table(:)
+     real(8), allocatable :: atom_percent(:)
+     integer :: sab_table
   end type Material
 
-  type :: ExtSource
+!=====================================================================
+! EXTSOURCE describes an external source of neutrons for a
+! fixed-source problem or for the starting source in a criticality
+! problem
+!=====================================================================
+
+  type ExtSource
      integer :: type                    ! type of source, e.g. 'box' or 'cell'
      real(8), allocatable :: values(:)  ! values for particular source type
   end type ExtSource
 
+!=====================================================================
+! ACEREACTION contains the cross-section and secondary energy and
+! angle distributions for a single reaction in a continuous-energy
+! ACE-format table
+!=====================================================================
+
+  type AceReaction
+     integer :: MT
+     real(8) :: Q_value
+     real(8) :: TY
+     integer :: energy_index
+     real(8), allocatable :: sigma(:)
+     real(8), allocatable :: ang_cos(:,:)
+     real(8), allocatable :: ang_pdf(:,:)
+     real(8), allocatable :: ang_cdf(:,:)
+  end type AceReaction
+
+!=====================================================================
+! ACECONTINUOUS contains all the data for an ACE-format
+! continuous-energy cross section. The ACE format (A Compact ENDF
+! format) is used in MCNP and several other Monte Carlo codes.
+!=====================================================================
+
+  type AceContinuous
+     character(20) :: name
+     real(8) :: awr
+     real(8) :: temp
+     real(8), allocatable :: energy(:)
+     real(8), allocatable :: sigma_t(:)
+     real(8), allocatable :: sigma_a(:)
+     real(8), allocatable :: sigma_el(:)
+     real(8), allocatable :: heating(:)
+
+     integer :: nu_p_type
+     real(8), allocatable :: nu_p_energy(:)
+     real(8), allocatable :: nu_p_value(:)
+     
+     integer :: nu_t_type
+     real(8), allocatable :: nu_t_energy(:)
+     real(8), allocatable :: nu_t_value(:)
+
+     real(8), allocatable :: nu_d_energy(:)
+     real(8), allocatable :: nu_d_value(:)
+     real(8), allocatable :: nu_d_precursor_const(:,:)
+     real(8), allocatable :: nu_d_precursor_energy(:,:)
+     real(8), allocatable :: nu_d_precursor_prob(:,:)
+     type(AceReaction), pointer :: reactions(:)
+
+  end type AceContinuous
+
+!=====================================================================
+! ACETHERMAL contains S(a,b) data for thermal neutron scattering,
+! typically off of light isotopes such as water, graphite, Be, etc
+!=====================================================================
+     
+  type AceThermal
+     character(20) :: name
+     real(8) :: awr
+     real(8) :: temperature
+     real(8), allocatable :: inelastic_e_in(:)
+     real(8), allocatable :: inelastic_sigma(:) 
+     real(8), allocatable :: inelastic_e_out(:,:)
+     real(8), allocatable :: inelastic_mu_out(:,:)
+     real(8), allocatable :: elastic_e_in(:)
+     real(8), allocatable :: elastic_P(:)
+     real(8), allocatable :: elastic_mu_out(:)
+  end type AceThermal
+
+!=====================================================================
+! LISTDATA Data stored in a linked list. In this case, we store the
+! (key,value) pair for a dictionary. Note that we need to store the
+! key in addition to the value for collision resolution.
+!=====================================================================
+
+  ! Key length for dictionary
+  integer, parameter :: DICT_KEY_LENGTH = 20
+
+  type ListData
+     character(len=DICT_KEY_LENGTH) :: key
+     integer                        :: value
+  end type ListData
+
+!=====================================================================
+! LINKEDLIST stores a simple linked list
+!=====================================================================
+
+  type LinkedList
+     type(LinkedList), pointer :: next
+     type(ListData)            :: data
+  end type LinkedList
+
+!=====================================================================
+! LINKEDLISTGRID stores a sorted list of energies for the unionized
+! energy grid as a linked list
+!=====================================================================
+
+  type LinkedListGrid
+     type(LinkedListGrid), pointer :: next
+     real(8)                   :: energy
+  end type LinkedListGrid
+
+!=====================================================================
+! HASHLIST - Since it's not possible to directly do an array of
+! pointers, this derived type provides a pointer
+!=====================================================================
+
+  type HashList
+     type(LinkedList), pointer :: list
+  end type HashList
+
+!=====================================================================
+! DICTIONARY provides a dictionary data structure of (key,value) pairs
+!=====================================================================
+
+  type Dictionary
+     type(HashList), pointer :: table(:)
+  end type Dictionary
+
+!=====================================================================
+! XSDATA contains data read in from a SERPENT xsdata file
+!=====================================================================
+
+  type xsData
+     character(10) :: alias
+     character(10) :: id
+     integer :: type
+     integer :: zaid
+     integer :: isomeric
+     real(8) :: awr
+     real(8) :: temp
+     integer :: binary
+     character(100) :: path
+  end type xsData
+
 end module types
diff --git a/src/unittest.f90 b/src/unittest.f90
new file mode 100644
index 0000000000..07c4426c89
--- /dev/null
+++ b/src/unittest.f90
@@ -0,0 +1,49 @@
+program unittest
+
+  use global
+  use energy_grid, only: add_grid_points
+
+  call test_energy_grid()
+
+contains
+
+!=====================================================================
+! TEST_ENERGY_GRID is a test of the energy grid unionization
+! subroutine on a simplified set of data
+!=====================================================================
+
+  subroutine test_energy_grid()
+
+    type(LinkedListGrid), pointer :: list => null()
+    type(LinkedListGrid), pointer :: current => null()
+    real(8), allocatable :: energy(:)
+    integer :: n, i
+    integer :: n_list
+
+    n = 10
+    allocate(energy(n))
+    do i = 1,n
+       energy(i) = 1.1248_8*i
+    end do
+
+    ! create list
+    allocate(list)
+    current => list
+    current%energy = 4.23_8
+    allocate(current%next)
+    current => current%next
+    current%energy = 7.23_8
+
+    call add_grid_points(list, energy, n)
+
+    current => list
+    do while(associated(current))
+       print *,current%energy
+       current => current%next
+    end do
+
+    deallocate(energy)
+
+  end subroutine test_energy_grid
+
+end program unittest