A multisided C2 B-spline patch over extraordinary vertices in quadrilateral meshes

We propose a generalised B-spline construction that extends uniform bicubic B-splines to multisided regions spanned over extraordinary vertices in quadrilateral meshes. We show how the structure of the generalised B´ezier patch introduced by V´arady et al. can be adjusted to work with B-spline basis functions. We create ribbon surfaces based on B-splines using special basis functions. The resulting multisided surfaces are C 2 continuous internally and connect with G 2 continuity to adjacent regular and other multisided B-splines patches. We visually assess the quality of these surfaces and compare them to Catmull-Clark limit surfaces on several challenging geometrical conﬁgurations.


Introduction
Constructing arbitrary topology surfaces based on B-splines has been a long-standing challenge in geometric design and computer graphics, and recently also in (isogeometric) analysis. Extraordinary vertices and faces introduce problems in the definition of the basis functions which are typically only well defined over regular control mesh regions. Many solutions have been proposed over time that try to create smooth surfaces over these extraordinary regions of the mesh, including subdivision surfaces.
Many of the commonly used subdivision surfaces generalise B-spline surfaces governed by quad-dominant meshes by creating a smooth surface by an iterative mesh refinement process. This process extends standard uniform B-spline subdivision masks/stencils also to faces with any number of sides (extraordinary faces) and vertices with valency other than four (extraordinary vertices). The process reproduces B-spline surfaces in regular regions, but there exists no such (closed-form) representation in extraordinary regions. It is possible to evaluate the surface at arbitrary parametric positions [1], but in many cases this is not efficient and/or restricts the topology of the control mesh. Moreover, subdivision surfaces are typically only G 1 continuous at extraordinary vertices and introduce other undesirable artefacts in the surrounding area, especially over high valency regions. Therefore, many subdivision techniques strive to increase the (visual) continuity and fidelity at said places of the surface. However, a simple technique that preserves both good shape and offers continuity higher than G 1 has not been found yet [2].
Another approach is to construct smooth tensor-product (macro-)patches, one per each quad incident with an extraordinary vertex, that connect smoothly to the surrounding regular regions and to each other. This includes for instance the construction of [3] and the recent series of articles by Karčiauskas and Peters; see e.g. [4,5] and the references cited therein. In contrast, multisided patches strive to cover extraordinary regions by a single patch.
The generalised Bézier patch [6] is a generalisation of tensor product Bézier patches to faces with an arbitrary number of sides. The patches are able to interpolate positions and derivatives and can be smoothly joined to other multisided and quadrilateral Bézier patches. Their structure is very simple as it is a surface defined as a combination of ribbons, partial Bézier surfaces. This simple structure of the generalised Bézier patch motivated us to find a way to adjust this structure to make it also applicable to the B-spline case over quad meshes of arbitrary manifold topology.
We propose a multisided B-spline patch that generalises uniform bicubic B-splines to faces with an arbitrary number of sides. The new construction is based on the generalised Bézier patch, in that it builds on many of the same elements: ribbon surfaces, generalised barycentric coordinates, and suitable blending functions.
The rest of the paper is structured as follows. In Section 2 we give an overview of related work. Then in Section 3 we recall both generalised barycentric coordinates and generalised Bézier patches. In Section 4 we introduce our multisided B-spline patch construction and afterwards, in Section 5, we show several interesting results. We discuss the proposed method in Section 6 and conclude the paper in Section 7.

Related Work
There exist many approaches that choose to use quadrilateral elements instead of multisided patches at extraordinary points. Catmull-Clark limit surfaces can be approximated using Gregory patches situated around extraordinary vertices [7] to efficiently create a global G 1 surface. For higher order continuity, there exist several options. These options typically vary in the number of elements, degree of the elements, and the visual quality obtained in extraordinary regions. For example, Loop & Schaefer [8] create n patches of degree bi-7, whereas Karčiuaskas & Peters [4] create patches of degree bi-6. We do not consider such constructions; we rather aim for a single multisided patch that smoothly connects to adjacent regular regions.
Manifold-based constructions [9,10] generalise B-splines of arbitrary degree to meshes with arbitrary topology based on the theory of manifolds of Grimm & Hughes [11]. A manifold on an irregular mesh is defined by associating a special chart to each face (or even edge and vertex), and these overlapping charts are combined using the partition of unity method, resulting in a surface of theoretically arbitrary smoothness. This approach is not suitable for completing B-spline complexes defined over regular regions without modifying them.
Loop & DeRose [12] propose a multisided generalisation of B-spline surfaces based on S-patches [13]. The C 1 data of surrounding regular regions is used to construct a G 1 S-patch of high degree that fits within bi-quadratic and bi-cubic B-spline surfaces. The approach suffers from the large number of control points that need to be determined for each S-patch, due to their high degree, and it does not directly apply at extraordinary vertices.
The Gregory corner interpolator patch [14] can be fitted into a regular B-spline mesh and can interpolate boundary curve positions and first and second derivatives by combining n bi-linearly parameterised surfaces that lie in the corners of the patches. The patches have been used in combination with several initial steps of Catmull-Clark subdivision to create arbitrary topology cubic B-spline surfaces [15], where the extraordinary regions are filled with Gregory patches. Another, similar approach [16] uses multisided Zheng & Ball patches [17] to create multisided patches over extraordinary regions that fit within bi-quadratic B-spline surfaces. However, these patches offer only a limited number of sides.
Multisided (generalised) Bézier patches [18], which are the most relevant construction for us, are recalled in detail in Section 3.2.

Preliminaries
We introduce several preliminary topics on which our new construction is built. As with most multisided patch approaches, we also use generalised barycentric coordinates to parametrise the patches, and so discuss these first.

Generalised Barycentric Coordinates
Generalised barycentric coordinates provide a coordinate system in which any point p on a planar polygon can be expressed as a weighted combination of the polygon's vertices. There are many types of these coordinates and an extensive overview of generalised barycentric coordinates can be found in [19]. We recall the main properties of the coordinates.
With a cyclic ordering of the vertices of a polygon v i , i = 1, . . . , n, the barycentric coordinate functions φ i , i = 1, . . . , n partition unity n i=1 φ i = 1, are non-negative φ i ≥ 0, ∀i on the polygon, and satisfy the barycentric property n i=1 φ i (p)v i = p. Linear interpolation is achieved on the boundary for smooth coordinates, and at the vertices the Lagrange property is satisfied φ i (v j ) = δ i j , where δ i j is the Kronecker delta. For triangles the coordinates reduce to the unique barycentric coordinates. For higher valencies this uniqueness property no longer holds.
To parametrise a polygon in R 3 we use a regular parametrisation domain, a regular polygon of valency n, which we parametrise by Wachspress coordinates [20,21]. Although many other options exist, Wachspress coordiantes provide, arguably, the simplest and most efficient choice due to their closed rational form.

Multisided Bézier Construction
The new multisided B-spline construction has a similar structure to the generalized Bézier patch [18]. The generalised Bézier patch is a multisided combination of n Bézier ribbons, which are partial Bézier surfaces. Each ribbon R i has local parameters constructed from the generalised barycentric coordinate functions corresponding to its end vertices v i−1 and v i . The local parameters are defined in terms of the generalised barycentric coordinates φ = [φ 1 , . . . , φ n ] as the side parameter These parameters locally behave like bi-linear coordinates. The side parameter s i increases linearly on side Γ i and is 0 on Γ i−1 and 1 on Γ i+1 . Similarly, the distance parameter h i is 0 on the whole of Γ i and increases linearly on sides Γ i−1 and Γ i+1 , and into the polygon. The side and distance parameters handle the parameterisation of the individual ribbons R i (s i , h i ).
Each patch is a combination of n ribbons that need to be blended together. For this purpose, special blending functions are attached to the control points that are common to more than one ribbon. We can decompose the blending functions into the left (α d i ) and right (β d i ) counterpart for each ribbon: Depending on the number of derivatives that need to be reproduced, the degree q of the terms in the blending functions is set accordingly. For G 1 continuity for now, we set q = 1. The blending functions assure that on each edge Γ i only the relevant basis functions of the associated ribbon contribute to the surface. Then the full definition of the ribbons is  For a G 1 construction, we need to ensure that each ribbon consists out of two rows of control points. This is also true for the ribbons on the adjacent sides Γ i−1 and Γ i . Therefore, in the corners many control points are shared between ribbons. These blocks of control points need to be blended using the blending functions.
The patch is defined to be the sum of the n ribbon surfaces, but this leads to a patch that is not necessarily affine invariant as the basis functions multiplied by the blending functions do not sum to unity. To create a patch that is affine invariant, there are two choices. Either the patch is normalised by the sum of the basis functions or a central control point C is created to which the deficient (or excess) weight is assigned where s i and h i stand for s i (φ) and h i (φ), respectively, for brevity. The control structure of the generalised Bézier patch is simple and the number of control points increases linearly with the number of sides (in contrast to the S-patch setting).
Although Bézier curves and surfaces have a strong relation with B-splines, a straightforward extension of the generalised Bézier patch to generalised B-spline patches is not easily made. This can be attributed to the difference in basis functions, topological restrictions, and mainly to required internal smoothness.

Multisided B-spline Construction
For the multisided B-spline construction we look for a generalisation of uniform cubic B-splines to multisided regions, in our case to the one-ring neighbourhoods of extraordinary vertices. These extraordinary vertices must be, as in the case of most similar constructions including [15], isolated, meaning that no two extraordinary vertices are adjacent to each other. However, to be able to process meshes which have extraordinary faces and (nonisolated) vertices, we use up to two Catmull-Clark subdivision steps until all extraordinary vertices have been correctly isolated.
For each extraordinary vertex with valency n, a single multisided patch is created, covering its one-ring neighbourhood. Like the generalised Bézier patch, each is a combination of n ribbons. In our case, the ribbons are composite B-spline surfaces with special basis functions. A schematic view of a multisided cubic B-spline surface is shown in Figure 2, left. The regular (bi-cubic) patches defined by the regular part of the control mesh leave a multisided hole for which each side is adjacent to two regular regions (so not a single Bézier patch). Therefore, the multisided face that we create must connect smoothly to both of these patches on each of its sides, and be itself smooth. For the creation of these special ribbons we first need to define appropriate basis functions.

Extended Cubic Basis Functions
Simply interchanging the Bernstein basis functions for Bspline basis functions leads to problems. As the ribbons are now a composite surface along the patch boundary, they must also be extended smoothly towards the inside of the patch. If this is not done, C 2 continuity will not be preserved towards the centre of the patch (along the spokes of the extraordinary vertex). Moreover, unlike the Bernstein polynomials, some of the standard Bspline basis functions have non-zero values and derivatives at the 'wrong' knots, as we now need to work with two knot intervals (over [0, 2] with three knots), and not just one interval [0, 1] as in the Bézier setting. Therefore, we need to create special basis functions that are C 2 everywhere inside the patch and reproduce the B-spline basis functions at certain values within the extended interval and vanish at others, as follows.
We define the functions E 3 j (u) over the interval [0, 2], to be 3  Secondly, we want to guarantee that the ribbon of side Γ i will not contribute on the distant sides Γ i+2 , · · · , Γ i−2 . This translates into the requirement that its positional contribution should vanish on the 'far side', i.e., E 3 j (2) = 0. Likewise, for G 1 continuity the first derivative should also vanish, i.e., E 3 j (u)| u=2 = 0. Finally, if G 2 continuity is desired, then the second derivative at u = 2 should vanish, too.
Unfortunately, merely extending the standard cubic basis functions to have an extended support only partly satisfies the conditions, as only some of the endpoint conditions are met. As the functions E j are fine on the interval [0, 1], we only need to find suitable extensions to [1,2] that assure the listed conditions for G 2 continuity. However, having this piece-wise representation of E j adds additional constraints. As we do not want to have lower than C 2 continuity at any point on the interval, care must be taken to create a C 2 continuous join at u = 1.
We now detail the construction of each of the three basis func-tions; see Figure 3. The first and second basis functions are simply E 3 0 (u) ≡ N 3 0 (u), and E 3 1 (u) ≡ N 3 1 (u) with u ∈ [0, 2]. As these functions both vanish at u = 2 up to second derivative, they satisfy all properties that we need to guarantee smoothness of the ribbons. The third basis function is a piece-wise function, but the second polynomial piece on the interval [1,2] is quintic: It is impossible to create a cubic function that reproduces up to the second derivative that of N 3 2 (u) at u = 1 and also has vanishing derivatives up to order two at u = 2. Therefore, a quintic function is needed to join smoothly to N 3 2 (u) on [1,2]. The quintic function connects to the cubic B-spline function with C 2 continuity and at the endpoint of the interval it vanishes along with its first and second derivative. The coefficients for the quintic function were obtained by solving the system where c i are the coefficients for a univariate quintic Bézier polynomial over [0, 1]. We only need to determine the first three coefficients c 0 , c 1 , c 2 , as the other three coefficients, c 3 , c 4 , c 5 , can simply be set to 0 to let all of the necessary derivatives vanish. Finally, the piece is shifted from the interval [0, 1] to [1,2].
Having defined the extended basis functions, we are now ready to define the B-spline ribbons.

B-spline Ribbons and Patch
The definition of the multisided B-spline patch is very similar to that of the generalised Bézier patch. The multisided B-spline is a combination of n B-spline ribbons. Again, the local ribbon parameters, s i and h i , are constructed from generalised barycentric coordinates and are used to parametrise each ribbon. However, this time they are scaled so that they span the interval [0, 2], ready to be used with the extended basis functions. Specifically, Along the direction ofs i normal cubic B-spline basis functions are used, and alongh i the extended basis functions E 3 j are employed. The distribution of the control point blending functions is also similar to the generalised Bézier patch. All control points of the ribbon require at least one blending function. The control points in the middle row of control points require a combination of two blending functions; see Figure 2, right. These control points need to be blended with respect to the adjacent B-spline ribbons on Γ i−1 and Γ i+1 . Both of these B-spline ribbons need at least three rows of control points to be able to join with curvature continuity to adjacent patches. This causes the blending functions of this row of control points to overlap as there are only 5 columns of control points in each ribbon. We can create suitable blending functions by simply multiplying the two respective blending functions together. This completes the B-spline ribbon definition: with v ∈ [0, 2]. Here, N 3 j are the uniform cubic B-spline basis functions and E 3 k are extended cubic basis functions as defined above. The distribution of the blending functions is governed by The ribbon blending functions are of the familiar form, but instead the functions α i and β i are cubed. This is different from the blending functions of the (enhanced) generalised Bézier patch [6], where squared terms provide already sufficient properties to obtain curvature continuity.
The cubic B-spline ribbon is a piece-wise surface, but is still C 2 continuous in the parametric direction along u and also at the join along u = 1. The blending functions introduce singularities in the corners of the patches, due to the denominator being zero there. The singularities can be easily removed as the evaluated position is simply the limit position of the control vertex at this corner. This position can be found by applying the regular limit stencil to the one-ring of vertices around the corner vertices, for example b i 11 and b i 31 , surrounding the extraordinary vertex.
Again, the ribbons can be combined into the final surface: Similarly to the generalised Bézier patch, this patch can be made affine invariant either by normalisation by dividing it by the sum of basis and blending functions: or by introduction of an additional point that takes care of the excess or deficient weight in the basis functions. The structure of the generalised B-spline patch is simple and requires only the vertices of the B-spline control net. This is in contrast to other methods like [16] that first need to extract the Bézier control points. The surfaces are affine invariant, preserve the convex hull property, and edits to the control mesh have local influence. The internal continuity of the patch is C 2 as the ribbons themselves are piece-wise surfaces. It remains to show that the patch joins smoothly with the surrounding bi-cubic patches in regular regions.

Continuity
To show that our patches can maintain continuity conditions with adjacent patches, either regular or multisided, we need to show that on each side Γ i only the ribbon R i contributes towards the final patch there. We do this by examining the contributions of the basis functions attached to the control points of the ribbons. By construction, the basis functions of the distant sides R i−2 , . . . , R i+2 do not contribute at all to R i due to the sodesigned extended cubic basis functions. The only two ribbons that might contribute are R i−1 and R i+1 . Let us examine the situation of control points in the left corner of ribbon R i . Here, a control point contribution is influenced by two basis functions: as expected. Here, all the blending functions β i−1 associated to all non-zero basis functions of R i−1 are zero, and thus so is the positional contribution of R i−1 on Γ i . This shows that only R i contributes to the position of the patch on Γ i , proving C 0 continuity.
We now show that also G 1 continuity is maintained. We differentiate the involved basis functions in an arbitrary direction not along s i ; this derivative is for simplicity denoted by ∂. For a control point in the left corner of R i we have where we make use of the convention ∂(B(s, 0)) = ∂(B(s, h))| h=0 . Since we are evaluating on Γ i , we have that both β i−1 and ∂(β i−1 ) vanish. In addition ∂(α i ) = 0. Thus we are left with ∂(B 3 j−1,k+1 (s i ,h i )), which means that the derivative across Γ i is determined solely by the derivatives of the basis functions of R i , as we needed to show for G 1 continuity.
Finally, the curvature continuity of the patches is established by showing that the second derivative of the patch on Γ i is solely determined by R i . Again, the construction of the extended basis functions allows us to say that they do not influence the curvature on the distant sides as their second derivatives vanish there. What remains to be shown is that only the second derivative of the basis functions of R i contributes to the curvature of the patch at Γ i . On this side the following conditions hold for the blending functions: β i−1 = ∂(β i−1 ) = ∂ 2 (β i−1 ) = 0. In addition we know that ∂(α i ) = ∂ 2 (α i ) = 0. Using this we can show for a control point on the left side of R i that as required. This shows that also the second derivative on Γ i is solely determined by the differential terms of R i and our patch reproduces all derivatives of the B-spline ribbon there up to order 2. Similar derivations hold for the control points in the right corner of the ribbon. This is also true for the centre column, where the two blending functions are combined; the chain rule applied to their product again ensures that the correct terms vanish at Γ i . This completes the G 2 proof.

Basis Functions
All basis functions of the generalised B-spline patches do not sum to unity and either weight deficiency or weight excess occurs towards the centre of the patch. This was also observed in the generalised Bézier patch [6]. Although this does not compromise the continuity and the affine invariance of the patch, as it can be restored by either normalisation or using an additional central control point to assign the deficient/excess weight to (see Section 3.2), it does however influence the shape of the resulting surface.
In these areas the surface decreases in curvature and tends to show flat spots. This becomes especially prevalent in the triangular case and in the octagonal case. We did not investigate patches of higher valencies as these rarely occur in practice, but as seen in Figure 4 it can be shown that increasing the valency will increase the deficient weight towards the centre of the patch.

Results
In this section we examine the visual quality of the proposed multisided construction, focusing especially on challenging configurations.

Challenging Cases
In Figure 5, we show several challenging geometries with the resulting multisided B-spline patches surrounded by one ring of regular bi-cubic patches. Each one of the geometries requires a single Catmull-Clark subdivision step to be able to generate a multisided patch. As can be seen from the reflection line renderings, the multisided patches join smoothly to the adjacent regular patches. Internally the patches are smooth and show no unexpected changes in the reflection lines. Naturally, when comparing to the Catmull-Clark limit surface there are differences in the non-regular regions. The Catmull-Clark limit surfaces show pinching artefacts around the extraordinary vertices.

Arbitrary Topology Meshes
The patches can be easily integrated into surfaces with multiple extraordinary vertices and faces. Figure 6 shows several such examples. Our patches span over the regions with extraordinary vertices. As can be seen, when multisided patches are adjacent to one another they are still able to join smoothly (with G 2 continuity). Figure 5: Several different mesh geometries with an extraordinary vertex for which a multisided generalised B-spline patch is generated. From top to bottom: Control mesh, shaded multisided B-spline patch and surrounding regular patches, shaded Catmull-Clark limit surface, reflection lines on the multisided B-spline patch and surrounding regular patches, reflection lines on the Catmull-Clark limit surface, and mean curvature of the multisided B-spline patch (from red via green to blue). Observe that our construction does not suffer from pinching artefacts near the extraordinary points.

Weight Deficiency
The effect of the weight deficiency/excess is particularly noticeable in high valency patches oriented in a 'parabolic' manner. This can be seen in the bottom row of Figure 6, where the octagonal patch in the centre of the surface appears to flatten off.
Also the triangular case is affected by this problem. In Figure 7, the effect of weight deficiency is shown on the resulting surface. Here the normalised version of the patches S N (φ) are used. Despite having a nice smooth transition along the patch boundaries, a flat spot is visible towards the centre of the patches. It is possible to change the centre of the surface by using S C (φ), the version of the patch which adds an additional control point to the surface. This takes careful tuning however, and a general best position for this central control point is not readily available, although the Catmull-Clark limit position of the extraordinary vertex may provide a useful hint. In the results presented in this paper, we take the average p of the one-ring neighbourhood of the extraordinary vertex e, and then displace the extraordinary vertex along the difference. Namely, we set C := e + 3 4 (e − p).

Discussion
Our extended basis functions E 3 i (u) present just one of the possible choices for these functions. As mentioned in Section 4.1, as long as the functions vanish at the end of the interval (u = 2) and they connect with sufficient smoothness at u = 1, they are valid functions. A piece-wise cubic representation of E 3 2 certainly exists, but that would require more than two pieces and thus introduce an extra (unwanted) knot. We have chosen for the representation that stays the closest to the original B-spline basis functions, because of their good shape properties. Furthermore, by tuning the functions carefully, it might be possible to minimise the weight deficiency/excess towards the centre of the patches, whilst still meeting the aforementioned continuity conditions.
As shown in Section 5.3, the weight deficiency/excess creates noticeable flat spots in high valency patches. We have not investigated means to reduce this shape defect. Varady et al. [6] suggest to distribute the missing weight on the inner ring of control points, but this only works properly when the inner ring of control points does not contribute to the continuity conditions at the boundary. In contrast, in our B-spline based construction each control point is needed in at least two separate ribbons at a time. Therefore, the weight can only be assigned to an additional central control point, or potentially several of these. As shown in Section 5.3, it is possible to influence the shape at the centre of the patch. With this it is possible to interpolate the limit position of the extraordinary vertex as obtained through standard Catmull-Clark subdivision, if so desired.
We use a regular parameterisation domain to retrieve generalised barycentric coordinates efficiently so that we can evaluate the patches using GPU tessellation [22]. We do not follow the exact same tessellation strategies featured in [22], but instead we choose to subdivide the polygonal domain similarly to a Catmull-Clark refinement step. This is because this way we can easily match tessellation levels with adjacent regular quads to achieve crack-free tessellation and thus artefact-free renderings. As the multisided B-spline patches have a piece-wise ribbon representation, it is easier this way to ensure equal level of detail along the edges of patches, avoiding creating gaps or overlaps in the tessellation process.
This proves to be a very efficient means to evaluate the multisided patches; we can easily reach interactive rates of performance. For instance the mesh featured in the third row of Figure 6 can be rendered in 3 milliseconds using an NVIDIA Titan V graphics card at the maximum tessellation level of 64 × 64. Providing equivalent level of detail using ordinary Catmull-Clark subdivision would require 9 levels of subdivision.
Our patch currently only works in combination with Catmull-Clark subdivision and tensor-product uniform cubic B-splines. However, it is definitely possible to vary the degree of the patches such that they can be used in combination with other bi-degree tensor product B-spline surfaces, or potentially also box-splines (such as in Loop subdivision surfaces), and also for the case of extraordinary faces. For each degree d, special extended basis functions E d j (u) need to be found, now potentially over more knot spans, such that their properties still guarantee the desired smoothness, typically C d−1 . With increased degree comes increased separation of extraordinary vertices and also the number of piece-wise surfaces in each B-spline ribbon. This requires a separate and careful investigation.

Conclusion
We have presented a multisided B-spline construction based on the generalised Bézier patch. The new multisided patches are constructed out of B-spline ribbons that use mostly the regular B-spline basis functions, but need an extra adjustment. These special basis functions are required to be able to guarantee G 2 smoothness with adjacent patches and C 2 smoothness within the patches. The proposed structure of the patches is simple and they can be efficiently rendered using modern graphics hardware.
In our future work, we will focus on generalising our multisided patches to other degrees (especially two, four and five) and levels of continuity so that they can be used in combination with other B-spline or box-spline based subdivision schemes other than that of Catmull and Clark.