First order expansion of the equation

Substituting in (1), we get

${\displaystyle (\mathbf {K} _{0}+\delta \mathbf {K} )(\mathbf {x} _{0i}+\delta \mathbf {x} _{i})=\left(\lambda _{0i}+\delta \lambda _{i}\right)\left(\mathbf {M} _{0}+\delta \mathbf {M} \right)\left(\mathbf {x} _{0i}+\delta \mathbf {x} _{i}\right),}$

which expands to

${\displaystyle {\begin{aligned}\mathbf {K} _{0}\mathbf {x} _{0i}&+\delta \mathbf {K} \mathbf {x} _{0i}+\mathbf {K} _{0}\delta \mathbf {x} _{i}+\delta \mathbf {K} \delta \mathbf {x} _{i}=\\[6pt]&\lambda _{0i}\mathbf {M} _{0}\mathbf {x} _{0i}+\lambda _{0i}\mathbf {M} _{0}\delta \mathbf {x} _{i}+\lambda _{0i}\delta \mathbf {M} \mathbf {x} _{0i}+\delta \lambda _{i}\mathbf {M} _{0}\mathbf {x} _{0i}+\\&\quad \lambda _{0i}\delta \mathbf {M} \delta \mathbf {x} _{i}+\delta \lambda _{i}\delta \mathbf {M} \mathbf {x} _{0i}+\delta \lambda _{i}\mathbf {M} _{0}\delta \mathbf {x} _{i}+\delta \lambda _{i}\delta \mathbf {M} \delta \mathbf {x} _{i}.\end{aligned}}}$

Canceling from (0) (${\displaystyle \mathbf {K} _{0}\mathbf {x} _{0i}=\lambda _{0i}\mathbf {M} _{0}\mathbf {x} _{0i}}$) leaves

${\displaystyle {\begin{aligned}\delta \mathbf {K} \mathbf {x} _{0i}+&\mathbf {K} _{0}\delta \mathbf {x} _{i}+\delta \mathbf {K} \delta \mathbf {x} _{i}=\lambda _{0i}\mathbf {M} _{0}\delta \mathbf {x} _{i}+\lambda _{0i}\delta \mathbf {M} \mathbf {x} _{0i}+\delta \lambda _{i}\mathbf {M} _{0}\mathbf {x} _{0i}+\\&\lambda _{0i}\delta \mathbf {M} \delta \mathbf {x} _{i}+\delta \lambda _{i}\delta \mathbf {M} \mathbf {x} _{0i}+\delta \lambda _{i}\mathbf {M} _{0}\delta \mathbf {x} _{i}+\delta \lambda _{i}\delta \mathbf {M} \delta \mathbf {x} _{i}.\end{aligned}}}$

Removing the higher-order terms, this simplifies to

${\displaystyle \mathbf {K} _{0}\delta \mathbf {x} _{i}+\delta \mathbf {K} \mathbf {x} _{0i}=\lambda _{0i}\mathbf {M} _{0}\delta \mathbf {x} _{i}+\lambda _{0i}\delta \mathbf {M} \mathrm {x} _{0i}+\delta \lambda _{i}\mathbf {M} _{0}\mathbf {x} _{0i}.\qquad (3)}$

In other words, ${\displaystyle \delta \lambda _{i}}$ no longer denotes the exact variation of the eigenvalue but its first order approximation.

As the matrix is symmetric, the unperturbed eigenvectors are ${\displaystyle M}$ orthogonal and so we use them as a basis for the perturbed eigenvectors. That is, we want to construct

${\displaystyle \delta \mathbf {x} _{i}=\sum _{j=1}^{N}\varepsilon _{ij}\mathbf {x} _{0j}\qquad (4)\quad }$ with ${\displaystyle \varepsilon _{ij}=\mathbf {x} _{0j}^{T}M\delta \mathbf {x} _{i}}$,

where the ε_ij are small constants that are to be determined.

In the same way, substituting in (2), and removing higher order terms, we get ${\displaystyle \delta \mathbf {x} _{j}\mathbf {M} _{0}\mathbf {x} _{0i}+\mathbf {x} _{0j}\mathbf {M} _{0}\delta \mathbf {x} _{i}+\mathbf {x} _{0j}\delta \mathbf {M} _{0}\mathbf {x} _{0i}=0\quad {(5)}}$

The derivation can go on with two forks.

Perturbation of an implicit function.

In the next paragraph, we shall use the Implicit function theorem (Statement of the theorem ); we notice that for a continuously differentiable function ${\displaystyle f:\mathbb {R} ^{n+m}\to \mathbb {R} ^{m},\;f:(x,y)\mapsto f(x,y)}$, with an invertible Jacobian matrix ${\displaystyle J_{f,b}(x_{0},y_{0})}$, from a point ${\displaystyle (x_{0},y_{0})}$ solution of ${\displaystyle f(x_{0},y_{0})=0}$, we get solutions of ${\displaystyle f(x,y)=0}$ with ${\displaystyle x}$ close to ${\displaystyle x_{0}}$ in the form ${\displaystyle y=g(x)}$ where ${\displaystyle g}$ is a continuously differentiable function ; moreover the Jacobian marix of ${\displaystyle g}$ is provided by the linear system

${\displaystyle J_{f,y}(x,g(x))J_{g,x}(x)+J_{f,x}(x,g(x))=0\quad (6)}$.

As soon as the hypothesis of the theorem is satisfied, the Jacobian matrix of ${\displaystyle g}$ may be computed with a first order expansion of ${\displaystyle f(x_{0}+\delta x,y_{0}+\delta y)=0}$, we get

${\displaystyle J_{f,x}(x,g(x))\delta x+J_{f,y}(x,g(x))\delta y=0}$; as ${\displaystyle \delta y=J_{g,x}(x)\delta x}$, it is equivalent to equation ${\displaystyle (6)}$.

Eigenvalue perturbation: a theoretical basis.

We use the previous paragraph (Perturbation of an implicit function) with somewhat different notations suited to eigenvalue perturbation; we introduce ${\displaystyle {\tilde {f}}:\mathbb {R} ^{2n^{2}}\times \mathbb {R} ^{n+1}\to \mathbb {R} ^{n+1}}$, with

• ${\displaystyle {\tilde {f}}(K,M,\lambda ,x)={\binom {f(K,M,\lambda ,x)}{f_{n+1}(x)}}}$ with

${\displaystyle f(K,M,\lambda ,x)=Kx-\lambda x,f_{n+1}(M,x)=x^{T}Mx-1}$. In order to use the Implicit function theorem, we study the invertibility of the Jacobian ${\displaystyle J_{{\tilde {f}};\lambda ,x}(K,M;\lambda _{0i},x_{0i})}$ with

${\displaystyle J_{{\tilde {f}};\lambda ,x}(K,M;\lambda _{i},x_{i})(\delta \lambda ,\delta x)={\binom {-Mx_{i}}{0}}\delta \lambda +{\binom {K-\lambda M}{2x_{i}^{T}M}}\delta x_{i}}$. Indeed, the solution of

${\displaystyle J_{{\tilde {f}};\lambda _{0i},x_{0i}}(K,M;\lambda _{0i},x_{0i})(\delta \lambda _{i},\delta x_{i})=}$${\displaystyle {\binom {y}{y_{n+1}}}}$ may be derived with computations similar to the derivation of the expansion.

$${\displaystyle \delta \lambda _{i}=-x_{0i}^{T}y,\;{\text{ and }}(\lambda _{0i}-\lambda _{0j})x_{0j}^{T}M\delta x_{i}=x_{j}^{T}y,j=1,\dots ,n,j\neq i\;;}$$

${\displaystyle {\text{ or }}x_{0j}^{T}M\delta x_{i}=x_{j}^{T}y/(\lambda _{0i}-\lambda _{0j}),{\text{ and }}\;2x_{0i}^{T}M\delta x_{i}=y_{n+1}}$

When ${\displaystyle \lambda _{i}}$ is a simple eigenvalue, as the eigenvectors ${\displaystyle x_{0j},j=1,\dots ,n}$ form an orthonormal basis, for any right-hand side, we have obtained one solution therefore, the Jacobian is invertible.

The implicit function theorem provides a continuously differentiable function ${\displaystyle (K,M)\mapsto (\lambda _{i}(K,M),x_{i}(K,M))}$ hence the expansion with little o notation: ${\displaystyle \lambda _{i}=\lambda _{0i}+\delta \lambda _{i}+o(\|\delta K\|+\|\delta M\|)}$ ${\displaystyle x_{i}=x_{0i}+\delta x_{i}+o(\|\delta K\|+\|\delta M\|)}$. with

${\displaystyle \delta \lambda _{i}=\mathbf {x} _{0i}^{T}\delta \mathbf {K} \mathbf {x} _{0i}-\lambda _{0i}\mathbf {x} _{0i}^{T}\delta \mathbf {M} \mathrm {x} _{0i};}$ ${\displaystyle \delta x_{i}=\mathbf {x} _{0j}^{T}\mathbf {M} _{0}\delta \mathbf {x} _{i}\mathbf {x} _{0j}{\text{ with}}}$${\displaystyle \mathbf {x} _{0j}^{T}\mathbf {M} _{0}\delta \mathbf {x} _{i}={\frac {-\mathbf {x} _{0j}^{T}\delta \mathbf {K} \mathbf {x} _{0i}+\lambda _{0i}\mathbf {x} _{0j}^{T}\delta \mathbf {M} \mathrm {x} _{0i}}{(\lambda _{0j}-\lambda _{0i})}},i=1,\dots n;j=1,\dots n;j\neq i.}$ This is the first order expansion of the perturbed eigenvalues and eigenvectors. which is proved.

The results

This means it is possible to efficiently do a sensitivity analysis on λ_i as a function of changes in the entries of the matrices. (Recall that the matrices are symmetric and so changing K_kℓ will also change K_ℓk, hence the (2 − δ_kℓ) term.)

${\displaystyle {\begin{aligned}{\frac {\partial \lambda _{i}}{\partial \mathbf {K} _{(k\ell )}}}&={\frac {\partial }{\partial \mathbf {K} _{(k\ell )}}}\left(\lambda _{0i}+\mathbf {x} _{0i}^{\top }\left(\delta \mathbf {K} -\lambda _{0i}\delta \mathbf {M} \right)\mathbf {x} _{0i}\right)=x_{0i(k)}x_{0i(\ell )}\left(2-\delta _{k\ell }\right)\\{\frac {\partial \lambda _{i}}{\partial \mathbf {M} _{(k\ell )}}}&={\frac {\partial }{\partial \mathbf {M} _{(k\ell )}}}\left(\lambda _{0i}+\mathbf {x} _{0i}^{\top }\left(\delta \mathbf {K} -\lambda _{0i}\delta \mathbf {M} \right)\mathbf {x} _{0i}\right)=-\lambda _{i}x_{0i(k)}x_{0i(\ell )}\left(2-\delta _{k\ell }\right).\end{aligned}}}$

Similarly

${\displaystyle {\begin{aligned}{\frac {\partial \mathbf {x} _{i}}{\partial \mathbf {K} _{(k\ell )}}}&=\sum _{j=1 \atop j\neq i}^{N}{\frac {x_{0j(k)}x_{0i(\ell )}\left(2-\delta _{k\ell }\right)}{\lambda _{0i}-\lambda _{0j}}}\mathbf {x} _{0j}\\{\frac {\partial \mathbf {x} _{i}}{\partial \mathbf {M} _{(k\ell )}}}&=-\mathbf {x} _{0i}{\frac {x_{0i(k)}x_{0i(\ell )}}{2}}(2-\delta _{k\ell })-\sum _{j=1 \atop j\neq i}^{N}{\frac {\lambda _{0i}x_{0j(k)}x_{0i(\ell )}}{\lambda _{0i}-\lambda _{0j}}}\mathbf {x} _{0j}\left(2-\delta _{k\ell }\right).\end{aligned}}}$

Eigenvalue sensitivity, a small example

A simple case is ${\displaystyle K={\begin{bmatrix}2&b\\b&0\end{bmatrix}}}$; however you can compute eigenvalues and eigenvectors with the help of online tools such as (see introduction in Wikipedia WIMS) or using Sage SageMath. You get the smallest eigenvalue ${\displaystyle \lambda =-\left[{\sqrt {b^{2}+1}}+1\right]}$ and an explicit computation ${\displaystyle {\frac {\partial \lambda }{\partial b}}={\frac {-x}{\sqrt {x^{2}+1}}}}$; more over, an associated eigenvector is ${\displaystyle {\tilde {x}}_{0}=[x,-({\sqrt {x^{2}+1}}+1))]^{T}}$; it is not an unitary vector; so ${\displaystyle x_{01}x_{02}={\tilde {x}}_{01}{\tilde {x}}_{02}/\|{\tilde {x}}_{0}\|^{2}}$; we get ${\displaystyle \|{\tilde {x}}_{0}\|^{2}=2{\sqrt {x^{2}+1}}({\sqrt {x^{2}+1}}+1)}$ and ${\displaystyle {\tilde {x}}_{01}{\tilde {x}}_{02}=-x({\sqrt {x^{2}+1}}+1)}$ ; hence ${\displaystyle x_{01}x_{02}=-{\frac {x}{2{\sqrt {x^{2}+1}}}}}$; for this example , we have checked that ${\displaystyle {\frac {\partial \lambda }{\partial b}}=2x_{01}x_{02}}$ or ${\displaystyle \delta \lambda =2x_{01}x_{02}\delta b}$.

Example 1

Consider the following matrix ${\displaystyle B(\epsilon )=\epsilon {\begin{bmatrix}\cos(2/\epsilon )&,\sin(2/\epsilon )\\\sin(2/\epsilon )&,s\cos(2/\epsilon )\end{bmatrix}}}$ and ${\displaystyle A(\epsilon )=I-e^{-1/\epsilon ^{2}}B;}$ ${\displaystyle A(0)=I.}$ For ${\displaystyle \epsilon \neq 0}$, the matrix ${\displaystyle A(\epsilon )}$ has eigenvectors ${\displaystyle \Phi ^{1}=[\cos(1/\epsilon ),-\sin(1/\epsilon )]^{T};\Phi ^{2}=[\sin(1/\epsilon ),-\cos(1/\epsilon )]^{T}}$ belonging to eigenvalues ${\displaystyle \lambda _{1}=1-e^{-1/\epsilon ^{2})},\lambda _{2}=1+e^{-1/\epsilon ^{2})}}$. Since ${\displaystyle \lambda _{1}\neq \lambda _{2}}$ for ${\displaystyle \epsilon \neq 0}$ if ${\displaystyle u^{j}(\epsilon ),j=1,2,}$ are any normalized eigenvectors belonging to ${\displaystyle \lambda _{j}(\epsilon ),j=1,2}$ respectively then ${\displaystyle u^{j}=e^{\alpha _{j}(\epsilon )}\Phi ^{j}(\epsilon )}$ where ${\displaystyle \alpha _{j},j=1,2}$ are real for ${\displaystyle \epsilon \neq 0.}$ It is obviously impossible to define ${\displaystyle \alpha _{1}(\epsilon )}$ , say, in such a way that ${\displaystyle u^{1}(\epsilon )}$ tends to a limit as ${\displaystyle \epsilon \rightarrow 0,}$ because ${\displaystyle |u^{1}(\epsilon )|=|\cos(1/\epsilon )|}$ has no limit as ${\displaystyle \epsilon \rightarrow 0.}$

Note in this example that ${\displaystyle A_{jk}(\epsilon )}$ is not only continuous but also has continuous derivatives of all orders. Rellich draws the following important consequence. << Since in general the individual eigenvectors do not depend continuously on the perturbation parameter even though the operator ${\displaystyle A(\epsilon )}$ does, it is necessary to work, not with an eigenvector, but rather with the space spanned by all the eigenvectors belonging to the same eigenvalue. >>

Example 2

This example is less nasty that the previous one. Suppose ${\displaystyle [K_{0}]}$ is the 2x2 identity matrix, any vector is an eigenvector; then ${\displaystyle u_{0}=[1,1]^{T}/{\sqrt {2}}}$ is one possible eigenvector. But if one makes a small perturbation, such as

${\displaystyle [K]=[K_{0}]+{\begin{bmatrix}\epsilon &0\\0&0\end{bmatrix}}}$

Then the eigenvectors are ${\displaystyle v_{1}=[1,0]^{T}}$ and ${\displaystyle v_{2}=[0,1]^{T}}$; they are constant with respect to ${\displaystyle \epsilon }$ so that ${\displaystyle \|u_{0}-v_{1}\|}$ is constant and does not go to zero.

Books

• Ren-Cang Li (2014). "Matrix Perturbation Theory". In Hogben, Leslie (ed.). Handbook of linear algebra (Second ed.). ISBN 978-1466507289.
• Rellich, F., & Berkowitz, J. (1969). Perturbation theory of eigenvalue problems. CRC Press.{{cite book}}: CS1 maint: multiple names: authors list (link).
• Bhatia, R. (1987). Perturbation bounds for matrix eigenvalues. SIAM.

Report

• Rellich, Franz (1954). Perturbation theory of eigenvalue problems. New-York: Courant Institute of Mathematical Sciences, New-York University.

Journal papers

• Simon, B. (1982). Large orders and summability of eigenvalue perturbation theory: a mathematical overview. International Journal of Quantum Chemistry, 21(1), 3-25.
• Crandall, M. G., & Rabinowitz, P. H. (1973). Bifurcation, perturbation of simple eigenvalues, and linearized stability. Archive for Rational Mechanics and Analysis, 52(2), 161-180.
• Stewart, G. W. (1973). Error and perturbation bounds for subspaces associated with certain eigenvalue problems. SIAM review, 15(4), 727-764.
• Löwdin, P. O. (1962). Studies in perturbation theory. IV. Solution of eigenvalue problem by projection operator formalism. Journal of Mathematical Physics, 3(5), 969-982.