PaperlayerView exact revision on arXiv

Source-linked AI summary

A Finite Element Method for the Surface Biharmonic Problem without Stabilization or Piola Transforms

Jiachen Gu, Xia Ji, Shuo Zhang

arXiv:2610.23492v1math.NAcs.NA

TL;DR

The surface biharmonic problem calls for a stable discretization that avoids Piola transforms and stabilization while retaining the Laplacian form. The paper constructs a scalar-assembly REQ method and proves weak discrete Miranda–Talenti estimates. It establishes well-posedness and the stated convergence rates on sufficiently refined interpolating surface meshes.

  • Problem

    On folded meshes, face gradients belong to different tangent planes, and Piola coupling alone does not provide the Miranda–Talenti control needed to stabilize the pure Laplacian form.

    Source: §1 · pp. 2–3

  • Method

    The paper constructs a reduced enriched quartic element with scalar assembly and an unstabilized elementwise Laplacian form, then bounds vertex residuals by separating planar and angle-defect contributions.

    Source: §1 · pp. 3–4

  • Results

    The method is well posed and, under the stated regularity assumptions, achieves first-order energy-norm and second-order broken H1 and L2 convergence.

    Source: §1 · p. 4

  • Takeaways & Limitations

    The estimates yield uniform full broken H2 control on sufficiently refined interpolating surface meshes in the stated mesh families.

    Source: §1 · p. 4

    The analysis assumes uniformly shape-regular, quasi-uniform closed mesh families with uniformly bounded vertex valences.

    Source: §2.5 · p. 12

Abstract

from arXiv · show

For the biharmonic problem on closed smooth surfaces, we construct a nonconforming finite element method that requires neither stabilization nor the Piola transform. The reduced enriched quartic (REQ) element relies solely on vertex values and scalar edge moments. The same local construction applies to every triangle of a piecewise planar surface. The bilinear form consists of elementwise Laplacian inner products and permits standard finite element assembly. By separating the planar vertex residual from the angle defect contribution, we establish weak discrete Miranda--Talenti-type estimates on both planar meshes and closed piecewise planar surfaces. These estimates guarantee uniform stability on sufficiently refined interpolating surface meshes. We prove the well-posedness of the method, first-order convergence in the energy norm, and second-order convergence in the broken $H^1$ and $L^2$ norms. Numerical tests on three surfaces show the predicted convergence rates and examine the discrete Miranda--Talenti constants and angle defect scaling.

1. Introduction

The paper develops an unstabilized, scalar-assembly finite element method for the surface biharmonic problem, avoiding Piola transforms. A weak discrete Miranda–Talenti estimate supports its stability and convergence analysis.

1. Introduction

  • The Laplacian form directly discretizes the surface problem without an explicit curvature coefficient, but its stability requires control of the full broken Hessian by the piecewise Laplacian.The paper motivates this form by the Bochner–Grisvard identity, which relates Hessian and Laplacian energies through a Gaussian-curvature term.Source: §1 · p. 2
  • On folded meshes, adjacent face gradients lie in different tangent planes, complicating extensions of planar elements with vertex derivative continuity.Prior surface NZT work used Piola coupling and conormal-derivative jump stabilization; a stabilized Morley method instead allowed scalar assembly but used a different energy.Source: §1 · p. 3
  • The REQ element combines scalar assembly with the unstabilized Laplacian form, using local constructions on piecewise planar meshes without Piola transforms.Its 12-dimensional space contains all quadratic polynomials and uses vertex values, one function moment per edge, and two conormal derivative moments per edge.Source: §1 · pp. 3–4
  • A weak discrete Miranda–Talenti estimate follows by splitting vertex residuals into a planar contribution and an angle-defect correction.The planar part is bounded by patchwise Laplacian energy, while the angle-defect correction is controlled by a lower-order gradient term.Source: §1 · p. 4
  • Under the stated regularity assumptions, the method is well posed and converges first order in the energy norm and second order in the broken H1 and L2 norms.On the considered shape-regular, quasi-uniform interpolating surface mesh families, the discrete Miranda–Talenti estimate gives uniform full broken H2 control for sufficiently small h.Source: §1 · p. 4
  • Numerical tests on a sphere, torus, and implicit surface illustrate the predicted convergence rates and examine Miranda–Talenti constants and angle-defect scaling.The paper also studies perturbations of the mesh families.Source: §1 · p. 4

2. The REQ Element and Discrete Miranda–Talenti-Type Estimates

The section develops the REQ element and analyzes its broken-Laplacian energy through vertex residuals to establish discrete Miranda–Talenti-type estimates on planar meshes and, under stated mesh-family assumptions, on closed piecewise planar surfaces. It also completes the proof of the stated inequalities for these estimates.

2.1. Piecewise planar triangulations and broken differential operators.

  • The signed angle defect is the closure defect of a vertex fan after rigid unfolding; zero defect means intrinsic flatness, not necessarily coplanarity of adjacent triangles.The vertex patches need not be convex.Source: §2.1 · p. 5

2.2. Reduced enriched quartic (REQ) element.

  • The 12-dimensional local REQ space has unisolvent degrees of freedom, contains P2(K), reproduces quadratics by interpolation, and has only harmonic quadratics in its harmonic kernel.The construction reduces an enriched quartic space using edge constraints; the interpolation map is defined for H3(K) functions.Source: §2.2 · p. 6

2.3. Global REQ space.

  • REQ assembly uses scalar vertex and edge data and elementwise planar contributions, without identifying tangent planes or vertex derivatives.Each triangle uses the same local construction, with stiffness contributions computed from its planar geometry.Source: §2.3 · pp. 7–8
  • On a connected closed piecewise planar triangulation, the broken Laplacian form is positive definite on the zero-mean REQ space, without a mesh-size restriction.Its kernel consists of global constants.Source: §2.3 · p. 8

2.4. Vertex residuals and rigid unfolding.

  • Rigid unfolding decomposes edge residuals into vertex contributions, and the vertex-patch residuals admit bounds for interior vertices and planar homogeneous-boundary patches.For surface fans, the final-edge contribution includes a holonomy correction; boundary patches have no closure rotation.Source: §2.4 · pp. 9, 11

2.5. Discrete Miranda–Talenti-type estimates.

  • On planar meshes, the discrete Miranda–Talenti-type estimate holds for functions in Vh,0; on closed surface mesh families, a uniform estimate holds for sufficiently small h under the stated assumptions.The closed-surface result assumes uniform Poincaré control and a uniformly bounded discrete curvature scale.Source: §2.5 · p. 12
  • The proof establishes the first inequality in (2.30), derives the reverse from the pointwise relation |∆v| ≤ 2|D2v|, and obtains (2.31) from (2.35).The reverse inequality follows directly from the pointwise relation; the final estimate follows from (2.35).Source: §2.5 · p. 13

3. Discretization and Error Analysis

This section develops geometric, interpolation, stability, and consistency estimates for the surface discretization and uses them to formulate energy-norm and broken H1 error estimates. It relates approximation and transport bounds to first-order energy convergence and second-order L2 convergence.

3.1. Continuous problem and surface approximation.

  • The continuous Laplacian controls the covariant Hessian with a geometry-dependent constant, while the angle-defect analysis bounds its scaling and establishes convergence under mesh refinement.The angle argument uses curvature-dependent local angle changes and a uniformly bounded number of incident triangles.Source: §3.1 · pp. 14, 16

3.2. Finite element discretization.

  • The discrete problem has a unique solution on each connected closed triangulation with nondegenerate triangles and is uniformly stable for sufficiently small h under the geometric assumptions.The uniform estimate is independent of h.Source: §3.2 · p. 16

3.3. Global surface interpolation and approximation.

  • For w in H3(γ), the broken local interpolant satisfies ∥we − Πhwe∥m,h ≤ Cγ h^(3−m)∥w∥3,γ for m = 0, 1, 2; the global operator enforces compatibility across edges.The global construction uses shared vertex values and compatible function and conormal-derivative moments.Source: §3.3 · p. 17

3.4. Conforming companion operator and consistency estimates.

  • The mean-preserving conforming companion differs from each discrete function by O(h^2) in L2 and scaled H1 norms, and its H1 norm is controlled by the discrete Laplacian.Interpolants also satisfy a second-order bilinear consistency bound in h.Source: §3.4 · pp. 18, 20

3.5. Energy-norm and broken H1 error estimates.

  • The analysis establishes energy-norm and broken H1 error estimates for the discrete solution under the assumptions of Theorem 3.8.The supplied theorem passages identify the estimates but do not show their rates.Source: §3.5 · p. 21
  • First-order energy convergence combines quadratic-reproduction approximation error with Laplacian-transport error, while the full broken H1 estimate yields second-order L2 convergence.The energy estimate also contains the direct Laplacian transport error.Source: §3.5 · p. 22
  • The natural energy norm is ∥∆hvh∥0,h, and broken H2 control means no jump or stabilization seminorm enters it.The broken H2 control is provided by estimate (2.31).Source: §3.5 · p. 22
  • Geometric perturbations beyond direct Laplacian transport are order h2, arising from area and tangential-gradient transport and reflected in the second-order consistency bound.Source: §3.5 · p. 22
  • Higher-order surface approximations can reduce geometric errors, but third-order L2 estimates require additional approximation or dual-consistency arguments.Source: §3.5 · p. 22

4. Numerical Experiments

Experiments on the sphere, torus, and implicit surface test convergence, while additional mesh studies examine discrete Miranda–Talenti constants and angle-defect scaling. Across the examples, measured errors broadly follow the predicted rates, with finite-mesh diagnostics consistent with the stated estimates.

4.1. Example 1: Unit sphere.

  • On refined sphere meshes, E0 and E1 converge at second order and E∆ at first order with both mesh-size measures.The numerical solution approximates the cubic spherical harmonic; Table 1 reports errors and experimental rates.Source: §4.1 · p. 23, Table 1 · p. 24

4.2. Example 2: Torus.

  • On finer torus meshes, E0 and E1 are broadly consistent with second-order convergence and E∆ with first-order convergence.Coarser meshes show preasymptotic behavior; the final E0 fluctuation may relate to local differences between independently generated meshes.Source: §4.2 · pp. 23–24, Table 2 · p. 24

4.3. Example 3: General implicit surface.

  • On the implicit surface, E0 and E⋆_1 show overall second-order convergence, while E∆ shows first-order convergence.Fluctuations in rmax may relate to a few exceptionally large elements in independent non-nested meshes; the overall error decay remains consistent with the analysis.Source: Table 3 · p. 25, §4.3 · p. 25

4.4. Discrete Miranda–Talenti constants.

  • The discrete Miranda–Talenti study uses five meshes per example and three projected vertex-perturbation meshes per original mesh, preserving connectivity and degrees of freedom.All 60 meshes pass the geometric check with minimum angles exceeding 18°.Source: §4.4 · pp. 25–26
  • Across the tested refinement range, CMT,h lies between 1.15326 and 1.28695, and perturbations change original-mesh values by less than 0.21%.Seven finest-mesh computations fail the accuracy checks and are omitted; the reported constants are consistent with the discrete estimate over the tested range.Source: Table 4 · p. 26, §4.4 · p. 26, Figure 7 · p. 27

4.5. Angle defect scaling.

  • Normalized angle defects stay below 0.81 on the sphere, 4.0 on the torus, and 20.1 on the implicit surface; absolute defects decrease on all three.Original implicit-surface normalized values rise from 8.10 to 20.03 without reaching a plateau; the finite-mesh tests do not establish uniform bounds under further refinement.Source: Figure 8 · p. 27, §4.5 · p. 27, Table 5 · p. 28

5. Conclusions and Outlook

The REQ method combines scalar vertex values and edge moments with an unstabilized Laplacian scheme; residual estimates yield stability and the stated convergence rates. Numerical tests examine those rates, Miranda–Talenti constants, and angle-defect scaling.

Conclusions and Outlook

  • The 12-dimensional REQ element uses vertex values and scalar edge moments, allowing standard assembly without transporting vertex gradients.Separating planar vertex residuals from angle-defect corrections gives a weak discrete Miranda–Talenti-type estimate and uniform broken H2 stability on sufficiently refined interpolating meshes.Source: §5 · p. 28
  • The method converges at first order in the energy norm and second order in the broken H1 and L2 norms.The analysis combines stability with approximation and consistency estimates.Source: §5 · p. 28
  • The global discrete harmonic kernel contains only constants on connected closed meshes with non-degenerate triangles.Source: §5 · p. 28
  • Three numerical examples exhibit the predicted convergence rates, while tested perturbations produce little variation in the Miranda–Talenti constants.Absolute angle defects decrease under refinement, but normalized defects on the implicit surface have not reached a plateau.Source: §5 · p. 28
  • Further work targets residual analysis on locally refined meshes and surfaces with boundary or evolving geometry, as well as effective preconditioners.Source: §5 · p. 28

A.1. Elimination of the enriched bubble function coefficients.

This appendix establishes that the enriched bubble coefficients must vanish on every non-degenerate triangle, using polynomial decomposition and an algebraic rank argument.

A.1. Elimination of the enriched bubble function coefficients

  • The local function is decomposed into a quartic polynomial and three enriched bubble terms, with barycentric coordinates used for the calculation.Source: §A.1 · p. 28
  • In barycentric coordinates, the Laplacian has coefficients g22, g23, and g33 satisfying positive-definiteness conditions on a non-degenerate triangle.Source: §A.1 · p. 29
  • Restricting the polynomial relations to two edges and the third edge parameterization yields algebraic conditions for the bubble coefficients.Source: §A.1 · p. 29
  • Eliminating the metric ratios leaves one candidate compatible with positive definiteness, but an additional minor is nonzero there.The resulting matrix has full column rank on every non-degenerate triangle, so all three enriched bubble coefficients vanish.Source: §A.1 · p. 29

A.2. Elimination of cubic and quartic harmonic parts.

The appendix eliminates the cubic and quartic harmonic parts by expressing them in complex coordinates and applying the element’s reduction constraints along edges.

A.2. Elimination of cubic and quartic harmonic parts

  • After placing a vertex at the origin and identifying the element plane with the complex plane, the harmonic cubic and quartic parts have real-part forms.Source: §A.2 · p. 29
  • The degree-at-most-two part can be subtracted under the reduction constraints, and an edge residual is defined from endpoint values and an edge moment.Source: §A.2 · p. 29
  • The reduction constraints yield an edge relation involving the quartic normal derivative, which is combined with differentiation of the harmonic forms along the conormal.Source: §A.2 · pp. 29–30
  • The constraints imply α + 2βτ = 0; applying this relation on two distinct edges from one vertex forces β and α to vanish.Source: §A.2 · p. 30

Appendix B. Technical Details for the Vertex-Residual Estimates

The vertex-residual estimate is developed by aligning neighboring element planes around an interior edge and fixing a common edge orientation.

Appendix B. Technical Details for the Vertex-Residual Estimates

  • For an interior edge shared by two faces, the estimate setup rigidly aligns their planes and uses a common edge orientation.Source: Appendix B · p. 30

B.1. Complete cancellation of interior edge integrals.

Interior-edge integrals cancel by separating polynomial and enriched contributions; cyclicity extends the enriched cancellation across all three edges.

B.1. Complete cancellation of interior edge integrals.

  • The polynomial contribution vanishes because its second tangential derivative lies in P1(e), orthogonal to the edge jump.The polynomial trace belongs to P3(e), so its second tangential derivative lies in P1(e).Source: §B.1 · p. 30
  • The enriched edge jump has profile χ(t), while only one bubble contributes conormal derivatives, with profile ψ(t).The two bubbles containing λk have zero conormal derivatives on the shared edge.Source: §B.1 · p. 30
  • Cyclic application across k = 1, 2, 3 makes the enriched contribution vanish on every edge, proving the cancellation identity.Combining the polynomial and enriched parts proves (B.1).Source: §B.1 · p. 30

B.2. Uniformity of local vertex patch estimates.

The local patch estimates are made uniform by normalizing admissible patches and using compactness, then transferring the bound back under a uniform valence bound.

B.2. Uniformity of local vertex patch estimates.

  • After fixing patch valence, normalized admissible configurations form a compact set with uniformly non-degenerate triangles.Closed planar patches have pm+1 = p1, while unfolded surface fans need not be closed.Source: §B.2 · p. 30
  • Continuity on the compact configuration set gives a uniform bound for the symmetric bilinear form Bflat_z.The constraint matrices and bilinear-form matrix depend continuously on the configuration.Source: §B.2 · pp. 30–31
  • A compactness contradiction establishes a positive lower bound on the form for vectors orthogonal to the harmonic kernel.A hypothetical normalized sequence with vanishing form converges to a vector both in and orthogonal to the kernel, contradicting unit norm.Source: §B.2 · p. 31
  • For arbitrary vectors, subtracting their kernel projection reduces the estimate to the orthogonal component without changing the cyclic determinant sum.The harmonic projection lies in the kernel, and the determinant sum is unchanged by adding a common harmonic quadratic polynomial.Source: §B.2 · p. 31
  • A uniform valence upper bound permits only finitely many patch valences, so undoing normalization proves (2.25) with a uniform constant.Normalization preserves local constraints and moment conditions while leaving the relevant ratio invariant.Source: §B.2 · pp. 30–31
Loading 2610.23492v1…