# Verified Error Bounds for Real Eigenvalues of Real Symmetric and Persymmetric Matrices

Authors: Zhe Li; Xueqing Wang
Journal: Mathematical Problems in Engineering
 (2020)

---

## Abstract

This paper mainly investigates the verification of real eigenvalues of the real symmetric and persymmetric matrices. For a real symmetric or persymmetric matrix, we useeig code in Matlab to obtain its real eigenvalues on the basis of numerical computation and provide an algorithm to compute verified error bound such that there exists a perturbation matrix of the same type within the computed error bound whose exact real eigenvalues are the computed real eigenvalues.


---

## Body

## 1. Introduction

The eigenvalues of a matrix are one of the important tools to solve the complex mathematical problems in the fields of image processing, quantum mechanics, chemistry, and so on [1–3]. The matrices involved in practical problems often have algebraic structure. The preservation of matrix algebraic structure is helpful to keep the physical background of the matrix eigenvalue problem [4, 5]. For a matrix with special algebraic structure, the computation of its eigenvalues is very important in practical problems.Many scholars have carried out a lot of work to compute the eigenvalues of matrices with special algebraic structure. Using QR decomposition, Bunse-Gerstner et al. [6] designed a stable algorithm to compute the eigenvalues of structured matrices. Higham and Higham [7, 8] estimated the backward error bound and the condition number for computing the generalized eigenvalue and generalized eigenvector of structured matrix. Based on interval calculation, Rump [9] provided the sensitivity analysis of eigenvalues of structured matrix under structured perturbation and proved that, for circulant matrix, Toeplitz matrix, symmetric Toeplitz matrix, and symmetric matrix, the structured condition number is equal to the unstructured condition number under normwise perturbations. Byers and Kressner [10] used the condition number of invariant subspace of the structured matrix to compute the eigenvalues of the structured matrix. Shiozaki [11] proposed an effective algorithm to compute the structured eigenvalues in the field of quantum chemistry. Alon et al. [12] studied the concentration of the maximum eigenvalue of random symmetric matrix whose diagonal and upper diagonal entries are independent real random variables. By Gaussian probability density function with the same mean and variance, Edwards and Jones [13] proposed a straightforward method to analyse the characteristic spectrum of the large symmetric matrix. Given a symmetric matrix whose entries depend on a parameter, Hiriart-Urruty and Ye [14] investigated the first-order sensitivity of all the eigenvalues. Hladík et al [15] considered the eigenvalue problem about the symmetric matrix with perturbed interval entries. Hernandez et al. [16] proposed a greedy algorithm to compute the selected eigenpairs of a large sparse symmetric matrix by exploiting localization features of the eigenvector. Reid [17] showed some useful eigenvalue and eigenvector properties of symmetric and persymmetric matrices.In this paper, we use Rump interval method and Kantorovich theorem to compute the verified error bounds of real eigenvalues of the given real symmetric and persymmetric matrices, such that there exists a perturbation matrix of the same type within the computed error bound whose exact real eigenvalues are the computed real eigenvalues. To be precise, we transform the verification of real eigenvalues of real symmetric and persymmetric matrices into the verification of a root of a nonlinear system. We utilize the Rump interval method to compute the constants appearing in Kantorovich theorem about the nonlinear system and then use Kantorovich theorem to compute verified error bound of the zero vector as an approximate solution.The paper is organized as follows. Section2 is devoted as a preparation of this paper. The main theory and algorithm are, respectively, given in Section 3. Section 4 gives some examples to demonstrate the performance of our algorithm.


## 2. Notation and Preliminaries

Letℕ and ℝ denote the set of natural numbers and real numbers, respectively. For a matrix A∈ℝm×n, A: is a vector obtained by reshaping all elements of A into a single column vector, Ai1:i2, designates a submatrix of A by selecting from the i1th to i2th rows, and A:j1:j2 represents a submatrix of A by selecting from the j1th to j2th columns. Let Om,n denote an m×n zero matrix and In be the identity matrix of order n. For an m×n matrix A, let nullA denote the nullspace of A.Definition 1 (see [18]).
Given a matrixA∈ℝm×n with m≥n, the corank of A is defined by corankA=n−rankA. For a threshold δ>0, if the singular values σ1A,…,σnA of matrix A satisfy that(1)σ1A≥⋯≥σn−qA>δ≥σn−q+1A≥⋯≥σnA,then we say that A has numerical δ-corank q, denoted by corankδA=q.
LetIℝ represent the set of all intervals. A matrix and vector with interval entries are, respectively, called interval vector and interval matrix. Given an interval matrix A∈Iℝn×n if an arbitrary real matrix A satisfying A∈A is nonsingular, then we call the interval matrix A a nonsingular interval matrix. Rump [19] developed INTLAB toolbox in Matlab for interval arithmetic. For a nonlinear system if the Jacobian matrix of the system is Lipschitz continuous on some domain, Kantorovich [20] established Kantorovich theorem, which gives a sufficient condition to judge whether the Newton iterative method converges or not by the information of the initial approximate value on some domain.Lemma 1 (see [21]).
GivenA,R∈ℝn×n if the spectral radius of the matrix I−RA is less than 1, then A is nonsingular.Theorem 1 (see [22]).
Ifverifylss function in INTLAB runs successfully for a given interval matrix A∈Iℝn×n and interval right-hand side vector b∈Iℝn, then the computed interval X⊂Iℝn satisfies the following condition:(2)∑A,b≔x∈ℝn:Ax=b,∀A∈A,b∈b⊆X.Theorem 2 (see [20]).
Letf:ℝn⟶ℝn with f=f1,…,fnT, where f1,…,fn are continuously differentiable functions and f′x denotes the Jacobian matrix of the system fx=0. Suppose x˜∈ℝn is an approximate solution satisfying the condition that f′x˜ is invertible. Let B be a constant such that f′x˜−1≤B and κ a Lipschitz constant such that(3)f′u−f′v≤κu−v,u,v∈Ω,where Ω is a sufficiently large region containing x˜ and η is a constant such that(4)f′x˜−1fx˜≤η.
Ifh=Bκη≤1/2 and U¯x˜,ρ=x:x−x˜≤ρ⊂Ω for(5)ρ=1−1−2hhη,then there exists x^∈U¯x˜,ρ such that fx^=0.Remark 1.
As pointed by Rall [23], one may take the region Ω in Theorem 2 as U¯x˜,2η. If κ is the Lipschitz constant for this region, then U¯x˜,ρ⊂U¯x˜,2η if and only if h≤1/2.


### 2.1. Main Result

For a square matrixA if AT=A, then A is called a symmetric matrix. For a matrix A∈ℝn×n with entries ai,j, 1≤i,j≤n, define a subtranspose matrix As with entries an+1−j,n+1−i, 1≤i,j≤n. If As=A, then A is called a persymmetric matrix.LetAsym and Apersym, respectively, denote the n×n symmetric and persymmetric matrices. Let Esym and Epersym, respectively, represent the corresponding perturbation matrices, namely,(6)Esym=ε1,1ε1,2ε1,3…ε1,nε1,2ε2,1ε2,2…ε2,n−1ε1,3ε2,2ε3,1…ε3,n−2⋮⋮⋮⋱⋮ε1,nε2,n−1ε3,n−2…εn,1,Epersym=ε1,1ε1,2ε1,3…ε1,nε2,1ε2,2ε2,3…ε1,n−1ε3,1ε3,2ε3,3…ε1,n−2⋮⋮⋮⋱⋮εn,1εn−1,1εn−2,1…ε1,1.Letλ˜1sym,λ˜2sym,…,λ˜ksym be the set of all distinct eigenvalues of Asym computed by eig code in MATLAB. For each s=1,2,…,k, the singular value decomposition of Asym−λ˜ssymIn is(7)Asym−λ˜ssymIn=Uλ˜ssymΣλ˜ssymVλ˜ssymT=u1λ˜ssym,…,unλ˜ssymdiagσ1λ˜ssym,…,σnλ˜ssymv1λ˜ssym,…,vnλ˜ssymT.Assume that, for eachs=1,2,…,k, corankδAsym−λ˜ssymIn=qs, where δ is a positive real number that is very close to zero.Ifλ^1sym,λ^2sym,…,λ^ksym are all exact eigenvalues of Asym, then, for each s=1,2,…,k, the vectors(8)v1λ^ssym,v2λ^ssym,…,vn−qsλ^ssym,un−qs+1λ^ssym,un−qs+2λ^ssym,…,unλ^ssym,are linearly independent.Therefore, we can make the following reasonable assumption.Assumption 1.
For eachs=1,2,…,k, the vectors(9)v1λ˜ssym,v2λ˜ssym,…,vn−qsλ˜ssym,un−qs+1λ˜ssym,un−qs+2λ˜ssym,…,unλ˜ssym,are linearly independent.
Fors=1,2,…,k, define(10)Cssymεsym=Asym+Esym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs,where(11)εsym=ε1,1,ε1,2,…,ε1,n,ε2,1,…,ε2,n−1,…,εn−1,1,εn−1,2,εn,1T.Lemma 2.
For eachs=1,2,…,k if corankδAsym−λ˜ssymIn=qs, then the matrix(12)Asym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs,is nonsingular.Proof.
(13)Asym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs=Uλ˜ssymOqs,qsOqs,nIqsΣ1:n−qs,1:n−qsλ˜ssymOn−qs,qsOn−qs,qsOqs,n−qsΣn−qs+1:n,n−qs+1:nλ˜ssymIqsOqs,n−qsIqsOqs,qsV:,1:n−qsλ˜ssymTOn−qs,qsU:,n−qs+1:nλ˜ssymTOqs,qsΣn−qs+1:n,n−qs+1:nλ˜ssymIqs.
This completes the proof.For a persymmetric matrixApersym, let λ˜1persym,λ˜2persym,…,λ˜kpersym be the set of all distinct real eigenvalues of Apersym computed by eig code in MATLAB. Let J be the n×n inverse identity matrix; then, JApersym+Epersym−λ˜spersymIn is a symmetric matrix. For each s=1,2,…,k, the singular value decomposition of JApersym−λ˜spersymIn is(14)JApersym−λ˜spersymIn=Uλ˜spersymΣλ˜spersymVλ˜spersymT=u1λ˜spersym,…,unλ˜spersymdiagσ1λ˜spersym,…,σnλ˜spersymv1λ˜spersym,…,vnλ˜spersymT.Assume that, for eachs=1,2,…,k, corankδJApersym−λ˜spersymIn=qs, and then define(15)Cspersymεpersym=JApersym+Epersym−λ˜spersymInU:,n−qs+1:nλ˜spersymU:,n−qs+1:nλ˜spersymTOqs,qs,where(16)εpersym=εn,1,εn−1,1,…,ε1,1,εn−1,2,…,ε1,2,…,ε2,n−1,ε1,n−1,ε1,nT.For simplifying, we writeεstruct to stand for both εsym and εpersym and we write Csstruct to denote both Cssym and Cspersym. Then, define the constant K by(17)K=min1nn+qsCsstruct0−1∞,s=1,2,…,k.Lemma 3.
Ifεstruct∞<K, then, for each s=1,2,…,k, the matrix Csstructεstruct is nonsingular.Proof.
Ifεstruct∞<K, for each s=1,2,…,k,(18)In+1−Csstruct0−1Csstructεstruct2≤n+qsIn+qs−Csstruct0−1Csstructεstruct∞≤n+qsCsstruct0−1∞Csstruct0−Csstructεstruct∞<nn+qsεstruct∞Csstruct0−1∞<1.
From Lemma1, we can deduce that, for each s=1,2,…,k, the matrix Csstructεstruct is nonsingular.
For eachs=1,2,…,k, let Wsεstruct and Fsεstruct be, respectively, the first n-rows and the last qs-rows of the solution of the following linear system:(19)CsstructεstructWsεstructFsεstruct=On,qsIqs.
According to Cramer’s rule, we have the following easy lemma.Lemma 4.
For eachs=1,2,…,k, the matrix Fsεstruct is a symmetric matrix.Proposition 1.
Suppose thatε˜struct∞<K. If Fsε˜struct=Oqs,qs holds for s=1,2,…,k, then(20)corankAstruct+E˜struct−λ˜sstructIn=qs.Proof.
Firstly, consider the case for the symmetric matrix. Ifε˜sym∞<K, by Lemma 3, we know that, for each s=1,2,…,k, the matrix Cssymε˜sym is nonsingular. Next assume that s is an arbitrary fixed integer from 1 to k. If Fsε˜sym=Oqs,qs, then(21)Asym+Esym−λ˜ssymInWsε˜sym=On,qs,U:,n−qs+1:nλ˜ssymTWsε˜sym=Iqs.
Notice that the columns ofWsε˜sym are linearly independent. Hence, corankAsym+E˜sym−λ˜ssymIn≥qs. Assume that corankAsym+E˜sym−λ˜ssymIn>qs; then, there exists a nonzero vector βs∈nullAsym+E˜sym−λ˜ssymIn such that the matrix Wsε˜symβs is full rank, and then, for a nonzero vector bs∈ℝqs+1, the following equality(22)Asym+Esym−λ˜ssymInU:,n−qs+1:nλ˜ssymTWsε˜symβsbs=0,holds, which leads to a conditions. Consequently, corankAsym+E˜sym−λ˜ssymIn=qs.
Consider the case of persymmetric matrix. According to the above conclusion, we know thatcorankJApersym+E˜persym−λ˜spersymIn=qs. Apparently, corankApersym+E˜persym−λ˜spersymIn=qs.Lemma 5.
The following nonlinear system,(23)Gεstruct=F1εstruct1:qs,1F2εstruct2:qs,2⋮Fkεstructqs:qs,qs=0,is an underdetermined nonlinear system.Proof.
System (23) consists of nn+1/2 variables, and ∑i=1sqiqi+1/2 equations. Since(24)q12+q22+⋯+qs2+q1+q2+⋯+qs2≤q12+q22+⋯+qs2+n2=q1+q2+⋯+qs2−2∑1≤i≤j≤sqiqj+n2≤n2+n2−∑1≤i≤j≤sqiqj≤nn+12,then the proof is completed.
In the following, assume thats is an arbitrary fixed integer from 1 to k. By Lemma 3, there exists a neighborhood with the center of the origin such that, for an arbitrary εstruct belonging to this neighborhood, the matrix Csstructεstruct is nonsingular. Moreover, in this neighborhood, each entry of the solution of system (19) has the partial derivative with respect to each variable.
Concerning symmetric matrix, for each pairi,j with 1≤i≤n, 1≤j≤n−i+1, differentiating both sides of system (19) with respect to εi,j gives the following system:(25)Asym+Esym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs∂Wsε∂εi,j∂Fsε∂εi,j=t−Oj+i−2,qsWsεi,:On+qs−j−i+1,qs−Oi−1,qsWsεj+i−1,:On+qs−i,qs,where(26)t=12,j=1,1,2≤j≤n−i+1.
For pairsi1,j1 and i2,j2 with 1≤i1≤n, 1≤j1≤n−i1+1, 1≤i2≤n, and 1≤j2≤n−i2+1, differentiating both sides of system (25) with respect to εi1,j1 and εi2,j2, produces the following system:(27)Asym+Esym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs∂2Wsε∂εi1,j1∂εi2,j2∂2Fsε∂εi1,j1∂εi2,j2=t1−Oj1+i1−2,qs∂Wsε∂εi2,j2i1,:On+qs−j1−i1+1,qs−Oi1−1,qs∂Wsε∂εi2,j2j1+i1−1,:On+qs−i1,qs+t2−Oj2+i2−2,qs∂Wsε∂εi1,j1i2,:On+qs−j2−i2+1,qs−Oi2−1,qs∂Wsε∂εi1,j1j2+i2−1,:On+qs−i2,qs,where(28)t1=12,j1=1,1,2≤j1≤n−i1+1,t2=12,j2=1,1,2≤j2≤n−i2+1.
Considering the persymmetric matrix, for each pairi,j with 1≤i≤n and 1≤j≤n−i+1 differentiating both sides of system (19) with respect to εi,j yields the following system:(29)JApersym+Epersym−λ˜spersymInU:,n−qs+1:nλ˜spersymU:,n−qs+1:nλ˜spersymTOqs,qs∂Wsε∂εi,j∂Fsε∂εi,j=t−Oj−1,qsWsεn−i+1,:On+qs−j,qs−On−i,qsWsεj,:Oqs+i−1,qs,where(30)t=12,j=n−i+1,1,1≤j≤n−i.
For pairsi1,j1 and i2,j2 with 1≤i1≤n, 1≤j1≤n−i1+1, 1≤i2≤n, and 1≤j2≤n−i2+1, differentiating both sides of system (29) with respect to εi1,j1 and εi2,j2, provides the following system:(31)JApersym+Epersym−λ˜spersymInU:,n−qs+1:nλ˜spersymU:,n−qs+1:nλ˜spersymTOqs,qs∂2Wsε∂εi1,j1∂εi2,j2∂2Fsε∂εi1,j1∂εi2,j2=t1−Oj1−1,qs∂Wsε∂εi2,j2n−i1+1,:On+qs−j1,qs−On−i1,qs∂Wsε∂εi2,j2j1,:Oqs+i1−1,qs+t2−Oj2−1,qs∂Wsε∂εi1,j1n−i2+1,:On+qs−j2,qs−On−i2,qs∂Wsε∂εi1,j1j2,:Oqs+i2−1,qs,where(32)t1=12,j1=n−i1+1,1,1≤j1≤n−i1,t2=12,j2=n−i2+1,1,1≤j2≤n−i2.
Solving (25) and (27) or (29) and (31) at εstruct=0 can obtain Jacobi matrix G′0 and Hessian matrix G″0.Assumption 2.
The Jacobi matrixG′0 is full rank.
Suppose that, for an index setℐ=i1,j1,i2,j2,…,im,jm, the following matrix,(33)∂G0∂εi1,j1,∂G0∂εi2,j2,…,∂G0∂εim,jm,is nonsingular, and then define the following nonlinear system,(34)G˜εi1,j1,εi2,j2,…,εim,jm=0,by(35)G˜εi1,j1,εi2,j2,…,εim,jm=Gεstruct|εi,j=0fori,j∉ℐ.Remark 2.
By Theorem2, we compute verified error bound when the zero vector is the solution of the nonlinear system (34).
For system (34), define the constants B and η in Theorem 2 by(36)B=G˜′0−1∞,η=G˜′0−1G˜0∞.
Define anm-dimensional interval vector Ω=Ωi1,j1,Ωi2,j2,…,Ωim,jmT with each entry as interval −2η,2η. Define the corresponding interval perturbation matrix Ω˜ by setting Ω˜i,j=Ωi,j for i,j∈ℐ and Ω˜i,j=0 for i,j∉ℐ.
Byverifylss function, we solve interval linear systems (19), (25), (27) or (19), (29), (31) at E=Ω˜, and then we can obtain an interval tensor Hs,t:1≤s≤t≤m that satisfies(37)Hs,t⊃∂2GΩ∂εis,js∂εit,jt,1≤s≤t≤m.
Define the Lipschitz constantκ in Theorem 2 by(38)κ=max1≤s≤t≤mmaxH∞:∀H∈Hs,t.


## 2.1. Main Result

For a square matrixA if AT=A, then A is called a symmetric matrix. For a matrix A∈ℝn×n with entries ai,j, 1≤i,j≤n, define a subtranspose matrix As with entries an+1−j,n+1−i, 1≤i,j≤n. If As=A, then A is called a persymmetric matrix.LetAsym and Apersym, respectively, denote the n×n symmetric and persymmetric matrices. Let Esym and Epersym, respectively, represent the corresponding perturbation matrices, namely,(6)Esym=ε1,1ε1,2ε1,3…ε1,nε1,2ε2,1ε2,2…ε2,n−1ε1,3ε2,2ε3,1…ε3,n−2⋮⋮⋮⋱⋮ε1,nε2,n−1ε3,n−2…εn,1,Epersym=ε1,1ε1,2ε1,3…ε1,nε2,1ε2,2ε2,3…ε1,n−1ε3,1ε3,2ε3,3…ε1,n−2⋮⋮⋮⋱⋮εn,1εn−1,1εn−2,1…ε1,1.Letλ˜1sym,λ˜2sym,…,λ˜ksym be the set of all distinct eigenvalues of Asym computed by eig code in MATLAB. For each s=1,2,…,k, the singular value decomposition of Asym−λ˜ssymIn is(7)Asym−λ˜ssymIn=Uλ˜ssymΣλ˜ssymVλ˜ssymT=u1λ˜ssym,…,unλ˜ssymdiagσ1λ˜ssym,…,σnλ˜ssymv1λ˜ssym,…,vnλ˜ssymT.Assume that, for eachs=1,2,…,k, corankδAsym−λ˜ssymIn=qs, where δ is a positive real number that is very close to zero.Ifλ^1sym,λ^2sym,…,λ^ksym are all exact eigenvalues of Asym, then, for each s=1,2,…,k, the vectors(8)v1λ^ssym,v2λ^ssym,…,vn−qsλ^ssym,un−qs+1λ^ssym,un−qs+2λ^ssym,…,unλ^ssym,are linearly independent.Therefore, we can make the following reasonable assumption.Assumption 1.
For eachs=1,2,…,k, the vectors(9)v1λ˜ssym,v2λ˜ssym,…,vn−qsλ˜ssym,un−qs+1λ˜ssym,un−qs+2λ˜ssym,…,unλ˜ssym,are linearly independent.
Fors=1,2,…,k, define(10)Cssymεsym=Asym+Esym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs,where(11)εsym=ε1,1,ε1,2,…,ε1,n,ε2,1,…,ε2,n−1,…,εn−1,1,εn−1,2,εn,1T.Lemma 2.
For eachs=1,2,…,k if corankδAsym−λ˜ssymIn=qs, then the matrix(12)Asym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs,is nonsingular.Proof.
(13)Asym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs=Uλ˜ssymOqs,qsOqs,nIqsΣ1:n−qs,1:n−qsλ˜ssymOn−qs,qsOn−qs,qsOqs,n−qsΣn−qs+1:n,n−qs+1:nλ˜ssymIqsOqs,n−qsIqsOqs,qsV:,1:n−qsλ˜ssymTOn−qs,qsU:,n−qs+1:nλ˜ssymTOqs,qsΣn−qs+1:n,n−qs+1:nλ˜ssymIqs.
This completes the proof.For a persymmetric matrixApersym, let λ˜1persym,λ˜2persym,…,λ˜kpersym be the set of all distinct real eigenvalues of Apersym computed by eig code in MATLAB. Let J be the n×n inverse identity matrix; then, JApersym+Epersym−λ˜spersymIn is a symmetric matrix. For each s=1,2,…,k, the singular value decomposition of JApersym−λ˜spersymIn is(14)JApersym−λ˜spersymIn=Uλ˜spersymΣλ˜spersymVλ˜spersymT=u1λ˜spersym,…,unλ˜spersymdiagσ1λ˜spersym,…,σnλ˜spersymv1λ˜spersym,…,vnλ˜spersymT.Assume that, for eachs=1,2,…,k, corankδJApersym−λ˜spersymIn=qs, and then define(15)Cspersymεpersym=JApersym+Epersym−λ˜spersymInU:,n−qs+1:nλ˜spersymU:,n−qs+1:nλ˜spersymTOqs,qs,where(16)εpersym=εn,1,εn−1,1,…,ε1,1,εn−1,2,…,ε1,2,…,ε2,n−1,ε1,n−1,ε1,nT.For simplifying, we writeεstruct to stand for both εsym and εpersym and we write Csstruct to denote both Cssym and Cspersym. Then, define the constant K by(17)K=min1nn+qsCsstruct0−1∞,s=1,2,…,k.Lemma 3.
Ifεstruct∞<K, then, for each s=1,2,…,k, the matrix Csstructεstruct is nonsingular.Proof.
Ifεstruct∞<K, for each s=1,2,…,k,(18)In+1−Csstruct0−1Csstructεstruct2≤n+qsIn+qs−Csstruct0−1Csstructεstruct∞≤n+qsCsstruct0−1∞Csstruct0−Csstructεstruct∞<nn+qsεstruct∞Csstruct0−1∞<1.
From Lemma1, we can deduce that, for each s=1,2,…,k, the matrix Csstructεstruct is nonsingular.
For eachs=1,2,…,k, let Wsεstruct and Fsεstruct be, respectively, the first n-rows and the last qs-rows of the solution of the following linear system:(19)CsstructεstructWsεstructFsεstruct=On,qsIqs.
According to Cramer’s rule, we have the following easy lemma.Lemma 4.
For eachs=1,2,…,k, the matrix Fsεstruct is a symmetric matrix.Proposition 1.
Suppose thatε˜struct∞<K. If Fsε˜struct=Oqs,qs holds for s=1,2,…,k, then(20)corankAstruct+E˜struct−λ˜sstructIn=qs.Proof.
Firstly, consider the case for the symmetric matrix. Ifε˜sym∞<K, by Lemma 3, we know that, for each s=1,2,…,k, the matrix Cssymε˜sym is nonsingular. Next assume that s is an arbitrary fixed integer from 1 to k. If Fsε˜sym=Oqs,qs, then(21)Asym+Esym−λ˜ssymInWsε˜sym=On,qs,U:,n−qs+1:nλ˜ssymTWsε˜sym=Iqs.
Notice that the columns ofWsε˜sym are linearly independent. Hence, corankAsym+E˜sym−λ˜ssymIn≥qs. Assume that corankAsym+E˜sym−λ˜ssymIn>qs; then, there exists a nonzero vector βs∈nullAsym+E˜sym−λ˜ssymIn such that the matrix Wsε˜symβs is full rank, and then, for a nonzero vector bs∈ℝqs+1, the following equality(22)Asym+Esym−λ˜ssymInU:,n−qs+1:nλ˜ssymTWsε˜symβsbs=0,holds, which leads to a conditions. Consequently, corankAsym+E˜sym−λ˜ssymIn=qs.
Consider the case of persymmetric matrix. According to the above conclusion, we know thatcorankJApersym+E˜persym−λ˜spersymIn=qs. Apparently, corankApersym+E˜persym−λ˜spersymIn=qs.Lemma 5.
The following nonlinear system,(23)Gεstruct=F1εstruct1:qs,1F2εstruct2:qs,2⋮Fkεstructqs:qs,qs=0,is an underdetermined nonlinear system.Proof.
System (23) consists of nn+1/2 variables, and ∑i=1sqiqi+1/2 equations. Since(24)q12+q22+⋯+qs2+q1+q2+⋯+qs2≤q12+q22+⋯+qs2+n2=q1+q2+⋯+qs2−2∑1≤i≤j≤sqiqj+n2≤n2+n2−∑1≤i≤j≤sqiqj≤nn+12,then the proof is completed.
In the following, assume thats is an arbitrary fixed integer from 1 to k. By Lemma 3, there exists a neighborhood with the center of the origin such that, for an arbitrary εstruct belonging to this neighborhood, the matrix Csstructεstruct is nonsingular. Moreover, in this neighborhood, each entry of the solution of system (19) has the partial derivative with respect to each variable.
Concerning symmetric matrix, for each pairi,j with 1≤i≤n, 1≤j≤n−i+1, differentiating both sides of system (19) with respect to εi,j gives the following system:(25)Asym+Esym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs∂Wsε∂εi,j∂Fsε∂εi,j=t−Oj+i−2,qsWsεi,:On+qs−j−i+1,qs−Oi−1,qsWsεj+i−1,:On+qs−i,qs,where(26)t=12,j=1,1,2≤j≤n−i+1.
For pairsi1,j1 and i2,j2 with 1≤i1≤n, 1≤j1≤n−i1+1, 1≤i2≤n, and 1≤j2≤n−i2+1, differentiating both sides of system (25) with respect to εi1,j1 and εi2,j2, produces the following system:(27)Asym+Esym−λ˜ssymInU:,n−qs+1:nλ˜ssymU:,n−qs+1:nλ˜ssymTOqs,qs∂2Wsε∂εi1,j1∂εi2,j2∂2Fsε∂εi1,j1∂εi2,j2=t1−Oj1+i1−2,qs∂Wsε∂εi2,j2i1,:On+qs−j1−i1+1,qs−Oi1−1,qs∂Wsε∂εi2,j2j1+i1−1,:On+qs−i1,qs+t2−Oj2+i2−2,qs∂Wsε∂εi1,j1i2,:On+qs−j2−i2+1,qs−Oi2−1,qs∂Wsε∂εi1,j1j2+i2−1,:On+qs−i2,qs,where(28)t1=12,j1=1,1,2≤j1≤n−i1+1,t2=12,j2=1,1,2≤j2≤n−i2+1.
Considering the persymmetric matrix, for each pairi,j with 1≤i≤n and 1≤j≤n−i+1 differentiating both sides of system (19) with respect to εi,j yields the following system:(29)JApersym+Epersym−λ˜spersymInU:,n−qs+1:nλ˜spersymU:,n−qs+1:nλ˜spersymTOqs,qs∂Wsε∂εi,j∂Fsε∂εi,j=t−Oj−1,qsWsεn−i+1,:On+qs−j,qs−On−i,qsWsεj,:Oqs+i−1,qs,where(30)t=12,j=n−i+1,1,1≤j≤n−i.
For pairsi1,j1 and i2,j2 with 1≤i1≤n, 1≤j1≤n−i1+1, 1≤i2≤n, and 1≤j2≤n−i2+1, differentiating both sides of system (29) with respect to εi1,j1 and εi2,j2, provides the following system:(31)JApersym+Epersym−λ˜spersymInU:,n−qs+1:nλ˜spersymU:,n−qs+1:nλ˜spersymTOqs,qs∂2Wsε∂εi1,j1∂εi2,j2∂2Fsε∂εi1,j1∂εi2,j2=t1−Oj1−1,qs∂Wsε∂εi2,j2n−i1+1,:On+qs−j1,qs−On−i1,qs∂Wsε∂εi2,j2j1,:Oqs+i1−1,qs+t2−Oj2−1,qs∂Wsε∂εi1,j1n−i2+1,:On+qs−j2,qs−On−i2,qs∂Wsε∂εi1,j1j2,:Oqs+i2−1,qs,where(32)t1=12,j1=n−i1+1,1,1≤j1≤n−i1,t2=12,j2=n−i2+1,1,1≤j2≤n−i2.
Solving (25) and (27) or (29) and (31) at εstruct=0 can obtain Jacobi matrix G′0 and Hessian matrix G″0.Assumption 2.
The Jacobi matrixG′0 is full rank.
Suppose that, for an index setℐ=i1,j1,i2,j2,…,im,jm, the following matrix,(33)∂G0∂εi1,j1,∂G0∂εi2,j2,…,∂G0∂εim,jm,is nonsingular, and then define the following nonlinear system,(34)G˜εi1,j1,εi2,j2,…,εim,jm=0,by(35)G˜εi1,j1,εi2,j2,…,εim,jm=Gεstruct|εi,j=0fori,j∉ℐ.Remark 2.
By Theorem2, we compute verified error bound when the zero vector is the solution of the nonlinear system (34).
For system (34), define the constants B and η in Theorem 2 by(36)B=G˜′0−1∞,η=G˜′0−1G˜0∞.
Define anm-dimensional interval vector Ω=Ωi1,j1,Ωi2,j2,…,Ωim,jmT with each entry as interval −2η,2η. Define the corresponding interval perturbation matrix Ω˜ by setting Ω˜i,j=Ωi,j for i,j∈ℐ and Ω˜i,j=0 for i,j∉ℐ.
Byverifylss function, we solve interval linear systems (19), (25), (27) or (19), (29), (31) at E=Ω˜, and then we can obtain an interval tensor Hs,t:1≤s≤t≤m that satisfies(37)Hs,t⊃∂2GΩ∂εis,js∂εit,jt,1≤s≤t≤m.
Define the Lipschitz constantκ in Theorem 2 by(38)κ=max1≤s≤t≤mmaxH∞:∀H∈Hs,t.


## 3. Main Algorithm

We propose Algorithm1 for computing verified error bounds such that a slightly perturbed matrix is guaranteed to possess an eigenvalue of geometric multiplicity q within the computed bounds if the algorithm is successful.Theorem 3.
If Algorithm1 is successful, then there exists a perturbation matrix E˜struct of the same type, whose entries satisfy the condition(39)ε˜i,j≤ρ,i,j∈ℐ,ε˜i,j=0,i,j∉ℐ.λ˜1struct,λ˜2struct,…,λ˜kstruct are all exact real eigenvalues of matrix Astruct+E˜struct. Furthermore, for each s=1,2,…,k, qs is the geometric multiplicity of eigenvalue λ˜sstruct.Algorithm 1: VerifyEig.
InputAstruct: an n×n symmetric or persymmetric matrix; δ: the tolerance of numerical rank.Output {λ˜sstruct∈ℝ:s=1,2,…,k }, qs∈ℝ:s=1,2,…,k, ρ∈ℝ and an index set ℐ.Step 1 Use eig code to calculate all distinct real eigenvalues {λ˜sstruct∈ℝ:s=1,2,…,k} of Astruct.Step 2 For each s=1,2,…,k, compute the singular value decomposition of the matrix Asym−λ˜ssymIn or JApersym−λ˜spersymIn.Step 3 Solve (25) and (27) or (29) and (31) to obtain Jacobi matrix G′0.Step 4 If G′0 is full rank, choose the index set ℐ such that G˜′0 is nonsingular.Step 5 Compute the constants B,η by (36).Step 6 Compute the constant κ by (38).Step 7 If h=κBη≤1/2, then compute ρ by (5).Step 8 If ρ<K, return λ˜sstruct:s=1,2,…,k, qs:s=1,2,…,k, ρ and ℐ.Proof.
If Algorithm1 is successful, then by Theorem 2, we know that a perturbation matrix E˜struct whose entries satisfy condition (39), such that G˜ε˜i1,j1,ε˜i2,j2,…,ε˜im,jm=0. Thus, there exists an interval matrix E˜struct whose entries satisfy condition (39), such that Gε˜struct=0. Finally, it follows by Proposition 1 that, for each s=1,2,…,k, corankAstruct+E˜struct−λ˜sstructIn=qs.


## 4. Examples

In this section, we show the performance of Algorithm1. The following experiments are carried out in Matlab R2012a with INTLAB V5 under Windows 7.Example 1.
Given a symmetric matrix(40)A=1001000011000001100001001000000010200000200000201,then by Algorithm 1, we get q1=q4=1,q2=2,q3=3, λ˜1=−1,λ˜2=0,λ˜3=2,λ˜4=3, and ρ=0.Example 2.
For a symmetric matrix,(41)A=1101000011100000011000001001000000001002000001200000021000002001,applying Algorithm 1 yields q1=q6=2,q2=q3=q4=q5=1, λ˜1=−1.0000,λ˜2=−0.6180,λ˜3=0.3820,λ˜4=1.6180,λ˜5=2.6180,λ˜6=3.0000, and ρ=5.5359e−16.Example 3.
Given a symmetric matrix(42)A=100110000011100000011000000110100000100010000000001001000000110000000110000001001,then the computed results of Algorithm 1 are q1=q3=q5=1,q2=q4=3, λ˜1=−0.7321,λ˜2=0,λ˜3=1.0000,λ˜4=2.0000,λ˜5=2.7321, and ρ=5.1876e−16.Example 4.
For ann×n symmetric matrix A, whose entries are uniformly distributed in the interval 0,1, Table 1 shows the radius ρ computed by Algorithm 1 for different n.Table 1
The calculation results ofρ in Example 4.
nρ41.0257e − 1553.1887e − 1465.7511e − 1574.4464e − 1581.3765e − 1591.9717e − 14102.1013e − 15201.2099e − 14Example 5.
Given a persymmetric matrix(43)A=1000101020000000100010001000000010010000200100001,Algorithm 1 outputs q1=1,q2=3, λ˜1=2,λ˜2=1, and ρ=0.Example 6.
Given a persymmetric matrix(44)A=1000000202100000000020000001020000001000101000100000002010100001,then by Algorithm 1, we obtain q1=q5=2,q2=q3=q4=1, λ˜1=1.0000,λ˜2=0,λ˜3=2.4142,λ˜4=−0.4142,λ˜5=2.0000, and ρ=1.2527e−15.Example 7.
Given a persymmetric matrix(45)A=100000012020000001000002000000100200000000000000001000200000000000000020102000001,then by Algorithm 1, we obtain q1=q2=q4=q5=1,q3=2, λ˜1=2.0000,λ˜2=1.0000,λ˜3=0,λ˜4=2.4142,λ˜5=−0.4142, and ρ=1.9386e−15.


## 5. Conclusion and Future Work

This paper demonstrated how to compute the validated and narrow error bounds for the symmetric and persymmetric matrices, such that there exists a perturbation matrix of the same type within the computed error bound whose exact real eigenvalues are the computed real eigenvalues. We will make further efforts to extend the results for other types of structured matrices.



---

*Source: 9240852-2020-12-21.xml*