Show the the Dimension of the Space S of Continuous Piecewise Quadratic Polynomials

Piecewise Polynomial

Spline Basics

Carl de Boor , in Handbook of Computer Aided Geometric Design, 2002

6.1 PIECEWISE POLYNOMIALS

A piecewise polynomial of order k with break sequence ξ (necessarily strictly increasing) is, by definition, any function f that, on each of the half-open intervals [ξ _j ‥ ξ _j+1), agrees with some polynomial of degree < k. The term 'order' used here is not standard but handy.

Note that this definition makes a piecewise polynomial function right-continuous, meaning that, for any x, f(x) = f(x+):= lim _h↓0 f(x + h). This choice is arbitrary, but has become standard. Keep in mind that, at its break ξ _j , the piecewise polynomial function f has, in effect, two values, namely its limit from the left, f(ξ _j -), and its limit from the right, f(ξ _j +) = f(ξ _j ).

The set of all piecewise polynomial functions of order k with break sequence ξ is denoted here

${\displaystyle \prod {}_{<k,\xi .}}$

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780444511041500071

B-Spline Approximation and the de Boor Algorithm

Ron Goldman , in Pyramid Algorithms, 2003

7.7.1 Elementary Properties of the B-Spline Basis Functions

In this section we shall study the elementary properties of the B-spline basis functions. These characteristics of the B-splines {N_{k, n} (t)} both mirror and are mirrored in the elementary features of B-spline curves derived in Section 7.4. Below we list these features and then derive each of these properties from the corresponding properties of B-spline curves.

1.

Piecewise polynomial

2.

Continuity of order C^n−μ at knots of multiplicity μ

3.

Compact support

4.

Partition of unity

5.

Nonnegativity

6.

Spline basis

7.

Unimodality

▪

Piecewise polynomial. From the up recurrence we know that the B-spline basis functions are B-spline curves. Therefore, the B-splines must be piecewise polynomials.

▪

Continuity. Again since by construction the B-spline basis functions are B-spline curves, the B-splines must have continuity of order n − μ at knots of multiplicity μ.

▪

Compact support. By the de Boor algorithm, the only B-splines that are nonzero over the parameter interval [t_k , t _k+1] are N_{k-n, n} (t), …, N_{k, n} (t). Hence the B-spline N_{k, n} (t) is nonzero only for values of t in the parameter interval [t_k , t _k+n+1]—that is, support{N_{k, n} (t)} = [t_k , t _k+n+1]. Therefore, from now on, whenever we want to make explicit the knots on which N_{k, n} (t) depends, we shall write N_{k, n} (t|t_k , …, t _k+n+1). By Equation (7.7), the compact support of the B-splines is equivalent to the local control property for the control points of B-spline curves.

▪

Partition of unity. The B-splines form a partition of unity. This result can be proved from the down recurrence (Equation (7.8)) by induction on n. This property can also be derived from the de Boor algorithm by setting P_k = 1 for all k and observing that since at every stage of the algorithm we are taking affine combinations of the nodes the value at every interior node is also equal to one. Hence the value at any apex must be one. Therefore,

$1=S(t)={{\displaystyle \sum }}_k{N}_{k,n}(t)P_k={\Sigma }_k{N}_{k,n}(t).$

The partition of unity property of the B-spline basis functions is equivalent to the affine invariance of B-spline curves.

▪

Nonnegativity. Recall that for any parameter interval the labels along the edges of the de Boor algorithm are nonnegative. Since the B-spline N_{k, n} (t) represents the sum over all paths from the kth position at the base to the various apexes of the de Boor triangles, it follows that the B-splines too are nonnegative. The partition of unity and nonnegativity of the B-spline basis functions are equivalent to the affine invariance and the convex hull properties of B-spline curves.

▪

Spline basis. To prove that the B-splines {N_{k, n} (t)} with knots {t_j } form a basis for the space of all splines S(t) with knots {t_j }, we need to show that the B-splines span this space and are linearly independent. But by Theorem 7.3, every spline is a B-spline curve; that is, every spline S(t) with knots {t_j } can be generated from the de Boor algorithm for some set of control points {P_k }. Therefore, by Equation (7.7), $S(t)={\Sigma }_k{N}_{k,n}(t)P_k$ , so the B-splines {N_{k, n} (t)} do indeed span the space of all splines with knots {t_j }. To prove that the B-splines are linearly independent, we must show that if ${\Sigma }_k{c}_k{N}_{k,n}(t)\equiv 0$ , then C_k = 0 for all k. Let's restrict our attention to the parameter interval [t_i , t _i+1]. Over this interval N_{i-n, n} (t), …, N_{i, n} (t) are the only nonzero B-splines, so over this interval

${\Sigma }_k{c}_k{N}_{k,n}(t)=\underset{h=0}{\overset{n}{\Sigma }}{c}_{i-n+h}{N}_{i-n+h,n}(t)\cdot$

Moreover, over the interval [t_i , t _i+1], the B-splines N_{i-n, n} (t), …, N_{i, n} (t) are polynomials, and by Section 7.1 these polynomials are just the progressive basis functions $b_0^n(t),\mathrm{\dots },b_n^n(t)$ , which form a polynomial basis. Therefore, if

$0\equiv {\displaystyle \sum _{h=0}^n{c}_{i-n+h}{N}_{i-n+h,n}(t)}={\displaystyle \sum _{h=0}^n{c}_{i-n+h}{b}_h^n(t),}$

then C _i-n+h = 0 for all h. Hence the B-splines are indeed linearly independent. The linear independence of the B-splines is equivalent to the nondegeneracy of B-spline curves.

▪

Unimodality. Recall that a function is said to be unimodal if it has only one local maximum. The B-splines {N_{k, n} (t)} are unimodal in t. To understand why, consider the graph of the function N_{k, n} (t)—that is, the curve S(t) = (t, N_{k, n} (t)). The function F(t) = t is a polynomial and hence certainly a spline (see Section 7.5, Exercise 2). Since the B-splines form a basis for the space of all splines, there must be constants {c_j } such that $t={{\displaystyle \sum }}_j{c}_j{N}_{j,n}(t)$ . (We shall derive explicit expressions for the constants {c_j } in Section 7.7.2, but for now all we need to know is that such constants exist.) Therefore, S(t) = Σ _j (C_j , δ_{j, k} )N_{j, n} (t). Thus the control points for the graph of N_{k, n} (t) all lie along the t-axis except for the point at (C_k , 1). Therefore, the control polygon for the graph of N_{k, n} (t) has only one local maximum (see Figure 7.31). But the graph of N_{k, n} (t) is a B-spline curve, and by Theorem 7.4, B-spline curves are variation diminishing. Therefore, the graph of N_{k, n} (t) can oscillate no more than its control polygon. Hence N_{k, n} (t) has only one local maximum.

Figure 7.31. Graph of the cubic B-spline N _{0, 3}(t) (light) together with its control polygon (dark).

Exercises

1.

Let S(t) be a spline of degree n with knots {t_k } whose support lies in [t_n , t _k+n+1].

a.

Prove that there is a constant c such that S(t) = cN_{k, n} (t).

b.

Conclude from part (a) that the B-splines {N_{k, n} (t)} have minimal support. That is, show that if S(t) is a spline of degree n with knots {t_k } whose support lies in a closed subinterval of [t_n , t _k+n+1], then S(t) is identically zero.

2.

Prove that Sign $Alternations\{{\Sigma }_k{c}_k{N}_{k,n}(t)\}\le Sign$ Alternations {c_k }.

3.

Let τ _k+i = at _k+i + b, i = 0, …, n + 1, for some fixed constants a > 0 and b. Show that

$N_k^n(at+b|{\tau }_k,\dots ,{\tau }_{k+n+1})=N_k^n(t|t_k,\dots ,t_{k+n+1}).$

Compare this result to Section 7.4, Exercise 2.

4.

Show, by example, that the B-splines {N_{k, n} (t)} are not necessarily unimodal in k for n ≥ 6.

5.

The B-splines {N_{k, n} (t)} are called the de Boor normalized B-splines. The Schoenberg normalized B-splines {M_{k, n} (t)} are defined by setting

$M_{k,n}(t)=\frac{N_{k,n}(t)}{t_{k+n+1}-t_k}=\frac{N_{k,n}(t)}{Support}$

a.

Prove that $M_{k,n}(t)=\frac{t-t_k}{t_{k+n+1}-t_k}{M}_{k,n-1}(t)+\frac{t_{k+n+1}-t}{t_{k+n+1}-t_k}{M}_{k+1,n-1(t)}$

b.

Conclude from part (a) that the Schoenberg normalized B-splines {M_{k, n} (t)} are unimodal in k.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9781558603547500088

Splines on Surfaces

Marian Neamtu , in Handbook of Computer Aided Geometric Design, 2002

9.2 SCALAR SPLINES ON SMOOTH SURFACES

Having looked at examples of splines on specific manifolds, let us now consider the case of general surfaces. We will first define the context of our setting. In this section the surface Swill be C ^∞-smooth in the usual sense of differential geometry. This is because our main interest will be in smooth splines on S, which would be an awkward requirement if Sitself were not sufficiently smooth. Of course, for all practical purposes the infinite differentiability condition can be relaxed to match the order of smoothness of f. Also, we will restrict ourselves here to the case of scalar-valued functions f(i.e., k = 1 in (9.2)).

Splines are usually obtained by dividing up the domain of their definition, S in our case, into disjoint subsets. The most familiar means of partitioning a given domain is to triangulate it. If the domain is planar, the resulting spline space is the well-known space of piecewise polynomials on triangulations. While there exist a great variety of bivariate splines, corresponding to many types of grids and polygonal partitions, piecewise polynomials on planar triangulations are the most universal since all other types can be viewed and represented as splines on triangulations. Therefore, to address the problem on general surfaces S, we will also assume that the sought-for spline space will correspond to a triangulation of S.

As usual, a triangulation of Sis a collection Δ of geodesic triangles on S, whose interiors are disjoint and such that the union of all triangles in Δ is S, i.e., S = ∪ _T∈Δ T. Here, a geodesic triangle T is a closed subset of S, homeomorphic to a planar triangle, whose boundary consists of three geodesic segments on Sconnecting a triple of points in T, called the vertices of T. These three geodesics are the edges of T. The vertices of T will be denoted by V(T) and the edges by E(T). The edges will be indexed by V(T) such that e _v ∈ E(T) will denote the edge of T opposite to v∈ V(T).

Perhaps the most popular way of representing splines on triangulations is by means of the Bernstein-Bézier formalism. Bernstein-Bézier techniques are extremely useful tools for constructing piecewise polynomial surfaces over triangulated planar domains. As is demonstrated in many chapters of the book (see chapter 4), they play an important role in CAGD, data fitting, computer vision, and elsewhere (see [22],[39]). To define analogs of Bernstein-Bézier techniques associated with S, we need to find finite-dimensional spaces of functions that play the role of ordinary polynomials in the plane. This in turn leads us to the question of defining an appropriate counterpart of the space of linear functions, since then higher-degree "polynomials" will be obtained as products of the "linear" functions. To be able to express such generalized polynomials in a Bernstein-Bézier-like form, in every triangle T ∈ Δ, we also need to find analogs of the well-known barycentric coordinates. Such generalized coordinates will be denoted as

$b_v^T:T\to R,v\in V\left(T\right),T\in \Delta .$

This reduces the problem of an appropriate definition of "piecewise polynomials" on Δ to finding a reasonable way to define the functions b ^T v. To answer this query, let us list some of the main properties associated with the classical Bernstein-Bézier formalism:

(1)

Non-negativity and partition of unity of Bernstein polynomials associated with any given triangle;

(2)

The convex hull property;

(3)

Reproduction of algebraic polynomials;

(4)

Affine invariance i.e., Bézier coefficients will not change after transforming the Euclidean plane by an affine transformation;

(5)

The possibility of obtaining arbitrarily smooth piecewise polynomials on any triangulation;

(6)

Smoothness between adjacent Bézier triangles can be expressed in terms of local Bézier coefficients, corresponding to these triangles.

While this may seem somewhat surprising, in the context of splines on general surfaces the most fundamental of the above items turn out to be properties (5) and (6). That is, spline spaces on general surfaces should be such that they can be used to build functions that are arbitrarily globally smooth (provided the degree is sufficiently high) and yet these spaces should be flexible enough so that it is possible to construct local methods of reconstruction. Another important property of the classical piecewise polynomials, that should be maintained in the general case, is that they are generated in a simple way from linear functions on each triangle. In our setting, the analog of the space of linear functions will be the three-dimensional space L ^T , which is the linear span of the barycentric coordinates, i.e., L ^T := span{b ^T _v , v∈ V(T)].

It is not difficult to see that the spaces L ^T cannot be quite arbitrary. In particular, the barycentric coordinates should interpolate at the vertices, and in fact

(9.8) $b_v^T\left(V\right)=1and{b}_v^T\left(e_v\right)=0,v\in V\left(T\right).$

This implies that L ^T restricted to any of the three edges of T is two dimensional, rather than three dimensional. In addition, a requirement for global continuity (C ⁰) of a spline function in the "linear" spline space

$S_1^{-1}\left(\Delta \right):=\left\{f:S\to R,f{\backslash }_T\in L^T,T\in \Delta \right\},$

is that the spaces L ^T |_eand LT′|_e, corresponding to any pair of neighboring triangles T and T′, sharing a common edge e = T ∩ T′, should be identical, i.e.,

(9.9) $L^T{|}_e=L^{T^t}{|}_e.$

Otherwise it would be impossible to join two neighboring triangular patches in this (or any induced higher-degree) spline space continuously, let alone smoothly.

Assuming that all spaces L ^T , T ∈ Δ, satisfy (9.8) and (9.9), we can define the space of continuous splines of degree n on Sas

$S_n^0\left(\Delta \right):=\left\{f\in C^o\left(S\right),f{|}_T\in p_n^T,T\in \Delta \right\},$

where p ^T _n is the space of "polynomials" of degree n on T ∈ Δ, defined as

$p_n^T:=\{f={\displaystyle \sum _{\left|i\right|=n}{c}_i^T{B}_i^n,c_i^T\in R}\},$

and where B ^n,T _i are the "Bernstein polynomials" or B-polynomials

$B_i^{n,T}\left(s\right):=\frac{n!}{{\displaystyle \prod {}_{v\in v\left(T\right)}\left(i_v!\right)}}{\displaystyle \prod _{v\in V\left(T\right)}{\left(b_v^T\left(s\right)\right)}^{iv},s\in T},$

for every multi-index i = (i _v ) _v∈V(T)with $\left|i\right|:={\Sigma }_{v\in V\left(T\right)}{i}_v=n.$ . Note that the B-polynomials reduce to the ordinary (algebraic) polynomials if b ^T _v are the usual barycentric coordinates associated with a planar triangle T. Also note that

${\displaystyle \sum _{\left|i\right|=n}{B}_i^{n,T}={\left({\displaystyle \sum _{v\in V\left(T\right)}{b}_v^T}\right)}^n,}$

and that the B-polynomials are linearly independent (hence form a basis for P ^T _n ), as a consequence of (9.8).

Do conditions (9.8) and (9.9) guarantee that the spaces S ⁰ _n (Δ), n>1, will automatically contain nontrivial smooth functions, C ₁say? It should not come as a surprise that the answer is negative. This can be seen in the planar case S = R ² if we choose a set of non-standard functions b ^T _v . Intuitively, the reason for this is that conditions (9.8) and (9.9) do not enforce compatibility of derivatives of these functions across the edges of the triangulation.

Example 3

In Example 2, the sphere Swas parametrized using the spherical coordinates. An alternative is to parametrize Sby the standard octahedron with vertices v ₁ = (0,0,1), v ₂ = (1,0,0), v ₃ = (0,1,0), v ₄ = (−1,0,0), v ₅ = (0, −1, 0), v ₆ = (0, 0, −1) (or by any polyhedron inscribed in S). The vertices of the octahedron give rise to a triangulation Δ of the sphere, consisting of eight spherical triangles. Consider two adjacent triangles in Δ, say T and T′, determined by their vertices V(T) = {v ₁, v ₃, v ₂} and V(T′) = {v ₁, v ₃, v ₄}, respectively. Let us now define the functions b ^T _v , b ^T′ _v for the two triangles as the usual barycentric coordinates associated with the planar triangles with vertices V(T) and V(T′). One can check that this leads to

$b_{v1}^T\left(s\right)=\frac{z}{x+y+z},b_{v3}^T\left(s\right)=\frac{y}{x+y+z},b_{v2}^T\left(s\right)\frac{x}{x+y+z},$

where s = (x, y, z) ∈ T and

$b_{v1}^{T\text{'}}\left(s\right)=\frac{z}{-x+y+z},b_{v3}^{T\text{'}}\left(s\right)=\frac{y}{-x+y+z},b_{v2}^{T\text{'}}\left(s\right)\frac{-x}{-x+y+z},$

where s = (x, y, z) ∈ T′. The reader is invited to verify that (9.8) and (9.9) are satisfied in this case. However, given a function f ∈ P ^T ₂, it may be impossible to find an f′ ∈ P ^T′ ₂such that the two functions join smoothly (C ¹) across the common edge v ₁ v ₃(for example, take f = (b ^T _v1)²).

The following proposition gives a necessary condition for a smooth join between neighboring triangles [54].

Proposition 1

Let T, T′ ∈ Δ be two adjacent triangles on Ssuch that T ∪ T′ is homeomorphic to a disk in R ². Let b ^T _v ∈ C ^∞(T), v∈ V(T), and b ^T′ _v ∈ C ^∞(T′), v∈ V{T′). Suppose that for every n and every f ∈ P ^T _n there exists an f′ ∈ P ^T′ _n such that f and f′ join with C ⁿ⁻¹ continuity along the common edge T ∩ T′. Then each b ^T _v ∈ T, can be extended as a C ^∞ function on T ∪ T′ which, when restricted to T′, belongs to L ^T′ . Equivalently, there exists a three-dimensional space L of C ^∞ functions on T ∪ T′ such that L| _T = L ^T an d L| _T′ = L ^T′ .

A consequence of this result is that the choice of the functions b ^T _v is greatly restricted. In particular, the proposition says that all barycentric coordinates b ^T _v corresponding to any given triangle T, can be smoothly extended across edges of T to neighboring triangles and hence all spaces L ^T must locally belong to a single three-dimensional space of "linear functions" L. The following example shows that such a space L might not contain globally continuous functions, i.e., for some surfaces one may be able to define L only locally.

Example 4

Let Sbe the circular cylinder parametrized by

$\left(\sin \varphi ,\cos \varphi ,z\right)\in R^3,\varphi \in R,z\in R,$

and let

$L:=span\left\{1,\varphi ,z\right\}.$

Hence, the point on Swhose parameters are (ϕ, z) is assigned the value a + bϕ + cz, where a, b, c are real coefficients. Note that L is not defined globally since there is no globally continuous function in it. However, on every "strip"

$\{\left(\varphi ,z\right)\varphi \in \left(\alpha ,\beta \right),z\in R\},$

where β– α≤ 2π, the space L is well defined and generated by C ^∞ functions.

If Δ is a geodesic triangulation of Sconsisting of "small" triangles (i.e., triangles for which the values of ϕare in an interval of length not exceeding pI), then it is not difficult to see that for every such triangle T, dim(L ^T | _e ) = 2, e ∈ E(T), where L ^T := L| _T . This follows from the observation that every f ∈ L vanishes along geodesics. Namely, if f(ϕ, z) = a + bϕ + cz, where b ∞ 0, then f(ϕ, z) = 0 along the curve

$\left\{(\sin \left(\left(a+cz\right)/b\right),\cos \left(\left(a+cz\right)/b\right),z),z\in R\right\},$

which is a helix, hence a geodesic. Otherwise, if b = 0, z is a constant, which also corresponds to a geodesic on S.

The above discussion shows that for every sufficiently small geodesic triangle T ∈ Δ, one can define cylindrical barycentric coordinates. It turns out that many properties of these coordinates are similar to the properties of their planar counterparts. In particular, the cylindrical barycentric coordinates make it possible to define a spline space, which can be used to construct smooth splines on S. Figure 9.5shows an example of a cylindrical triangulation Δ and Figure 9.6shows a C ¹ smooth spline in the space S ⁰ ₅(Δ), corresponding to this triangulation.

Figure 9.5. Triangulation of a cylinder.

Figure 9.6. A C ¹ smooth quintic spline on the cylinder.

Example 4 raises the question whether a space L, or a collection of spaces L ^T satisfying the conditions implied by Proposition 1, always exists on any Sand, if so, how can one find such spaces. A simple approach to constructing L is to take advantage of our knowledge of S, i.e., the implicit assumption here that we can evaluate Sat any point. In particular, we can define L as the span of the three functions x(s), y(s), z(s), the Cartesian coordinates of s i.e., s = (x(s), y(s), z(s)), s∈ S. This works well in the important special case of a sphere in R ³.

Example 5.

The problem of constructing spherical analogs of Bézier triangles has an interesting history. It has received ample attention in the spline literature for the obvious reason that splines on a sphere have many potential applications in geosciences, including meteorology, geophysics, and geodesy. Researchers had been searching for many years for appropriate spherical analogs of Bézier methods, but were hampered by the difficulty of defining spherical barycentric coordinates. For example, several candidates for such coordinates have been introduced in [13],[14],[41], but they lacked many key properties of the planar coordinates. As it turned out, the mentioned attempts were destined to be unsuccessful for the simple reason that they all insisted on the partition of unity property, which is well known for the classical barycentric coordinates. Namely, it was shown in [14] that spherical barycentric coordinates that sum to one and satisfy a list of other reasonable assumptions, do not exist. This negative result provided an explanation why the various earlier generalizations were unsuccessful in building smooth splines on spherical triangulations.

A breakthrough in this development came when the authors of [1] realized, by studying identities of spherical trigonometry, that there is in fact a natural way of defining barycentric coordinates for spherical triangles. Let T be a spherical triangle with vertices v ₁, v ₂, v ₃∈ S. Thus the edges of T are the three shortest geodesies connecting the vertices. One can define the spherical barycentric coordinates of a point s∈ T as

$b_{vi}^T\left(s\right):=\frac{\sin {\delta }_i}{\sin \left({\delta }_i+{\gamma }_i\right)},i=1,2,3,$

where Δ _i is the geodesic distance (measured along the great circle passing through s and V _i ) of u _i and s, and γ _i is the geodesic distance of s and V _i (see Figure 9.7).

Figure 9.7. Definition of spherical barycentric coordinates.

Despite the inevitable fact that the above coordinates do not add up to one, it has been shown in [1] that they resemble the standard barycentric coordinates in many respects. First, the barycentric coordinates are infinitely smooth functions on T and they satisfy (9.8). The space L ^T , the span of the three coordinates, reduces to dimension two along the edges of T and in fact along every great circle intersecting T. More precisely, the restriction of L ^T to any such great circle can be shown to be the linear span of {sin α, cos α}, where αis the arclength distance measured along this circle. Another important consequence of the above definition is that the spaces L ^T and L ^T′ , associated with neighboring triangles T and T′, satisfy (9.9). Furthermore, the properties of C7 are consistent with those specified in Proposition 1. In particular, this space can be extended to a three-dimensional space L of infinitely differentiate functions over all of S. Conversely, for any triangle T the space L ^T is just the restriction of L to T.

The space L has many interesting properties that indicate that spherical barycentric coordinates, as defined here, are unique. More precisely, L is the only three-dimensional space of functions on the sphere Sthat is rotationally invariant and such that its dimension is reduced along great circles. Moreover, L, can be easily described using the idea suggested earlier. Namely, L is precisely the span of the Cartesian coordinates x(s), y(s), z(s), viewed as functions on S. This is equivalent to saying that L is the space of spherical harmonics of degree one [49].

The intimate connection between spherical harmonics and spherical barycentric coordinates is not unexpected. In retrospect, it seems quite obvious that these functions should have been considered early on as the most natural candidates for "linear functions" on the sphere. As a matter of fact, such functions, along with the corresponding barycentric coordinates defined above, had been considered some 150 years before the paper [1] appeared. The authors of that paper found out only later that their idea of defining barycentric coordinates was not new after all, since the same definition had already been given by Möbius [47].

Spherical barycentric coordinates give rise to Bernstein-Bézier-type methods, with immediate applications to a variety of problems on the sphere. Indeed, because of the close analogy with standard Bernstein-Bézier techniques, virtually all of the classical methods for piecewise polynomials on planar triangulations can be carried over to the spherical setting, and indeed to any setting where barycentric coordinates are available. A detailed treatment of some of these methods has been given in [ 2],[3]. The spherical Bernstein-Bézier methods are also of interest in the design of surfaces, especially star-like surfaces, even though some of the geometric properties of planar Bézier methods are missing on the sphere, such as the convex hull property.

Typical scattered data interpolation/approximation methods on the sphere start with a triangulation of the sphere, for example the so-called Delaunay triangulation. We refer the reader to [69], for a survey on triangulations, and also to [11],[39], for a discussion of triangulation methods on general surfaces. Methods for triangulating scattered data points on the sphere are discussed in [41],[58],[66]. Figures 9.9and 9.10show wire plots of smooth C ¹ quadratic and cubic spherical splines, respectively, corresponding to (a part of) the triangulation in Figure 9.8. Details about various data-fitting methods that lead to such spherical splines are discussed in [3].

Figure 9.9. A C ¹ smooth quadratic spherical spline.

Figure 9.10. A C ¹ smooth cubic spherical spline.

Figure 9.8. A triangulation of the sphere.

Besides the sphere, the suggestion to use Cartesian coordinates to define the space L makes sense, as long as the triangulation Δ of S consists of triangles whose edges reduce the dimension of L. This is equivalent to the condition that the edges of Δ are planar and in the same plane with the origin (0, 0, 0). This is acceptable for certain surfaces (for example, for star-like surfaces [3]), but is more restrictive in cases in which the assumption of coplanarity would result in severely distorted triangles (relative to the geodesic ones).

To cope with this difficulty, a possibility is to make a coordinate transformation to minimize the distortion of the triangles, e.g., one could shift the origin of the Cartesian system. More generally, one could in principle take "any" three-dimensional space K of smooth functions in R ³ and define L as the restriction of H to S. In the spherical case, this would lead to L = H|s, where H:= span{x, y, z}. In the planar case S = R ², we can think of L as H|s, where H:= span{l, x, y}. In this way one can also interpret the construction of L for the cylinder. In particular, we can take for any fixed α∈ R,

$H:=span\left\{\arctan \left(\frac{x\cos \alpha -y\sin \alpha }{x\sin \alpha +y\cos \alpha }\right),z,1\right\},$

in which case L is the restriction of H to the strip

$\left\{\left(\varphi ,z\right),\alpha -\pi /2<\varphi <\alpha +\pi /2,z\in R\right\}.$

On the other hand, it is not clear how to choose H for a general surface Sso as to obtain "least distorted" triangulations.

Above, we have seen two examples of non-planar surfaces for which it is possible to construct meaningful analogs of spline spaces. The characteristic property shared by both types of splines, as well as by the classical bivariate splines on planar triangulations, is that they can be generated by a three-dimensional space L of "linear functions" on S.

This space is such that

•

Functions in L vanish along geodesics on S i.e., for every nontrivial f ∈ L, the set C:= {s∈ S, f(s) = 0} is a geodesic. Conversely, for every geodesic Con S, there exists a nonzero function f ∈ L vanishing on C.

•

L is invariant under isometric transformations of S. Thus, if I: S→ Sis an isometry (e.g., a rotation, if Sis the sphere), then f o I ∈ L for all f ∈ L.

It is clear that these two properties imply that L is quite special and that for a given surface Sone cannot expect to have more than one such space. In fact, in the three mentioned cases, L is uniquely determined by the above two properties. It is an intriguing open question whether there always exists a space satisfying at least the first property. This in turn is closely related to the central issue of this section of whether one can find analogs of spline spaces on general surfaces. Using the well-known Beltrami's Theorem about the existence of local geodesic mappings [18], it can be shown that L can be found for surfaces of constant (Gaussian) curvature [54]. However, it is not known if one can go beyond such surfaces.

It should be said that, for all practical purposes, we need not require that L strictly satisfy the mentioned conditions. The lack of a theoretical proof of the existence of L should not prevent us from being able to establish useful spline spaces on S. For example, for a given fixed geodesic triangulation Δ, one could use a space L which has a reduced dimension along all of its edges, but which does not necessarily have the property that all functions in L vanish along geodesics. Such space L will still allow the construction of barycentric coordinates for all triangles of Δ, hence also the construction of a spline space corresponding to Δ.

Another possibility is to relax the assumption that the triangles in Δ are strictly geodesic. In fact, parametric surfaces Scomposed of triangular Bézier patches are already equipped with a natural triangulation Δ in which every triangle corresponds to a Bézier patch. Such a triangulation is in general not geodesic. In this situation there is a particularly simple way to choose the barycentric functions b ^T _v . Recall that Example 3 gives reasons why it is in general not a good idea to use the standard barycentric coordinates corresponding to the triangular facets of the piecewise linear interpolant to the vertices of Δ. This is because this will in general not permit us to design smooth splines over S. However, it turns out that this approach will work if Sis a composite surface consisting of triangular polynomial patches. Namely, the fact that neighboring triangular patches are joined smoothly guarantees that the ordinary planar barycentric coordinates corresponding to these triangles are compatible, in the sense of Proposition 1. Thus in any given triangle these coordinates can be extended to the neighboring triangles as smooth functions (where the degree of smoothness will depend on the smoothness of S).

We have seen that the success or failure of constructing spline spaces on general surfaces hinges upon the existence of barycentric coordinates. It is also clear from our examples, that even for simple surfaces such constructions may be far from trivial. Still, the discussed framework of splines on surfaces offers many benefits, compared to other existing reconstruction methods. The main argument supporting this claim is that the Bernstein-Bézier formalism for splines on surfaces is the same for all surfaces. That is, it is essentially irrelevant whether we work on the sphere, in the plane, or on any other surface. As a consequence, we can use the same algorithms for splines on all surfaces, as long as we use a correct procedure to compute the barycentric coordinates (which do depend on the type of the surface).

Example 6.

To illustrate the above point, suppose that T, T′ ∈ Δ are two adjacent triangles on S, with vertices V(T) = {v ₁, v ₃, v ₂} and V(T′) = (v ₁, v ₂, v ₃). Let f and f′ be two "polynomials" of degree n of the form

$f={\displaystyle \sum _{\left|i\right|=n}{c}_i^T{B}_i^{n,T},f\text{'}{\displaystyle \sum _{\left|i\right|=n}{c}_i^{T\text{'}}{B}_i^{n,T\text{'}}},}$

where c ^T _i , c ^T′ _i ∈ R are given Bézier coefficients. Then f and f′ join continuously along the edge V ₁ v ₃if

$c_{i1,i_3,0}^{T\text{'}}=c_{i_1,i_3,0}^T{i}_1+i_3\mathrm{=n},$

and they join with C ¹ continuity if and only if

$c_{i_1,i_3,1}^{T\text{'}}=c_{i_1+1{,}_{i_{3,}0}}^T{b}_{v1}^T\left(v_4\right)+c_{i_1,i_3+1,0}^T{b}_{v_3}^T\left(v_4\right)+c_{i_1,i_3,1}^T{b}_{v_2}^T\left(v_4\right),i_1+i_3=n-1,$

where we used the abbreviations i _j := i _{v j} , j = 1, 2, 3. Moreover, the values b ^T _{v i} (v ₄) are obtained by first extending b ^T _{v j} . as smooth functions to T ∪ ^T′ , and then evaluating them at v ₄. This extension is possible as long as the barycentric coordinates are compatible in the sense of Proposition 1.

The above conditions for a smooth join between two generalized Bézier triangles can be immediately seen to be formally identical to the corresponding conditions in the planar setting. This explains why we have not addressed in this section the actual reconstruction problem (e.g., interpolation or approximation) using splines on surfaces. The reason is that the same methods known in the plane can be transformed almost automatically to any surface. Various extensions of such planar methods to the sphere are discussed in [3] and it is indeed clear from that paper that, except for a few details, the planar methods carry over to the spherical setting, and indeed to the setting of any smooth surface S.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780444511041500101

Bézier Triangles

Gerald Farin , in Curves and Surfaces for CAGD (Fifth Edition), 2002

17.8 Nonparametric Patches

In an analogy to the univariate case, we may write the function

(17.26) $z=\begin{array}{l}{\displaystyle \sum _{\left|\text{i}\right|=n}{b}_{\text{i}}{B}_{\text{i}}^n\left(\text{u}\right)}\hfill \end{array}$

as a surface

$\left[\begin{array}{l}u\\ v\\ w\\ z\end{array}\right]=\begin{array}{l}{\displaystyle \sum _{\left|\text{i}\right|=n}\left[\begin{array}{l}i/n\\ j/n\\ k/n\\ b_i\end{array}\right]B_{\text{i}}^n\left(\text{u}\right).}\hfill \end{array}$

Thus the abscissa values of the control polygon of a nonparametric patch are given by the triples i/n, as illustrated in Figure 17.14. The last equation holds because of the linear precision property of the Bernstein polynomials $B_{\text{i}}^n$ ,

Figure 17.14. Nonparametric patches: the abscissas of the control net are the n-partition points of the domain triangle.

$u=\begin{array}{l}{\displaystyle \sum _{\left|\text{i}\right|=n}\frac{i}{n}{B}_{\text{i}}^n\left(\text{u}\right),}\hfill \end{array}$

and analogous formulas for v and w. The proof is by degree elevation from 1 to n of the linear function u. Example 17.4 shows a nonparametric patch.

Nonparametric Bézier triangles play an important role in the investigation of spaces of piecewise polynomials, as studied in approximation theory. Their use has facilitated the investigation of one of the main open questions in that field: what is the dimension of those function spaces? (See, for instance, Alfeld and Schumaker [7].) They have also been useful in defining nonparametric piecewise polynomial interpolants; see, for example, Barnhill and Farin [29], Farin [194], Petersen [471], or Sablonnière [525].

Example 17.4

A nonparametric quadratic patch.

The bivariate function z = x ² + y ² may be written as a quadratic nonparametric Bézier patch over the triangle (0, 0), (2, 0), (0, 2). Its coefficients are:

$\begin{array}{lll}\left[\begin{array}{l}0\\ 2\\ 4\end{array}\right]\hfill & \hfill & \hfill \\ \left[\begin{array}{l}0\\ 1\\ 0\end{array}\right]\hfill & \left[\begin{array}{l}1\\ 1\\ 0\end{array}\right]\hfill & \hfill \\ \left[\begin{array}{l}0\\ 0\\ 0\end{array}\right]\hfill & \left[\begin{array}{l}1\\ 0\\ 0\end{array}\right]\hfill & \left[\begin{array}{l}2\\ 0\\ 4\end{array}\right]\hfill \end{array}.$

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B978155860737850017X

B-Spline Curves

Gerald Farin , in Curves and Surfaces for CAGD (Fifth Edition), 2002

8.8 B-Splines

Consider a knot sequence u ₀, …, u_M and the set of piecewise polynomials of degree n defined over it, where each function in that set is n – r_i times continuously differentiable at knot u_i . All these piecewise polynomials form a linear space, with dimension

(8.13) $\text{dim}=\left(n+1\right)+{\displaystyle \sum _{i=1}^{M-1}{r}_i.}$

For a proof, suppose we want to construct an element of our piecewise polynomial linear space. The number of independent constraints that we can impose on an arbitrary element, or its number of degrees of freedom, is equal to the dimension of the considered linear space. We may start by completely specifying the first polynomial segment, defined over [u ₀, u ₁]; we can do this in n + 1 ways, which is the number of coefficients that we can specify for a polynomial of degree n. The next polynomial segment, defined over [u ₁ u ₂], must agree with the first segment in position and n – r ₁ derivatives at u ₁, thus leaving only r ₁ coefficients to be chosen for the second segment. Continuing further, we obtain (8.13).

We are interested in B-spline curves that are piecewise polynomials over the special knot sequence [u_{n –1}, u_L ]. The dimension of the linear space that they form is L + 1, which also happens to be the number of B-spline vertices for a curve in this space. If we can define L + 1 linearly independent piecewise polynomials in our linear function space, we have found a basis for this space. We proceed as follows.

Define functions $N_i^n\left(u\right)$ , called B-splines by defining their de Boor ordinates to satisfy d_i = 1 and d_j = 0 for all j ≠ i. The $N_i^n\left(u\right)$ are clearly elements of the linear space formed by all piecewise polynomials over [u_{n– 1}, u_L ]. They have local support:

$N_i^n\left(u\right)\ne 0\text{only}\text{if}u\in \left[u_{i-1},u_{i+n}\right].$

This follows because knot insertion, and hence the de Boor algorithm, is a local operation; if a new knot is inserted, only those Greville abscissae that are "close" will be affected.

B-splines also have minimal support: if a piecewise polynomial with the same smoothness properties over the same knot vector has less support than $N_i^n$ , it must be the zero function. All piecewise polynomials defined over [u_{i– 1}, u_{i +n} ], the support region of $N_i^n$ , are elements of a function space of dimension 2n + 1, according to (8.13). A support region that is one interval "shorter" defines a function space of dimension 2n. The requirement of vanishing n – r_{i –1} derivatives at u_i–1 and of vanishing n – r_{i +n} derivatives at u_{i +n} imposes 2n conditions on any element in the linear space of functions over [u_{i– 1}, u_{i +n–1}] The additional requirement of assuming a nonzero value at some point in the support region raises the number of independent constraints to 2n + 1, too many to be satisfied by an element of the function space with dimension 2n.

Another important property of the $N_i^n$ is their linear independence. To demonstrate this independence, we must verify that

(8.14) ${\displaystyle \sum _{j=0}^L{c}_j{N}_j^n\left(u\right)\equiv 0}$

implies c_j = 0 for all j. It is sufficient to concentrate on one interval [u_I , u_{I +1}] with u_I < u_{I +1}. Because of the local support property of B-splines, (8.14) reduces to

$\begin{array}{cc}{\displaystyle \sum _{j=I-n+1}^{I+1}{c}_j{N}_j^n\left(u\right)\equiv 0}& \text{for}\text{u}\in \left[u_I,u_{I+1}\right].\end{array}$

We have completed our proof if we can show that the linear space of piecewise polynomials defined over [u_I–n , u_{I +n+1}] does not contain a nonzero element that vanishes over [u_I , u_{I +1}]. Such a piecewise polynomial cannot exist: it would have to be a nonzero local support function over [u_{I +1},u_{I +n+1}]. The existence of such a function would contradict the fact that B-splines are of minimal local support.

Because the B-splines $N_i^n$ are linearly independent, every piecewise polynomial s over [u_{n– 1}, u_L ] may be written uniquely in the form

(8.15) $s\left(u\right)={\displaystyle \sum _{j=0}^L{\text{d}}_j{N}_j^n\left(u\right).}$

The B-splines thus form a basis for this space. This reveals the origin of their name, which is short for Basis splines. Figure 8.17 gives examples of some cubic B-splines.

Figure 8.17. B-splines: some cubic examples.

If we set all d_i = 1 in (8.15), the function s(u) will be identically equal to 1, thus asserting that B-splines form a partition of unity.

B-spline curves are simply the parametric equivalent of (8.15):

$x\left(u\right)={\displaystyle \sum _{j=0}^L{\text{d}}_j{N}_j^n\left(u\right).}$

Just as the de Casteljau algorithm for Bézier curves is related to the recursion of Bernstein polynomials, the de Boor algorithm yields a recursion for B-splines. It is given by

(8.16) $N_I^n\left(u\right)=\frac{u-u_{l-1}}{u_{l+n-1}-u_{l-1}}{N}_l^{n-1}\left(u\right)+\frac{u_{l+n}-u}{u_{l+n}-u_l}{N}_{l+1}^{n-1}\left(u\right),$

with the "anchor" for the recursion being given by

(8.17) $N_i^0\left(u\right)=\{\begin{array}{ll}1\hfill & \text{if}{u}_{i-1}\le u<u_i\hfill \\ 0\hfill & \text{else}\hfill \end{array}.$

Its proof relates the local recursion (8.10) to the global indexing scheme. An example is shown in Figure 8.18.

Figure 8.18. The B-spline recursion: top, two linear B-splines yield a quadratic one; bottom, two quadratic B-splines yield a cubic one.

Equation (8.16) is due to L. Mansfield, C. de Boor, and M. Cox; see de Boor [137] and Cox [129]. For an illustration of (8.16), see Figure 8.18. This formula shows that a B-spline of degree n is a strictly convex combination of two lower-degree ones; it is therefore a very stable formula from a numerical viewpoint. If B-spline curves must be evaluated repeatedly at the same parameter values u_k , it is a good idea to compute the values for $N_i^n\left(u_k\right)$ using (8.16) and then to store them.

A comment on end knot multiplicities: the widespread data format IGES uses two additional knots at the ends of the knot sequence; in our terms, it adds knots u_{– 1} and u_{L +2n–1}. The reason is that formulas like (8.16) seemingly require the presence of these knots. Since they are multiplied only by zero factors, their values have no influence on any computation. There is no reason to store completely inconsequential data, and hence the "leaner" notation of this chapter.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9781558607378500089

Interrogation of Subdivision Surfaces

Malcolm Sabin , in Handbook of Computer Aided Geometric Design, 2002

13.2 HISTORICAL BACKGROUND

One of the first ways in which subdivision impinged on the CAGD community was through graphics. The hidden surface problem spawned a large number of different approaches during the early 1970s, and one of the novelties then was the idea of Warnock[19], whose algorithm was essentially 'if the picture is simple, draw it; if not, use this same approach for a quarter at a time.'. This was an image space algorithm and the subdivision concept was soon applied in object space algorithms and then in parameter space algorithms for surfaces. [18],?] This was stimulated by the development of knot insertion algorithms, notably the Oslo algorithm [6] and the Boehm algorithm [2].

The motivation for these was based on two observations:-

1.

That for any piecewise polynomial curve, typically represented by the combination of B-spline basis functions with point-valued coefficients, knots can be inserted at places where there are no actual discontinuities of derivative, and that if enough knots are inserted at the same place, the original curve falls into two pieces, each represented by its own subset of control points.

2.

That such a curve always lies inside the convex hull of its control points.

If we wish to determine such things as the intersection points of curves, or the intersections of curves or surfaces with planes or with each other, these two properties give an elegant recursive approach.

As a simple example, if we wish to calculate the points of intersection of a curve with a plane, we have a recursive algorithm exactly echoing the Warnock structure

'If the curve is simple get the intersection directly; if not, apply this same algorithm to two halves of it, taking the union of any points which are found in the halves'

There is one addition, which makes it much more efficient than just dividing the curve into a large number of simple pieces ab initio and considering them one at a time. If, at each level, we test whether the convex hull crosses the plane, we can discover whether the curve lies entirely on one side of the plane, in which case there is no intersection, and no further subdivision need take place.

Such an algorithm promises both efficiency and robustness.

Algorithms in this style are capable of supporting the wide range of interrogations which are necessary in industrial use of surfaces represented on the computer; they give the very high levels of robustness which are essential if the surface handling is embedded inside a solid modelling context, and with the dramatic growth of real memories in the last decade or two are potentially the fastest available approach.

There is, however, a significant amount of sophistication which can and must be added to the basic idea in order to turn it into a widely applicable tool, and this chapter deals with that detailed technology by taking the simple approach and exploring how it needs to be tuned.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780444511041500149

Multiscale Wavelet Methods for Partial Differential Equations

Christoph Schwab , Tobias von Petersdorff , in Wavelet Analysis and Its Applications, 1997

Abstract

We analyze multiscale Galerkin methods for strongly elliptic boundary integral equations of order zero on closed surfaces in ℝ³ . Piecewise polynomial, discontinuous multiwavelet bases of any polynomial degree are constructed explicitly. We show that optimal convergence rates in the boundary energy norm and in certain negative norms can be achieved with "compressed" stiffness matrices containing O(N (log N)²) nonvanishing entries, where N denotes the number of degrees of freedom on the boundary manifold. We analyze a quadrature scheme giving rise to fully discrete methods. We show that the fully discrete scheme preserves the asymptotic accuracy of the Galerkin scheme with exact integration and without compression. The overall computational complexity of our algorithm is O(N (log N)⁴) kernel evaluations. The implications of the results for the numerical solution of elliptic boundary value problems in or exterior to bounded, three-dimensional domains are discussed.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/S1874608X9780009X

Handbook of Numerical Methods for Hyperbolic Problems

J. Qiu , Q. Zhang , in Handbook of Numerical Analysis, 2016

4.3 Scalar Equation with Discontinuous Initial Solution

It is well known that the numerical oscillation will happen when the initial solution contains a discontinuous point and the piecewise polynomials k ≥ 1. Consider the linear constant hyperbolic equation, namely f(u) = βu in (1). The detailed analysis shows that the pollution region around the characteristic line across the discontinuous point is only restricted in a narrow zone, and the optimal error estimate outside the pollution region is also preserved. The next theorem (Zhang and Shu, 2014) stated this double-optimal result for RKDG(3,3,k) method with arbitrary k ≥ 0.

Theorem 3

Assume T = Mτ for simplicity, and k ≥ 1. Under the standard CFL condition, namely, τ ≤ λh with λ being suitably small, there is the optimal error estimate for RKDG(3,3,k) method to solve the linear constant hyperbolic equation

$\parallel u(\cdot ,T)-u_h^M{\parallel }_{L^2(\mathbb{R}\setminus {\mathcal{R}}_T)}\le C_1(h^{k+1}+{\tau }^3),$

out of the pollution region

${\mathcal{R}}_T=\left(\beta T-C_2\sqrt{T\beta }{h}^{1/2}log\frac{1}{h},\beta T+C_3\sqrt{T\beta }{h}^{1/2}log\frac{1}{h}\right).$

The above bounding constants are all independent of h,   τ and λ⁻¹ .

This theorem is proved by energy analysis with two special weight functions, in order to detect the left side and right side of pollution region. The analysis is long and complex, which involves many technical points, for example, the superconvergence results, the generalized slope function and the highest frequency component, as well as the complex treatment of those errors coming from the Runge–Kutta time marching.

The above result is independent on λ ⁻¹ and hence also holds for the semidiscrete DG method. The numerical results given in Zhang and Shu (2014) verify the sharpness of the above results.

Remark 4

Similar works have been carried out by many authors. For example, Johnson and Pitkäranta (1986) considered the space-time DG method and proved that the pollution region at the crosswind direction has the width of $\mathcal{O}(h^{1/2}log(1/h))$ , and Cockburn and Guzman (2008) considered RKDG(2,2,1) method and proved similar results.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/S1570865916300011

Multiscale Wavelet Methods for Partial Differential Equations

Panayot S. Vassilevski , Junping Wang , in Wavelet Analysis and Its Applications, 1997

9.1 Galerkin discretization

The Galerkin method for the approximation of u is based on the variational problem (9.3). Let V  = V_h be a C ⁰ -conforming finite element space of piecewise polynomials corresponding to a quasiuniform triangulation T_h of Ω. The Galerkin approximation is a function u_h ϵ V_h such that

(9.1.1) $\begin{array}{l}{b}_{є}\left(u_h,\psi \right)\equiv є\left(\nabla u_h,\nabla \psi \right)+\left(\psi ,\underset{¯}{b}\cdot \nabla u\right)+\delta \left(\underset{¯}{b}\cdot \nabla u,\underset{¯}{b}\cdot \nabla \psi \right)\hfill \\ -є\delta {\displaystyle \sum _{T\mathrm{є}{T}_h}}{\displaystyle \underset{T}{\int }}\Delta u_h\left(\underset{¯}{b}\cdot \nabla \psi \right)dx\hfill \\ =\left(f_{\delta },\psi \right)\hfill \end{array}$

for all ψ ϵ V_h. For continuous piecewise linear functions, one has ∆u_h   =   0 on each element. It follows that the discrete problem seeks u_h ∈ V_h such that

(9.1.2) $\begin{array}{l}\mathrm{ϵ}\left(\nabla u_h,\nabla \psi \right)+\left(\psi ,\underset{¯}{b}\cdot \nabla u_h\right)+\mathrm{\delta }\left(\underset{¯}{b}\cdot \nabla u_h,\underset{¯}{b}\cdot \nabla \psi \right)=\left(f_{\delta },\psi \right)\hfill \\ \hfill \end{array}$

for all ψ ϵ V_h .

The convection term is assumed to satisfy

(9.1.3) $\nabla \cdot \underset{¯}{b}\le 0\mathrm{in}\Omega \text{.}$

For a convergence analysis of the streamline diffusion finite element approximation u_h , we refer to [23] and [2].

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/S1874608X97800040

24th European Symposium on Computer Aided Process Engineering

Rubén Ruiz-Femenia , ... Ignacio E. Grossmann , in Computer Aided Chemical Engineering, 2014

2.2 Discretization using orthogonal collocation

We transform the disjunctive multistage problem into a discretized GDP problem by orthogonal discretization, a simultaneous method that fully discretizes the DAE system by approximating the control and state variables as piecewise polynomials functions over finite elements ( Kameswaran and Biegler, 2006). Figure 2 shows how the time horizon is discretized. Accordingly, at each collocation point the state variable is represented by:

Figure 2. Discretization scheme used to apply the orthogonal collocation method.

(1) $x_{sik}=x_{si}^0+h_{si}{\displaystyle \sum _{j=1}^K{\Omega }_j\left(T_k\right){\dot{x}}_{sij},s=1,\dots ,S,}i=1,\dots ,I,k=1,\dots ,K$

where $x_{si}^0$ is the value of the state variable at the beginning of element i in stage s, ${\dot{x}}_{sij}$ is the value of its first derivative in element i at the collocation point k in stage s, h_is is the length of element i, in stage s, Ω _j (T_k ) is the interpolation polynomial of order K for collocation point j, and T_k , is the non-dimensional time coordinate. We enforce continuity of the state variable across finite element boundaries in each stage by $x_{si}^0=x_{s,i-1,K}$ for all s  =   1,…,S, i  =   2,…,I. Additional stage transition conditions map the differential state variable values across the stage boundaries:

(2) $x_{s1}^0-x_{s-1,I,K}=0,s=2,\dots ,S$

For each stage, the collocation method requires the time to be discretized over each finite element at the selected collocation points t_sik :

(3) $t_{sik}=t_{si}^0+h_{si}{T}_k,s=1,\dots ,S,i=1,\dots ,I,k=1,\dots ,K$

where $t_{si}^0$ is the value of the time at the beginning of element i in stage s. Time continuity between stages and between elements within a stage is also enforced by the following constraints:

(4) $\begin{array}{l}{t}_{s,1}^0=t_{s-1,I,K},s=2,\dots ,S\\ t_{si}^0=t_{s,i-1,K},s=1,\dots ,S,i=2,\dots ,I\end{array}$

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780444634566500570

clancynohnerhed93.blogspot.com

Source: https://www.sciencedirect.com/topics/computer-science/piecewise-polynomial

Share this post