# A Globally Convergent Inexact Newton-Like Cayley Transform Method for Inverse Eigenvalue Problems

Authors: Yonghui Ling; Xiubin Xu
Journal: Journal of Applied Mathematics
 (2013)

---

## Abstract

We propose an inexact Newton method for solving inverse eigenvalue problems (IEP). This method is globalized by employing the classical backtracking techniques. A global convergence analysis of this method is provided and the R-order convergence property is proved under some mild assumptions. Numerical examples demonstrate that the proposed method is very effective in solving the IEP with distinct eigenvalues.


---

## Body

## 1. Introduction

In the present paper, we consider inverse eigenvalue problems (IEP) which are defined as follows. Letc=(c1,c2,…,cn)⊤∈ℝnand let {Ai}i=0nbe a sequence of real symmetricn×nmatrices. Define
(1)A(c)=A0+∑i=1n‍ciAi
and denote its eigenvalues by{λi(c)}i=1nwith the increasing orderλ1(c)≤λ2(c)≤⋯≤λn(c). Givennreal numbers{λi*}i=1nwhich are arranged in increasing orderλ1*≤λ2*≤⋯≤λn*, the IEP is to find a vectorc*∈ℝnsuch that
(2)λi(c*)=λi*,i=1,2,…,n.
Such vectorc*is called a solution of the IEP. This type of inverse problem arises in a variety of applications, for instance, the inverse Toeplitz eigenvalue problem [1, 2], inverse Sturm-Liouville’s problem, inverse vibrating string problem, and the pole assignment problem; see [3–5] and the references therein for more details on these applications.Definef:ℝn→ℝnby
(3)f(c)=(λ1(c)-λ1*,λ2(c)-λ2*,…,λn(c)-λn*)⊤.
Then, solving IEP (2) is equivalent to solving the nonlinear equationf(c)=0onℝn. It is clear thatc*is a solution of the IEP if and only ifc*is a solution of the equationf(c)=0. Based on this equivalence, Newton’s method can be applied to the IEP, and it converges quadratically [6]. As it is known, each iteration of Newton’s method involves solving a complete eigenproblem for the matrixA(c). To overcome this drawback, different Newton-like methods have been proposed and studied [7, 8]. To alleviate the over-solving problem, Bai et al. presented in [9] an inexact Cayley transform method for solving the nonlinear systemf(c)=0. To avoid the computation of the approximate Jacobian equations, Shen and Li proposed in [10, 11] Ulm-like methods for solving the IEPs. However, all these numerical methods for solving the IEP converge only locally.In this paper, we study the numerical methods with global convergence property for solving the IEP. Since the IEP is a nonlinear equation, we review some classical work on solving the general nonlinear equationf(c)=0. Among the inexact Newton-type methods where a line search procedure is used, an inexact Newton backtracking method was proposed in [12]. It performed backtracking along the inexact Newton step, and computational results on a large set of test problems have shown its robustness and efficiency [13, 14].The purpose of the present paper is motivated by the inexact Newton backtracking method proposed in [12], to propose an inexact Newton-type method which combines with the Cayley transform method for solving the IEP. In the backtracking procedure, we use the Rayleigh quotient instead of the classical merit function and therefore reduce the computational cost. Under the classical assumption, which is also used in [7, 9, 10], that the given eigenvalues are distinct and the Jacobian matrix J(c*) is invertible, we show that this method is globally convergent. Some numerical examples are reported to illustrate the effectiveness of the proposed method with distinct eigenvalues.The paper is organized as follows. In Section2, a global inexact Newton-type algorithm is proposed. The global convergence analysis is given in Section 3. And finally in Section 4, some numerical examples are given to confirm the numerical effectiveness and the good performance of our algorithm.


## 2. A Globally Inexact Newton-Like Cayley Transform Method

In this section, we present our algorithm. Letℝn×ndenote the set of all realn×nmatrices. Let ∥·∥ and ∥·∥F denote the 2-norm and the Frobenius norm inℝn, respectively. The induced 2-norm inℝn×nis also denoted by ∥·∥; that is,
(4)∥A∥∶=supx∈ℝn,x≠0∥Ax∥∥x∥,A∈ℝn×n.
Then, we have∥A∥≤∥A∥Ffor any A∈ℝn×n. For x∈ℝn and a positive numberr, we useB(x,r)to stand for the open ball with radius r and center x. Letλ1(c)≤λ2(c)≤⋯≤λn(c)be the eigenvalues of matrixA(c) and let{qi(c)}i=1nbe the normalized eigenvectors corresponding to{λi(c)}i=1n. DefineJ(c)=([J(c)]ij)by
(5)[J(c)]ij∶=qi(c)⊤Ajqi(c),1≤i,j≤n.
Let{λi*}i=1nbe given withλ1*≤λ2*≤⋯≤λn*and writeλ*=(λ1*,λ2*,…,λn*)⊤.The Cayley transform method for computing approximately the eigenproblem of the matrixA(c)was proposed in [6] and was applied in [9, 10]. We now recall this method and then apply it to our algorithm. Suppose thatc*is a solution of the IEP. Then, there exists an orthogonal matrixQ*such that(6)Q*⊤A(c*)Q*=diag{λ1*,λ2*,…,λn*}≜Λ*.
Assume thatc, ρ, andPare the current approximations of c*, λ*, and Q*, respectively. Define eY:=P⊤Q*, where Y is a skew-symmetric matrix. Then, (6) can be rewritten as
(7)P⊤A(c*)P=eYΛ*e-YP⊤A(c*)P=(I+Y+12Y2+⋯)Λ*(I-Y+12Y2+⋯)P⊤A(c*)P=Λ*+YΛ*-Λ*Y+O(∥Y∥2).
Based on (7), we define the new approximationcnewofc*by neglecting the second-order terms inY:
(8)P⊤A(cnew)P=Λ*+YΛ*-Λ*Y.
By equating the diagonal elements in (8), we have
(9)λi*=pi(c)⊤A(cnew)pi(c),i=1,2,…,n,
where{pi(c)}i=1nare the column vectors ofP. Thus, once we getcnewby solving the Jacobian equation, we can obtainYby equating the off-diagonal elements in (8); that is,
(10)[Y]ij∶=pi(c)⊤A(cnew)pj(c)λj*-λi*,1≤i≠j≤n.In order to update the new approximationPnewofQ*, we construct an orthogonal matrixUusing Cayley's transform
(11)U∶=(I+12Y)(I-12Y)-1,
and setPnew=PU; that is, we can obtainPnewby solving
(12)(I+12Y)Pnew=(I-12Y)P.
Finally, the new approximations of eigenvalues can be obtained by
(13)ρinew∶=pi(cnew)⊤A(cnew)pi(cnew),i=1,2,…,n,
where{pi(cnew)}i=1nare the column vectors ofPnew.Note that (12) can be computed as follows. ComputeH:=(I-(1/2)Y)P⊤and lethibe theith column ofHat first. Then, solvewi, i=1,2,…,n, iteratively from thenlinear systems:
(14)(I+12Y)wi=hi,i=1,2,…,n.
Finally, setPnew:=[w1,…,wn]⊤. SincePis an orthogonal matrix andYis skew-symmetric matrix, we see thatPnewmust be orthogonal. To maintain the orthogonality ofPnew, (14) cannot be solved inexactly. One could expect that it requires only few iterations to solve each system of (14). This is due to the fact that as{ck}converges toc*,∥Yk∥converges to zero; see [6,(3.64)]. Consequently, the coefficient matrix on the left-hand side of (14) approaches the identity matrix in the limit.For solving the general nonlinear equationf(x)=0, linesearch techniques [15] are often used to enlarge the convergence basin of a locally convergent method. They are based on a globally convergent method for a problem of the formminx∈ℝnM(x), whereMis an appropriately chosen merit function whose global minimum is a zero off. In these cases, for a given directions∈ℝn, we have the iteration form xk+1=xk+αs, whereα∈(0,1]is such thatM(xk+αs)<M(xk). The existence of such anαis ensured if there exists anα0>0such thatM(xk+αs)<M(xk)for allα<α0.In typical linesearch strategies, the step lengthαis chosen by using so-called backtracking approach. Among the backtracking method, inexact Newton backtracking methods (INB) [12] is a globally convergent process where thekth iteration of an inexact Newton method is embedded in a backtracking strategy. The merit functionMof INB usually usedM:=∥f∥; see, for example, [12–14, 16]. Thanks to (3), this will involve computingλi(ck)ofA(ck)which are costly to compute. Our intention here is to replace them by the Rayleigh quotient (see (13)). In Section 3 we will show that this replacement retain superlinear and global convergence.The details of our algorithm for solving the IEP are specified as Algorithm1. In Step 7, the following sufficient decrease in the merit function∥ρ(c)-λ*∥based on the Rayleigh quotient is provided:
(15)∥ρ(ck+Δck)-λ*∥≤(1-ξ(1-ηk))∥ρ(ck)-λ*∥.
The while loop in Step 7 is also called backtracking loop below. Note that if we set
(16)η¯k1∶=∥ρ(ck)-λ*∥∥λ*∥,η¯k2∶=∥ρ(ck)-λ*∥∥ρ(ck-1)-λ*∥,
when ck is away from an accumulation point c*, thenη¯k1>1orη¯k2>1may be possible. Accordingly, the forcing termη¯kgiven in (20) guarantees thatη¯k∈[0,1)for eachk.


## 3. Convergence Analysis

In this section, we analyze the global behavior of Algorithm1. We will show that if the given eigenvalues are distinct and if there exists an accumulation pointc*of{ck}such that the Jacobian matrixJ(c*)is invertible, then the iterations are guaranteed to remain near c* and ρ(c*)-λ*=0, ck→c*as k→∞. Furthermore, forksufficiently large, we have the equalityck+1=ck+Δc-k. Thus, we obtain that the ultimate rate of convergence isβwhich depends on the choices of theη¯kgiven in (20).Algorithm 1 (Inexact Newton-Like Backtracking Cayley Transform Method for IEP).
For anyc0∈ℝn, ηmax∈[0,1), ξ∈(0,1), 0<θmin<θmax<1. Compute the orthonormal eigenvectors {qi(c0)}i=1n and the eigenvalues {λi(c0)}i=1n of A(c0). Let
(17)P0:=[p1(c0),p2(c0),…,pn(c0)]P0=[q1(c0),q2(c0),…,qn(c0)],ρ(c0):=[ρ1(c0),ρ2(c0),…,ρn(c0)]⊤(c0)ρ=[λ1(c0),λ2(c0),…,λn(c0)]⊤.
Fork=0,1,2,…until convergence do the following.
Step
1. Form[Jk]ij=pi(ck)⊤Ajpi(ck)for1≤i,j≤n.
Step
2. SolveΔc-kinexactly from the approximate Jacobian equation:
(18)JkΔc-k+ρ(ck)-λ*=0
such that
(19)∥JkΔc-k+ρ(ck)-λ*∥≤η¯k∥ρ(ck)-λ*∥,
whereη¯0∈(0,1)and
(20)η¯k:=min{∥ρ(ck)-λ*∥β∥λ*∥β,∥ρ(ck)-λ*∥β∥ρ(ck-1)-λ*∥β,ηmax},mmmmmmmmmmmmmm..β∈(1,2],k=1,2,….
Step
3. SetΔck=Δc-kandηk=η¯k.
Step
4. Form the skew-symmetric matrixY¯kby (10) withcnew=ck+Δck.
Step
5. Compute matrix P¯k:=[p1(ck+Δck),…,pn(ck+Δck)] by solving (12).
Step
6. Compute ρ(ck+Δck)=(ρ1(ck+Δck),…,ρn(ck+Δck))⊤ by (13).
Step
7. While ∥ρ(ck+Δck)-λ*∥>(1-ξ(1-ηk))∥ρ(ck)-λ*∥ do the following:
(21)chooseθ∈[θmin,θmax],thenupdateΔck⟵θΔck,ηk⟵1-θ(1-ηk).
Step
8. Setck+1=ck+Δck. As Steps 4–6, compute, respectively, the new approximations Pk+1:=[p1(ck+1),…,pn(ck+1)]andρ(ck+1):=(ρ1(ck+1),…,ρn(ck+1))⊤.It is worth noting that if{ck}has no accumulation point, {ck}has one or more accumulation points and the Jacobian matrix is singular at each of them, or the vectorΔc¯kcomputed by solving the Jacobian equation (18) is such thatΔc-k=0, then our algorithm fails.It is clear that ifρ(ck)-λ*→0ask→∞andc*is an accumulation point of{ck}, then ρ(c*)-λ*=0. Let{Pk}be generated by Algorithm 1 (see Step 7) and defineEk:=Pk-Q*for eachk=0,1,…. The following lemma is taken from [17, Lemma 2].Lemma 2 (see [17]).
For anyc, c¯∈ℝn, one has
(22)∥A(c)-A(c¯)∥≤L∥c-c¯∥,
whereL:=(∑i=1n‍∥Ai∥2)1/2.Based on Lemma2, the following lemma is a straightforward application of [9, Lemma 4].Lemma 3.
Assume thatc*is an accumulation point of{ck}. Let the given eigenvalues{λi*}i=1nbe distinct. Then
(23)∥ρ(ck)-λ*∥≤nL∥ck-c*∥+μ2∥Ek∥,
for anyk∈ℕ, whereμ2:=2n·∥A(c*)∥.As shown in [18, Theorem 2.3], in the case when the given eigenvalues{λi*}i=1nare distinct, the eigenvalues ofA(c)are distinct too for any pointcin some neighborhood ofc*. It follows that the functionf(·)defined in (3) is analytic in the same neighborhood. However, ifcis not near the solutionc*, the analyticity of the functionf(·)cannot be guaranteed.For any symmetric matrixX∈ℝn×n, setΛ(X):=diag(λ1(X),…,λn(X)), whereλi(X), i=1,2,…,n are the eigenvalues ofX. As proved by D. Sun and J. Sun [19, Theorem 4.7], Λ(·) is a strongly semismooth function. Based on this result, for anyc∈ℝnandɛ>0, there existsδ>0sufficiently small such that
(24)∥ρ(c¯)-ρ(c)-J(c)(c¯-c)∥≤ɛ∥c¯-c∥2,foranyc¯∈B(c,δ).The following lemma says that the backtracking loop in Step 7 of Algorithm1 terminates after a finite number of steps.Lemma 4.
There existsηkmin∈[0,1)such that, for anyηk∈[ηkmin,1), there isΔcksatisfying
(25)∥ρ(ck+Δck)-λ*∥≤[1-ξ(1-ηk)]∥ρ(ck)-λ*∥.Proof.
By using the strong semismoothness of all eigenvalues of a real symmetric matrix [19], for any givenɛ>0, there existsδ>0sufficiently small such that
(26)∥ρ(ck+Δc)-ρ(ck)-JkΔc∥≤ɛ∥Δc∥
whenever∥Δc∥≤δ. Choose
(27)ɛ=(1-ξ)(1-η¯k)∥ρ(ck)-λ*∥∥Δc¯k∥,
and set
(28)ηkmin∶=max{η¯k,1-(1-η¯k)δ∥Δc¯k∥}.
For anyηk∈[ηkmin,1), letΔck:=((1-ηk)/(1-η¯k))Δc¯k. Then, by the definition ofηkmingiven in (28), we have
(29)∥Δck∥=1-ηk1-η¯k∥Δc¯k∥≤1-ηkmin1-η¯k∥Δc¯k∥≤(1-η¯k)δ∥Δc¯k∥·11-η¯k∥Δc¯k∥=δ.
On the other hand, by (19), one gets that
(30)∥ρ(ck)-λ*+JkΔck∥≤ηk-η¯k1-η¯k∥ρ(ck)-λ*∥+1-ηk1-η¯k∥ρ(ck)-λ*+JkΔc¯k∥≤ηk-η¯k1-η¯k∥ρ(ck)-λ*∥+1-ηk1-η¯k·η¯k∥ρ(ck)-λ*∥=ηk∥ρ(ck)-λ*∥.
This together with (26) and (27) yields that
(31)∥ρ(ck+Δck)-λ*∥≤∥ρ(ck+Δck)-ρ(ck)-JkΔck∥+∥ρ(ck)-λ*+JkΔck∥≤ɛ∥Δck∥+ηk∥ρ(ck)-λ*∥=ɛ1-ηk1-η¯k∥Δc¯k∥+ηk∥ρ(ck)-λ*∥=(1-ξ)(1-ηk)∥ρ(ck)-λ*∥+ηk∥ρ(ck)-λ*∥=[1-ξ(1-ηk)]∥ρ(ck)-λ*∥,
and the proof is completed.Next, we give sufficient conditions for Algorithm1 not to break down in the backtracking loop in Step7.Lemma 5.
Ifρ(ck)-λ*≠0and there existsΓ>0for which
(32)∥Δck∥≤Γ(1-ηk)∥ρ(ck)-λ*∥,
then the backtracking loop terminates.Proof.
For constantΓin (32) and the givenξ∈(0,1), choosingɛ=(1-ξ)/Γ, there existsδ>0sufficiently small such that (26) holds whenever∥Δc∥≤δ. We chooseηk∈[η¯k,1]satisfying
(33)1-ηk<δΓ∥ρ(ck)-λ*∥.
It follows that
(34)∥Δck∥≤Γ(1-ηk)∥ρ(ck)-λ*∥<δ,
which givesck+Δck∈B(ck,δ). Thus, we have
(35)∥ρ(c+Δck)-ρ(ck)-JkΔck∥≤1-ξΓ∥Δck∥.
This together with (19) gives that
(36)∥ρ(ck+Δck)-λ*∥≤∥ρ(ck)-λ*+JkΔck∥+∥ρ(ck+Δck)-ρ(ck)-JkΔck∥≤ηk∥ρ(ck)-λ*∥+(1-ξ)(1-ηk)∥ρ(ck)-λ*∥=(1-ξ(1-ηk))∥ρ(ck)-λ*∥.
It follows that the backtracking loop terminates. This completes the proof.The next lemma gives condition under which (32) is satisfied.Lemma 6.
Assume that the Jacobian matrixJkis invertible and setMk:=∥Jk-1∥. Then there existsΓsuch that (32) holds.Proof.
By using condition (19), one has that
(37)∥Δck∥≤∥Jk-1∥∥JkΔck∥≤Mk(∥ρ(ck)-λ*∥+∥ρ(ck)-λ*+JkΔck∥)≤Mk(1+ηk)∥ρ(ck)-λ*∥=Mk1+ηk1-ηk(1-ηk)∥ρ(ck)-λ*∥≤Mk1+ηmax1-ηmax(1-ηk)∥ρ(ck)-λ*∥.
We finish the proof by takingΓ:=Mk(1+ηmax)/(1-ηmax).Lemmas5 and 6 yield the result below.Corollary 7.
Assume thatρ(ck)-λ*≠0andJkare invertible. SetMk:=∥Jk-1∥and Γ:=Mk(1+ηmax)/ (1-ηmax). Then, the backtracking loop in Step 7 of Algorithm 1 terminates with
(38)1-ηk≥min{1-η¯k,δθminΓ∥ρ(ck)-λ*∥},
for anyδ>0small enough such that (26) holds whenever∥Δc∥≤δ.Proof.
Suppose thatηkis the final value determined by the while-loop. Ifηk=η¯k, then (38) is trivial. Assume thatηk≠η¯k; that is, the body of the while-loop has been executed at least once. Denoting the penultimate value byηk-, then it follows from (33) that1-ηk-≥δ/Γ∥ρ(ck)-λ*∥. Thus, we have
(39)1-ηk=θ(1-ηk-)≥θmin(1-ηk-)≥δθminΓ∥ρ(ck)-λ*∥.
This completes the proof.Lemma 8.
Assume thatc*is an accumulation point of{ck}such that there exists a constantΓindependent ofkfor which
(40)∥Δck∥≤Γ(1-ηk)∥ρ(ck)-λ*∥,
wheneverckis sufficiently nearc*andkis sufficient large. Thenck→c*ask→∞.Proof.
Suppose thatck↛c*. Then, there exists δ>0 sufficiently small such that there are infinitely many k for which ck∉B(c*,δ) and (40) hold wheneverck∈B(c*,δ)forksufficiently large.
Sincec*is an accumulation point of{ck}, there exists subsequence {cki}⊂{ck} such that cki∈B(c*,δ/i) for i sufficiently large. Choose ℓi>0 satisfyingki+ℓi<ki+1andcki+ℓi∉B(c*,δ). It follows that
(41)∥ρ(cki+j-1)-λ*∥≤1Γ(1-ηki+j-1)(∥ρ(cki+j-1)-λ*∥-∥ρ(cki+j)-λ*∥),mmmmmmmmmmmmmmmmmmmmm1m,j=1,…,ℓi.
Then, by (40), we have, forisufficiently large,
(42)δ2≤∥cki+ℓi-cki∥≤∑k=kiki+ℓi-1‍∥Δck∥≤∑k=kiki+ℓi-1‍Γ(1-ηk)∥ρ(ck)-λ*∥≤Γt(∥ρ(cki)-λ*∥-∥ρ(cki+1)-λ*∥).
Note thatcki→c*asi→∞. It follows that∥ρ(cki)-λ*∥-∥ρ(cki+1)-λ*∥→0, which is a contradiction. Therefore,ck→c*ask→∞.Lemma 9.
Assume thatc* is an accumulation point of {ck} such that J(c*) is invertible. Set M=∥J(c*)-1∥ and let δ1 be such that 0<δ1≤1/2μ1M, where μ1:=2n·max1≤j≤n∥Aj∥. Suppose that ∥Ek∥≤δ1 for k sufficiently large. Then, for all k sufficiently large,
(43)∥Jk-1∥≤∥J(c*)-1∥1-μ1δ1∥J(c*)-1∥≤2M.Proof.
By the definitions of[Jk]ijand[J(c*)]ij, for allksufficiently large,
(44)|[Jk]ij-[J(c*)]ij|≤|pi(ck)⊤Ajpi(ck)-pi(ck)⊤Ajqi(c*)|+|pi(ck)⊤Ajqi(c*)-qi(c*)⊤Ajqi(c*)|≤∥pi(ck)⊤∥∥Aj∥∥pi(ck)-qi(c*)∥+∥pi(ck)⊤-qi(c*)⊤∥∥Aj∥∥qi(c*)∥=2∥Aj∥∥pi(ck)-qi(c*)∥,1≤i,j≤n.
Then, we have, for allksufficiently large,
(45)∥Jk-J(c*)∥≤∥Jk-J(c*)∥F≤2n·max1≤j≤n∥Aj∥·max1≤i≤n∥pi(ck)-qi(c*)∥.
Noting that pi(ck)-qi(c*)is theith column ofEk, then∥pi(ck)-qi(c*)∥≤∥Ek∥fori=1,…,n. So, for allksufficiently large,
(46)∥J(c*)-1∥∥Jk-J(c*)∥≤μ1δ1M≤12<1.
It follows from Banach lemma thatJkis invertible and∥Jk-1∥≤2Mfor allksufficiently large. This completes the proof.Lemmas6, 8, and 9 yield the result below.Corollary 10.
Assume thatc* is an accumulation point of {ck} such that J(c*) is invertible. Set M:=∥J(c*)-1∥ and let δ1 be determined by Lemma 9. Assume that ∥Ek∥≤δ1 for k sufficiently large. Then ck→c* as k→∞.Proof.
By Lemma9, we knowJkis invertible and∥Jk-1∥≤2Mfor allksufficiently large. From the proof of Lemma 6, (40) holds for the constantΓ=2M(1+ηmax)/(1-ηmax)independent ofk. Therefore, ck→c* as k→∞ follows from Lemma 8.Corollary 11.
Assume thatc*is an accumulation point of{ck}such thatJ(c*)is invertible. Then, ρ(ck)-λ*→0ask→∞. Moreover, for allksufficiently large, one hasηk=η¯k.Proof.
IfΔckis computed in the backtracking loop, the backtracking terminates withηksuch that (38) holds. Sinceck→c*by Corollary 10, we haveck∈B(c*,δ)for allksufficiently large. Thus, the series ∑k=0∞‍(1-ηk) is divergent. Then, we have
(47)∏1≤i<k‍[1-ξ(1-ηi)]≤∏1≤i<k‍e-ξ(1-ηi)∏1≤i<k‍[1-ξ(1-ηi)]=e-ξ·∑i=1k-1(1-ηi)⟶0,k⟶∞,
and so
(48)∥ρ(ck)-λ*∥≤[1-ξ(1-ηk-1)]∥ρ(ck-1)-λ*∥≤∥ρ(c0)-λ*∥∏0≤i<k‍[1-ξ(1-ηi)]⟶0,k⟶∞.
This together with (38) yields thatηk=η¯kfor allksufficiently large.Lemma 12.
Assume thatc*is an accumulation point of{ck}such thatJ(c*)is invertible. Then, there existsδ2>0such that
(49)∥ρ(ck)-λ*-Jk(ck-c*)∥≤μ3∥Ek∥F2
wheneverck∈B(c*,δ2)andksufficiently large, whereμ3:=2max1≤i≤n|λi*|.Proof.
SetHk:=Q(c*)⊤Pk-I. Then, we have, for allksufficiently large,
(50)Pk⊤A(c*)Pk=(I+Hk)⊤Λ*(I+Hk)=Λ*+Λ*Hk+Hk⊤Λ*+Hk⊤Λ*Hk.
It follows that, for allksufficiently large,
(51)pi(ck)⊤A(c*)pi(ck)=λi*+2λi*[Hk]ii+∑j=1n‍λj*[Hk]ij2,i=1,2,…,n.
Set
(52)vi:=2λi*[Hk]ii+∑j=1n‍λj*[Hk]ij2,i=1,2,…,n,
and writev=(v1,v2,…,vn)⊤. It follows thatJkc*=λ*-ak+v, where
(53)[ak]i=pi(ck)⊤A0pi(ck),i=1,2,…,n.
In view of thatJkck=ρ(ck)-akfor anyk∈ℕ, one has that, for allksufficiently large,
(54)v=ρ(ck)-λ*-Jk(ck-c*).
Since
(55)I+HkHk⊤+HkHk⊤=(I+Hk)(I+Hk)⊤=Q(c*)⊤PkPk⊤Q(c*)=I,
we have
(56)[Hk]ii=-12∑j=1n‍[Hk]ij2,i=1,2,…,n.
It follows from (52) that, for allksufficiently large,
(57)vi2=(2λi*[Hk]ii+∑j=1n‍λj*[Hk]ij2)2≤2[4(λi*)2[Hk]ii2+(∑j=1n‍λj*[Hk]ij2)2].
Combining (56) with (57), one has that, for allksufficiently large,
(58)∑i=1n‍vi2≤2[∑i=1n‍4(λi*)2[Hk]ii2+∑i=1n‍(∑j=1n‍λj*[Hk]ij2)2]∑i=1n‍vi2=2[∑i=1n‍(λi*)2(∑j=1n‍[Hk]ij2)2+∑i=1n‍(∑j=1n‍λj*[Hk]ij2)2]∑i=1n‍vi2≤2[max1≤i≤n|λi*|2·∑i=1n‍(∑j=1n‍[Hk]ij2)2z∑i=1n‍vi2zzzz+max1≤i≤n|λi*|2·∑i=1n‍(∑j=1n‍[Hk]ij2)2]∑i=1n‍vi2≤4max1≤i≤n|λi*|2·(∑i=1n‍∑j=1n‍[Hk]ij2)2∑i=1n‍vi2=4max1≤i≤n|λi*|2·∥Hk∥F4,
which gives
(59)∥v∥≤2max1≤i≤n|λi*|·∥Hk∥F2.
Note that,
(60)∥Hk∥F=∥Q(c*)⊤Pk-I∥F=∥Pk-Q(c*)∥F=∥Ek∥F.
Therefore, we obtain from (54) that (49) holds for allksufficiently large. This completes the proof.Lemma 13.
Assume thatc*is an accumulation point of{ck}such thatJ(c*)is invertible. Then, there existδ3>0,δ4∈(0,δ1]sufficiently small such that∥Ek∥≤δ4and
(61)∥ck+1-c*∥≤2M(μ3∥Ek∥F2+(2μ3)β∥ρ(ck)-λ*∥β)
wheneverck∈B(c*,δ3)andksufficiently large, whereδ1is determined by Lemma 9.Proof.
Thanks to Corollary11, we haveηk=η¯kandΔck=Δc¯kfor allksufficiently large. Thus, the Jacobian equation (18) is equivalent to, for allksufficiently large,
(62)Jkck+1-λ*+ak=0.
Assume that the residual of this approximate Jacobian equation is defined byrk; that is, for allksufficiently large,
(63)rk=Jkck+1-λ*+ak.
This together withJkc*-λ*+ak=vgives, for allksufficiently large,
(64)Jk(c*-ck+1)=v-rk.
By (19), Lemmas 9 and 12, we obtain, for allksufficiently large,
(65)∥ck+1-c*∥≤∥Jk-1∥(∥v∥+∥rk∥)∥ck+1-c*∥≤2M(μ3∥Ek∥F2+∥ρ(ck)-λ*∥β+1∥λ*∥β).
It follows from Lemma 3 that
(66)∥ρ(ck)-λ*∥≤nL∥ck-c*∥+μ2∥Ek∥,∀k∈ℕ.
Thus, we can chooseδ3>0and0<δ4≤δ1sufficiently small such that
(67)nLδ3+μ2δ4≤1,
whenever∥ck-c*∥≤δ3and∥Ek∥F≤δ4. Therefore, combining this with the definition ofμ3given in Lemma 12, (61) follows.Lemma 14 (see [6]).
There exist two positive numbersδ5 and ω1 such that, for any orthogonal matrix P with ∥P-Q(c*)∥<δ5, the skew-symmetric X defined by eX:=P⊤Q(c*) satisfies ∥X∥≤ω1∥P-Q(c*)∥.Based on Lemma14, by using the similar arguments in the proof of [10, Lemma 5], we can obtain the following lemma. If we write eXk:=Pk⊤Q*, then there exists C>0 such that ∥Yk∥≤C(∥ck+1-c*∥+∥Xk∥).Lemma 15.
Suppose that the given eigenvalues{λi*}i=1n are distinct and the Jacobian matrix J(c*) is invertible. Then, there exist ω2>0 and 0<δ6<min{δ5,1/(1+ω1)C} such that, forksufficiently large, if ∥ck+1-c*∥≤δ6 and ∥Ek∥≤δ6, then
(68)∥Ek+1∥≤ω2(∥ck+1-c*∥+∥Ek∥2),
whereω1is determined by Lemma 14.In order to prove our global convergence result for Algorithm1, we introduce some notations. Letζ1:=2(nL+μ2)/μ3andζ2:=2M(μ3+ζ1β). Set
(69)τ∶=min{1,(2M(μ3+ζ1β))-(1/(β-1))mmmm1m×(ω2(ζ2+1))-(1/(β-1))},(70)δ∶=min{τ,δ2,δ3,δ4,δ6}.
Our main global convergence result is as follows.Theorem 16.
Assume that{ck} is generated by Algorithm 1. Suppose that the given eigenvalues {λi*}i=1n are distinct and c* is an accumulation point of {ck} such that J(c*) is invertible. Then, ρ(ck)-λ*→0 and ck→c* as k→∞. Moreover, the convergence is of R-order β.Proof.
It follows immediately from Corollaries10 and 11 that ρ(ck)-λ*→0 and ck→c* as k→∞. For the τ given in (69), there exists k0 sufficiently large such that ∥ck0-c*∥≤τ and ∥Ek0∥≤τ. Set γ:=δ/τ, where δ and τ are given in (70) and (69), respectively. Then, γ≤1. We will show that, for all k≥k0 sufficiently large,
(71)∥ck-c*∥≤τ·γβk-k0,(72)∥Ek∥≤τ·γβk-k0.
Suppose that (71) and (72) hold for some k=ℓ≥k0. Consider the case k=ℓ+1. Thanks to Lemma 3, we have, for all ℓ≥k0,
(73)∥ρ(cℓ)-λ*∥≤nL∥cℓ-c*∥+μ2∥Eℓ∥≤nLτ·γβℓ-k0+μ2τ·γβℓ-k0=(nL+μ2)τ·γβℓ-k0.
Then, by using Lemma 13, one has that, for all ℓ≥k0,
(74)∥cℓ+1-c*∥≤2M(μ3∥Eℓ∥F2+(2μ3)β∥ρ(cℓ)-λ*∥)≤2Mμ3τ2(γβℓ-k0)2+2M(2μ3)β(nL+μ3)βτβγβℓ-k0+1≤2M[μ3+(2(nL+μ2)μ3)β]τβγβℓ-k0+1≤τγβℓ-k0+1,
where the last inequality follows from the definition ofτin (69). By Lemma 15, we have, for all ℓ≥k0,
(75)∥Eℓ+1∥≤ω2(∥cℓ+1-c*∥+∥Eℓ∥2)∥Eℓ+1∥≤2Mω2(μ3+ζ1β)τβγβℓ-k0+1+ω2τ2(γβℓ-k0)2∥Eℓ+1∥≤ω(2M(μ3+ζ1β)+1)τβ·γβℓ-k0+1∥Eℓ+1∥=ω2(ζ2+1)τβ·γβℓ-k0+1≤τ·γβℓ-k0+1.
Therefore, we conclude that (71) and (72) hold for allk≥k0. Moreover, we see from (71) thatckconverges toc*with R-orderβ. This completes the proof.If{ck}generated by Algorithm 1 converges to a solution at which the Jacobian matrix is invertible, then the ultimate rate of convergence is governed by the choices of theηk(k=0,1,…)as in the local theory of [9].


## 4. Numerical Examples

In this section, we illustrate the effectiveness of Algorithm1 in solving IEP on three examples. The tests were carried out in MATLAB 7.0 running on a PC Intel Pentium P6200 of 2.13 GHz CPU.The given parameters used in our algorithm wereη0=0.5, ηmax=0.9, ξ=10-4, θmin=0.1, and θmax=0.9. In the while loop, we choose θ∈[θmin,θmax] to minimize ∥ρ(ck+θΔck)-λ*∥ if 80 iterations of the backtracking loop fail to produce the sufficient decrease in ∥ρ(c)-λ*∥.Linear systems (14) and (19) are all solved iteratively by the QMR method [20] using the MATLAB qmr function. In order to guarantee the orthogonality of Pnew in (14), this system is solved up to machine precision eps (≈2.2×10-16). The inner loop stopping tolerance for (18) is given by (20). The stopping criterion of the outer iteration in our algorithm is
(76)∥Pk⊤A(ck)Pk-Λ*∥F≤10-10.Example 1.
This is an inverse Toeplitz eigenvalue problem (see [2] for more details on this inverse problem) with distinct eigenvalues. The basis matrices {Ai}i=15 are given as follows:
(77)A1=[1000001000001000001000001],A2=[0100010100010100010100010],A3=[0010000010100010100000100],A4=[0001000001000001000001000],A5=[0000100000000000000010000].
The given real eigenvalues and a solution, respectively, are
(78)λ*=(-5.2361,-1.5876,-0.7639,-0.5555,18.1431)⊤,c*=(2,3,4,5,6)⊤.
In Table 1, we report our numerical results for various starting points:
(79)(a)c0=(1,2,3,4,5)⊤,(b)c0=(1,5,10,15,20)⊤,(c)c0=(11,12,13,14,15)⊤,(d)c0=(21,38,46,63,81)⊤,(e)c0=(101,112,123,134,145)⊤,
where c0, errs, ite. and c* stand for the starting point, the error value of the left-hand side of (76) for the last three iterates of the algorithm, the number of outer iteration, and the accumulation point corresponding to the starting point.Table 1
Numerical results for Example1.
c
0
β
=
1.5
β
=
1.8
β
=
2.0
(a)
errs
1.0659
e
-
02
1.0659
e
-
02
1.0659
e
-
02
1.1323
e
-
05
1.1323
e
-
05
1.1323
e
-
05
8.0567
e
-
12
8.0567
e
-
12
8.0567
e
-
12
ite.
6
6
6
c
*
(
2.0000,3.2926,3.4471,4.9014,6.5529
)
⊤
(b)
errs
3.8846
e
-
03
3.8846
e
-
03
3.8846
e
-
03
1.2491
e
-
06
1.2491
e
-
06
1.2491
e
-
06
7.9990
e
-
14
7.9990
e
-
14
7.9990
e
-
14
ite.
9
9
9
c
*
(
2,3
,
4,5
,
6
)
⊤
(c)
errs
5.0389
e
-
05
5.0389
e
-
05
5.0389
e
-
05
3.5838
e
-
10
3.5838
e
-
10
3.5838
e
-
10
5.0337
e
-
15
5.0337
e
-
15
5.0337
e
-
15
ite.
13
13
13
c
*
(
2.0000,3.2926,3.4471,4.9014,6.5529
)
⊤
(d)
errs
1.3806
e
-
04
1.3806
e
-
04
1.3806
e
-
04
6.4118
e
-
09
6.4118
e
-
09
6.4118
e
-
09
9.8493
e
-
15
9.8493
e
-
15
9.8493
e
-
15
ite.
75
75
75
c
*
(
2.0000,3.2926,3.4471,4.9014,6.5529
)
⊤
(e)
errs
1.2028
e
-
04
7.7624
e
-
03
7.7624
e
-
03
3.8608
e
-
10
6.1243
e
-
06
6.1243
e
-
06
5.3680
e
-
15
2.1634
e
-
12
2.1634
e
-
12
ite.
13
5
5
c
*
(
2,3
,
4,5
,
6
)
⊤Example 2.
This is a Toeplitz-plus-Hankel inverse eigenvalue problem (see [1] for more details on this inverse problem) with distinct eigenvalues. The basis matrices{Ai}i=17are given as follows:
(80)A1=[-1000000010000000100000001000000010000000100000001],A2=[0-100000-101000001010000010100000101000001010000010],A3=[00-100000-201000-10001000100010001000100010000000100],A4=[000-100000-201000-200010-1000001010000000100000001000],A5=[0000-100000-201000-200010-200000-100000001000000010000],A6=[00000-100000-201000-200000-200000-200000-10000000100000],A7=[000000-100000-200000-200000-200000-200000-200000-1000000].
The given real eigenvalues are
(81)λ*=(-35.4513,-13.6805,-9.5675,mm,m-8.5489,8.7666,11.8220,26596)⊤,
andc*=(2,3,4,5,6,7,8)⊤is a solution. In Table 2, we report our numerical results for various starting points:
(82)(a)c0=(1,2,3,4,5,6,7)⊤,(b)c0=(1,3,5,7,9,11,13)⊤,(c)c0=(11,13,15,17,19,21,23)⊤,(d)c0=(50,52,56,58,62,65,68)⊤,(e)c0=(101,102,103,104,106,108,110)⊤.Table 2
Numerical results for Example2.
c
0
β
=
1.5
β
=
1.8
β
=
2.0
(a)
errs
8.2730
e
-
03
8.2730
e
-
03
8.2730
e
-
03
2.2866
e
-
05
2.2866
e
-
05
2.2866
e
-
05
5.5845
e
-
11
5.5845
e
-
11
5.5845
e
-
11
ite.
6
6
6
c
*
(
2,3
,
4,5
,
6,7
,
8
)
⊤
(b)
errs
8.2129
e
-
03
4.0190
e
-
03
1.7235
e
-
04
2.0027
e
-
05
1.5988
e
-
06
3.7839
e
-
09
5.0828
e
-
11
2.8276
e
-
13
1.9386
e
-
14
ite.
7
7
7
c
*
(
2.8703,2.0639,4.9434,4.1053,6.7530,6.3287,8.4794
)
⊤
(c)
errs
3.4941
e
-
03
3.4941
e
-
03
3.4941
e
-
03
2.5608
e
-
06
2.5608
e
-
06
2.5608
e
-
06
6.7550
e
-
13
6.7550
e
-
13
6.7550
e
-
13
ite.
11
11
11
c
*
(
2.8703,2.0639,4.9434,4.1053,6.7530,6.3287,8.4794
)
⊤
(d)
errs
2.9914
e
-
04
2.9914
e
-
04
2.9914
e
-
04
8.2616
e
-
09
8.2616
e
-
09
8.2616
e
-
09
1.8946
e
-
14
1.8946
e
-
14
1.8946
e
-
14
ite.
13
13
13
c
*
(
1.1787
,
-
0.0035,2.0401,1.6976,5.3794,6.2068,8.5273
)
⊤
(e)
errs
2.1008
e
-
04
2.1008
e
-
04
2.1008
e
-
04
1.5703
e
-
08
1.5703
e
-
08
1.5703
e
-
08
2.2331
e
-
14
2.2331
e
-
14
2.2331
e
-
14
ite.
13
13
13
c
*
(
2.8703,2.0639,4.9434,4.1053,6.7530,6.3287,8.4794
)
⊤Example 3 (see [10]).
GivenB=[bij]8×8:=I+VV⊤, where
(83)V=[1-1-3-5-611-2-5-171-1-1518111201-12011130-12.50.20.30.50.62-0.20.30.50.8]8×5,
the basis matrices {Ak}k=18 are defined from B as follows:
(84)Ak=bkkekek⊤+∑j=1k-1‍bkj(ekek⊤+ejek⊤),k=1,2,…,8,
whereekis thekth column of the identity matrixI. The given real eigenvalues and a solution, respectively, are
(85)λ*=(1.0e+04)×(-1.1672,-0.3139,-0.0140,-0.0008,mmmmmmmmm1m0.0003,0.0075,0.0327,1.7372)⊤,
andc*=(2,3,4,5,6,7,8,9)⊤. Table 3 shows the numerical results for the following various starting points:
(86)(a)c0=(1,1,1,1,1,1,1,1)⊤,(b)c0=(2.5,2.5,2.5,2.5,2.5,2.5,2.5,2.5)⊤,(c)c0=(1,2,3,4,5,6,7,8)⊤,(d)c0=(2,4,6,8,10,12,14,16)⊤,
where
(87)c1*=(2,3,4,5,6,7,8,9)⊤,c2*=(-2.7065,-7.0690,8.4130,1.1038,-3.6637,mmm8.8346,5.5078,9.1273)⊤,c3*=(1.2877,2.4677,4.5379,-6.1037,1.7930,mmm1.1093,1.1024,5.3320)⊤
are the three accumulation points of this problem.Table 3
Numerical results for Example3.
c
0
β
=
1.5
β
=
1.8
β
=
2.0
(a)
errs
5.2594
e
-
03
2.2949
e
-
04
6.3777
e
-
02
1.5831
e
-
08
2.1853
e
-
10
5.2251
e
-
05
9.8210
e
-
12
5.3338
e
-
12
1.1737
e
-
11
ite.
17
17
17
c
*
c
2
*
c
3
*
c
2
*
(b)
errs
2.3747
e
-
02
2.4283
e
-
02
5.3326
e
-
02
6.6052
e
-
07
6.8845
e
-
07
3.5462
e
-
05
9.2871
e
-
12
7.8077
e
-
12
1.2096
e
-
11
ite.
13
14
14
c
*
c
2
*
c
2
*
c
2
*
(c)
errs
5.0900
e
-
04
5.0902
e
-
04
5.0902
e
-
04
5.5760
e
-
09
5.5819
e
-
09
5.5827
e
-
09
1.4600
e
-
11
1.3617
e
-
11
1.7533
e
-
11
ite.
10
10
10
c
*
c
1
*
c
1
*
c
1
*
(d)
errs
1.7014
e
-
02
1.7014
e
-
02
1.7014
e
-
02
6.9575
e
-
07
6.9575
e
-
07
6.9575
e
-
07
1.5224
e
-
11
1.7256
e
-
11
6.2077
e
-
11
ite.
9
9
9
c
*
c
1
*
c
1
*
c
1
*We observe from Tables1–3 that our algorithm is convergent for different starting points. We also see that our algorithm converges to a solution of the IEP, which is not necessarily equal to the original one. An interesting question is to consider the performance of the algorithm when the starting pointc0is chosen as random vector, which needs future study.



---

*Source: 630618-2013-12-05.xml*