Efficient Approaches for Verifying the Existence and Bound of Inverse of Linear Operators in Hilbert Spaces

This paper describes some numerical verification procedures to prove the invertibility of a linear operator in Hilbert spaces and to compute a bound on the norm of its inverse. These approaches improve on previous procedures that use an orthogonal projection of the Hilbert space and its a priori error estimations. Several verified examples which confirm the effectiveness of the new procedures are presented.


Introduction
Let X , Y be complex Hilbert spaces endowed with inner products ( u, v ) X , ( u, v ) Y and norms u X = √ ( u, u ) X , u Y = √ ( u, u ) Y , respectively, and D(A) be a complex Banach space satisfying D(A) ⊂ X ⊂ Y . We assume that the embedding D(A) → X is compact and that there exists an embedding constant C p > 0 for X → Y such that u Y ≤ C p u X , ∀u ∈ X . (1) Here, a concrete value of C p satisfying (1) should be estimated. We consider linear operators A : D(A) → Y , Q : X → Y , and The aim of this paper is to present computer-assisted procedures for proving the invertibility of L defined by (2) and for computing a constant M > 0 satisfying in a mathematically rigorous sense. The constant M in (3) represents an upper bound on the operator norm L −1 L(Y ,X ) . For example, when one aims at proving the existence of a solution of a nonlinear functional equation with a mathematically rigorous error bound by using Newton-type or Newton-Kantorovich-type arguments (e.g., [14,17]), the operator L is the linearization of a given nonlinear problem. Therefore, the invertibility of L and bounding the norm for L −1 play essential roles in most computer-assisted rigorous approaches to nonlinear problems, and it is also desirable to obtain M in (3) as small as possible [7, 13-15, 17, 25]. We note that determining the existence and norm bounds for L −1 is also applied to eigenvalue exclosure for self-adjoint or non-self-adjoint eigenvalue problems in Hilbert spaces [11,23]. Concerning nonlinear problems in Hilbert or Banach spaces, we are able to refer to powerful computer-assisted proofs for the ordinary/partial differential equations: Gameiro-Lessard [1], Arioli-Koch [2], Day-Lessard-Mischaikow [4], Figueras-Gameiro-Lessard-de la Llave [5], Nagatou-Plum-McKenna [12], Hungria-Lessard-Mireles James [8], van den Berg-Williams' [3], Oishi [16], and references therein.
We have previously proposed two types of numerical verification approaches [10,24] to assure the invertibility of L defined by (2) and to bound M in (3). These are based on orthogonal projections to Galerkin approximations with constructive a priori error estimations and can be applied to the case in which the operator L is non-self-adjoint.
The first approach in [10] transforms the problem L ψ = φ for φ ∈ Y into an equivalent fixed-point problem on X , constructs a validated bound M, and verifies the invertibility of L . Although numerous computer-assisted proofs have demonstrated the effectiveness of this approach [13,14,23], it has a restriction such that M in (3) does not converge to the exact operator norm L −1 L(Y ,X ) even though the dimension of the Galerkin space increases.
In [24], some of the authors considered another estimation for M in (3). This second approach avoided the fixed-point formulation and showed that M is expected to converge, as the the dimension of the Galerkin space increases, to its exact operator norm of L −1 under suitable assumptions. However, the invertibility criterion of the second approach in [24] is sometimes more difficult compared to that of the first approach, in which case it requires a greater computational cost [10,23].
The aim of the present paper is to propose some improved procedures relative to those in [10,24] obtaining (3), as well as new invertibility criteria for L . The essential principle is, in fact, the same as that in [21] for second-order linear elliptic operators. However, the proposed procedures in this paper are a generalization of [21], and this generalization allows an extension to operators in Hilbert spaces other than second-order elliptic operators (e.g., see Sect. 5.2). These procedures are based on a projection and constructive a priori error estimations in Hilbert spaces as well as a computable upper bound M A > 0 satisfying which assures the invertibility of L under some assumptions (see the proof of Theorem 3). Although, in the ideal situation (h → 0), our proposed approaches have the same bound M as in [10] or [24], some verification examples reveal that the proposed procedures have better bounds than the approach in [10] or [24]. We also note that the estimation (4) and the guaranteed upper bound M A should be helpful for invertibility confirmation such as the eigenvalue excluding approach in [11] or the Newton-Kantorovich-type verification approach in [17].
Recently, Sekine-Nakao-Oishi [20] introduced the following 2 × 2 block operator matrix on and proposed an interesting verification technique by effective use of D. The method is based on the fact that, under the condition that the finite-dimensional operator T : X h → X h is invertible, if its Schur complement defined by S = W − ZT −1 Y : X * → X * is also invertible, then D itself is reversible. By using this principle, the solution of L ψ = φ for φ ∈ Y can be enclosed without estimating the norm of L −1 directly, in a way similar to the Gaussian elimination method for simultaneous linear equations. In [20], they also presented a formulation of the numerical verification of the solution for nonlinear elliptic problems according to this methodology, which can be considered as an alternative approach to the method described in this paper. In the future, it will be important to compare the strengths and weaknesses of such verification methods with those described here.
The remainder of the present paper is organized as follows. The next section is devoted to the description of assumptions on the given linear operators and some finite-dimensional approximation subspaces with related constants. We summarize our previous results based on [10,24] on the invertibility of L and M. Section 3 proposes two constructive norm estimations of (4), as well as the invertibility criteria of L . In Sect. 4 we propose a new bound on the norm of the inverse of L under the invertibility condition described in the previous section and present various considerations. Several verification examples of the proposed procedures comparing them with the previous approaches are reported on in the final section.

Galerkin Approximation, Related Constants, and Known Results
This section describes assumptions on the linear operator L defined by (2) and introduces a finite-dimensional approximation subspace with related constants. Assume that the operator A has the following properties A1 and A2.

A2 It holds that
A typical example of A1 is when A is the Laplacian and A2 is derived from partial integration (see Sect. 5 or [14,Chapter 4]). From A1, A2, and (1), we have We define a finite-dimensional approximation subspace X h ⊂ X dependent on the parameter h > 0. For example, in the case of a partial differential equation problem, X h is a finite element subspace with mesh size h. We also define the orthogonal projection P h : X → X h and the projection P 0 : Y → X h by respectively. Because X h is a closed subspace of X , each element u ∈ X can be uniquely decomposed as u = u h + u * for u h ∈ X h and u * ∈ X * := (I − P h )X . Now we assume that P h , A, and Q have the properties A3, A4, and A5.
A3 There exists C(h) > 0 satisfying C(h) → 0 as h → 0 and A4 There exist τ 1 > 0 and τ 2 > 0 such that A5 There exists τ 3 > 0 satisfying Assumption A3 corresponds to error estimation of the orthogonal projection P h defined by (7). We emphasize that the estimation (9) is indispensable in our argument, and the compactness of the embedding D(A) → X plays an essential role in obtaining the constant C(h) with the desired properties. Assumptions A4 and A5 indicate detailed information about the boundedness of the operator Q. We note that concrete values of the constants C(h) and τ i (i = 1, 2, 3) have to be known and have to be evaluated rigorously, and the constants τ i depend on Q and h. Concrete examples of C(h) and τ i are shown in Sect. 5. We define a basis function of is positive definite and can be decomposed as where H denotes conjugate transposition. For example, L i is taken to be the Cholesky factor of A i ; i.e., L i is a lower triangular matrix. By using L i , for u h ∈ X h and u = [u i ] ∈ C N , it holds that We also assume that G has an inverse and we define positive values ρ,ρ by respectively. In actual computation, ρ andρ are obtained as upper bounds of matrix 2-norms, and evaluation of them, including a proof of the invertibility of G, can be reduced to the verified computation of the maximum singular value of a matrix [19]. We next present two previous computer-assisted proofs of the invertibility of L , as well as the computable bound of M in (3).

Theorem 1 [10, Theorem 3] If
then L defined by (2) has an inverse and M > 0 for (3) is obtained by Remark 1 Theorem 1 is based on a fixed-point problem on X and it is an improved version of that appearing in [23, Theorem 5.1]. In Theorem 1, assuming that κ and C(h) converge to 0 as h → 0 implies that M for (18) converges to C p ρ. Generally, M/(C p ρ) → 1 does not indicate convergence to the exact operator norm L −1 L(Y ,X ) .

Remark 2
The matrix G is the Galerkin approximation of the operator L andρ in (16) reflects an approximation of (20) is expected to converge to the exact operator norm of L −1 , as h → 0. However, it has been reported that sometimes the criterion κ < 1 is harder to satisfy than κ < 1 for fixed h, experimentally [14].

Invertibility Verification of L and Constructive Norm Bounds for A
In this section, we construct a computable upper bound M A > 0 satisfying (4). The inequality (4) and A1 give that L is one-to-one. Then, by using the compactness of an operator A −1 Q on X and the Fredholm alternative, the invertibility of L is assured (see the proof of Theorem 3, [10, proof of Theorem 3], [23, proof of Theorem 5.1]). Let us define an Hermitian and positive semidefinite N × N matrix E by and the related matrix 2-norms by respectively. Note that if the matrix E given by (21) can be decomposed as E = L 3 L H 3 , then ρ L 1 and ρ L 2 can be written as We also note that the matrix E appears in the actual computations since it derives the quadrature form of estimation such that Using ρ L 1 , ρ L 2 ,ρ, and ρ, we have the following norm estimations.

Theorem 3 If it holds that
Proof Since each u ∈ D(A) can be uniquely decomposed as u = u h + u * with u h = P h u and u * = (I − P h )u, setting f = L u we have Using A3, A4, (24), (6), and A5, we have the following estimation: Therefore, under the assumption κ L 1 < 1, the estimation (4) holds with (31). In particular, and A1 show that L is one-to-one. Furthermore, for any given φ ∈ Y , the problem is equivalent to Since A −1 Q : X → X is compact, the Fredholm alternative holds for problem (34), whereby L being one-to-one implies that (34), and hence (33), is uniquely solvable. Therefore, L is bijective and hence invertible.
We show the second approach for (4), which assures the invertibility of L . (4) is given by

Theorem 4 If it holds that
Proof Since each u ∈ D(A) can be uniquely decomposed as u = u h + u * , where u h = P h u and u * = (I − P h )u, by setting f = L u, we obtain the same inequality (32) as in the case of Theorem 3. By using A3, A4, and (25), we get the following: Then under the assumption κ L 2 < 1, the estimation (4) holds with (36).

Upper Bound M > 0
This section presents some formulas for upper bound M satisfying (3) which are expected to be an improvement on previous results (18) and (20).
Applying the invertibility theorems for L described in Sect. 3, we have the following theorem.

Theorem 5 Assume that there exists a constant M
Proof For each f ∈ Y , we set u = L −1 f ∈ D(A) and represent u as u = u h + u * with u h = P h u ∈ X h and u * = (I − P h )u ∈ X * . For the finite-dimensional term, inequalities (26), (6), and A5 imply that and inequalities (27) and A4 imply that For the infinite-dimensional term, from A3 and (4), we have Therefore, by substituting (41) into (39) and (40), it follows that respectively, and then u 2 X = u h 2 X + u * 2 X implies the desired conclusions. Remark 4 When h → 0, the estimation in (37) shows that M/(ρ C p ) → 1 and the estimation in (38) shows that M/ρ → 1. Therefore, in the ideal situation (h → 0), we have the same bound M, namely, (37) versus (18) and (38) versus (20).

Remark 5
In order to obtain bound M A > 0 in (4) by Theorems 3 or 4, the additional computational cost of the 2-norm of ρ L 1 or ρ L 2 is required. According to our numerical examples in the next section, the additional costs for ρ L 1 or ρ L 2 could be seen as approximately the same order of magnitude for ρ orρ.

Verification Examples
This section reports on several computer-assisted results for the estimations of M in (3) proposed by Theorems 1, 2 and 5, as well as the invertibility verification of L .
All computations were carried out on the Fujitsu PRIMERGY CX2570 M4; Intel Xeon Gold 6140 (Skylake-SP); 2.3 GHz (Turbo 3.7 GHz) by using INTerval LABoratory Version 11, a toolbox in MATLAB R2019a (9.7.0.1261785) 64-bit (glnxa64) developed by Rump [18] for self-validating algorithms. Therefore, all numerical values in these tables are verified data in the sense of mathematically strict rounding error control.

Second-order Differential Operators
For a unit square region = (0, 1) × (0, 1) and some integer m, let H m ( ) denote the complex L 2 -Sobolev space of order m on . We define the Hilbert space H 1 0 ( ) := {u ∈ H 1 ( ) | u = 0 on ∂ } with inner product ( ∇u, ∇v ) L 2 ( ) and norm u H 1 0 ( ) := ∇u L 2 ( ) , where ( u, v ) L 2 ( ) denotes the L 2 -inner product on . We consider a twodimensional second-order non-self-adjoint operator respectively. We are able to set It is well known that A1 holds [6], and A2 is obtained by partial integration. We form a linear and uniform triangular mesh on with the side length of each finite element being h > 0 and set X h as a classical P1 finite element subspace of H 1 0 ( ). Then, for A3, P h is now the usual H 1 0 -projection, and it is well known that (7) holds with C(h) = 0.493h and C p = 1/(π √ 2) [14].
Concerning A4 and A5, we can take Note that when b is differentiable, we can derive more accurate estimates for τ i (i = 1, 2, 3) [21]. For more general setting of the second-order operator see, e.g., [14].

Example 1
Our first verification example is a non-self-adjoint case of  Table 3, and the symbol "-" indicates that the invertibility criterion does not hold. The value M A adopted is the smaller of the results from (31) and (36). The computational cost is concentrated in the singular value computations; elapsed times for obtaining ρ,ρ, ρ L 2 , and ρ L 1 are 19.08, 16.05, 17.12 and 17.44 s, respectively, for 1/h = 40. In Table 1, (37) gives the best result. This may be because the value of C p ρ is close toρ which is the approximation of the exact operator norm L −1 L(L 2 ,H 1 0 ) . On the other hand, the difference between (37) and (38) tends to get smaller as the value of 1/h increases.

Example 2
The next example is the case of a linearized equation of semilinear PDEs:  Values in bold are the smallest M Equation (46) comes from a reaction-diffusion equation and, by the Newton-Raphson method with usual floating-point arithmetic, we have (at least) three finite element approximate solutions denoted by u h . We refer to these approximate solution as "1st," "2nd," and "3rd" Fig. 1 shows the shape of u h . The linearized operator L at u h is defined by which is the case of b = 0. Tables 3, 4    The computational cost is concentrated in the singular value computations; elapsed times for obtaining ρ,ρ, ρ L 2 , and ρ L 1 are 5.81, 3.82, 4.29 and 4.46 s, respectively, for 1/h = 40.

Example 3
We take Re = 5776 and a = 1.019. Tables 6, 7 and 8 show the verification results of the operator (48) for μ = −100 + 1552.59i, μ = −200 + 1552.59i, and μ = −500 + 1552.59i, respectively. Let us remark that, especially for h smaller than 1/400, the round-off error for ρ and ρ L 1 tends to be a serious problem (see also [23]). For each case, the effectiveness of our proposed procedure and the validity of Theorem 5 are also shown. In both Tables 6, 7 and 8, it can be seen that, as 1/h increases, (38) gives the best results and tends to approach the value of ρ, which is an approximation of the operator norm L −1 L(L 2 ,H 2

Remark 6
In the verification examples treated in the present paper, we focus on linear operators only. We emphasize again that also in the context of computer-assisted proofs for nonlinear equations including ordinary/partial differential equations, the verification of the invertibility of L and the computation of a norm bound on L −1 plays an essential role, and our proposed approaches are expected to provide accurate and efficient enclosure results for solutions of nonlinear problems. We will report on concrete computer-assisted proofs of nonlinear equations using our proposed method in subsequent papers.

Conclusion
We proposed some improvements on the computer-assisted procedure for verifying the invertibility of a linear operator in a Hilbert space and the computation of a norm bound of its inverse.
Although an additional computation of a 2-norm bound is required, verification examples confirm that the upper bound given by (37) improves on (18) and the upper bound given by (38) improves on (20) except for the case shown in Table 3. The invertibility criteria in Theorems 3 and 4 and the guaranteed upper bound M A should be helpful for invertibility confirmation such as eigenvalue exclosure for self-adjoint or nonself-adjoint eigenvalue problems in Hilbert spaces in [11,23]. We conclude that an appropriate procedure can be selected taking into consideration the given problem or the computational cost, for example.