Eigenvalue perturbation

Setting of perturbation for a generalized eigenvalue problem

Suppose the solutions to the generalized eigenvalue problem are known to be

${\displaystyle {\mathbf{K}}_0{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{M}}_0{\mathbf{x}}_{0i}.(0)}$

where ${\displaystyle {\mathbf{K}}_0}$ and ${\displaystyle {\mathbf{M}}_0}$ are matrices. That is, we know the eigenvalues λ_0i and eigenvectors x_0i for i = 1, ..., N. An important note is that the eigenvalues are required to be distinct.

In order to perturb the matrices, one must find the eigenvalues and eigenvectors of

${\displaystyle \mathbf{K}{\mathbf{x}}_i={\lambda }_i\mathbf{M}{\mathbf{x}}_i(1)}$

where

${\displaystyle \begin{array}{cc}\mathbf{K}& ={\mathbf{K}}_0+\delta \mathbf{K}\\ \mathbf{M}& ={\mathbf{M}}_0+\delta \mathbf{M}\end{array}}$

with the perturbations ${\displaystyle \delta \mathbf{K}}$ and ${\displaystyle \delta \mathbf{M}}$ much smaller than ${\displaystyle \mathbf{K}}$ and ${\displaystyle \mathbf{M}}$ respectively. Then the new eigenvalues and eigenvectors are expected to be similar to the original, plus small perturbations:

${\displaystyle \begin{array}{cc}{\lambda }_i& ={\lambda }_{0i}+\delta {\lambda }_i\\ {\mathbf{x}}_i& ={\mathbf{x}}_{0i}+\delta {\mathbf{x}}_i\end{array}}$

Steps

Under the assumption that the matrices are symmetric, positive definite, and assume the eigenvectors are scaled such that

${\displaystyle {\mathbf{x}}_{0j}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}={\delta }_{ij},}$${\displaystyle {\mathbf{x}}_i^T\mathbf{M}{\mathbf{x}}_j={\delta }_{ij}(2)}$

where δ_ij is the Kronecker delta. Now the equation to be solved is

${\displaystyle \mathbf{K}{\mathbf{x}}_i-{\lambda }_i\mathbf{M}{\mathbf{x}}_i=0.}$

In this article, the study is restricted to first order perturbation.

First order expansion of the equation

Substituting in (1) results in

${\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{array}{cc}{\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=\\ & {\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}+\\ & {\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{array}}$

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

${\displaystyle \begin{array}{cc}\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{array}}$

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}.(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 can be used as a basis for the perturbed eigenvectors. This is the same as

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

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

In the same way, substituting in (2), and removing higher order terms, ${\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(5)}$

The derivation is then split into two paths.

First path: get first eigenvalue perturbation

Eigenvalue perturbation

Starting with (3)${\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};}$which is then left multiplied with ${\displaystyle {\mathbf{x}}_{0i}^T}$ along with using (2) as well as its first order variation (5); one gets

${\displaystyle {\mathbf{x}}_{0i}^T\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^T\delta \mathbf{M}{\mathrm{x}}_{0i}+\delta {\lambda }_i}$

or

${\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}}$

This is the first order perturbation of the generalized Rayleigh quotient with fixed ${\displaystyle x_{0i}}$: ${\displaystyle R(K,M;x_{0i})=x_{0i}^{T}K{x}_{0i}/x_{0i}^{T}M{x}_{0i},\text{ with }{x}_{0i}^{T}M{x}_{0i}=1}$

Moreover, for ${\displaystyle M=I}$, the formula ${\displaystyle \delta {\lambda }_i=x_{0i}^T\delta K{x}_{0i}}$ should be compared with Bauer-Fike theorem which provides a bound for eigenvalue perturbation.

Eigenvector perturbation

One then left multiplies (3) with ${\displaystyle x_{0j}^T}$ for ${\displaystyle j\ne i}$ and get

${\displaystyle {\mathbf{x}}_{0j}^T{\mathbf{K}}_0\delta {\mathbf{x}}_i+{\mathbf{x}}_{0j}^T\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0j}^T{\mathbf{M}}_0\delta {\mathbf{x}}_i+{\lambda }_{0i}{\mathbf{x}}_{0j}^T\delta \mathbf{M}{\mathrm{x}}_{0i}+\delta {\lambda }_i{\mathbf{x}}_{0j}^T{\mathbf{M}}_0{\mathbf{x}}_{0i}.}$

Recalling that ${\displaystyle {\mathbf{x}}_{0j}^{T}K={\lambda }_{0j}{\mathbf{x}}_{0j}^{T}M\text{ and }{\mathbf{x}}_{0j}^T{\mathbf{M}}_0{\mathbf{x}}_{0i}=0,}$ for ${\displaystyle j\ne i}$, one may substitute for

${\displaystyle {\lambda }_{0j}{\mathbf{x}}_{0j}^T{\mathbf{M}}_0\delta {\mathbf{x}}_i+{\mathbf{x}}_{0j}^T\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0j}^T{\mathbf{M}}_0\delta {\mathbf{x}}_i+{\lambda }_{0i}{\mathbf{x}}_{0j}^T\delta \mathbf{M}{\mathrm{x}}_{0i}.}$

or

${\displaystyle ({\lambda }_{0j}-{\lambda }_{0i}){\mathbf{x}}_{0j}^T{\mathbf{M}}_0\delta {\mathbf{x}}_i+{\mathbf{x}}_{0j}^T\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0j}^T\delta \mathbf{M}{\mathrm{x}}_{0i}.}$

As the eigenvalues are assumed to be simple, for ${\displaystyle j\ne i}$

${\displaystyle {ϵ}_{ij}={\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\ne i.}$

Moreover (5) (the first order variation of (2) ) yields ${\displaystyle 2{ϵ}_{ii}=2{\mathbf{x}}_{0i}^T{\mathbf{M}}_0\delta x_i=-{\mathbf{x}}_{0i}^T\delta M{\mathbf{x}}_{0i}.}$ All of the components of ${\displaystyle \delta x_i}$ have now been obtained.

Second path: Straightforward manipulations

Substituting (4) into (3) and rearranging gives

${\displaystyle \begin{array}{cccc}{\mathbf{K}}_0\sum _{j=1}^N{\epsilon }_{ij}{\mathbf{x}}_{0j}+\delta \mathbf{K}{\mathbf{x}}_{0i}& ={\lambda }_{0i}{\mathbf{M}}_0\sum _{j=1}^N{\epsilon }_{ij}{\mathbf{x}}_{0j}+{\lambda }_{0i}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{M}}_0{\mathbf{x}}_{0i}& & (5)\\ \sum _{j=1}^N{\epsilon }_{ij}{\mathbf{K}}_0{\mathbf{x}}_{0j}+\delta \mathbf{K}{\mathbf{x}}_{0i}& ={\lambda }_{0i}{\mathbf{M}}_0\sum _{j=1}^N{\epsilon }_{ij}{\mathbf{x}}_{0j}+{\lambda }_{0i}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{M}}_0{\mathbf{x}}_{0i}& & \\ (\text{applying }{\mathbf{K}}_0\text{ to the sum})\\ \sum _{j=1}^N{\epsilon }_{ij}{\lambda }_{0j}{\mathbf{M}}_0{\mathbf{x}}_{0j}+\delta \mathbf{K}{\mathbf{x}}_{0i}& ={\lambda }_{0i}{\mathbf{M}}_0\sum _{j=1}^N{\epsilon }_{ij}{\mathbf{x}}_{0j}+{\lambda }_{0i}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{M}}_0{\mathbf{x}}_{0i}& & (\text{using Eq. }(1))\end{array}}$

Because the eigenvectors are M₀-orthogonal when M₀ is positive definite, one can remove the summations by left-multiplying by ${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}}$:

${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}{\epsilon }_{ii}{\lambda }_{0i}{\mathbf{M}}_0{\mathbf{x}}_{0i}+{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\epsilon }_{ii}{\mathbf{x}}_{0i}+{\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}.}$

By use of equation (1) again:

${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{K}}_0{\epsilon }_{ii}{\mathbf{x}}_{0i}+{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\epsilon }_{ii}{\mathbf{x}}_{0i}+{\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}.(6)}$

The two terms containing ε_ii are equal because left-multiplying (1) by ${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}}$ gives

${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{K}}_0{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}.}$

Canceling those terms in (6) leaves

${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}.}$

Rearranging gives

${\displaystyle \delta {\lambda }_i=\frac{{\mathbf{x}}_{0i}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}}{{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}}}$

But by (2), this denominator is equal to 1. Thus

${\displaystyle \delta {\lambda }_i={\mathbf{x}}_{0i}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}.}$

Then, as ${\displaystyle {\lambda }_i\ne {\lambda }_k}$ for ${\displaystyle i\ne k}$ (this is the assumption of simple eigenvalues) by left-multiplying equation (5) by ${\displaystyle {\mathbf{x}}_{0k}^{\mathrm{\top }}}$:

${\displaystyle {\epsilon }_{ik}=\frac{{\mathbf{x}}_{0k}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}}{{\lambda }_{0i}-{\lambda }_{0k}},i\ne k.}$

Or by changing the name of the indices:

${\displaystyle {\epsilon }_{ij}=\frac{{\mathbf{x}}_{0j}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}}{{\lambda }_{0i}-{\lambda }_{0j}},i\ne j.}$

To find ε_ii, using the fact that

${\displaystyle {\mathbf{x}}_i^{\mathrm{\top }}\mathbf{M}{\mathbf{x}}_i=1}$

implies:

${\displaystyle {\epsilon }_{ii}=-\frac{1}{2}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}.}$

Summary of the first order perturbation result

In the case where all the matrices are Hermitian positive definite and all the eigenvalues are distinct,

${\displaystyle \begin{array}{cc}{\lambda }_i& ={\lambda }_{0i}+{\mathbf{x}}_{0i}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}\\ {\mathbf{x}}_i& ={\mathbf{x}}_{0i}\left(1-\frac{1}{2}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}\right)+\sum _{\genfrac{}{}{0ex}{}{j=1}{j\ne i}}^N\frac{{\mathbf{x}}_{0j}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}}{{\lambda }_{0i}-{\lambda }_{0j}}{\mathbf{x}}_{0j}\end{array}}$

for infinitesimal ${\displaystyle \delta \mathbf{K}}$ and ${\displaystyle \delta \mathbf{M}}$ (the higher order terms in (3) being neglected).

A proof that higher order terms may be neglect may be derived using the implicit function theorem.

Theoretical derivation

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

<span class="mwe-math-element mwe-math-element-inline" typeof="mw:Extension/math" about="#mwt108" data-mw='{"name":"math","attrs":{},"body":{"extsrc":" J_{f,y}(x,g(x)) J_{g,x}(x)+J_{f,x}(x,g(x))=0 \\quad (6) "}}' id="mwAQI"><span class="mwe-math-mathml-a11y" style="display: none;"><math xmlns="http://www.w3.org/1998/Math/MathML" class="mathjax_ignore" alttext="{\displaystyle J_{f,y}(x,g(x))J_{g,x}(x)+J_{f,x}(x,g(x))=0\quad (6)}"> <semantics> <mrow class="MJX-TeXAtom-ORD"> <mstyle displaystyle="true" scriptlevel="0"> <msub> <mi>J</mi> <mrow class="MJX-TeXAtom-ORD"> <mi>f</mi> <mo>,</mo> <mi>y</mi> </mrow> </msub> <mo stretchy="false">(</mo> <mi>x</mi> <mo>,</mo> <mi>g</mi> <mo stretchy="false">(</mo> <mi>x</mi> <mo stretchy="false">)</mo> <mo stretchy="false">)</mo> <msub> <mi>J</mi> <mrow class="MJX-TeXAtom-ORD"> <mi>g</mi> <mo>,</mo> <mi>x</mi> </mrow> </msub> <mo stretchy="false">(</mo> <mi>x</mi> <mo stretchy="false">)</mo> <mo>+</mo> <msub> <mi>J</mi> <mrow class="MJX-TeXAtom-ORD"> <mi>f</mi> <mo>,</mo> <mi>x</mi> </mrow> </msub> <mo stretchy="false">(</mo> <mi>x</mi> <mo>,</mo> <mi>g</mi> <mo stretchy="false">(</mo> <mi>x</mi> <mo stretchy="false">)</mo> <mo stretchy="false">)</mo> <mo>=</mo> <mn>0</mn> <mspace width="1em"></mspace> <mo stretchy="false">(</mo> <mn>6</mn> <mo stretchy="false">)</mo> </mstyle> </mrow> <annotation encoding="application/x-tex">{\displaystyle J_{f,y}(x,g(x))J_{g,x}(x)+J_{f,x}(x,g(x))=0\quad (6)}</annotation> </semantics> </math></span><img src="https://wikimedia.org/api/rest_v1/media/math/render/svg/9e86fa6289ecb94fde259da91675bc60ea0b11f0" class="mwe-math-fallback-image-inline mw-invert skin-invert" aria-hidden="true" style="vertical-align: -1.005ex; width:43.625ex; height:3.009ex;" alt="{\displaystyle J_{f,y}(x,g(x))J_{g,x}(x)+J_{f,x}(x,g(x))=0\quad (6)}"/></span>.

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)=(\genfrac{}{}{0ex}{}{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)=(\genfrac{}{}{0ex}{}{-M{x}_i}{0})\delta \lambda +(\genfrac{}{}{0ex}{}{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 (\genfrac{}{}{0ex}{}{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\ne 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(\Vert \delta K\Vert +\Vert \delta M\Vert )}$ ${\displaystyle x_i=x_{0i}+\delta x_i+o(\Vert \delta K\Vert +\Vert \delta M\Vert )}$. 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\ne i.}$ This is the first order expansion of the perturbed eigenvalues and eigenvectors. which is proved.

Results of sensitivity analysis with respect to the entries of the matrices

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

Why generalized eigenvalues?

In the entry applications of eigenvalues and eigenvectors we find numerous scientific fields in which eigenvalues are used to obtain solutions. Generalized eigenvalue problems are less widespread but are a key in the study of vibrations. They are useful when we use the Galerkin method or Rayleigh-Ritz method to find approximate solutions of partial differential equations modeling vibrations of structures such as strings and plates; the paper of Courant (1943) is fundamental. The Finite element method is a widespread particular case.

In classical mechanics, generalized eigenvalues may crop up when we look for vibrations of multiple degrees of freedom systems close to equilibrium; the kinetic energy provides the mass matrix ${\displaystyle M}$, the potential strain energy provides the rigidity matrix ${\displaystyle K}$. For further details, see the first section of this article of Weinstein (1941, in French)

With both methods, we obtain a system of differential equations or Matrix differential equation ${\displaystyle M\ddot{x}+B\dot{x}+Kx=0}$ with the mass matrix ${\displaystyle M}$ , the damping matrix ${\displaystyle B}$ and the rigidity matrix ${\displaystyle K}$. If we neglect the damping effect, we use ${\displaystyle B=0}$, we can look for a solution of the following form ${\displaystyle x=e^{i\omega t}u}$; we obtain that ${\displaystyle u}$ and ${\displaystyle {\omega }^2}$are solution of the generalized eigenvalue problem ${\displaystyle -{\omega }^{2}Mu+Ku=0}$

Setting of perturbation for a generalized eigenvalue problem

Suppose we have solutions to the generalized eigenvalue problem,

${\displaystyle {\mathbf{K}}_0{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{M}}_0{\mathbf{x}}_{0i}.(0)}$

where ${\displaystyle {\mathbf{K}}_0}$ and ${\displaystyle {\mathbf{M}}_0}$ are matrices. That is, we know the eigenvalues λ_0i and eigenvectors x_0i for i = 1, ..., N. It is also required that the eigenvalues are distinct.

Now suppose we want to change the matrices by a small amount. That is, we want to find the eigenvalues and eigenvectors of

${\displaystyle \mathbf{K}{\mathbf{x}}_i={\lambda }_i\mathbf{M}{\mathbf{x}}_i(1)}$

where

${\displaystyle \begin{array}{cc}\mathbf{K}& ={\mathbf{K}}_0+\delta \mathbf{K}\\ \mathbf{M}& ={\mathbf{M}}_0+\delta \mathbf{M}\end{array}}$

with the perturbations ${\displaystyle \delta \mathbf{K}}$ and ${\displaystyle \delta \mathbf{M}}$ much smaller than ${\displaystyle \mathbf{K}}$ and ${\displaystyle \mathbf{M}}$ respectively. Then we expect the new eigenvalues and eigenvectors to be similar to the original, plus small perturbations:

${\displaystyle \begin{array}{cc}{\lambda }_i& ={\lambda }_{0i}+\delta {\lambda }_i\\ {\mathbf{x}}_i& ={\mathbf{x}}_{0i}+\delta {\mathbf{x}}_i\end{array}}$

Steps

We assume that the matrices are symmetric and positive definite, and assume we have scaled the eigenvectors such that

${\displaystyle {\mathbf{x}}_{0j}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}={\delta }_{ij},}$${\displaystyle {\mathbf{x}}_i^T\mathbf{M}{\mathbf{x}}_j={\delta }_{ij}(2)}$

where δ_ij is the Kronecker delta. Now we want to solve the equation

${\displaystyle \mathbf{K}{\mathbf{x}}_i-{\lambda }_i\mathbf{M}{\mathbf{x}}_i=0.}$

In this article we restrict the study to first order perturbation.

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{array}{cc}{\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=\\ & {\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}+\\ & {\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{array}}$

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

${\displaystyle \begin{array}{cc}\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{array}}$

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}.(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{\epsilon }_{ij}{\mathbf{x}}_{0j}(4)}$ with ${\displaystyle {\epsilon }_{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(5)}$

The derivation can go on with two forks.

First fork: get first eigenvalue perturbation

Eigenvalue perturbation

We start with (3)${\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};}$

we left multiply with ${\displaystyle {\mathbf{x}}_{0i}^T}$ and use (2) as well as its first order variation (5); we get

${\displaystyle {\mathbf{x}}_{0i}^T\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^T\delta \mathbf{M}{\mathrm{x}}_{0i}+\delta {\lambda }_i}$

or

${\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}}$

We notice that it is the first order perturbation of the generalized Rayleigh quotient with fixed ${\displaystyle x_{0i}}$: ${\displaystyle R(K,M;x_{0i})=x_{0i}^{T}K{x}_{0i}/x_{0i}^{T}M{x}_{0i},\text{ with }{x}_{0i}^{T}M{x}_{0i}=1}$

Moreover, for ${\displaystyle M=I}$, the formula ${\displaystyle \delta {\lambda }_i=x_{0i}^T\delta K{x}_{0i}}$ should be compared with Bauer-Fike theorem which provides a bound for eigenvalue perturbation.

Eigenvector perturbation

We left multiply (3) with ${\displaystyle x_{0j}^T}$ for ${\displaystyle j\ne i}$ and get

${\displaystyle {\mathbf{x}}_{0j}^T{\mathbf{K}}_0\delta {\mathbf{x}}_i+{\mathbf{x}}_{0j}^T\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0j}^T{\mathbf{M}}_0\delta {\mathbf{x}}_i+{\lambda }_{0i}{\mathbf{x}}_{0j}^T\delta \mathbf{M}{\mathrm{x}}_{0i}+\delta {\lambda }_i{\mathbf{x}}_{0j}^T{\mathbf{M}}_0{\mathbf{x}}_{0i}.}$

We use ${\displaystyle {\mathbf{x}}_{0j}^{T}K={\lambda }_{0j}{\mathbf{x}}_{0j}^{T}M\text{ and }{\mathbf{x}}_{0j}^T{\mathbf{M}}_0{\mathbf{x}}_{0i}=0,}$ for ${\displaystyle j\ne i}$.

${\displaystyle {\lambda }_{0j}{\mathbf{x}}_{0j}^T{\mathbf{M}}_0\delta {\mathbf{x}}_i+{\mathbf{x}}_{0j}^T\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0j}^T{\mathbf{M}}_0\delta {\mathbf{x}}_i+{\lambda }_{0i}{\mathbf{x}}_{0j}^T\delta \mathbf{M}{\mathrm{x}}_{0i}.}$

or

${\displaystyle ({\lambda }_{0j}-{\lambda }_{0i}){\mathbf{x}}_{0j}^T{\mathbf{M}}_0\delta {\mathbf{x}}_i+{\mathbf{x}}_{0j}^T\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0j}^T\delta \mathbf{M}{\mathrm{x}}_{0i}.}$

As the eigenvalues are assumed to be simple, for ${\displaystyle j\ne i}$

${\displaystyle {ϵ}_{ij}={\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\ne i.}$

Moreover (5) (the first order variation of (2) ) yields ${\displaystyle 2{ϵ}_{ii}=2{\mathbf{x}}_{0i}^T{\mathbf{M}}_0\delta x_i=-{\mathbf{x}}_{0i}^T\delta M{\mathbf{x}}_{0i}.}$ We have obtained all the components of ${\displaystyle \delta x_i}$ .

Second fork: Straightforward manipulations

Substituting (4) into (3) and rearranging gives

${\displaystyle \begin{array}{cccc}{\mathbf{K}}_0\sum _{j=1}^N{\epsilon }_{ij}{\mathbf{x}}_{0j}+\delta \mathbf{K}{\mathbf{x}}_{0i}& ={\lambda }_{0i}{\mathbf{M}}_0\sum _{j=1}^N{\epsilon }_{ij}{\mathbf{x}}_{0j}+{\lambda }_{0i}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{M}}_0{\mathbf{x}}_{0i}& & (5)\\ \sum _{j=1}^N{\epsilon }_{ij}{\mathbf{K}}_0{\mathbf{x}}_{0j}+\delta \mathbf{K}{\mathbf{x}}_{0i}& ={\lambda }_{0i}{\mathbf{M}}_0\sum _{j=1}^N{\epsilon }_{ij}{\mathbf{x}}_{0j}+{\lambda }_{0i}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{M}}_0{\mathbf{x}}_{0i}& & \\ (\text{applying }{\mathbf{K}}_0\text{ to the sum})\\ \sum _{j=1}^N{\epsilon }_{ij}{\lambda }_{0j}{\mathbf{M}}_0{\mathbf{x}}_{0j}+\delta \mathbf{K}{\mathbf{x}}_{0i}& ={\lambda }_{0i}{\mathbf{M}}_0\sum _{j=1}^N{\epsilon }_{ij}{\mathbf{x}}_{0j}+{\lambda }_{0i}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{M}}_0{\mathbf{x}}_{0i}& & (\text{using Eq. }(1))\end{array}}$

Because the eigenvectors are M₀-orthogonal when M₀ is positive definite, we can remove the summations by left-multiplying by ${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}}$:

${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}{\epsilon }_{ii}{\lambda }_{0i}{\mathbf{M}}_0{\mathbf{x}}_{0i}+{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\epsilon }_{ii}{\mathbf{x}}_{0i}+{\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}.}$

By use of equation (1) again:

${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{K}}_0{\epsilon }_{ii}{\mathbf{x}}_{0i}+{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\epsilon }_{ii}{\mathbf{x}}_{0i}+{\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}.(6)}$

The two terms containing ε_ii are equal because left-multiplying (1) by ${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}}$ gives

${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{K}}_0{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}.}$

Canceling those terms in (6) leaves

${\displaystyle {\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{K}{\mathbf{x}}_{0i}={\lambda }_{0i}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}+\delta {\lambda }_i{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}.}$

Rearranging gives

${\displaystyle \delta {\lambda }_i=\frac{{\mathbf{x}}_{0i}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}}{{\mathbf{x}}_{0i}^{\mathrm{\top }}{\mathbf{M}}_0{\mathbf{x}}_{0i}}}$

But by (2), this denominator is equal to 1. Thus

${\displaystyle \delta {\lambda }_i={\mathbf{x}}_{0i}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}.}$

Then, as ${\displaystyle {\lambda }_i\ne {\lambda }_k}$ for ${\displaystyle i\ne k}$ (assumption simple eigenvalues) by left-multiplying equation (5) by ${\displaystyle {\mathbf{x}}_{0k}^{\mathrm{\top }}}$:

${\displaystyle {\epsilon }_{ik}=\frac{{\mathbf{x}}_{0k}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}}{{\lambda }_{0i}-{\lambda }_{0k}},i\ne k.}$

Or by changing the name of the indices:

${\displaystyle {\epsilon }_{ij}=\frac{{\mathbf{x}}_{0j}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}}{{\lambda }_{0i}-{\lambda }_{0j}},i\ne j.}$

To find ε_ii, use the fact that:

${\displaystyle {\mathbf{x}}_i^{\mathrm{\top }}\mathbf{M}{\mathbf{x}}_i=1}$

implies:

${\displaystyle {\epsilon }_{ii}=-\frac{1}{2}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}.}$

Summary of the first order perturbation result

In the case where all the matrices are Hermitian positive definite and all the eigenvalues are distinct,

${\displaystyle \begin{array}{cc}{\lambda }_i& ={\lambda }_{0i}+{\mathbf{x}}_{0i}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}\\ {\mathbf{x}}_i& ={\mathbf{x}}_{0i}\left(1-\frac{1}{2}{\mathbf{x}}_{0i}^{\mathrm{\top }}\delta \mathbf{M}{\mathbf{x}}_{0i}\right)+\sum _{\genfrac{}{}{0ex}{}{j=1}{j\ne i}}^N\frac{{\mathbf{x}}_{0j}^{\mathrm{\top }}\left(\delta \mathbf{K}-{\lambda }_{0i}\delta \mathbf{M}\right){\mathbf{x}}_{0i}}{{\lambda }_{0i}-{\lambda }_{0j}}{\mathbf{x}}_{0j}\end{array}}$

for infinitesimal ${\displaystyle \delta \mathbf{K}}$ and ${\displaystyle \delta \mathbf{M}}$ (the higher order terms in (3) being neglected).

So far, we have not proved that these higher order terms may be neglected. This point may be derived using the implicit function theorem; in next section, we summarize the use of this theorem in order to obtain a first order expansion.

Theoretical derivation

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

<span class="mwe-math-element mwe-math-element-inline" typeof="mw:Extension/math" about="#mwt232" data-mw='{"name":"math","attrs":{},"body":{"extsrc":" J_{f,y}(x,g(x)) J_{g,x}(x)+J_{f,x}(x,g(x))=0 \\quad (6) "}}' id="mwAhw"><span class="mwe-math-mathml-a11y" style="display: none;"><math xmlns="http://www.w3.org/1998/Math/MathML" class="mathjax_ignore" alttext="{\displaystyle J_{f,y}(x,g(x))J_{g,x}(x)+J_{f,x}(x,g(x))=0\quad (6)}"> <semantics> <mrow class="MJX-TeXAtom-ORD"> <mstyle displaystyle="true" scriptlevel="0"> <msub> <mi>J</mi> <mrow class="MJX-TeXAtom-ORD"> <mi>f</mi> <mo>,</mo> <mi>y</mi> </mrow> </msub> <mo stretchy="false">(</mo> <mi>x</mi> <mo>,</mo> <mi>g</mi> <mo stretchy="false">(</mo> <mi>x</mi> <mo stretchy="false">)</mo> <mo stretchy="false">)</mo> <msub> <mi>J</mi> <mrow class="MJX-TeXAtom-ORD"> <mi>g</mi> <mo>,</mo> <mi>x</mi> </mrow> </msub> <mo stretchy="false">(</mo> <mi>x</mi> <mo stretchy="false">)</mo> <mo>+</mo> <msub> <mi>J</mi> <mrow class="MJX-TeXAtom-ORD"> <mi>f</mi> <mo>,</mo> <mi>x</mi> </mrow> </msub> <mo stretchy="false">(</mo> <mi>x</mi> <mo>,</mo> <mi>g</mi> <mo stretchy="false">(</mo> <mi>x</mi> <mo stretchy="false">)</mo> <mo stretchy="false">)</mo> <mo>=</mo> <mn>0</mn> <mspace width="1em"></mspace> <mo stretchy="false">(</mo> <mn>6</mn> <mo stretchy="false">)</mo> </mstyle> </mrow> <annotation encoding="application/x-tex">{\displaystyle J_{f,y}(x,g(x))J_{g,x}(x)+J_{f,x}(x,g(x))=0\quad (6)}</annotation> </semantics> </math></span><img src="https://wikimedia.org/api/rest_v1/media/math/render/svg/9e86fa6289ecb94fde259da91675bc60ea0b11f0" class="mwe-math-fallback-image-inline mw-invert skin-invert" aria-hidden="true" style="vertical-align: -1.005ex; width:43.625ex; height:3.009ex;" alt="{\displaystyle J_{f,y}(x,g(x))J_{g,x}(x)+J_{f,x}(x,g(x))=0\quad (6)}"/></span>.

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)=(\genfrac{}{}{0ex}{}{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)=(\genfrac{}{}{0ex}{}{-M{x}_i}{0})\delta \lambda +(\genfrac{}{}{0ex}{}{K-\lambda M}{2x_i^{T}M})\delta x_i}$. Indeed, the solution of