diff --git a/src/orange/OrangeTypes.cc b/src/orange/OrangeTypes.cc index 0aaf28421c..fa2fc25385 100644 --- a/src/orange/OrangeTypes.cc +++ b/src/orange/OrangeTypes.cc @@ -125,6 +125,7 @@ char const* to_cstring(SurfaceType value) "kz", "sq", "gq", + "inv", // clang-format on }; return to_cstring_impl(value); diff --git a/src/orange/OrangeTypes.hh b/src/orange/OrangeTypes.hh index 8da6a57cfb..a070bf1180 100644 --- a/src/orange/OrangeTypes.hh +++ b/src/orange/OrangeTypes.hh @@ -102,8 +102,13 @@ enum class Sense : bool /*! * Enumeration for mapping surface classes to integers. * - * These are ordered by number of coefficients needed for their representation: - * 1 for `p.|sc|c.c`, 3 for `c.`, 4 for `[ps]|k.`, 7 for `sq`, and 10 for `gq`. + * These are ordered roughly by complexity. The storage requirement for + * corresponding surfaces are: + * - 1 for `p.|sc|c.c`, + * - 3 for `c.`, + * - 4 for `[ps]|k.`, + * - 7 for `sq`, and + * - 10 for `gq`. * * See \c orange/surf/SurfaceTypeTraits.hh for how these map to classes. */ @@ -126,6 +131,7 @@ enum class SurfaceType : unsigned char kz, //!< Cone parallel to Z axis sq, //!< Simple quadric gq, //!< General quadric + inv, //!< Involute size_ //!< Sentinel value for number of surface types }; @@ -211,6 +217,16 @@ enum class BoundaryResult : bool exiting = true }; +//---------------------------------------------------------------------------// +/*! + * Chirality of a twirly object (currently only Involute). + */ +enum class Chirality : bool +{ + left, //!< Sinistral, spiraling counterclockwise + right, //!< Dextral, spiraling clockwise +}; + //---------------------------------------------------------------------------// /*! * Volume logic encoding. diff --git a/src/orange/detail/SurfacesRecordBuilder.cc b/src/orange/detail/SurfacesRecordBuilder.cc index 17b29aa4a4..bd8004b5cd 100644 --- a/src/orange/detail/SurfacesRecordBuilder.cc +++ b/src/orange/detail/SurfacesRecordBuilder.cc @@ -39,6 +39,13 @@ auto SurfacesRecordBuilder::operator()(VecSurface const& surfaces) -> result_typ // Functor to save the surface type and data, and the data offset auto emplace_surface = [this](auto&& s) { + if constexpr (std::remove_reference_t::surface_type() + == SurfaceType::inv) + { + // See discussion on + // https://github.com/celeritas-project/celeritas/pull/1342 + CELER_NOT_IMPLEMENTED("runtime involute support"); + } types_.push_back(s.surface_type()); auto data = s.data(); auto real_range = reals_.insert_back(data.begin(), data.end()); diff --git a/src/orange/orangeinp/IntersectRegion.cc b/src/orange/orangeinp/IntersectRegion.cc index bd63cd8ebc..15e79b4ae7 100644 --- a/src/orange/orangeinp/IntersectRegion.cc +++ b/src/orange/orangeinp/IntersectRegion.cc @@ -17,9 +17,11 @@ #include "corecel/math/SoftEqual.hh" #include "geocel/BoundingBox.hh" #include "geocel/Types.hh" +#include "orange/OrangeTypes.hh" #include "orange/orangeinp/detail/PolygonUtils.hh" #include "orange/surf/ConeAligned.hh" #include "orange/surf/CylCentered.hh" +#include "orange/surf/Involute.hh" #include "orange/surf/PlaneAligned.hh" #include "orange/surf/SimpleQuadric.hh" #include "orange/surf/SphereCentered.hh" @@ -647,6 +649,82 @@ void InfWedge::output(JsonPimpl* j) const to_json_pimpl(j, *this); } +//---------------------------------------------------------------------------// +// Involute +//---------------------------------------------------------------------------// +/*! + * Construct with prarameters and half height. + */ +Involute::Involute(Real3 const& radii, + Real2 const& displacement, + Chirality sign, + real_type halfheight) + : radii_(radii), a_(displacement), t_bounds_(), sign_(sign), hh_{halfheight} +{ + CELER_VALIDATE(radii_[0] > 0, + << "nonpositive involute radius: " << radii_[0]); + CELER_VALIDATE(radii_[1] > radii_[0], + << "inner cylinder radius " << radii_[1] + << " is not greater than involute radius " << radii_[0]); + CELER_VALIDATE(radii_[2] > radii_[1], + << "outer cylinder radius " << radii_[2] + << " is not greater than inner cyl radius " << radii_[1]); + + CELER_VALIDATE(a_[1] > a_[0], + << "nonpositive delta displacment: " << a_[1] - a_[0]); + CELER_VALIDATE(hh_ > 0, << "nonpositive half-height: " << hh_); + + for (auto i : range(2)) + { + t_bounds_[i] = std::sqrt( + clamp_to_nonneg(ipow<2>(radii_[i + 1] / radii_[0]) - 1)); + } + auto outer_isect = t_bounds_[0] + 2 * constants::pi - (a_[1] - a_[0]); + CELER_VALIDATE(t_bounds_[1] < outer_isect, + << "radial bounds result in angular overlap: " + << outer_isect - t_bounds_[1]); +} + +//---------------------------------------------------------------------------// +/*! + * Build surfaces. + */ +void Involute::build(IntersectSurfaceBuilder& insert_surface) const +{ + using InvSurf = ::celeritas::Involute; + + insert_surface(Sense::outside, PlaneZ{-hh_}); + insert_surface(Sense::inside, PlaneZ{hh_}); + insert_surface(Sense::outside, CCylZ{radii_[1]}); + insert_surface(Sense::inside, CCylZ{radii_[2]}); + // Make an inside and outside involute + Real2 const xy{0, 0}; + auto sense = (sign_ == Chirality::right ? Sense::outside : Sense::inside); + static char const* names[] = {"invl", "invr"}; + + for (auto i : range(2)) + { + insert_surface(sense, + InvSurf{xy, + radii_[0], + eumod(a_[i], 2 * constants::pi), + sign_, + t_bounds_[0], + t_bounds_[1] + a_[1] - a_[0]}, + std::string{names[i]}); + sense = flip_sense(sense); + } +} + +//---------------------------------------------------------------------------// +/*! + * Write output to the given JSON object. + */ +void Involute::output(JsonPimpl* j) const +{ + to_json_pimpl(j, *this); +} + //---------------------------------------------------------------------------// // PARALLELEPIPED //---------------------------------------------------------------------------// @@ -687,18 +765,18 @@ void Parallelepiped::build(IntersectSurfaceBuilder& insert_surface) const constexpr auto Y = to_int(Axis::y); constexpr auto Z = to_int(Axis::z); - // cache trigonometric values + // Cache trigonometric values real_type sinth, costh, sinphi, cosphi, sinal, cosal; sincos(theta_, &sinth, &costh); sincos(phi_, &sinphi, &cosphi); sincos(alpha_, &sinal, &cosal); - // base vectors + // Base vectors auto a = hpr_[X] * Real3{1, 0, 0}; auto b = hpr_[Y] * Real3{sinal, cosal, 0}; auto c = hpr_[Z] * Real3{sinth * cosphi, sinth * sinphi, costh}; - // positioning the planes + // Position the planes auto xnorm = make_unit_vector(cross_product(b, c)); auto ynorm = make_unit_vector(cross_product(c, a)); auto xoffset = dot_product(a, xnorm); diff --git a/src/orange/orangeinp/IntersectRegion.hh b/src/orange/orangeinp/IntersectRegion.hh index 79062225a4..a72106d914 100644 --- a/src/orange/orangeinp/IntersectRegion.hh +++ b/src/orange/orangeinp/IntersectRegion.hh @@ -317,6 +317,53 @@ class InfWedge final : public IntersectRegionInterface Turn interior_; }; +//---------------------------------------------------------------------------// +/*! + * An involute "blade" centered on the origin. + * + * This is the intersection of two parallel involutes with a cylindrical shell. + * The three radii, which must be in ascending order, are that of the involute, + * the inner cylinder, and the outer cylinder. + * + * The "chirality" of the involute is viewed from the \em +z axis looking down: + * whether it spirals to the right or left. + */ +class Involute final : public IntersectRegionInterface +{ + public: + // Construct with radius + explicit Involute(Real3 const& radii, + Real2 const& displacement, + Chirality chirality, + real_type halfheight); + + // Build surfaces + void build(IntersectSurfaceBuilder&) const final; + + // Output to JSON + void output(JsonPimpl*) const final; + + //// ACCESSORS //// + + //! Radii: Rdius of involute, minimum radius, maximum radius + Real3 const& radii() const { return radii_; } + //! Displacement angle + Real2 const& displacement_angle() const { return a_; } + //! Angular bounds of involute + Real2 const& t_bounds() const { return t_bounds_; } + //! Chirality of involute: turning left or right + Chirality chirality() const { return sign_; } + //! Halfheight + real_type halfheight() const { return hh_; } + + private: + Real3 radii_; + Real2 a_; + Real2 t_bounds_; + Chirality sign_; + real_type hh_; +}; + //---------------------------------------------------------------------------// /*! * A general parallelepiped centered on the origin. diff --git a/src/orange/orangeinp/IntersectSurfaceBuilder.cc b/src/orange/orangeinp/IntersectSurfaceBuilder.cc index ec0225123c..b116e10cce 100644 --- a/src/orange/orangeinp/IntersectSurfaceBuilder.cc +++ b/src/orange/orangeinp/IntersectSurfaceBuilder.cc @@ -245,6 +245,7 @@ CSB_INSTANTIATE(ConeAligned); CSB_INSTANTIATE(ConeAligned); CSB_INSTANTIATE(SimpleQuadric); CSB_INSTANTIATE(GeneralQuadric); +CSB_INSTANTIATE(Involute); #undef CSB_INSTANTIATE //! \endcond diff --git a/src/orange/orangeinp/ObjectIO.json.cc b/src/orange/orangeinp/ObjectIO.json.cc index 028315e1bd..a3d87e00c5 100644 --- a/src/orange/orangeinp/ObjectIO.json.cc +++ b/src/orange/orangeinp/ObjectIO.json.cc @@ -209,6 +209,15 @@ void to_json(nlohmann::json& j, Sphere const& cr) { j = {{"_type", "sphere"}, SIO_ATTR_PAIR(cr, radius)}; } +void to_json(nlohmann::json& j, Involute const& cr) +{ + j = {{"_type", "involute"}, + SIO_ATTR_PAIR(cr, radii), + SIO_ATTR_PAIR(cr, displacement_angle), + SIO_ATTR_PAIR(cr, t_bounds), + SIO_ATTR_PAIR(cr, chirality), + SIO_ATTR_PAIR(cr, halfheight)}; +} //!@} //---------------------------------------------------------------------------// diff --git a/src/orange/orangeinp/ObjectIO.json.hh b/src/orange/orangeinp/ObjectIO.json.hh index 5a7654a9a6..e1a9e0191b 100644 --- a/src/orange/orangeinp/ObjectIO.json.hh +++ b/src/orange/orangeinp/ObjectIO.json.hh @@ -43,6 +43,7 @@ class InfWedge; class Parallelepiped; class Prism; class Sphere; +class Involute; //---------------------------------------------------------------------------// @@ -74,6 +75,7 @@ void to_json(nlohmann::json& j, InfWedge const& cr); void to_json(nlohmann::json& j, Parallelepiped const& cr); void to_json(nlohmann::json& j, Prism const& cr); void to_json(nlohmann::json& j, Sphere const& cr); +void to_json(nlohmann::json& j, Involute const& cr); //---------------------------------------------------------------------------// } // namespace orangeinp diff --git a/src/orange/orangeinp/detail/LocalSurfaceInserter.cc b/src/orange/orangeinp/detail/LocalSurfaceInserter.cc index 0594915021..a8848c7554 100644 --- a/src/orange/orangeinp/detail/LocalSurfaceInserter.cc +++ b/src/orange/orangeinp/detail/LocalSurfaceInserter.cc @@ -203,6 +203,7 @@ LSI_INSTANTIATE(ConeAligned); LSI_INSTANTIATE(ConeAligned); LSI_INSTANTIATE(SimpleQuadric); LSI_INSTANTIATE(GeneralQuadric); +LSI_INSTANTIATE(Involute); #undef LSI_INSTANTIATE //! \endcond diff --git a/src/orange/orangeinp/detail/SurfaceHashPoint.hh b/src/orange/orangeinp/detail/SurfaceHashPoint.hh index 8012136c94..a1c225d8fc 100644 --- a/src/orange/orangeinp/detail/SurfaceHashPoint.hh +++ b/src/orange/orangeinp/detail/SurfaceHashPoint.hh @@ -66,6 +66,11 @@ struct SurfaceHashPoint return norm(s.origin()); } + real_type operator()(Involute const& s) const + { + return s.displacement_angle(); + } + real_type operator()(SimpleQuadric const& s) const { return std::sqrt(s.zeroth()); diff --git a/src/orange/surf/FaceNamer.cc b/src/orange/surf/FaceNamer.cc index 6f94f25d3a..00b08b592f 100644 --- a/src/orange/surf/FaceNamer.cc +++ b/src/orange/surf/FaceNamer.cc @@ -111,6 +111,15 @@ std::string FaceNamer::Impl::operator()(Sphere const&) const return "s"; } +//---------------------------------------------------------------------------// +/*! + * Construct a name for an involute. + */ +std::string FaceNamer::Impl::operator()(Involute const&) const +{ + return "inv"; +} + //---------------------------------------------------------------------------// /*! * Construct a name for a cone. diff --git a/src/orange/surf/FaceNamer.hh b/src/orange/surf/FaceNamer.hh index 02169281d1..194f0bae5b 100644 --- a/src/orange/surf/FaceNamer.hh +++ b/src/orange/surf/FaceNamer.hh @@ -87,6 +87,8 @@ class FaceNamer std::string operator()(SimpleQuadric const&) const; std::string operator()(GeneralQuadric const&) const; + + std::string operator()(Involute const&) const; }; }; diff --git a/src/orange/surf/Involute.cc b/src/orange/surf/Involute.cc index 5e86b31e5e..7b9185961c 100644 --- a/src/orange/surf/Involute.cc +++ b/src/orange/surf/Involute.cc @@ -23,7 +23,7 @@ namespace celeritas Involute::Involute(Real2 const& origin, real_type radius, real_type displacement, - Sign sign, + Chirality sign, real_type tmin, real_type tmax) : origin_(origin), r_b_(radius), a_(displacement), tmin_(tmin), tmax_(tmax) @@ -33,7 +33,7 @@ Involute::Involute(Real2 const& origin, CELER_EXPECT(tmin_ >= 0); CELER_EXPECT(tmax_ > tmin_ && tmax_ < 2 * constants::pi + tmin_); - if (sign) + if (sign == Chirality::right) { a_ = constants::pi - a_; r_b_ = -r_b_; diff --git a/src/orange/surf/Involute.hh b/src/orange/surf/Involute.hh index ccff066b9f..70206142b2 100644 --- a/src/orange/surf/Involute.hh +++ b/src/orange/surf/Involute.hh @@ -59,17 +59,17 @@ class Involute //! \name Type aliases using Intersections = Array; using StorageSpan = Span; - using Sign = detail::InvoluteSolver::Sign; using Real2 = Array; //@} public: //// CLASS ATTRIBUTES //// - //! Surface type identifier - static SurfaceType surface_type() + // Surface type identifier + static CELER_CONSTEXPR_FUNCTION SurfaceType surface_type() + { - CELER_NOT_IMPLEMENTED("runtime involute"); + return SurfaceType::inv; } //! Safety cannot be trivially calculated @@ -81,7 +81,7 @@ class Involute explicit Involute(Real2 const& origin, real_type radius, real_type displacement, - Sign sign, + Chirality sign, real_type tmin, real_type tmax); @@ -98,10 +98,10 @@ class Involute CELER_FUNCTION real_type r_b() const { return std::fabs(r_b_); } //! Displacement angle - CELER_FUNCTION real_type a() const { return a_; } + CELER_FUNCTION real_type displacement_angle() const { return a_; } // Orientation of the involute curve - inline CELER_FUNCTION Sign sign() const; + inline CELER_FUNCTION Chirality sign() const; //! Get bounds of the involute CELER_FUNCTION real_type tmin() const { return tmin_; } @@ -155,9 +155,9 @@ CELER_FUNCTION Involute::Involute(Span data) /*! * Orientation of the involute curve. */ -CELER_FUNCTION auto Involute::sign() const -> Sign +CELER_FUNCTION auto Involute::sign() const -> Chirality { - return r_b_ > 0 ? Sign::counterclockwise : Sign::clockwise; + return r_b_ > 0 ? Chirality::left : Chirality::right; } //---------------------------------------------------------------------------// @@ -200,7 +200,7 @@ CELER_FUNCTION SignedSense Involute::calc_sense(Real3 const& pos) const xy[0] = pos[0] - origin_[0]; xy[1] = pos[1] - origin_[1]; - if (this->sign() == Sign::clockwise) + if (this->sign() == Chirality::right) { xy[0] = negate(xy[0]); } @@ -300,7 +300,7 @@ CELER_FORCEINLINE_FUNCTION Real3 Involute::calc_normal(Real3 const& pos) const + a_; Real3 normal_ = {std::sin(angle), -std::cos(angle), 0}; - if (this->sign() == Sign::clockwise) + if (this->sign() == Chirality::right) { normal_[0] = negate(normal_[0]); } diff --git a/src/orange/surf/SoftSurfaceEqual.cc b/src/orange/surf/SoftSurfaceEqual.cc index cc0fec4b11..d2b3025622 100644 --- a/src/orange/surf/SoftSurfaceEqual.cc +++ b/src/orange/surf/SoftSurfaceEqual.cc @@ -199,6 +199,20 @@ bool SoftSurfaceEqual::operator()(GeneralQuadric const& a, && soft_eq_(a.zeroth(), b.zeroth()); } +//---------------------------------------------------------------------------// +/*! + * Compare two centered involutes for near equality. + */ +bool SoftSurfaceEqual::operator()(Involute const& a, Involute const& b) const +{ + return this->soft_eq_(a.r_b(), b.r_b()) + && this->soft_eq_(a.displacement_angle(), b.displacement_angle()) + && a.sign() == b.sign() && this->soft_eq_(a.tmin(), b.tmin()) + && this->soft_eq_(a.tmax(), b.tmax()) + && this->soft_eq_distance({a.origin()[0], a.origin()[1], 0}, + {b.origin()[0], b.origin()[1], 0}); +} + //---------------------------------------------------------------------------// // PRIVATE HELPER FUNCTIONS //---------------------------------------------------------------------------// diff --git a/src/orange/surf/SoftSurfaceEqual.hh b/src/orange/surf/SoftSurfaceEqual.hh index cc4abd8247..300e2951f9 100644 --- a/src/orange/surf/SoftSurfaceEqual.hh +++ b/src/orange/surf/SoftSurfaceEqual.hh @@ -73,6 +73,8 @@ class SoftSurfaceEqual bool operator()(GeneralQuadric const&, GeneralQuadric const&) const; + bool operator()(Involute const&, Involute const&) const; + private: SoftEqual<> soft_eq_; diff --git a/src/orange/surf/SurfaceClipper.cc b/src/orange/surf/SurfaceClipper.cc index 3727bcbd97..51fa9c2984 100644 --- a/src/orange/surf/SurfaceClipper.cc +++ b/src/orange/surf/SurfaceClipper.cc @@ -168,6 +168,16 @@ void SurfaceClipper::operator()(GeneralQuadric const&) const *int_ = BoundingBox{}; } +//---------------------------------------------------------------------------// +/*! + * Clip the bounding boxes to an involute. + */ +void SurfaceClipper::operator()(Involute const&) const +{ + // We no longer can guarantee any point being inside the shape; reset it + *int_ = BoundingBox{}; +} + //---------------------------------------------------------------------------// /*! * Clip a variant surface. diff --git a/src/orange/surf/SurfaceClipper.hh b/src/orange/surf/SurfaceClipper.hh index 745bebb88d..208ced60ff 100644 --- a/src/orange/surf/SurfaceClipper.hh +++ b/src/orange/surf/SurfaceClipper.hh @@ -70,6 +70,8 @@ class SurfaceClipper void operator()(GeneralQuadric const&) const; + void operator()(Involute const&) const; + // Apply to a surface with unknown type void operator()(VariantSurface const& surf) const; diff --git a/src/orange/surf/SurfaceFwd.hh b/src/orange/surf/SurfaceFwd.hh index 4b9327d01e..719b2058b1 100644 --- a/src/orange/surf/SurfaceFwd.hh +++ b/src/orange/surf/SurfaceFwd.hh @@ -20,13 +20,13 @@ class CylAligned; template class CylCentered; class GeneralQuadric; +class Involute; class Plane; template class PlaneAligned; class SimpleQuadric; class Sphere; class SphereCentered; -class Involute; //---------------------------------------------------------------------------// } // namespace celeritas diff --git a/src/orange/surf/SurfaceIO.cc b/src/orange/surf/SurfaceIO.cc index f3b56611ba..65d110fefb 100644 --- a/src/orange/surf/SurfaceIO.cc +++ b/src/orange/surf/SurfaceIO.cc @@ -11,6 +11,7 @@ #include #include "corecel/cont/ArrayIO.hh" +#include "orange/OrangeTypes.hh" #include "ConeAligned.hh" // IWYU pragma: associated #include "CylAligned.hh" // IWYU pragma: associated @@ -71,6 +72,22 @@ std::ostream& operator<<(std::ostream& os, GeneralQuadric const& s) return os; } +//---------------------------------------------------------------------------// +std::ostream& operator<<(std::ostream& os, Involute const& s) +{ + std::string sign{"ccw"}; + real_type a = s.displacement_angle(); + if (s.sign() == Chirality::right) + { + sign = "cw"; + a = constants::pi - a; + } + + return os << "Involute " << sign << ": r=" << s.r_b() << ", a=" << a + << ", t={" << s.tmin() << ',' << s.tmax() + << "} at x=" << s.origin()[0] << ", y=" << s.origin()[1]; +} + //---------------------------------------------------------------------------// std::ostream& operator<<(std::ostream& os, Plane const& s) { @@ -109,15 +126,6 @@ std::ostream& operator<<(std::ostream& os, SphereCentered const& s) return os; } -//---------------------------------------------------------------------------// -std::ostream& operator<<(std::ostream& os, Involute const& s) -{ - os << "Involute: r, a, sign, tmin, tmax =" << s.r_b() << ' ' << s.a() - << ' ' << s.sign() << ' ' << s.tmin() << ' ' << s.tmax() << " at " - << s.origin(); - return os; -} - //---------------------------------------------------------------------------// #undef ORANGE_INSTANTIATE_SHAPE_STREAM } // namespace celeritas diff --git a/src/orange/surf/SurfaceIO.hh b/src/orange/surf/SurfaceIO.hh index fa4428012b..3b89b492b3 100644 --- a/src/orange/surf/SurfaceIO.hh +++ b/src/orange/surf/SurfaceIO.hh @@ -29,6 +29,8 @@ std::ostream& operator<<(std::ostream&, CylCentered const&); std::ostream& operator<<(std::ostream&, GeneralQuadric const&); +std::ostream& operator<<(std::ostream&, Involute const&); + std::ostream& operator<<(std::ostream&, Plane const&); template @@ -39,8 +41,6 @@ std::ostream& operator<<(std::ostream&, SimpleQuadric const&); std::ostream& operator<<(std::ostream&, Sphere const&); std::ostream& operator<<(std::ostream&, SphereCentered const&); - -std::ostream& operator<<(std::ostream&, Involute const&); //!@} //---------------------------------------------------------------------------// } // namespace celeritas diff --git a/src/orange/surf/SurfaceTypeTraits.hh b/src/orange/surf/SurfaceTypeTraits.hh index 309270d688..9c44568be8 100644 --- a/src/orange/surf/SurfaceTypeTraits.hh +++ b/src/orange/surf/SurfaceTypeTraits.hh @@ -51,6 +51,7 @@ ORANGE_SURFACE_TRAITS(ky, ConeAligned); ORANGE_SURFACE_TRAITS(kz, ConeAligned); ORANGE_SURFACE_TRAITS(sq, SimpleQuadric); ORANGE_SURFACE_TRAITS(gq, GeneralQuadric); +ORANGE_SURFACE_TRAITS(inv, Involute); // clang-format on #undef ORANGE_SURFACE_TRAITS @@ -90,6 +91,10 @@ visit_surface_type(F&& func, SurfaceType st) ORANGE_ST_VISIT_CASE(kz); ORANGE_ST_VISIT_CASE(sq); ORANGE_ST_VISIT_CASE(gq); + case SurfaceType::inv: + // Prevented by input reader: see + // orange/detail/SurfacesRecordBuilder.cc + CELER_ASSERT_UNREACHABLE(); default: CELER_ASSERT_UNREACHABLE(); } diff --git a/src/orange/surf/detail/AllSurfaces.hh b/src/orange/surf/detail/AllSurfaces.hh index d651e96b32..458e4bd01e 100644 --- a/src/orange/surf/detail/AllSurfaces.hh +++ b/src/orange/surf/detail/AllSurfaces.hh @@ -12,6 +12,7 @@ #include "../CylAligned.hh" #include "../CylCentered.hh" #include "../GeneralQuadric.hh" +#include "../Involute.hh" #include "../Plane.hh" #include "../PlaneAligned.hh" #include "../SimpleQuadric.hh" diff --git a/src/orange/surf/detail/InvoluteSolver.hh b/src/orange/surf/detail/InvoluteSolver.hh index e08febf67f..052918db72 100644 --- a/src/orange/surf/detail/InvoluteSolver.hh +++ b/src/orange/surf/detail/InvoluteSolver.hh @@ -48,21 +48,18 @@ class InvoluteSolver using Intersections = Array; using SurfaceState = celeritas::SurfaceState; using Real2 = Array; - - //! Enum defining chirality of involute - enum Sign : bool - { - counterclockwise = 0, - clockwise, //!< Apply symmetry when solving - }; + //@} static inline CELER_FUNCTION real_type line_angle_param(real_type u, real_type v); public: // Construct Involute from parameters - inline CELER_FUNCTION InvoluteSolver( - real_type r_b, real_type a, Sign sign, real_type tmin, real_type tmax); + inline CELER_FUNCTION InvoluteSolver(real_type r_b, + real_type a, + Chirality sign, + real_type tmin, + real_type tmax); // Solve fully general case inline CELER_FUNCTION Intersections operator()( @@ -86,7 +83,7 @@ class InvoluteSolver //! Get involute parameters CELER_FUNCTION real_type r_b() const { return r_b_; } CELER_FUNCTION real_type a() const { return a_; } - CELER_FUNCTION Sign sign() const { return sign_; } + CELER_FUNCTION Chirality sign() const { return sign_; } //! Get bounds of the involute CELER_FUNCTION real_type tmin() const { return tmin_; } @@ -97,7 +94,7 @@ class InvoluteSolver // Involute parameters real_type r_b_; real_type a_; - Sign sign_; + Chirality sign_; // Bounds real_type tmin_; @@ -111,7 +108,7 @@ class InvoluteSolver * Construct from involute parameters. */ CELER_FUNCTION InvoluteSolver::InvoluteSolver( - real_type r_b, real_type a, Sign sign, real_type tmin, real_type tmax) + real_type r_b, real_type a, Chirality sign, real_type tmin, real_type tmax) : r_b_(r_b), a_(a), sign_(sign), tmin_(tmin), tmax_(tmax) { CELER_EXPECT(r_b > 0); @@ -122,14 +119,15 @@ CELER_FUNCTION InvoluteSolver::InvoluteSolver( } //---------------------------------------------------------------------------// - /*! * Find all roots for involute surfaces that are within the bounds and result - * in positive distances. Performed by doing a Regular Falsi Iteration on the + * in positive distances. + * + * Performed by doing a Regula Falsi Iteration on the * root function, \f[ * f(t) = r_b * [v{cos(a+t) + tsin(a+t)} + u{sin(a+t) - tcos(a+t)}] + xv - yu * \f] - * where the Regular Falsi Iteration is given by: \f[ + * where the Regula Falsi Iteration is given by: \f[ * t_c = [t_a*f(t_b) - t_b*f(t_a)] / [f(t_b) - f(t_a)] * \f] * where \em t_c replaces the bound with the same sign (e.g. \em t_a and \em @@ -152,7 +150,7 @@ InvoluteSolver::operator()(Real3 const& pos, real_type v = dir[1]; // Mirror systemm for counterclockwise involutes - if (sign_) + if (sign_ == Chirality::right) { x = -x; u = -u; @@ -221,7 +219,7 @@ InvoluteSolver::operator()(Real3 const& pos, // Only iterate when roots have different signs if (signum(ft_lower) != signum(ft_upper)) { - // Regular Falsi Iteration: Sometimes will slowly converge + // Regula Falsi Iteration real_type t_gamma = find_root_between(t_lower, t_upper); // Convert root to distance and store if positive and in interval diff --git a/src/orange/surf/detail/SurfaceTransformer.cc b/src/orange/surf/detail/SurfaceTransformer.cc index f0222b8364..4ee8ca9f2c 100644 --- a/src/orange/surf/detail/SurfaceTransformer.cc +++ b/src/orange/surf/detail/SurfaceTransformer.cc @@ -16,6 +16,7 @@ #include "../CylAligned.hh" #include "../CylCentered.hh" #include "../GeneralQuadric.hh" +#include "../Involute.hh" #include "../Plane.hh" #include "../PlaneAligned.hh" #include "../SimpleQuadric.hh" @@ -216,6 +217,15 @@ GeneralQuadric SurfaceTransformer::operator()(GeneralQuadric const& other) const qprime[0][0]); } +//---------------------------------------------------------------------------// +/*! + * Transform an Involute. + */ +Involute SurfaceTransformer::operator()(Involute const&) const +{ + CELER_NOT_IMPLEMENTED("transformed involutes"); +} + //---------------------------------------------------------------------------// } // namespace detail } // namespace celeritas diff --git a/src/orange/surf/detail/SurfaceTransformer.hh b/src/orange/surf/detail/SurfaceTransformer.hh index 9b62a899a7..97b7323cda 100644 --- a/src/orange/surf/detail/SurfaceTransformer.hh +++ b/src/orange/surf/detail/SurfaceTransformer.hh @@ -53,6 +53,8 @@ class SurfaceTransformer GeneralQuadric operator()(GeneralQuadric const&) const; + Involute operator()(Involute const&) const; + private: Transformation tr_; }; diff --git a/src/orange/surf/detail/SurfaceTranslator.cc b/src/orange/surf/detail/SurfaceTranslator.cc index abf0b94ff7..d5e2f78aed 100644 --- a/src/orange/surf/detail/SurfaceTranslator.cc +++ b/src/orange/surf/detail/SurfaceTranslator.cc @@ -7,6 +7,7 @@ //---------------------------------------------------------------------------// #include "SurfaceTranslator.hh" +#include "corecel/Constants.hh" #include "corecel/math/ArrayOperators.hh" #include "corecel/math/ArrayUtils.hh" #include "orange/MatrixUtils.hh" @@ -178,6 +179,24 @@ GeneralQuadric SurfaceTranslator::operator()(GeneralQuadric const& other) const second, make_array(other.cross()), real_type(2) * newfirst, newzeroth}; } +//---------------------------------------------------------------------------// +/*! + * Construct a translated Involute. + */ +Involute SurfaceTranslator::operator()(Involute const& other) const +{ + using constants::pi; + Real3 origin = tr_.transform_up({other.origin()[0], other.origin()[1], 0}); + + Involute invo{{origin[0], origin[1]}, + other.r_b(), + other.displacement_angle(), + other.sign(), + other.tmin(), + other.tmax()}; + return invo; +} + //---------------------------------------------------------------------------// } // namespace detail } // namespace celeritas diff --git a/src/orange/surf/detail/SurfaceTranslator.hh b/src/orange/surf/detail/SurfaceTranslator.hh index 6348bbfa45..625acdf26e 100644 --- a/src/orange/surf/detail/SurfaceTranslator.hh +++ b/src/orange/surf/detail/SurfaceTranslator.hh @@ -54,6 +54,8 @@ class SurfaceTranslator GeneralQuadric operator()(GeneralQuadric const&) const; + Involute operator()(Involute const&) const; + private: Translation tr_; }; diff --git a/test/orange/OrangeJson.test.cc b/test/orange/OrangeJson.test.cc index 32eea8e111..12e69405a0 100644 --- a/test/orange/OrangeJson.test.cc +++ b/test/orange/OrangeJson.test.cc @@ -5,6 +5,7 @@ //---------------------------------------------------------------------------// //! \file orange/OrangeJson.test.cc //---------------------------------------------------------------------------// +#include #include #include @@ -847,6 +848,130 @@ TEST_F(InputBuilderTest, DISABLED_universe_union_boundary) } } +//---------------------------------------------------------------------------// +// TODO: see celeritas-project/celeritas#1342 +TEST_F(InputBuilderTest, DISABLED_involute) +{ + { + SCOPED_TRACE("involute"); + auto result = this->track({0, 2.1, 0}, {1, 0, 0}); + static char const* const expected_volumes[] + = {"channel", "blade", "rest", "shell"}; + EXPECT_VEC_EQ(expected_volumes, result.volumes); + // Float and double produce different results + if constexpr (CELERITAS_REAL_TYPE == CELERITAS_REAL_TYPE_DOUBLE) + { + static real_type const expected_distances[] = { + 0.531630657850489, + 0.660884355302089, + 2.21189389295726, + 1.13321161570553, + }; + EXPECT_VEC_SOFT_EQ(expected_distances, result.distances); + } + else + { + static real_type const expected_distances[] + = {0.531625, 0.66089, 2.21189389295726, 1.13321161570553}; + EXPECT_VEC_SOFT_EQ(expected_distances, result.distances); + } + } + + { + SCOPED_TRACE("involute"); + auto result = this->track({0, 3.5, 0}, {0, -1, 0}); + static char const* const expected_volumes[] + = {"rest", "blade", "channel", "center", "rest", "shell"}; + EXPECT_VEC_EQ(expected_volumes, result.volumes); + static real_type const expected_distances[] = { + 0.528306129329604, 0.530097149547556, 0.441596721122841, 4, 2, 1}; + EXPECT_VEC_SOFT_EQ(expected_distances, result.distances); + } +} + +TEST_F(InputBuilderTest, DISABLED_involute_cw) +{ + { + SCOPED_TRACE("involute"); + auto result = this->track({-2, -1.5, 0}, {1, 0, 0}); + static char const* const expected_volumes[] + = {"rest", "center", "rest", "blade", "rest", "shell"}; + EXPECT_VEC_EQ(expected_volumes, result.volumes); + static real_type const expected_distances[] = { + 0.6771243444677, + 2.6457513110646, + 0.17135844958615, + 0.50955324797963, + 1.7043118904498, + 1.0615967635369, + }; + EXPECT_VEC_SOFT_EQ(expected_distances, result.distances); + } +} + +TEST_F(InputBuilderTest, DISABLED_involute_fuel) +{ + { + SCOPED_TRACE("involute"); + auto result = this->track({1.75, 3.5, 0}, {0, -1, 0}); + static char const* const expected_volumes[] = { + "clad2", + "fuel2", + "clad2", + "rest2", + "middle", + "rest1", + "clad1", + "fuel1", + "clad1", + "rest1", + "middle", + "rest2", + "shell", + }; + // Float and double produce different results + if constexpr (CELERITAS_REAL_TYPE == CELERITAS_REAL_TYPE_DOUBLE) + { + static real_type const expected_distances[] = { + 0.12694500489541, + 0.16591646127344, + 0.46300647647806, + 0.30743347115084, + 0.65134147906653, + 2.1115149489926, + 0.30217303181663, + 0.70489828892381, + 0.39023902667286, + 0.06188891786551, + 0.65134147906653, + 1.1601750562823, + 1.0868748563143, + }; + EXPECT_VEC_SOFT_EQ(expected_distances, result.distances); + } + else + { + static real_type const expected_distances[] = { + 0.12694500489541, + 0.16591646127344, + 0.46300647647806, + 0.30743347115084, + 0.65134147906653, + 2.1115149489926, + 0.30217303181663, + 0.70489828892381, + 0.390229, + 0.0618988, + 0.65134147906653, + 1.1601750562823, + 1.0868748563143, + }; + EXPECT_VEC_SOFT_EQ(expected_distances, result.distances); + } + EXPECT_VEC_EQ(expected_volumes, result.volumes); + } +} + //---------------------------------------------------------------------------// } // namespace test } // namespace celeritas diff --git a/test/orange/data/inputbuilder-involute-cw.org.json b/test/orange/data/inputbuilder-involute-cw.org.json new file mode 100644 index 0000000000..4545fe8c2c --- /dev/null +++ b/test/orange/data/inputbuilder-involute-cw.org.json @@ -0,0 +1,183 @@ +{ +"_format": "ORANGE", +"_version": 0, +"tol": { +"abs": 1e-05, +"rel": 1e-05 +}, +"universes": [ +{ +"_type": "unit", +"bbox": [ +[ +-5.0, +-5.0, +-1.0 +], +[ +5.0, +5.0, +1.0 +] +], +"md": { +"name": "involute" +}, +"surface_labels": [ +"blade@mz", +"blade@pz", +"bound@cz", +"blade@cz", +"blade@cz", +"blade@invl", +"blade@invr" +], +"surfaces": { +"data": [ +-1.0, +1.0, +25.0, +4.0, +16.0, +0.0, +0.0, +-1.0, +3.141592653589793, +1.7320508075688772, +4.36517666724533, +0.0, +0.0, +-1.0, +2.64939933255188, +1.7320508075688772, +4.36517666724533 +], +"sizes": [ +1, +1, +1, +1, +1, +6, +6 +], +"types": [ +"pz", +"pz", +"czc", +"czc", +"czc", +"inv", +"inv" +] +}, +"volume_labels": [ +"[EXTERIOR]@involute", +"center", +"blade", +"rest", +"shell" +], +"volumes": [ +{ +"faces": [ +0, +1, +2 +], +"flags": 1, +"logic": "0 1 ~ & 2 ~ & ~" +}, +{ +"bbox": [ +[ +-2.0, +-2.0, +-1.0 +], +[ +2.0, +2.0, +1.0 +] +], +"faces": [ +0, +1, +3 +], +"logic": "0 1 ~ & 2 ~ &" +}, +{ +"bbox": [ +[ +-4.0, +-4.0, +-1.0 +], +[ +4.0, +4.0, +1.0 +] +], +"faces": [ +0, +1, +3, +4, +5, +6 +], +"logic": "0 1 ~ & 2 & 3 ~ & 4 & 5 ~ &" +}, +{ +"bbox": [ +[ +-4.0, +-4.0, +-1.0 +], +[ +4.0, +4.0, +1.0 +] +], +"faces": [ +0, +1, +3, +4, +5, +6 +], +"flags": 1, +"logic": "0 1 ~ & 3 ~ & 0 1 ~ & 2 ~ & ~ & 0 1 ~ & 2 & 3 ~ & 4 & 5 ~ & ~ &" +}, +{ +"bbox": [ +[ +-5.0, +-5.0, +-1.0 +], +[ +5.0, +5.0, +1.0 +] +], +"faces": [ +0, +1, +2, +4 +], +"flags": 1, +"logic": "0 1 ~ & 2 ~ & 0 1 ~ & 3 ~ & ~ &" +} +] +} +] +} \ No newline at end of file diff --git a/test/orange/data/inputbuilder-involute-fuel.org.json b/test/orange/data/inputbuilder-involute-fuel.org.json new file mode 100644 index 0000000000..dfbc1777e2 --- /dev/null +++ b/test/orange/data/inputbuilder-involute-fuel.org.json @@ -0,0 +1,391 @@ +{ +"_format": "ORANGE", +"_version": 0, +"tol": { +"abs": 1e-05, +"rel": 1e-05 +}, +"universes": [ +{ +"_type": "unit", +"bbox": [ +[ +-5.0, +-5.0, +-1.0 +], +[ +5.0, +5.0, +1.0 +] +], +"md": { +"name": "involute" +}, +"surface_labels": [ +"blade1@mz", +"blade1@pz", +"bound@cz", +"blade1@cz", +"fuel1@cz", +"fuel1@cz", +"fuel1@invl", +"fuel1@invr", +"blade1@cz", +"blade1@invl", +"blade1@invr", +"blade2@cz", +"fuel2@cz", +"fuel2@cz", +"fuel2@invl", +"fuel2@invr", +"blade2@cz", +"blade2@invl", +"blade2@invr" +], +"surfaces": { +"data": [ +-1.0, +1.0, +25.0, +2.25, +3.24, +4.840000000000001, +0.0, +0.0, +-1.0, +3.0473448739820994, +1.4966629547095767, +2.0852555003701343, +0.0, +0.0, +-1.0, +2.9216811678385075, +1.4966629547095767, +2.0852555003701343, +6.25, +0.0, +0.0, +-1.0, +3.141592653589793, +1.118033988749895, +2.6054471128368992, +0.0, +0.0, +-1.0, +2.827433388230814, +1.118033988749895, +2.6054471128368992, +9.0, +10.240000000000002, +14.44, +0.0, +0.0, +2.0, +0.4084070449666731, +1.2489995996796799, +1.741213148283943, +0.0, +0.0, +2.0, +0.5340707511102649, +1.2489995996796799, +1.741213148283943, +16.0, +0.0, +0.0, +2.0, +0.3141592653589793, +1.118033988749895, +2.0462100729278565, +0.0, +0.0, +2.0, +0.6283185307179586, +1.118033988749895, +2.0462100729278565 +], +"sizes": [ +1, +1, +1, +1, +1, +1, +6, +6, +1, +6, +6, +1, +1, +1, +6, +6, +1, +6, +6 +], +"types": [ +"pz", +"pz", +"czc", +"czc", +"czc", +"czc", +"inv", +"inv", +"czc", +"inv", +"inv", +"czc", +"czc", +"czc", +"inv", +"inv", +"czc", +"inv", +"inv" +] +}, +"volume_labels": [ +"[EXTERIOR]@involute", +"center", +"fuel1", +"clad1", +"rest1", +"middle", +"fuel2", +"clad2", +"rest2", +"shell" +], +"volumes": [ +{ +"faces": [ +0, +1, +2 +], +"flags": 1, +"logic": "0 1 ~ & 2 ~ & ~" +}, +{ +"bbox": [ +[ +-1.5, +-1.5, +-1.0 +], +[ +1.5, +1.5, +1.0 +] +], +"faces": [ +0, +1, +3 +], +"logic": "0 1 ~ & 2 ~ &" +}, +{ +"bbox": [ +[ +-2.2, +-2.2, +-1.0 +], +[ +2.2, +2.2, +1.0 +] +], +"faces": [ +0, +1, +4, +5, +6, +7 +], +"logic": "0 1 ~ & 2 & 3 ~ & 4 & 5 ~ &" +}, +{ +"bbox": [ +[ +-2.5, +-2.5, +-1.0 +], +[ +2.5, +2.5, +1.0 +] +], +"faces": [ +0, +1, +3, +4, +5, +6, +7, +8, +9, +10 +], +"flags": 1, +"logic": "0 1 ~ & 2 & 7 ~ & 8 & 9 ~ & 0 1 ~ & 3 & 4 ~ & 5 & 6 ~ & ~ &" +}, +{ +"bbox": [ +[ +-2.5, +-2.5, +-1.0 +], +[ +2.5, +2.5, +1.0 +] +], +"faces": [ +0, +1, +3, +8, +9, +10 +], +"flags": 1, +"logic": "0 1 ~ & 3 ~ & 0 1 ~ & 2 & 3 ~ & 4 & 5 ~ & ~ & 0 1 ~ & 2 ~ & ~ &" +}, +{ +"bbox": [ +[ +-3.0, +-3.0, +-1.0 +], +[ +3.0, +3.0, +1.0 +] +], +"faces": [ +0, +1, +8, +11 +], +"flags": 1, +"logic": "0 1 ~ & 3 ~ & 0 1 ~ & 2 ~ & ~ &" +}, +{ +"bbox": [ +[ +-3.8, +-3.8, +-1.0 +], +[ +3.8, +3.8, +1.0 +] +], +"faces": [ +0, +1, +12, +13, +14, +15 +], +"logic": "0 1 ~ & 2 & 3 ~ & 4 ~ & 5 &" +}, +{ +"bbox": [ +[ +-4.0, +-4.0, +-1.0 +], +[ +4.0, +4.0, +1.0 +] +], +"faces": [ +0, +1, +11, +12, +13, +14, +15, +16, +17, +18 +], +"flags": 1, +"logic": "0 1 ~ & 2 & 7 ~ & 8 ~ & 9 & 0 1 ~ & 3 & 4 ~ & 5 ~ & 6 & ~ &" +}, +{ +"bbox": [ +[ +-4.0, +-4.0, +-1.0 +], +[ +4.0, +4.0, +1.0 +] +], +"faces": [ +0, +1, +11, +16, +17, +18 +], +"flags": 1, +"logic": "0 1 ~ & 3 ~ & 0 1 ~ & 2 & 3 ~ & 4 ~ & 5 & ~ & 0 1 ~ & 2 ~ & ~ &" +}, +{ +"bbox": [ +[ +-5.0, +-5.0, +-1.0 +], +[ +5.0, +5.0, +1.0 +] +], +"faces": [ +0, +1, +2, +16 +], +"flags": 1, +"logic": "0 1 ~ & 2 ~ & 0 1 ~ & 3 ~ & ~ &" +} +] +} +] +} \ No newline at end of file diff --git a/test/orange/data/inputbuilder-involute.org.json b/test/orange/data/inputbuilder-involute.org.json new file mode 100644 index 0000000000..41ddd61556 --- /dev/null +++ b/test/orange/data/inputbuilder-involute.org.json @@ -0,0 +1,217 @@ +{ +"_format": "ORANGE", +"_version": 0, +"tol": { +"abs": 1e-05, +"rel": 1e-05 +}, +"universes": [ +{ +"_type": "unit", +"bbox": [ +[ +-5.0, +-5.0, +-1.0 +], +[ +5.0, +5.0, +1.0 +] +], +"md": { +"name": "involute" +}, +"surface_labels": [ +"blade@mz", +"blade@pz", +"bound@cz", +"blade@cz", +"blade@cz", +"blade@invl", +"blade@invr", +"channel@invr" +], +"surfaces": { +"data": [ +-1.0, +1.0, +25.0, +4.0, +16.0, +0.0, +0.0, +1.0, +0.0, +1.7320508075688772, +4.36517666724533, +0.0, +0.0, +1.0, +0.4921933210379129, +1.7320508075688772, +4.36517666724533, +0.0, +0.0, +1.0, +0.9843866420758258, +1.7320508075688772, +4.36517666724533 +], +"sizes": [ +1, +1, +1, +1, +1, +6, +6, +6 +], +"types": [ +"pz", +"pz", +"czc", +"czc", +"czc", +"inv", +"inv", +"inv" +] +}, +"volume_labels": [ +"[EXTERIOR]@involute", +"center", +"blade", +"channel", +"rest", +"shell" +], +"volumes": [ +{ +"faces": [ +0, +1, +2 +], +"flags": 1, +"logic": "0 1 ~ & 2 ~ & ~" +}, +{ +"bbox": [ +[ +-2.0, +-2.0, +-1.0 +], +[ +2.0, +2.0, +1.0 +] +], +"faces": [ +0, +1, +3 +], +"logic": "0 1 ~ & 2 ~ &" +}, +{ +"bbox": [ +[ +-4.0, +-4.0, +-1.0 +], +[ +4.0, +4.0, +1.0 +] +], +"faces": [ +0, +1, +3, +4, +5, +6 +], +"logic": "0 1 ~ & 2 & 3 ~ & 4 ~ & 5 &" +}, +{ +"bbox": [ +[ +-4.0, +-4.0, +-1.0 +], +[ +4.0, +4.0, +1.0 +] +], +"faces": [ +0, +1, +3, +4, +6, +7 +], +"logic": "0 1 ~ & 2 & 3 ~ & 4 ~ & 5 &" +}, +{ +"bbox": [ +[ +-4.0, +-4.0, +-1.0 +], +[ +4.0, +4.0, +1.0 +] +], +"faces": [ +0, +1, +3, +4, +5, +6, +7 +], +"flags": 1, +"logic": "0 1 ~ & 3 ~ & 0 1 ~ & 2 ~ & ~ & 0 1 ~ & 2 & 3 ~ & 4 ~ & 5 & ~ & 0 1 ~ & 2 & 3 ~ & 5 ~ & 6 & ~ &" +}, +{ +"bbox": [ +[ +-5.0, +-5.0, +-1.0 +], +[ +5.0, +5.0, +1.0 +] +], +"faces": [ +0, +1, +2, +4 +], +"flags": 1, +"logic": "0 1 ~ & 2 ~ & 0 1 ~ & 3 ~ & ~ &" +} +] +} +] +} \ No newline at end of file diff --git a/test/orange/orangeinp/IntersectRegion.test.cc b/test/orange/orangeinp/IntersectRegion.test.cc index 7c69e44e28..605eef7a72 100644 --- a/test/orange/orangeinp/IntersectRegion.test.cc +++ b/test/orange/orangeinp/IntersectRegion.test.cc @@ -9,6 +9,7 @@ #include "orange/BoundingBoxUtils.hh" #include "orange/MatrixUtils.hh" +#include "orange/OrangeTypes.hh" #include "orange/orangeinp/CsgTreeUtils.hh" #include "orange/orangeinp/IntersectSurfaceBuilder.hh" #include "orange/orangeinp/detail/CsgUnitBuilder.hh" @@ -1301,6 +1302,203 @@ TEST_F(InfWedgeTest, half_turn) } } +//---------------------------------------------------------------------------// +// Involute +//---------------------------------------------------------------------------// +class InvoluteTest : public IntersectRegionTest +{ + public: + inline static constexpr auto ccw = Chirality::left; + inline static constexpr auto cw = Chirality::right; +}; + +TEST_F(InvoluteTest, single) +{ + { + // involute + auto result = this->test( + "invo", + Involute({1.0, 2.0, 4.0}, {0, 0.15667 * constants::pi}, cw, 1.0)); + + static char const expected_node[] = "all(+0, -1, +2, -3, +4, -5)"; + + EXPECT_EQ(expected_node, result.node); + EXPECT_VEC_SOFT_EQ((Real3{-4, -4, -1}), result.exterior.lower()); + EXPECT_VEC_SOFT_EQ((Real3{4, 4, 1}), result.exterior.upper()); + } + + static char const* const expected_surfaces[] = { + "Plane: z=-1", + "Plane: z=1", + "Cyl z: r=2", + "Cyl z: r=4", + "Involute cw: r=1, a=0, t={1.7321,4.3652} at x=0, y=0", + "Involute cw: r=1, a=0.49219, t={1.7321,4.3652} at x=0, y=0", + }; + EXPECT_VEC_EQ(expected_surfaces, surface_strings(this->unit())); + + auto node_strings = md_strings(this->unit()); + static char const* const expected_node_strings[] = { + "", + "", + "invo@mz", + "invo@pz", + "", + "invo@cz", + "invo@cz", + "", + "invo@invl", + "invo@invr", + "", + "", + }; + EXPECT_VEC_EQ(expected_node_strings, node_strings); +} + +// Counterclockwise adjacent involutes +TEST_F(InvoluteTest, two_ccw) +{ + { + // involute + auto result = this->test( + "top", + Involute({1.0, 2.0, 4.0}, {0, 0.15667 * constants::pi}, ccw, 1.0)); + + static char const expected_node[] = "all(+0, -1, +2, -3, -4, +5)"; + + EXPECT_EQ(expected_node, result.node); + EXPECT_VEC_SOFT_EQ((Real3{-4, -4, -1}), result.exterior.lower()); + EXPECT_VEC_SOFT_EQ((Real3{4, 4, 1}), result.exterior.upper()); + } + { + // bottom + auto result = this->test( + "bottom", + Involute({1.0, 2.0, 4.0}, + {0.15667 * constants::pi, 0.31334 * constants::pi}, + ccw, + 1.0)); + + // Float and double produce different results + if constexpr (CELERITAS_REAL_TYPE == CELERITAS_REAL_TYPE_DOUBLE) + { + static char const expected_node[] = "all(+0, -1, +2, -3, -5, +6)"; + EXPECT_EQ(expected_node, result.node); + } + else + { + static char const expected_node[] = "all(+0, -1, +2, -3, -5, +7)"; + EXPECT_EQ(expected_node, result.node); + } + + EXPECT_VEC_SOFT_EQ((Real3{-4, -4, -1}), result.exterior.lower()); + EXPECT_VEC_SOFT_EQ((Real3{4, 4, 1}), result.exterior.upper()); + } + + static char const* const expected_surfaces[] = { + "Plane: z=-1", + "Plane: z=1", + "Cyl z: r=2", + "Cyl z: r=4", + "Involute ccw: r=1, a=0, t={1.7321,4.3652} at x=0, y=0", + "Involute ccw: r=1, a=0.49219, t={1.7321,4.3652} at x=0, y=0", + "Involute ccw: r=1, a=0.98439, t={1.7321,4.3652} at x=0, y=0", + }; + EXPECT_VEC_EQ(expected_surfaces, surface_strings(this->unit())); + + auto node_strings = md_strings(this->unit()); + static char const* const expected_node_strings[] = { + "", + "", + "bottom@mz,top@mz", + "bottom@pz,top@pz", + "", + "bottom@cz,top@cz", + "bottom@cz,top@cz", + "", + "top@invl", + "", + "bottom@invl,top@invr", + "", + "", + "bottom@invr", + "", + }; + EXPECT_VEC_EQ(expected_node_strings, node_strings); +} + +// Clockwise varient of previous +TEST_F(InvoluteTest, two_cw) +{ + { + // involute + auto result = this->test( + "top", + Involute({1.0, 2.0, 4.0}, {0, 0.15667 * constants::pi}, cw, 1.0)); + + static char const expected_node[] = "all(+0, -1, +2, -3, +4, -5)"; + + EXPECT_EQ(expected_node, result.node); + EXPECT_VEC_SOFT_EQ((Real3{-4, -4, -1}), result.exterior.lower()); + EXPECT_VEC_SOFT_EQ((Real3{4, 4, 1}), result.exterior.upper()); + } + { + // bottom + auto result = this->test( + "bottom", + Involute({1.0, 2.0, 4.0}, + {0.15667 * constants::pi, 0.31334 * constants::pi}, + cw, + 1.0)); + + // Float and double produce different results + if constexpr (CELERITAS_REAL_TYPE == CELERITAS_REAL_TYPE_DOUBLE) + { + static char const expected_node[] = "all(+0, -1, +2, -3, +5, -6)"; + EXPECT_EQ(expected_node, result.node); + } + else + { + static char const expected_node[] = "all(+0, -1, +2, -3, +5, -7)"; + EXPECT_EQ(expected_node, result.node); + } + + EXPECT_VEC_SOFT_EQ((Real3{-4, -4, -1}), result.exterior.lower()); + EXPECT_VEC_SOFT_EQ((Real3{4, 4, 1}), result.exterior.upper()); + } + + static char const* const expected_surfaces[] = { + "Plane: z=-1", + "Plane: z=1", + "Cyl z: r=2", + "Cyl z: r=4", + "Involute cw: r=1, a=0, t={1.7321,4.3652} at x=0, y=0", + "Involute cw: r=1, a=0.49219, t={1.7321,4.3652} at x=0, y=0", + "Involute cw: r=1, a=0.98439, t={1.7321,4.3652} at x=0, y=0", + }; + EXPECT_VEC_EQ(expected_surfaces, surface_strings(this->unit())); + + auto node_strings = md_strings(this->unit()); + static char const* const expected_node_strings[] = { + "", + "", + "bottom@mz,top@mz", + "bottom@pz,top@pz", + "", + "bottom@cz,top@cz", + "bottom@cz,top@cz", + "", + "top@invl", + "bottom@invl,top@invr", + "", + "", + "bottom@invr", + "", + "", + }; + EXPECT_VEC_EQ(expected_node_strings, node_strings); +} + //---------------------------------------------------------------------------// // PARALLELEPIPED //---------------------------------------------------------------------------// diff --git a/test/orange/orangeinp/UnitProto.test.cc b/test/orange/orangeinp/UnitProto.test.cc index 4e9ffd7c63..fa628741c1 100644 --- a/test/orange/orangeinp/UnitProto.test.cc +++ b/test/orange/orangeinp/UnitProto.test.cc @@ -16,6 +16,7 @@ #include "corecel/math/ArrayOperators.hh" #include "corecel/math/ArrayUtils.hh" #include "orange/OrangeInputIO.json.hh" +#include "orange/OrangeTypes.hh" #include "orange/orangeinp/CsgObject.hh" #include "orange/orangeinp/InputBuilder.hh" #include "orange/orangeinp/Shape.hh" @@ -36,6 +37,10 @@ namespace test //---------------------------------------------------------------------------// using SPConstObject = std::shared_ptr; using SPConstProto = std::shared_ptr; +//! Enum defining chirality +using Sign = celeritas::Chirality; +Sign ccw = celeritas::Chirality::left; +Sign cw = celeritas::Chirality::right; //---------------------------------------------------------------------------// // Construction helper functions @@ -78,6 +83,16 @@ SPConstObject make_box(std::string&& label, Real3 const& lo, Real3 const& hi) return result; } +SPConstObject make_inv(std::string&& label, + Real3 const& radii, + Real2 const& displacement, + Sign sign, + real_type halfheight) +{ + return make_shape( + std::move(label), radii, displacement, sign, halfheight); +} + SPConstProto make_daughter(std::string label) { UnitProto::Input inp; @@ -742,6 +757,157 @@ TEST_F(InputBuilderTest, universe_union_boundary) this->run_test(*outer); } +/*! + * Generate input for a universe with two involutes and two cylinders + */ +TEST_F(InputBuilderTest, involute) +{ + auto involute = std::make_shared([] { + auto invo1 = make_inv( + "blade", {1.0, 2.0, 4.0}, {0, 0.15667 * constants::pi}, ccw, 1.0); + auto invo2 + = make_inv("channel", + {1.0, 2.0, 4.0}, + {0.15667 * constants::pi, 0.31334 * constants::pi}, + ccw, + 1.0); + auto cyl = make_cyl("bound", 5.0, 1.0); + auto system = make_cyl("system", 4.0, 1.0); + auto inner = make_cyl("center", 2.0, 1.0); + UnitProto::Input inp; + inp.boundary.interior = cyl; + inp.boundary.zorder = ZOrder::media; + inp.label = "involute"; + + inp.materials.push_back(make_material(SPConstObject{inner}, 1)); + inp.materials.push_back(make_material(SPConstObject{invo1}, 2)); + inp.materials.push_back(make_material(SPConstObject{invo2}, 3)); + inp.materials.push_back( + make_material(make_rdv("rest", + {{Sense::inside, system}, + {Sense::outside, inner}, + {Sense::outside, invo1}, + {Sense::outside, invo2}}), + 5)); + inp.materials.push_back( + make_material(make_rdv("shell", + {{Sense::inside, inp.boundary.interior}, + {Sense::outside, system}}), + 5)); + + return inp; + }()); + + this->run_test(*involute); +} + +/*! + * Involute Blade + */ +TEST_F(InputBuilderTest, involute_cw) +{ + auto involute = std::make_shared([] { + auto invo1 = make_inv( + "blade", {1.0, 2.0, 4.0}, {0, 0.15667 * constants::pi}, cw, 1.0); + auto cyl = make_cyl("bound", 5.0, 1.0); + auto system = make_cyl("system", 4.0, 1.0); + auto inner = make_cyl("center", 2.0, 1.0); + UnitProto::Input inp; + inp.boundary.interior = cyl; + inp.boundary.zorder = ZOrder::media; + inp.label = "involute"; + + inp.materials.push_back(make_material(SPConstObject{inner}, 1)); + inp.materials.push_back(make_material(SPConstObject{invo1}, 2)); + inp.materials.push_back( + make_material(make_rdv("rest", + {{Sense::inside, system}, + {Sense::outside, inner}, + {Sense::outside, invo1}}), + 4)); + inp.materials.push_back( + make_material(make_rdv("shell", + {{Sense::inside, inp.boundary.interior}, + {Sense::outside, system}}), + 5)); + + return inp; + }()); + + this->run_test(*involute); +} + +/*! + * Clockwise and Counterclockwise fuel blade + */ +TEST_F(InputBuilderTest, involute_fuel) +{ + auto involute = std::make_shared([] { + auto inner1 = make_cyl("center", 1.5, 1.0); + auto cyl = make_cyl("bound", 5.0, 1.0); + auto invo1 = make_inv( + "blade1", {1.0, 1.5, 2.5}, {0, 0.1 * constants::pi}, cw, 1.0); + auto invo2 = make_inv("fuel1", + {1.0, 1.8, 2.2}, + {0.03 * constants::pi, 0.07 * constants::pi}, + cw, + 1.0); + auto outer1 = make_cyl("middle_1", 2.5, 1.0); + auto inner2 = make_cyl("middle_2", 3.0, 1.0); + auto invo3 = make_inv("blade2", + {2.0, 3.0, 4.0}, + {0.1 * constants::pi, 0.2 * constants::pi}, + ccw, + 1.0); + auto invo4 = make_inv("fuel2", + {2.0, 3.2, 3.8}, + {0.13 * constants::pi, 0.17 * constants::pi}, + ccw, + 1.0); + auto outer2 = make_cyl("outer", 4.0, 1.0); + + UnitProto::Input inp; + inp.boundary.interior = cyl; + inp.boundary.zorder = ZOrder::media; + inp.label = "involute"; + + inp.materials.push_back(make_material(SPConstObject{inner1}, 1)); + inp.materials.push_back(make_material(SPConstObject{invo2}, 2)); + inp.materials.push_back(make_material( + make_rdv("clad1", {{Sense::inside, invo1}, {Sense::outside, invo2}}), + 3)); + inp.materials.push_back( + make_material(make_rdv("rest1", + {{Sense::inside, outer1}, + {Sense::outside, invo1}, + {Sense::outside, inner1}}), + 4)); + inp.materials.push_back(make_material( + make_rdv("middle", + {{Sense::inside, inner2}, {Sense::outside, outer1}}), + 5)); + inp.materials.push_back(make_material(SPConstObject{invo4}, 6)); + inp.materials.push_back(make_material( + make_rdv("clad2", {{Sense::inside, invo3}, {Sense::outside, invo4}}), + 7)); + inp.materials.push_back( + make_material(make_rdv("rest2", + {{Sense::inside, outer2}, + {Sense::outside, invo3}, + {Sense::outside, inner2}}), + 8)); + inp.materials.push_back( + make_material(make_rdv("shell", + {{Sense::inside, inp.boundary.interior}, + {Sense::outside, outer2}}), + 9)); + + return inp; + }()); + + this->run_test(*involute); +} + //---------------------------------------------------------------------------// } // namespace test } // namespace orangeinp diff --git a/test/orange/surf/Involute.test.cc b/test/orange/surf/Involute.test.cc index 69566aa389..e42ea08d33 100644 --- a/test/orange/surf/Involute.test.cc +++ b/test/orange/surf/Involute.test.cc @@ -23,19 +23,18 @@ namespace test //---------------------------------------------------------------------------// using constants::pi; -using Sign = Involute::Sign; +using Sign = Chirality; +Sign ccw = Chirality::left; +Sign cw = Chirality::right; using Real2 = Involute::Real2; -constexpr auto ccw = Sign::counterclockwise; -constexpr auto cw = Sign::clockwise; - //---------------------------------------------------------------------------// TEST(InvoluteTest, construction) { auto check_props = [](Involute const& invo) { EXPECT_VEC_SOFT_EQ((Real2{1, 0}), invo.origin()); EXPECT_SOFT_EQ(2.0, invo.r_b()); - EXPECT_SOFT_EQ(pi - 0.2, invo.a()); + EXPECT_SOFT_EQ(pi - 0.2, invo.displacement_angle()); EXPECT_EQ(invo.sign(), cw); EXPECT_SOFT_EQ(1.0, invo.tmin()); EXPECT_SOFT_EQ(3.0, invo.tmax()); diff --git a/test/orange/surf/LocalSurfaceVisitor.test.cc b/test/orange/surf/LocalSurfaceVisitor.test.cc index f5881baa64..f9e5605e3f 100644 --- a/test/orange/surf/LocalSurfaceVisitor.test.cc +++ b/test/orange/surf/LocalSurfaceVisitor.test.cc @@ -173,6 +173,17 @@ TEST_F(SurfaceActionTest, surface_traits_visitor) // Check that all surface types can be visited and are consistent for (auto st : range(SurfaceType::size_)) { + if (st == SurfaceType::inv) + { + if (CELERITAS_DEBUG) + { + // TODO: see celeritas-project/celeritas#1342 + EXPECT_THROW(visit_surface_type(get_surface_type, st), + DebugError); + } + continue; + } + SurfaceType actual_st = visit_surface_type(get_surface_type, st); EXPECT_EQ(st, actual_st); } diff --git a/test/orange/surf/SoftSurfaceEqual.test.cc b/test/orange/surf/SoftSurfaceEqual.test.cc index 5fb7d314e3..30694d355d 100644 --- a/test/orange/surf/SoftSurfaceEqual.test.cc +++ b/test/orange/surf/SoftSurfaceEqual.test.cc @@ -192,6 +192,55 @@ TEST_F(SoftSurfaceEqualTest, general_quadric) softeq_(ref, SurfaceTranslator(Translation{{large, 0, 0}})(ref))); } +TEST_F(SoftSurfaceEqualTest, involute) +{ + using Sign = Chirality; + Sign ccw = Chirality::left; + Sign cw = Chirality::right; + + Involute ref_ccw{{1.0, 0.0}, 1.0, 2.0, ccw, 1.0, 2.0}; + Involute ref_cw{{1.0, 0.0}, 1.0, 2.0, cw, 1.0, 2.0}; + + // Counterclockwise + EXPECT_TRUE(softeq_( + ref_ccw, Involute{{1.0, 0.0}, 1.0 + small, 2.0, ccw, 1.0, 2.0})); + EXPECT_FALSE(softeq_( + ref_ccw, Involute{{1.0, 0.0}, 1.0 + large, 2.0, ccw, 1.0, 2.0})); + + EXPECT_TRUE(softeq_( + ref_ccw, Involute{{1.0, 0.0}, 1.0, 2.0 + small, ccw, 1.0, 2.0})); + EXPECT_FALSE(softeq_( + ref_ccw, Involute{{1.0, 0.0}, 1.0, 2.0 + large, ccw, 1.0, 2.0})); + + EXPECT_TRUE(softeq_( + ref_ccw, Involute{{1.0, 0.0}, 1.0, 2.0, ccw, 1.0 + small, 2.0})); + EXPECT_FALSE(softeq_( + ref_ccw, Involute{{1.0, 0.0}, 1.0, 2.0, ccw, 1.0 + large, 2.0})); + + EXPECT_TRUE(softeq_( + ref_ccw, Involute{{1.0, 0.0}, 1.0, 2.0, ccw, 1.0, 2.0 + small})); + EXPECT_FALSE(softeq_( + ref_ccw, Involute{{1.0, 0.0}, 1.0, 2.0, ccw, 1.0, 2.0 + large})); + + EXPECT_TRUE(softeq_( + ref_ccw, Involute{{1.0 + small, 0.0}, 1.0, 2.0, ccw, 1.0, 2.0})); + EXPECT_FALSE(softeq_( + ref_ccw, Involute{{1.0 + large, 0.0}, 1.0, 2.0, ccw, 1.0, 2.0})); + + EXPECT_FALSE(softeq_(ref_ccw, ref_cw)); + + // Clockwise + EXPECT_TRUE( + softeq_(ref_cw, Involute{{1.0, 0.0}, 1.0 + small, 2.0, cw, 1.0, 2.0})); + EXPECT_FALSE( + softeq_(ref_cw, Involute{{1.0, 0.0}, 1.0 + large, 2.0, cw, 1.0, 2.0})); + + EXPECT_TRUE(softeq_( + ref_cw, Involute{{1.0, 0.0}, 1.0, 2.0 + small, cw, 1.0 + small, 2.0})); + EXPECT_FALSE(softeq_( + ref_cw, Involute{{1.0, 0.0}, 1.0, 2.0 + large, cw, 1.0 + large, 2.0})); +} + //---------------------------------------------------------------------------// } // namespace test } // namespace celeritas diff --git a/test/orange/surf/detail/InvoluteSolver.test.cc b/test/orange/surf/detail/InvoluteSolver.test.cc index 5b984dd21c..e7ae9d9cc2 100644 --- a/test/orange/surf/detail/InvoluteSolver.test.cc +++ b/test/orange/surf/detail/InvoluteSolver.test.cc @@ -19,7 +19,9 @@ namespace test { //---------------------------------------------------------------------------// using constants::pi; -using Sign = InvoluteSolver::Sign; +using Sign = Chirality; +Sign ccw = Chirality::left; +Sign cw = Chirality::right; using SurfaceSate = celeritas::SurfaceState; TEST(SolveSurface, no_roots) @@ -30,7 +32,7 @@ TEST(SolveSurface, no_roots) { real_type r_b = 3.0; real_type a = 0; - auto sign = InvoluteSolver::counterclockwise; + auto sign = ccw; real_type x = 0; real_type y = -2; @@ -55,7 +57,7 @@ TEST(SolveSurface, no_roots) { real_type r_b = 0.75; real_type a = 0; - auto sign = InvoluteSolver::counterclockwise; + auto sign = ccw; real_type x = -7; real_type y = -1; @@ -80,7 +82,7 @@ TEST(SolveSurface, no_roots) { real_type r_b = 0.75; real_type a = 0; - auto sign = InvoluteSolver::counterclockwise; + auto sign = ccw; real_type x = -7; real_type y = -1; @@ -105,7 +107,7 @@ TEST(SolveSurface, no_roots) { real_type r_b = 1.1; real_type a = 0.5 * pi; - auto sign = InvoluteSolver::clockwise; + auto sign = cw; real_type x = -0.2; real_type y = 1.1; @@ -130,7 +132,7 @@ TEST(SolveSurface, no_roots) { real_type r_b = 3.0; real_type a = pi; - auto sign = InvoluteSolver::counterclockwise; + auto sign = ccw; real_type x = -4.101853006408607; real_type y = -5.443541628262038; @@ -157,7 +159,7 @@ TEST(SolveSurface, one_root) { real_type r_b = 1.0; real_type a = 0; - auto sign = InvoluteSolver::counterclockwise; + auto sign = ccw; real_type x = 0; real_type y = 0; @@ -182,7 +184,7 @@ TEST(SolveSurface, one_root) { real_type r_b = 1.5; real_type a = 0; - auto sign = InvoluteSolver::counterclockwise; + auto sign = ccw; real_type x = -1.5; real_type y = 1.0; @@ -207,7 +209,7 @@ TEST(SolveSurface, one_root) { real_type r_b = 0.5; real_type a = 0.6 * pi; - auto sign = InvoluteSolver::clockwise; + auto sign = cw; real_type x = -4.0; real_type y = 2.0; @@ -233,7 +235,7 @@ TEST(SolveSurface, one_root) { real_type r_b = 1.1; real_type a = 1.5 * pi; - auto sign = InvoluteSolver::counterclockwise; + auto sign = ccw; real_type x = 0.0058102462574510716; real_type y = -1.1342955336941216; @@ -260,7 +262,7 @@ TEST(SolveSurface, two_roots) { real_type r_b = 1.1; real_type a = 0.5 * pi; - auto sign = InvoluteSolver::clockwise; + auto sign = cw; real_type x = -0.2; real_type y = 1.1; @@ -294,7 +296,7 @@ TEST(SolveSurface, two_roots) { real_type r_b = 1.1; real_type a = 1.5 * pi; - auto sign = InvoluteSolver::clockwise; + auto sign = cw; real_type x = -0.0001; real_type y = -1.11; @@ -332,7 +334,7 @@ TEST(SolveSurface, three_roots) { real_type r_b = 1.1; real_type a = 1.5 * pi; - auto sign = InvoluteSolver::counterclockwise; + auto sign = ccw; real_type x = -6.8653259986571326; real_type y = -0.30468105643505367; @@ -357,7 +359,7 @@ TEST(SolveSurface, tangents) // Direction (0,1) real_type r_b = 1.0; real_type a = 0; - auto sign = InvoluteSolver::counterclockwise; + auto sign = ccw; real_type tmin = 0.33 * pi; real_type tmax = 0.67 * pi; InvoluteSolver solve(r_b, a, sign, tmin, tmax); @@ -633,7 +635,7 @@ TEST(Components, calc_dist) { real_type r_b = 1.1; real_type a = 0.5 * pi; - auto sign = InvoluteSolver::clockwise; + auto sign = cw; real_type x = -0.2; real_type y = 1.1; diff --git a/test/orange/surf/detail/SurfaceTranslator.test.cc b/test/orange/surf/detail/SurfaceTranslator.test.cc index d0e61dcc26..d97cb33460 100644 --- a/test/orange/surf/detail/SurfaceTranslator.test.cc +++ b/test/orange/surf/detail/SurfaceTranslator.test.cc @@ -120,6 +120,35 @@ TEST_F(SurfaceTranslatorTest, general_quadric) EXPECT_SOFT_EQ(7.0, distances[1]); } +TEST_F(SurfaceTranslatorTest, involute) +{ + using Real2 = Involute::Real2; + using Sign = Chirality; + Sign ccw = Chirality::left; + Sign cw = Chirality::right; + // See Involute.tst.cc + // Cocklwise involute + { + auto invo = translate(Involute{{1, 0}, 2.0, 0.2, cw, 1.0, 3.0}); + EXPECT_VEC_SOFT_EQ((Real2{3, 3}), invo.origin()); + EXPECT_SOFT_EQ(2.0, invo.r_b()); + EXPECT_SOFT_EQ(0.2, invo.displacement_angle()); + EXPECT_TRUE(cw == invo.sign()); + EXPECT_SOFT_EQ(1.0, invo.tmin()); + EXPECT_SOFT_EQ(3.0, invo.tmax()); + } + // Counterclockwise involute + { + auto invo = translate(Involute{{1, 0}, 2.0, 0.2, ccw, 1.0, 3.0}); + EXPECT_VEC_SOFT_EQ((Real2{3, 3}), invo.origin()); + EXPECT_SOFT_EQ(2.0, invo.r_b()); + EXPECT_SOFT_EQ(0.2, invo.displacement_angle()); + EXPECT_TRUE(ccw == invo.sign()); + EXPECT_SOFT_EQ(1.0, invo.tmin()); + EXPECT_SOFT_EQ(3.0, invo.tmax()); + } +} + //---------------------------------------------------------------------------// } // namespace test } // namespace detail