SplineUtils Namespace Reference

Utilities useful for splines in general Reference for B-Spline implementation: https://mat.fsv.cvut.cz/gcg/sbornik/prochazkova.pdf. More...

Functions
std::vector< Point >	circularControlPoints (const libMesh::Point &start_point, const libMesh::Point &end_point, const libMesh::RealVectorValue &parallel_direction, const unsigned int num_cps)
	Creates control points for the special case of parallel lines using an interpolating circle. More...
std::vector< Point >	bSplineControlPoints (const libMesh::Point &start_point, const libMesh::Point &end_point, const libMesh::RealVectorValue &start_direction, const libMesh::RealVectorValue &end_direction, const unsigned int cps_per_half, const libMesh::Real sharpness)
	Creates control points for an open uniform BSpline. More...
std::vector< Point >	controlPointsAlongLine (const libMesh::Point &start_point, const libMesh::Point &end_point, const libMesh::RealVectorValue &direction_vector, const libMesh::Real sharpness, const unsigned int num_cps, const bool reverse_order=false)
	Creates control points along an extrapolated line up to a certain intercept. More...
libMesh::Point	makeControlPoint (const libMesh::Point &start_point, const libMesh::Point &end_point, const libMesh::RealVectorValue &direction_vector, const libMesh::Real sharpness)
	Creates a control point. More...
std::pair< Point, Point >	closestPoints (const libMesh::Point &point_1, const libMesh::Point &point_2, const libMesh::RealVectorValue &direction_1, const libMesh::RealVectorValue &direction_2)
	Determines the two points defining the shortest line segment between two lines in space. More...

Detailed Description

Utilities useful for splines in general Reference for B-Spline implementation: https://mat.fsv.cvut.cz/gcg/sbornik/prochazkova.pdf.

Function Documentation

bSplineControlPoints()

std::vector< Point > SplineUtils::bSplineControlPoints	(	const libMesh::Point &	start_point,
		const libMesh::Point &	end_point,
		const libMesh::RealVectorValue &	start_direction,
		const libMesh::RealVectorValue &	end_direction,
		const unsigned int	cps_per_half,
		const libMesh::Real	sharpness
	)

Creates control points for an open uniform BSpline.

Parameters

start_point	start point for spline
end_point	end point for spline
start_direction	incoming direction (derivative)
end_direction	outgoing direction (derivative)
cps_per_half	number of control points per half of the spline – will be along extrapolated line from the point and its direction
sharpness	number in [0,1] that determines the sharpness of the curve

Definition at line 67 of file SplineUtils.C.

Referenced by Moose::BSpline::createControlPoints(), and BSplineCurveGenerator::generate().

{
  // check that start_point is different from end_point
  if ((start_point - end_point).norm_sq() < libMesh::TOLERANCE * libMesh::TOLERANCE)
    mooseError("Start and end points must be different.");
  // check that neither direction is the zero vector
  if (start_direction.norm_sq() < libMesh::TOLERANCE * libMesh::TOLERANCE ||
      end_direction.norm_sq() < libMesh::TOLERANCE * libMesh::TOLERANCE)
    mooseError("Direction vectors must have a nonzero magnitude.");
  mooseAssert(sharpness <= 1.0 && sharpness >= 0, "Sharpness must be in [0,1]!");

  // check if directions are parallel, use a special case if so
  const bool parallel = ((start_direction.cross(end_direction)).norm() < libMesh::TOLERANCE);
  if (parallel)
  {
    // We should not revert to a circle for this case
    if (MooseUtils::absoluteFuzzyEqual(start_direction.unit(), -end_direction.unit()))
      mooseError("Start and end directions are opposite, this is not currently supported");
    mooseWarning("Directions are parallel! Attempting to use circular control points...");
    unsigned int num_cps = 2 * cps_per_half + 1;
    if (num_cps < 30)
      mooseWarning("Number of control points required for circular control points is much greater "
                   "than for BSplines. Current num_cps: " +
                   std::to_string(num_cps));
    return SplineUtils::circularControlPoints(start_point, end_point, start_direction, num_cps);
  }

  // find closest points --> these will be identical if the extrapolated lines intersect
  const auto closest_points_vector =
      SplineUtils::closestPoints(start_point, end_point, start_direction, end_direction);
  const auto & closest_point_1 = closest_points_vector.first;
  const auto & closest_point_2 = closest_points_vector.second;

  std::vector<Point> first_half = SplineUtils::controlPointsAlongLine(
      start_point, closest_point_1, start_direction, sharpness, cps_per_half);
  std::vector<Point> second_half = SplineUtils::controlPointsAlongLine(
      end_point, closest_point_2, end_direction, sharpness, cps_per_half, true);

  std::vector<Point> control_points;

  // put it all together
  for (const auto i : index_range(first_half))
    control_points.push_back(first_half[i]);
  for (const auto i : index_range(second_half))
    control_points.push_back(second_half[i]);

  return control_points;
}

circularControlPoints()

std::vector< Point > SplineUtils::circularControlPoints	(	const libMesh::Point &	start_point,
		const libMesh::Point &	end_point,
		const libMesh::RealVectorValue &	parallel_direction,
		const unsigned int	num_cps
	)

Creates control points for the special case of parallel lines using an interpolating circle.

Parameters

start_point	start point for spline
end_point	end point for spline
parallel_direction	direction of both curves to be connected. must be the same for both
num_cps	total number of control points to be created

Definition at line 19 of file SplineUtils.C.

Referenced by bSplineControlPoints().

{
  // establish the direction vector between the two points
  const auto orthogonal_direction = start_point - end_point;

  // circular control points will only honor the derivatives if the points are colinear along the
  // orthogonal direction
  if (!(MooseUtils::absoluteFuzzyEqual(parallel_direction.unit() * orthogonal_direction.unit(), 0)))
    mooseError("Parallel lines cannot be interpolated using points along a circle! Direction "
               "vectors are incompatible with point locations.\nStart point: ",
               start_point,
               "\nEnd point: ",
               end_point,
               "\nDirection: ",
               parallel_direction);

  // define the distance between the two points
  const auto distance_between = orthogonal_direction.norm();

  // connecting circle is assumed to be between each point along the orthogonal_direction vector
  const auto radius = distance_between / 2.0;

  // normalize vectors
  const auto unit_normal = orthogonal_direction / distance_between;
  const auto unit_parallel = parallel_direction / parallel_direction.norm();

  // calculate circle center
  const auto center_point = end_point + radius * unit_normal;

  // initialize vector of control points
  std::vector<Point> control_points;

  // loop over the number of control points to generate the control points
  mooseAssert(num_cps > 1, "This routine does not support a single control point");
  for (const auto i : make_range(num_cps))
  {
    const auto t = (Real)i / (Real)(num_cps - 1);
    control_points.push_back(center_point + radius * (std::cos(t * M_PI) * unit_normal +
                                                      std::sin(t * M_PI) * unit_parallel));
  }

  return control_points;
}

closestPoints()

std::pair< Point, Point > SplineUtils::closestPoints	(	const libMesh::Point &	point_1,
		const libMesh::Point &	point_2,
		const libMesh::RealVectorValue &	direction_1,
		const libMesh::RealVectorValue &	direction_2
	)

Determines the two points defining the shortest line segment between two lines in space.

If the lines intersect, the two points returned will be identical.

Parameters

point_1	starting point of the first line
point_2	starting point of the second line
direction_1	direction of the first line
direction_2	direction of the second line

Definition at line 184 of file SplineUtils.C.

Referenced by bSplineControlPoints().

{
  const auto n_vec = direction_1.cross(direction_2);

  // check if n_vec is the zero vector (indicates parallel lines)
  if (n_vec.norm_sq() < libMesh::TOLERANCE * libMesh::TOLERANCE)
    mooseError("Lines are parallel! Infinitely many closest points exist.");

  const auto n1_vec = direction_1.cross(n_vec);
  const auto n2_vec = direction_2.cross(n_vec);

  // calculate closest points
  // take absolute value to ensure specified directions are honored
  libMesh::Real point_a_dir_scalar =
      std::abs((point_2 - point_1) * n2_vec / (direction_1 * n2_vec));
  libMesh::Real point_b_dir_scalar =
      std::abs((point_1 - point_2) * n1_vec / (direction_2 * n1_vec));
  libMesh::Point point_a = point_1 + point_a_dir_scalar * direction_1;
  libMesh::Point point_b = point_2 + point_b_dir_scalar * direction_2;

  // points will be returned as a vector containing two points
  const auto closest_points = std::make_pair(point_a, point_b);

  return closest_points;
}

controlPointsAlongLine()

std::vector< Point > SplineUtils::controlPointsAlongLine	(	const libMesh::Point &	start_point,
		const libMesh::Point &	end_point,
		const libMesh::RealVectorValue &	direction_vector,
		const libMesh::Real	sharpness,
		const unsigned int	num_cps,
		const bool	reverse_order = false
	)

Creates control points along an extrapolated line up to a certain intercept.

Parameters

start_point	start point for the line
end_point	end point of the line
direction_vector	direction of the line
sharpness	number in [0,1] that determines the sharpness of the future curve. In this context, it will determine how close to the end point the control points are placed.
reverse_order	boolean to return the vector of control points in the reverse order

Definition at line 122 of file SplineUtils.C.

Referenced by bSplineControlPoints().

{
  // initialize a vector of control points with the starting point inserted
  std::vector<Point> control_points;

  mooseAssert(num_cps != 0, "Number of control points cannot be zero!");
  if (num_cps == 1)
    control_points.push_back(start_point);

  // create the control points by varying the sharpness linearly
  else
  {
    for (const auto i : libMesh::make_range(num_cps))
    {
      const auto point_sharpness = (Real)i / (Real)(num_cps - 1) * sharpness;
      control_points.push_back(
          SplineUtils::makeControlPoint(start_point, end_point, direction_vector, point_sharpness));
    }

    // reverse_order parameter only set to true when required
    if (reverse_order)
      std::reverse(control_points.begin(), control_points.end());
  }
  return control_points;
}

makeControlPoint()

libMesh::Point SplineUtils::makeControlPoint	(	const libMesh::Point &	start_point,
		const libMesh::Point &	end_point,
		const libMesh::RealVectorValue &	direction_vector,
		const libMesh::Real	sharpness
	)

Creates a control point.

Subroutine of controlPointsAlongLine

Parameters

start_point	start point for the line
end_point	end point of the line
direction_vector	direction of the line
sharpness	number in [0,1] that determines the sharpness of the future curve. In this context, it will determine how close to the end point the control point is placed.

Definition at line 154 of file SplineUtils.C.

Referenced by controlPointsAlongLine().

{
  // check that sharpness value is ok
  mooseAssert(sharpness <= 1.0 && sharpness >= 0, "Sharpness must be in [0,1]!");

  // normalize direction vector
  const auto unit_direction = direction_vector.unit();

  // calculate the distance between points
  const auto distance = (end_point - start_point).norm();
  const auto distance_sq = distance * distance;

  // initialize return quantity and attempt to calculate control vector
  libMesh::Point control_point = start_point + sharpness * distance * unit_direction;

  // check if we went the wrong way, further from the end_point
  if ((end_point - control_point).norm_sq() > distance_sq)
  {
    control_point = start_point - sharpness * distance * unit_direction;
    mooseAssert((end_point - control_point).norm_sq() <= distance_sq,
                "Also going further from end point");
  }

  return control_point;
}