sample_external_source now checks to make sure the sampled xyz exists within the problem domain. If it does not, the source site is resampled.

src/DEPENDENCIES

@@ -329,6 +329,7 @@ solver_interface.o: vector_header.o
 source.o: bank_header.o
 source.o: constants.o
 source.o: error.o
+source.o: geometry.o
 source.o: geometry_header.o
 source.o: global.o
 source.o: math.o

src/source.F90

@@ -3,6 +3,7 @@ module source
   use bank_header,     only: Bank
   use constants
   use error,           only: fatal_error
+  use geometry,        only: find_cell
   use geometry_header, only: BASE_UNIVERSE
   use global
   use math,            only: maxwell_spectrum, watt_spectrum
@@ -73,6 +74,8 @@ contains
     real(8) :: p_max(3)   ! maximum coordinates of source
     real(8) :: a          ! Arbitrary parameter 'a'
     real(8) :: b          ! Arbitrary parameter 'b'
+    logical :: found      ! Does the source particle exist within geometry?
+    type(Particle) :: p         ! Temporary particle for using find_cell
 
     ! Set weight to one by default
     site % wgt = ONE
@@ -80,11 +83,24 @@ contains
     ! Sample position
     select case (external_source % type_space)
     case (SRC_SPACE_BOX)
-      ! Coordinates sampled uniformly over a box
-      p_min = external_source % params_space(1:3)
-      p_max = external_source % params_space(4:6)
-      r = (/ (prn(), i = 1,3) /)
-      site % xyz = p_min + r*(p_max - p_min)
+      ! Set particle defaults
+      call p % initialize()
+      ! Repeat sampling source location until a good site has been found
+      found = .false.
+      do while (.not.found)
+        ! Coordinates sampled uniformly over a box
+        p_min = external_source % params_space(1:3)
+        p_max = external_source % params_space(4:6)
+        r = (/ (prn(), i = 1,3) /)
+        site % xyz = p_min + r*(p_max - p_min)
+
+        ! Fill p with needed data
+        p % coord0 % xyz = site % xyz
+
+        ! Now search to see if location exists in geometry
+        call find_cell(p, found)
+
+      end do
 
     case (SRC_SPACE_POINT)
       ! Point source
@@ -146,7 +162,7 @@ contains
   end subroutine sample_external_source
 
 !===============================================================================
-! GET_SOURCE_PARTICLE returns the next source particle 
+! GET_SOURCE_PARTICLE returns the next source particle
 !===============================================================================
 
   subroutine get_source_particle(p, index_source)