Dispersion-minimizing quadrature rules for $C^1$ quadratic isogeometric analysis

TL;DR

This paper introduces optimized two-point quadrature rules for $C^1$ quadratic isogeometric analysis that significantly reduce dispersion errors and computational costs in wave propagation and structural vibration simulations.

Contribution

It presents novel dispersion-minimizing quadrature rules requiring only two points per element, simplifying integration and improving accuracy in isogeometric analysis.

Findings

Two-point quadrature rules reduce dispersion error effectively.

Proposed rules outperform standard quadrature in numerical tests.

Cost savings due to fewer quadrature points.

Abstract

We develop quadrature rules for the isogeometric analysis of wave propagation and structural vibrations that minimize the discrete dispersion error of the approximation. The rules are optimal in the sense that they only require two quadrature points per element to minimize the dispersion error [1], and they are equivalent to the optimized blending rules we recently described. Our approach further simplifies the numerical integration: instead of blending two three-point standard quadrature rules, we construct directly a single two-point quadrature rule that reduces the dispersion error to the same order for uniform meshes with periodic boundary conditions. Also, we present a 2.5-point rule for both uniform and non-uniform meshes with arbitrary boundary conditions. Consequently, we reduce the computational cost by using the proposed quadrature rules. Various numerical examples demonstrate the performance of these quadrature rules.

Equations

$$\Delta u+{\omega }^{2}u$$

$$u$$

$$B(u,v)=0\forall v\in H_0^1(\Omega ),$$

$$B(w,v)=a(w,v)-{\omega }^{2}b(w,v)$$

$$B(u_h,v_h)=0\forall v_h\in V_h.$$

$${\int }_{\hat{K}}\hat{f}(\hat{\mathbf{x}})\text{d}\hat{\mathbf{x}}\approx \sum _{l=1}^{N_q}{\hat{\varpi }}_l\hat{f}(\widehat{n_l}),$$

$${\int }_{K}f(\mathbf{x})\text{d}\mathbf{x}\approx \sum _{l=1}^{N_q}{\varpi }_{l,K}f(n_{l,K}),$$

$${\tilde{B}}_h(u_h,v_h)={\tilde{a}}_h({\tilde{u}}_h,v_h)-{\omega }^2{\tilde{b}}_h({\tilde{u}}_h,v_h)\forall v_h\in V_h,$$

$${\tilde{a}}_h(w,v)=\sum _{K\in {\mathcal{T}}_h}\sum _{l=1}^{N_q}{\varpi }_{l,K}^{(1)}\nabla w(n_{l,K}^{(1)})\cdot \nabla v(n_{l,K}^{(1)}),$$

$${\tilde{b}}_h(w,v)=\sum _{K\in {\mathcal{T}}_h}\sum _{l=1}^{N_q}{\varpi }_{l,K}^{(2)}w(n_{l,K}^{(2)})v(n_{l,K}^{(2)}),$$

$$(\mathbf{K}-{\omega }^2\mathbf{M})\mathbf{U}=0$$

$$U(x)=\sum _{a\in \mathbb{Z}}{U}^a{\varphi }_a(x)$$

$${\tilde{B}}_h(U,v_h)=0.$$

$$\displaystyle\big{(}K_{2}-\Lambda^{2}M_{2}\big{)}(U^{j-2}+U^{j+2})-\big{(}K_{1}+\Lambda^{2}M_{1}\big{)}(U^{j-1}+U^{j+1})+\big{(}K_{0}-\Lambda^{2}M_{0}\big{)}U^{j}=0,$$

$$K_0$$

$$K_1$$

$$K_2$$

$$M_0$$

$$M_1$$

$$M_2$$

$$\displaystyle 2\big{(}K_{2}-\Lambda^{2}M_{2}\big{)}\cos(2\mu^{(2)}_{NQ_{2}}h)-2\big{(}K_{1}+\Lambda^{2}M_{1}\big{)}\cos(\mu^{(2)}_{NQ_{2}}h)+\big{(}K_{0}-\Lambda^{2}M_{0}\big{)}=0,$$

$${\mu }_{N{Q}_2}^{(2)}h$$

$$T_3$$

$$T_5$$

$$\displaystyle\qquad\qquad\qquad\qquad\quad+(\varpi_{2}+\varpi_{1})\big{(}(180n_{1}^{4}-360n_{1}^{3}+120n_{1}^{2}+60n_{1}-17)\varpi_{1}$$

$$\displaystyle\qquad\qquad\qquad\qquad\quad+(180n_{2}^{4}-360n_{2}^{3}+120n_{2}^{2}+60n_{2}-17)\varpi_{2}\big{)}\Big{)}.$$

$$T_3$$

$$T_5$$

$$3n_1(1-n_1{)}^2{\varpi }_1+3n_2(1-n_2{)}^2{\varpi }_2$$

$$3n_1^2(1-n_1){\varpi }_1+3n_2^2(1-n_2){\varpi }_2$$

$$n_1^3{\varpi }_1+n_2^3{\varpi }_2+(1-n_1{)}^3{\varpi }_1+(1-n_2{)}^3{\varpi }_2$$

$$n_1$$

$$n_2$$

$${\varpi }_1$$

$${\varpi }_2$$

$$n_1$$

$$n_2$$

$${\varpi }_1$$

$${\varpi }_2$$

$$n_1$$

Full text

Keywords:

isogeometric analysis; quadrature rule; dispersion analysis; spectrum analysis

1 Introduction

Quadrature rules play an important role in the implementation of various numerical methods for solving partial differential equations. Fewer quadrature points result in a lower computational cost, however, the reduction of the quadrature points should not reduce the quality of the approximation. The design of efficient quadrature rules for isogeometric analysis (see Hughes et al. [18, 9, 12, 11]) is of interest as the continuity properties of the spline basis functions may require fewer quadrature points. The quadrature rules should preserve the optimal convergence of the numerical approximation to the exact solution. Traditionally, Gauss rules for discontinuous polynomial spaces are used, however, these choices are far from being optimal in general [4].

The construction of efficient quadrature rules for isogeometric analysis was initially considered by Hughes et al. in [21] in 2010. Taking advantage of the smoothness of the basis functions across element boundaries, a half-point rule that is independent of the polynomial order of the basis functions was developed. The new rule has advantages when compared to the traditional ones. The rule is optimal as it exactly integrates the spline basis functions with the minimum of number of quadrature sampling points. The rule is designed for uniform univariate splines and is Gaussian, that is, optimal in the sense of the minimum number of quadrature points. However, the rule is exact only for infinite domains or for the spline spaces that have a special structure at the boundaries of finite domains. To make the rule exact for a general spline space over finite domains, additional quadrature points are introduced at the boundary elements, resulting in nearly-optimal quadrature rules [5]. These non-Gaussian rules come as solutions of non-linear, possibly ill-conditioned, systems and possess both positive and negative weights.

Other works in this direction are reported in [27, 17, 5, 6, 7]. Optimal and reduced quadrature rules for tensor product and hierarchically refined splines for isogeometric analysis were developed in [17]. Gaussian rules for spline spaces of various degrees and continuities were derived in [6, 7]. Using the homotopy continuation argument [9], Gaussian rules can be derived by continuously modifying the spline space (knot vector) and by tracing numerically the rule, which is given by solving a certain algebraic system. These rules guarantee exactness of the integration up to machine precision, and the property of being Gaussian also directly implies that all weights are positive [24].

The study of dispersion error minimization for isogeometric analysis was initially studied numerically in Puzyrev et al. [26] and analytically in Calo et al. [10]. For general dispersion analysis of isogeometric discretizations, we refer the readers to [20, 19] and the references therein. Particularly, in Hughes et al. [20], a duality principle between the dispersion analysis and the spectral analysis was established and the analysis unified.

The study of dispersion analysis of the finite element method has a rich literature; see for example Thomson and Pinsky[28, 29], Ihlenburg and Babuska [22], Ainsworth [1, 2, 3], and others [14, 15, 16]. Thomson and Pinsky studied the dispersive effects of the finite element methods with Legendre, spectral, and Fourier local approximation basis for the Helmholtz equation in [28]. They found that the choice of basis functions had a negligible effect on the dispersion errors. This is due to the low continuity ($C^{0}$ continuity) of the basis functions. Hughes et al. [20] showed that the dispersion error of the isogeometric analysis with high continuity (up to $C^{p-1}$ for $p$-th order basis function) on the basis functions is smaller than that of the lower continuity finite element counterparts.

The $2p$-optimal convergence rate of the dispersion error for the $p$-th order standard finite elements was established in [1]. In 2009, Ainsworth and Wajid [2] extended this analysis to arbitrary spectral element methods. Based on Marfurt’s conjecture [23] that the most promising and efficient method for computing wave propagation is to blend the finite element method with the spectral element method with appropriate weights, Ainsworth and Wajid beautifully established the optimal blending of these two methods in [3]. A superconvergence (order $2p+2$ for $p$-th order polynomial approximation) result was obtained for arbitrary order of polynomial approximation, which includes the fourth order superconvergence result obtained by a modified integration rule for linear finite elements in [13].

To the best of our knowledge, this is the first paper studying the design of optimal quadrature rules which minimize the dispersion errors of the isogeometric analysis for the wave propagation and structural vibration problems. The dispersion error-minimizing quadratures, that combine Gauss-Legendre and Gauss-Lobatto rules proposed in [26, 10], are not efficient as the two traditional quadrature rules are used for each integration evaluation. Herein, we design quadrature rules that minimize the dispersion error and minimize the number of quadrature points. A rule that has minimal number of evaluation points per element (two in the case of a uniform mesh with periodic boundary conditions) is the solution of a non-linear system of algebraic equations which, due to the low polynomial degree, admits a closed form formula. We also design a quadrature rule that minimizes dispersion for the larger $C^{0}$ quadratic space, which leads to a quadrature rule that uses 2.5 points per element as it exactly integrates discontinuous cubic functions on the mesh. This rule is effective for non-uniform meshes and arbitrary boundary conditions.

The rest of this paper is organized as follows. Section 2 describes the isogeometric discretization of an eigenvalue problem. In Section 3, we present the constraints minimizing both dispersion error and the number of quadrature points and set up the equations for the quadrature weights and points. Both two-point and 2.5-point rules are considered here. Section 4 studies numerical examples to demonstrate the performance of both the two-point and 2.5-point rules. Concluding remarks are given in Section 5.

2 Model problem and its discretization

This section follows closely [10]. In order to illustrate the main ideas, we consider stationary waves as described by the Helmholtz equation

[TABLE]

where $\Omega=[0,1]^{d}\subset\mathbb{R}^{d},d=1,2,3,\Delta=\nabla^{2}$ is the Laplacian and $\omega=\omega_{f}/c$ with $\omega_{f}$ being the frequency of a particular sinusoidal wave and $c$ being the speed of sound of the medium.

For an open set $S\subset\mathbb{R}^{d}$ with Lipschitz boundary, we denote by $W^{m,p}(S)$ a Sobolev space equipped with the norm $\|\cdot\|_{m,p,S}$ and the semi-norm $|\cdot|_{m,p,S}$ where $m$ is the weak derivative order and $p$ corresponds to the $p$ in $L^{p}$ space. We use standard notation. If $p=2$, we omit $p$ and utilize $H^{m}(S)$ for Hilbert spaces and $H^{m}_{0}(S)$ for Hilbert spaces with functions vanishing at the boundary for $m>0$. The variational formulation of (2.1) is to find $u\in H^{1}_{0}(\Omega)$ such that

[TABLE]

where

[TABLE]

with $a(w,v)=(\nabla w,\nabla v)$ and $b(w,v)=(w,v)$. Here $(\cdot,\cdot)$ is the $L^{2}$ inner product.

Let $\mathcal{T}_{h}$ be a discretization of the bounded open domain $\Omega$ and we denote each element as $K$ such that $\bar{\Omega}=\cup_{K\in\mathcal{T}_{h}}K$. Let $h=\max_{K\in\mathcal{T}_{h}}\text{diameter}(K)$. The Galerkin-type numerical methods seek $u_{h}\in V_{h}$ such that

[TABLE]

Different trial spaces $V_{h}$ lead to different numerical methods. We focus on isogeometric analysis in this work. We denote $\phi_{a}=\phi_{a}(x)$ the B-spline basis functions we use in isogeometric analysis. Then $V_{h}=\text{span}\{\phi_{a}\}$.

In practice, the integrals involved in $a(u_{h},v_{h})$ and $b(u_{h},v_{h})$ are evaluated numerically, that is, approximated by quadrature rules. On a reference element $\hat{K}$, a quadrature rule is of the form

[TABLE]

where $\hat{\varpi}_{l}$ are the weights, $\hat{n_{l}}$ are the nodes, and $N_{q}$ is the number of nodes. For each element $K$, we assume that there is an invertible affine map $\sigma$ such that $K=\sigma(\hat{K})$, which leads to the correspondence between the functions on $K$ and $\hat{K}$. Assuming $J_{K}$ is the corresponding Jacobian of the mapping, (2.5) induces a quadrature rule over the element $K$ given by

[TABLE]

where $\varpi_{l,K}=\text{det}(J_{K})\hat{\varpi}_{l}$ and $n_{l,K}=\sigma(\hat{n}_{l})$.

Applying quadrature rules to (2.4), we have the approximate form

[TABLE]

where

[TABLE]

and

[TABLE]

where $\{\varpi_{l,K}^{(1)},n_{l,K}^{(1)}\}$ and $\{\varpi_{l,K}^{(2)},n_{l,K}^{(2)}\}$ specify two (possibly different) quadrature rules.

With quadrature rules, if we substitute the basis functions of $V_{h}$ into (2.7), this leads to the linear algebra problem

[TABLE]

where $\mathbf{K}$ and $\mathbf{M}$ are the global stiffness and mass matrices with entries $\mathbf{K}_{ab}=\tilde{a}_{h}(\phi_{a},\phi_{b})$, $\mathbf{M}_{ab}=\tilde{b}_{h}(\phi_{a},\phi_{b}),$ and $\mathbf{U}$ is the unknown vector. This system (2.10) admits a nontrivial solution in the view of [1].

3 Quadratures that minimize the dispersion error

For simplicity, we consider a one-dimensional problem where $\Omega=\mathbb{R}$. We describe the general framework for all orders of isogeometric analysis and then focus on the $C^{1}$ quadratic case as follows.

We denote by $C^{k}_{p}$ the space of piecewise polynomials of order $p$ and continuity $k$. Let $\phi_{a}=\phi_{a}(x)$ be a B-spline basis function of isogeometric analysis of order $p$ with maximum continuity $C^{p-1}$, thus $\phi_{a}\in C^{p-1}_{p}$. The functions for the integration corresponding to the stiffness $\tilde{a}_{h}(\phi_{a},\phi_{b})$ are in the space $C^{p-2}_{2p-2}$ while the functions for the integration corresponding to the mass $\tilde{b}_{h}(\phi_{a},\phi_{b})$ are in the space $C^{p-1}_{2p}$. Thus, to integrate both stiffness and mass matrices exactly, one needs quadrature rules which integrate all the functions in the space $C^{p-2}_{2p}$ exactly.

An $m$ point Gauss-Legendre quadrature rule, denoted by $G_{m}$, integrates the space $C^{-1}_{2m-1}$ exactly, while an $m$ point Gauss-Lobatto quadrature rule, denoted by $L_{m}$, integrates the space $C^{-1}_{2m-3}$ exactly. For $p$-th order isogeometric elements, to integrate both the stiffness and mass matrices exactly, the rule $G_{p+1}$ is enough since $C^{p-2}_{2p}$ is a subset of $C^{-1}_{2p+1}$ but $G_{p}$ is not enough. Thus, an optimized rule (minimized number of quadrature points) can be developed by considering the $p-2$ order of continuity; see [6].

In the view of analysis in [10] in a one-dimensional setting, the stiffness matrix is integrated exactly if the space of $C_{2p-1}^{p-2}$ is fully integrated for $p$-th order isogeometric elements while the mass matrix can be under-integrated to minimize the dispersion errors. Our aim is to develop a quadrature rule that minimizes the dispersion errors.

Below, we focus on isogeometric analysis with $C^{1}$ quadratic B-spline basis functions. Let $\phi_{a}=\phi_{a}(x)$ be a $C^{1}$ quadratic B-spline basis function. We seek an approximation of the form

[TABLE]

satisfying

[TABLE]

The quadrature rules we use to integrate (2.8) and (2.9) are the same. In the one-dimension case, we seek a quadrature rule which integrates the stiffness matrix exactly. We denote our new quadrature rule $NQ_{2}$ for $C^{1}$ quadratic isogeometric elements with nodes $n_{1},n_{2}$ in the reference interval $[0,1]$ and weights $\varpi_{1},\varpi_{2}$.

We apply this quadrature rule $NQ_{2}$ to (3.2) to obtain the following equation for the value $U^{j}$ of the approximation at node $x_{j}=jh,j\in\mathbb{Z}$

[TABLE]

where $\Lambda=\omega h$ and

[TABLE]

We assume that the equation admits nontrivial Bloch wave [25] solutions in the form $U^{j}=e^{ij\mu^{(2)}_{Q}h}$, where $i^{2}=-1$ and the subindex $Q$ denotes the corresponding numerical quadrature, then (3.3) simplifies to

[TABLE]

which is known as the discrete dispersion relation for the discrete method with a particular quadrature rule. Solving (3.5) for $\mu^{(2)}_{NQ_{2}}h$ and writing the expression as a series in $\Lambda$, we obtain

[TABLE]

where

[TABLE]

We seek a quadrature rule that reduces the dispersion error as much as possible. Since we consider a two-point rule, we have four degrees of freedom (two nodes and weights). To integrate exactly the $C_{3}^{0}$ space, the rule has to integrate the basis which, taking into account a repetitive pattern on uniform elements, involves three basis functions. This leads to the following algebraic system

[TABLE]

The system is built over the unit interval. The first two equations correspond to the dispersion error-minimization, while the last three represent the exactness of the rule on the $C_{3}^{0}$ spline space. There are five equations and four unknowns. However, symbolic calculations show that the equation $T_{3}=0$ is a redundant equation and this system has the following equivalent solutions:

Solution 1

[TABLE]

Solution 2

[TABLE]

Solution 3

[TABLE]

Solution 4

[TABLE]

These two-point quadrature rules lead to the discrete dispersion relation

[TABLE]

which gives the optimal dispersion error

[TABLE]

Remark 3.1**.**

The last three equations in (3.8) define a one-parameter family of quadrature rules that integrates exactly the $C_{3}^{0}$ space. Among these rules, we seek those that further reduce the dispersion error. Such rules result in an approximation error of order seven, see (3.14).

Remark 3.2**.**

Taking the difference between $\mu^{(2)}_{NQ_{2}}h$ and $\Lambda$ gives an error representation of the dispersion error, which is of order seven. Symbolic calculations show that a three-point (or more points) rule does not increase the convergence order in the dispersion error. The extra degrees of freedom obtained by utilizing more quadrature points can be used to integrate a larger space than the polynomial space $C^{0}_{3}$. We present this alternative in the following subsection. The error with an order of seven as in (3.14) is the minimized dispersion error for $C^{1}$ quadratic isogeometric analysis.

Remark 3.3**.**

This optimized dispersion error expansion coincides with the expansion derived from the blending schemes proposed in [10]. Further calculation shows that both the blending schemes and these new two-point quadrature rules give the same discrete dispersion relation.

3.1 A 2.5-point rule

The two-point quadrature rules (3.9) to (3.12) integrate the space $C^{0}_{3}$ exactly. Alternatively, one can derive a 2.5-point rule which integrates the space $C_{3}^{-1}$ exactly as well as minimizes the dispersion errors. We denote $G_{2.5}$ a three-point quadrature rule with one of the points fixed at one of the element boundaries, that might be shared with another element. That is, this class of quadratures generalizes the Gauss-Radau family of quadratures. We denote this quadrature rule $G_{2.5}$ for $C^{1}$ quadratic isogeometric elements with nodes $n_{1},n_{2},n_{3}$ in the reference interval $[0,1]$ and weights $\varpi_{1},\varpi_{2},\varpi_{3}$. By setting $n_{3}=1$, we say this is a 2.5-point rule per element as the mapping is continuous across the element interface, thus the evaluation at the interface does not need to be recomputed from the neighboring element. By the same process we described to solve (3.8), we obtain the following quadrature rule.

[TABLE]

Alternatively, one can fix a quadrature point at $n_{3}=0$ and derive an alternative, but equivalent rule.

Remark 3.4**.**

For multidimensional case, we assume that a tensor product grid is placed on the domain $\Omega$. Then, we conclude that the above derivations are independent of the spatial dimension and the same rule remains valid for each dimension; more details are referred to [10]. Moreover, a duality between dispersion analysis and spectrum analysis in error expansion form was established in [10].

3.2 Rule near the boundary elements

Naturally, since the 2.5-point rule integrates the $C^{-1}_{3}$ space exactly, it can be applied to both the interior and boundary elements over the domain. This is different for the two-point rule as it takes advantage of the $C^{0}$ continuity of the integrand across element interfaces. This continuity assumption is lost at the boundary elements. Figure 1 (b) shows that the derivatives of the first (last) two basis function do not vanish at $x=0$ ($x=1$) associated with the left (right) boundary element, as is the case of all other element interfaces in the domain.

Special treatment near these boundary elements is required. Simply, we propose to apply the 2.5-point rule at the boundary elements. This can be done in the implementation efficiently.

4 Numerical experiments

In this section, we present the numerical tests of the problem (2.1) in one and two dimensions (denoted with 1D, 2D) with uniform meshes. The comparisons of the numerical eigenvalues and eigenfunctions obtained from isogeometric analysis with those from finite elements are studied in [12, 20, 19, 26] and significant advantages of isogeometric elements over finite elements are shown. In this section, we show the numerical results of isogeometric elements obtained by the proposed new rules as well as the comparisons with those obtained by the optimal-blending rules.

The exact eigenvalues and eigenfunctions of the one dimensional problem (2.1) are

[TABLE]

respectively, while those of the two dimension problems are

[TABLE]

respectively. We sort the approximate eigenvalues both in one and two dimension in the ascending order. In the following figures, we present the eigenvalue (EV) errors as well as the eigenfunction (EF) errors in both $L^{2}$-norm and energy norm.

There are mainly three different optimally-blended rules proposed in [10, 26]: three-point Gauss-Legendre rule with three-point Gauss-Lobatto rule; two-point Gauss-Legendre rule with three-point Gauss-Lobatto rule, and three-point Gauss-Legendre rule with two-point Gauss-Legendre rule. For the comparison with our new rules, we choose the last one as it requires fewer evaluation points, and we denote this one as the blending rule for the tests.

In the one dimensional case, we consider the isogeometric elements with fully-integrated inner products and quadrature-rule approximated inner products. These quadrature rules include the two-point rule, 2.5-point rule, and the blending rule.

Figure 2 shows the relative eigenvalue errors $\dfrac{\lambda_{j}^{h}-\lambda_{j}}{\lambda_{j}}$, the $L^{2}$ eigenfunction errors $\|u_{j}^{h}-u_{j}\|_{0}$, and the scaled energy-norm errors $\dfrac{\|u_{j}^{h}-u_{j}\|_{E}}{\sqrt{\lambda_{j}}}$ in both linear and logarithmic scales when full integration and the two-point, 2.5-point, and the blending rules are applied for the inner product for $C^{1}$ quadratic isogeometric elements. We scale the energy-norm error in the view of the generalized Pythagorean eigenvalue error theorem [26, 10]. In the implementation, for this example, the full integration is realized by applying the three-point Gauss rule. In Figure 2, we observe that two extra orders of convergence in the eigenvalue errors when the two-point, 2.5-point, and the blending rules are applied. Also, we observe that the two-point rule, 2.5-point rule, and the optimal-blending rule lead to the same results. This verifies numerically that these newly-developed rules are equivalent to the optimally-blended rules.

A similar behavior is observed on non-uniform meshes as shown in Figure 3 and in 2D as shown in Figure 4. Figure 3 shows an equivalence between the 2.5-point and the blending rules for $C^{1}$ quadratic isogeometric analysis on non-uniform meshes, precisely the stretching meshes with stretching factors 1.05 and 1.07. A simple 2D test example on a uniform $64\times 64$ mesh is shown in Figure 4. We observe that both two-point and 2.5-point rule lead to the same results as those from optimally-blended rules. For more numerical results in 2D, we refer to the paper [26].

5 Concluding remarks

We introduce new quadrature rules that minimize the dispersion error and possess the optimal septic convergence order. The optimal rule requires two nodes per element and arises from an algebraic system that admits a closed-form solution. This rule is valid for uniform grids with periodic boundary conditions. We also introduce a 2.5-rule which exactly integrates discontinuous, cubic functions on general grids. The optimal rule combined with the generalized 2.5-point Gauss-Radau type of rule on the boundary elements remains dispersion-minimizing on finite domains with arbitrary boundary conditions. Moreover, compared with the optimal blending schemes proposed in [10], our approach further reduces the number of quadrature points, which brings a significant computational speed-up to the application problems such as the wave propagation or structural vibrations, particularly in three dimensions. Future work in this direction includes further studies on non-uniform meshes as well as extension to higher order isogeometric analysis.

6 Acknowledgments

This publication was made possible by the CSIRO Professorial Chair in Computational Geoscience of Curtin University, additional support was provided by the National Priorities Research Program grant 7-1482-1-278 from the Qatar National Research Fund (a member of The Qatar Foundation), and by the European Union’s Horizon 2020 Research and Innovation Program of the Marie Skłodowska-Curie grant agreement No. 644202. The Spring 2016 Trimester on ”Numerical methods for PDEs”, organized with the collaboration of the Centre Emile Borel at the Institut Henri Poincare in Paris supported VMC’s visit to IHP in October, 2016.