ISOGEOMETRIC COLLOCATION METHODS
F.AURICCHIO *,{,L.BEIR ~A O DA VEIGA y ,T.J.R.HUGHES z ,A.REALI *,jj and
G.SANGALLI x
*Structural Mechanics Department,
University of Pavia
y
Mathematics Department \F.Enriques",
University of Milan z Institute
for Computational Engineering and Sciences,University of Texas at Austin
x
Mathematics Department,University of Pavia
{
auricchio@unipv.it jj alessandro.reali@unipv.it y
lourenco.beirao@unimi.it z
hughes@https://www.doczj.com/doc/d73386925.html, x
giancarlo.sangalli@unipv.it
Received 22October 2009Revised 12February 2010Communicated by F.Brezzi
We initiate the study of collocation methods for NURBS-based isogeometric analysis.The idea is to connect the superior accuracy and smoothness of NURBS basis functions with the low computational cost of collocation.We develop a one-dimensional theoretical analysis,and perform numerical tests in one,two and three dimensions.The numerical results obtained conˉrm theoretical results and illustrate the potential of the methodology.
Keywords :Isogeometric analysis;collocation methods;non-uniform rational B-splines.AMS Subject Classiˉcation:22E46,53C35,57S20
1.Introduction
Isogeometric analysis is a computational mechanics technology based on functions used to represent geometry.1,4à8,10à12,19,20,27,30The idea is to build a geometry model and,rather than develop a ˉnite element model approximating the geometry,to directly use the functions describing the geometry in analysis.
In Computer-Aided Design (CAD)systems used in engineering,Non-Uniform Rational B-Splines (NURBS)are the dominant technology.When a NURBS model is constructed,the basis functions used to deˉne the geometry can be systematically enriched by h -,p -,or k -reˉnement (i.e.smooth order elevation;see Cottrell et al.12)
Mathematical Models and Methods in Applied Sciences Vol.20,No.11(2010)2075à
2107
#c World Scientiˉc Publishing Company DOI:10.1142/S0218202510004878
2075
2076 F.Auricchio et al.
without altering the geometry or its parametrization.This means that adaptive mesh reˉnement techniques can be utilized without a link to the CAD database,in contrast withˉnite element methods.This appears to be a distinct advantage of isogeometric analysis overˉnite element analysis.In addition,on a per-degree-of-freedom basis, isogeometric analysis exhibits superior accuracy and robustness compared withˉnite element analysis.1,8,9,20This is particularly true when a k-reˉned basis is adopted.In this case,the upper part of the discrete spectrum is better behaved,20resulting in better conditioned discrete systems.Therefore,isogeometric analysis appears to o?er several important advantages over classicalˉnite element analysis.
However,the use of high-degree basis function,as in k-reˉnement,raises the issue of an e±cient implementation.Until now,isogeometric analysis has been imple-mented by Galerkin formulations.In this case,the e±ciency issue is related to the numerical quadrature rules which are adopted when assembling the system of equations.While elementwise Gauss quadrature is optimal in theˉnite element context,it is sub-optimal when k-reˉned isogeometric discretizations of Galerkin type are considered.More e±cient special rules have been derived for isogeometric analysis by Hughes et al.21Taking inspiration from that work,we initiate here the study of even more e±cient collocation-based methods.
The structure of the paper is the following.In Sec.2,we present a brief review of NURBS and isogeometric analysis,and we introduce the proposed collocation scheme.In Sec.3,we develop a one-dimensional theoretical analysis of the method, which serves the dual purpose of providing theoretical background and guiding the selection of collocation points.In Sec.4,we present numerical tests on simple elliptic problems,in one,two and three dimensions.We test the accuracy of the method, study the behavior of the discrete eigenspectrum and discuss the performance of the scheme with respect to the choice of collocation points.Finally,in Sec.5,we present an analytical cost comparison between the isogeometric Galerkin method and the isogeometric collocation scheme.We draw conclusions in Sec.6.
2.NURBS-Based Isogeometric Analysis
Non-Uniform Rational B-Splines(NURBS)are a standard tool for describing and modeling curves and surfaces in computer-aided design and computer graphics(see Piegl and Tiller25and Rogers28for an extensive description of these functions and their properties).In this work,we use NURBS as an analysis tool,as proposed by Hughes et al.19The aim of this section is to present a short description of B-splines and NURBS,followed by a simple discussion on the basics of isogeometric analysis and by an introduction to the proposed collocation method.
2.1.B-splines and NURBS
B-splines in the plane are piecewise polynomial curves composed of linear combinations of B-spline basis functions.The coe±cients(B i)are points in the plane,referred to as control points.
A knot vector is a set of non-decreasing real numbers representing coordinates in the parametric space of the curve
f 1?0;...; n tp t1?1
g ;
e2:1T
where p is the order of the B-spline and n is the number of basis functions (and control points)necessary to describe it.The interval ? 1; n tp t1 is called a patch .A knot vector is said to be uniform if its knots are uniformly-spaced and non-uniform otherwise;it is said to be open if its ˉrst and last knots have multiplicity p t1.In what follows,we always employ open knot vectors.Basis functions formed from open knot vectors are interpolatory at the ends of the parametric interval ?0;1 but are not,in general,interpolatory at interior knots.
Given a knot vector,univariate B-spline basis functions are deˉned recursively starting with p ?0(piecewise constants)
N i ;0e T?1if i < i t1
0otherwise :(
e2:2T
For p >1:
N i ;p e T? à i i tp à i N i ;p à1e Tt i tp t1à i tp t1à i t1N i t1;p à1e Tif i < i tp t1
0otherwise ;
8
><>:e2:3T
where,in (2.3),we adopt the convention 0=0?0.
In Fig.1we present an example consisting of n ?9cubic basis functions generated from the open knot vector f 0;0;0;0;1=6;1=3;1=2;2=3;5=6;1;1;1;1g .
If internal knots are not repeated,B-spline basis functions are C p à1-continuous.If a knot has multiplicity k ,the basis is C p àk -continuous at that knot.In particular,when a knot has multiplicity p ,the basis is C 0and interpolates the control point at
0.20.4
0.60.81
00.20.40.60.81ξ
N i ,3
Fig. 1.Cubic basis functions formed from the open knot vector f 0;0;0;0;1=6;1=3;1=2;2=3;5=6;1;1;1;1g :
Isogeometric Collocation Methods 2077
that location.We deˉne
^S
n
?span f N i;pe T;i?1;...;n g:e2:4TBy means of tensor products,a multi-dimensional B-spline region can be con-structed.We discuss here the case of a two-dimensional region,the higher-dimensional case being analogous.Consider the knot vectors f 1?0;...; ntpt1?1g and f 1?0;...; mtqt1?1g,and an n?m net of control points B i;j.One-dimensional basis functions N i;p and M j;q(with i?1;...;n and j?1;...;m)of order p and q, respectively,are deˉned from the knot vectors,and the B-spline region is the image of the map S:?0;1 ??0;1 !" given by
Se ; T?
X n
i?1X m
j?1
N i;pe TM j;qe TB i;j:e2:5T
The two-dimensional parametric space is the domain?0;1 ??0;1 .Observe that the two knot vectors f 1?0;...; ntpt1?1g and f 1?0;...; mtqt1?1g generate a mesh of rectangular elements in the parametric space in a natural way.Analogous to (2.4),we deˉne
^S
nm
?span f N i;pe TM j;qe T;i?1;...;n;j?1;...;m g:e2:6TIn general,a rational B-spline in R d is the projection onto d-dimensional physical space of a polynomial B-spline deˉned in(dt1)-dimensional homogeneous coordi-nate space.For a complete discussion see the book by Farin17and references therein. In this way,a great variety of geometrical entities can be constructed and,in par-ticular,all conic sections in physical space can be obtained exactly.The projective transformation of a B-spline curve yields a rational polynomial curve.Note that when we refer to the\degree"or\order"of a NURBS curve,we mean the degree or order, respectively,of the polynomial curve from which the rational curve was generated.
To obtain a NURBS curve in R2,we start from a set B!i2R3(i?1;...;n)of control points(\projective points")for a B-spline curve in R3with knot vector?. Then the control points for the NURBS curve are
?B i k??B!i k
!i
;k?1;2;e2:7T
where?B i k is the k th component of the vector B i and!i??B!i 3is referred to as the i th weight.The NURBS basis functions of order p are then deˉned as
R p
i e T?
N i;pe T!i
P
n
^i?1
N^i;p^i
:e2:8T
The NURBS curve is deˉned by
Ce T?
X n
i?1R p ie TB i:e2:9T
2078 F.Auricchio et al.
Analogously to B-splines,NURBS basis functions on the two-dimensional para-metric space^ ??0;1 ??0;1 are deˉned as
R p;q i;je ; T?
N i;pe TM j;qe T!i;j
P
n
^i?1
P
m
^j?1
N^i;pe TM^j;qe T!^i;^j
;e2:10T
where!i;j?eB!i;jT3.Observe that the continuity and support of NURBS basis functions are the same as for B-splines.
NURBS regions,similarly to B-spline regions,are deˉned in terms of the basis functions(2.10).In particular,we assume from now on that the physical domain is a NURBS region associated with the n?m net of control points B i;j,and we introduce the geometrical map F:^ !" given by
Fe ; T?
X n
i?1X m
j?1
R p;q
i;j
e ; TB i;j:e2:11T
For the sake of simplicity,the extension to three(or more)dimensions is not discussed here(but it is completely analogous).
2.2.Basics of isogeometric analysis
In one-space dimension,let F be the parametrization of the physical interval?a;b , which is assumed piecewise smooth with piecewise smooth inverse.Following the isoparametric approach,the space of NURBS functions on?a;b is deˉned as the span of the push-forward of the basis functions(2.8)
V n?span f R p
i
Fà1;i?1;...;n g:e2:12T
Similarly,in two dimensions,let F be the parametrization of the physical domain ,as given in(2.11).We assume that the Jacobian of F and its inverse are non-singular.The image of the elements in the parametric space are elements in the physical space.The physical mesh is therefore made of elements Fee i; it1T?e j; jt1TT,with i?1;...;ntp;j?1;...;mtq.We denote by h the mesh-size, that is,the maximum diameter of the elements of the physical mesh.
Again,the space of NURBS functions on is deˉned as the span of the push-forward of the basis functions(2.10)
V nm?span f R p;q
i;j
Fà1;with i?1;...;n;j?1;...;m g:e2:13T
The numerical solution of a boundary-value problem is sought in the space V n or V nm,endowed with suitable boundary conditions.
The interested reader mayˉnd more details on isogeometric analysis as well as many interesting applications in a number of recently published papers.6,11,12,10,19,27
Isogeometric Collocation Methods2079
2.3.Isogeometric collocation method
In the ˉrst part of this section we present the concept of isogeometric collocation method in a general setting,while in the second part we give a more detailed description for a particular case.
Consider the following boundary-value problem
D u ?f in ;
G u ?g on @ ;(
e2:14T
where the solution u : !R ,D represents a scalar di?erential operator,G u :@ !R r is a vector operator (r 2N t)representing boundary conditions,and f ;g are given data.Note that both operators D and G may depend on the position x 2 [@ .
We introduce a ˉnite set of collocation points f ^ #g #2I in ^ ??0;1 ??0;1 ,divi-ded into two distinct sets I D and I G ,such that I ?I D [I G ,and
#2I G )^ #2@^
;#I D tr á#I G ?dim eV nm T;
where #indicates the cardinality of a set.Note that the collocation points ^ #are
typically built in a tensor product fashion.The choice of the collocation points is crucial for the stability and good behavior of the discrete problem.Finally let #?F e^ #T,for all #2I .
Then,the isogeometric collocation approximation of (2.14)reads:ˉnd u nm 2V nm such that
D u nm e #T?f e #T8#2I D ;
G u nm e #T?g e #T8#2I G :(
e2:15T
As a particular example,we consider the case where the collocation points are chosen using the Greville abscissae 14,15of the knot vectors.To simplify the pre-sentation,we restrict our description to the case of a second-order operator,still indicated by D ,with boundary conditions represented by the scalar operator G u :@ !R .
Let "
i ,i ?1;...;n ,be the Greville abscissae related to the knot vector f 1;...; n tp t1g :
"
i ? i t1t i t2tááát i tp p
:e2:16T
Analogously,let "
j ,j ?1;...;m ,be the Greville abscissae related to the knot vector f 1;...; m tq t1g .It is easy to see that "
1?" 1?0," n ?" m ?1,while all the remaining points are in e0;1T.We deˉne the collocation points ij 2 by the tensor product structure
ij ?F e^ ij T;
^ ij ?e"
i ;" j T2^ ;for i ?1;...;n ;j ?1;...;m .
2080 F.Auricchio et al.
Then,the isogeometric collocation problem with Greville abscissae reads:ˉnd u nm 2V nm such that D u nm e ij T?f e ij Ti ?2;...;n à1;j ?2;...;m à1;
G u nm e ij T?g e ij Tei ;j T2f 1;n g ?f 1;...;m g [f 1;...;n g ?f 1;m g :(
e2:17T
We wish to remark that collocating at Greville abscissae as in the example above is just one of the possible options.In the following,we will also consider a couple of other interesting choices (see Sec.3.3).3.Theoretical Analysis
This section is devoted to the theoretical analysis of the NURBS-based collocation approach.We ˉrst present an abstract framework,which is based on the results of De Boor and Swartz.13Within this framework,we then present a theoretical analysis of the one-dimensional model problem.Finally,we discuss the di±culties that are encountered in the theoretical analysis of the multi-dimensional case.
Given any A ;B 2R ,we write A .B whenever A CB with the positive constant C depending only on the problem data.We deˉne the symbol &analogously,and write A ’B whenever both A .B and A &B hold,in general with di?erent constants.3.1.The abstract problem and its discretization
Taking the initial steps from De Boor and Swartz,13the main idea of the proof is to rewrite the di?erential operator (of order n )D as a perturbation of the n th derivative operator.As it will be clearer in the following,this can be achieved using a Green function.13
Given a Banach space W ,endowed with the norm jj ájj W ,a linear operator T :W !W ,and datum 2W ,we consider the following abstract problem:Find y 2W such that
eI tT Ty ? :
e3:1T
The following assumption of well-posedness is made.
Assumption 3.1.There exists one and only one solution to problem (3.1)and
jjeI tT Tz jj W ’jj z jj W
8z 2W :
e3:2T
Observe that (3.2)implies the continuity of T .In consideration of the eventual discretization of (3.1),we also require the next assumption.
Assumption 3.2.T is compact,and in particular there exists a compact subspace W tof W such that T is continuous from W to W t,that is,jj Tz jj W t.jj z jj W 8z 2W .Given an n -dimensional space W n &W ,and a linear projector Q n :W !W n ,we consider the following discretization of (3.1):Find y n 2W n such that
Q n eI tT Ty n ?Q n :
e3:3TIsogeometric Collocation Methods 2081
2082 F.Auricchio et al.
Observe that Q neItTTy n?y ntQ n Ty n.In the following sections the operator Q n will denote interpolation.We make the following two assumptions on Q n and W n. Assumption3.3.There exists a constant C s>0,independent of n,such that
jj Q n z jj W C s jj z jj W,for all z2W.
Assumption3.4.With the notation of Assumption3.2,
jj zàz n jj W c n jj z jj Wt8z2Wt;e3:4T
min
z n2W n
and the positive constants c n!0as n!1.
With these assumptions we have the following result.
Lemma3.1.There exists"n2N,depending only on the operator T and the problem datum,such that the operatoreItQ n TT:W!W is an isomorphism for all n!"n. Moreover,for n!"n
jjeItQ n TTz jj W’jj z jj W8z2W:e3:5TProof.Let z2W.Adding and subtracting any z n2W n,using that Q n z n?z n and the continuity of Q n gives
jjeIàQ nTz jj W jj zàz n jj Wtjj z nàQ n z jj W
?jj zàz n jj Wtjj Q nezàz nTjj W
e1tC sTjj zàz n jj W:e3:6TSelecting z n to be the minimizer in(3.4),the above bound gives
jjeIàQ nTz jj We1tC sTc n jj z jj Wt:e3:7TAdding and subtracting Tz,using the triangle inequality,and applying(3.7)yields jjeItQ n TTz jj W jjeItTTz jj WtjjeIàQ nTTz jj W
jjeItTTz jj Wte1tC sTc n jj Tz jj Wt:e3:8TBound(3.8)combined with(3.2)and Assumption3.2gives
jjeItQ n TTz jj W.jj z jj Wte1tC sTc n jj z jj W;e3:9Twhich,due to Assumption3.4,gives
jjeItQ n TTz jj W.jj z jj W:e3:10TVery similar steps as the ones above give
jjeItQ n TTz jj W&jj z jj Wàe1tC sTc n jj z jj W;e3:11Twhich for n su±ciently large yields
jjeItQ n TTz jj W&jj z jj W:e3:12TBound(3.5)is proved by(3.10)and(3.12).
Therefore,the operator eI tQ n T Tis one-to-one and its restriction eI tQ n T T:W n !W n is an isomorphism.The following corollary is immediate.Corollary 3.1.For n !"n
,Problem (3.3)has a unique solution .We can now derive the main result of this section.
Proposition 3.1.Let Assumptions 3.1à3.4hold .Let y be the solution of (3.1)and
y n the solution of (3.3).Then
jj y ày n jj W .e1tC s Tmin z n 2W n
jj y àz n jj W
for all n su±ciently large .
Proof.Applying the operator Q n to both sides of (3.1)and subtracting (3.3),it follows that
Q n eI tT Tey ày n T?0:
e3:13T
First using Lemma 3.1,then applying (3.13)we infer
jj y ày n jj W .jjeI tQ n T Tey ày n Tjj W ?jjeI àQ n Tey ày n Tjj W :
e3:14T
From (3.14),noting that Q n y n ?y n and using (3.6),we get
jj y ày n jj W .jjeI àQ n Ty jj W e1tC s Tmin z n 2W n
jj y àz n jj W :
e3:15T
3.2.Theoretical analysis in one dimension
In this section we restrict our attention to the one-dimensional case of a second-order di?erential equation and prove that,if the collocation points are chosen suitably,the collocation method converges with optimal rate.3.2.1.The continuous problem
Let f ,a 0,a 1,be real functions in C 0?a ;b ,with a
u 00
ex Tta 1ex Tu 0ex Tta 0ex Tu ex T?f ex T8x 2ea ;b T
BC i eu T?g i i ?1;2;(e3:16T
where u 0;u 00represent the ˉrst and second derivatives of u ,respectively.We will
indicate the derivative operator of order i also as D i ,i 2N .We assume that (3.16)has one and only one solution u ,and that the boundary condition operators BC i are linearly independent on Ker eD 2T,that is,on the space of linear functions.
Under this assumption,the Green's function H 0:?a ;b ??a ;b !R is well deˉned for the problem D 2v ? with datum 2C 0?a ;b and the homogeneous boundary
Isogeometric Collocation Methods 2083
conditions BC ievT?0,i?1;2.As a consequence,all elements of the space
C2BC?f v2C2?a;b j BC ievT?g i for i?1;2g
can be written uniquely as
vexT?u0exTtZ
b
a
H0ex;tTD2vetTd t;e3:17T
where u0is the unique linear function on?a;b which satisˉes BC ieu0T?g i,i?1;2.
Introducing the partial derivatives of H0
H kex;tT?
@
@x
k
H0ex;tT8ex;tT2?a;b ??a;b ;k?1;2;
and deˉning y?D2u,using(3.17)one gets
D k uexT?D k u0exTtZ
b
a
H kex;tTyetTd t8x2?a;b ;k?0;1;2:e3:18T
Introducing the operator T:C0?a;b !C0?a;b as
TzexT?
X1
i?0a iexT
Z
b
a
H iex;tTzetTd t8z2C0?a;b ;x2?a;b ;e3:19T
and deˉning2C0?a;b as
exT?fexTà
X1
i?0
a iexTD i u0exT8x2?a;
b ;e3:20T
the problem(3.16)can then be rewritten in the abstract framework of the previous section:substituting(3.18)into(3.16),and using the deˉnitions(3.19)à(3.20),we obtain the problem ofˉnding y2C0?a;b such that
eItTTy? :e3:21T
Indeed,the new problem(3.21)is by construction equivalent to Problem(3.16)up to the transformations
u?u0tZ
b
a
H0eá;tTyetTd t;y?D2u:
For the analysis of(3.16)we select the C0-topology.We will denote by jjájj C0?jjájj L1the usual L1-norm and by
jáj W i;1?jj D iájj L1;
jjájj W k;1?
X k
i?0jj D iájj L1;
2084 F.Auricchio et al.
the usual Sobolev seminorms and norms.Having the correspondence
W ?C 0?a ;b ?f D 2v j v 2C 2?a ;b ;v ea T?v eb T?0g ;
e3:22T
Assumption 3.1is equivalent to the well-posedness of (3.16):given f 2C 0?a ;b there exists a unique solution u 2C 2?a ;b to (3.16),with continuous dependence.This is assumed to hold true.Furthermore,Assumption 3.2is veriˉed,as shown in the next lemma.
Lemma 3.2.The operator T :C 0?a ;b !C 0?a ;b in (3.19)is compact .In particular it holds
jj Tz jj W 1;1?a ;b .jj z jj C 0?a ;b
8z 2C 0?a ;b :
Proof.Given any z 2C 0?a ;b ,let v 2C 2?a ;b be given by
v ?
Z b
a
H 0eá;t Tz et Td t ;so that D 2v ?z .Note that,due to the homogeneous boundary conditions satisˉed
by v ,
jj D i v jj L 1?a ;b .jj D 2v jj L 1?a ;b jj D 2v jj C 0?a ;b
i ?0;1;2:
e3:23T
We then easily have,by deˉnition of T ,v and using (3.23)
jj Tz jj W 1;1?a ;b ?X 1
i ?0a i
D i v
W
1;1?a ;b
.
X 1i ?0
jj D i v jj W 1;1?a ;b
.jj D 2
v jj C 0?a ;b ?jj z jj C 0?a ;b :
e3:24T
Finally,the compactness of the operator follows from the compact inclusion of W 1;1?a ;b into L 1?a ;b .
3.3.The discretized di?erential problem
Given n 2N ,let V n t2&C 2?a ;b be a NURBS space of dimension n t2on the
interval ?a ;b ,associated with a spline space ^S
n t2&C 2?0;1 on the parametric interval ?0;1 .Recalling the standard assumptions on the one-dimensional geo-metrical map F of Sec.2.2,we now consider DF >0on the parametric domain ?0;1 .Given for all n 2N , 1< 2<ááá< n assigned collocation points in ?a ;b ,we con-sider the following discrete problem:Find u n 2V n t2such that
u 00
n e j Tta 1e j Tu 0n e j Tta 0e j Tu n e j T?f e j Tj ?1;...;n
BC i u n ?g i
i ?1;2:(e3:25T
Isogeometric Collocation Methods 2085
Notice that in one dimension the boundary conditions are imposed exactly by the collocation scheme.Following the notation of the abstract setting presented in Sec.3.1, we now denote as W n the space
W n?f D2v n j v n2V nt2;veaT?vebT?0g;e3:26Twhich is n-dimensional.
Remark 3.1.Since the isoparametric paradigm(2.13)is applied for the con-struction of V nt2,it follows that the space of global linear polynomials is contained in V nt2.As this subspace is the kernel of the operator D2,the space W n can equivalently be written as
W n?f D2v n j v n2V nt2g:
For the continuous di?erential problem it is easy to see that the above discrete problem is equivalent toˉnding y n2W n such that
?eItTTy n e jT?e jTj?1;...;ne3:27Tup to the transformations
u n?u0tZ
b
a
H0eá;tTy netTd t;y n?D2u n;
where u0is,as in(3.17),the linear lifting of the boundary conditions.
Problem(3.25)is what is solved in practice,while problem(3.27)is what we use for the theoretical analysis of the method.Denoting by Q n:C0?a;b !W n the interpolation operator
eQ nezTTe jT?ze jT8z2C0?a;b ;j?1;...;n;
problem(3.27)ˉts within the framework of Sec.3.1.Then,the optimality of the collocation scheme(3.25)follows from Proposition3.1,if Assumptions3.3and3.4 hold true.The approximation properties of NURBS spaces with respect to the mesh-size h have been studied,in the L2setting,by Bazilevs et al.6Assumption3.4follows from the analogous result in the L1setting,which is stated in the next proposition. Theorem3.1.There exists a constant C a>0,dependent on the geometry map F,
weight!,and on the degree p of the NURBS space V n,such that
min
w n2W n
jj wàw n jj L1?a;b C a h l jj w jj W l;1?a;b ;8w2W l;1?a;b ;e3:28Twith0l pà1.
Proof.A simple extension to the L1setting of the results in Bazilevs et al.6for the L2framework gives
min
v n2V n;0
jj v00àv00n jj L1?a;b Ch l jj v jj W lt2;1?a;b 8v2W lt2;1?a;b ;veaT?vebT?0;
e3:29T2086 F.Auricchio et al.
where V n;0is the subspace of all functions in V nt2which are zero at a and b.Because of(3.22)and(3.26),any w2W l;1?a;b can be written as v00with v2W lt2;1?a;b and veaT?vebT?0.By deˉnition,any w n2W n can be written as v00n with v n2V n;0. Therefore,ˉrst using the above observations and(3.29),and then,since v is null at the boundaries,a Poincar e inequality,we obtain
min w n2W n jj wàw n jj L1?a;b ?min
v n2V n;0
jj v00àv00n jj L1?a;b
Ch l jj v jj W lt2;1?a;b
Ch l jj w jj W l;1?a;b :e3:30T
Theorem 3.1implies Assumption 3.4.The following theorem gives su±cient conditions to guarantee the stability of Q n,for n large enough,i.e.Assumption3.3.
Theorem 3.2.Given^C s>0,there exists"n2N such that,for any n>"n the following holds.Let
^T
n
?f D2s n j s n2^S nt2g&C0?0;1 ;
which is n-dimensional,and let^Q n:C0?0;1 !^T n be the interpolation operator on n points^ 1<^ 2<ááá<^ n in?0;1 .Let the stability condition jj^Q n w jj L1 ^C
s
jj w jj L1;8w2C0?0;1 ,hold.Then the interpolation operator Q n:C0?a;b !W n associated with the points j?Fe^ jTis stable,that is jj Q n w jj L1C s jj w jj L1;8w2 C0?a;b .The constants"n and C s depend only on^C s and the geometry para-metrization F:?0;1 !?a;b .
Proof.Let v?R
b
a
H0eá;tTwetTd t.Note that by deˉnition veaT?vebT?0and
w?D2v.Then,the interpolation problem on?a;b ,that is,the problem ofˉnding w n?Q n w2W n such that
w ne jT?we jT;1j n;
is equivalent toˉnding v n2V n such that
eD2v nTe jT?eD2vTe jT;1j n;v neaT?v nebT?0;e3:31Twith w n?D2v n.
Consider the function of^x2?0;1
s ne^xT?v neFe^xTT!e^xT;e3:32Twhere!e^xTis the denominator in the deˉnition(2.8)of the NURBS,i.e.
!e^xT?
X n
^i?1
N^i;pe^xT!^i:e3:33TObserve that by deˉnition s n2^S nt2,and that s ne0T?s ne1T?0.
Isogeometric Collocation Methods2087
Now we want to write the interpolation problem(3.31)in terms of s n.Denoting by G?Fà1the inverse of the geometric parametrization F,we have for all x2?a;b
v nexT?s neGexTT
!eGexTT
:e3:34T
In the following we denote by^D i the i-derivative with respect to^https://www.doczj.com/doc/d73386925.html,puting the derivatives in(3.32)and using a Poincar e inequality yields
jj^D2s n jj L1?0;1 .jj v n jj W2;1ea;bT.jj D2v n jj L1ea;bT;e3:35Twhere the hidden constant in(3.35)depends on F and!.
Substituting(3.32)in(3.31),using change of variable and the chain rule for derivatives,in(3.34),we get
X2
i?0
b ie^ jTe^D i s nTe^ jT?eD2vTe jT;1j ne3:36Twith suitable factors b ie^xT.In particular,the leading coe±cient is
b2e^xT?eDGeFe^xTTT2
!e^xT
:
Due to the assumption that^DFe^xT>0for all x2?0;1 ,the inverse mapping satisˉes DGexT>0for all x2?a;b .Therefore b2e^xT!c>08^x2?0;1 ,and we can divide (3.36)by it and obtain
e^D2s nTe^ jTt
X1
i?0c ie^ jTe^D i s nTe^ jT?
1
b2e^ jT
eD2v FTe^ jT;j?1;2:e3:37T
The above problem can be written in the form
^Q n eIt^TTy n?^Q n
eD2vT F
b2
;e3:38T
with^T:C0?0;1 !C0?0;1 a compact operator and^D2s n?y n,s ne0T?s ne1T?0, by following the steps of Sec.3.2.1.
By deˉnition and by the same derivative transformation used above,we have
jj w n jj L1ea;bT?jj D2v n jj L1ea;bT
’^D2s nt
X1
i?0c ie;^D i s nT
L1e0;1T
?jjeIt^TTy n jj L1e0;1T:e3:39TA triangle inequality and the fact that^Q n y n?y n,give
jjeIt^TTy n jj L1e0;1TjjeIt^Q n^TTy n jj L1e0;1TtjjeIà^Q nT^T y n jj L1e0;1T:e3:40T2088 F.Auricchio et al.
Equation(3.38)and the stability of the spline interpolant^Q n yield immediately
jjeIt^Q n^TTy n jj L1e0;1T?^Q n D 2v F b2
L1e0;1T
.D 2v F b2
L1e0;1T
.jj D2v jj L1ea;bT?jj w jj L1ea;bT:e3:41TUsing the L1stability of the interpolation operator^Q n and standard approximation properties for splines14it follows easily
jjeIà^Q nTjj L1e0;1T.c n jj jj W1;1e0;1T82W1;1e0;1T;
with lim n!1c n?0.As a consequence,the second term in(3.40)can be bounded using the compactness of^T and recalling(3.35).One gets
jjeIà^Q nT^T y n jj L1e0;1T.c n jj^T y n jj W1;1e0;1T.c n jj y n jj L1e0;1T
?c n jj^D2s n jj L1e0;1T.c n jj D2v n jj L1ea;bT
?c n jj w n jj L1ea;bT:e3:42TCombining(3.39)with(3.40)à(3.42)it follows that,for n su±ciently large,
jj w n jj L1ea;bT.jj w jj L1ea;bT:e3:43T
Combining the results of this section with Proposition3.1yields
Theorem3.3.Given a NURBS space V nt2on?a;b and a set of collocation points f 1;...; n g2?a;b ,let the hypotheses of Theorem3.2on the collocation points and the spline space^S nt2be satisˉed.Then
jj D2euàu nTjj L1?a;b .h l jj u jj W lt2;1?a;b ;e3:44Twhere u is the solution of the continuous problem(3.16)and u n is the solution of the collocation problem(3.25).
The NURBS collocation method is therefore optimal in the W2;1seminorm.The convergence in the full W2;1norm can be recovered using a Poincar e inequality
jj uàu n jj W2;1?a;b .h l jj u jj W lt2;1?a;b :e3:45TTherefore,in order to guarantee the optimal behavior of the collocation scheme (3.25),the collocation points can be selected according to Theorem3.2.Notice that the space^T n,which is made of second derivatives of the splines in^S nt2,is still a space of splines.In particular,^T n is the space of splines of degree pà2on the knot vector obtained by removing theˉrst two and last two knots from the knot vector of^S nt2. The optimal selection of points for interpolation of one-dimensional splines is addressed in various papers.The only choice which is proved to be stable for any mesh and degree is the one proposed by Demko.15These points are referred to as
Isogeometric Collocation Methods2089
2090 F.Auricchio et al.
Demko abscissae.They are the extrema of the Chebyshev splines,i.e.the splines that have the maximum number of oscillations,for which the extrema take the values?1. The Demko abscissae are obtained by an iterative algorithm(see,e.g.chbpnt in MATLAB and also the book by de Boor14).A di?erent approach which is proposed in the engineering literature23is to collocate at the Greville abscissae.These are obtained as knots averages14,15(see(2.16))and Greville abscissae interpolation is proved to be stable up to degree3,while there are examples24of instability for degrees higher than19on particular non-uniform meshes.a In Sec.4.1,we numerically test both possibilities,i.e.we consider:
.collocation at the images(though the geometry map F)of the Demko abscissae of the space^T n(we refer to these as\second-derivative Demko abscissae");
.collocation at the images of the Greville abscissae of the space^S nt2,after removing theˉrst and last knots(we refer to these as simply\Greville abscissae").
In both cases,boundary conditions will be imposed strongly on the discrete solution, according to formulation(3.25).It will be shown that theˉrst choice(Demko abscissae)is superior,even though it represents a very small advantage over the second choice in many cases of practical interest.
In Fig.2we compare a sample of Demko abscissae for the space^T n,Greville abscissae for^T n,and Greville abscissae for^S nt2after removing theˉrst and last knots.We consider a uniform open knot vector with polynomial degree p?3in the ˉrst frame,with p?7in the second frame,and a non-uniform open knot vector with knot spans forming a geometric sequence with p?7in the last frame.Coordinates of collocation points are reported in Tables1à3.
3.4.Remarks on the theoretical analysis in higher dimensions
The extension to the multi-dimensional case of the analysis of Sec.3.2,based on the abstract framework of Sec.3.1,poses various di±culties.One of the main issues is that the space consisting of the Laplacian of multivariate splines,or more generally of second-order derivatives of multivariate splines,is not a complete space of splines of a given degree.Therefore,stable interpolation in this space is di±cult to determine. Thus,we only consider numerical experiments in Sec.4.
Other approaches for the theoretical analysis of the collocation method with spline spaces in higher dimensions are given,for particular cases,in Arnold et al.3and Prenter.26
4.Numerical Tests
In this section we present some numerical results for the collocation approach pro-posed.Weˉrst study a one-dimensional second-order di?erential equation with a To be precise,we refer to cases where the lengths of consecutive knot spans are increased by a constant scaling factor,forming a geometric sequence.
?0.2
?0.2?0.2
0.2
0.4
0.6
0.8
1
1.2
00.511.5200.20.40.60.81 1.2
00.511.520
0.2
0.4
0.6
0.8
1
1.2
https://www.doczj.com/doc/d73386925.html,parison of Demko abscissae for ^T
n (ˉrst line of each frame),Greville abscissae for ^T n (second line),and Greville abscissae for ^S n t2after removing the ˉrst and last knots (third line).Uniform open knot
vector with p ?3(top frame),with p ?7(middle frame),and non-uniform open knot vector with knot spans forming a geometric sequence with p ?7(bottom frame).
Table 1.Coordinates of collocation points in Fig.2(upper frame):uniform mesh,p ?3.
Demko abscissae for ^T n 0
0.25000.50000.7500 1.0000Greville abscissae for ^T n 0
0.25000.50000.7500 1.0000Greville abscissae for ^S
n t20.0833
0.2500
0.5000
0.7500
0.9167
Isogeometric Collocation Methods 2091
Dirichlet and Neumann boundary conditions,for linear and nonlinear parametriza-tions,either using Greville or second-derivative Demko abscissae (see Sec.3.3).We also explore the 1D spectral approximation properties of the method.We then move to higher dimensions,focusing on an elliptic model problem deˉned on a quarter of an annulus in 2D and on a cube in 3D.
4.1.One-dimensional source problem with Dirichlet
boundary conditions
The ˉrst test case we study is the following source problem deˉned on the domain ?0;1 :
àu 00
tu 0tu ?e1t4 2Tsin e2 x Tà2 cos e2 x T;8x 2e0;1T;
u e0T?u e1T?0; e4:1T
which admits the exact solution:
u ?sin e2 x T:
e4:2T
This problem is numerically solved using the collocation method outlined in the previous sections,considering a linear parametrization (F ?Id )and employing both Greville and second-derivative Demko abscissae.In Figs.3à5,we report log-scale plots of the relative errors for di?erent degrees of approximations in L 1-,W 1;1-and W 2;1-norms,respectively,for both choices of abscissae,recalling that
j áj W i ;1?jj D i ájj L 1
and
jj ájj W k ;1?
X k i ?0
jj D i ájj L 1:
Denoting by p the degree of the approximation,the ˉgures show that in the ˉrst two
norms an order of convergence p is attained for even degrees,while an order p à1is attained for odd degrees.Instead in the W 2;1-norm,we observe the expected optimal
Table 2.Coordinates of collocation points in Fig.2(middle frame):uniform mesh,p ?7.
Demko abscissae for ^T n 00.04290.15690.31440.50000.68560.84310.9571 1.0000Greville abscissae for ^T
n 00.05000.15000.30000.50000.70000.85000.9500 1.0000Greville abscissae for ^S
n t20.03570.10710.21430.35710.50000.64290.78570.89290.9643Table 3.Coordinates of collocation points in Fig.2(bottom frame):knot spans forming a geometric
sequence with p ?7.Demko abscissae for ^T n 00.00720.03360.10220.24380.45680.71340.9208 1.0000
Greville abscissae for ^T n 00.00500.02500.09000.29000.49000.68500.8650 1.0000Greville abscissae for ^S
n t20.00360.01790.06430.20710.35000.49290.63570.77500.9036
2092 F.Auricchio et al.
??????log 10
( n )
l o g 10
( ||u
e x
u h ||L ∞ / ||u e x ||L ∞)
12
?????log 10
( n )
l o g 10( ||u
e x
u h ||L / ||u e x ||
L
∞
∞??)
Fig.3.1D source problem with Dirichlet boundary conditions (linear parametrization),using Greville (left)and second-derivative Demko (right)abscissae.Relative error in L 1-norm for di?erent degrees of approximation.
??????log 10
( n )
l o g 10
( |u
e x
u h |
W
1,∞
/ |u e x |
W
1,∞
)
??????log 10
( n )
l o g 10( |u
e x
u h |W 1,∞ / |u e x |W 1,∞)
??Fig.4.1D source problem with Dirichlet boundary conditions (linear parametrization),using Greville (left)and second-derivative Demko (right)abscissae.Relative error in W 1;1-norm for di?erent degrees of approximation.
??????log 10( n )
l o g 10( |u
e x
u h |W 2,∞ / |u e x |W 2,∞)
??????log 10( n )
l o g 10( |u
e x
u h |
W 2,∞
?? / |u e x |W 2,∞
)
Fig.5.1D source problem with Dirichlet boundary conditions (linear parametrization),using Greville (left)and second-derivative Demko (right)abscissae.Relative error in W 2;1-norm for di?erent degrees of approximation.
Isogeometric Collocation Methods
2093
order of convergence,i.e.p à1,for all approximation degrees,in agreement with the predictions of the theory (see Sec.3.3).Here,and in the following examples,the use of second-derivative Demko abscissae gives the same convergence rates as Greville,but a very slight di?erence in the magnitude of errors is observed for certain degrees.We have not included results for p ?2,since this case is not covered by the theory.However,the results obtained are in complete agreement with the ones observed for higher degrees (i.e.order p for the L 1-and W 1;1-norms and p à1for the W 2;1-norm).
As proven earlier,the optimal theoretical convergence rates are also expected when nonlinear parametrizations are employed,provided that the geometry map is kept ˉxed during the reˉnement process.In order to test this,we refer to a particular nonlinear mesh (an illustration of which is given in Fig.6)and we solve the source problem with Dirichlet boundary conditions,using standard h -reˉnement procedures such that the geometry map is unchanged.The results are reported in Figs.7à9where it is possible to see the same orders of convergence observed in the case of linear parametrization.In particular,Fig.9shows how the W 2;1-norm errors converge with the expected optimal orders (i.e.p à1).
Finally,we consider a problem which shows a signiˉcant advantage in using second-derivative Demko abscissae,i.e.collocation points that guarantee optimal interpolation on the second derivatives of the discrete space.The di?erential equation is now
u 00
ex T?f ex T8x 2e0;1T;
u e0T?u e1T?0;(e4:3T
and,therefore,the collocation scheme is in fact reduced to a simple interpolation of second derivatives.We compare with the other two choices,i.e.collocation at Greville
0.1
0.2
0.3
0.4
0.50.6
0.7
0.8
0.9
1
00.10.20.30.40.50.60.70.80.91ξ
x
Fig.6.Nonlinear geometry map.
2094 F.Auricchio et al.