AXOM
Axom provides a robust, flexible software infrastructure for the development of multi-physics applications and computational tools.
axom::primal::NURBSPatch< T, NDIMS > Class Template Reference

Represents a NURBS patch defined by a 2D array of control points. More...

#include </home/docs/checkouts/readthedocs.org/user_builds/axom/checkouts/develop/src/axom/primal/geometry/NURBSPatch.hpp>

Public Types

using VectorType = Vector< T, NDIMS >
 
using PlaneType = Plane< T, NDIMS >
 
using CoordsVec = axom::Array< PointType, 1 >
 
using CoordsMat = axom::Array< PointType, 2 >
 
using WeightsVec = axom::Array< T, 1 >
 
using WeightsMat = axom::Array< T, 2 >
 
using KnotVectorType = KnotVector< T >
 
using BoundingBoxType = BoundingBox< T, NDIMS >
 
using OrientedBoundingBoxType = OrientedBoundingBox< T, NDIMS >
 
using NURBSCurveType = primal::NURBSCurve< T, NDIMS >
 
using TrimmingCurveType = primal::NURBSCurve< T, 2 >
 
using TrimmingCurveVec = axom::Array< TrimmingCurveType >
 
using ParameterPointType = Point< T, 2 >
 
using ParameterBoundingBoxType = BoundingBox< T, 2 >
 

Public Member Functions

 AXOM_STATIC_ASSERT_MSG ((NDIMS==1)||(NDIMS==2)||(NDIMS==3), "A NURBS Patch object may be defined in 1-, 2-, or 3-D")
 
 AXOM_STATIC_ASSERT_MSG (std::is_arithmetic< T >::value, "A NURBS Patch must be defined using an arithmetic type")
 
Query/modify patch properties (degree, rationality, ...)
void setParameters (int npts_u, int npts_v, int deg_u, int deg_v)
 Reset the degree and resize arrays of points (and weights) More...
 
void setDegree_u (int deg)
 Reset the knot vector in u. More...
 
void setDegree_v (int deg)
 Reset the knot vector in v. More...
 
void setDegree (int deg_u, int deg_v)
 Reset the knot vector and increase the number of control points. More...
 
int getDegree_u () const
 Returns the degree of the NURBS Patch on the first axis. More...
 
int getDegree_v () const
 Returns the degree of the NURBS Patch on the second axis. More...
 
int getOrder_u () const
 Returns the order (degree + 1) of the NURBS Patch on the first axis. More...
 
int getOrder_v () const
 Returns the order of the NURBS Patch on the second axis. More...
 
void setNumControlPoints (int npts_u, int npts_v)
 Set the number control points in u. More...
 
int getNumControlPoints_u () const
 Returns the number of control points in the NURBS Patch on the first axis. More...
 
int getNumControlPoints_v () const
 Returns the number of control points in the NURBS Patch on the second axis. More...
 
void setNumControlPoints_u (int npts)
 Set the number control points in u. More...
 
void setNumControlPoints_v (int npts)
 Set the number control points in v. More...
 
void clear ()
 Clears the list of control points, make nonrational. More...
 
bool isRational () const
 Use array size as flag for rationality. More...
 
void makeRational ()
 Make trivially rational. If already rational, do nothing. More...
 
void makeNonrational ()
 Make nonrational by shrinking array of weights. More...
 
bool isValidNURBS () const
 Function to check if the NURBS surface is valid. More...
 
Query/modify patch's geometry (control points, weights, bounding box, ...)
PointTypeoperator() (int ui, int vi)
 Retrieves the control point at index (idx_p, idx_q) More...
 
const PointTypeoperator() (int ui, int vi) const
 Retrieves the vector of control points at index idx. More...
 
CoordsMatgetControlPoints ()
 Returns a reference to the NURBS patch's control points. More...
 
const CoordsMatgetControlPoints () const
 Returns a reference to the NURBS patch's control points. More...
 
const T & getWeight (int ui, int vi) const
 Get a specific weight. More...
 
void setWeight (int ui, int vi, T weight)
 Set the weight at a specific index. More...
 
WeightsMatgetWeights ()
 Returns a reference to the NURBS patch's weights. More...
 
const WeightsMatgetWeights () const
 Returns a const reference to the NURBS patch's weights. More...
 
BoundingBoxType boundingBox () const
 Returns an axis-aligned bounding box containing the patch. More...
 
OrientedBoundingBoxType orientedBoundingBox () const
 Returns an oriented bounding box containing the patch. More...
 
Methods for (untrimmed) Patch Parameterization.

These methods operate on the parameterization of the patch leaving the geometry unchanged.

void setKnot_u (int idx, T knot)
 Set the knot value in the u vector at a specific index. More...
 
void setKnot_v (int idx, T knot)
 Set the knot value in the v vector at a specific index. More...
 
void setKnots_u (const axom::Array< T > &knots, int degree)
 Set the u knot vector by an axom::Array. More...
 
void setKnots_v (const axom::Array< T > &knots, int degree)
 Set the v knot vector by an axom::Array. More...
 
void setKnots_u (const KnotVectorType &knotVector)
 Set the u knot vector by a KnotVector object. More...
 
void setKnots_v (const KnotVectorType &knotVector)
 Set the v knot vector by a KnotVector object. More...
 
KnotVectorTypegetKnots_u ()
 Return a reference to the KnotVector instance on the first axis. More...
 
const KnotVectorTypegetKnots_u () const
 Return a const reference to the KnotVector instance on the first axis. More...
 
getMinKnot_u () const
 Get the minimum knot value in the u-axis. More...
 
getMaxKnot_u () const
 Get the maximum knot value in the u-axis. More...
 
getParameterSpaceDiagonal () const
 Get the length of the parameter space bounding box. More...
 
KnotVectorTypegetKnots_v ()
 Return a reference to the KnotVector instance on the second axis. More...
 
const KnotVectorTypegetKnots_v () const
 Return a const reference to the KnotVector instance on the second axis. More...
 
getMinKnot_v () const
 Get the minimum knot value in the v-axis. More...
 
getMaxKnot_v () const
 Get the maximum knot value in the v-axis. More...
 
int getNumKnots_u () const
 Return the length of the knot vector on the first axis. More...
 
int getNumKnots_v () const
 Return the length of the knot vector on the second axis. More...
 
axom::IndexType insertKnot_u (T u, int target_multiplicity=1)
 Insert a knot to the u knot vector to have the given multiplicity. More...
 
axom::IndexType insertKnot_v (T v, int target_multiplicity=1)
 Insert a knot to the v knot vector to have the given multiplicity. More...
 
void reverseTrimmingCurves ()
 Reverses all trimming curves in the patch. More...
 
void reverseOrientation (int axis)
 Reverses the order of one direction of the NURBS patch's control points and weights. More...
 
void reverseOrientation_u ()
 Reverses the order of the control points, weights, and knots on the first axis. More...
 
void reverseOrientation_v ()
 Reverses the order of the control points, weights, and knots on the second axis. More...
 
void swapAxes ()
 Swap the axes such that s(u, v) becomes s(v, u) More...
 
void normalize ()
 Normalize the knot vectors to the span [0, 1]. More...
 
void normalize_u ()
 Normalize the knot vector in u to the span [0, 1]. More...
 
void normalize_v ()
 Normalize the knot vector in v to the span [0, 1]. More...
 
void normalizeBySpan ()
 Normalize to the span [0, N] x [0, M] where N and M are the number of spans in u and v. More...
 
void rescale (T a, T b)
 Rescale both knot vectors to the span of [a, b]. More...
 
void rescale_u (T a, T b)
 Rescale the knot vector in u to the span of [a, b]. More...
 
void rescale_v (T a, T b)
 Rescale the knot vector in v to the span of [a, b]. More...
 
bool isValidParameter_u (T u, T EPS=1e-8) const
 Function to check if the u parameter is within the knot span. More...
 
bool isValidParameter_v (T v, T EPS=1e-8) const
 Function to check if the v parameter is within the knot span. More...
 
bool isValidInteriorParameter_u (T t) const
 Checks if given u parameter is interior to the knot span. More...
 
bool isValidInteriorParameter_v (T t) const
 Checks if given v parameter is interior to the knot span. More...
 
void scaleParameterSpace (double scaleFactor, bool removeTrimmingCurves=false)
 Scale the parameter space of the NURBS patch geometry linearly (by tangents) in all directions. More...
 
void expandParameterSpace (double expansionAmount_u, double expansionAmount_v, bool removeTrimmingCurves=false)
 Expand the parameter space of the NURBS patch geometry linearly (by tangents) in all directions by a fixed amount. More...
 
Functions to evaluate (untrimmed) patch and its derivatives and normals along points or lines

These methods operate only on the geometry of the patch as defined by the control points and weights. They do not modify nor interact with the trimming curves.

PointType evaluate (T u, T v) const
 Evaluate the NURBS patch geometry at a particular parameter value t. More...
 
NURBSCurveType isocurve (T uv, int axis) const
 Returns a NURBS patch isocurve for a fixed parameter value of u or v. More...
 
NURBSCurveType isocurve_u (T u) const
 Returns a NURBS patch isocurve with a fixed value of u. More...
 
NURBSCurveType isocurve_v (T v) const
 Returns a NURBS patch isocurve with a fixed value of v. More...
 
void evaluateDerivatives (T u, T v, int d, axom::Array< VectorType, 2 > &ders, bool evalByTotalOrder=false) const
 Evaluate the NURBS patch geometry derivatives up to order d at parameter u, v. More...
 
void evaluateFirstDerivatives (T u, T v, PointType &eval, VectorType &Du, VectorType &Dv) const
 Evaluates all first derivatives of the NURBS patch geometry at (u, v) More...
 
void evaluateLinearDerivatives (T u, T v, PointType &eval, VectorType &Du, VectorType &Dv, VectorType &DuDv) const
 Evaluates all linear derivatives of the NURBS patch geometry at (u, v) More...
 
void evaluateSecondDerivatives (T u, T v, PointType &eval, VectorType &Du, VectorType &Dv, VectorType &DuDu, VectorType &DvDv, VectorType &DuDv) const
 Evaluates all second derivatives of the NURBS patch geometry at (u, v) More...
 
VectorType du (T u, T v) const
 Computes a tangent in u of the NURBS patch geometry at (u, v) More...
 
VectorType dv (T u, T v) const
 Computes a tangent in v of the NURBS patch geometry at (u, v) More...
 
VectorType dudu (T u, T v) const
 Computes the second derivative in u of the NURBS patch geometry patch at (u, v) More...
 
VectorType dvdv (T u, T v) const
 Computes the second derivative in v of the NURBS patch geometry patch at (u, v) More...
 
VectorType dudv (T u, T v) const
 Computes the mixed second derivative in u and v of the NURBS patch geometry at (u, v) More...
 
VectorType dvdu (T u, T v) const
 Computes the mixed second derivative in u and v of the NURBS patch geometry at (u, v) More...
 
VectorType normal (T u, T v) const
 Computes the normal vector to the NURBS patch geometry at (u, v) More...
 
VectorType calculateUntrimmedPatchNormal (int npts=20) const
 Calculate the average normal for the (untrimmed) patch. More...
 
Methods for (trimmed) Surface Geometry.

These methods operate on the visible geometry of the NURBS surface as defined by the patch geometry and the trimming curves which define visibility.

const TrimmingCurveVecgetTrimmingCurves () const
 Get array of trimming curves. More...
 
TrimmingCurveVecgetTrimmingCurves ()
 Get mutable array of trimming curves. More...
 
const TrimmingCurveTypegetTrimmingCurve (int idx) const
 Get a trimming curve by index. More...
 
void addTrimmingCurve (const TrimmingCurveType &curve)
 Add a trimming curve. More...
 
void addTrimmingCurves (const TrimmingCurveVec &curves)
 Add array of trimming curves. More...
 
void setTrimmingCurves (const TrimmingCurveVec &curves)
 Set the array of trimming curves. More...
 
void clearTrimmingCurves ()
 Clear trimming curves, but DON'T mark as untrimmed. More...
 
int getNumTrimmingCurves () const
 Get number of trimming curves. More...
 
bool isTriviallyTrimmed (double tol=1e-8) const
 Predicate to check if the patch is "trivially trimmed" in parameter space. More...
 
bool isTrimmed () const
 use boolean flag for trimmed-ness More...
 
bool isInvisible () const
 Predicate to check if the patch is entirely invisible due to trimming state. More...
 
void markAsTrimmed ()
 Mark as trimmed. More...
 
void makeUntrimmed ()
 Delete all trimming curves. More...
 
void makeTriviallyTrimmed ()
 Make trivially trimmed by adding trimming curves at each boundary. More...
 
bool isVisible (T u, T v) const
 Check if a parameter point is visible on the NURBS patch via a trim test. More...
 
void clip (T min_u, T max_u, T min_v, T max_v, bool normalizeParameters=false)
 Restrict the edges of a NURBS surface to the given parameter values, if necessary. More...
 
void clipToCurves (double padding=1e-5)
 Clip the edges of a NURBS surface to the AABB of the trimming curves. More...
 
Functions dealing with (untrimmed) patch subdivision
axom::Array< BezierPatch< T, NDIMS > > extractBezier () const
 Splits the NURBS patch geometry (at each internal knot) into several Bezier patches. More...
 
axom::Array< NURBSPatch< T, NDIMS > > extractTrimmedBezier () const
 Splits the NURBS patch geometry into several one-span trimmed NURBS surfaces. More...
 
VectorType calculateTrimmedPatchNormal (int npts=20, bool useBezierExtraction=true) const
 Calculate the average normal for the trimmed patch. More...
 
std::pair< VectorType, double > calculateTrimmedPatchNormalArea (int npts=20, bool useBezierExtraction=true) const
 Calculate the unsigned surface area and integrated normal for the trimmed patch. More...
 
NURBSPatch cleanedTrimmedRepresentation (bool normalize_by_span=true, bool ensure_trimmed=true, bool clip_to_curves=true) const
 Return a "clean" trimmed representation suitable for moment and GWN calculation. More...
 
template<int ORDER, int NVALS = 4 + (ORDER >= 0 ? 3 : 0) + (ORDER >= 1 ? 9 : 0) + (ORDER >= 2 ? 27 : 0)>
primal::Vector< T, NVALS > calculateSurfaceMoments (int npts=20, bool useBezierExtraction=false) const
 Calculate the surface moments for the trimmed patch for GWN evaluation. More...
 

Constructors for NURBSPatch

The NURBSPatch class provides a variety of constructors to support flexible initialization. The default constructor creates an empty, invalid patch. Other constructors allow you to specify degrees, control point counts, knot vectors, and weights, either using C-style arrays or axom::Array containers, in both 1D and 2D forms. You can also construct a NURBSPatch from a BezierPatch.

Depending on the constructor used, the resulting instance will have:

  • Control points and weights arrays sized according to the provided parameters.
  • Knot vectors initialized either as uniform (for degree/size-based constructors) or from provided arrays/objects.
  • Rationality determined by the presence of weights (patches are rational if weights are provided, nonrational otherwise).
  • Trimming curves are always empty and the patch is untrimmed by default.

For 1D arrays, the mapping of control points and weights to the patch is lexicographical, i.e.

 pts[0]                 -> nodes[0, 0],      ..., pts[npts_v]        -> nodes[0, npts_v]
 pts[npts_v+1]          -> nodes[1, 0],      ..., pts[2*npts_v]      -> nodes[1, npts_v]
                                             ...
 pts[npts_u*(npts_v-1)] -> nodes[npts_u, 0], ..., pts[npts_u*npts_v] -> nodes[npts_u, npts_v]

All constructors ensure that the patch is internally consistent and valid, provided the input parameters are valid.

using PointType = Point< T, NDIMS >
 
 NURBSPatch ()
 Default constructor for an empty (invalid) NURBS patch. More...
 
 NURBSPatch (int deg_u, int deg_v)
 Constructor for a simple NURBS surface that reserves space for the minimum (sensible) number of points for the given degrees. More...
 
 NURBSPatch (int npts_u, int npts_v, int deg_u, int deg_v)
 Constructor for a simple NURBS surface that reserves space for npts_u * npts_v control points. More...
 
 NURBSPatch (const BezierPatch< T, NDIMS > &bezierPatch)
 Constructor for a NURBS surface from a Bezier surface. More...
 
 NURBSPatch (ArrayView< const PointType, 2 > controlPoints, ArrayView< const T, 2 > weights, const KnotVectorType &knotVector_u, const KnotVectorType &knotVector_v)
 Constructor for a NURBSPatch from 2D ArrayViews of control points and weights and KnotVectors for the u- and v- directions. More...
 
 NURBSPatch (ArrayView< const PointType, 2 > controlPoints, const KnotVectorType &knotVector_u, const KnotVectorType &knotVector_v)
 Constructor for a NURBSPatch from 2D ArrayViews of control points and KnotVectors for the u- and v- directions. More...
 
 NURBSPatch (ArrayView< PointType, 2 > controlPoints, const KnotVectorType &knotVector_u, const KnotVectorType &knotVector_v)
 
 NURBSPatch (ArrayView< PointType, 2 > controlPoints, ArrayView< T, 2 > weights, const KnotVectorType &knotVector_u, const KnotVectorType &knotVector_v)
 
 NURBSPatch (const PointType *pts, int npts_u, int npts_v, int deg_u, int deg_v)
 Constructor for a NURBS Patch from an array of coordinates and degrees. More...
 
 NURBSPatch (const PointType *pts, const T *weights, int npts_u, int npts_v, int deg_u, int deg_v)
 Constructor for a NURBS Patch from arrays of coordinates and weights. More...
 
 NURBSPatch (const CoordsVec &pts, int npts_u, int npts_v, int deg_u, int deg_v)
 Constructor for a NURBS Patch from 1D arrays of coordinates and degrees. More...
 
 NURBSPatch (const CoordsVec &pts, const WeightsVec &weights, int npts_u, int npts_v, int deg_u, int deg_v)
 Constructor for a NURBS Patch from 1D arrays of coordinates and weights. More...
 
 NURBSPatch (const CoordsMat &pts, int deg_u, int deg_v)
 Constructor for a NURBS Patch from 2D arrays of coordinates and degrees. More...
 
 NURBSPatch (const CoordsMat &pts, const WeightsMat &weights, int deg_u, int deg_v)
 Constructor for a NURBS Patch from 2D arrays of coordinates and weights. More...
 
 NURBSPatch (const PointType *pts, int npts_u, int npts_v, const T *knots_u, int nkts_u, const T *knots_v, int nkts_v)
 Constructor for a NURBS Patch from C-style arrays of coordinates and knot vectors. More...
 
 NURBSPatch (const PointType *pts, const T *weights, int npts_u, int npts_v, const T *knots_u, int nkts_u, const T *knots_v, int nkts_v)
 Constructor for a NURBS Patch from C-style arrays of coordinates and weights. More...
 
 NURBSPatch (const CoordsVec &pts, int npts_u, int npts_v, const axom::Array< T > &knots_u, const axom::Array< T > &knots_v)
 Constructor for a NURBS Patch from 1D axom::Array arrays of coordinates and knots. More...
 
 NURBSPatch (const CoordsVec &pts, const WeightsVec &weights, int npts_u, int npts_v, const axom::Array< T > &knots_u, const axom::Array< T > &knots_v)
 Constructor for a NURBS Patch from 1D axom::Array arrays of coordinates, weights, and knots. More...
 
 NURBSPatch (const CoordsVec &pts, int npts_u, int npts_v, const KnotVectorType &knotvec_u, const KnotVectorType &knotvec_v)
 Constructor for a NURBS Patch from 1D axom::Array arrays of coordinates and KnotVectors. More...
 
 NURBSPatch (const CoordsVec &pts, const WeightsVec &weights, int npts_u, int npts_v, const KnotVectorType &knotvec_u, const KnotVectorType &knotvec_v)
 Constructor for a NURBS Patch from 1D axom::Array arrays of coordinates, weights, and KnotVectors. More...
 
 NURBSPatch (const CoordsMat &pts, const axom::Array< T > &knots_u, const axom::Array< T > &knots_v)
 Constructor for a NURBS Patch from 2D axom::Array array of coordinates and array of knots. More...
 
 NURBSPatch (const CoordsMat &pts, const WeightsMat &weights, const axom::Array< T > &knots_u, const axom::Array< T > &knots_v)
 Constructor for a NURBS Patch from 2D axom::Array array of coordinates, weights, and array of knots. More...
 
 NURBSPatch (const CoordsMat &pts, const KnotVectorType &knotvec_u, const KnotVectorType &knotvec_v)
 Constructor for a NURBS Patch from 1D axom::Array array of coordinates and KnotVector objects. More...
 
 NURBSPatch (const CoordsMat &pts, const WeightsMat &weights, const KnotVectorType &knotvec_u, const KnotVectorType &knotvec_v)
 Constructor for a NURBS Patch from 2D axom::Array array of coordinates, weights, and KnotVector objects. More...
 

Functions dealing with trimmed patch subdivision

bool split (T u, T v, NURBSPatch &p1, NURBSPatch &p2, NURBSPatch &p3, NURBSPatch &p4, bool normalizeParameters=false) const
 Splits a NURBS surface into four NURBS patches. More...
 
bool split_u (T u, NURBSPatch &p1, NURBSPatch &p2, bool normalizeParameters=false) const
 Split the NURBS surface in two along the u direction. More...
 
bool split_v (T v, NURBSPatch &p1, NURBSPatch &p2, bool normalizeParameters=false) const
 Split the NURBS surface in two along the v direction. More...
 
void nearBisectOnLongerAxis (NURBSPatch &p1, NURBSPatch &p2) const
 Split the patch in the longer parametric direction, determined via maximum control-net polyline length. More...
 
void diskSplit (T u, T v, T r, NURBSPatch &the_disk, NURBSPatch &the_rest, bool clipDisk=true) const
 For a disk of radius r and center (u, v), split a NURBS surface into the portion inside/outside. More...
 
void diskSplit (T u, T v, T r, NURBSPatch &the_disk, NURBSPatch &the_rest, bool &isDiskInside, bool &isDiskOutside, bool ignoreInteriorDisk, double disk_padding) const
 Split a NURBS surface into two by cutting out a disk of radius r centered at (u, v) More...
 
std::ostream & print (std::ostream &os) const
 Simple formatted print of a NURBS Patch instance. More...
 
bool operator== (const NURBSPatch< T, NDIMS > &lhs, const NURBSPatch< T, NDIMS > &rhs)
 Equality operator for NURBS patches. More...
 
bool operator!= (const NURBSPatch< T, NDIMS > &lhs, const NURBSPatch< T, NDIMS > &rhs)
 Inequality operator for NURBS patches. More...
 

Detailed Description

template<typename T, int NDIMS>
class axom::primal::NURBSPatch< T, NDIMS >

Represents a NURBS patch defined by a 2D array of control points.

Template Parameters
Tthe coordinate type, e.g., double, float, etc.

A NURBS patch has degrees p and q with knot vectors of length r+1 and s+1 respectively. There is a control net of (n + 1) * (m + 1) points with r+1 = n+p+2 and s+1 = m+q+2. Optionally has an equal number of weights for rational patches.

The patch must be open (clamped on all boundaries) and continuous (unless p = 0 or q = 0)

Nonrational NURBS patches are identified by an empty weights array.

A NURBS surface is identified as trimmed with an internal flag, as an empty trimming curve vector indicates a patch with no visible surface.

Member Typedef Documentation

◆ PointType

template<typename T , int NDIMS>
Overload for non const axom::primal::NURBSPatch< T, NDIMS >::PointType

This is an overloaded member function, provided for convenience. It differs from the above function only in what argument(s) it accepts.

◆ VectorType

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::VectorType = Vector<T, NDIMS>

◆ PlaneType

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::PlaneType = Plane<T, NDIMS>

◆ CoordsVec

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::CoordsVec = axom::Array<PointType, 1>

◆ CoordsMat

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::CoordsMat = axom::Array<PointType, 2>

◆ WeightsVec

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::WeightsVec = axom::Array<T, 1>

◆ WeightsMat

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::WeightsMat = axom::Array<T, 2>

◆ KnotVectorType

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::KnotVectorType = KnotVector<T>

◆ BoundingBoxType

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::BoundingBoxType = BoundingBox<T, NDIMS>

◆ OrientedBoundingBoxType

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::OrientedBoundingBoxType = OrientedBoundingBox<T, NDIMS>

◆ NURBSCurveType

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::NURBSCurveType = primal::NURBSCurve<T, NDIMS>

◆ TrimmingCurveType

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::TrimmingCurveType = primal::NURBSCurve<T, 2>

◆ TrimmingCurveVec

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::TrimmingCurveVec = axom::Array<TrimmingCurveType>

◆ ParameterPointType

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::ParameterPointType = Point<T, 2>

◆ ParameterBoundingBoxType

template<typename T , int NDIMS>
using axom::primal::NURBSPatch< T, NDIMS >::ParameterBoundingBoxType = BoundingBox<T, 2>

Constructor & Destructor Documentation

◆ NURBSPatch() [1/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( )
inline

Default constructor for an empty (invalid) NURBS patch.

Note
An empty NURBS patch is not valid

◆ NURBSPatch() [2/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( int  deg_u,
int  deg_v 
)
inline

Constructor for a simple NURBS surface that reserves space for the minimum (sensible) number of points for the given degrees.

Parameters
[in]deg_u,deg_vThe patch's degrees on the first and second axis
Precondition
deg_u, deg_v both greater than or equal to 0, or both -1

◆ NURBSPatch() [3/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( int  npts_u,
int  npts_v,
int  deg_u,
int  deg_v 
)
inline

Constructor for a simple NURBS surface that reserves space for npts_u * npts_v control points.

Parameters
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]deg_u,deg_vThe patch's degrees on the first and second axis
Precondition
Requires npts_d > deg_d and deg_d >= 0 for d = u, v

References axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), and SLIC_ASSERT.

◆ NURBSPatch() [4/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const BezierPatch< T, NDIMS > &  bezierPatch)
inlineexplicit

Constructor for a NURBS surface from a Bezier surface.

Parameters
[in]bezierPatchthe Bezier patch to convert to a NURBS patch

◆ NURBSPatch() [5/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( ArrayView< const PointType, 2 >  controlPoints,
ArrayView< const T, 2 >  weights,
const KnotVectorType knotVector_u,
const KnotVectorType knotVector_v 
)
inline

Constructor for a NURBSPatch from 2D ArrayViews of control points and weights and KnotVectors for the u- and v- directions.

Parameters
[in]controlPoints2D ArrayView of control points
[in]weights2D ArrayView of weights
[in]knotVector_u,knotVector_vThe knot vectors in the u- and v- directions
Precondition
If controlPoints is not empty, its sizes must match the number of control points implied by the knot vectors: knotVector_u.getNumControlPoints() * knotVector_v.getNumberControlPoints()
If weights is not empty, its dimensions must match that of the controlPoints
The KnotVector degrees must be valid: knotVector_u.getDegree() >= -1 and knotVector_v.getDegree() >= -1.
If the degrees are not both -1 (i.e., not empty), then: They must both be non-negative and the number of control points must be greature than the degree for both axes.

References AXOM_MAYBE_UNUSED, axom::ArrayView< T, DIM, SPACE >::data(), axom::ArrayView< T, DIM, SPACE >::empty(), axom::primal::KnotVector< T >::getDegree(), axom::primal::KnotVector< T >::getNumControlPoints(), axom::primal::NURBSPatch< T, NDIMS >::isValidNURBS(), axom::utilities::max(), axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), axom::ArrayBase< T, DIM, ArrayType >::shape(), and SLIC_ASSERT.

◆ NURBSPatch() [6/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( ArrayView< const PointType, 2 >  controlPoints,
const KnotVectorType knotVector_u,
const KnotVectorType knotVector_v 
)
inline

Constructor for a NURBSPatch from 2D ArrayViews of control points and KnotVectors for the u- and v- directions.

Parameters
[in]controlPoints2D ArrayView of control points
[in]weights2D ArrayView of weights
[in]knotVector_u,knotVector_vThe knot vectors in the u- and v- directions

◆ NURBSPatch() [7/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( ArrayView< PointType, 2 >  controlPoints,
const KnotVectorType knotVector_u,
const KnotVectorType knotVector_v 
)
inline

◆ NURBSPatch() [8/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( ArrayView< PointType, 2 >  controlPoints,
ArrayView< T, 2 >  weights,
const KnotVectorType knotVector_u,
const KnotVectorType knotVector_v 
)
inline

◆ NURBSPatch() [9/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const PointType pts,
int  npts_u,
int  npts_v,
int  deg_u,
int  deg_v 
)
inline

Constructor for a NURBS Patch from an array of coordinates and degrees.

Parameters
[in]ptsA 1D C-style array of npts_u*npts_v control points
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]deg_u,deg_vThe patch's degree on the first and second axis
Precondition
Requires that npts_d >= deg_d + 1 and deg_d >= 0 for d = u, v

◆ NURBSPatch() [10/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const PointType pts,
const T *  weights,
int  npts_u,
int  npts_v,
int  deg_u,
int  deg_v 
)
inline

Constructor for a NURBS Patch from arrays of coordinates and weights.

Parameters
[in]ptsA 1D C-style array of (ord_u+1)*(ord_v+1) control points
[in]weightsA 1D C-style array of (ord_u+1)*(ord_v+1) positive weights
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]deg_u,deg_vThe patch's degree on the first and second axis
Precondition
Requires that npts_d >= deg_d + 1 and deg_d >= 0 for d = u, v

◆ NURBSPatch() [11/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsVec pts,
int  npts_u,
int  npts_v,
int  deg_u,
int  deg_v 
)
inline

Constructor for a NURBS Patch from 1D arrays of coordinates and degrees.

Parameters
[in]ptsA 1D axom::Array of npts_u*npts_v control points
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]deg_u,deg_vThe patch's degree on the first and second axis
Precondition
Requires that npts_d >= deg_d + 1 and deg_d >= 0 for d = u, v

◆ NURBSPatch() [12/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsVec pts,
const WeightsVec weights,
int  npts_u,
int  npts_v,
int  deg_u,
int  deg_v 
)
inline

Constructor for a NURBS Patch from 1D arrays of coordinates and weights.

Parameters
[in]ptsA 1D axom::Array of (ord_u+1)*(ord_v+1) control points
[in]weightsA 1D axom::Array of (ord_u+1)*(ord_v+1) positive weights
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]deg_u,deg_vThe patch's degree on the first and second axis
Precondition
Requires that npts_d >= deg_d + 1 and deg_d >= 0 for d = u, v

◆ NURBSPatch() [13/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsMat pts,
int  deg_u,
int  deg_v 
)
inline

Constructor for a NURBS Patch from 2D arrays of coordinates and degrees.

Parameters
[in]ptsA 2D axom::Array of (npts_u, npts_v) control points
[in]deg_u,deg_vThe patch's degree on the first and second axis
Precondition
Requires that npts_d >= deg_d + 1 and deg_d >= 0 for d = u, v

◆ NURBSPatch() [14/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsMat pts,
const WeightsMat weights,
int  deg_u,
int  deg_v 
)
inline

Constructor for a NURBS Patch from 2D arrays of coordinates and weights.

Parameters
[in]ptsA 2D axom::Array of (ord_u+1, ord_v+1) control points
[in]weightsA 2D axom::Array of (ord_u+1, ord_v+1) positive weights
[in]deg_u,deg_vThe patch's degree on the first and second axis
Precondition
Requires that npts_d >= deg_d + 1 and deg_d >= 0 for d = u, v

◆ NURBSPatch() [15/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const PointType pts,
int  npts_u,
int  npts_v,
const T *  knots_u,
int  nkts_u,
const T *  knots_v,
int  nkts_v 
)
inline

Constructor for a NURBS Patch from C-style arrays of coordinates and knot vectors.

Parameters
[in]ptsA 1D C-style array of npts_u*npts_v control points
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]knots_uA 1D C-style array of npts_u + deg_u + 1 knots
[in]nkts_uThe number of knots in the u direction
[in]knots_vA 1D C-style array of npts_v + deg_v + 1 knots
[in]nkts_vThe number of knots in the v direction

For clamped and continuous patch axes, npts and the knot vector uniquely determine the degree

Precondition
Requires valid pointers and knot vectors

◆ NURBSPatch() [16/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const PointType pts,
const T *  weights,
int  npts_u,
int  npts_v,
const T *  knots_u,
int  nkts_u,
const T *  knots_v,
int  nkts_v 
)
inline

Constructor for a NURBS Patch from C-style arrays of coordinates and weights.

Parameters
[in]ptsA 1D C-style array of npts_u*npts_v control points
[in]weightsA 1D C-style array of npts_u*npts_v positive weights
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]knots_uA 1D C-style array of npts_u + deg_u + 1 knots
[in]nkts_uThe number of knots in the u direction
[in]knots_vA 1D C-style array of npts_v + deg_v + 1 knots
[in]nkts_vThe number of knots in the v direction

For clamped and continuous patch axes, npts and the knot vector uniquely determine the degree

Precondition
Requires valid pointers and knot vectors

◆ NURBSPatch() [17/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsVec pts,
int  npts_u,
int  npts_v,
const axom::Array< T > &  knots_u,
const axom::Array< T > &  knots_v 
)
inline

Constructor for a NURBS Patch from 1D axom::Array arrays of coordinates and knots.

Parameters
[in]ptsA 1D axom::Array of npts_u*npts_v control points
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]knots_uAn axom::Array of npts_u + deg_u + 1 knots
[in]knots_vAn axom::Array of npts_v + deg_v + 1 knots

For clamped and continuous patch axes, npts and the knot vector uniquely determine the degree

Precondition
Requires a valid knot vector and npts_d > deg_d

◆ NURBSPatch() [18/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsVec pts,
const WeightsVec weights,
int  npts_u,
int  npts_v,
const axom::Array< T > &  knots_u,
const axom::Array< T > &  knots_v 
)
inline

Constructor for a NURBS Patch from 1D axom::Array arrays of coordinates, weights, and knots.

Parameters
[in]ptsA 1D axom::Array of npts_u*npts_v control points
[in]weightsA 1D axom::Array of npts_u*npts_v positive weights
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]knots_uAn axom::Array of npts_u + deg_u + 1 knots
[in]knots_vAn axom::Array of npts_v + deg_v + 1 knots

For clamped and continuous patch axes, npts and the knot vector uniquely determine the degree

Precondition
Requires a valid knot vector and npts_d > deg_d

◆ NURBSPatch() [19/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsVec pts,
int  npts_u,
int  npts_v,
const KnotVectorType knotvec_u,
const KnotVectorType knotvec_v 
)
inline

Constructor for a NURBS Patch from 1D axom::Array arrays of coordinates and KnotVectors.

Parameters
[in]ptsA 1D axom::Array of npts_u*npts_v control points
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]knotvec_u,knotvec_vKnotVector objects for the first and second axis

For clamped and continuous patch axes, npts and the knot vectoruniquely determine the degree

Precondition
Requires a valid knot vector and npts_d > deg_d

◆ NURBSPatch() [20/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsVec pts,
const WeightsVec weights,
int  npts_u,
int  npts_v,
const KnotVectorType knotvec_u,
const KnotVectorType knotvec_v 
)
inline

Constructor for a NURBS Patch from 1D axom::Array arrays of coordinates, weights, and KnotVectors.

Parameters
[in]ptsA 1D axom::Array of npts_u*npts_v control points
[in]weightsA 1D axom::Array of npts_u*npts_v positive weights
[in]npts_u,npts_vThe number of control points on the first and second axis
[in]knotvec_u,knotvec_vKnotVector objects for the first and second axis

For clamped and continuous patch axes, npts and the knot vector uniquely determine the degree

Precondition
Requires a valid knot vector and npts_d > deg_d

◆ NURBSPatch() [21/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsMat pts,
const axom::Array< T > &  knots_u,
const axom::Array< T > &  knots_v 
)
inline

Constructor for a NURBS Patch from 2D axom::Array array of coordinates and array of knots.

Parameters
[in]ptsA 2D axom::Array of (npts_u, npts_v) control points
[in]knots_uAn axom::Array of npts_u + deg_u + 1 knots
[in]knots_vAn axom::Array of npts_v + deg_v + 1 knots

For clamped and continuous patch axes, npts and the knot vector uniquely determine the degree

Precondition
Requires a valid knot vector and npts_d > deg_d

◆ NURBSPatch() [22/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsMat pts,
const WeightsMat weights,
const axom::Array< T > &  knots_u,
const axom::Array< T > &  knots_v 
)
inline

Constructor for a NURBS Patch from 2D axom::Array array of coordinates, weights, and array of knots.

Parameters
[in]ptsA 2D axom::Array of (ord_u+1, ord_v+1) control points
[in]weightsA 2D axom::Array of (ord_u+1, ord_v+1) positive weights
[in]knots_uAn axom::Array of npts_u + deg_u + 1 knots
[in]knots_vAn axom::Array of npts_v + deg_v + 1 knots

For clamped and continuous patch axes, npts and the knot vector uniquely determine the degree

Precondition
Requires a valid knot vector and npts_d > deg_d

◆ NURBSPatch() [23/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsMat pts,
const KnotVectorType knotvec_u,
const KnotVectorType knotvec_v 
)
inline

Constructor for a NURBS Patch from 1D axom::Array array of coordinates and KnotVector objects.

Parameters
[in]ptsA 2D axom::Array of (ord_u+1, ord_v+1) control points
[in]knotvec_u,knotvec_vKnotVector objects for the first and second axis

For clamped and continuous patch axes, npts and the knot vector uniquely determine the degree

Precondition
Requires a valid knot vector and npts_d > deg_d

◆ NURBSPatch() [24/24]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::NURBSPatch ( const CoordsMat pts,
const WeightsMat weights,
const KnotVectorType knotvec_u,
const KnotVectorType knotvec_v 
)
inline

Constructor for a NURBS Patch from 2D axom::Array array of coordinates, weights, and KnotVector objects.

Parameters
[in]ptsA 2D axom::Array of (ord_u+1, ord_v+1) control points
[in]weightsA 2D axom::Array of (ord_u+1, ord_v+1) positive weights
[in]knotvec_u,knotvec_vKnotVector objects for the first and second axis

For clamped and continuous patch axes, npts and the knot vector uniquely determine the degree

Precondition
Requires a valid knot vector and npts_d > deg_d

Member Function Documentation

◆ AXOM_STATIC_ASSERT_MSG() [1/2]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::AXOM_STATIC_ASSERT_MSG ( (NDIMS==1)||(NDIMS==2)||(NDIMS==3)  ,
"A NURBS Patch object may be defined in 1-  ,
2-  ,
or 3-D"   
)

◆ AXOM_STATIC_ASSERT_MSG() [2/2]

template<typename T , int NDIMS>
axom::primal::NURBSPatch< T, NDIMS >::AXOM_STATIC_ASSERT_MSG ( std::is_arithmetic< T >::value  ,
"A NURBS Patch must be defined using an arithmetic type"   
)

◆ setParameters()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setParameters ( int  npts_u,
int  npts_v,
int  deg_u,
int  deg_v 
)
inline

Reset the degree and resize arrays of points (and weights)

Parameters
[in]npts_u,npts_vThe target number of control points on the first and second axis
[in]deg_u,deg_vThe target degrees on the first and second axis
Warning
This method will replace existing knot vectors with a uniform one.

References axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::primal::NURBSPatch< T, NDIMS >::makeNonrational(), axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), and SLIC_ASSERT.

◆ setDegree_u()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setDegree_u ( int  deg)
inline

Reset the knot vector in u.

Parameters
[in]degThe target degree
Warning
This method does NOT change the existing control points, i.e. does not perform degree elevation/reduction. Will replace existing knot vector with a uniform one.
Precondition
Requires deg_u < npts_u and deg >= 0

References axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_u(), axom::primal::KnotVector< T >::makeUniform(), and SLIC_ASSERT.

◆ setDegree_v()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setDegree_v ( int  deg)
inline

Reset the knot vector in v.

Parameters
[in]degThe target degree
Warning
This method does NOT change the existing control points, i.e. does not perform degree elevation/reduction. Will replace existing knot vector with a uniform one.
Precondition
Requires deg_v < npts_v and deg >= 0

References axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_v(), axom::primal::KnotVector< T >::makeUniform(), and SLIC_ASSERT.

◆ setDegree()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setDegree ( int  deg_u,
int  deg_v 
)
inline

Reset the knot vector and increase the number of control points.

Parameters
[in]deg_uThe target degree in u
[in]deg_vThe target degree in v
Warning
This method does NOT change the existing control points, i.e. is not performing degree elevation or reduction.
Precondition
Requires deg_u < npts_u and deg_v < npts_v

References axom::primal::NURBSPatch< T, NDIMS >::setDegree_u(), and axom::primal::NURBSPatch< T, NDIMS >::setDegree_v().

◆ getDegree_u()

template<typename T , int NDIMS>
int axom::primal::NURBSPatch< T, NDIMS >::getDegree_u ( ) const
inline

Returns the degree of the NURBS Patch on the first axis.

References axom::primal::KnotVector< T >::getDegree().

◆ getDegree_v()

template<typename T , int NDIMS>
int axom::primal::NURBSPatch< T, NDIMS >::getDegree_v ( ) const
inline

Returns the degree of the NURBS Patch on the second axis.

References axom::primal::KnotVector< T >::getDegree().

◆ getOrder_u()

template<typename T , int NDIMS>
int axom::primal::NURBSPatch< T, NDIMS >::getOrder_u ( ) const
inline

Returns the order (degree + 1) of the NURBS Patch on the first axis.

References axom::primal::KnotVector< T >::getDegree().

◆ getOrder_v()

template<typename T , int NDIMS>
int axom::primal::NURBSPatch< T, NDIMS >::getOrder_v ( ) const
inline

Returns the order of the NURBS Patch on the second axis.

References axom::primal::KnotVector< T >::getDegree().

◆ setNumControlPoints()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setNumControlPoints ( int  npts_u,
int  npts_v 
)
inline

Set the number control points in u.

Parameters
[in]nptsThe target number of control points
Warning
This method does NOT maintain the patch shape, i.e. is not performing knot insertion/removal. Will replace existing knot vectots with uniform ones.

References axom::primal::NURBSPatch< T, NDIMS >::getDegree_u(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_v(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::primal::KnotVector< T >::makeUniform(), axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), and SLIC_ASSERT.

◆ getNumControlPoints_u()

template<typename T , int NDIMS>
int axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_u ( ) const
inline

Returns the number of control points in the NURBS Patch on the first axis.

References axom::ArrayBase< T, DIM, ArrayType >::shape().

◆ getNumControlPoints_v()

template<typename T , int NDIMS>
int axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_v ( ) const
inline

Returns the number of control points in the NURBS Patch on the second axis.

References axom::ArrayBase< T, DIM, ArrayType >::shape().

◆ setNumControlPoints_u()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setNumControlPoints_u ( int  npts)
inline

Set the number control points in u.

Parameters
[in]nptsThe target number of control points
Warning
This method does NOT maintain the patch shape, i.e. is not performing knot insertion/removal.

References axom::primal::NURBSPatch< T, NDIMS >::getDegree_u(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_v(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::primal::KnotVector< T >::makeUniform(), axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), and SLIC_ASSERT.

◆ setNumControlPoints_v()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setNumControlPoints_v ( int  npts)
inline

Set the number control points in v.

Parameters
[in]nptsThe target number of control points
Warning
This method does NOT maintain the patch shape, i.e. is not performing knot insertion/removal.

References axom::primal::NURBSPatch< T, NDIMS >::getDegree_v(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_u(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::primal::KnotVector< T >::makeUniform(), axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), and SLIC_ASSERT.

◆ clear()

◆ isRational()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::isRational ( ) const
inline

Use array size as flag for rationality.

References axom::Array< T, DIM, SPACE, StoragePolicy >::empty().

◆ makeRational()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::makeRational ( )
inline

◆ makeNonrational()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::makeNonrational ( )
inline

Make nonrational by shrinking array of weights.

References axom::Array< T, DIM, SPACE, StoragePolicy >::clear().

◆ isValidNURBS()

◆ operator()() [1/2]

template<typename T , int NDIMS>
PointType& axom::primal::NURBSPatch< T, NDIMS >::operator() ( int  ui,
int  vi 
)
inline

Retrieves the control point at index (idx_p, idx_q)

◆ operator()() [2/2]

template<typename T , int NDIMS>
const PointType& axom::primal::NURBSPatch< T, NDIMS >::operator() ( int  ui,
int  vi 
) const
inline

Retrieves the vector of control points at index idx.

◆ getControlPoints() [1/2]

template<typename T , int NDIMS>
CoordsMat& axom::primal::NURBSPatch< T, NDIMS >::getControlPoints ( )
inline

Returns a reference to the NURBS patch's control points.

◆ getControlPoints() [2/2]

template<typename T , int NDIMS>
const CoordsMat& axom::primal::NURBSPatch< T, NDIMS >::getControlPoints ( ) const
inline

Returns a reference to the NURBS patch's control points.

◆ getWeight()

template<typename T , int NDIMS>
const T& axom::primal::NURBSPatch< T, NDIMS >::getWeight ( int  ui,
int  vi 
) const
inline

Get a specific weight.

Parameters
[in]uiThe index of the weight on the first axis
[in]viThe index of the weight on the second axis
Precondition
Requires that the surface be rational

References axom::primal::NURBSPatch< T, NDIMS >::isRational(), and SLIC_ASSERT.

◆ setWeight()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setWeight ( int  ui,
int  vi,
weight 
)
inline

Set the weight at a specific index.

Parameters
[in]uiThe index of the weight in on the first axis
[in]viThe index of the weight in on the second axis
[in]weightThe updated value of the weight
Precondition
Requires that the surface be rational
Requires that the weight be positive

References axom::primal::NURBSPatch< T, NDIMS >::isRational(), and SLIC_ASSERT.

◆ getWeights() [1/2]

template<typename T , int NDIMS>
WeightsMat& axom::primal::NURBSPatch< T, NDIMS >::getWeights ( )
inline

Returns a reference to the NURBS patch's weights.

◆ getWeights() [2/2]

template<typename T , int NDIMS>
const WeightsMat& axom::primal::NURBSPatch< T, NDIMS >::getWeights ( ) const
inline

Returns a const reference to the NURBS patch's weights.

◆ boundingBox()

template<typename T , int NDIMS>
BoundingBoxType axom::primal::NURBSPatch< T, NDIMS >::boundingBox ( ) const
inline

Returns an axis-aligned bounding box containing the patch.

References axom::Array< T, DIM, SPACE, StoragePolicy >::data(), and axom::Array< T, DIM, SPACE, StoragePolicy >::size().

◆ orientedBoundingBox()

template<typename T , int NDIMS>
OrientedBoundingBoxType axom::primal::NURBSPatch< T, NDIMS >::orientedBoundingBox ( ) const
inline

Returns an oriented bounding box containing the patch.

References axom::Array< T, DIM, SPACE, StoragePolicy >::data(), and axom::Array< T, DIM, SPACE, StoragePolicy >::size().

◆ setKnot_u()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setKnot_u ( int  idx,
knot 
)
inline

Set the knot value in the u vector at a specific index.

Parameters
[in]idxThe index of the knot
[in]knotThe updated value of the knot

◆ setKnot_v()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setKnot_v ( int  idx,
knot 
)
inline

Set the knot value in the v vector at a specific index.

Parameters
[in]idxThe index of the knot
[in]knotThe updated value of the knot

◆ setKnots_u() [1/2]

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setKnots_u ( const axom::Array< T > &  knots,
int  degree 
)
inline

Set the u knot vector by an axom::Array.

Parameters
[in]knotsThe new knot vector

◆ setKnots_v() [1/2]

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setKnots_v ( const axom::Array< T > &  knots,
int  degree 
)
inline

Set the v knot vector by an axom::Array.

Parameters
[in]knotsThe new knot vector

◆ setKnots_u() [2/2]

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setKnots_u ( const KnotVectorType knotVector)
inline

Set the u knot vector by a KnotVector object.

Parameters
[in]knotVectorThe new knot vector

◆ setKnots_v() [2/2]

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setKnots_v ( const KnotVectorType knotVector)
inline

Set the v knot vector by a KnotVector object.

Parameters
[in]knotVectorThe new knot vector

◆ getKnots_u() [1/2]

template<typename T , int NDIMS>
KnotVectorType& axom::primal::NURBSPatch< T, NDIMS >::getKnots_u ( )
inline

Return a reference to the KnotVector instance on the first axis.

◆ getKnots_u() [2/2]

template<typename T , int NDIMS>
const KnotVectorType& axom::primal::NURBSPatch< T, NDIMS >::getKnots_u ( ) const
inline

Return a const reference to the KnotVector instance on the first axis.

◆ getMinKnot_u()

template<typename T , int NDIMS>
T axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u ( ) const
inline

Get the minimum knot value in the u-axis.

References axom::primal::KnotVector< T >::getMinKnot().

◆ getMaxKnot_u()

template<typename T , int NDIMS>
T axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u ( ) const
inline

Get the maximum knot value in the u-axis.

References axom::primal::KnotVector< T >::getMaxKnot().

◆ getParameterSpaceDiagonal()

template<typename T , int NDIMS>
T axom::primal::NURBSPatch< T, NDIMS >::getParameterSpaceDiagonal ( ) const
inline

◆ getKnots_v() [1/2]

template<typename T , int NDIMS>
KnotVectorType& axom::primal::NURBSPatch< T, NDIMS >::getKnots_v ( )
inline

Return a reference to the KnotVector instance on the second axis.

◆ getKnots_v() [2/2]

template<typename T , int NDIMS>
const KnotVectorType& axom::primal::NURBSPatch< T, NDIMS >::getKnots_v ( ) const
inline

Return a const reference to the KnotVector instance on the second axis.

◆ getMinKnot_v()

template<typename T , int NDIMS>
T axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v ( ) const
inline

Get the minimum knot value in the v-axis.

References axom::primal::KnotVector< T >::getMinKnot().

◆ getMaxKnot_v()

template<typename T , int NDIMS>
T axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v ( ) const
inline

Get the maximum knot value in the v-axis.

References axom::primal::KnotVector< T >::getMaxKnot().

◆ getNumKnots_u()

template<typename T , int NDIMS>
int axom::primal::NURBSPatch< T, NDIMS >::getNumKnots_u ( ) const
inline

Return the length of the knot vector on the first axis.

References axom::primal::KnotVector< T >::getNumKnots().

◆ getNumKnots_v()

template<typename T , int NDIMS>
int axom::primal::NURBSPatch< T, NDIMS >::getNumKnots_v ( ) const
inline

Return the length of the knot vector on the second axis.

References axom::primal::KnotVector< T >::getNumKnots().

◆ insertKnot_u()

template<typename T , int NDIMS>
axom::IndexType axom::primal::NURBSPatch< T, NDIMS >::insertKnot_u ( u,
int  target_multiplicity = 1 
)
inline

Insert a knot to the u knot vector to have the given multiplicity.

Parameters
[in]uThe parameter value of the knot to insert
[in]target_multiplicityThe multiplicity of the knot to insert
Returns
The index of the new knot

Algorithm A5.3 on p. 155 of "The NURBS Book"

Note
If the knot is already present, it will be inserted up to the given multiplicity, or the maximum permitted by the degree
Precondition
Requires u in the span of the knots (up to a small tolerance)
Note
If u is outside the knot span up this tolerance, it is clamped to the span
Returns
The (maximum) index of the new knot

References axom::utilities::clampVal(), axom::primal::KnotVector< T >::findSpan(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_u(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_u(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_v(), axom::primal::KnotVector< T >::insertKnotBySpan(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::primal::NURBSPatch< T, NDIMS >::isValidParameter_u(), axom::utilities::min(), axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), SLIC_ASSERT, and SLIC_ASSERT_MSG.

◆ insertKnot_v()

template<typename T , int NDIMS>
axom::IndexType axom::primal::NURBSPatch< T, NDIMS >::insertKnot_v ( v,
int  target_multiplicity = 1 
)
inline

Insert a knot to the v knot vector to have the given multiplicity.

Parameters
[in]vThe parameter value of the knot to insert
[in]target_multiplicityThe multiplicity of the knot to insert
Returns
The index of the new knot

Algorithm A5.3 on p. 155 of "The NURBS Book"

Note
If the knot is already present, it will be inserted up to the given multiplicity, or the maximum permitted by the degree
Precondition
Requires v in the span of the knots (up to a small tolerance)
Note
If v is outside the knot span up this tolerance, it is clamped to the span
Returns
The (maximum) index of the new knot

References axom::utilities::clampVal(), axom::primal::KnotVector< T >::findSpan(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_v(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_u(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_v(), axom::primal::KnotVector< T >::insertKnotBySpan(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::primal::NURBSPatch< T, NDIMS >::isValidParameter_v(), axom::utilities::min(), axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), SLIC_ASSERT, and SLIC_ASSERT_MSG.

◆ reverseTrimmingCurves()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::reverseTrimmingCurves ( )
inline

Reverses all trimming curves in the patch.

Warning
Trimming curves should be oriented CCW in parameter space, and this method may make them CW

◆ reverseOrientation()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::reverseOrientation ( int  axis)
inline

Reverses the order of one direction of the NURBS patch's control points and weights.

This method does not affect the position of the patch in space, or its trimming curves, but it does reverse the patch's normal vectors.

Parameters
[in]axisorientation of patch. 0 to reverse in u, 1 for reverse in v

References axom::primal::NURBSPatch< T, NDIMS >::reverseOrientation_u(), and axom::primal::NURBSPatch< T, NDIMS >::reverseOrientation_v().

◆ reverseOrientation_u()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::reverseOrientation_u ( )
inline

◆ reverseOrientation_v()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::reverseOrientation_v ( )
inline

◆ swapAxes()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::swapAxes ( )
inline

Swap the axes such that s(u, v) becomes s(v, u)

This method does not affect the position of the patch in space, or its trimming curves.

References axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::ArrayBase< T, DIM, ArrayType >::shape(), and axom::utilities::swap().

◆ normalize()

◆ normalize_u()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::normalize_u ( )
inline

◆ normalize_v()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::normalize_v ( )
inline

◆ normalizeBySpan()

◆ rescale()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::rescale ( a,
b 
)
inline

Rescale both knot vectors to the span of [a, b].

Parameters
[in]aThe lower bound of the new knot vector
[in]bThe upper bound of the new knot vector
Precondition
Requires a < b

References axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), axom::primal::KnotVector< T >::rescale(), and SLIC_ASSERT.

◆ rescale_u()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::rescale_u ( a,
b 
)
inline

Rescale the knot vector in u to the span of [a, b].

Parameters
[in]aThe lower bound of the new knot vector
[in]bThe upper bound of the new knot vector
Precondition
Requires a < b

References axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::KnotVector< T >::rescale(), and SLIC_ASSERT.

◆ rescale_v()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::rescale_v ( a,
b 
)
inline

Rescale the knot vector in v to the span of [a, b].

Parameters
[in]aThe lower bound of the new knot vector
[in]bThe upper bound of the new knot vector
Precondition
Requires a < b

References axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), axom::primal::KnotVector< T >::rescale(), and SLIC_ASSERT.

◆ isValidParameter_u()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::isValidParameter_u ( u,
EPS = 1e-8 
) const
inline

Function to check if the u parameter is within the knot span.

References axom::primal::KnotVector< T >::getNumKnots().

◆ isValidParameter_v()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::isValidParameter_v ( v,
EPS = 1e-8 
) const
inline

Function to check if the v parameter is within the knot span.

References axom::primal::KnotVector< T >::getNumKnots().

◆ isValidInteriorParameter_u()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::isValidInteriorParameter_u ( t) const
inline

Checks if given u parameter is interior to the knot span.

References axom::primal::KnotVector< T >::isValidInteriorParameter().

◆ isValidInteriorParameter_v()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::isValidInteriorParameter_v ( t) const
inline

Checks if given v parameter is interior to the knot span.

References axom::primal::KnotVector< T >::isValidInteriorParameter().

◆ scaleParameterSpace()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::scaleParameterSpace ( double  scaleFactor,
bool  removeTrimmingCurves = false 
)
inline

Scale the parameter space of the NURBS patch geometry linearly (by tangents) in all directions.

Parameters
[in]scaleFactorThe multiplicative factor to expand each knot vector by
[in]removeTrimmingCurvesIf true, the resulting patch has no trimming curves

Algorithm from Wolters, Hans J., "Extensions: Extrapolation Methods for CAD", 1999

Note
This function only affects the geometry of the untrimmed NURBS patch, and does not affect any existing trimming curves (unless explicitly removed)
Postcondition
If removeTrimmingCurves is false, the resulting patch will be trimmed.
Warning
Method becomes numerically unstable for large values of scaleFactor, or for rational patches with a large range of weights.

References axom::primal::NURBSPatch< T, NDIMS >::expandParameterSpace(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), SLIC_ASSERT, and SLIC_WARNING_IF.

◆ expandParameterSpace()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::expandParameterSpace ( double  expansionAmount_u,
double  expansionAmount_v,
bool  removeTrimmingCurves = false 
)
inline

Expand the parameter space of the NURBS patch geometry linearly (by tangents) in all directions by a fixed amount.

Parameters
[in]expansionAmount_uThe absolute additive amount by which the u is expanded
[in]expansionAmount_vThe absolute additive amount by which the v is expanded
[in]removeTrimmingCurvesIf true, the resulting patch has no trimming curves

Algorithm from Wolters, Hans J., "Extensions: Extrapolation Methods for CAD", 1999

Note
This function only affects the geometry of the untrimmed NURBS patch, and does not affect any existing trimming curves (unless explicitly removed)
Postcondition
If removeTrimmingCurves is false, the resulting patch will be trimmed.
Warning
Method becomes numerically unstable for values of expansionAmount that are a large fraction of the existing parameter size, or for rational patches with a large range of weights.

References axom::primal::Vector< T, NDIMS >::array(), axom::Array< T, DIM, SPACE, StoragePolicy >::clear(), axom::Array< T, DIM, SPACE, StoragePolicy >::fill(), axom::primal::KnotVector< T >::getDegree(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_u(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_v(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_u(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_v(), axom::primal::KnotVector< T >::getNumKnots(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::primal::NURBSPatch< T, NDIMS >::isTrimmed(), axom::primal::NURBSPatch< T, NDIMS >::makeTriviallyTrimmed(), axom::utilities::min(), axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), SLIC_ASSERT, and SLIC_WARNING_IF.

◆ evaluate()

template<typename T , int NDIMS>
PointType axom::primal::NURBSPatch< T, NDIMS >::evaluate ( u,
v 
) const
inline

Evaluate the NURBS patch geometry at a particular parameter value t.

Parameters
[in]uThe parameter value on the first axis
[in]vThe parameter value on the second axis

Adapted from Algorithm A3.5 on page 103 of "The NURBS Book"

Precondition
Requires u, v in the span of each knot vector (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::KnotVector< T >::calculateBasisFunctionsBySpan(), axom::utilities::clampVal(), axom::primal::KnotVector< T >::findSpan(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_u(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_v(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), and axom::primal::Point< T, NDIMS >::zero().

◆ isocurve()

template<typename T , int NDIMS>
NURBSCurveType axom::primal::NURBSPatch< T, NDIMS >::isocurve ( uv,
int  axis 
) const
inline

Returns a NURBS patch isocurve for a fixed parameter value of u or v.

Parameters
[in]uvparameter value at which to construct the isocurve
[in]axisorientation of curve. 0 for fixed u, 1 for fixed v
Returns
c The isocurve C(v) = S(u, v) for fixed u or C(u) = S(u, v) for fixed v
Precondition
Requires uv be in the span of the relevant knot vector

References axom::primal::NURBSPatch< T, NDIMS >::isocurve_u(), axom::primal::NURBSPatch< T, NDIMS >::isocurve_v(), and SLIC_ASSERT.

◆ isocurve_u()

◆ isocurve_v()

◆ evaluateDerivatives()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::evaluateDerivatives ( u,
v,
int  d,
axom::Array< VectorType, 2 > &  ders,
bool  evalByTotalOrder = false 
) const
inline

Evaluate the NURBS patch geometry derivatives up to order d at parameter u, v.

Parameters
[in]uThe parameter value on the first axis
[in]vThe parameter value on the second axis
[in]ordThe maximum (or total) order of evaluated derivatives to evaluate
[out]dersA matrix of size d+1 x d+1 containing the derivatives
[in]evalByTotalOrderIf true, only evaluate derivatives up to total order d instead of maximum component order d.

ders[i][j] is the derivative of S with respect to u i times and v j times. For consistency, ders[0][0] contains the evaluation point stored as a vector

If evalByTotalOrder is false (default behavior), then ders contains for d == 3

   | eval ,  S_u  ,   S_uu |

ders = | S_v , S_uv , S_uuv | | S_vv , S_vvu , S_uuvv |

If evalByTotalOrder is true, then ders[i][j] == 0 whenever i + j > d

   | eval , S_u  , S_uu |

ders = | S_v , S_uv , 0 | | S_vv , 0 , 0 |

Implementation adapted from Algorithm A3.6 on p. 111 of "The NURBS Book". Rational derivatives from Algorithm A4.4 on p. 137 of "The NURBS Book".

Precondition
Requires u, v be in the span of the knots (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::utilities::binomialCoefficient(), axom::utilities::clampVal(), axom::primal::KnotVector< T >::derivativeBasisFunctionsBySpan(), axom::primal::NURBSPatch< T, NDIMS >::du(), axom::primal::NURBSPatch< T, NDIMS >::dv(), axom::Array< T, DIM, SPACE, StoragePolicy >::fill(), axom::primal::KnotVector< T >::findSpan(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_u(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_v(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::utilities::min(), axom::Array< T, DIM, SPACE, StoragePolicy >::resize(), and axom::primal::Point< T, NDIMS >::zero().

◆ evaluateFirstDerivatives()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::evaluateFirstDerivatives ( u,
v,
PointType eval,
VectorType Du,
VectorType Dv 
) const
inline

Evaluates all first derivatives of the NURBS patch geometry at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
[out]evalThe point value of the NURBS patch at (u, v)
[out]DuThe vector value of S_u(u, v)
[out]DvThe vector value of S_v(u, v)

Uses a specialized formula that is faster than evaluateDerivatives, avoiding repeated dynamic allocations in the generic implementation. This is performance-critical for 3D winding-number evaluations

Precondition
We require evaluation of the patch at u and v between 0 and 1

References axom::utilities::clampVal(), axom::primal::KnotVector< T >::derivativeBasisFunctionsBySpan(), axom::primal::NURBSPatch< T, NDIMS >::du(), axom::primal::NURBSPatch< T, NDIMS >::dv(), axom::primal::KnotVector< T >::findSpan(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_u(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_v(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::utilities::min(), and axom::Array< T, DIM, SPACE, StoragePolicy >::resize().

◆ evaluateLinearDerivatives()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::evaluateLinearDerivatives ( u,
v,
PointType eval,
VectorType Du,
VectorType Dv,
VectorType DuDv 
) const
inline

Evaluates all linear derivatives of the NURBS patch geometry at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
[out]evalThe point value of the NURBS patch at (u, v)
[out]DuThe vector value of S_u(u, v)
[out]DvThe vector value of S_v(u, v)
[out]DuDvThe vector value of S_uv(u, v) == S_vu(u, v)
Precondition
We require evaluation of the patch at u and v (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::NURBSPatch< T, NDIMS >::evaluateDerivatives().

◆ evaluateSecondDerivatives()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::evaluateSecondDerivatives ( u,
v,
PointType eval,
VectorType Du,
VectorType Dv,
VectorType DuDu,
VectorType DvDv,
VectorType DuDv 
) const
inline

Evaluates all second derivatives of the NURBS patch geometry at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
[out]evalThe point value of the NURBS patch at (u, v)
[out]DuThe vector value of S_u(u, v)
[out]DvThe vector value of S_v(u, v)
[out]DuDuThe vector value of S_uu(u, v)
[out]DvDvThe vector value of S_vv(u, v)
[out]DuDvThe vector value of S_uu(u, v)
Precondition
We require evaluation of the patch at u and v (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::NURBSPatch< T, NDIMS >::evaluateDerivatives().

◆ du()

template<typename T , int NDIMS>
VectorType axom::primal::NURBSPatch< T, NDIMS >::du ( u,
v 
) const
inline

Computes a tangent in u of the NURBS patch geometry at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
Precondition
We require evaluation of the patch at u and v (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::NURBSPatch< T, NDIMS >::evaluateDerivatives().

◆ dv()

template<typename T , int NDIMS>
VectorType axom::primal::NURBSPatch< T, NDIMS >::dv ( u,
v 
) const
inline

Computes a tangent in v of the NURBS patch geometry at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
Precondition
We require evaluation of the patch at u and v (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::NURBSPatch< T, NDIMS >::evaluateDerivatives().

◆ dudu()

template<typename T , int NDIMS>
VectorType axom::primal::NURBSPatch< T, NDIMS >::dudu ( u,
v 
) const
inline

Computes the second derivative in u of the NURBS patch geometry patch at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
Precondition
We require evaluation of the patch at u and v (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::NURBSPatch< T, NDIMS >::evaluateDerivatives().

◆ dvdv()

template<typename T , int NDIMS>
VectorType axom::primal::NURBSPatch< T, NDIMS >::dvdv ( u,
v 
) const
inline

Computes the second derivative in v of the NURBS patch geometry patch at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
Precondition
We require evaluation of the patch at u and v (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::NURBSPatch< T, NDIMS >::evaluateDerivatives().

◆ dudv()

template<typename T , int NDIMS>
VectorType axom::primal::NURBSPatch< T, NDIMS >::dudv ( u,
v 
) const
inline

Computes the mixed second derivative in u and v of the NURBS patch geometry at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
Precondition
We require evaluation of the patch at u and v (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::NURBSPatch< T, NDIMS >::evaluateDerivatives().

◆ dvdu()

template<typename T , int NDIMS>
VectorType axom::primal::NURBSPatch< T, NDIMS >::dvdu ( u,
v 
) const
inline

Computes the mixed second derivative in u and v of the NURBS patch geometry at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
Precondition
We require evaluation of the patch at u and v (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::NURBSPatch< T, NDIMS >::dudv().

◆ normal()

template<typename T , int NDIMS>
VectorType axom::primal::NURBSPatch< T, NDIMS >::normal ( u,
v 
) const
inline

Computes the normal vector to the NURBS patch geometry at (u, v)

Parameters
[in]uParameter value at which to evaluate on the first axis
[in]vParameter value at which to evaluate on the second axis
Precondition
We require evaluation of the patch at u and v (up to a small tolerance)
Note
If u/v is outside the knot span up this tolerance, it is clamped to the span

References axom::primal::Vector< T, NDIMS >::cross_product(), and axom::primal::NURBSPatch< T, NDIMS >::evaluateFirstDerivatives().

◆ calculateUntrimmedPatchNormal()

template<typename T , int NDIMS>
VectorType axom::primal::NURBSPatch< T, NDIMS >::calculateUntrimmedPatchNormal ( int  npts = 20) const
inline

Calculate the average normal for the (untrimmed) patch.

Parameters
[in]nptsThe number of quadrature nodes used in each component integral

Algorithm from "Mean normal vector to a surface bounded by Bézier curves" by Kenji Ueda, 1996

Projects the 4 boundary curves of the patch along each coordinate axis, then computes the 2D area of that projection to get the corresponding component of the average surface normal.

Evaluates the integral with Gauss-Legendre quadrature on each boundary curve

Returns
The calculated mean surface normal

References axom::primal::evaluate_area_integral(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_u(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_v(), axom::primal::NURBSPatch< T, NDIMS >::getKnots_u(), axom::primal::NURBSPatch< T, NDIMS >::getKnots_v(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_u(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_v(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), and SLIC_ASSERT.

◆ getTrimmingCurves() [1/2]

template<typename T , int NDIMS>
const TrimmingCurveVec& axom::primal::NURBSPatch< T, NDIMS >::getTrimmingCurves ( ) const
inline

Get array of trimming curves.

◆ getTrimmingCurves() [2/2]

template<typename T , int NDIMS>
TrimmingCurveVec& axom::primal::NURBSPatch< T, NDIMS >::getTrimmingCurves ( )
inline

Get mutable array of trimming curves.

◆ getTrimmingCurve()

template<typename T , int NDIMS>
const TrimmingCurveType& axom::primal::NURBSPatch< T, NDIMS >::getTrimmingCurve ( int  idx) const
inline

Get a trimming curve by index.

References axom::Array< T, DIM, SPACE, StoragePolicy >::size(), and SLIC_ASSERT.

◆ addTrimmingCurve()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::addTrimmingCurve ( const TrimmingCurveType curve)
inline

◆ addTrimmingCurves()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::addTrimmingCurves ( const TrimmingCurveVec curves)
inline

◆ setTrimmingCurves()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::setTrimmingCurves ( const TrimmingCurveVec curves)
inline

Set the array of trimming curves.

◆ clearTrimmingCurves()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::clearTrimmingCurves ( )
inline

Clear trimming curves, but DON'T mark as untrimmed.

References axom::Array< T, DIM, SPACE, StoragePolicy >::clear().

◆ getNumTrimmingCurves()

template<typename T , int NDIMS>
int axom::primal::NURBSPatch< T, NDIMS >::getNumTrimmingCurves ( ) const
inline

Get number of trimming curves.

References axom::Array< T, DIM, SPACE, StoragePolicy >::size().

◆ isTriviallyTrimmed()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::isTriviallyTrimmed ( double  tol = 1e-8) const
inline

Predicate to check if the patch is "trivially trimmed" in parameter space.

A patch is considered trivially trimmed if it is marked as trimmed and has exactly four trimming curves that form an axis-aligned rectangle on the patch's parametric boundary:

  • Two curves are horizontal (v=min_v and v=max_v) and have opposite directions
  • Two curves are vertical (u=min_u and u=max_u) and have opposite directions
  • Each curve is approximately linear in (u,v) space
  • Curve endpoints match the patch's (min_u,max_u,min_v,max_v) corner coordinates
Parameters
[in]tolThreshold for squared distance used by NURBSCurve::isLinear and for the axis-alignment check on the curve endpoints.

References axom::primal::NURBSPatch< T, NDIMS >::du(), axom::primal::NURBSPatch< T, NDIMS >::dv(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getNumTrimmingCurves(), and axom::primal::NURBSPatch< T, NDIMS >::isTrimmed().

◆ isTrimmed()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::isTrimmed ( ) const
inline

use boolean flag for trimmed-ness

◆ isInvisible()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::isInvisible ( ) const
inline

Predicate to check if the patch is entirely invisible due to trimming state.

A patch is considered invisible when it is marked as trimmed, but has no trimming curves. This represents a trimmed patch whose visible region is empty.

References axom::primal::NURBSPatch< T, NDIMS >::getNumTrimmingCurves(), and axom::primal::NURBSPatch< T, NDIMS >::isTrimmed().

◆ markAsTrimmed()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::markAsTrimmed ( )
inline

Mark as trimmed.

◆ makeUntrimmed()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::makeUntrimmed ( )
inline

Delete all trimming curves.

References axom::Array< T, DIM, SPACE, StoragePolicy >::clear().

◆ makeTriviallyTrimmed()

◆ isVisible()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::isVisible ( u,
v 
) const
inline

Check if a parameter point is visible on the NURBS patch via a trim test.

Parameters
[in]uThe parameter value on the first axis
[in]vThe parameter value on the second axis

Checks for containment of the parameter point in the collection of trimming curves via an even-odd rule.

If the collection of trimming curves does not form closed loops, then the (now fractional) generalized winding number is rounded to the nearest integer before the even-odd rule is applied.

References axom::primal::NURBSPatch< T, NDIMS >::isTrimmed(), and axom::primal::KnotVector< T >::isValidParameter().

◆ clip()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::clip ( min_u,
max_u,
min_v,
max_v,
bool  normalizeParameters = false 
)
inline

Restrict the edges of a NURBS surface to the given parameter values, if necessary.

Parameters
[in]min_uThe minimum value of the u parameter
[in]max_uThe maximum value of the u parameter
[in]min_vThe minimum value of the v parameter
[in]max_vThe maximum value of the v parameter
[in]normalizeIf true, normalize the patch to the range [0, 1]^2

If a given parameter is less/greater than the minimum/maximum knot value, then the patch is not changed along that direction and axis

Precondition
Requires that min_u < max_u and min_v < max_v
Postcondition
A patch whose parameter space is a subset of [min_u, max_u] x [min_v, max_v]
See also
NURBSPatch::split()

References axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMaxKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_u(), axom::primal::NURBSPatch< T, NDIMS >::getMinKnot_v(), axom::primal::NURBSPatch< T, NDIMS >::normalize(), SLIC_ASSERT, axom::primal::NURBSPatch< T, NDIMS >::split_u(), and axom::primal::NURBSPatch< T, NDIMS >::split_v().

◆ clipToCurves()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::clipToCurves ( double  padding = 1e-5)
inline

Clip the edges of a NURBS surface to the AABB of the trimming curves.

Parameters
[in]paddingThe amount to be left on each side of the AABB after clipping
See also
NURBSPatch::clip()

References axom::primal::BoundingBox< T, NDIMS >::addBox(), axom::primal::BoundingBox< T, NDIMS >::getMax(), and axom::primal::BoundingBox< T, NDIMS >::getMin().

◆ extractBezier()

template<typename T , int NDIMS>
axom::Array<BezierPatch<T, NDIMS> > axom::primal::NURBSPatch< T, NDIMS >::extractBezier ( ) const
inline

Splits the NURBS patch geometry (at each internal knot) into several Bezier patches.

If either degree_u or degree_v is zero, the resulting Bezier patches along that axis will be disconnected and order 0

This method ignores any trimming curves in the patch, and returns all extracted patches of the untrimmed patch.

Algorithm A5.7 on p. 177 of "The NURBS Book"

Returns
An array of Bezier patches ordered lexicographically (in v, then u)

References axom::primal::NURBSPatch< T, NDIMS >::getDegree_u(), axom::primal::NURBSPatch< T, NDIMS >::getDegree_v(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_u(), axom::primal::NURBSPatch< T, NDIMS >::getNumControlPoints_v(), axom::primal::KnotVector< T >::getNumKnotSpans(), axom::primal::NURBSPatch< T, NDIMS >::getWeight(), axom::primal::NURBSPatch< T, NDIMS >::isRational(), axom::utilities::max(), and axom::Array< T, DIM, SPACE, StoragePolicy >::size().

◆ extractTrimmedBezier()

template<typename T , int NDIMS>
axom::Array<NURBSPatch<T, NDIMS> > axom::primal::NURBSPatch< T, NDIMS >::extractTrimmedBezier ( ) const
inline

Splits the NURBS patch geometry into several one-span trimmed NURBS surfaces.

Returns
An array of Bezier patches ordered lexicographically (in v, then u)
See also
extractBezier

References axom::primal::NURBSPatch< T, NDIMS >::getKnots_u(), axom::primal::NURBSPatch< T, NDIMS >::getKnots_v(), axom::primal::KnotVector< T >::getUniqueKnots(), and axom::Array< T, DIM, SPACE, StoragePolicy >::size().

◆ calculateTrimmedPatchNormal()

template<typename T , int NDIMS>
VectorType axom::primal::NURBSPatch< T, NDIMS >::calculateTrimmedPatchNormal ( int  npts = 20,
bool  useBezierExtraction = true 
) const
inline

Calculate the average normal for the trimmed patch.

Parameters
[in]nptsThe number of quadrature nodes used in each component integral

Decomposes the NURBS surface into trimmed Bezier components (to ensure differentiability of the integrand) and evaluates the integral numerically on each component using trimming curves

Evaluates the integral with Gauss-Legendre quadrature on each boundary curve

Returns
The calculated mean surface normal

References axom::primal::evaluate_area_integral(), axom::primal::NURBSPatch< T, NDIMS >::extractTrimmedBezier(), axom::primal::NURBSPatch< T, NDIMS >::getTrimmingCurves(), axom::primal::NURBSPatch< T, NDIMS >::normal(), and SLIC_ASSERT.

◆ calculateTrimmedPatchNormalArea()

template<typename T , int NDIMS>
std::pair<VectorType, double> axom::primal::NURBSPatch< T, NDIMS >::calculateTrimmedPatchNormalArea ( int  npts = 20,
bool  useBezierExtraction = true 
) const
inline

Calculate the unsigned surface area and integrated normal for the trimmed patch.

Parameters
[in]nptsThe number of quadrature nodes used in each component integral
[in]useBezierExtractionSet whether to do Bezier extraction on input for exponential convergence of general NURBS input

Decomposes the surface into trimmed Bezier components and evaluates the area and integrated normal numerically using trimming curves. Avoids the redundant geometry processing of computing each value separately

References axom::primal::area(), axom::primal::evaluate_area_integral(), axom::primal::NURBSPatch< T, NDIMS >::normal(), and SLIC_ASSERT.

◆ cleanedTrimmedRepresentation()

template<typename T , int NDIMS>
NURBSPatch axom::primal::NURBSPatch< T, NDIMS >::cleanedTrimmedRepresentation ( bool  normalize_by_span = true,
bool  ensure_trimmed = true,
bool  clip_to_curves = true 
) const
inline

Return a "clean" trimmed representation suitable for moment and GWN calculation.

Parameters
[in]normalize_by_spanNormalize the patch parameter space
[in]ensure_trimmedIf the patch isn't trimmed, make trivially trimmed
[in]clip_to_curvesClip patch parameter space to AABB of the trimming curves

References axom::primal::NURBSPatch< T, NDIMS >::clipToCurves(), axom::primal::NURBSPatch< T, NDIMS >::getNumTrimmingCurves(), axom::primal::NURBSPatch< T, NDIMS >::isTrimmed(), axom::primal::NURBSPatch< T, NDIMS >::makeTriviallyTrimmed(), and axom::primal::NURBSPatch< T, NDIMS >::normalizeBySpan().

◆ calculateSurfaceMoments()

template<typename T , int NDIMS>
template<int ORDER, int NVALS = 4 + (ORDER >= 0 ? 3 : 0) + (ORDER >= 1 ? 9 : 0) + (ORDER >= 2 ? 27 : 0)>
primal::Vector<T, NVALS> axom::primal::NURBSPatch< T, NDIMS >::calculateSurfaceMoments ( int  npts = 20,
bool  useBezierExtraction = false 
) const
inline

Calculate the surface moments for the trimmed patch for GWN evaluation.

Parameters
[in]nptsThe number of quadrature nodes used in each component integral
[in]useBezierExtractionSet whether to do Bezier extraction on input for exponential convergence of general NURBS input

Decomposes the surface into trimmed Bezier components and evaluates up to second order moments without recomputing quadrature nodes, avoiding the redundant processing of evaluating each separately.

References axom::primal::Vector< T, NDIMS >::cross_product(), axom::primal::evaluate_area_integral(), and axom::primal::Vector< T, NDIMS >::norm().

◆ split()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::split ( u,
v,
NURBSPatch< T, NDIMS > &  p1,
NURBSPatch< T, NDIMS > &  p2,
NURBSPatch< T, NDIMS > &  p3,
NURBSPatch< T, NDIMS > &  p4,
bool  normalizeParameters = false 
) const
inline

Splits a NURBS surface into four NURBS patches.

Parameters
[in]uparameter value at which to bisect on the first axis
[in]vparameter value at which to bisect on the second axis
[out]p1First output NURBS patch
[out]p2Second output NURBS patch
[out]p3Third output NURBS patch
[out]p4Fourth output NURBS patch
[in]normalizeParametersIf true, normalize the knot span of the output patches to the range [0, 1]^2

v = v_max

| | | | p3 | p4 | | | | -----—(u,v)------— | | | | p1 | p2 | | | | -------------------— u = u_max u/v_min

Note
If u/v is not strictly interior to the knot span, will return an invalid NURBS for the invalid portion and the original surface for the rest
Returns
True if and only if the patch was split (i.e., u, v is in the knot span)

References axom::primal::NURBSPatch< T, NDIMS >::normalize(), and axom::primal::NURBSPatch< T, NDIMS >::split_v().

◆ split_u()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::split_u ( u,
NURBSPatch< T, NDIMS > &  p1,
NURBSPatch< T, NDIMS > &  p2,
bool  normalizeParameters = false 
) const
inline

Split the NURBS surface in two along the u direction.

Note
If u is not strictly interior to the knot span, will return an invalid NURBS for the invalid portion and the original surface for the rest
Returns
True if and only if the patch was split (i.e., u is in the knot span)

References axom::primal::NURBSPatch< T, NDIMS >::getTrimmingCurves(), and axom::primal::NURBSPatch< T, NDIMS >::normalize_u().

◆ split_v()

template<typename T , int NDIMS>
bool axom::primal::NURBSPatch< T, NDIMS >::split_v ( v,
NURBSPatch< T, NDIMS > &  p1,
NURBSPatch< T, NDIMS > &  p2,
bool  normalizeParameters = false 
) const
inline

Split the NURBS surface in two along the v direction.

Note
If v is not strictly interior to the knot span, will return an invalid NURBS for the invalid portion and the original surface for the rest
Returns
True if and only if the patch was split (i.e., v is in the knot span)

References axom::primal::NURBSPatch< T, NDIMS >::getTrimmingCurves(), and axom::primal::NURBSPatch< T, NDIMS >::normalize_v().

◆ nearBisectOnLongerAxis()

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::nearBisectOnLongerAxis ( NURBSPatch< T, NDIMS > &  p1,
NURBSPatch< T, NDIMS > &  p2 
) const
inline

Split the patch in the longer parametric direction, determined via maximum control-net polyline length.

We measure the maximum polyline length of the control net in each parametric direction,

  • u-length: max over segments ||P(i+1,j) - P(i,j)||
  • v-length: max over segments ||P(i,j+1) - P(i,j)|| and split along the direction with the larger maximum. This approximates the geometric stretch along that axis.
Note
This heuristic considers only control points (and ignores weights for rational patches).

References axom::utilities::max().

◆ diskSplit() [1/2]

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::diskSplit ( u,
v,
r,
NURBSPatch< T, NDIMS > &  the_disk,
NURBSPatch< T, NDIMS > &  the_rest,
bool  clipDisk = true 
) const
inline

For a disk of radius r and center (u, v), split a NURBS surface into the portion inside/outside.

Parameters
[in]uThe x-coordinate of the center of the disk
[in]vThe y-coordinate of the center of the disk
[in]rThe radius of the disk
[out]the_diskThe NURBS surface inside the disk
[out]the_restThe NURBS surface outside the disk
[in]clipDiskIf true, the returned disk is clipped to the disk boundary

◆ diskSplit() [2/2]

template<typename T , int NDIMS>
void axom::primal::NURBSPatch< T, NDIMS >::diskSplit ( u,
v,
r,
NURBSPatch< T, NDIMS > &  the_disk,
NURBSPatch< T, NDIMS > &  the_rest,
bool &  isDiskInside,
bool &  isDiskOutside,
bool  ignoreInteriorDisk,
double  disk_padding 
) const
inline

Split a NURBS surface into two by cutting out a disk of radius r centered at (u, v)

Parameters
[in]uThe x-coordinate of the center of the disk
[in]vThe y-coordinate of the center of the disk
[in]rThe radius of the disk
[out]the_diskThe NURBS surface inside the disk
[out]the_restThe NURBS surface outside the disk
[out]isDiskInsideTrue if the disk is entirely inside the trimming curves
[out]isDiskOutsideTrue if the disk is entirely outside the trimming curves
[in]ignoreInteriorDiskIf true, don't perform subdivision if disk is entirely inside the trimming curves
[in]disk_paddingHow much space is retained around the disk in each direction after clipping
Note
Function arguments suited for use in GWN evaluation

References axom::primal::NURBSPatch< T, NDIMS >::addTrimmingCurve(), axom::Array< T, DIM, SPACE, StoragePolicy >::begin(), axom::Array< T, DIM, SPACE, StoragePolicy >::clear(), axom::Array< T, DIM, SPACE, StoragePolicy >::end(), axom::primal::NURBSCurve< T, NDIMS >::getMaxKnot(), axom::primal::NURBSCurve< T, NDIMS >::getMinKnot(), axom::primal::NURBSPatch< T, NDIMS >::makeTriviallyTrimmed(), axom::primal::NURBSPatch< T, NDIMS >::markAsTrimmed(), axom::Array< T, DIM, SPACE, StoragePolicy >::push_back(), axom::primal::NURBSCurve< T, NDIMS >::reverseOrientation(), axom::Array< T, DIM, SPACE, StoragePolicy >::size(), axom::sort(), and axom::primal::NURBSCurve< T, NDIMS >::split().

◆ print()

template<typename T , int NDIMS>
std::ostream& axom::primal::NURBSPatch< T, NDIMS >::print ( std::ostream &  os) const
inline

Simple formatted print of a NURBS Patch instance.

Parameters
osThe output stream to write to
Returns
A reference to the modified ostream

Friends And Related Function Documentation

◆ operator==

template<typename T , int NDIMS>
bool operator== ( const NURBSPatch< T, NDIMS > &  lhs,
const NURBSPatch< T, NDIMS > &  rhs 
)
friend

Equality operator for NURBS patches.

Parameters
[in]lhsThe left-hand side NURBS patch
[in]rhsThe right-hand side NURBS patch
Returns
True if the two patches are equal, false otherwise

◆ operator!=

template<typename T , int NDIMS>
bool operator!= ( const NURBSPatch< T, NDIMS > &  lhs,
const NURBSPatch< T, NDIMS > &  rhs 
)
friend

Inequality operator for NURBS patches.

Parameters
[in]lhsThe left-hand side NURBS patch
[in]rhsThe right-hand side NURBS patch
Returns
True if the two patches are not equal, false otherwise

The documentation for this class was generated from the following file: