B-splines

Functions and classes for B-spline basis functions.

class pyiga.bspline.BSplineFunc(kvs, coeffs)

Any function that is given in terms of a tensor product B-spline basis with coefficients.

Parameters:
  • kvs (seq) – tuple of d KnotVector.
  • coeffs (ndarray) – coefficient array

kvs represents a tensor product B-spline basis, where the i-th KnotVector describes the B-spline basis in the i-th coordinate direction.

coeffs is the array of coefficients with respect to this tensor product basis. The length of its first d axes must match the number of degrees of freedom in the corresponding KnotVector. Trailing axes, if any, determine the output dimension of the function. If there are no trailing dimensions or only a single one of length 1, the function is scalar-valued.

For convenience, if coeffs is a vector, it is reshaped to the proper size for the tensor product basis. The result is a scalar-valued function.

kvs

the knot vectors representing the tensor product basis

Type:seq
coeffs

the coefficients for the function or geometry

Type:ndarray
sdim

dimension of the parameter domain

Type:int
dim

dimension of the output of the function

Type:int
apply_matrix(A)

Apply a matrix to each control point of this function 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.

as_nurbs()

Return a NURBS version of this function with constant weights equal to 1.

as_vector()

Convert a scalar function to a 1D vector function.

boundary(bdspec)

Return one side of the boundary as a BSplineFunc.

Parameters:bdspec – the side of the boundary to return; see compute_dirichlet_bc()
Returns:BSplineFunc – 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)
copy()

Return a copy of this geometry.

cylinderize(z0=0.0, z1=1.0, support=(0.0, 1.0))

Create a patch with one additional space dimension by linearly extruding along the new axis from z0 to z1.

By default, the new knot vector will be defined over the interval (0, 1). A different interval can be specified through the support parameter.

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.

grid_eval(gridaxes)

Evaluate the function on a tensor product grid.

Parameters:gridaxes (seq) – list of 1D vectors describing the tensor product grid.

Note

The gridaxes should be given in reverse order, i.e., the x axis last.

Returns:ndarray – array of function values; shape corresponds to input grid.
grid_hessian(gridaxes)

Evaluate the Hessian matrix of a scalar or vector function on a tensor product grid.

Parameters:gridaxes (seq) – list of 1D vectors describing the tensor product grid.

Note

The gridaxes should be given in reverse order, i.e., the x axis last.

Returns:ndarray – array of the components of the Hessian; symmetric part only, linearized. I.e., in 2D, it contains per grid point a 3-component vector corresponding to the derivatives (d_xx, d_xy, d_yy), and in 3D, a 6-component vector with the derivatives (d_xx, d_xy, d_xz, d_yy, d_yz, d_zz). If the input function is vector-valued, one such Hessian vector is computed per component of the function.

Thus, the output is an array of shape grid_shape x num_comp x num_hess, where grid_shape is the shape of the tensor grid described by the gridaxes, num_comp is the number of components of the function, and num_hess is the number of entries in the symmetric part of the Hessian as described above. The axis corresponding to num_comp is elided if the input function is scalar.

grid_jacobian(gridaxes)

Evaluate the Jacobian on a tensor product grid.

Parameters:gridaxes (seq) – list of 1D vectors describing the tensor product grid.

Note

The gridaxes should be given in reverse order, i.e., the x axis last.

Returns:ndarray – array of Jacobians (dim × sdim); shape corresponds to input grid. For scalar functions, the output is a vector of length sdim (the gradient) per grid point.
is_scalar()

Returns True if the function is scalar-valued.

is_vector()

Returns True if the function is vector-valued.

perturb(noise)

Create a copy of this function where all coefficients are randomly perturbed by noise of the given magnitude.

pointwise_eval(points)

Evaluate the B-spline 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 B-spline 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.
rotate_2d(angle)

Rotate a geometry with dim = 2 by the given angle and return the result.

scale(factor)

Scale all control points either by a scalar factor or componentwise by a vector 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.

transformed_jacobian(geo)

Create a function which evaluates the physical (transformed) gradient of the current function after a geometry transform.

translate(offset)

Return a version of this geometry translated by the specified offset.

class pyiga.bspline.KnotVector(knots, p)

Represents an open B-spline knot vector together with a spline degree.

Parameters:
  • knots (ndarray) – the 1D knot vector. Should be an open knot vector, i.e., the first and last knot should be repeated p+1 times. Interior knots may be single or repeated up to p times.
  • p (int) – the spline degree.

This class is commonly used to represent the B-spline basis over the given knot vector with the given spline degree. The B-splines are normalized in the sense that they satisfy a partition of unity property.

Tensor product B-spline bases are typically represented simply as tuples of the univariate knot spans. E.g., (kv1, kv2) would represent the tensor product basis formed from the two B-spline bases over the KnotVector instances kv1 and kv2.

A more convenient way to create knot vectors is the make_knots() function.

kv

vector of knots

Type:ndarray
p

spline degree

Type:int
copy()

Return a copy of this knot vector.

findspan(u)

Returns an index i such that kv[i] <= u < kv[i+1] (except for the boundary, where u <= kv[m-p] is allowed) and p <= i < len(kv) - 1 - p

first_active(k)

Index of first active basis function in interval (kv[k], kv[k+1])

first_active_at(u)

Index of first active basis function in the interval which contains u.

greville()

Compute Gréville abscissae for this knot vector

mesh

Return the mesh, i.e., the vector of unique knots in the knot vector.

mesh_span_indices()

Return an array of indices i such that kv[i] != kv[i+1], i.e., the indices of the nonempty spans. Return value has length self.numspans.

mesh_support_idx(j)

Return the first and last mesh index of the support of i

mesh_support_idx_all()

Compute an integer array of size N × 2, where N = self.numdofs, which contains for each B-spline the result of mesh_support_idx().

meshsize_avg()

Compute average length of the knot spans of this knot vector

numdofs

Number of basis functions in a B-spline basis defined over this knot vector

numspans

Number of nontrivial intervals in the knot vector

refine(new_knots=None)

Return the refinement of this knot vector by inserting new_knots, or performing uniform refinement if none are given.

support(j=None)

Support of the knot vector or, if j is passed, of the j-th B-spline

support_idx(j)

Knot indices of support of j-th B-spline

class pyiga.bspline.PhysicalGradientFunc(func, geo)

A class for function objects which evaluate physical (transformed) gradients of scalar functions with geometry transforms.

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.bspline.active_ev(knotvec, u)

Evaluate all active B-spline basis functions at the points u.

Returns an array of shape (p+1, u.size) if u is an array.

pyiga.bspline.collocation(kv, nodes)

Compute collocation matrix for B-spline basis at the given interpolation nodes.

Parameters:
  • kv (KnotVector) – the B-spline knot vector
  • nodes (array) – array of nodes at which to evaluate the B-splines
Returns:

A Scipy CSR matrix with shape (len(nodes), kv.numdofs) whose entry at (i,j) is the value of the j-th B-spline evaluated at nodes[i].

pyiga.bspline.collocation_derivs(kv, nodes, derivs=1)

Compute collocation matrix and derivative collocation matrices for B-spline basis at the given interpolation nodes.

Returns a list of derivs+1 sparse CSR matrices with shape (nodes.size, kv.numdofs).

pyiga.bspline.collocation_derivs_info(kv, nodes, derivs=1)

Similar to collocation_info(), but the second array also contains coefficients for computing derivatives up to order derivs. It has shape (derivs + 1) x len(nodes) x (p + 1).

Corresponds to a row-wise representation of the matrices computed by collocation_derivs().

pyiga.bspline.collocation_info(kv, nodes)

Return two arrays: one containing the index of the first active B-spline per evaluation node, and one containing, per node, the coefficient vector of length p+1 for the linear combination of basis functions which yields the point evaluation at that node.

Corresponds to a row-wise representation of the collocation matrix (see collocation()).

pyiga.bspline.deriv(knotvec, coeffs, deriv, u)

Evaluate a derivative of the spline with given B-spline coefficients at all points u.

Parameters:
  • knotvec (KnotVector) – B-spline knot vector
  • coeffs (ndarray) – 1D array of coefficients, length knotvec.numdofs
  • deriv (int) – which derivative to evaluate
  • u (ndarray) – 1D array of evaluation points
Returns:

ndarray – the vector of function derivatives

pyiga.bspline.ev(knotvec, coeffs, u)

Evaluate a spline with given B-spline coefficients at all points u.

Parameters:
  • knotvec (KnotVector) – B-spline knot vector
  • coeffs (ndarray) – 1D array of coefficients, length knotvec.numdofs
  • u (ndarray) – 1D array of evaluation points
Returns:

ndarray – the vector of function values

pyiga.bspline.interpolate(kv, func, nodes=None)

Interpolate function in B-spline basis at given nodes (or Gréville abscissae by default)

pyiga.bspline.knot_insertion(kv, u)

Return a sparse matrix of size (n+1) x n, with n = kv.numdofs´, which maps coefficients from `kv to a new knot vector obtained by inserting the new knot u into kv.

pyiga.bspline.load_vector(kv, f)

Compute the load vector (L_2 inner products of basis functions with f).

pyiga.bspline.make_knots(p, a, b, n, mult=1)

Create an open knot vector of degree p over an interval (a,b) with n knot spans.

This automatically repeats the first and last knots p+1 times in order to create an open knot vector. Interior knots are single by default, i.e., have maximum continuity.

Parameters:
  • p (int) – the spline degree
  • a (float) – the starting point of the interval
  • b (float) – the end point of the interval
  • n (int) – the number of knot spans to divide the interval into
  • mult (int) – the multiplicity of interior knots
Returns:

KnotVector – the new knot vector

pyiga.bspline.numdofs(kvs)

Convenience function which returns the number of dofs in a single knot vector or in a tensor product space represented by a tuple of knot vectors.

pyiga.bspline.project_L2(kv, f)

Compute the B-spline coefficients for the L_2-projection of f.

pyiga.bspline.prolongation(kv1, kv2)

Compute prolongation matrix between B-spline bases.

Given two B-spline bases, where the first spans a subspace of the second one, compute the matrix which maps spline coefficients from the first basis to the coefficients of the same function in the second basis.

Parameters:
  • kv1 (KnotVector) – source B-spline basis knot vector
  • kv2 (KnotVector) – target B-spline basis knot vector
Returns:

csr_matrix – sparse matrix which prolongs coefficients from kv1 to kv2

pyiga.bspline.single_ev(knotvec, i, u)

Evaluate i-th B-spline at all points u

pyiga.bspline.tp_bsp_eval_pointwise(kvs, coeffs, points)

Evaluate a tensor-product B-spline function at an unstructured list of points.

Parameters:
  • kvs – tuple of KnotVector instances representing a tensor-product B-spline basis
  • coeffs (ndarray) – coefficient array; see BSplineFunc for details
  • points – an array or sequence such that points[i] is an array containing the coordinates for dimension i, where i = 0, …, len(kvs) - 1 (in xyz order). All arrays must have the same shape.
Returns:

An ndarray containing the function values of the spline function at the points.

pyiga.bspline.tp_bsp_eval_with_jac_pointwise(kvs, coeffs, points)

Evaluate the values and Jacobians of a tensor-product B-spline function at an unstructured list of points.

Parameters:
  • kvs – tuple of KnotVector instances representing a tensor-product B-spline basis
  • coeffs (ndarray) – coefficient array; see BSplineFunc for details
  • points – an array or sequence such that points[i] is an array containing the coordinates for dimension i, where i = 0, …, len(kvs) - 1 (in xyz order). All arrays must have the same shape.
Returns:

A pair of `ndarray`s – one for the values and one for the Jacobians.

pyiga.bspline.tp_bsp_jac_pointwise(kvs, coeffs, points)

Evaluate the Jacobian of a tensor-product B-spline function at an unstructured list of points.

Parameters:
  • kvs – tuple of KnotVector instances representing a tensor-product B-spline basis
  • coeffs (ndarray) – coefficient array; see BSplineFunc for details
  • points – an array or sequence such that points[i] is an array containing the coordinates for dimension i, where i = 0, …, len(kvs) - 1 (in xyz order). All arrays must have the same shape.
Returns:

An ndarray containing the Jacobians of the spline function at the points.