primal: Add method for checking if a sphere contains a point

[adayton1]

The current method is to check if the signed distance is less than (or equal to?) zero. But there's probably a faster way avoiding square roots.

[kennyweiss]

We have something very similar here in primal::in_sphere(BoundingBox<2>, Sphere<2>), which uses the squared distance to each of the (four) vertices. It should be easy to add an overload primal::in_sphere(Point, Sphere) in 2D and 3D.

See:

axom/src/axom/primal/operators/in_sphere.hpp

Lines 132 to 178 in 3077763

	/*!
	* \brief Tests whether a bounding box lies inside a 2D sphere
	*
	* \param [in] bb the bounding box
	* \param [in] circle the sphere
	*/
	template <typename T>
	inline bool in_sphere(const BoundingBox<T, 2>& bb, const Sphere<T, 2>& circle)
	{
	// Check if any corner of the bounding box is outside the sphere.
	// This version requires a lot of multiplications, but no square roots
	// and a higher likelihood of early returns.
	auto the_max = bb.getMax();
	auto the_min = bb.getMin();
	auto radius = circle.getRadius();
	auto center = circle.getCenter();
	if((center[0] - the_min[0]) * (center[0] - the_min[0]) +
	(center[1] - the_min[1]) * (center[1] - the_min[1]) >
	radius * radius)
	{
	return false;
	}
	if((center[0] - the_max[0]) * (center[0] - the_max[0]) +
	(center[1] - the_min[1]) * (center[1] - the_min[1]) >
	radius * radius)
	{
	return false;
	}
	if((center[0] - the_min[0]) * (center[0] - the_min[0]) +
	(center[1] - the_max[1]) * (center[1] - the_max[1]) >
	radius * radius)
	{
	return false;
	}
	if((center[0] - the_max[0]) * (center[0] - the_max[0]) +
	(center[1] - the_max[1]) * (center[1] - the_max[1]) >
	radius * radius)
	{
	return false;
	}
	return true;
	}

[kennyweiss]

We could also add something like primal::Sphere::computeSquaredDistance(Point) and use it in primal::Sphere::getOrientation(Point)

axom/src/axom/primal/geometry/Sphere.hpp

Lines 122 to 168 in 3077763

	/*!
	* \brief Computes the signed distance of a point to the Sphere's boundary.
	*
	* \param [in] q The test point
	* \return d the computed signed distance of the point \a q to the sphere.
	*
	* \note The signed distance of a point \a q is:
	* <ul>
	* <li> negative inside the sphere </li>
	* <li> positive outside the sphere </li>
	* <li> zero on the boundary </li>
	* </ul>
	*/
	AXOM_HOST_DEVICE inline T computeSignedDistance(const PointType& q) const
	{
	return VectorType(m_center, q).norm() - m_radius;
	}
	/*!
	* \brief Computes the orientation of a point with respect to the Sphere.
	*
	* \param [in] q The test point
	* \param [in] TOL user-supplied tolerance. Optional. Default is 1.e-9.
	* \return orient the orientation of \a q with respect to the sphere.
	*
	* \note This method returns one of the following values:
	* <ul>
	* <li> <b>ON_BOUNDARY</b> : if \a q is on the sphere's boundary </li>
	* <li> <b>ON_POSITIVE_SIDE</b> : if \a q is outside the sphere </li>
	* <li> <b>ON_NEGATIVE_SIDE</b> : if \a q is inside the sphere </li>
	* </ul>
	*
	* \see OrientationResult for the list of possible return values.
	*
	*/
	AXOM_HOST_DEVICE
	inline int getOrientation(const PointType& q, double TOL = 1.e-9) const
	{
	const T signed_distance = this->computeSignedDistance(q);
	if(axom::utilities::isNearlyEqual(signed_distance, 0., TOL))
	{
	return primal::ON_BOUNDARY;
	}
	return (signed_distance < T {0}) ? primal::ON_NEGATIVE_SIDE : primal::ON_POSITIVE_SIDE;
	}