DeepAI

Inverse Inequality Estimates with Symbolic Computation

In the convergence analysis of numerical methods for solving partial differential equations (such as finite element methods) one arrives at certain generalized eigenvalue problems, whose maximal eigenvalues need to be estimated as accurately as possible. We apply symbolic computation methods to the situation of square elements and are able to improve the previously known upper bound, given in "p- and hp-finite element methods" (Schwab, 1998), by a factor of 8. More precisely, we try to evaluate the corresponding determinant using the holonomic ansatz, which is a powerful tool for dealing with determinants, proposed by Zeilberger in 2007. However, it turns out that this method does not succeed on the problem at hand. As a solution we present a variation of the original holonomic ansatz that is applicable to a larger class of determinants, including the one we are dealing with here. We obtain an explicit closed form for the determinant, whose special form enables us to derive new and tight upper resp. lower bounds on the maximal eigenvalue, as well as its asymptotic behaviour.

• 15 publications
• 2 publications
• 2 publications
01/24/2020

Guaranteed Lower Eigenvalue Bound of Steklov Operator with Conforming Finite Element Methods

For the eigenvalue problem of the Steklov differential operator, by foll...
01/15/2020

Finite Element Approximation of Transmission Eigenvalues for Anisotropic Media

The transmission eigenvalue problem arises from the inverse scattering t...
04/22/2022

Convergence of Lagrange Finite Element Methods for Maxwell Eigenvalue Problem in 3D

We prove convergence of the Maxwell eigenvalue problem using quadratic o...
03/03/2019

Optimal quadratic element on rectangular grids for H^1 problems

In this paper, a piecewise quadratic finite element method on rectangula...
04/09/2020

Finite Element Approximation of the Modified Maxwell's Stekloff Eigenvalues

The modified Maxwell's Stekloff eigenvalue problem arises recently from ...
08/21/2022

Efficient numerical method for reliable upper and lower bounds on homogenized parameters

A numerical procedure providing guaranteed two-sided bounds on the effec...
05/21/2020

A cookbook for finite element methods for nonlocal problems, including quadrature rules and approximate Euclidean balls

The implementation of finite element methods (FEMs) for nonlocal models ...

1. Introduction

Interdisciplinary collaborations between different areas of mathematics can be hard work because of different terminology and the difficulty of recognizing the applicability of the methods from one field in the other field. In Linz there is an almost-20-year tradition of bringing together researchers from numerical mathematics and symbolic computation [20, 21, 15, 1], which at the beginning faced exactly these kinds of problems. Additionally, there could be the risk that the results are only interesting for one community and not rewarded by the other one. Fortunately, this didn’t happen in our case: in the current work, we use and invent tools at the frontier of symbolic computation to solve a problem that arose at the frontier of numerical analysis research. Hence, this work improves the knowledge and tools for both communities.

Inverse inequalities of the form

 (1) ||vn||X(Ω) ≤c1(h,n)||vn||Y(Ω)for all vn∈Vn, (2) ||vn||Z(∂Ω) ≤c2(h,n)||vn||Y(Ω)for all vn∈Vn

play an important role in the analysis and design of numerical methods for partial differential equations [4, 22, 2, 6] and in the construction of efficient solvers for the arising linear systems of those methods [10, 23]. Bounds for the constants of the type (1)–(2) have been studied for example in [22, 25, 24, 9, 7], where the asymptotic behaviour with respect to and is covered but usually these constants are over estimated. In many numerical methods a precise knowledge of these constants is required, which motivates this work where we use and present tools from symbolic computation to derive precise estimates.

Here , is a bounded and open set with sufficiently smooth boundary , describing a finite element with diameter which is used in the numerical method (intervals for , often triangles or quadrilaterals for , usually tetrahedra or hexahedra for , …). Let be some infinite-dimensional space of functions defined on such that the solution of the PDE is an element of . With we denote a family of finite-dimensional (usually closed) subspaces of whose dimension depends on ; the desired solution of the PDE is approximated by an element of . Moreover we have some given norms , and which are induced by certain inner products , and , which are used in the analysis of the numerical methods. In general the constants and of (1) and (2) depend on the diameter  and on the parameter  reflecting the dimension of the space . The dependence with respect to the diameter  is obtained by transforming Equations (1) and (2) to a reference domain , i.e.

 (3) ||^vn||X(^Ω) ≤^c1(n)||^vn||Y(^Ω)for all ^vn∈^Vn, (4) ||^vn||Z(∂^Ω) ≤^c2(n)||^vn||Y(^Ω)for all ^vn∈^Vn

and applying a scaling argument [22, 25, 2]. The more challenging problem is to find precise estimates for the constants and with respect to the parameter . The best possible constants by definition are given by

 (5) ^c1(n) =sup^vn∈^Vn||^vn||X(^Ω)||^vn||Y(^Ω)= ⎷sup^vn∈^Vn(^vn,^vn)X(^Ω)(^vn,^vn)Y(^Ω), (6) ^c2(n)

Introducing for the basis functions , i.e., , we further obtain

 (^c1(n))2=supv–n∈Rn(Knv–n,v–n)ℓ2(Mnv–n,v–n)ℓ2and(^c2(n))2=supv–n∈Rn(Lnv–n,v–n)ℓ2(Mnv–n,v–n)ℓ2,

with the symmetric and positive (semi-) definite matrices

 Kn(i,j):=(φj,φi)X(^Ω), Mn(i,j):=(φj,φi)Y(^Ω), and Ln(i,j):=(φj,φi)Z(^∂Ω)

for . Note that the rows of the matrices are naturally defined via the second argument of the inner product and the columns are given naturally via the first argument of the inner product. Hence, the constants and are given by the largest eigenvalues of the generalized eigenvalue problems

 (7) Knx––n=λnMnx––nandLnx––n=μnMnx––n,

i.e.

 (^c1(n))2=λnand(^c2(n))2=μn.

In this work we want to solve the problems (3) and (4), i.e., determine estimates for and , for the reference domain with

 (u,v)X(^Ω) (u,v)Y(^Ω) =∫^Ωu(x,y)v(x,y)dxdy,

for , where is the space of polynomials of degree less than , i.e.

 ^Vn={xiyj:0≤i,j

In Section 2 we state the problem in detail and derive its formulation as a generalized eigenvalue problem of the form (7). The difficulty is now to find an accurate estimate for the largest eigenvalues and for general parameter . Of course one can compute the eigenvalues exactly for a given fixed parameter , which is done for example in [19]. But to derive their exact values or precise estimates for a general parameter  one needs techniques from symbolic computation. In Section 3 we use the HolonomicFunctions package [11, 12] to prove a closed-form representation of the characteristic polynomial of our eigenvalue problem, in the spirit of the holonomic ansatz [27] for evaluating determinants. The holonomic ansatz is a very powerful method (it was the key to the “holy grail of enumerative combinatorics” [14]), and a very flexible one, too [16]. For our purposes here we had to adapt this algorithm, which led to a new variant that is applicable to a larger class of determinants. In Section 4, our representation of the characteristic polynomial is used to derive and prove the estimates from below and above

 14√n(n−1)(n+1)(n+2)≤^c1(n)≤12√2√n(n−1)(n+1)(n+2)

for the constant  in the inequality

 ||∂^x^un||L2(^Ω)≤^c1(n)||^un||L2(^Ω)for all ^un∈^Vn.

In Lemmas 4.5 and 4.8 we give much sharper estimates for ; as an application, they allow to tune the parameters of numerical methods precisely. In Section 5 the same representation is used to investigate the asymptotic behaviour of the eigenvalues; on the way we discover some interesting connections to the Taylor expansions of trigonometric functions. As an encore, we deal with the second inequality (4) in Section 6. It turns out that it is considerably simpler and we are able to derive the exact value of .

Throughout the paper, we employ the following notation: denotes the Pochhammer symbol, also known as rising factorial, defined for all nonnegative integers  by

 (a)n:=a⋅(a+1)⋯(a+n−1) for n>0and(a)0:=1.

We use for the floor function, and for the ceiling function, i.e., the largest integer below , resp. the smallest integer above . For a polynomial  we refer to the degree of  with respect to the variable  by . By we denote the Kronecker delta symbol, i.e., if and . If is the matrix , then we use for its determinant the short-hand notation . The determinant of the matrix is defined to be .

2. The Maximal Eigenvalue Problem

Let the reference domain be defined by , the open square of size  centered around the origin. For and , we define and to be the unique integers in satisfying . In other words,

 χn(k):=⌊k−1n⌋andρn(k):=k−1modn.

For the rest of this section we fix and write shortly and . By employing the standard monomial basis , we obtain the matrix with entries defined by

 (8) mi,j:=∫^Ωφiφjdxdt(1≤i,j≤n2)

and the matrix with entries

 (9) ki,j:=∫^Ω(∂xφi)(∂xφj)dxdt(1≤i,j≤n2).

Since and are just monomials, these integrals can be evaluated in a straightforward manner:

 mi,j =∫1−1(∫1−1xρ(i)tχ(i)xρ(j)tχ(j)dx)dt =∫1−11−(−1)ρ(i)+ρ(j)+1ρ(i)+ρ(j)+1tχ(i)+χ(j)dt =1−(−1)ρ(i)+ρ(j)+1ρ(i)+ρ(j)+1⋅1−(−1)χ(i)+χ(j)+1χ(i)+χ(j)+1.

Similarly

 ki,j =∫1−1(∫1−1ρ(i)ρ(j)xρ(i)−1tχ(i)xρ(j)−1tχ(j)dx)dt =∫1−1ρ(i)ρ(j)1−(−1)ρ(i)+ρ(j)−1ρ(i)+ρ(j)−1tχ(i)+χ(j)dt =ρ(i)ρ(j)1−(−1)ρ(i)+ρ(j)−1ρ(i)+ρ(j)−1⋅1−(−1)χ(i)+χ(j)+1χ(i)+χ(j)+1,

where we assumed that ; otherwise the integral equals to .

We are interested in computing the maximal such that . In the following we derive an equivalent formulation of this problem that involves smaller matrices. For this purpose let

 ai,j:=1−(−1)i+j−1i+j−1andbi,j:=(i−1)(j−1)1−(−1)i+j−3i+j−3,

such that the matrix entries and can be written as

 mi,j =aχ(i)+1,χ(j)+1⋅aρ(i)+1,ρ(j)+1 ki,j =aχ(i)+1,χ(j)+1⋅bρ(i)+1,ρ(j)+1.

This shows that the matrices and can be written as Kronecker products:

 Mn=An⊗AnandKn=An⊗Bn.

These representations as Kronecker products are quite natural since the used basis functions and the reference domain

itself have tensor product structure. In particular we then obtain,

 det(Kn−λnMn)=det(An⊗(Bn−λnAn))=det(An)ndet(Bn−λnAn)n.

So the problem is equivalent to computing the maximal such that

 det(Bn−λnAn)=0.

3. Determinant Evaluation

According to the previous discussion, we are now interested in evaluating the determinant

 det(Bn−λAn)=det1≤i,j≤n((1−(−1)i+j−1)((i−1)(j−1)i+j−3−λi+j−1))

for symbolic ; the desired maximal eigenvalue  is then just the largest root of the obtained polynomial. We see that the matrix has zeros at all positions for which

is an odd integer. By applying the permutation

to the rows and to the columns of the matrix, we decompose it into block form and obtain

 det(Bn−λAn)=2n∣∣ ∣∣A(0)⌊n/2⌋00A(1)⌈n/2⌉∣∣ ∣∣=2ndet(A(0)⌊n/2⌋)⋅det(A(1)⌈n/2⌉)

where the subscripts indicate the dimensions of the square matrices and , whose entries are independent of the dimension and given by

 (10) a(0)i,j :=(2i−1)(2j−1)2i+2j−3−λ2i+2j−1, (11) a(1)i,j :=4(i−1)(j−1)2i+2j−5−λ2i+2j−3.

Hence the matrices and start as follows:

 A(0) =⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝1−λ31−λ51−λ71−λ9⋯1−λ595−λ7157−λ973−λ11⋯1−λ7157−λ9259−λ113511−λ13⋯1−λ973−λ113511−λ134913−λ15⋯⋮⋮⋮⋮⋱⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠, A(1) =⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝−λ−λ3−λ5−λ7⋯−λ343−λ585−λ7127−λ9⋯−λ585−λ7167−λ983−λ11⋯−λ7127−λ983−λ113611−λ13⋯⋮⋮⋮⋮⋱⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠.
Theorem 3.1.

Let and be defined as in (10) and (11), then the following identities hold for all nonnegative integers :

 detA(0)n =det1≤i,j≤na(0)i,j=(−1)nh(0)n⋅F2n(λ), detA(1)n =det1≤i,j≤na(1)i,j=(−1)nh(1)n⋅λF2n−1(λ),

where

 (12) Fn(λ) :=ν∑j=0(−4)j−ν(2ν−2j+1)n(2j−2ν+n)!λjwithν=ν(n):=⌊n2⌋, (13) h(ℓ)n :=12nn∏i=1((i−1)!)2(i−ℓ+12)n.
Corollary 3.2.

For all nonnegative integers  we have

 det(Bn−λAn)=(−2)nh(0)⌊n/2⌋h(1)⌈n/2⌉λFn−1(λ)Fn(λ).

The key ingredient for the proof of Theorem 3.1 is the following lemma which shows that the quantities and defined there are basically the entries of the last column of the inverses of and , respectively.

Lemma 3.3.

With

 p(0)n,j =22n+2j−3(32)2n−1(n+12)j−1(n−1)!(2j−1)!n−1∑m=02n−2m−2∑k=0(−1)j+m(2m+1)2kλm4m+kk!(2m+k−n−j+2)!, p(1)n,j =4j−n(4n−3)!(n−12)j−1(2n−2)!(n−1)!(2j−2)!n−1∑m=02n−2m−2∑k=0(−1)j+m(2m)2kλm4m+kk!(2m+k−n−j+2)!,

the following identities hold for all nonnegative integers  and for :

 n∑j=1a(0)i,jp(0)n,j =δi,nF2n(λ), n∑j=1a(1)i,jp(1)n,j =δi,nλF2n−1(λ).
Proof.

These identities can be proven routinely using the holonomic systems approach [26]. We have carried out the necessary calculations using the HolonomicFunctions package [11, 12]. The results are documented in the supplementary electronic material [13].

First we derive, using holonomic closure properties and creative telescoping, a (left Gröbner) basis for the set of recurrence equations that satisfies. Again applying closure properties (in this case for multiplication) one obtains recurrences for the product , and by creative telescoping, for its definite sum, which we denote by (it is the left-hand side of the first identity). Here we face the problem of poles inside the summation range that are introduced by the certificate of the telescopic relation. We solve this issue by constructing a different certificate, free of the problematic denominators, using an ansatz reminiscent of the polynomial ansatz [11, Sec. 3.4]. The recurrences for have the following form (some polynomial coefficients are omitted for space reasons):

 16n4(n+1)2(2n+1)4(2n+3)2(4n+1)(i−n+1)2(2i+2n+3)2 ×(4i2+2i+λ−4n2−2n)2li,n+2=(⋯)li+1,n+(⋯)li,n+1+(⋯)li,n, 2n(2n+1)(i−n+1)(2i+2n+3)(4i2+2i+λ−4n2−2n)li+1,n+1= (⋯)li+1,n+(⋯)li,n+1+(⋯)li,n, 2(n−1)n(2n−1)(2n+1)(4n+1)2(4n+3)(i−n+1)2(i−n+2)(2i+2n+3) ×(4i2+2i+λ−4n2−2n)li+2,n=(⋯)li+1,n+(⋯)li,n+1+(⋯)li,n.

From their support and their leading coefficients it becomes clear that when we want to use them to compute for all , then we have to give the initial conditions , , , and . By verifying that they all equal  we have shown that the first identity holds for .

For we can construct, by holonomic substitution, a univariate recurrence satisfied by . It turns out that the corresponding operator is a left multiple of the second-order operator that annihilates . Also in this case, the proof can be completed by checking a few initial conditions. The proof of the second identity is established in an analogous way. ∎

Lemma 3.4.

The following determinant evaluations hold for all nonnegative integers :

 det1≤i,j≤n(12i+2j−1) =12nn∏i=1((i−1)!)2(i+12)n(=h(0)n), det1≤i,j≤n(12i+2j−3) =12nn∏i=1((i−1)!)2(i−12)n(=h(1)n).
Proof.

These determinants are special cases of Cauchy’s classic double alternant [3]

 det1≤i,j≤n(1xi+yj)=∏1≤i

where are indeterminates; see also [17, Thm. 12, Eq. (5.5)] and [18, Thm. 15] for a proof by factor exhaustion. In order to obtain the first assertion, we specialize and for , and obtain:

 det1≤i,j≤n(12i+2j−1) =∏1≤i

The second assertion is derived in a completely analogous way. Note also that these two determinants can be proven routinely using the holonomic ansatz [27]. ∎

Proof of Theorem 3.1.

Lemma 3.3

shows that the vector

is, up to a scalar multiple, the -th column of for . Since the entries of this vector (and of course those of the matrices itself) are polynomials in , this shows that and that . Note that both polynomials and have degree  in . Next we argue that also the determinants of and have degree  in , which is the maximal possible—taking into account that the matrix entries are linear polynomials in . Observe that the matrix entries in Lemma 3.4 are precisely . Thus Lemma 3.4 implies that converges to a nonzero constant (only depending on ) as goes to infinity. Hence for , which means that the two determinants are now determined up to a multiplicative constant not depending on . By noting that the polynomials are monic and that the expressions given in Lemma 3.4 are, up to sign, the leading coefficients of and , respectively, the assertion of the theorem is proven. ∎

Note that our proof of the determinant evaluations in Theorem 3.1 is very reminiscent of Zeilberger’s holonomic ansatz [27]. In fact, the only difference is that we chose to normalize the vector in a different way as Zeilberger would do it: while he suggests the normalization , we normalize such that and is a monic polynomial with . (The same discussion applies to , of course.)

In the original formulation of the holonomic ansatz, i.e., with the normalization , the final result in the case of success is a holonomic recurrence, i.e., a linear recurrence with polynomial coefficients, for . However, this ansatz is not at all guaranteed to succeed: even if the matrix entries are holonomic, this doesn’t mean that the sequence of quotients of consecutive determinants is a holonomic sequence. The determinant of is such an example: the polynomials satisfy the second-order recurrence

 (4n+3)F2n+4(λ)+(4n+5)(16n2+40n−2λ+21)F2n+2(λ)+(4n+7)λ2F2n(λ)=0,

which means that (most likely) the quotient doesn’t satisfy a holonomic recurrence of any order. (We have strong evidence that this quotient is non-holonomic, but we haven’t tried to prove this rigorously.) Provided that this is true, the original holonomic ansatz must fail.

Thanks to the additional parameter that appears polynomially in the matrix entries, we can identify the determinant of in the denominators of the inverse matrix. Thus a natural normalization of the vector  would be such that . In that case, the final result would be a holonomic recurrence for ; hence this variant is applicable when the determinant itself is a holonomic sequence in . Unfortunately, that’s not the case for the matrix because of the non-holonomic prefactor . This explains why we had to choose yet another normalization, in order to separate the holonomic and the non-holonomic part of the determinant. For each part then we had to prove a different determinant evaluation: for the holonomic “polynomial part” this was done in Lemma 3.3, for the non-holonomic “constant part” in Lemma 3.4. It is not unlikely that there are many more examples of determinants where the original holonomic ansatz fails, but where the modifications described here lead to success.

At the end of this section we want to briefly discuss an alternative way to derive the polynomials . In our above considerations we started with the monomial basis when formulating the eigenvalue problem. Alternatively, one could employ the Legendre basis leading to the following determinant:

 Dn=det1≤i,j≤n(∫1−1P′i(x)P′j(x)dx−λ∫1−1Pi(x)Pj(x)dx)

(note that only the matrix entries on the main diagonal depend on ). Indeed any basis for the space can be used for the computation of the eigenvalues given in (7). So by construction, this determinant leads to the same family of polynomials , and in fact we have that does not depend on . Doing the same block decomposition as before, we obtain the two families of matrices

 (2m(2m+1)−δi,j2λ4i+1)1≤i,j≤nand(2m(2m−1)−δi,j2λ4i−1)1≤i,j≤n

where stands for , whose determinants are given by

 (−1)n2n(54)nF2n+1(λ)% resp.(−1)n2n(34)nF2n(λ).

Note that these determinants are “nicer” than the ones we considered above, because their leading coefficients form holonomic sequences (actually they are hypergeometric). So it seems that we should have started with this formulation. But there is also a drawback: the matrix entries are defined in terms of , which on the one hand yields nicely structured matrices (constant along “hooks”, with a perturbation on the diagonal) such as

 ⎛⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜ ⎜⎝2−2x32222…212−2x7121212…21230−2x113030…2123056−2x1556…212305690−2x19…⋮⋮⋮⋮⋮⋱⎞⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟ ⎟⎠,

but on the other hand requires case distinctions that make the proofs of the relevant identities (the analog of Lemma 3.3) more complicated.

4. Upper and Lower Bounds on the Maximal Root of Fn(λ)

In this section we give lower and upper bounds on the maximal root of . Recall that we are interested in the maximal root of

 det(Bn−λAn)=cnλFn(λ)Fn−1(λ).

We will prove that the maximal root of is equal to the maximal root of . We prove this in Lemma 4.6 which is based on Lemmas 4.2, 4.3, and 4.5, which are technical in nature. A lower and an upper bound on the maximal root of are given in Lemma 4.5. A better upper bound is given in Lemma 4.8. These two lemmas are based on Lemma 4.3. Recall the definition of .

Definition 4.1.

To simplify notation in this section, we introduce the polynomials

 (14) fj(n):=(n−2j+1)4j4j(2j)!,

which correspond (up to sign) to the coefficients of :

 Fn(λ)=ν(n)∑j=0(−1)jfj(n)λν(n)−j=