Geometry¶
Classes and functions for creating and manipulating tensor product B-spline and NURBS patches.
See Geometry manipulation in pyiga for examples on how to create custom geometries.
-
class
pyiga.geometry.
NurbsFunc
(kvs, coeffs, weights, premultiplied=False)¶ Any function that is given in terms of a tensor product NURBS basis with coefficients and weights.
Parameters: - kvs (seq) – tuple of d
KnotVector
s. - coeffs (ndarray) – coefficient array; see
BSplineFunc
for format. The constructor may modify coeffs during premultiplication! - weights (ndarray) – coefficients for weight function in the same format as coeffs. If weights=None is passed, the weights are assumed to be given as the last vector component of coeffs instead.
- premultiplied (bool) – pass True if the coefficients are already premultiplied by the weights.
-
kvs
¶ the knot vectors representing the tensor product basis
Type: seq
-
coeffs
¶ the premultiplied coefficients for the function, including the weights in the last vector component
Type: ndarray
-
sdim
¶ dimension of the parameter domain
Type: int
-
dim
¶ dimension of the output of the function
Type: int
The evaluation functions have the same prototypes and behavior as those in
BSplineFunc
.-
apply_matrix
(A)¶ Apply a matrix to each control point of this function, leave the weights unchanged, and return the result.
A should either be a single matrix or an array of matrices, one for each control point. Standard numpy broadcasting rules apply.
-
boundary
(bdspec)¶ Return one side of the boundary as a
NurbsFunc
.Parameters: bdspec – the side of the boundary to return; see compute_dirichlet_bc()
Returns: NurbsFunc
– representation of the boundary side; has sdim reduced by 1 and the same dim as this function
-
bounding_box
(grid=1)¶ Compute a bounding box for the image of this geometry.
By default, only the corners are taken into account. By choosing grid > 1, a finer grid can be used (for non-convex geometries).
Returns: a tuple of (lower,upper) limits per dimension (in XY order)
-
coeffs_weights
()¶ Return the non-premultiplied coefficients and weights as a pair of arrays.
-
copy
()¶ Return a copy of this geometry.
-
eval
(*x)¶ Evaluate the function at a single point of the parameter domain.
Parameters: *x – the point at which to evaluate the function, in xyz order
-
find_inverse
(x, tol=1e-08)¶ Find the coordinates in the parameter domain which correspond to the physical point x.
-
is_scalar
()¶ Returns True if the function is scalar-valued.
-
is_vector
()¶ Returns True if the function is vector-valued.
-
pointwise_eval
(points)¶ Evaluate the NURBS function at an unstructured list of points.
Parameters: points – an array or sequence such that points[i] is an array containing the coordinates for dimension i, where i = 0, …, sdim - 1 (in xyz order). All arrays must have the same shape. Returns: An ndarray containing the function values at the points.
-
pointwise_jacobian
(points)¶ Evaluate the Jacobian of the NURBS function at an unstructured list of points.
Parameters: points – an array or sequence such that points[i] is an array containing the coordinates for dimension i, where i = 0, …, sdim - 1 (in xyz order). All arrays must have the same shape. Returns: An ndarray containing the Jacobian matrices at the points, i.e., a matrix of size dim x sdim per evaluation point.
-
scale
(factor)¶ Scale all control points either by a scalar factor or componentwise by a vector, leave the weights unchanged, and return the resulting new function.
-
support
¶ Return a sequence of pairs (lower,upper), one per source dimension, which describe the extent of the support in the parameter space.
-
translate
(offset)¶ Return a version of this geometry translated by the specified offset.
- kvs (seq) – tuple of d
-
class
pyiga.geometry.
UserFunction
(f, support, dim=None, jac=None)¶ A function (supporting the same basic protocol as
BSplineFunc
) which is given in terms of a user-defined callable.Parameters: - f (callable) – a function of d variables; may be scalar or vector-valued
- support – a sequence of d pairs of the form (lower,upper) describing
the support of the function (see
BSplineFunc.support
) - dim (int) – the dimension of the function output; by default, is automatically determined by calling f
- jac (callable) – optionally, a function evaluating the Jacobian matrix of the function
The
sdim
attribute (seeBSplineFunc.sdim
) is determined from the length of support.-
boundary
(bdspec)¶ Return one side of the boundary as a
UserFunction
.Parameters: bdspec – the side of the boundary to return; see compute_dirichlet_bc()
Returns: UserFunction
– representation of the boundary side; has sdim reduced by 1 and the same dim as this function
-
bounding_box
(grid=1)¶ Compute a bounding box for the image of this geometry.
By default, only the corners are taken into account. By choosing grid > 1, a finer grid can be used (for non-convex geometries).
Returns: a tuple of (lower,upper) limits per dimension (in XY order)
-
find_inverse
(x, tol=1e-08)¶ Find the coordinates in the parameter domain which correspond to the physical point x.
-
is_scalar
()¶ Returns True if the function is scalar-valued.
-
is_vector
()¶ Returns True if the function is vector-valued.
-
pyiga.geometry.
bspline_quarter_annulus
(r1=1.0, r2=2.0)¶ A B-spline approximation of a quarter annulus in the first quadrant.
Parameters: - r1 (float) – inner radius
- r2 (float) – outer radius
Returns: BSplineFunc
2D geometry
-
pyiga.geometry.
circle
(r=1.0)¶ Construct a circle with radius r using NURBS.
-
pyiga.geometry.
circular_arc
(alpha, r=1.0)¶ Construct a circular arc with angle alpha and radius r.
The arc is centered at the origin, starts on the positive x axis and travels in counterclockwise direction.
-
pyiga.geometry.
circular_arc_3pt
(alpha, r=1.0)¶ Construct a circular arc with angle alpha and radius r using 3 control points.
The angle alpha must be between 0 and pi.
-
pyiga.geometry.
circular_arc_5pt
(alpha, r=1.0)¶ Construct a circular arc with angle alpha and radius r using 5 control points.
-
pyiga.geometry.
circular_arc_7pt
(alpha, r=1.0)¶ Construct a circular arc with angle alpha and radius r using 7 control points.
-
pyiga.geometry.
disk
(r=1.0)¶ A NURBS representation of a circular disk.
The parametrization has four boundary singularities where the determinant of the Jacobian becomes 0, at the bottom, top, left and right edges.
Parameters: r (float) – radius Returns: NurbsFunc
2D geometry
-
pyiga.geometry.
identity
(extents)¶ Identity mapping (using linear splines) over a d-dimensional box given by extents as a list of (min,max) pairs or of
KnotVector
.Returns: BSplineFunc
geometry
-
pyiga.geometry.
line_segment
(x0, x1, support=(0.0, 1.0), intervals=1)¶ Return a
BSplineFunc
which describes the line between the vectors x0 and x1.If specified, support describes the interval in which the function is supported; by default, it is the interval (0,1).
If specified, intervals is the number of intervals in the underlying linear spline space. By default, the minimal spline space with 2 dofs is used.
-
pyiga.geometry.
outer_product
(G1, G2)¶ Compute the outer product of two
BSplineFunc
orNurbsFunc
geometries. This means that given two input functions\[G_1(y), G_2(x),\]it returns a new function
\[G(x,y) = G_1(y) G_2(x),\]where the multiplication is componentwise in the case of vector functions. The returned function is a
NurbsFunc
if either input is and aBSplineFunc
otherwise. It has source dimension (BSplineFunc.sdim
) equal to the sum of the source dimensions of the input functions.G1 and G2 should have the same image dimension (
BSplineFunc.dim
), and the output then has the same as well. However, broadcasting according to standard Numpy rules is permissible; e.g., one function can be vector-valued and the other scalar-valued.The coefficients of the result are the pointwise products of the coefficients of the input functions over a new tensor product spline space.
-
pyiga.geometry.
outer_sum
(G1, G2)¶ Compute the outer sum of two
BSplineFunc
orNurbsFunc
geometries. This means that given two input functions\[G_1(y), G_2(x),\]it returns a new function
\[G(x,y) = G_1(y) + G_2(x).\]The returned function is a
NurbsFunc
if either input is and aBSplineFunc
otherwise. It has source dimension (BSplineFunc.sdim
) equal to the sum of the source dimensions of the input functions.G1 and G2 should have the same image dimension (
BSplineFunc.dim
), and the output then has the same as well. However, broadcasting according to standard Numpy rules is permissible; e.g., one function can be vector-valued and the other scalar-valued.The coefficients of the result are the pointwise sums of the coefficients of the input functions over a new tensor product spline space.
-
pyiga.geometry.
perturbed_square
(num_intervals=5, noise=0.02)¶ Randomly perturbed unit square.
Unit square with given number of intervals per direction; the control points are perturbed randomly according to the given noise level.
Returns: BSplineFunc
2D geometry
-
pyiga.geometry.
quarter_annulus
(r1=1.0, r2=2.0)¶ A NURBS representation of a quarter annulus in the first quadrant. The ‘bottom’ and ‘top’ boundaries (with respect to the reference domain) lie on the x and y axis, respectively.
Parameters: - r1 (float) – inner radius
- r2 (float) – outer radius
Returns: NurbsFunc
2D geometry
-
pyiga.geometry.
semicircle
(r=1.0)¶ Construct a semicircle in the upper half-plane with radius r using NURBS.
-
pyiga.geometry.
tensor_product
(G1, G2, *Gs)¶ Compute the tensor product of two or more
BSplineFunc
orNurbsFunc
functions. This means that given two input functions\[G_1(y), G_2(x),\]it returns a new function
\[G(x,y) = G_2(x) \times G_1(y),\]where \(\times\) means that vectors are joined together. The resulting
BSplineFunc
orNurbsFunc
has source dimension (BSplineFunc.sdim
) equal to the sum of the source dimensions of the input functions, and target dimension (BSplineFunc.dim
) equal to the sum of the target dimensions of the input functions.
-
pyiga.geometry.
twisted_box
()¶ A 3D volume that resembles a box with its right face twisted and bent upwards.
Corresponds to gismo data file twistedFlatQuarterAnnulus.xml.
Returns: BSplineFunc
3D geometry
-
pyiga.geometry.
unit_cube
(dim=3, num_intervals=1)¶ The dim-dimensional unit cube with num_intervals intervals per coordinate direction.
Returns: BSplineFunc
geometry
-
pyiga.geometry.
unit_square
(num_intervals=1)¶ Unit square with given number of intervals per direction.
Returns: BSplineFunc
2D geometry