Tweak ORANGE construction (#1178)

* Do not provide default tolerance for ORANGE input
* Take ORANGE input by rvalue to capture metadata
* Change default tracking tolerance so that the square doesn't round off for O(1)
* Force inline
* Fix json comparer output for real numbers
* Rearrange geometry test harnesses
* Fix latin

src/orange/OrangeInput.hh

@@ -128,7 +128,7 @@ struct OrangeInput
     Tolerance<> tol;
 
     //! Whether the unit definition is valid
-    explicit operator bool() const { return !universes.empty(); }
+    explicit operator bool() const { return !universes.empty() && tol; }
 };
 
 //---------------------------------------------------------------------------//

src/orange/OrangeInputIO.json.cc

@@ -550,6 +550,13 @@ void from_json(nlohmann::json const& j, OrangeInput& value)
     {
         j.at("tol").get_to(value.tol);
     }
+    else
+    {
+        value.tol = Tolerance<>::from_default();
+        CELER_LOG(debug) << "No input tolerance provided: setting default "
+                            "tolerance";
+    }
+    CELER_ENSURE(value);
 }
 
 //---------------------------------------------------------------------------//

src/orange/OrangeParams.cc

@@ -9,6 +9,7 @@
 
 #include <fstream>
 #include <initializer_list>
+#include <limits>
 #include <numeric>
 #include <utility>
 #include <vector>
@@ -114,16 +115,19 @@ OrangeParams::OrangeParams(G4VPhysicalVolume const*)
  *
  * Volume and surface labels must be unique for the time being.
  */
-OrangeParams::OrangeParams(OrangeInput input)
+OrangeParams::OrangeParams(OrangeInput&& input)
 {
     CELER_VALIDATE(input, << "input geometry is incomplete");
     CELER_VALIDATE(!input.universes.empty(),
                    << "input geometry has no universes");
 
-    if (!input.tol)
-    {
-        input.tol = Tolerance<>::from_default();
-    }
+    // Save global bounding box
+    bbox_ = [&input] {
+        auto& global = input.universes[orange_global_universe.unchecked_get()];
+        CELER_VALIDATE(std::holds_alternative<UnitInput>(global),
+                       << "global universe is not a SimpleUnit");
+        return std::get<UnitInput>(global).bbox;
+    }();
 
     // Create host data for construction, setting tolerances first
     HostVal<OrangeParamsData> host_data;
@@ -142,9 +146,9 @@ OrangeParams::OrangeParams(OrangeInput input)
             detail::UnitInserter{&insert_universe_base, &host_data},
             detail::RectArrayInserter{&insert_universe_base, &host_data}};
 
-        for (auto const& u : input.universes)
+        for (auto&& u : input.universes)
         {
-            std::visit(insert_universe, u);
+            std::visit(insert_universe, std::move(u));
         }
 
         univ_labels_ = LabelIdMultiMap<UniverseId>{std::move(universe_labels)};
@@ -161,17 +165,13 @@ OrangeParams::OrangeParams(OrangeInput input)
               [](SimpleUnitRecord const& su) { return su.simple_safety; })
           && host_data.rect_arrays.empty();
 
-    CELER_VALIDATE(std::holds_alternative<UnitInput>(input.universes.front()),
-                   << "global universe is not a SimpleUnit");
-    bbox_ = std::get<UnitInput>(input.universes.front()).bbox;
-
     // Update scalars *after* loading all units
     CELER_VALIDATE(host_data.scalars.max_logic_depth
                        < detail::LogicStack::max_stack_depth(),
                    << "input geometry has at least one volume with a "
                       "logic depth of"
                    << host_data.scalars.max_logic_depth
-                   << " (surfaces are nested too deeply); but the logic "
+                   << " (a volume's CSG tree is too deep); but the logic "
                       "stack is limited to a depth of "
                    << detail::LogicStack::max_stack_depth());
 

src/orange/OrangeParams.hh

@@ -45,7 +45,7 @@ class OrangeParams final : public GeoParamsSurfaceInterface,
     explicit OrangeParams(G4VPhysicalVolume const*);
 
     // ADVANCED usage: construct from explicit host data
-    explicit OrangeParams(OrangeInput input);
+    explicit OrangeParams(OrangeInput&& input);
 
     //! Whether safety distance calculations are accurate and precise
     bool supports_safety() const final { return supports_safety_; }

src/orange/OrangeTypes.cc

@@ -16,11 +16,31 @@
 
 namespace celeritas
 {
+namespace
+{
+//---------------------------------------------------------------------------//
+/*!
+ * Whether the square of the tolerance rounds to zero for O(1) operations.
+ *
+ * This is an important criterion for construction and tracking operations that
+ * involve \c sqrt, for example:
+ * - Rotation matrix simplification
+ * - Shape simplification
+ * - Tracking to quadric surfaces
+ */
+template<class T>
+constexpr bool is_tol_sq_nonnegligible(T value)
+{
+    return T{1} - ipow<2>(value) != T{1};
+}
+
+}  // namespace
+
 //---------------------------------------------------------------------------//
 /*!
  * Use a relative error of \f$ \sqrt(\epsilon_\textrm{machine}) \f$ .
  *
- * Technically we're rounding the machine epsilon to a nearby power of 10. We
+ * Technically we're rounding the machine epsilon to a nearby value. We
  * could use numeric_limits<real_type>::epsilon instead.
  */
 template<class T>
@@ -29,14 +49,16 @@ Tolerance<T> Tolerance<T>::from_default(real_type length)
     constexpr real_type sqrt_emach = [] {
         if constexpr (std::is_same_v<real_type, double>)
         {
-            return 1.e-8;
+            // std::sqrt(LimitsT::epsilon()) = 1.4901161193847656e-08
+            return 1.5e-8;
         }
         else if constexpr (std::is_same_v<real_type, float>)
         {
-            return 5.e-3f;
+            // std::sqrt(LimitsT::epsilon()) = 0.00034526698f
+            return 3e-4f;
         }
     }();
-    static_assert(real_type{1} - ipow<2>(sqrt_emach) != real_type{1},
+    static_assert(is_tol_sq_nonnegligible(sqrt_emach),
                   "default tolerance is too low");
 
     return Tolerance<T>::from_relative(sqrt_emach, length);
@@ -66,6 +88,7 @@ Tolerance<T> Tolerance<T>::from_relative(real_type rel, real_type length)
     CELER_VALIDATE(length > 0,
                    << "length scale " << length
                    << " is invalid [must be positive]");
+
     Tolerance<T> result;
     result.rel = rel;
     result.abs = rel * length;

src/orange/OrangeTypes.hh

@@ -18,7 +18,6 @@
 #include "corecel/Types.hh"
 #include "corecel/cont/Array.hh"
 #include "corecel/math/NumericLimits.hh"
-#include "geocel/Types.hh"
 #include "geocel/Types.hh"  // IWYU pragma: export
 
 namespace celeritas
@@ -289,6 +288,8 @@ struct Daughter
  * \note For historical reasons, the absolute tolerance used by \c SoftEqual
  * defaults to 1/100 of the relative tolerance, whereas with \c Tolerance the
  * equivalent behavior is setting a length scale of 0.01.
+ *
+ * \todo Move this to a separate file.
  */
 template<class T = ::celeritas::real_type>
 struct Tolerance

src/orange/detail/UnitInserter.cc

@@ -121,7 +121,7 @@ struct NumIntersectionGetter
 
 //---------------------------------------------------------------------------//
 //! Construct surface labels, empty if needed
-std::vector<Label> make_surface_labels(UnitInput const& inp)
+std::vector<Label> make_surface_labels(UnitInput& inp)
 {
     CELER_EXPECT(inp.surface_labels.empty()
                  || inp.surface_labels.size() == inp.surfaces.size());
@@ -131,24 +131,25 @@ std::vector<Label> make_surface_labels(UnitInput const& inp)
 
     for (auto i : range(inp.surface_labels.size()))
     {
-        Label surface_label = inp.surface_labels[i];
+        Label surface_label = std::move(inp.surface_labels[i]);
         if (surface_label.ext.empty())
         {
             surface_label.ext = inp.label.name;
         }
         result[i] = std::move(surface_label);
     }
+    inp.surface_labels.clear();
     return result;
 }
 
 //---------------------------------------------------------------------------//
 //! Construct volume labels from the input volumes
-std::vector<Label> make_volume_labels(UnitInput const& inp)
+std::vector<Label> make_volume_labels(UnitInput& inp)
 {
     std::vector<Label> result;
     for (auto const& v : inp.volumes)
     {
-        Label vl = v.label;
+        Label vl = std::move(v.label);
         if (vl.ext.empty())
         {
             vl.ext = inp.label.name;
@@ -196,7 +197,7 @@ UnitInserter::UnitInserter(UniverseInserter* insert_universe, Data* orange_data)
 /*!
  * Create a simple unit and return its ID.
  */
-UniverseId UnitInserter::operator()(UnitInput const& inp)
+UniverseId UnitInserter::operator()(UnitInput&& inp)
 {
     CELER_VALIDATE(inp,
                    << "simple unit '" << inp.label
@@ -296,10 +297,12 @@ UniverseId UnitInserter::operator()(UnitInput const& inp)
 
     CELER_ASSERT(unit);
     simple_units_.push_back(unit);
+    auto surf_labels = make_surface_labels(inp);
+    auto vol_labels = make_volume_labels(inp);
     return (*insert_universe_)(UniverseType::simple,
-                               inp.label,
-                               make_surface_labels(inp),
-                               make_volume_labels(inp));
+                               std::move(inp.label),
+                               std::move(surf_labels),
+                               std::move(vol_labels));
 }
 
 //---------------------------------------------------------------------------//

src/orange/detail/UnitInserter.hh

@@ -42,7 +42,7 @@ class UnitInserter
     UnitInserter(UniverseInserter* insert_universe, Data* orange_data);
 
     // Create a simple unit and store in in OrangeParamsData
-    UniverseId operator()(UnitInput const& inp);
+    UniverseId operator()(UnitInput&& inp);
 
   private:
     Data* orange_data_{nullptr};

src/orange/surf/LocalSurfaceVisitor.hh

@@ -152,9 +152,8 @@ CELER_FUNCTION T LocalSurfaceVisitor::make_surface(LocalSurfaceId id) const
  * collection.
  */
 template<class T, class U>
-CELER_FUNCTION T LocalSurfaceVisitor::get_item(Items<T> const& items,
-                                               ItemRange<T> const& offsets,
-                                               ItemId<U> item)
+CELER_FORCEINLINE_FUNCTION T LocalSurfaceVisitor::get_item(
+    Items<T> const& items, ItemRange<T> const& offsets, ItemId<U> item)
 {
     CELER_EXPECT(*offsets.end() <= items.size());
     CELER_EXPECT(item < offsets.size());

test/celeritas/ext/GeantVolumeMapper.test.cc

@@ -200,6 +200,7 @@ void NestedTest::build_orange()
 
     OrangeInput input;
     input.universes.push_back(std::move(ui));
+    input.tol = Tolerance<>::from_default();
     auto geo = std::make_shared<OrangeParams>(std::move(input));
 #if CELERITAS_CORE_GEO == CELERITAS_CORE_GEO_ORANGE
     geo_params_ = std::move(geo);

test/geocel/g4/GeantGeo.test.cc

@@ -45,36 +45,12 @@ class GeantGeoTest : public GeantGeoTestBase
     virtual SpanStringView expected_log_levels() const { return {}; }
 };
 
+//---------------------------------------------------------------------------//
 class FourLevelsTest : public GeantGeoTest
 {
     std::string geometry_basename() const override { return "four-levels"; }
 };
 
-class SolidsTest : public GeantGeoTest
-{
-    std::string geometry_basename() const override { return "solids"; }
-
-    SpanStringView expected_log_levels() const final
-    {
-        static std::string_view const levels[] = {"error"};
-        return make_span(levels);
-    }
-};
-
-class CmseTest : public GeantGeoTest
-{
-  public:
-    std::string geometry_basename() const override { return "cmse"; }
-};
-
-class ZnenvTest : public GeantGeoTest
-{
-  public:
-    std::string geometry_basename() const override { return "znenv"; }
-};
-
-//---------------------------------------------------------------------------//
-
 TEST_F(FourLevelsTest, accessors)
 {
     auto const& geom = *this->geometry();
@@ -351,6 +327,16 @@ TEST_F(FourLevelsTest, safety)
 }
 
 //---------------------------------------------------------------------------//
+class SolidsTest : public GeantGeoTest
+{
+    std::string geometry_basename() const override { return "solids"; }
+
+    SpanStringView expected_log_levels() const final
+    {
+        static std::string_view const levels[] = {"error"};
+        return make_span(levels);
+    }
+};
 
 TEST_F(SolidsTest, accessors)
 {
@@ -605,6 +591,11 @@ TEST_F(SolidsTest, reflected_vol)
 }
 
 //---------------------------------------------------------------------------//
+class CmseTest : public GeantGeoTest
+{
+  public:
+    std::string geometry_basename() const override { return "cmse"; }
+};
 
 TEST_F(CmseTest, trace)
 {
@@ -678,6 +669,11 @@ TEST_F(CmseTest, trace)
 }
 
 //---------------------------------------------------------------------------//
+class ZnenvTest : public GeantGeoTest
+{
+  public:
+    std::string geometry_basename() const override { return "znenv"; }
+};
 
 TEST_F(ZnenvTest, trace)
 {

test/geocel/vg/Vecgeom.test.cc

@@ -1197,18 +1197,6 @@ class SolidsGeantTest : public VecgeomGeantTestBase
 
 //---------------------------------------------------------------------------//
 
-#define ZnenvGeantTest TEST_IF_CELERITAS_GEANT(ZnenvGeantTest)
-class ZnenvGeantTest : public VecgeomGeantTestBase
-{
-  public:
-    SPConstGeo build_geometry() final
-    {
-        return this->load_g4_gdml("znenv.gdml");
-    }
-};
-
-//---------------------------------------------------------------------------//
-
 TEST_F(SolidsGeantTest, DISABLED_dump)
 {
     this->geometry();
@@ -1472,6 +1460,18 @@ TEST_F(SolidsGeantTest, reflected_vol)
 
 //---------------------------------------------------------------------------//
 
+#define ZnenvGeantTest TEST_IF_CELERITAS_GEANT(ZnenvGeantTest)
+class ZnenvGeantTest : public VecgeomGeantTestBase
+{
+  public:
+    SPConstGeo build_geometry() final
+    {
+        return this->load_g4_gdml("znenv.gdml");
+    }
+};
+
+//---------------------------------------------------------------------------//
+
 TEST_F(ZnenvGeantTest, trace)
 {
     // NOTE: This tests the capability of the G4PVDivision conversion based on