From 63dd99c089ed922628600f3dafdac729a608d3e0 Mon Sep 17 00:00:00 2001
From: Patrick Shriwise <pshriwise@gmail.com>
Date: Wed, 3 Apr 2019 13:51:35 -0500
Subject: [PATCH] Addressing @smharper's comments.

---
 include/openmc/geometry.h |  4 +++-
 src/geometry.cpp          |  6 +-----
 src/lattice.cpp           | 28 ++++++++++++----------------
 3 files changed, 16 insertions(+), 22 deletions(-)

diff --git a/include/openmc/geometry.h b/include/openmc/geometry.h
index 39f500b847..8e6e8a82ef 100644
--- a/include/openmc/geometry.h
+++ b/include/openmc/geometry.h
@@ -39,7 +39,9 @@ struct BoundaryInfo {
 //! Check two distances by coincidence tolerance
 //==============================================================================
 
-bool coincident(double d1, double d2);
+inline bool coincident(double d1, double d2) {
+  return std::abs(d1 - d2) < FP_COINCIDENT;
+}
 
 //==============================================================================
 //! Check for overlapping cells at a particle's position.
diff --git a/src/geometry.cpp b/src/geometry.cpp
index 225a4970f4..d70e1ed0bf 100644
--- a/src/geometry.cpp
+++ b/src/geometry.cpp
@@ -32,10 +32,6 @@ std::vector<int64_t> overlap_check_count;
 // Non-member functions
 //==============================================================================
 
-bool coincident(double d1, double d2) {
-  return std::abs(d1 - d2) < FP_COINCIDENT;
-}
-
 bool check_cell_overlap(Particle* p)
 {
   int n_coord = p->n_coord_;
@@ -400,7 +396,7 @@ BoundaryInfo distance_to_boundary(Particle* p)
     // a higher level then we need to make sure that the higher level boundary
     // is selected.  This logic must consider floating point precision.
     double& d = info.distance;
-    if (d_surf < d_lat && !coincident(d_surf, d_lat)) {
+    if (d_surf < d_lat - FP_COINCIDENT) {
       if (d == INFINITY || (d - d_surf)/d >= FP_REL_PRECISION) {
         d = d_surf;
 
diff --git a/src/lattice.cpp b/src/lattice.cpp
index 87ec3c83fd..dc1b1d8fe8 100644
--- a/src/lattice.cpp
+++ b/src/lattice.cpp
@@ -711,16 +711,23 @@ const
 std::array<int, 3>
 HexLattice::get_indices(Position r, Direction u) const
 {
-  // result variables
-  int ix = 0;
-  int ia = 0;
-  int iz = 0;
-
   // Offset the xyz by the lattice center.
   Position r_o {r.x - center_.x, r.y - center_.y, r.z};
   if (is_3d_) {r_o.z -= center_.z;}
 
+  // Convert coordinates into skewed bases.  The (x, alpha) basis is used to
+  // find the index of the global coordinates to within 4 cells.
+  double alpha = r_o.y - r_o.x / std::sqrt(3.0);
+  int ix = std::floor(r_o.x / (0.5*std::sqrt(3.0) * pitch_[0]));
+  int ia = std::floor(alpha / pitch_[0]);
+
+  // Add offset to indices (the center cell is (i_x, i_alpha) = (0, 0) but
+  // the array is offset so that the indices never go below 0).
+  ix += n_rings_-1;
+  ia += n_rings_-1;
+
   // Index the z direction, accounting for coincidence
+  int iz{};
   if (is_3d_) {
     double iz_ {r_o.z / pitch_[1] + 0.5 * n_axial_};
     long iz_close {std::lround(iz_)};
@@ -731,17 +738,6 @@ HexLattice::get_indices(Position r, Direction u) const
     }
   }
 
-  // Convert coordinates into skewed bases.  The (x, alpha) basis is used to
-  // find the index of the global coordinates to within 4 cells.
-  double alpha = r_o.y - r_o.x / std::sqrt(3.0);
-  ix = std::floor(r_o.x / (0.5*std::sqrt(3.0) * pitch_[0]));
-  ia = std::floor(alpha / pitch_[0]);
-
-  // Add offset to indices (the center cell is (i_x, i_alpha) = (0, 0) but
-  // the array is offset so that the indices never go below 0).
-  ix += n_rings_-1;
-  ia += n_rings_-1;
-
   // Calculate the (squared) distance between the particle and the centers of
   // the four possible cells.  Regular hexagonal tiles form a Voronoi
   // tessellation so the xyz should be in the hexagonal cell that it is closest