PaperlayerView exact revision on arXiv

Source-linked AI summary

Fast Isogeometric Boundary Element Method based on Independent Field Approximation

Benjamin Marussig, Jürgen Zechner, Gernot Beer, Thomas-Peter Fries

arXiv:1406.0306v3math.NA

TL;DR

The paper addresses the computational cost and discretization rigidity of isogeometric BEM for elasticity. It introduces independent geometry, displacement, and traction approximations with H-matrix system compression. Numerical tests show optimal convergence and accuracy together with reduced storage and computational effort relative to conventional isogeometric BEM formulations.

  • Problem

    Isogeometric BEM for elasticity faces fully populated system matrices and motivates methods that reduce computational effort while retaining accurate CAD-based analysis.

  • Method

    The method uses independent NURBS discretizations for geometry, displacement, and traction together with hierarchical matrices for system-matrix approximation.

  • Results

    Numerical tests in two and three dimensions confirm optimal convergence and accuracy, with potential reductions in storage requirements and computational effort compared with conventional isogeometric BEM.

  • Takeaways & Limitations

    Independent field discretization avoids redundant refinement information while preserving exact geometry and enabling flexible treatment of boundary fields.

  • Takeaways & Limitations

    The optimal position of the discontinuous collocation points remains to be investigated beyond the examples shown.

Abstract

from arXiv · show

An isogeometric boundary element method for problems in elasticity is presented, which is based on an independent approximation for the geometry, traction and displacement field. This enables a flexible choice of refinement strategies, permits an efficient evaluation of geometry related information, a mixed collocation scheme which deals with discontinuous tractions along non-smooth boundaries and a significant reduction of the right hand side of the system of equations for common boundary conditions. All these benefits are achieved without any loss of accuracy compared to conventional isogeometric formulations. The system matrices are approximated by means of hierarchical matrices to reduce the computational complexity for large scale analysis. For the required geometrical bisection of the domain, a strategy for the evaluation of bounding boxes containing the supports of NURBS basis functions is presented. The versatility and accuracy of the proposed methodology is demonstrated by convergence studies showing optimal rates and real world examples in two and three dimensions.

1. Introduction

Isogeometric BEM combines CAD-based boundary representations with boundary-integral analysis, avoiding difficult volume meshing while retaining exact geometry. The paper introduces an independent discretization for geometry, displacement, and traction, combined with H-matrices for efficient elasticity analysis.

  • Motivation: CAD-based BEM avoids the difficult generation of volume discretizations required by isogeometric finite element analysis.BEM also uses boundary representations, allowing CAD models to be used directly for simulation.
  • Computational challenge: Because BEM system matrices are fully populated, fast methods such as H-matrices are needed to reduce computational effort.The introduction identifies H-matrices alongside FMM, wavelets, and fast Fourier methods as complexity-reduction strategies.
  • Related work: Prior fast BEM work with Bézier, B-spline, and NURBS functions has addressed only a limited set of applications and formulations.The cited examples include multilevel FMM, wavelets for trimmed NURBS, and customized fast techniques for Laplace problems.
  • Proposed formulation: The proposed formulation independently discretizes geometry, displacement, and traction, enabling discontinuous tractions on non-smooth boundaries while preserving continuous displacements.Independent refinement also supports finer resolution for selected fields and reduces geometry-related evaluation and right-hand-side storage requirements.
  • Proposed formulation: H-matrices are applied to the system matrices, with NURBS-specific construction modifications described alongside the formulation.The method retains exact CAD geometry despite violating the conventional isoparametric concept.

2. Boundary Integral Equation

The boundary integral formulation converts the Lamé-Navier elasticity problem into relations among boundary displacements and tractions. Somigliana’s identity uses displacement and traction fundamental solutions, with singular integrals handled through boundary operators and principal-value treatment.

  • Governing equation: The Lamé-Navier equation governs the displacement field of the elastic body, with λ and µ as material constants.The closed boundary is Γ, and its outward normal defines the boundary operators.
  • Boundary operators: The trace operator maps domain displacements to boundary displacements, while the conormal derivative maps displacements to surface tractions.These boundary quantities provide the data used in the boundary-value problem and its integral formulation.
  • Boundary conditions: The elasticity problem prescribes displacements on Dirichlet boundaries and tractions on Neumann boundaries, partitioning Γ into ΓD and ΓN.The prescribed data are denoted gD for displacements and gN for surface tractions.
  • Boundary integral equation: Somigliana’s identity relates boundary displacement and traction data through fundamental solutions for displacements and tractions.The trace operator shifts the evaluation points to the boundary, while U and T represent the corresponding fundamental solutions.
  • Singular operators: The single-layer operator is weakly singular, whereas the double-layer operator is strongly singular and requires a Cauchy principal-value interpretation.The limiting process contributes an integral-free term, with c = 1/2 on smooth surfaces.

3. Discretisation

The discretisation independently approximates geometry, displacement, and traction fields, enabling field-specific refinement while preserving the exact NURBS geometry. Anchors and mixed collocation accommodate continuous displacements and discontinuous tractions, while reducing right-hand-side cost for common boundary conditions.

  • 3.2. Subparametric Patches: Independent approximations for geometry, displacement, and traction fields define subparametric patches with more unknown-field parameters than geometry parameters.This generalises the conventional isoparametric approach and permits refinement tailored to each field.
  • 3.2.1. Geometry Representation: The geometry remains at its initial CAD-based representation, so the computational boundary equals the exact design-model boundary without geometry approximation error.Unnecessary geometry refinement is therefore avoided throughout the analysis.
  • 3.2.1. Geometry Representation: Geometry evaluation has quadratic complexity in the order of the basis functions, making order elevation of the geometry basis computationally undesirable.The increased effort particularly affects tangent-vector, Gram-determinant, and outward-normal calculations.
  • 3.2.3. Individual Refinement: Displacement basis functions remain C0-continuous, whereas traction basis functions are C−1-continuous to represent traction jumps at corners and edges.The traction basis is constructed by knot insertion from the initial basis, while field refinements are applied only when needed.
  • 3.3. Collocation: Anchors based on Greville abscissae, with an offset for discontinuous functions, determine collocation points and enable mixed collocation on different boundary parts.Displacement anchors are used on ΓN, while discontinuous traction anchors are used on ΓD.
  • 3.3. Collocation: The right-hand-side computational effort can be O(nm) with m ≪ n, because m depends on nonhomogeneous prescribed data and homogeneous-condition entries are omitted.Without this reduction, system-matrix storage is O(n^2).

4. Hierarchical Matrices

The method reduces the quadratic storage and matrix-vector costs of isogeometric BEM by combining H-matrix approximation with ACA and geometry-aware clustering of NURBS supports.

  • 4. Hierarchical Matrices: O(n^2) storage and matrix-vector effort are reduced by applying hierarchical matrices to the fully populated system matrices.H-matrices provide asymptotic O(nlogn) storage for fully populated matrices.
  • 4.1. Matrix Approximation: ACA constructs low-rank matrix blocks from selected rows and columns, with adaptive rank controlled by the approximation quality ϵH.Pivot points guide the selected rows and columns, while the rank k is determined adaptively.
  • 4.1. Matrix Approximation: For elasticity, decomposing Cauchy data by spatial direction preserves suitable block structure for ACA when the fundamental solution is tensor-valued.Without decomposition, singular values may decay less favorably and increase rank and calculation time.
  • 4.2. Geometrical Bisection: Geometrical clustering labels collocation points and NURBS anchors, then recursively bisects bounding boxes along their largest extension until the minimum leaf size is reached.The resulting cluster trees support the partition into near-field full blocks and far-field blocks eligible for approximation.
  • 4.2. Geometrical Bisection: Bounding boxes for NURBS basis supports are obtained from Bézier segments and the convex hull of their control grids, enabling geometry-aware bisection.The construction is illustrated for a cubic NURBS curve with separate traction and displacement supports.
  • 4.2. Geometrical Bisection: The H-matrix partition classifies sufficiently separated row-column clusters as far field and recursively refines the remaining near-field blocks.The admissibility factor η controls how many blocks become low-rank candidates.

5. Numerical Results

Numerical tests assess integration efficiency, operator accuracy, convergence, matrix compression, and storage savings across two- and three-dimensional problems. The proposed subparametric formulation preserves accuracy while reducing geometry-evaluation and system-storage costs.

  • Single Patch Integration: For geometry order pg = 5, superparametric integration takes almost twice as long as integration at the lowest geometry order.The number of quadrature points stays constant, so the added cost is attributed to more expensive geometry evaluations; Kτ is affected slightly more because it requires the outward normal.
  • Single Patch Integration: An order elevation by 2 yields approximately 25% and 35% speedups for Vτ and Kτ, respectively, over the isogeometric formulation.The speedup comes from refining Cauchy-data bases while retaining the original geometry representation; Kτ benefits additionally from avoiding refined normal evaluations.
  • Operator accuracy: Optimal convergence is demonstrated for the discrete single- and double-layer boundary integral operators.The operators are evaluated with an integration error of ϵQ = 10^-11, and the reported error decreases optimally.
  • Tunnel Excavation in Two Dimensions: Smooth NURBS approximations are more accurate than Bézier segments when relative L2-error is plotted against total degrees of freedom n.Bézier segments can appear better against mesh parameter h because each element contains more degrees of freedom, but NURBS are superior with respect to n.
  • Three-dimensional tests: All torus tests achieve optimal convergence rates, while H-matrix compression does not affect convergence before the prescribed approximation error ϵH is reached.For the exterior Neumann problem, the relative error remains below the compression threshold across the tested orders.
  • Storage and compression: Subparametric patches preserve result quality and make the finest system matrix approximately 100 times smaller than the isoparametric matrix.The storage reduction increases with the number of unknowns because redundant refinement of the geometry representation is avoided.

6. Conclusion

The paper introduces a subparametric collocation BEM that independently discretizes geometry and Cauchy data while using H-matrices for compression. Numerical tests report optimal accuracy and reduced storage and computational effort, with collocation-point placement remaining an open issue.

  • Method: The method combines NURBS collocation, independent geometry and Cauchy-data descriptions, and H-matrices in a subparametric BEM.The formulation avoids strict isoparametric refinement by assigning separate representations to geometry, displacement, and traction fields.
  • Efficiency: Independent field refinement removes redundant geometry control parameters and reduces right-hand-side effort while H-matrices reduce left-hand-side complexity to O(nlogn).The separation of field discretizations is especially useful when known boundary data are defined over whole patches.
  • Results: Numerical tests in two and three dimensions confirm optimal convergence and accuracy together with reduced storage and computational effort.The reported tests compare the method with conventional isogeometric BEM formulations.
  • Approximation choices: B-spline and NURBS Cauchy-data approximations produce almost identical results, whereas non-smooth Bézier segments reduce accuracy with respect to degrees of freedom.The Bézier comparison concerns the approximation itself and is distinct from Bézier extraction.
  • Limitation: The optimal placement of discontinuous collocation points remains to be investigated beyond the examples presented.The selected points are reported as suitable for the demonstrated examples, but their optimal position is unresolved.

Appendix A. Basis Function Derivatives

Appendix A presents derivative expressions for B-spline and NURBS basis functions, including bivariate derivatives in a parametric direction.

  • Basis Function Derivatives: The appendix introduces the computation of first derivatives for B-spline basis functions.
  • Basis Function Derivatives: It gives the corresponding first-derivative expression for NURBS basis functions using weighted basis-function terms.
  • Basis Function Derivatives: For bivariate NURBS bases, the derivative in one parametric direction is expressed using the correlated parametric coordinates and orders.

Appendix B. Number of Elementary Operations

Appendix B frames operation counting for basis-function evaluation in terms of elementary floating-point operations.

  • Number of Elementary Operations: Elementary operations are treated as floating-point operations such as multiplication or division.
  • Number of Elementary Operations: The appendix derives the number of elementary operations required to evaluate B-spline functions.
  • Number of Elementary Operations: The same operation-counting analysis covers evaluation of NURBS functions.

Appendix B.1. Basis Functions and Derivatives

The appendix details operation counts for evaluating univariate and bivariate B-spline/NURBS basis functions and their derivatives, including reuse of intermediate results. It also records a correction to the original NURBS derivative definition.

  • Univariate basis functions: Reusing intermediate results reduces the derivative evaluation cost to 2 operations for each N′i,p.The non-zero basis functions of order p−1 are computed first, then the recurrence is applied to all non-zero terms.
  • Correction: The original article contains an error in the definition of NURBS derivatives, corrected in equations (A.2) and (A.3).The correction concerns the derivative definitions used in the appendix.
  • Derivatives: First derivatives require additional operations for basis derivatives, weighted products, tensor products, and rational normalization in both parametric directions.The operation accounting separately includes derivative construction in the r and s directions.
  • Bivariate basis functions: Bivariate basis evaluation combines univariate basis computations, tensor products, weighting, and rational normalization.The listed operations include tensor products of basis functions and their weighted forms, followed by numerator–denominator division for each NURBS function.

Appendix B.2. Curve and Surface Evaluations

Curve and surface evaluation combines basis-function operations with mapping operations to obtain points and tangents. The resulting operation counts are summarized for geometry-related information.

  • Geometry evaluation: Mapping operations must be added to basis-function and derivative costs to evaluate points of B-spline and NURBS objects.For constant order p in all parametric directions, the appendix gives additional mapping operations depending on the parametric dimension d.
  • Geometry evaluation: Surface tangents require applying the mapping twice, once for each intrinsic direction.The operations associated with the mapping in equation (B.11) are therefore doubled for surface tangents.
  • Operation summary: Table B.1 summarizes elementary-operation counts for evaluating geometry-related information, including points and tangents.For surfaces, the order is equal in each direction and tangents are calculated in both intrinsic directions.

Appendix C. Numerical Integration

Numerical integration is performed element-wise over NURBS patches after subdividing knot-span elements into integration regions. Accurate integration is important because the integrands combine rational geometry, basis functions, and the fundamental solution.

  • Element-wise integration: System-matrix coefficients are evaluated numerically by integrating each NURBS patch element-wise over subdivided integration regions.Elements are defined by non-zero knot spans in the Cauchy-data description, and each matrix entry is expressed as a sum of regional integrals.
  • Integration regions: The integration procedure therefore evaluates contributions over subdivided regions rather than treating the entire NURBS patch as one integration domain.The subdivision produces regions denoted ˆτ⊡ for assembling the system-matrix entries.
  • Integration accuracy: Accurate integration is crucial because the integrand contains the fundamental solution, NURBS basis functions, and coordinate transformations, all represented as rational functions.The exact geometry and higher-order basis functions can reduce the degrees of freedom needed relative to conventional BEM formulations.

Appendix C.1. Regular and Nearly Singular Integration

The method uses regular subdivision, singular subdivision, reference-element quadrature, and adaptive order increases to integrate regular and nearly singular regions. Nearly singularity is identified by comparing region size with distance to the target.

  • Regular and nearly singular subdivision: Regions are first subdivided so their Cartesian edge lengths are approximately equal, while nearly singular integrals are identified using diameter-to-distance comparisons.The subdivision strategy is designed for extensive three-dimensional NURBS patches and uses a heuristic criterion for near singularity.
  • Regular and nearly singular subdivision: Figure C.27 contrasts regular subdivision into integration regions with singular subdivision that isolates a regularisation area and separately treated regions.The singular case labels the regularisation area ˆτs and regions ˆτ1, ˆτ2, and ˆτ3.
  • Adaptive quadrature: Adaptive integration compares quadrature results at orders G0 and G1, increasing the order further when their variation exceeds the allowed error ϵQ.If the error criterion is not satisfied, the procedure restarts with a higher quadrature order.

Appendix C.2. Singular Integration

Singular collocation-point integrals are regularized locally with specialized transformations, while the remaining integration region uses the regular geometry-driven subdivision strategy.

  • Singular integrals are isolated within a rectangular neighborhood around the collocation point and treated with regularisation techniques.The surrounding boundary segment is subdivided before the singular region is processed.
  • Duffy transformation treats the weakly singular integral V, while the Guiggiani–Gigante procedure handles strongly singular entries in K.
  • The leftover integration area is subdivided using equilibrated aspect ratios and a geometry-driven hierarchical criterion.This follows the strategy used for regular integration.
Loading 1406.0306v3…