It is well known that the Newton and the Chebyshev methods for nonlinear systems require the solving of a linear system at each iteration step. In this note we shall study two modified methods which avoid the solving of the linear systems by using the Schultz method to approximate the inverses of the Fréchet derivatives. At the same time we shall use the particularities of the nonlinear systems arising from eigenproblems, since the Fréchet derivatives of order higher than two are the null multilinear operators. Some numerical examples will be provided in the end of this note.

Denote $V={𝕂}^n$ and let $A=\left(a_{ij}\right)\in {𝕂}^{n\times n}$, where $𝕂=ℝ$ or $ℂ$. We recall that the scalar $\lambda \in 𝕂$ is an eigenvalue of $A$ if there exists $v\in V$, $v\ne 0$ such that

The vector $v$ is called an eigenvector corresponding to the eigenvalue $\lambda$. Since for an eigenvalue $\lambda$ the eigenpair $(v,\lambda )$ is not uniquely determined, it is necessary to impose a supplementary condition. Different Newton-type methods were studied in the papers [1]-[4], [6], [9]-[11], [16]-[18], [20], [21]. It is worth mentioning that the Rayleigh quotient method is equivalent with a certain Newton method (”the scaled Newton method” [10]).

We shall consider a ”norming” function $G:V\to 𝕂,$ $G\left(0\right)\ne 1$ and besides (1) the equation

The function $G$ may be chosen in different ways (see [16] and [3]):

We shall consider in the theoretical results hereafter the choice II.

Let $X=V\times 𝕂(={𝕂}^{n+1})$ and for $x=\left(\genfrac{}{}{0.0pt}{}{v}{\lambda }\right)\in X$ take

where the norm on $V$ is one of the usual norms.

Consider the system

with the mapping $F$ given by

Denoting $v=(x^{\left(1\right)},x^{\left(2\right)},\mathrm{\dots },x^{\left(n\right)})$ and $\lambda =x^{\left(n+1\right)}$ then the system (2) can be written explicitly

It can be easily seen that the Fréchet derivatives of $F$ are given by the following relations:

for all $x_0=\left(\genfrac{}{}{0.0pt}{}{v_0}{{\lambda }_0}\right)$, $h=\left(\genfrac{}{}{0.0pt}{}{u}{\alpha }\right)$, $k=\left(\genfrac{}{}{0.0pt}{}{w}{\beta }\right)\in X$.

Since the Fréchet derivative of order $2$ does not depend on $x_0$ it is obvious that the derivatives of order higher than $2$ an the null multilinear operators, so for any fixed $x_0\in X$

It is easy to verify that when we use the $\max$ norm on $V$ and $G$ is given by choice II

The following result concerning the invertibility of $F^{\prime }$ at a solution hold.

Proof. The corresponding result for the choice I of $G$ was proved by Yamamoto [21], for which

The stated affirmation follows immediately observing that the two matrices differ by a nonzero factor in the last row. $\mathrm{\square }$

We shall study first the convergence of the sequences ${\left(x_k\right)}_{k\ge 0}\subset X$ and $\left({\mathrm{\Gamma }}_k\right)\subset ℒ\left(X\right)\left(={𝕂}^{(n+1)\times (n+1)}\right)$ generated by the following Newton-type process applied to the nonlinear system (2), initially proposed by Ul’m [19] and studied by Diaconu and Pǎvǎloiu [7]:

$x_0\in X$ and ${\mathrm{\Gamma }}_0\in ℒ\left(X\right)$ being given.

We shall need the following preliminary result.

Denoting ${\overline{B}}_r\left(x_0\right)=\{x\in X:\Vert x-x_0\Vert \le r\}$ we can state the following result.

Proof. It can be easily seen that in our hypotheses the derivatives $F^{\prime }\left(x\right)$ are invertible for all $x\in {\overline{B}}_r\left(x_0\right)$.

Using the inequality

and applying the Banach lemma we get

Taking into account b) and c) it follows that

which together with relation (4) imply

i.e. $x_1\in {\overline{B}}_r\left(x_0\right)$.

Denote ${\rho }_1=\frac{100}{81}{\beta }^2\Vert F\left(x_1\right)\Vert$ and ${\delta }_1=\Vert I-F^{\prime }\left(x_1\right){\mathrm{\Gamma }}_1\Vert$. If we take $x=x_1$ in (3) then an elementary reasoning shows that

whence, by lemma 2 it follows that $\max \{{\rho }_1,{\delta }_1\}\le \frac{1}{9}{d}^2$.

It can be easily proved by induction that the following relations hold for $k=0,1,..$

• •

  $x_k\in {\overline{B}}_r\left(x_0\right);$

• •

  ${\delta }_k:=\Vert I-F^{\prime }\left(x_k\right){\mathrm{\Gamma }}_k\Vert \le \frac{1}{9}{d}^{2^k};$

• •

  ${\rho }_k:=\frac{100}{81}{\beta }^2\Vert F\left(x_k\right)\Vert \le \frac{1}{9}{d}^{2^k};$

• •

  $\Vert x_{k+1}-x_k\Vert \le \frac{d^{2^k}}{10\beta }$.

From the above properties it results that the sequence ${\left(x_k\right)}_{k\ge 0}$ is Cauchy and therefore there exists $x^{\ast }\in {\overline{B}}_r\left(x_0\right)$ such that $x^{\ast }=\underset{k\to \mathrm{\infty }}{lim}{x}_k$. The last inequality above implies that for all $m\in$ $\mathbb{N}$

which leads to the first estimation from the enounce.

The convergence of the sequence ${\left({\mathrm{\Gamma }}_k\right)}_{k\ge 0}$ is infered from the inequalities

which lead to the second stated estimation. $\mathrm{\square }$

As we can see, the Newton-type method (4) has the $r$-convergence order at least 2. The conditions from the above theorem assure that the eigenvalue ${\lambda }^{\ast }={lim}_{k\to \mathrm{\infty }}{x}_k^{\left(n+1\right)}$ is simple, according to lemma 1.

We shall consider now the following sequences given by the Chebyshev-type method, initially proposed by Diaconu [8]

where $x_0\in X$ and $B_0\in ℒ\left(X\right)$ are given.

For the study of the above method we need the following auxilliary result.

As for the previous method, we shall consider the elements $x_0\in X$ and the ball ${\overline{B}}_r\left(x_0\right)$.

Proof. From a), b) and the Banach lemma it easily follows that for any $x\in {\overline{B}}_r\left(x_0\right)$ the Jacobian of $F$is invertible and

For the norms of $B_0$ and $C_0$, taking into account the hypotheses, we get

and

so $\max \{\Vert B_0\Vert ,\Vert C_0\Vert \}\le a$.

From (6) we have that

whence, taking into account d) it follows that $x_1\in {\overline{B}}_r\left(x_0\right)$.

Further, by the identity (3) and by (6) one obtains

whence

Denoting ${\rho }_1=a^2\Vert F\left(x_1\right)\Vert$ and taking into account the inequality

it follows

From the third relation of (6) we get

i.e.,

By lemma 4 the above inequalities imply that $\max \{{\rho }_1,{\delta }_1\}$ $\le \frac{1}{7}{d}^3$.

Assume now that the following properties hold:

• •

  $x_0,x_1,..,x_k\in {\overline{B}}_r\left(x_0\right);$

• •

  ${\rho }_i:=a^2\Vert F\left(x_i\right)\Vert \le \frac{1}{7}{d}^{3^i}$ and ${\delta }_i:=\parallel I-F^{\prime }\left(x_i\right)B_i\parallel \le \frac{1}{7}{d}^{3^i},i=0,..,k$.

It easily follows that $\max \{\Vert B_k\Vert ,\Vert C_k\Vert \}\le a$ and

From the above formula it follows that $x_{k+1}\in {\overline{B}}_r\left(x_0\right)$:

Denoting ${\rho }_{k+1}=a^2\Vert F\left(x_{k+1}\right)\Vert$ and ${\delta }_{k+1}=\Vert I-F^{\prime }\left(x_{k+1}\right)B_{k+1}\Vert$, the following relations are obtained in the same manner as for ${\rho }_1$ and ${\delta }_1$:

whence, by lemma 4, we get that $\max \{{\rho }_{k+1},{\delta }_{k+1}\}\le \frac{1}{7}{d}^{3^{k+1}}$ and the induction is proved.

We will show now that ${\left(x_k\right)}_{k\ge 0}$ is a Cauchy sequence. Indeed,

for all $k,m\in$ $\mathbb{N}$, which implies that ${\left(x_k\right)}_{k\ge 0}$ converges. Denoting $x^{\ast }=\underset{k\to \mathrm{\infty }}{lim}{x}_k$ we obtain

The convergence of ${\left(B_k\right)}_{k\ge 0}$ is obtained from the third relation of (6):

Denoting $B^{\ast }=lim{B}_k$ it easily follows that $B^{\ast }=F^{\prime }{\left(x^{\ast }\right)}^{-1}$ and that

The proof is completed. $\mathrm{\square }$

We shall consider two test matrices¹¹1These matrices are available from MatrixMarket at the following address: http://math.nist.gov/MatrixMarket/. in order to study the behavior of the considered methods. The programs were written in Matlab²²2MATLAB is a registered trademark of the MathWorks, Inc. and were run on a PC.

Pores1 matrix. This matrix arise from oil reservoir simulation. It is real, unsymmetric, of dimension 30 and has 20 real eigenvalues. We have chosen to study the largest eigenvalue ${\lambda }^{\ast }=-1.8363e+1$. The initial approximation was taken ${\lambda }_0={\lambda }^{\ast }+0.5$; for the initial vector $v_0$ we perturbed the solution $v^{\ast }$ (computed by Matlab and then properly scaled to fulfill the norming equation) with random vectors having the components uniformly distributed on (-$\epsilon$,$\epsilon$), $\epsilon =0.2$. The following results are typical for the runs made (we have considered here the same vector $\epsilon$ for the four initial approximations).

Fidap002 matrix. This real symmetric matrix of dimension $n=441$ arise from finite element modeling. Its eigenvalues are all simple and range from $-7\cdot {10}^8$ to $3\cdot {10}^6$. We have chosen to study the smallest eigenvalue, which is well separated. The initial approximation was taken ${\lambda }_0={\lambda }^{\ast }+{10}^3=-6.9996\cdot {10}^8+1000$; for the initial vector $v_0$ we perturbed the solution $v^{\ast }$ with random vectors having the components uniformly distributed on (-$\epsilon$,$\epsilon$), $\epsilon =0.1$. The following results are typical for the runs made (we have considered a common vector $\epsilon$).

If $\epsilon$ was increased to $0.15$ then the Newton-type method with choice I has not converged for some of the initial approximations, even if we took ${\lambda }^{\ast }$ for ${\lambda }_0$. The Newton-type method with choice II and the Chebyshev-type method converged as above. The explanation seem to reside in the fact that the eigenvector $v^{\ast }$ has a larger norm with choice II than with choice I, the relative error of $v^{\ast }+\epsilon$ with the second choice being much smaller than of $v^{\ast }+\epsilon$ with the first choice.