Gershgorin and Rayleigh Bounds on the Eigenvalues of the Finite-Element Global Matrices via Optimal Similarity Transformations

Isaac Fried, Roberto Riganti, Chen Yu

Department of Mathematics, Boston University, Boston, MA, USA.

DOI: 10.4236/am.2020.119060                

Abstract

The large finite element global stiffness matrix is an algebraic, discreet, even-order, differential operator of zero row sums. Direct application of the, practically convenient, readily applied, Gershgorin’s eigenvalue bounding theorem to this matrix inherently fails to foresee its positive definiteness, predictably, and routinely failing to produce a nontrivial lower bound on the least eigenvalue of this, theoretically assured to be positive definite, matrix. Considered here are practical methods for producing an optimal similarity transformation for the finite-elements global stiffness matrix, following which non trivial, realistic, lower bounds on the least eigenvalue can be located, then further improved. The technique is restricted here to the common case of a global stiffness matrix having only non-positive off-diagonal entries. For such a matrix application of the Gershgorin bounding method may be carried out by a mere matrix vector multiplication.

Keywords

Finite Elements, Global Stiffness Matrix, Gershgorin and Rayleigh Computed Upper and Lower Bounds on the Extremal Eigenvalues, Similarity Transformations

1. Introduction

Knowledge, even approximate, of the extremal eigenvalues of the, large, positive definite, global stiffness matrix generated by the finite-element method (see [1] ), is essential for assessing the correctness of the numerical solution of the global algebraic system of equations set up by this variational method. Other numerical procedures related to this large algebraic system, such as iterative solution methods, are also greatly affected by the condition number, the ratio of the largest to the smallest eigenvalue, of the matrix (see [2] ). In this note we, theoretically and numerically, examine various practical procedures for the numerical, a-posteriori, ever better, bounding of the eigenvalues of the large global stiffness matrices generated by the finite-element method, predominantly based on a preparatory, most favorable, similarity transformation. Special attention is paid in this note to the common class of matrices having non-positive off-diagonal entries (see [3] [4] [5] ). For such matrices, application of the Gershgorin bounding method may be most conveniently carried out by a mere matrix vector multiplication, which may be, in turn, very effectively carried out on the element level (see [6] ).

2. Gershgorin’s Eigenvalue Bounds

The utterly simple Gershgorin and Rayleigh eigenvalue bound theorems (see [7] [8] ) are of such fundamental importance in computational linear algebra that we find it good to fully repeat them here.

Throughout this paper we denote a number by a lower case Greek letter, a vector by a lower case Roman letter, and a matrix by an upper case Roman letter.

The eigenvalue spectrum $\lambda \left(A\right)$ of real symmetric matrix A is real. For such a matrix Gershgorin’s theorem assumes the simpler form.

Theorem (Gershgorin). Let $A=A\left(n\times n\right)=A^{\text{T}}$ be symmetric, and so of a real eigensystem. Then every eigenvalue of A lies in, or on the end points of, at least one of the intervals

$\left|\lambda -A_{ii}\right|\le \left|A_{i1}\right|+\left|A_{i2}\right|+\cdots +\left|A_{in}\right|,i=\mathrm{1,2,}\cdots \mathrm{,}n\mathrm{.}$ (1)

Namely,

$\begin{array}{l}\underset{i}{min}\left[A_{ii}-\left(\left|A_{i1}\right|+\left|A_{i2}\right|+\cdots +\left|A_{in}\right|\right)\right]\le \lambda \left(A\right)\\ \le \underset{i}{max}\left[A_{ii}+\left(\left|A_{i1}\right|+\left|A_{i2}\right|+\cdots +\left|A_{in}\right|\right)\right]\end{array}$ (2)

for every eigenvalue of A.

Proof. Let λ be an eigenvalue of A and $x=\left\{x_1{x}_2\cdots x_n\right\}\ne 0$ the corresponding eigenvector, so that $Ax=\lambda x$. Assume that for this λ the kth component of $x\mathrm{,}{x}_k$, is largest in magnitude and is normalized so that $\left|x_k\right|=1$ and $\left|x_i\right|\le 1$, $i\ne k$. The kth equation of $Ax=\lambda x$ is then

$A_{k1}{x}_1+A_{k2}{x}_2+\cdots +A_{kk}+\cdots +A_{kn}{x}_n=\lambda$ (3)

and

$\begin{array}{c}\left|\lambda -A_{kk}\right|=\left|A_{k1}{x}_1+A_{k2}{x}_2+\cdots +A_{kn}{x}_n\right|\\ \le \left|A_{k1}\right|\left|x_1\right|+\left|A_{k2}\right|\left|x_2\right|+\cdots +\left|A_{kn}\right|\left|x_n\right|\\ \le \left|A_{k1}\right|+\left|A_{k2}\right|+\cdots +\left|A_{kn}\right|\end{array}$ (4)

since $\left|x_i\right|<1$. We do not know what k is, but are sure that λ certainly lies in one of these intervals. □

Gershgorin’s theorem readily provides outside bounds on all eigenvalues of matrix A.

3. Positive-Definite Matrices with Non-Positive Off-Diagonal Entries

Positive definite and symmetric global stiffness matrices of a positive diagonal and non-negative off-diagonal entries are common in the finite elements, or finite differences, method. For such matrices the application of the Gershgorin theorem simplifies into a mere matrix-vector multiplication operation, efficiently carried out as a summation of such multiplications on the element level of the global mesh (see [6] ).

Say matrix A is symmetric, $A=A^{\text{T}}$ and positive (semi) definite, and of non-positive off-diagonal entries, then obviously

${\lambda }_1\left(A\right)\ge \underset{i}{\min }Ae,e_i=1\text{or}e=\left\{1,1,1,\cdots ,1\right\}$ (5)

where ${\lambda }_1\left(A\right)$ is the smallest eigenvalue of matrix A, and ${\lambda }_n\left(A\right)$ the largest. If matrix A is symmetric and positive definite, then ${\lambda }_1\left(A\right)>0$. For example

$Ae=\left[\begin{array}{ccc}1& -1& \\ -1& 2& -1\\ r& -1& 2\end{array}\right]\left[\begin{array}{c}1\\ 1\\ 1\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ 1\end{array}\right],\underset{i}{\min }\left[\begin{array}{c}0\\ 0\\ 1\end{array}\right]=0,0\le {\lambda }_1\left(A\right).$ (6)

Gershgorin’s theorem correctly predicts here that matrix A is, at least, positive semi-definite, but it fails to ascertain that it is actually positive definite, with a positive lowest eigenvalue. To obtain a more realistic, nontrivial, lower bound on the lowest eigenvalue of matrix A we resort to similarity transformations designed to rescue the lower bound from triviality.

4. Rayleigh Quotient

Let matrix $A=A^{\text{T}}$ be real and symmetric. If $x$ is an eigenvector of A for corresponding eigenvalue $\lambda$, $Ax=\lambda x$, then $x^{\text{T}}Ax=\lambda x^{\text{T}}x$, or $\lambda =x^{\text{T}}Ax/x^{\text{T}}x$. The rational function, for variable vector $x$

$R\left[x\right]=\frac{x^{\text{T}}Ax}{x^{\text{T}}x}\mathrm{,}x\ne 0$ (7)

is the Rayleigh quotient of matrix A. It has some interesting properties of great practical and theoretical reach.

First property: The quotient $R\left[x\right]$ produces very accurate eigenvalue approximations from even not so good eigenvector approximations. Indeed, let v be an eigenvector corresponding to some $\lambda$, and let

$x=v+ϵw,v^{\text{T}}v=1,w^{\text{T}}w=1,w^{\text{T}}v=v^{\text{T}}w=0$ (8)

where $ϵw$ is the error vector.

Then

$R\left[x\right]=\frac{x^{\text{T}}Ax}{x^{\text{T}}x}=\frac{{\left(v+ϵw\right)}^{\text{T}}A\left(v+ϵw\right)}{{\left(v+ϵw\right)}^{\text{T}}\left(v+ϵw\right)}=\frac{v^{\text{T}}Av+2ϵw^{\text{T}}Av+{ϵ}^2{w}^{\text{T}}Aw}{v^{\text{T}}v+2ϵw^{\text{T}}v+{ϵ}^2{w}^{\text{T}}w}$ (9)

or

$R\left[x\right]=\frac{\lambda +{ϵ}^2{w}^{\text{T}}Aw}{1+{ϵ}^2}=\lambda +{ϵ}^2\frac{w^{\text{T}}Aw-\lambda }{1+{ϵ}^2}$ (10)

since $w^{\text{T}}v=0$, implying that also $w^{\text{T}}Av=0$. Thus, here

$R\left[x\right]=\lambda +O\left({ϵ}^2\right)\mathrm{.}$ (11)

Rayleigh quotient $R\left[x\right]$ produces a very accurate, $O\left({ϵ}^2\right)$, approximation to an eigenvalue from an approximation x to the corresponding eigenvalue that is not very accurate, only $O\left(ϵ\right)$. Conversely, a very good approximation to an eigenvalue does not imply it coming from a very good approximation to an eigenvector.

Second property: If $A=A^{\text{T}}$, then for any vector $x\ne 0$

${\lambda }_1\left(A\right)\le R\left[x\right]=\frac{x^{\text{T}}Ax}{x^{\text{T}}x}\le {\lambda }_n\left(A\right)\text{or}{\lambda }_1\left(A\right)\le x^{\text{T}}Ax\le {\lambda }_n\left(A\right)\mathrm{,}{x}^{\text{T}}x=1.$ (12)

The Rayleigh quotient provides inner bounds on all eigenvalues of symmetric matrix A. Hence the Rayleigh bounds and the Gershgorin bounds complement each other.

Theorem (Rayleigh.) Let the eigenvalues of $A=A^{\text{T}}$ be arranged in the ascending order ${\lambda }_1\le {\lambda }_2\le \cdots \le {\lambda }_n$ , with orthogonal eigenvectors $x_1\mathrm{,}{x}_2\mathrm{,}\cdots \mathrm{,}{x}_n$ . Then

${\lambda }_{k+1}\le \frac{x^{\text{T}}Ax}{x^{\text{T}}x}\le {\lambda }_n\text{if}{x}^{\text{T}}{x}_1=x^{\text{T}}{x}_2=\cdots =x^{\text{T}}{x}_k=0,x\ne 0$ (13)

with the lower equality holding if and only if $x=x_{k+1}$, and the upper inequality holding if and only if $x=x_n$. Also

${\lambda }_1\le \frac{x^{\text{T}}Ax}{x^{\text{T}}x}\le {\lambda }_{n-k}\text{if}{x}^{\text{T}}{x}_n=x^{\text{T}}{x}_{n-1}=\cdots =x^{\text{T}}{x}_{n-k+1}=0,x\ne 0$ (14)

with the lower equality holding if and only if $x=x_1$, and the upper if and only if $x=x_{n-k}$. The two inequalities reduce to

${\lambda }_1\le \frac{x^{\text{T}}Ax}{x^{\text{T}}x}\le {\lambda }_n$ (15)

for arbitrary, nonzero, $x\in R^n$.

Proof. Vector $x\in R^n$, orthogonal to $x_1\mathrm{,}{x}_2\mathrm{,}\cdots \mathrm{,}{x}_k$ has the unique expansion

$x={\alpha }_{k+1}{x}_{k+1}+{\alpha }_{k+2}{x}_{k+2}+\cdots +{\alpha }_n{x}_n$ (16)

with which

$x^{\text{T}}Ax={\lambda }_{k+1}{\alpha }_{k+1}^2+{\lambda }_{k+2}{\alpha }_{k+2}^2+\cdots +{\lambda }_n{\alpha }_n^2.$ (17)

We normalize x by

$x^{\text{T}}x={\alpha }_{k+1}^2+{\alpha }_{k+2}^2+\cdots +{\alpha }_n^2=1$ (18)

and use this equation to eliminate ${\alpha }_{k+1}^2$ from $x^{\text{T}}Ax$ to be left with

$\lambda \left(x\right)=x^{\text{T}}Ax={\lambda }_{k+1}+{\alpha }_{k+2}^2\left({\lambda }_{k+2}-{\lambda }_{k+1}\right)+\cdots +{\alpha }_n^2\left({\lambda }_n-{\lambda }_{k+1}\right).$ (19)

By assumption ${\lambda }_j-{\lambda }_{k+1}\ge 0$ if $j>k+1$ and hence

$\lambda \left(x\right)={\lambda }_{k+1}+non-negativequantity$ (20)

or $\lambda \left(x\right)\ge {\lambda }_{k+1}$, with equality holding if and only if

${\alpha }_{k+2}^2\left({\lambda }_{k+2}-{\lambda }_{k+1}\right)+\cdots +{\alpha }_n^2\left({\lambda }_n-{\lambda }_{k+1}\right)=0.$ (21)

In case of distinct eigenvalues, ${\lambda }_j-{\lambda }_{k+1}\ne 0,j=k+2,\cdots ,n$, equality holds if and only if ${\alpha }_{k+2}={\alpha }_{k+3}=\cdots ={\alpha }_n=0$, and $\lambda \left(x_{k+1}\right)={\lambda }_{k+1}$. If eigenvalues repeat and ${\lambda }_j-{\lambda }_{k+1}=0$, then ${\alpha }_j$ need not be zero, but equality still holds if and only if x is in the invariant subspace spanned by the eigenvectors of ${\lambda }_{k+1}$.

To prove the upper bound we use

${\alpha }_n^2=1-{\alpha }_{k+1}^2-{\alpha }_{k+2}^2-\cdots -{\alpha }_{n-1}^2$ (22)

to eliminate it from $\lambda \left(x\right)$, so as to be left with

$\lambda \left(x\right)={\lambda }_n-{\alpha }_{k+1}^2\left({\lambda }_n-{\lambda }_{k+1}\right)-\cdots -{\alpha }_{n-1}^2\left({\lambda }_n-{\lambda }_{n-1}\right)$ (23)

and $\lambda \left(x\right)\le {\lambda }_n$ with equality holding if and only if $x=x_n$. The proof to the second part of the theorem is the same. □

5. Perron’s Theorem on Positive Matrices

Positive matrices are common in the finite-element method. The following is a symmetric version of Perron’s theorem on positive matrices (see [9] ).

Theorem (Perron). If A is a symmetric positive matrix, $A_{ij}>0$ , then the eigenvector corresponding to the largest (positive) eigenvalue of A is positive and unique.

Proof. If $x_n$ is a unit eigenvector corresponding to ${\lambda }_n$, and $x\ne x_n$ is such that $x^{\text{T}}x=1$, then

$x^{\text{T}}Ax<{\lambda }_n=\lambda \left(x_n\right)=x_n^{\text{T}}A{x}_n$ (24)

and ${\lambda }_n$ is certainly positive. Moreover, since $A_{ij}>0$ the components of $x_n$ cannot have different signs, for this would contradict the assumption that $x_n$ maximizes $\lambda \left(x\right)$. Say then that ${\left(x_n\right)}_i\ge 0$. But none of the ${\left(x_n\right)}_i$ components can be zero since $A{x}_n={\lambda }_n{x}_n$, and obviously $A{x}_n>0$. Hence $x_n>0$.

There can be no other positive vector orthogonal to $x_n$, and so the eigenvector, and also the largest eigenvalue ${\lambda }_n$, are unique. □

The following is another important and useful statement on positive matrices.

Theorem. Suppose that A has a positive inverse $A^{-1}$ . Let x be any vector satisfying $Ax-e=r$ , $e=\left\{1,1,\cdots ,1\right\}$, ${\Vert r\Vert }_{\infty }<1$ . Then

$\frac{{\Vert x\Vert }_{\infty }}{1+{\Vert r\Vert }_{\infty }}\le {\Vert A^{-1}\Vert }_{\infty }\le \frac{{\Vert x\Vert }_{\infty }}{1-{\Vert r\Vert }_{\infty }}\mathrm{.}$ (25)

Proof. Obviously $x=A^{-1}e+A^{-1}r$ so that

${\Vert x\Vert }_{\infty }\le {\Vert A^{-1}e\Vert }_{\infty }+{\Vert A^{-1}r\Vert }_{\infty }\le {\Vert A^{-1}\Vert }_{\infty }+{\Vert A^{-1}\Vert }_{\infty }{\Vert r\Vert }_{\infty }$ (26)

and

$\frac{{\Vert x\Vert }_{\infty }}{1+{\Vert r\Vert }_{\infty }}\le {\Vert A^{-1}\Vert }_{\infty }\mathrm{.}$ (27)

To prove the other bound write $x=A^{-1}e-\left(-A^{-1}r\right)$, observe that ${\Vert A^{-1}e\Vert }_{\infty }={\Vert A^{-1}\Vert }_{\infty }$, and have that

${\Vert x\Vert }_{\infty }\ge {\Vert A^{-1}e\Vert }_{\infty }-{\Vert A^{-1}r\Vert }_{\infty }\ge {\Vert A^{-1}\Vert }_{\infty }-{\Vert A^{-1}\Vert }_{\infty }{\Vert r\Vert }_{\infty }\mathrm{.}$ (28)

Hence, if ${\Vert r\Vert }_{\infty }<1$, then

$\frac{{\Vert x\Vert }_{\infty }}{1-{\Vert r\Vert }_{\infty }}\ge {\Vert A^{-1}\Vert }_{\infty }\mathrm{.}$ (29)

This completes the proof.

6. Similarity Transformation, Similar Matrices

We will make also great use of the following fact.

Theorem. If $\lambda$ is an eigenvalue of matrix A, then $\lambda$ is also an eigenvalue of similar matrix $P^{-1}AP$ . Similar matrices A and $P^{-1}AP$ have the same eigenvalues.

Proof. We have

$\begin{array}{l}Ax=\lambda x,P^{-1}Ax=\lambda P^{-1}x,x=Py,\\ P^{-1}APy=\lambda P^{-1}Py,P^{-1}APy=\lambda y,\end{array}$ (30)

and the result follows. □

In this note we are greatly interested in accurate eigenvectors.

To fix ideas consider the typical finite-elements, or finite-differences, matrix

$A=\left[\begin{array}{ccc}1& -1& \\ -1& 2& -1\\ & -1& 2\end{array}\right],A^{-1}=\left[\begin{array}{ccc}3& 2& 1\\ 2& 2& 1\\ 1& 1& 1\end{array}\right].$ (31)

The matrix A is symmetric and positive definite, ${\lambda }_1\left(A\right)>0$, with a corresponding eigenvector that is positive. The fact that $A^{-1}$ is completely positive is a manifestation of the physical fact that a point force applied to the system causes all free points of the system to move, and all in the same direction.

To save Gershgorin’s theorem from the trivial, and useless prediction ${\lambda }_1\left(A\right)\ge 0$, we propose to similarly transform matrix A with a positive diagonal matrix D, $D_{ii}>0$, such that

$\underset{i}{\min }{D}^{-1}ADe=\underset{i}{\min }{D}^{-1}AD\left[\begin{array}{c}1\\ 1\\ 1\end{array}\right]\le {\lambda }_1\left(A\right),e=\left[\begin{array}{c}1\\ 1\\ 1\end{array}\right].$ (32)

The fact that diagonal matrix D is positive, $D_{ii}>0$, guarantees that similar matrix $D^{-1}AD$ retains the entry sign pattern of original matrix A. Namely, ${\left(D^{-1}AD\right)}_{ii}>0$, ${\left(D^{-1}AD\right)}_{ij}\le 0$.

For example, with

$D=\left[\begin{array}{ccc}6& & \\ & 5& \\ & & 3\end{array}\right],D^{-1}ADe=\left[\begin{array}{c}1/6\\ 1/5\\ 1/3\end{array}\right],\text{and}{\lambda }_1\left(A\right)\ge \frac{1}{6}=\mathrm{0.167.}$ (33)

Now Gershgorin’s theorem correctly recognizes that matrix A is positive definite.

For a better bound, matrix D needs to raise the lowest entry and to lower the highest entry of $D^{-1}ADe$. Thus,

The optimal D is for which all entries of $D^{-1}ADe$ are equal, or

$D^{-1}ADe=\lambda e,A\left(De\right)=\lambda \left(De\right).$ (34)

Namely, the optimal D is such that De is the eigenvector corresponding to the lowest eigenvalue of matrix A. With such a better D, obtained here by inverse iterations on, the here positive, $A^{-1}$, $A_{ij}^{-1}>0$, we have

$\begin{array}{l}D=\left[\begin{array}{ccc}793& & \\ & 636& \\ & & 353\end{array}\right],D^{-1}ADe=\left[\begin{array}{c}157/793\\ 21/106\\ 70/353\end{array}\right],\text{and}\\ {\lambda }_1\left(A\right)\ge \frac{157}{793}=\mathrm{0.197982.}\end{array}$ (35)

Taking the diagonal entries of matrix D as vector x we obtain from Rayleigh’s quotient

$R\left[x\right]=\frac{x^{\text{T}}Ax}{x^{\text{T}}x}=0.198062\ge {\lambda }_1\left(A\right)$ (36)

or

$0.197982\le {\lambda }_1\left(A\right)\le 0.198062$ (37)

which are good enclosing bounds.

As for the highest eigenvalue of matrix A, ${\lambda }_3\left(A\right)$, we have from Gershgorin and Rayleigh that

$R\left[x\right]=\frac{x^{\text{T}}Ax}{x^{\text{T}}x}=3.24673\le {\lambda }_3\left(A\right)\le 4,x=\left\{19,-42,33\right\}$ (38)

with x obtained by the application of three steps of the power method to matrix A starting with $x_0=\left\{\mathrm{1,}-\mathrm{1,1}\right\}$ to have a reasonable approximate for the eigenvector corresponding to the highest eigenvalue of A.

Next we turn our attention to more realistic finite element examples.

7. Taut String (Cable)

The basic issues of this paper are best illustrated on some concrete cases. First we consider the one-dimensional problem of the taut string discretized by n linear finite elements of which the element stiffness matrix is

$k=\frac{1}{h}\left[\begin{array}{cc}1& -1\\ -1& 1\end{array}\right],{\lambda }_1\left(k\right)=0,{\lambda }_2\left(k\right)=\frac{2}{h}$ (39)

in which $h=1/n$ is the size of the element, n being the number of mesh subdivisions. Little element matrix k is positive semi-definite with a zero eigenvalue for eigenvector $x=\left\{\mathrm{1,1}\right\}$. Element stiffness matrix k is, moreover, such that $k_{ii}>0$, $k_{ij}\le 0$, and of this sign pattern is also any global matrix A assembled from it.

We assemble the element k matrices into the global matrix A, here over a mesh of five elements, impose the essential boundary condition of a fixed end point by deleting the first row and first column of the assembled global matrix, and are left with

$A=\left[\begin{array}{ccccc}2& -1& & & \\ -1& 2& -1& & \\ & -1& 2& -1& \\ & & -1& 2& -1\\ & & & -1& 1\end{array}\right],h=1,$

${\lambda }_1\left(A\right)=0.081014053,\text{and}{\lambda }_5\left(A\right)=\mathrm{3.6825.}$ (40)

Global stiffness matrix A is now positive definite, $\lambda \left(A\right)>0$. Matrix A is, moreover, tridiagonal with non-positive off-diagonal entries $A_{ij}\le 0$.

We know (see [10] [11] ) that as the size of matrix $A\left(n\times n\right)$ increases, its lowest eigenvalue decreases, actually $O\left(n^{-2}\right)$.

8. Taut String, Similarity Transformation Followed by a Gershgorin Bound

Being keenly interested in an actual, numerical lower bound on ${\lambda }_1\left(A\right)$ we naively apply Gershgorin’s theorem directly to global stiffness matrix A of Equation (40) and obtain from it the disappointing

$0\le \lambda \left(A\right)\le 4.$ (41)

The bounding fails to confirm that matrix A is positive definite, only that it is positive semi-definite. We remark that the lower Gershgorin bound is obtained here as

$\underset{i}{\min }{\left(Ae\right)}_i,e=\left\{1,1,1,\cdots ,1\right\}$ (42)

inasmuch as $A_{ii}>0$ and $A_{ij}\le 0$.

To avert a trivial Gershgorin bound we propose to similarly transform matrix A into P⁻¹AP with a positive diagonal matrix P to facilitate the inversion, and subsequent matrix-matrix multiplications, while still maintaining ${\left(P^{-1}AP\right)}_{ij}\le 0$.

For example, we take $P=\text{diag}\left[\mathrm{10,19,27,33,37}\right]$ and have

$\left(P^{-1}AP\right)e=\left\{\mathrm{0.100,0.053,0.074,0.061,0.108}\right\}\mathrm{,}$

$\underset{i}{\min }{\left(P^{-1}APe\right)}_i=\mathrm{0.053,}$ (43)

and

$0.053\le \lambda \left(A\right)\le 4$ (44)

triumphantly ascertaining that global stiffness matrix A is indeed positive definite.

We have already remarked that we desire matrix P to be such that all entries of $\left(P^{-1}AP\right)e$ are nearly equal, or that

$\left(P^{-1}AP\right)e=\lambda e,A\left(Pe\right)=\lambda \left(Pe\right)$ (45)

or that Pe is the (positive) eigenvector corresponding to ${\lambda }_1\left(A\right)$.

How to efficiently get good approximations to the first eigenvector of global stiffness matrix A is the subject of the rest of this note.

9. Taut String, Linearization of the Characteristic Polynomial

We know that as the size of matrix $A\left(n\times n\right)$ increases, its lowest eigenvalue decreases, tending to zero, actually $O\left(n^{-2}\right)$ (see again [10] [11] ). We start the linearization by writing

$Ax=\lambda x,x=\left\{1,\alpha ,\beta ,\gamma ,\delta \right\}$ (46)

and have by a successive linearization for $\lambda$ of each equation that, approximately

$1-15\lambda =0,\lambda =\frac{1}{15}=\mathrm{0.06667.}$ (47)

Then, linearly

$\alpha =2-\lambda ,\beta =3-4\lambda ,\gamma =4-10\lambda ,\delta =5-20\lambda ,$ (48)

and it follows that the eigenvector corresponding to the lowest eigenvalue of matrix A is approximately

$x=\left\{\mathrm{15,29,41,50,55}\right\}\mathrm{,}$ (49)

and, correspondingly

$R=R\left[x\right]=\frac{x^{\text{T}}Ax}{x^{\text{T}}x}=\mathrm{0.081117.}$ (50)

We put $P=\text{diagonal}\left[x\right]$ and have from

$G=G\left(x\right)=\underset{i}{\min }{P}^{-1}APe=0.06667$ (51)

the reasonable bounds

$G\le {\lambda }_1\left(A\right)\le R\text{or}0.06667\le {\lambda }_1\left(A\right)\le \mathrm{0.081117.}$ (52)

10. Taut String, Quadratization of the Characteristic Polynomial

Keeping quadratic terms of $\lambda$ in the characteristic equation of matrix $A=A\left(5\times 5\right)$ we have for vector x of Equation (46) the five quadratically curtailed equations of $Ax=\lambda x$ :

$\begin{array}{l}\alpha =2-\lambda \\ \beta =3-4\lambda +{\lambda }^2\\ \gamma =4-10\lambda +6{\lambda }^2\\ \delta =5-20+21{\lambda }^2\\ 35{\lambda }^2-15\lambda +1=0,\end{array}$ (53)

and

$\lambda =0.082578,x=\left\{1,1.917,2.677,3.215,3.492\right\}$ (54)

resulting in

$R=R\left[x\right]=0.0810117\text{and}G=\underset{i}{\min }{P}^{-1}APe=0.07919,$ (55)

and

$G\le {\lambda }_1\left(A\right)\le R\text{or}0.0792\le {\lambda }_1\left(A\right)\le \mathrm{0.081017.}$ (56)

The success of the linearization, or quadratization, of the characteristic equation hinges on the theoretically assured fact (see [10] [11] ) that ${\lambda }_1\left(A\right)\ll 1$, and that as the size of the finite-elements global stiffness matrix A increases ${\lambda }_1\left(A\right)$ further decreases.

11. Taut String, Linearization Following a Shift

To work with a matrix of a lesser minimal eigenvalue we turn to the shifted matrix

$B=A-RI,R=0.081117,{\lambda }_i\left(B\right)={\lambda }_i\left(A\right)-R$ (57)

from which we get the approximate linear characteristic polynomial

$9.85553\lambda +0.00102=0,$ (58)

then

$\lambda =-0.000102999,x=\left\{1,1.919,2.683,3.229,3.513\right\},$

$R=R\left[x\right]=0.081014052771,$ (59)

and

$G=\underset{i}{\min }{P}^{-1}APe=0.08101402206921$ (60)

with

$R-G=3.1\times {10}^{-8}.$ (61)

12. Taut String, Shifted Power Method

We continue with our quest for good approximations to the fundamental eigenvector of the finite elements global stiffness matrix A. We propose to first iterate for the highest eigenvalue and the corresponding eigenvector. We choose to use the power method since it requires but one vector-matrix multiplication per step. We start with

$x_0=\left\{1,-1,1,-1,1\right\},\text{to}\text{have}{x}_6=A{x}_5,R\left[x_6\right]=\mathrm{3.68242.}$ (62)

Then we undertake the shift

$B=R_{6}I-A,$ (63)

and apply to it the power method starting now with

$x_0=\left\{1,4,6,8,9\right\},x_1=A{x}_0$ (64)

to successively have

$\begin{array}{l}{R}_5=0.0821133,G_5=-0.00985\\ R_6=0.0817712,G_6=0.0108047\\ R_7=0.0815362,G_7=0.0255912\\ R_8=0.0813744,G_8=0.036619986\end{array}$ (65)

as upper and lower bounds on ${\lambda }_1\left(A\right)$.

13. Hanging String

The freely down hanging string (see [12] ) is of zero tension at the free lower end, and with the tension growing linearly towards the fixed upper hanging point, that carries the entire weight of the string.

The element stiffness matrix of the hanging string of size h is

$k_i=\frac{1}{h}\frac{2i-1}{2}\left[\begin{array}{cc}1& -1\\ -1& 1\end{array}\right]$ (66)

assembled into the global stiffness matrix

$A=\left[\begin{array}{ccccc}1& -1& & & \\ -1& 4& -3& & \\ & -3& 8& -5& \\ & & -5& 12& -7\\ & & & -7& 16\end{array}\right],A^{-1}=\frac{1}{315}\left[\begin{array}{ccccc}563& 248& 143& 80& 35\\ 248& 248& 143& 80& 35\\ 143& 143& 143& 80& 35\\ 80& 80& 80& 80& 35\\ 35& 35& 35& 35& 35\end{array}\right]$ (67)

including the essential boundary condition of zero movement of the hanging point. Global stiffness matrix A factors as

$2A=\left[\begin{array}{ccccc}1& & & & \\ -1& 1& & & \\ & -1& 1& & \\ & & -1& 1& \\ & & & -1& 1\end{array}\right]\left[\begin{array}{ccccc}1& & & & \\ & 3& & & \\ & & 5& & \\ & & & 7& \\ & & & & 9\end{array}\right]\left[\begin{array}{ccccc}1& -1& & & \\ & 1& -1& & \\ & & 1& -1& \\ & & & 1& -1\\ & & & & 1\end{array}\right]$ (68)

showing matrix A to be symmetric and positive definite. In fact, ${\lambda }_1\left(A\right)=0.398987923256$ and ${\lambda }_5\left(A\right)=22.007479$.

We are seeking good lower and upper bounds on ${\lambda }_1\left(A\right)>0$, the lowest eigenvalue of global stiffness A.

14. Hanging String, Linearization of the Characteristic Polynomial

Here we have

$x=\left\{1,\alpha ,\beta ,\gamma ,\delta \right\},\alpha =1-\lambda ,\beta =1-\frac{5}{3}\lambda ,\gamma =1-\frac{34}{15}\lambda ,\delta =1-\frac{298}{105}\lambda$ (69)

with the linearized characteristic polynomial for $\lambda$

$1069\lambda -315=0$ (70)

then

$\lambda =0.294668,x=\left\{1,0.705,0.509,0.332,0.164\right\},R=R\left[x\right]=0.42172494$ (71)

from which we get

$G=G\left(x\right)=\underset{i}{\min }{P}^{-1}APe=0.29467$ (72)

and then the reasonable bounds

$G\le {\lambda }_1\left(A\right)\le R,0.2947\le {\lambda }_1\left(A\right)\le 0.4217,{\lambda }_1\left(A\right)=\mathrm{0.398987923256.}$ (73)

15. Hanging String, Second Sweep for a Linearized Characteristic Polynomial

Now we undertake the shift

$B=IR-A,R=0.4217249428452726$ (74)

intended to further reduce ${\lambda }_1\left(B\right)$, and have for

$x=\left\{1,\alpha ,\beta ,\gamma ,\delta \right\},Bx=\lambda x$ (75)

a linearization resulting in

$\begin{array}{l}\alpha =0.578275+\lambda \\ \beta =0.356409+1.38552\lambda \\ \gamma =0.193228+1.57125\lambda \\ \delta =0.0650291+1.63685\lambda \end{array}$ (76)

and then the linearized characteristic polynomial

$0.339557-14.5657\lambda =0,\lambda =\mathrm{0.023312131.}$ (77)

From this we have

$x=\left\{1,0.601587,0.388708,0.229857,0.103188\right\}$ (78)

by which we compute for original matrix A

$R\left[x\right]=0.398988578,G\left[x\right]=0.3984128$ (79)

and $G\le {\lambda }_1\left(A\right)\le R$, translating into the good bounds

$0.39841\le {\lambda }_1\left(A\right)\le 0.3989886,{\lambda }_1\left(A\right)=\mathrm{0.398987923256.}$ (80)

16. Hanging String, Quadratization of the Characteristic Polynomial

Here, for vector x of Equation (75), and the quadratization of $Ax=\lambda x$ we have

$\alpha =1-\lambda ,\beta =1-\frac{5}{3}\lambda +\frac{1}{3}{\lambda }^2,\gamma =1-\frac{34}{15}\lambda +\frac{13}{15}{\lambda }^2,\delta =1-\frac{298}{105}\lambda +\frac{165}{105}{\lambda }^2,$ (81)

and, from the last equation of $Ax=\lambda x$, the quadratic approximate characteristic polynomial

$767{\lambda }^2-1069\lambda +315=0,\lambda =0.4231227,$ (82)

and have from this $\lambda$ that

$x=\left\{1,0.577,0.354,0.196,0.080\right\},R\left[x\right]=\mathrm{0.4022336291.}$ (83)

17. Hanging String, Power Method

The global stiffness matrix we take is still the modest $5\times 5$

$A=\left[\begin{array}{ccccc}1& -1& & & \\ -1& 4& -3& & \\ & -3& 8& -5& \\ & & -5& 12& -7\\ & & & -7& 16\end{array}\right],\begin{array}{ll}\hfill & {\lambda }_1\left(A\right)=0.398987923256\hfill \\ \hfill & {\lambda }_5\left(A\right)=22.007479\hfill \end{array}$ (84)

with corresponding eigenvectors

$x_1\left(A\right)=\left\{\mathrm{9.71716,5.84013,3.77107,2.22872,1}\right\}$

$x_5\left(A\right)=\left\{0.00252621,-0.0530693,0.317706,-0.858211,1\right\}.$ (85)

We start the power method with

$x_0=\left\{1,-1,1,-1,1\right\},R\left[x_0\right]=14.6$ (86)

to successively produce

$\begin{array}{l}R\left[x_1\right]=21.0301,R\left[x_2\right]=21.7884,R\left[x_3\right]=21.9497,\\ R\left[x_4\right]=21.9917,R\left[x_5\right]=22.0031,\end{array}$ (87)

then

$R\left[x_6\right]=22.0063$ for

$x_6=A{x}_5=\left\{0.00296949,-0.0578104,0.330741,-0.869775,1\right\}.$ (88)

It follows then, from Rayleigh and Gershgorin, that

$22.0063\le {\lambda }_6\left(A\right)\le 24,{\lambda }_6\left(A\right)=\mathrm{22.0075.}$ (89)

18. Hanging String, Shifted Power Method

Next we turn our attention to computing, an ever better, approximation to ${\lambda }_1\left(A\right)$, the lowest eigenvalue of global stiffness matrix A of the hanging string. For this we take the shifted matrix

$B=R\left[x_6\right]I-A,\lambda \left(B\right)=R\left[x_6\right]-\lambda \left(A\right).$ (90)

The power method applied to matrix B converges now to $x_1$, the eigenvector corresponding to the lowest eigenvalue of original matrix A.

We start with

$x_0=\left\{9,5,3,2,1\right\},R\left[x_0\right]=0.408333,G\left[x_0\right]=-0.333333$ (91)

and repeatedly compute

$x_1=B{x}_0=\left\{1,0.556701,0.345357,0.216497,0.103095\right\},$

$R\left[x_1\right]=0.403054,G\left[x_1\right]=0.02975,$

$x_2=B{x}_1=\left\{1,0.5593,0.351981,0.214014,0.0989981\right\},$

$R\left[x_2\right]=0.401485,G\left[x_2\right]=0.192842,$

$x_3=B{x}_2=\left\{1,0.562325,0.356026,0.213042,0.097039\right\},$

$R\left[x_3\right]=0.400867,G\left[x_3\right]=0.269716,$

$x_4=B{x}_3=\left\{1,0.565334,0.35880,0.212863,0.0961647\right\},$

$R\left[x_4\right]=0.400526,G\left[x_4\right]=0.306781,$

$x_5=B{x}_4=\left\{1,0.568152,0.360926,0.213109,0.095850\right\},$

$R\left[x_5\right]=0.400289,G\left[x_5\right]=0.325285,$

$x_6=B{x}_5=\left\{1,0.570725,0.362709,0.213587,0.0958295\right\},$

$R\left[x_6\right]=0.400105,G\left[x_6\right]=0.33514$

$x_7=B{x}_6=\left\{1,0.573053,0.364291,0.214189,0.0959672\right\},$

$R\left[x_7\right]=0.399954,G\left[x_7\right]=0.34100$

$x_8=B{x}_7=\left\{1,0.575153,0.365742,0.214857,0.0961906\right\},$

$R\left[x_8\right]=0.399827,G\left[x_8\right]=0.34503$

$x_9=B{x}_8=\left\{1,0.577051,0.367094,0.215554,0.0964599\right\},$

$R\left[x_9\right]=0.399719,G\left[x_9\right]=0.34823,$ (92)

and have that

$0.399719\le {\lambda }_1\left(A\right)\le 0.348229,{\lambda }_1\left(A\right)=0.398988,$ (93)

which we gladly accept.

19. The Four-Nodes Rectangular Membrane Element

Next we move on to two dimensional finite-element membrane problems, or the boundary value problem

$\frac{{\partial }^{2}u}{\partial x^2}+\frac{{\partial }^{2}u}{\partial y^2}+f\left(x\mathrm{,}y\right)=0\text{in}D$

$u=0\text{on}{S}_1,\text{and}\frac{\partial u}{\partial n}=0\text{on}{S}_2$ (94)

where n is an outwardly normal to boundary S encompassing domain D.

The element stiffness and mass matrices for the rectangular four-nodes membrane finite element of sides a and b are

$k=\frac{b}{6a}\left[\begin{array}{cccc}2& -2& 1& -1\\ -2& 2& -1& 1\\ 1& -1& 2& -2\\ -1& 1& -2& 2\end{array}\right]+\frac{a}{6b}\left[\begin{array}{cccc}2& 1& -2& -1\\ 1& 2& -1& -2\\ -2& -1& 2& 1\\ -1& -2& 1& 2\end{array}\right]$ (95)

and

$m=\frac{ab}{36}\left[\begin{array}{cccc}4& 2& 2& 1\\ 2& 4& 1& 2\\ 2& 1& 4& 2\\ 1& 2& 2& 4\end{array}\right]$ (96)

respectively. If $a=b=h$, then

$k=\frac{1}{6}\left[\begin{array}{cccc}4& -1& -1& -2\\ -1& 4& -2& -1\\ -1& -2& 4& -1\\ -2& -1& -1& 4\end{array}\right]\text{and}m=\frac{h^2}{36}\left[\begin{array}{cccc}4& 2& 2& 1\\ 2& 4& 1& 2\\ 2& 1& 4& 2\\ 1& 2& 2& 4\end{array}\right].$ (97)

with eigenvalues

$\lambda \left(k\right)=\frac{1}{6}\left\{0,4,6,6\right\},\lambda \left(m\right)=\frac{h^2}{36}\left\{1,3,3,9\right\}.$ (98)

Element stiffness matrix k of Equation (97) is of a positive diagonal and non-positive off-diagonal entries, and hence such is also any global matrix assembled from it.

Global stiffness matrix A for the modest $3\times 3$ mesh is

$A=\left[\begin{array}{ccccccccc}16& -2& 0& -2& -2& 0& 0& 0& 0\\ -2& 16& -2& -2& -2& -2& 0& 0& 0\\ 0& -2& 8& 0& -2& -1& 0& 0& 0\\ -2& -2& 0& 16& -2& 0& -2& -2& 0\\ -2& -2& -2& -2& 16& -2& -2& -2& -2\\ 0& -2& -1& 0& -2& 8& 0& -2& -1\\ 0& 0& 0& -2& -2& 0& 8& -1& 0\\ 0& 0& 0& -2& -2& -2& -1& 8& -1\\ 0& 0& 0& 0& -2& -1& 0& -1& 4\end{array}\right]$ (99)

with, as expected, $A_{ii}>0$ and $A_{ij}\le 0$. The eigenvalues and eigenvectors of global stiffness A are

$x_1=\left\{\mathrm{0.14854,0.2737,0.3541,0.2737,0.5209,0.6713,0.3541,0.6713,1}\right\}$

${\lambda }_1=1.61549$

$x_2=-\left\{0.1258,0.242,0.4179,0.242,0.2034,0.3329,0.418,0.3329,-1\right\}$

${\lambda }_2=5.07254$

$x_3=\left\{\mathrm{0,}-\mathrm{0.538359,}-\mathrm{1.97486,0.538359,0,}-\mathrm{1,1.97486,1,0}\right\}$

${\lambda }_3=6.94842$

$x_4=\left\{0.385,0.4198,2.934,0.4198,0.6593,-2.762,2.935,-2.762,1\right\}$

${\lambda }_4=8.20581$

$x_5=\left\{0,-0.119,0.5390,0.1197,0,-1,-0.5390,1,0\right\}$

${\lambda }_5=10.2995$

$x_6=\left\{-26.991,-19.691,13.26,-19.692,-12.05,7.953,13.26,7.953,1\right\}$

${\lambda }_6=12.189$

$x_7=\left\{-9.224,8.124,-0.246,8.124,-6.739,-0.292,-0.246,-0.292,1\right\}$

${\lambda }_7=18.0621$

$x_8=\left\{0,4.773,-0.7950,-4.773,0,-1,0.7950,1,0\right\}$

${\lambda }_8=18.752$

$x_9=\left\{6.105,-0.0677,1.487,-0.0677,-8.579,1.152,1.487,1.152,1\right\}$

${\lambda }_9=\mathrm{18.855.}$ (100)

20. Rectangular Membrane, Power Method

We prefer the power method as it requires only a vector-matrix multiplication, We start the method with the initial $x_0$

$x_0=\left\{6,0,2,0,-10,1,2,1,1\right\},R\left[x_0\right]=18.8027$ (101)

and keep multiplying, or raising to ever higher power, to have.

$x_1=A{x}_0=\left\{10,0.1724,3.017,0.1724,-16.03,1.983,3.017,1.983,1.897\right\}$

$R\left[x_1\right]=18.8418$

$x_2=A{x}_1=\left\{10,0.234,2.815,0.2342,-15.730,2.023,2.815,2.023,1.865\right\}$

$R\left[x_2\right]=18.8491$

$x_3=A{x}_2=\left\{10,0.266,2.703,0.266,-15.52,2.018,2.702,2.018,1.830\right\}$

$R\left[x_3\right]=18.8516$

$x_4=A{x}_3=\left\{10,0.280,2.638,0.280,-15.37,2.005,2.638,2.005,1.807\right\}$

$R\left[x_4\right]=18.8526$

$x_5=A{x}_4=\left\{10,0.283,2.5985,0.283,-15.25,1.991,2.599,1.99,1.791\right\}$

$R\left[x_5\right]=18.8531$

$x_6=A{x}_5=\left\{10,0.2791,2.573,0.2791,-15.16,1.98,2.573,1.980,1.779\right\}$

$R\left[x_6\right]=18.8533$

$x_7=A{x}_6=\left\{10,0.271,2.556,0.271,-15.09,1.971,2.556,1.971,1.769\right\}$

$R\left[x_7\right]=18.8535$

$x_8=\left\{10,0.2595,2.544,0.260,-15.025,1.964,2.544,1.964,1.76\right\}$

$R\left[x_8\right]=18.8537$ (102)

such that, by Rayleigh and Gershgorin

$18.8537\le {\lambda }_9\left[A\right]\le 32.$ (103)

21. Rectangular Membrane. Shifted Power Method

Now we go after the, more interesting, lowest eigenvalue of the global stiffness matrix A and replace it by the shifted matrix

$B=R\left[x_4\right]I-A.$ (104)

We start with $x_0$ and continue to have

$x_0=\left\{1,3,4,3,5,7,4,7,10\right\}$

$R\left[x_0\right]=1.64964,G\left[x_0\right]=-6$

$x_1=B{x}_0=\left\{1.441,2.815,3.85,2.815,5.348,6.954,3.85,6.95,10\right\}$

$R\left[x_1\right]=1.61878,G\left[x_1\right]=0.759863$

$x_2=B{x}_1=\left\{1.51,2.821,3.76,2.82,5.35,6.91,3.757,6.905,10\right\}$

$R\left[x_2\right]=1.61685,G\left[x_2\right]=1.39789$

$x_3=B{x}_2=\left\{1.519,2.816,3.700,2.816,5.329,6.868,3.700,6.868,10\right\}$

$R\left[x_3\right]=1.61624,G\left[x_3\right]=1.56061$

$x_4=B{x}_3=\left\{1.518,2.804,3.661,2.804,5.307,6.839,3.661,6.839,10\right\}$

$R\left[x_4\right]=1.61594,G\left[x_4\right]=1.57068$

$x_5=B{x}_4=\left\{1.514,2.793,3.634,2.793,5.289,6.815,3.634,6.816,10\right\}$

$R\left[x_5\right]=1.61577,G\left[x_5\right]=1.57916$

$x_6=B{x}_5=\left\{1.509,2.782,3.613,2.782,5.273,6.796,3.613,6.796,10\right\}$

$R\left[x_6\right]=1.61567,G\left[x_6\right]=1.58616$

$x_7=B{x}_6=\left\{1.505,2.773,3.598,2.773,5.260,6.780,3.598,6.78,10\right\}$

$R\left[x_7\right]=1.6156,G\left[x_7\right]=1.59187$

$x_8=B{x}_7=\left\{1.501,2.766,3.59,2.766,5.25,6.767,3.59,6.767,10\right\}$

$R\left[x_8\right]=1.61556,G\left[x_8\right]=1.59651,$ (105)

and

$1.59651\le {\lambda }_1\left(A\right)\le 1.61556,{\lambda }_1\left(A\right)=1.61549,$ (106)

which we appraise as quite good and tight numerical bounds.

22. Triangular Membrane, Triangular Finite Elements

The linear, first order, membrane finite element stiffness matrix k for a triangle of sides $L_1\mathrm{,}{L}_2\mathrm{,}{L}_3$ and area T is

$k=\frac{1}{4T}\left(L_1^2{C}_1+L_2^2{C}_2+L_3^2{C}_3\right)$ (107)

where $C_1\mathrm{,}{C}_2\mathrm{,}{C}_3$ are the constant matrices

$C_1=\frac{1}{2}\left[\begin{array}{ccc}2& -1& -1\\ -1& 0& 1\\ -1& 1& 0\end{array}\right],C_2=\frac{1}{2}\left[\begin{array}{ccc}0& -1& 1\\ -1& 2& -1\\ 1& -1& 0\end{array}\right],C_3=\frac{1}{2}\left[\begin{array}{ccc}0& 1& -1\\ 1& 0& -1\\ -1& -1& 2\end{array}\right].$ (108)

Otherwise, the element stiffness matrix is written, by the aid of the law of the cosines, as

$k=\frac{1}{4T}\left[\begin{array}{ccc}{L}_1^2& -L_1{L}_2\cos {\theta }_3& -L_1{L}_3\cos {\theta }_2\\ -L_1{L}_2\cos {\theta }_3& L_2^2& -L_2{L}_3\cos {\theta }_1\\ -L_1{L}_3\cos {\theta }_2& -L_2{L}_3\cos {\theta }_1& L_3^2\end{array}\right]$ (109)

where ${\theta }_i$ is the angle opposite $L_i$. If ${\theta }_i<\pi /2$, then $k_{ij}<0$, and so is any global stiffness matrix assembled from it.

We consider now an equilateral triangular membrane discretized by four equilateral triangular finite elements per side. The membrane is held fixed on one side, to produce the 10 × 10 symmetric and positive definite global stiffness matrix

$A=\left[\begin{array}{cccccccccc}6& -2& 0& 0& -1& 0& 0& 0& 0& 0\\ -2& 12& -2& 0& -2& -2& 0& 0& 0& 0\\ 0& -2& 12& -2& 0& -2& -2& 0& 0& 0\\ 0& 0& -2& 6& 0& 0& -1& 0& 0& 0\\ -1& -2& 0& 0& 6& -2& 0& -1& 0& 0\\ 0& -2& -2& 0& -2& 12& -2& -2& -2& 0\\ 0& 0& -2& -1& 0& -2& 6& 0& -1& 0\\ 0& 0& 0& 0& -1& -2& 0& 6& -2& -1\\ 0& 0& 0& 0& 0& -2& -1& -2& 6& -1\\ 0& 0& 0& 0& 0& 0& 0& -1& -1& 2\end{array}\right]$ (110)

with, as expected, $A_{ii}>0$ and $A_{ij}\le 0$. The ten computed eigenvalues of this matrix are:

$\lambda =\left\{0.61,2.29,4.06,4.74,6.67,6.98,8.45,10.55,14.82,14.83\right\}.$ (111)

23. Triangular Membrane, Power Method

We start our iterative quest for the highest eigenvalue of matrix A with $x_0$, and then continue with the power method to have

$x_0=\left\{1,-1,-1,1,-1,1,-1,-1,-1,1\right\}$

$x_1=A{x}_0=\left\{0.38,-0.5,-0.5,0.38,-0.25,1,-0.25,-0.25,-0.25,0.17\right\}$

$R\left[x_1\right]=13.9798$

$x_2=A{x}_1=\left\{0.22,-0.45,-0.45,0.22,-0.16,1,-0.16,-0.18,-0.18,0.05\right\}$

$R\left[x_2\right]=14.6874$

$x_3=A{x}_2=\left\{0.16,-0.44,-0.44,0.16,-0.14,1,-0.14,-0.17,-0.17,0.03\right\}$

$R\left[x_3\right]=14.8011$

$x_4=A{x}_3=\left\{0.13,-0.43,-0.43,0.13,-0.13,1,-0.13,-0.17,-0.17,0.03\right\}$

$R\left[x_4\right]=14.8234$

$x_5=A{x}_4=\left\{0.12,-0.42,-0.42,0.12,-0.13,1,-0.13,-0.17,-0.17,0.03\right\}$

$R\left[x_5\right]=14.8286$

$x_6=A{x}_5=\left\{0.11,-0.42,-0.42,0.11,-0.13,1,-0.13,-0.17,-0.17,0.03\right\}$

$R\left[x_6\right]=14.8301$

$x_7=A{x}_6=\left\{0.11,-0.41,-0.41,0.11,-0.13,1,-0.13,-0.18,-0.18,0.03\right\}$

$R\left[x_7\right]=14.8306$

$x_8=A{x}_7=\left\{0.11,-0.41,-0.41,0.11,-0.13,1,-0.13,-0.18,-0.18,0.03\right\}$

$R\left[x_8\right]=14.8308,$ (112)

and thus, by Rayleigh and Gershgorin

$14.8308\le {\lambda }_{10}\left(A\right)\le 24,{\lambda }_{10}\left(A\right)=\mathrm{14.831.}$ (113)

24. Triangular Membrane, Shifted Power Method

Next we turn our attention to the lowest eigenvalue ${\lambda }_1\left(A\right)$ of global stiffness matrix A of the triangular membrane. For this we apply the power method again, but this time to the shifted matrix

$B=R\left[x_8\right]I-A$ (114)

in which $R\left[x_8\right]$ is the rayleigh quotient for $x_8$ obtained by the power method in the previous section.

We start with an arbitrary $x_0$ and continue with the shifted power method to have

$x_0=\left\{1,2,2,1,4,5,4,7,7,10\right\}$

$R\left[x_0\right]=0.633962$

$x_1=B{x}_0=\left\{2.10,3.71,3.71,2.10,7.17,8.27,7.17,12.48,12.48,17.8\right\}$

$R\left[x_1\right]=0.619807$

$x_2=B{x}_1=\left\{3.12,4.99,4.99,3.12,9.58,11.0,9.58,16.63,16.63,23.83\right\}$

$R\left[x_2\right]=0.616615$

$x_3=B{x}_2=\left\{3.27,4.96,4.96,3.27,9.46,10.82,9.46,16.34,16.34,23.5\right\}$

$R\left[x_3\right]=0.615228,G\left[x_3\right]=0.0722$

$x_4=B{x}_3=\left\{3.36,4.94,4.94,3.36,9.37,10.69,9.37,16.1,16.1,23.3\right\}$

$R\left[x_4\right]=0.614575,G\left[x_4\right]=0.2632$

$x_5=B{x}_4=\left\{3.41,4.94,4.94,3.41,9.32,10.6,9.32,15.9,15.9,23.1\right\}$

$R\left[x_5\right]=0.614248,G\left[x_5\right]=0.3782$

$x_6=B{x}_5=\left\{3.45,4.94,4.94,3.45,9.29,10.6,9.29,15.9,15.9,23.0\right\}$

$R\left[x_6\right]=0.614075,G\left[x_6\right]=0.4505$

$x_7=B{x}_6=\left\{3.48,4.94,4.94,3.48,9.29,10.54,9.29,15.8,15.8,22.9\right\}$

$R\left[x_7\right]=0.613978,G\left[x_7\right]=0.4976$

$x_8=B{x}_7=\left\{3.51,4.95,4.95,3.51,9.29,10.53,9.29,15.8,15.8,22.9\right\}$

$R\left[x_8\right]=0.61392,G\left[x_8\right]=0.5291$ (115)

and

$0.5291\le {\lambda }_1\left(A\right)\le 0.61392,{\lambda }_1\left(A\right)=\mathrm{0.613799.}$ (116)

25. Conclusions

In this note we have considered the common case of the positive definite and symmetric global finite element stiffness matrices having non positive off-diagonal entries. The finite element global stiffness matrix readily becomes very large with an ever smaller least eigenvalue.

We have shown here how optimal similarity transformations of this matrix, requiring only matrix-vector multiplication, lead to eminently practical and reliable numerical iterative algorithms for tight realistic Rayleigh and Gershgorin bounds on the least eigenvalue of this large finite-elements global stiffness matrix.

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Strang, G. and Fix, G. (2008) An Analysis of the Finite Element Method. 2nd Edition, Wellesley Cambridge Press, Wellesley, Section 1.9.
[2]	Fried, I. (1973) Bounds on the Spectral and Maximum Norms of the Finite Element Stiffness, Flexibility and Mass Matrices. International Journal of Solids and Structures, 9, 1013-1034. https://doi.org/10.1016/0020-7683(73)90013-9
[3]	Fried, I. (1975) Monotone Finite Element Matrices and Their Computed Condition Numbers. Computers and Structures, 5, 317-319. https://doi.org/10.1016/0045-7949(75)90038-3
[4]	Mangasarian, O.L. (1968) Characterizations of Real Matrices of Monotone Kind. SIAM Review, 10, 439-441. https://doi.org/10.1137/1010095
[5]	Fujimoto, T. and Ranade, R. (2004) Two Characterizations of Inverse-Positive Matrices: The Hawkins-Simon Condition and the Le Chatelier-Braun Principle. Electronic Journal of Linear Algebra, 11, 59-65. https://www.scirp.org/journal/AM/ https://doi.org/10.13001/1081-3810.1122
[6]	Fried, I. (1970) A Computational Procedure for the Solution of Large Problems Arising from the Finite Element Method. International Journal for Numerical Methods in Engineering, 2, 477-494. https://doi.org/10.1002/nme.1620020403
[7]	Scott, D.S. (1985) On the Accuracy of the Gershgorin Circle Theorem for Bounding the Spread of a Real Symmetric Matrix. Linear Algebra and Applications, 65, 147-155. https://doi.org/10.1016/0024-3795(85)90093-X
[8]	Varga, R.S. (2004) Gershgorin and His Circles. Springer, New York. https://doi.org/10.1007/978-3-642-17798-9
[9]	MacCluer, C.R. (2000) The Many Proofs and Applications of Perron’s Theorem. SIAM Review, 42, 487-498. https://doi.org/10.1137/S0036144599359449
[10]	Fried, I. (2016) The Condition of the Least-Squares Finite Element Matrices. International Journal for Numerical Methods in Engineering, 106, 760-770. https://doi.org/10.1002/nme.5146
[11]	Fried, I. (1972) Condition of Finite Element Matrices Generated from Non-Uniform Meshes. AIAA Journal, 10, 219-221. https://doi.org/10.2514/3.6561
[12]	Verbin, Y. (2014) Boundary Conditions and Modes of the Vertically Hanging Chain. European Journal of Physics, 36, Article ID: 015005. https://doi.org/10.1088/0143-0807/36/1/015005