1 Introduction
There has recently been an increased academic and industrial interest in high-order finite elements due to their efficiency advantages over classical low-order elements, see e.g.[1, 2]. To achieve their full potential high-order methods require the curved domain boundaries to be approximated with non-planar elements. The intersection between the curved elements and a given line is required in a wide range of applications, including contact [3], mesh generation [4, 5] and immersed finite elements [6]. As known, intersection computations lead to an easy to formulate, but hard to solve, nonlinear root-finding problem. The prevalent technique in computational mechanics for solving such problems is Newton-Raphson iteration, which is in general not very robust.
In computer-aided geometric design (CAD) and manufacturing (CAM) intersection computation is a recurring task and, to this end, a number of ingenious methods have been developed [7]. Especially promising are the non-iterative methods with a origin in algebraic geometry, which are for the most part unknown in computational mechanics. Algebraic geometry deals with systems of polynomial equations and geometric objects defined by them and provides the most rigorous framework for intersection computations [8]. The specific technique considered in this paper is the implicisation technique proposed by Busé [9], which shares some commonalities with the method of moving lines/planes introduced by Sederberg and Chen [10]. Different from the original work, in this paper we derive the method from a purely linear algebra viewpoint. To follow the presented derivations it is sufficient to only know the notion of the algebraic degree of a curve or surface. The algebraic degree of a curve or surface is defined as its number of intersections with a line. Counting all intersections (real, complex, multiple ones and ones at infinity), the algebraic degree of a polynomial curve of degree
and the corresponding tensor-product surface are
and , respectively. As will become clear, adopting a linear algebra viewpoint has the added benefit that many linear algebra techniques, like preconditioning and various matrix decompositions, become readily available for intersection computations.In the following, we first discuss the intersection between a planar Lagrange curve and a line and provide an easy to follow illustrative example. This simple case is sufficient to introduce and discuss the key aspects of the proposed non-iterative technique. Its extension to the surface case is straightforward and is discussed last.
2 Intersection of lines with curves
2.1 Moving lines and intersections
Let be a planar parametric curve, with , of degree given either in Lagrange basis or power (monomial) basis with
(1) |
where are the nodal coordinates and are the coefficients in the power basis. As usual, the power basis contains the consecutive powers of from up to . The two basis are related by
(2) |
where is the Vandermonde matrix and is the parametric coordinate of the -th Lagrange node. Following a similar approach a curve given in any other polynomial basis can be re-expressed in the power basis.
To define a point on the curve as the intersection of several moving lines, or pencils of lines, consider
(3) |
where
is an auxiliary vector collecting the parameters of the line. For a fixed
equation (3) describes a line and the line moves with the parameter , see Figure (a)a . The three parameters of the line are assumed to be polynomial functions given by(4) |
The degree of the power basis has to be chosen sufficiently high in order to be able to compute all the intersection points (real, complex, multiple ones and ones at infinity). The number of intersection points is equivalent to the algebraic degree of the curve . A curve of degree has intersection points with a line. As will become clear, the number of intersections implies a constraint on the minimum possible degree for .
Next, the aim is to find several lines (3), or more specifically their coefficients in (4), with a common intersection point which is a point on the curve . It is required that each line satisfies at the common intersection point
(5) |
After introducing the definition of the curve (1) this yields
(6) |
The bracketed term can be expressed in a new power basis of dimension with
(7) |
where the matrix components contain the known coefficients . Equation (6) can now be rewritten as
(8) |
where the array contains the components of the yet unknown vectors sorted (by choice) in the following way
(9) |
It is required that (8) is always satisfied irrespective of , which is the case for the right null vectors of the matrix . The right null vectors are determined with a SVD, see e.g. [11], yielding the set of null vectors , where the index denotes the number of the null vector.
The number of null vectors of depends on the degrees and of the basis and , and the coefficients of the specific curve considered. For subsequent computations the number of null vectors must be more than the number of intersections of the curve with a line (or its algebraic degree). The non-square matrix has rows and columns. Hence, its number of right null vectors must be equal or greater than .^{1}^{1}1 The number of right null vectors is larger then, for instance, a quadratic curve is described with a cubic polynomial (i.e. but ). More precisely, if is the largest integer such that then the number of right null vectors of is exactly . In order to obtain the intersections it is necessary to have
(10) |
The set of null vector denoted with introduced in (3) yields a set of moving lines
(11) |
with a common intersection point on the curve . As an example, in Figure (a)a the description of a quadratic curve with by two moving lines and with is shown.
Next, the intersection of a given parametric line
(12) |
with the curve is considered, where and are two given vectors. We require that the determined set moving lines and intersect at the same point to be yet determined. According to (5) at the common intersection point it is required that
(13) |
These equations describe the intersection of each moving line with the given line and can be rewritten as
(14) |
As discussed the number of moving lines satisfies and all of them can be combined in one homogenous equation system
(15) |
with and representing the components of two matrices and . Both matrices have rows and or more columns. To obtain and that satisfy (15) the following generalised eigenvalue problem is considered
(16) |
That is, the eigenvalues are the parametric coordinates of the intersection points on the line
and the eigenvectors are (up to a multiplicative constant) the basis functions
evaluated at the intersection points . Unfortunately, the matrices and are not always square and computing the values that satisfy this equation requires non-standard linear algebra techniques. However, for the purposes of intersection computation it is sufficient to consider a square eigenvalue problem obtained from (16) by taking only some of its columns. The non-complex eigenvalues of this square eigenvalue problem contain all the intersection points between the given line and the curve. As will be discussed further below, some of these non-complex eigenvalues may not be actual intersection points, but they can easily be identified.A non-complex eigenvalue of (16), or its corresponding square eigenvalue problem, gives the potential intersection point . The respective unknown parameter value on the curve satisfies the equation , which is a hard to solve nonlinear problem. According to (15) and (16), however, the left null vector corresponding is proportional to the vector , which is exploited to determine . More precisely, if there is a single parameter value the ratio of any two consecutive components yields
(17) |
It is assumed here that the monomials in are sorted in increasing order so that the ratio of two consecutive entries is simply . If there are parameter values corresponding to the eigenvalue , then the corresponding left null vector space is given by a matrix with rows. Keeping the previous assumption on the ordering of the monomials in , we define the matrices by taking the columns to of . Then, the parameter values are obtained by solving a generalised eigenvalue problem
(18) |
2.2 Illustrative example
The intersection of a cubic Lagrange curve with a line depicted in Figure 6 is considered next. The cubic curve interpolating the points , , and is expressed in power basis with
(19) |
The intersections between the given line
(20) |
and the curve are sought.
As discussed, the degree of the auxiliary polynomial needs to satisfy . In this example choosing yields matrices and that have the dimensions and it is straightforward to compute the eigenvalues of the generalised eigenvalue problem (16). To illustrate the more challenging case with non-square matrices and , which turns out to be inevitable in the case of surfaces, we choose here such that
(21) |
The right null vectors that define the moving lines are obtained from (8) with and the components
(22) |
This matrix has the dimensions , its rank is seven and the dimension of its right null space is five. Each of the corresponding five independent null vectors defines one moving line. Figure (b)b depicts the first four moving lines , , and at the parameter values and .
Substituting the line (20) as in (15) gives the two non-square matrices (with four significant digits)
(23) |
and
(24) |
One approach to obtaining the generalised eigenvalues of matrices and is to use pencil reduction, see [12, 13, 14]
, which is not a widely used linear algebra operation and may introduce additional numerical issues because of several numerical rank estimations. Alternatively, the eigenvalue problems defined by square submatrices
and of the largest size, e.g. the first four columns of and , can be considered(25) |
Although there can only be three intersection points for a cubic curve this problem has four eigenvalues and eigenvectors. Notice that each column of , i.e , represents the intersection between the line and the moving line . Taking a different couple of submatrices and of largest size, like the last four columns of and , yields a different set of intersection points. Both sets of intersection points contain the three true intersection points in addition to one fictitious intersection point.
To demonstrate this observation, the intersection points with two different pairs of largest square submatrices and are computed. The eigenvalues for the problem defined by the first four columns of and are
(26) |
and the eigenvalues for the problem defined by the last four columns and are
(27) |
It is evident that for the first problem and for the second problem correspond to fictitious intersection points while the other three eigenvalues correspond to the true intersection points, see Figure 6. The coordinates of the three true intersection points are computed by introducing in the line equation (20) yielding
(28) |
The fictitious points can also be detected without computing several eigenvalue problems and comparing their eigenvalues. This is accomplished by determining the parametric coordinates of the intersection points on the curve . As indicated in (17) the parametric coordinates are computed using the eigenvectors. As an example, consider the eigenvalue problem (25) defined by first four colums of and . Its eigenvalues are given on (26) and the coordinates of the intersection points are
(29) | ||||
The respective parameters of the intersection points are
(30) |
A point is a true intersection point if and only if ; otherwise, it is a fictitious point. It can easily be found such that the fourth point is not an intersection point. Thus, in practice the intersection points are computed from a single couple of square matrices and .
3 Intersection of lines with surfaces
The extension of the introduced method to surfaces is straightforward. Let be a parametric surface, with , of bi-degree given either in Lagrange basis or power basis with
(31) |
where and are multi-indices, are the parametric surface coordinates, and and are the coefficients in the two basis. Usually, in finite element applications the degrees and of the surface are the same.
In line with the curve case, a point on the surface is defined as the intersection of several moving planes of the form
(32) |
where is an auxiliary vector collecting the parameters of the plane. Although is now a plane instead of a line, it is still denoted with the same symbol to keep the notation simple. The parameters are assumed to be of the following form
(33) |
The bi-degree of the power basis has to be sufficiently high to describe all the intersection points of the surface with a line.
The planes describing the surface have to satisfy
(34) |
which can be rearranged to
(35) |
with a new power basis of bi-degree and the array containing the sorted components of the vectors . The matrix has rows and columns. It has at least right null vectors, i.e. the difference between the number of columns and rows. The surface has the algebraic degree , which is equal to its number of intersections with a line. In order to obtain all the intersections the following condition has to be satisfied
(36) |
As in the curve case, the degree of along the direction can be chosen as yielding
(37) |
By symmetry, it is also valid to choose and . After the right null vectors of are computed the subsequent steps in computing the intersections are identical to the curve case.
Finally, for the sake of completeness, we mention that the introduced intersection algorithm also applies to triangular finite elements. If is a triangular parametric surface of degree then the degree of the auxiliary vector has to be chosen to satisfy ; see also [9, §3].
4 Conclusions
The introduced method is able to determine in one-shot all the intersection points between a line and a curve or surface. It can be applied to curves or surfaces given in any polynomial basis, like the Lagrange or Bernstein, after a straightforward conversion to the power basis. Aspects requiring further research include cases for which the eigenvalue problem (18) is degenerate and the preconditioning of the eigenvalue problem (25). In the accompanying implementation several choices have been considered for both.
References
- [1] Hughes TJR, Cottrell JA, Bazilevs Y. Isogeometric analysis: CAD, finite elements, NURBS, exact geometry and mesh refinement. Computer Methods in Applied Mechanics and Engineering 2005; 194:4135–4195.
- [2] Karniadakis G, Sherwin S. Spectral/hp element methods for computational fluid dynamics. Oxford University Press, 2013.
- [3] Wriggers P, Zavarise G. Computational contact mechanics. Encyclopedia of computational mechanics 2004; .
- [4] Xie ZQ, Sevilla R, Hassan O, Morgan K. The generation of arbitrary order curved meshes for 3d finite element analysis. Computational Mechanics 2013; 51:361–374.
- [5] Turner M, Peiró J, Moxey D. Curvilinear mesh generation using a variational framework. Computer-Aided Design 2018; 103:73–91.
- [6] Rüberg T, Cirak F. Subdivision-stabilised immersed b-spline finite elements for moving boundary flows. Computer Methods in Applied Mechanics and Engineering 2012; 209–212:266–283.
- [7] Patrikalakis NM, Maekawa T. Shape interrogation for computer aided design and manufacturing. Springer, 2009.
- [8] Sederberg TW. Computer Aided Geometric Design. Computer Aided Geometric Design Course Notes, Brigham Young University, 2012.
- [9] Busé L. Implicit matrix representations of rational Bézier curves and surfaces. Computer-Aided Design 2014; 46:14–24.
- [10] Sederberg TW, Chen F. Implicitization using moving curves and surfaces. SIGGRAPH 1995 Conference Proceedings, ACM, 1995; 301–308.
- [11] Strang G. Introduction to Linear Algebra. Fifth edn., Wellesley-Cambridge Press, 2016.
- [12] Xiao X, Sabin M, Cirak F. Interrogation of spline surfaces with application to isogeometric design and analysis of lattice-skin structures. Submitted for publication. 2018; .
- [13] Busé L, Luu Ba T. Matrix-based Implicit Representations of Rational Algebraic Curves and Applications. Computer Aided Geometric Design 2010; 27(9):681–699.
- [14] Beelen T, Van Dooren P. An improved algorithm for the computation of Kronecker’s canonical form of a singular pencil. Linear Algebra and its Applications 1988; 105:9–65.
Comments
There are no comments yet.