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

Added floating point equality check for lattice and parent boundaries.

Browse source
This commit is contained in:
Paul Romano 2011-11-30 23:11:25 -05:00
parent d8af37eef7
commit 1d81b40bde
2 changed files with 33 additions and 13 deletions
Download patch file Download diff file Expand all files Collapse all files

3
src/constants.f90
View file

@ -21,6 +21,9 @@ module constants
! Used for surface current tallies
real(8), parameter :: TINY_BIT = 1e-8
! User for precision in geometry
real(8), parameter :: FP_PRECISION = 1e-9_8
! Maximum number of collisions/crossings
integer, parameter :: MAX_EVENTS = 10000
integer, parameter :: MAX_SAMPLE = 100000

43
src/geometry.f90
View file

@ -550,10 +550,12 @@ contains
i_x = p % index_x
i_y = p % index_y
if (i_x < 1 .or. i_x > lat % n_x) then
message = "Reached edge of lattice."
message = "Reached edge of lattice (" // trim(int_to_str(i_x)) // &
"," // trim(int_to_str(i_y)) // ")."
call fatal_error()
elseif (i_y < 1 .or. i_y > lat % n_y) then
message = "Reached edge of lattice."
message = "Reached edge of lattice (" // trim(int_to_str(i_x)) // &
"," // trim(int_to_str(i_y)) // ")."
call fatal_error()
end if
@ -634,7 +636,9 @@ contains
v = p % uvw(2)
w = p % uvw(3)
! loop over all surfaces
! ==========================================================================
! FIND MINIMUM DISTANCE TO SURFACES IN CELL AND PARENT
dist = INFINITY
do i = 1, n_surf
if (i <= n1) then
@ -911,10 +915,11 @@ contains
dist = d
surf = -expression(i)
end if
end do
! Check lattice surfaces
! ==========================================================================
! FIND MINIMUM DISTANCE TO LATTICE SURFACES
in_lattice = .false.
if (p % lattice > 0) then
lat => lattices(p % lattice)
@ -922,9 +927,8 @@ contains
x = p % xyz_local(1)
y = p % xyz_local(2)
z = p % xyz_local(3)
x0 = lat % width_x * 0.5
y0 = lat % width_y * 0.5
x0 = lat % width_x * 0.5_8
y0 = lat % width_y * 0.5_8
! left and right sides
if (u > 0) then
@ -932,9 +936,19 @@ contains
else
d = -(x + x0)/u
end if
if (d < dist) then
dist = d
in_lattice = .true.
! If the lattice boundary is coincident with the parent cell boundary,
! we need to make sure that the lattice is not selected. This is
! complicated by the fact that floating point may determine that one
! is closer than the other (can't check direct equality). Thus, the
! logic here checks whether the relative difference is within floating
! point precision.
if (d < dist) then
if (abs(d - dist)/dist >= FP_PRECISION) then
dist = d
in_lattice = .true.
end if
end if
! top and bottom sides
@ -943,9 +957,12 @@ contains
else
d = -(y + y0)/v
end if
if (d < dist) then
dist = d
in_lattice = .true.
if (abs(d - dist)/dist >= FP_PRECISION) then
dist = d
in_lattice = .true.
end if
end if
elseif (lat % type == LATTICE_HEX) then