Matrices with special structures arise in numerous applications. In
some cases, such as quasidefinite matrices or their generalizations,
we can exploit this special structure. If the matrix H is quasidefinite,
we propose a new variant of the symmetric indefinite factorization. 
We show that linear system Hz = b, H quasidefinite
with a special structure, can be interpreted as an equilibrium system.
So, even if some blocks in H are ill--conditioned, the important part of solution vector z can be accurately computed. In the case of a
generalized quasidefinite matrix, we derive bounds on number of its
positive and negative eigenvalues

Sanja Singer

Singer Saša

Hrčak - Portal of scientific journals of Croatia

Mathematical Communications 4(1999), 19–25 19Symmetric indefinite factorization of quasidefinitematrices∗Sanja Singer† and Sasˇa Singer‡Abstract. Matrices with special structures arise in numerous ap-plications. In some cases, such as quasidefinite matrices or their gen-eralizations, we can exploit this special structure. If the matrix H isquasidefinite, we propose a new variant of the symmetric indefinite fac-torization. We show that linear system Hz = b, H quasidefinite with aspecial structure, can be interpreted as an equilibrium system. So, evenif some blocks in H are ill–conditioned, the important part of solutionvector z can be accurately computed. In the case of a generalized quasi-definite matrix, we derive bounds on number of its positive and negativeeigenvalues.Key words: quasidefinite matrices, inertia, special linear systems,accurate solutionAMS subject classifications: 15A23, 65F05, 65G061. IntroductionQuasidefinite matrices have the formH =[A BB∗ −D], (1)where A and D are symmetric (Hermitian) and positive definite. This class ofmatrices arises in optimization when barrier or interior–point methods are applied(see [6], [1]).For example, to solve the damped linear least squares problemminx‖b−Ax‖22 + δ2‖x‖22,∗This paper was presented at the Special Section: Appl. Math. Comput. at the 7th InternationalConference on OR, Rovinj, 1998.†Faculty of Mechanical Engineering and Naval Architecture, I. Lucˇic´a 5, HR-10 000 Zagreb,Croatia, e-mail: ssinger@math.hr‡Department of Mathematics, University of Zagreb, Bijenicˇka 30, HR-10 000 Zagreb, Croatia,e-mail: singer@math.hr20 Sanja Singer and Sasˇa Singerone needs to solve the linear system[A BBτ −D] [xy]=[b0](2)with A = D = δI, δ > 0.2. Factorization of quasidefinite matricesIt is a well known fact that a symmetric indefinite factorization of the matrix (1)can be computed using only 1 × 1 pivots. Here, we give a new, elementary andconstructive proof of this fact.Let A = G∗1G1 and D = G∗3G3 be Cholesky factorizations of A and D, respec-tively. The matrix H can be factorized asH =[G∗1 0T ∗ G∗2] [Ik 00 −I] [G1 T0 G2]=[G∗1G1 G∗1TT ∗G1 T∗T −G∗2G2](3)withT = G−∗1 B,G∗2G2 = D + B∗A−1B = G∗3G3 + T∗T. (4)From (4) we see that the matrix G2 can be computed as the upper triangularfactor from the QR factorization[TG3]= Q[G20].From (3) we can also conclude a somewhat surprising fact that the inertia ofH does not depend on B – the matrix H has exactly k positive and  negativeeigenvalues.To ensure a numerical stability we usually use Bunch–Kaufman or Bunch–Parlett pivoting. But, if the matrices A and D are both diagonal, this sparsitypattern cannot be preserved. In this case, we can use diagonal pivoting in blocks Aand D to reduce the fill–in of nonzero elements.Another reason why permutations should be avoided lies in parallel computation.An algorithm which includes pivoting or permutation processes is almost strictlysequential and cannot take advantage of parallel processing.But, here lies a danger. If pivots in A are too small compared with elements ofB, big elements can occur in T and G2. Note, if A and D are diagonal, computedelements of G1 and T have small relative errors. The diagonals of G1 and G3are computed using only square roots and off diagonal elements of G1 and G3 areexactly 0. Elements of T are computed from B and G1 using only one division foreach of them. For errors in computed matrices we have|(δG1)ii| ≤ ε|(G1)ii|,|(δG3)jj | ≤ ε|(G3)jj |, |(δT )ij | ≤ ε|(T )ij |, 1 ≤ i ≤ k, 1 ≤ j ≤ , (5)Factorization of quasidefinite matrices 21where ε is a rounding error of basic floating–point operations.What happens with the computed G2? Using Higham [2] Lemma 18.8., we mayinterpret the computed QR factorization of the matrixM =[T˜G˜3]:=[T + δTG3 + δG3](6)as an exact factorization of a slightly perturbed matrix M + δM . The backwarderror matrix δM satisfies|δM | ≤ f(m, ) εmG |M |, (7)where m = k +  and  are dimensions of M , f(m, ) is a moderate polynomial andG = m−1eeτ , with e = [1, 1, . . . , 1]τ .The following theorem can be proved.Theorem 1. Let H be a quasidefinite matrix with A and D diagonal. Let G2 in(3) be obtained by QR factorization. Suppose that in the floating point arithmeticthe symmetric indefinite factorization produces matrixR˜ =[G˜1 T˜0 G˜2].Then R˜ is an exact factor of the perturbed matrixH˜ := H + δH =[A + δA B + δBB∗ + δB∗ −D − δD]with|δA| ≤ εpA|δB| ≤ εp|B||δD| ≤ εpD + εm|M∗|G |M |≤ εpD + (1 + ε)2εm[ |B∗|G−11 G3 ] G [ G−11 |B|G3],where m = k + , εp := 2ε + ε2, εm := εf(m, )(2m+ εf(m, )m3).Proof. Multiplication of computed factors givesH˜ =[G˜∗1G˜1 G˜∗1T˜T˜ ∗G˜1 T˜∗T˜ − G˜∗2G˜2].From (5), it follows that |δG1| ≤ ε|G1|, andG˜∗1G˜1 = G∗1G1 + δG∗1G1 + G∗1δG1 + δG∗1δG1 := A + δA,or|δA| ≤ (2ε + ε2)A.22 Sanja Singer and Sasˇa SingerThe bound for |δB| can be obtained similarly, if we notice that |B| = G1|T |.Finally, from (6) it follows−D − δD = T˜ ∗T˜ − G˜∗2G˜2= T˜ ∗T˜ − (M∗ + δM∗)(M + δM)= −G˜∗3G˜3 − δM∗M −M∗δM − δM∗δM.By definition of G˜3, we haveG˜∗3G˜3 = D + δG∗3G3 + G∗3δG3 + δG∗3δG3,and|δD| ≤ |δG∗3G3 + G∗3δG3 + δG∗3δG3 + δM∗M + M∗δM + δM∗δM |≤ εpD + |δM∗| |M |+ |M∗| |δM |+ |δM∗| |δM |.The first inequality for |δD| follows from (7). The second inequality is a consequenceof|M | =[ |T + δT ||G3 + δG3|]≤ (1 + ε)[ |T ||G3|]= (1 + ε)[G−11 |B|G3].It is easy to see that the computed factorization of H is “good” if the elementsof |M∗|G |M | are comparable in magnitude with the elements of D.3. Inertia of generalized quasidefinite matricesHigham and Cheng [3] have studied the inertia of matrices arising in optimization.In view of this, it is an interesting question how the inertia of H depends on off–diagonal block B, if A and D from (1) are Hermitian, nonsingular, but possiblyindefinite.Definition 1. A nonsingular matrix H is generalized quasidefinite if the blocksA and D from (1) are nonsingular and indefinite.In this case we can give some bounds on the inertia of H . The inertia of H isan ordered triple (i+, i−, i0), where i+, i− and i0 are numbers of positive, negativeand zero eigenvalues of H , respectively.Theorem 2. Let H be a generalized quasidefinite matrix. Theninertia(H) = (i+(H), i−(H), 0)satisfiesmax{i+(A), − i+(D)} ≤ i+(H) ≤ min{ + i+(A), k + i−(D)},andi−(H) = k + − i+(H).Factorization of quasidefinite matrices 23Proof. If A and D are nonsingular, they have symmetric indefinite factorizationsA = P ∗1G∗1J1G1P1, D = P∗3G∗3J3G3P3,where, for i = 1, 3, Pi are permutation matrices, Ji are signature matrices and Giare block upper triangular matrices with diagonal blocks of order 1 or 2.Let T = J1G−∗1 P1B and S = −D − T ∗J1T . The matrix H can be written asH =[P ∗1G∗1 0T ∗ I] [J1 00 S] [G1P1 T0 I]=[P ∗1G1J1G1P1 P ∗1G∗1J1TT ∗J1G1P1 T∗J1T + S].Since H is nonsingular, we conclude that S is nonsingular too, and it can befactorized asS = P ∗2G∗2J2G2P2.Then−D = T ∗J1T + S,orP ∗2G∗2J2G2P2 = −P ∗3G∗3J3G3P3 − T ∗J1T= [ P ∗3G∗3 T ∗ ][ −J3 00 −J1] [G3P3T].Let inertia(A) = (i+(A), i−(A), 0) and inertia(D) = (i+(D), i−(D), 0). Nonsin-gularity of S implies that G2 can be obtained by the indefinite QR factorizationof [G3P3T]where the indefinite inner product is generated by the matrixJ :=[ −J3 00 −J1](see Singer [5]). The inertia of J is equal toinertia(J) = (i−(A) + i−(D), i+(A) + i+(D), 0).Construction of the indefinite factorization implies thatinertia(J2) = (i+(J2), − i+(J2), 0),and0 ≤ i+(J2) ≤ i−(A) + i−(D)0 ≤ − i+(J2) ≤ i+(A) + i+(D),ormax{0, − i+(A)− i+(D)} ≤ i+(J2) ≤ min{, i−(A) + i−(D)}.Finally, we conclude thatmax{i+(A), − i+(D)} ≤ i+(H) ≤ min{ + i+(A), k + i−(D)}.24 Sanja Singer and Sasˇa SingerIt is easy to show that the bounds from the previous theorem are nontrivial.Example 1. Let H be a generalized quasidefinite matrix with indefinite blocksA and D, such that inertia(A) = (1, 7, 0) and inertia(D) = (2, 4, 0). From theprevious theorem we havemax{1, 6− 2} ≤ i+(H) ≤ min{6 + 1, 8 + 4},or4 ≤ i+(H) ≤ 7, i−(H) = 14− i+(H).4. Accurate solution of special linear systemsFinally, we consider a special case of quasidefinite linear system (2), when A andD are diagonal matrices, with emphasis on accurate computation of the y part ofthe solution z.Such a system is a generalization of the so-called equilibrium system[Â B̂B̂τ 0] [x̂ŷ]=[b̂0], (8)with diagonal positive definite Â and full column rank B̂.This system can be written asÂx̂ + B̂ŷ = b̂B̂τ x̂ = 0.The first equation gives x̂ = Â−1 (̂b−B̂ŷ). Substituting this into the second equationgivesŷ = (B̂τ Â−1B̂)−1B̂τ Â−1b̂. (9)So, if we can bound ‖(B̂τ Â−1B̂)−1B̂τ Â−1‖ independent of Â, ŷ cannot be muchlarger than b̂, no matter how ill–conditioned Â is.Vavasis in [7] proves the following result (originally proved independently byStewart and Todd).Theorem 3. Let A denote the set of all k × k positive definite real diagonalmatrices. Let B̂ be a k×  real matrix of rank . Then there exist constants χ 0B andχ¯ 0B such that for any Â ∈ A(a) ‖(B̂τ Â−1B̂)−1B̂τ Â−1‖ ≤ χ 0B, and(b) ‖B̂(B̂τ Â−1B̂)−1B̂τ Â−1‖ ≤ χ¯ 0B.Here we assume that the matrix norm ‖ ‖ is induced by some vector norm. Our quasidefinite problem[A BBτ −D] [xy]=[b0], (10)Factorization of quasidefinite matrices 25can be written in a componentwise form asAx + By = bBτx−Dy = 0.Similarly, from the first equation we have x = A−1(b−By), andy = (BτA−1B + D)−1BτA−1b. (11)If ‖(BτA−1B + D)−1BτA−1‖ is bounded and depends only on B and D, ill–condition of A should not destroy the components of y. To show this, we transformthe system (10) into the form (8).Let D = ∆2 be the Cholesky factorization of D. DefineÂ =[A 00 I], B̂ =[B∆], b̂ =[b0]. (12)It is obvious that if A is positive definite, Â is positive definite too. Also, matrixB̂ has full column rank because ∆ is a nonsingular diagonal matrix.The connection between the solution ŷ of linear system (8) and the y–part of(10) is very simple.Theorem 4. If matrices Â and B̂ and vector b̂ are defined by (12) then ŷ = y.Proof. From (9), using (11) and (12), it followsŷ =([Bτ ∆] [ A−1 00 I] [B∆])−1 [Bτ ∆] [ A 00 I]−1 [b0]= (BτA−1B + D)−1[BτA−1 ∆] [ b0]= (BτA−1B + D)−1BτA−1b = y.Vavasis [7] and Hough and Vavasis [4] have developed algorithms for stablecomputation of ŷ component for the linear system (8). Theorem 4 shows that thesealgorithms are suitable for quasidefinite linear systems (10) as well.References[1] P.E.Gill, M.A. Saunders, J. R. Shinnerl, On the stability of Choleskyfactorization for symmetric quasidefinite systems, SIAM J. Matrix Anal. Appl.17(1996), 35–46.[2] N. J.Higham, Accuracy and Stability of Numerical Algorithms, SIAM, Philadel-phia, 1996.[3] N. J.Higham, S.H.Cheng, Modifying the inertia of matrices arising in opti-mization, Linear Algebra Appl. 275–276(1998), 261–279.26 Sanja Singer and Sasˇa Singer[4] P.D.Hough, S.A.Vavasis, Complete orthogonal decomposition for weightedleast squares, SIAM J. Matrix Anal. Appl. 18(1997), 369–392.[5] S. Singer, Indefinite QR Factorization and Its Applications, Ph.D. thesis, De-partment of Mathematics, University of Zagreb, 1997 (In Croatian).[6] R. J.Vanderbei, Symmetric quasidefinite matrices, SIAM J. Optim. 5(1995),100–113.[7] S.A.Vavasis, Stable numerical algorithms for equilibrium systems, SIAM J.Matrix Anal. Appl. 15(1994), 1108–1131.

Symmetric indefinite factorization of quasidefinite matrices

Abstract

Similar works

Full text

Available Versions

Hrčak - Portal of scientific journals of Croatia