# bezier.surface module¶

Helper for Bézier Surfaces / Triangles.

class bezier.surface.Surface(nodes, degree, base_x=0.0, base_y=0.0, width=1.0, _copy=True)

Bases: bezier._base.Base

Represents a Bézier surface.

We define a Bézier triangle as a mapping from the unit simplex in 2D (i.e. the unit triangle) onto a surface in an arbitrary dimension. We use barycentric coordinates

$\lambda_1 = 1 - s - t, \lambda_2 = s, \lambda_3 = t$

for points in the unit triangle $$\left\{(s, t) \mid 0 \leq s, t, s + t \leq 1\right\}$$:

As with curves, using these weights we get convex combinations of points $$v_{i, j, k}$$ in some vector space:

$B\left(\lambda_1, \lambda_2, \lambda_3\right) = \sum_{i + j + k = d} \binom{d}{i \, j \, k} \lambda_1^i \lambda_2^j \lambda_3^k \cdot v_{i, j, k}$

Note

We assume the nodes are ordered from left-to-right and from bottom-to-top. So for example, the linear triangle:

(0,0,1)

(1,0,0)  (0,1,0)


is ordered as

$\left[\begin{array}{c c c} v_{1,0,0} & v_{0,1,0} & v_{0,0,1} \end{array}\right]^T$

(0,0,2)

(1,0,1)  (0,1,1)

(2,0,0)  (1,1,0)  (0,2,0)


is ordered as

$\left[\begin{array}{c c c c c c} v_{2,0,0} & v_{1,1,0} & v_{0,2,0} & v_{1,0,1} & v_{0,1,1} & v_{0,0,2} \end{array}\right]^T$

the cubic triangle:

(0,0,3)

(1,0,2)  (0,1,2)

(2,0,1)  (1,1,1)  (0,2,1)

(3,0,0)  (2,1,0)  (1,2,0)  (0,3,0)


is ordered as

$\left[\begin{array}{c c c c c c c c c c} v_{3,0,0} & v_{2,1,0} & v_{1,2,0} & v_{0,3,0} & v_{2,0,1} & v_{1,1,1} & v_{0,2,1} & v_{1,0,2} & v_{0,1,2} & v_{0,0,3} \end{array}\right]^T$

and so on.

The index formula

$j + \frac{k}{2} \left(2 (i + j) + k + 3\right)$

can be used to map a triple $$(i, j, k)$$ onto the corresponding linear index, but it is not particularly insightful or useful.

>>> import bezier
>>> nodes = np.asfortranarray([
...     [0.0  , 0.0  ],
...     [0.5  , 0.0  ],
...     [1.0  , 0.25 ],
...     [0.125, 0.5  ],
...     [0.375, 0.375],
...     [0.25 , 1.0  ],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> surface
<Surface (degree=2, dimension=2)>

Parameters: nodes (numpy.ndarray) – The nodes in the surface. The rows represent each node while the columns are the dimension of the ambient space. degree (int) – The degree of the surface. This is assumed to correctly correspond to the number of nodes. Use from_nodes() if the degree has not yet been computed. base_x (Optional [ float ]) – The $$x$$-coordinate of the base vertex of the sub-triangle that this surface represents. See width() for more info. base_y (Optional [ float ]) – The $$y$$-coordinate of the base vertex of the sub-triangle that this surface represents. See width() for more info. width (Optional [ float ]) – The width of the sub-triangle that this surface represents. See width() for more info. _copy (bool) – Flag indicating if the nodes should be copied before being stored. Defaults to True since callers may freely mutate nodes after passing in.
area

float – The area of the current surface.

Raises: NotImplementedError – If the area isn’t already cached.
base_x

float – The x-coordinate of the base vertex.

See width() for more detail.

base_y

float – The y-coordinate of the base vertex.

See width() for more detail.

degree

int – The degree of the current shape.

dimension

int – The dimension that the shape lives in.

For example, if the shape lives in $$\mathbf{R}^3$$, then the dimension is 3.

edges

The edges of the surface.

>>> nodes = np.asfortranarray([
...     [0.0   ,  0.0   ],
...     [0.5   , -0.1875],
...     [1.0   ,  0.0   ],
...     [0.1875,  0.5   ],
...     [0.625 ,  0.625 ],
...     [0.0   ,  1.0   ],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> edge1, _, _ = surface.edges
>>> edge1
<Curve (degree=2, dimension=2)>
>>> edge1.nodes
array([[ 0.  ,  0.    ],
[ 0.5 , -0.1875],
[ 1.  ,  0.    ]])

Returns: The edges of the surface. Tuple [ Curve, Curve, Curve ]
elevate()

Return a degree-elevated version of the current surface.

Does this by converting the current nodes $$\left\{v_{i, j, k}\right\}_{i + j + k = d}$$ to new nodes $$\left\{w_{i, j, k}\right\}_{i + j + k = d + 1}$$. Does so by re-writing

$E\left(\lambda_1, \lambda_2, \lambda_3\right) = \left(\lambda_1 + \lambda_2 + \lambda_3\right) B\left(\lambda_1, \lambda_2, \lambda_3\right) = \sum_{i + j + k = d + 1} \binom{d + 1}{i \, j \, k} \lambda_1^i \lambda_2^j \lambda_3^k \cdot w_{i, j, k}$

In this form, we must have

\begin{split}\begin{align*} \binom{d + 1}{i \, j \, k} \cdot w_{i, j, k} &= \binom{d}{i - 1 \, j \, k} \cdot v_{i - 1, j, k} + \binom{d}{i \, j - 1 \, k} \cdot v_{i, j - 1, k} + \binom{d}{i \, j \, k - 1} \cdot v_{i, j, k - 1} \\ \Longleftrightarrow (d + 1) \cdot w_{i, j, k} &= i \cdot v_{i - 1, j, k} + j \cdot v_{i, j - 1, k} + k \cdot v_{i, j, k - 1} \end{align*}\end{split}

where we define, for example, $$v_{i, j, k - 1} = 0$$ if $$k = 0$$.

>>> nodes = np.asfortranarray([
...     [0.0 , 0.0 ],
...     [1.5 , 0.0 ],
...     [3.0 , 0.0 ],
...     [0.75, 1.5 ],
...     [2.25, 2.25],
...     [0.0 , 3.0 ],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> elevated = surface.elevate()
>>> elevated
<Surface (degree=3, dimension=2)>
>>> elevated.nodes
array([[ 0.  , 0.  ],
[ 1.  , 0.  ],
[ 2.  , 0.  ],
[ 3.  , 0.  ],
[ 0.5 , 1.  ],
[ 1.5 , 1.25],
[ 2.5 , 1.5 ],
[ 0.5 , 2.  ],
[ 1.5 , 2.5 ],
[ 0.  , 3.  ]])

Returns: The degree-elevated surface. Surface
evaluate_barycentric(lambda1, lambda2, lambda3, _verify=True)

Compute a point on the surface.

Evaluates $$B\left(\lambda_1, \lambda_2, \lambda_3\right)$$.

>>> nodes = np.asfortranarray([
...     [0.0  , 0.0  ],
...     [0.5  , 0.0  ],
...     [1.0  , 0.25 ],
...     [0.125, 0.5  ],
...     [0.375, 0.375],
...     [0.25 , 1.0  ],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> point = surface.evaluate_barycentric(0.125, 0.125, 0.75)
>>> point
array([[ 0.265625 , 0.73046875]])


However, this can’t be used for points outside the reference triangle:

>>> surface.evaluate_barycentric(-0.25, 0.75, 0.5)
Traceback (most recent call last):
...
ValueError: ('Parameters must be positive', -0.25, 0.75, 0.5)


or for non-Barycentric coordinates;

>>> surface.evaluate_barycentric(0.25, 0.25, 0.25)
Traceback (most recent call last):
...
ValueError: ('Values do not sum to 1', 0.25, 0.25, 0.25)


However, these “invalid” inputs can be used if _verify is False.

>>> surface.evaluate_barycentric(-0.25, 0.75, 0.5, _verify=False)
array([[ 0.6875 , 0.546875]])
>>> surface.evaluate_barycentric(0.25, 0.25, 0.25, _verify=False)
array([[ 0.203125, 0.1875 ]])

Parameters: lambda1 (float) – Parameter along the reference triangle. lambda2 (float) – Parameter along the reference triangle. lambda3 (float) – Parameter along the reference triangle. _verify (Optional [ bool ]) – Indicates if the barycentric coordinates should be verified as summing to one and all non-negative (i.e. verified as barycentric). Can either be used to evaluate at points outside the domain, or to save time when the caller already knows the input is verified. Defaults to True. The point on the surface (as a two dimensional NumPy array with a single row). numpy.ndarray ValueError – If the weights are not valid barycentric coordinates, i.e. they don’t sum to 1. (Won’t raise if _verify=False.) ValueError – If some weights are negative. (Won’t raise if _verify=False.)
evaluate_barycentric_multi(param_vals, _verify=True)

Compute multiple points on the surface.

Assumes param_vals has three columns of Barycentric coordinates. See evaluate_barycentric() for more details on how each row of parameter values is evaluated.

>>> nodes = np.asfortranarray([
...     [ 0. , 0.  ],
...     [ 1. , 0.75],
...     [ 2. , 1.  ],
...     [-1.5, 1.  ],
...     [-0.5, 1.5 ],
...     [-3. , 2.  ],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> surface
<Surface (degree=2, dimension=2)>
>>> param_vals = np.asfortranarray([
...     [0.   , 0.25, 0.75 ],
...     [1.   , 0.  , 0.   ],
...     [0.25 , 0.5 , 0.25 ],
...     [0.375, 0.25, 0.375],
... ])
>>> points = surface.evaluate_barycentric_multi(param_vals)
>>> points
array([[-1.75  , 1.75    ],
[ 0.    , 0.      ],
[ 0.25  , 1.0625  ],
[-0.625 , 1.046875]])


Note

There is also a Fortran implementation of this function, which will be used if it can be built.

Parameters: param_vals (numpy.ndarray) – Array of parameter values (as a Nx3 array). _verify (Optional [ bool ]) – Indicates if the coordinates should be verified. See evaluate_barycentric(). Defaults to True. Will also double check that param_vals is the right shape. The points on the surface. numpy.ndarray ValueError – If param_vals is not a 2D array and _verify=True.
evaluate_cartesian(s, t, _verify=True)

Compute a point on the surface.

Evaluates $$B\left(1 - s - t, s, t\right)$$ by calling evaluate_barycentric():

This method acts as a (partial) inverse to locate().

>>> nodes = np.asfortranarray([
...     [0.0 , 0.0  ],
...     [0.5 , 0.5  ],
...     [1.0 , 0.625],
...     [0.0 , 0.5  ],
...     [0.5 , 0.5  ],
...     [0.25, 1.0  ],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> point = surface.evaluate_cartesian(0.125, 0.375)
>>> point
array([[ 0.16015625, 0.44726562]])
>>> surface.evaluate_barycentric(0.5, 0.125, 0.375)
array([[ 0.16015625, 0.44726562]])

Parameters: s (float) – Parameter along the reference triangle. t (float) – Parameter along the reference triangle. _verify (Optional [ bool ]) – Indicates if the coordinates should be verified inside of the reference triangle. Defaults to True. The point on the surface (as a two dimensional NumPy array). numpy.ndarray
evaluate_cartesian_multi(param_vals, _verify=True)

Compute multiple points on the surface.

Assumes param_vals has two columns of Cartesian coordinates. See evaluate_cartesian() for more details on how each row of parameter values is evaluated.

>>> nodes = np.asfortranarray([
...     [ 0.0, 0.0],
...     [ 2.0, 1.0],
...     [-3.0, 2.0],
... ])
>>> surface = bezier.Surface(nodes, degree=1)
>>> surface
<Surface (degree=1, dimension=2)>
>>> param_vals = np.asfortranarray([
...     [0.0  , 0.0  ],
...     [0.125, 0.625],
...     [0.5  , 0.5  ],
... ])
>>> points = surface.evaluate_cartesian_multi(param_vals)
>>> points
array([[ 0.   , 0.   ],
[-1.625, 1.375],
[-0.5  , 1.5  ]])


Note

There is also a Fortran implementation of this function, which will be used if it can be built.

Parameters: param_vals (numpy.ndarray) – Array of parameter values (as a Nx2 array). _verify (Optional [ bool ]) – Indicates if the coordinates should be verified. See evaluate_cartesian(). Defaults to True. Will also double check that param_vals is the right shape. The points on the surface. numpy.ndarray ValueError – If param_vals is not a 2D array and _verify=True.
classmethod from_nodes(nodes, base_x=0.0, base_y=0.0, width=1.0, _copy=True)

Create a Surface from nodes.

Computes the degree based on the shape of nodes.

Parameters: nodes (numpy.ndarray) – The nodes in the surface. The rows represent each node while the columns are the dimension of the ambient space. base_x (Optional [ float ]) – The $$x$$-coordinate of the base vertex of the sub-triangle that this surface represents. base_y (Optional [ float ]) – The $$y$$-coordinate of the base vertex of the sub-triangle that this surface represents. width (Optional [ float ]) – The width of the sub-triangle that this surface represents. _copy (bool) – Flag indicating if the nodes should be copied before being stored. Defaults to True since callers may freely mutate nodes after passing in. The constructed surface. Surface
intersect(other, strategy=<IntersectionStrategy.geometric: 'geometric'>, _verify=True)

Find the common intersection with another surface.

Parameters: other (Surface) – Other surface to intersect with. strategy (Optional [ IntersectionStrategy ]) – The intersection algorithm to use. Defaults to geometric. _verify (Optional [ bool ]) – Indicates if extra caution should be used to verify assumptions about the algorithm as it proceeds. Can be disabled to speed up execution time. Defaults to True. List of intersections (possibly empty). TypeError – If other is not a surface. NotImplementedError – If at least one of the surfaces isn’t two-dimensional.
is_valid

bool – Flag indicating if the surface is “valid”.

Here, “valid” means there are no self-intersections or singularities.

This checks if the Jacobian of the map from the reference triangle is nonzero. For example, a linear “surface” with collinear points is invalid:

>>> nodes = np.asfortranarray([
...     [0.0, 0.0],
...     [1.0, 1.0],
...     [2.0, 2.0],
... ])
>>> surface = bezier.Surface(nodes, degree=1)
>>> surface.is_valid
False


while a quadratic surface with one straight side:

>>> nodes = np.asfortranarray([
...     [ 0.0  , 0.0  ],
...     [ 0.5  , 0.125],
...     [ 1.0  , 0.0  ],
...     [-0.125, 0.5  ],
...     [ 0.5  , 0.5  ],
...     [ 0.0  , 1.0  ],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> surface.is_valid
True


though not all higher degree surfaces are valid:

>>> nodes = np.asfortranarray([
...     [1.0, 0.0],
...     [0.0, 0.0],
...     [1.0, 1.0],
...     [0.0, 0.0],
...     [0.0, 0.0],
...     [0.0, 1.0],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> surface.is_valid
False

locate(point, _verify=True)

Find a point on the current surface.

Solves for $$s$$ and $$t$$ in $$B(s, t) = p$$.

This method acts as a (partial) inverse to evaluate_cartesian().

Note

A unique solution is only guaranteed if the current surface is valid. This code assumes a valid surface, but doesn’t check.

>>> nodes = np.asfortranarray([
...     [0.0 ,  0.0 ],
...     [0.5 , -0.25],
...     [1.0 ,  0.0 ],
...     [0.25,  0.5 ],
...     [0.75,  0.75],
...     [0.0 ,  1.0 ],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> point = np.asfortranarray([[0.59375, 0.25]])
>>> s, t = surface.locate(point)
>>> s
0.5
>>> t
0.25

Parameters: point (numpy.ndarray) – A (1xD) point on the surface, where $$D$$ is the dimension of the surface. _verify (Optional [ bool ]) – Indicates if extra caution should be used to verify assumptions about the inputs. Can be disabled to speed up execution time. Defaults to True. The $$s$$ and $$t$$ values corresponding to point or None if the point is not on the surface. NotImplementedError – If the surface isn’t in $$\mathbf{R}^2$$. ValueError – If the dimension of the point doesn’t match the dimension of the current surface.
nodes

numpy.ndarray – The nodes that define the current shape.

plot(pts_per_edge, color=None, ax=None, with_nodes=False)

Plot the current surface.

Parameters: pts_per_edge (int) – Number of points to plot per edge. color (Optional [ Tuple [ float, float, float ] ]) – Color as RGB profile. ax (Optional [ matplotlib.artist.Artist ]) – matplotlib axis object to add plot to. with_nodes (Optional [ bool ]) – Determines if the control points should be added to the plot. Off by default. The axis containing the plot. This may be a newly created axis. matplotlib.artist.Artist NotImplementedError – If the surface’s dimension is not 2.
subdivide()

Split the surface into four sub-surfaces.

Does so by taking the unit triangle (i.e. the domain of the surface) and splitting it into four sub-triangles

Then the surface is re-parameterized via the map to / from the given sub-triangles and the unit triangle.

For example, when a degree two surface is subdivided:

>>> nodes = np.asfortranarray([
...     [-1.0 , 0.0 ],
...     [ 0.5 , 0.5 ],
...     [ 2.0 , 0.0 ],
...     [ 0.25, 1.75],
...     [ 2.0 , 3.0 ],
...     [ 0.0 , 4.0 ],
... ])
>>> surface = bezier.Surface(nodes, degree=2)
>>> _, sub_surface_b, _, _ = surface.subdivide()
>>> sub_surface_b
<Surface (degree=2, dimension=2, base=(0.5, 0.5), width=-0.5)>
>>> sub_surface_b.nodes
array([[ 1.5   , 2.5   ],
[ 0.6875, 2.3125],
[-0.125 , 1.875 ],
[ 1.1875, 1.3125],
[ 0.4375, 1.3125],
[ 0.5   , 0.25  ]])

Returns: The lower left, central, lower right and upper left sub-surfaces (in that order). Tuple [ Surface, Surface, Surface, Surface ]
width

float – The “width” of the parameterized triangle.

When re-parameterizing (e.g. via subdivide()) we specialize the surface from the unit triangle to some sub-triangle. After doing this, we re-parameterize so that that sub-triangle is treated like the unit triangle.

To track which sub-triangle we are in during the subdivision process, we use the coordinates of the base vertex as well as the “width” of each leg.

>>> surface.base_x, surface.base_y
(0.0, 0.0)
>>> surface.width
1.0


Upon subdivision, the width halves (and potentially changes sign) and the vertex moves to one of four points:

>>> _, sub_surface_b, sub_surface_c, _ = surface.subdivide()
>>> sub_surface_b.base_x, sub_surface_b.base_y
(0.5, 0.5)
>>> sub_surface_b.width
-0.5
>>> sub_surface_c.base_x, sub_surface_c.base_y
(0.5, 0.0)
>>> sub_surface_c.width
0.5