1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230
|
<!DOCTYPE HTML PUBLIC "-//IETF//DTD HTML 2.0//EN">
<!--Converted with LaTeX2HTML 96.1-h (September 30, 1996) by Nikos Drakos (nikos@cbl.leeds.ac.uk), CBLU, University of Leeds -->
<HTML>
<HEAD>
<TITLE>Error Bounds for the Symmetric Eigenproblem </TITLE>
<META NAME="description" CONTENT="Error Bounds for the Symmetric Eigenproblem ">
<META NAME="keywords" CONTENT="slug">
<META NAME="resource-type" CONTENT="document">
<META NAME="distribution" CONTENT="global">
<LINK REL=STYLESHEET HREF="slug.css">
</HEAD>
<BODY LANG="EN" >
<A NAME="tex2html3974" HREF="node142.html"><IMG WIDTH=37 HEIGHT=24 ALIGN=BOTTOM ALT="next" SRC="http://www.netlib.org/utk/icons/next_motif.gif"></A> <A NAME="tex2html3972" HREF="node132.html"><IMG WIDTH=26 HEIGHT=24 ALIGN=BOTTOM ALT="up" SRC="http://www.netlib.org/utk/icons/up_motif.gif"></A> <A NAME="tex2html3966" HREF="node140.html"><IMG WIDTH=63 HEIGHT=24 ALIGN=BOTTOM ALT="previous" SRC="http://www.netlib.org/utk/icons/previous_motif.gif"></A> <A NAME="tex2html3976" HREF="node1.html"><IMG WIDTH=65 HEIGHT=24 ALIGN=BOTTOM ALT="contents" SRC="http://www.netlib.org/utk/icons/contents_motif.gif"></A> <A NAME="tex2html3977" HREF="node190.html"><IMG WIDTH=43 HEIGHT=24 ALIGN=BOTTOM ALT="index" SRC="http://www.netlib.org/utk/icons/index_motif.gif"></A> <BR>
<B> Next:</B> <A NAME="tex2html3975" HREF="node142.html">Error Bounds for the </A>
<B>Up:</B> <A NAME="tex2html3973" HREF="node132.html">Accuracy and Stability</A>
<B> Previous:</B> <A NAME="tex2html3967" HREF="node140.html">Error Bounds for Linear </A>
<BR> <P>
<H1><A NAME="SECTION04670000000000000000">Error Bounds for the Symmetric Eigenproblem </A></H1>
<A NAME="secsym"> </A>
<P>
The eigendecomposition<A NAME="5404"> </A> of
an <I>n</I>-by-<I>n</I> real symmetric matrix is the
factorization <IMG WIDTH=84 HEIGHT=15 ALIGN=BOTTOM ALT="tex2html_wrap_inline13904" SRC="img181.gif"> (<IMG WIDTH=86 HEIGHT=15 ALIGN=BOTTOM ALT="tex2html_wrap_inline18618" SRC="img634.gif"> in the complex Hermitian
case), where <I>Z</I> is orthogonal (unitary) and
<IMG WIDTH=154 HEIGHT=26 ALIGN=MIDDLE ALT="tex2html_wrap_inline18622" SRC="img635.gif"> is real and diagonal,
with <IMG WIDTH=142 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18624" SRC="img636.gif">.
The <IMG WIDTH=13 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18626" SRC="img637.gif"> are the <B>eigenvalues</B><A NAME="5406"> </A> of <I>A</I>, and the columns <IMG WIDTH=11 HEIGHT=17 ALIGN=MIDDLE ALT="tex2html_wrap_inline18630" SRC="img638.gif"> of
<I>Z</I> are the <B>eigenvectors</B><A NAME="5408"> </A>. This is also often written
<IMG WIDTH=77 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18634" SRC="img639.gif">. The eigendecomposition of a symmetric matrix is
computed
by the driver routines PxSYEV and PxSYEVX.
The complex counterparts of these routines, which compute the
eigendecomposition of complex Hermitian matrices, are
the driver routines PxHEEV and PxHEEVX
(see section <A HREF="node46.html#subsecdriveeig">3.2.3</A>).
<A NAME="5410"> </A><A NAME="5411"> </A>
<A NAME="5412"> </A><A NAME="5413"> </A>
<A NAME="5414"> </A><A NAME="5415"> </A>
<P>
The approximate error
bounds
<A NAME="5417"> </A>
<A NAME="5418"> </A>
<A NAME="5419"> </A>
for the computed eigenvalues
<IMG WIDTH=101 HEIGHT=33 ALIGN=MIDDLE ALT="tex2html_wrap_inline18636" SRC="img640.gif"> are
<BR><IMG WIDTH=326 HEIGHT=22 ALIGN=BOTTOM ALT="displaymath18596" SRC="img641.gif"><BR>
The approximate error bounds for the computed eigenvectors
<IMG WIDTH=11 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18638" SRC="img642.gif">, which bound the acute angles
between the computed eigenvectors and true
eigenvectors <IMG WIDTH=11 HEIGHT=17 ALIGN=MIDDLE ALT="tex2html_wrap_inline18630" SRC="img638.gif">, are
<A NAME="5425"> </A>
<A NAME="5426"> </A>
<BR><IMG WIDTH=331 HEIGHT=17 ALIGN=BOTTOM ALT="displaymath18597" SRC="img643.gif"><BR>
These bounds can be computed by the following code fragment:
<P>
<PRE> EPSMCH = PSLAMCH( ICTXT, 'E' )
* Compute eigenvalues and eigenvectors of A
* The eigenvalues are returned in W
* The eigenvector matrix Z overwrites A
CALL PSSYEV( 'V', UPLO, N, A, IA, JA, DESCA, W, Z, IZ, JZ,
$ DESCZ, WORK, LWORK, INFO )
IF( INFO.GT.0 ) THEN
PRINT *,'PSSYEV did not converge'
ELSE IF( N.GT.0 ) THEN
* Compute the norm of A
ANORM = MAX( ABS( W( 1 ) ), ABS( W( N ) ) )
EERRBD = EPSMCH * ANORM
* Compute reciprocal condition numbers for eigenvectors
CALL SDISNA( 'Eigenvectors', N, N, W, RCONDZ, INFO )
DO 10 I = 1, N
ZERRBD( I ) = EPSMCH * ( ANORM / RCONDZ( I ) )
10 CONTINUE
END IF</PRE>
<P>
For example, if <IMG WIDTH=315 HEIGHT=30 ALIGN=MIDDLE ALT="tex2html_wrap_inline18242" SRC="img571.gif"> and
<BR><BR>
<BR><IMG WIDTH=325 HEIGHT=67 ALIGN=BOTTOM ALT="displaymath18598" SRC="img644.gif"><BR>
then the eigenvalues, approximate error bounds, and true errors are
as follows.
<P>
<BR><IMG WIDTH=527 HEIGHT=90 ALIGN=BOTTOM ALT="tabular5436" SRC="img645.gif"><BR>
<P>
<B>Further Details: Error Bounds for the Symmetric Eigenproblem</B><A NAME="secDetailsym"> </A>
<P>
The usual error analysis of the
symmetric<A NAME="5459"> </A>
eigenproblem
using ScaLAPACK
driver <TT>PSSYEV</TT><A NAME="5461"> </A><A NAME="5462"> </A>
(see subsection <A HREF="node46.html#subsecdriveeig">3.2.3</A>)
is as follows [<A HREF="node189.html#parl80">101</A>]:
<P>
<BLOCKQUOTE> The computed eigendecomposition <IMG WIDTH=47 HEIGHT=16 ALIGN=BOTTOM ALT="tex2html_wrap_inline18680" SRC="img646.gif"> is nearly
the exact
eigendecomposition of <I>A</I>+<I>E</I>, namely,
<IMG WIDTH=232 HEIGHT=32 ALIGN=MIDDLE ALT="tex2html_wrap_inline18684" SRC="img647.gif">
is a true eigendecomposition so that <IMG WIDTH=54 HEIGHT=34 ALIGN=MIDDLE ALT="tex2html_wrap_inline18686" SRC="img648.gif"> is
orthogonal,
where <IMG WIDTH=146 HEIGHT=27 ALIGN=MIDDLE ALT="tex2html_wrap_inline18688" SRC="img649.gif"> and
<IMG WIDTH=107 HEIGHT=32 ALIGN=MIDDLE ALT="tex2html_wrap_inline18690" SRC="img650.gif">. Here <I>p</I>(<I>n</I>) is a modestly growing
function of <I>n</I>. We take <I>p</I>(<I>n</I>)=1 in the above code fragment.
Each computed eigenvalue <IMG WIDTH=13 HEIGHT=33 ALIGN=MIDDLE ALT="tex2html_wrap_inline18698" SRC="img651.gif"> differs from a
true <IMG WIDTH=13 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18626" SRC="img637.gif"> by at most
<BR><IMG WIDTH=387 HEIGHT=22 ALIGN=BOTTOM ALT="displaymath18599" SRC="img652.gif"><BR>
Thus, large eigenvalues (those near <IMG WIDTH=125 HEIGHT=27 ALIGN=MIDDLE ALT="tex2html_wrap_inline18702" SRC="img653.gif">)
are computed to high relative accuracy, <A NAME="5478"> </A> while small ones may not be.
<A NAME="5479"> </A>
<P>
The angular difference between the computed unit eigenvector
<IMG WIDTH=11 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18704" SRC="img654.gif"> and a true unit eigenvector <IMG WIDTH=11 HEIGHT=17 ALIGN=MIDDLE ALT="tex2html_wrap_inline18630" SRC="img638.gif"> satisfies the approximate bound
<A NAME="5480"> </A>
<BR><IMG WIDTH=378 HEIGHT=41 ALIGN=BOTTOM ALT="displaymath18600" SRC="img655.gif"><BR>
if <IMG WIDTH=41 HEIGHT=26 ALIGN=MIDDLE ALT="tex2html_wrap_inline18543" SRC="img626.gif"> is small enough.
Here, <IMG WIDTH=171 HEIGHT=28 ALIGN=MIDDLE ALT="tex2html_wrap_inline18712" SRC="img656.gif"> is the
<B>absolute gap</B><A NAME="5486"> </A><A NAME="5487"> </A>
between <IMG WIDTH=13 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18626" SRC="img637.gif"> and the nearest other eigenvalue. Thus, if <IMG WIDTH=13 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18626" SRC="img637.gif">
is close to other eigenvalues, its corresponding eigenvector <IMG WIDTH=11 HEIGHT=17 ALIGN=MIDDLE ALT="tex2html_wrap_inline18630" SRC="img638.gif">
may be inaccurate.
The gaps may be easily computed from the array of computed eigenvalues
by using subroutine <TT>SDISNA</TT><A NAME="5489"> </A><A NAME="5490"> </A>.
The gaps computed by <TT>SDISNA</TT> are ensured not to be so small as
to cause overflow when used as divisors.
<A NAME="5492"> </A>
<A NAME="5493"> </A>
<P>
Let <IMG WIDTH=11 HEIGHT=16 ALIGN=BOTTOM ALT="tex2html_wrap_inline17614" SRC="img477.gif"> be the invariant subspace spanned by a collection of eigenvectors
<IMG WIDTH=81 HEIGHT=27 ALIGN=MIDDLE ALT="tex2html_wrap_inline18722" SRC="img657.gif">, where <IMG WIDTH=12 HEIGHT=12 ALIGN=BOTTOM ALT="tex2html_wrap_inline18724" SRC="img658.gif"> is a subset of the
integers from 1 to <I>n</I>. Let <IMG WIDTH=11 HEIGHT=12 ALIGN=BOTTOM ALT="tex2html_wrap_inline17602" SRC="img475.gif"> be the corresponding true subspace. Then
<BR><IMG WIDTH=330 HEIGHT=41 ALIGN=BOTTOM ALT="displaymath18601" SRC="img659.gif"><BR>
<A NAME="5498"> </A>
where
<BR><IMG WIDTH=326 HEIGHT=35 ALIGN=BOTTOM ALT="displaymath18602" SRC="img660.gif"><BR>
is the absolute gap between the eigenvalues in <IMG WIDTH=12 HEIGHT=12 ALIGN=BOTTOM ALT="tex2html_wrap_inline18724" SRC="img658.gif"> and the nearest
other eigenvalue. Thus, a cluster<A NAME="5502"> </A> of
close eigenvalues that is
far away from any other eigenvalue may have a well-determined
invariant subspace <IMG WIDTH=11 HEIGHT=16 ALIGN=BOTTOM ALT="tex2html_wrap_inline17614" SRC="img477.gif"> even if its individual eigenvectors are
ill-conditioned.
</BLOCKQUOTE>
<P>
A small possibility exists that <TT>PSSYEV</TT> will fail to achieve the
above error bounds on a heterogeneous network of processors
<A NAME="5505"> </A> for reasons
discussed in section <A HREF="node134.html#sec_Hetero">6.2</A>. On a homogeneous network,
<TT>PSSYEV</TT> is as robust as the corresponding LAPACK routine <TT>SSYEV</TT>.
A future release will attempt to detect heterogeneity and warn the user
to use an alternative algorithm.
<P>
In contrast to LAPACK, where the same error analysis applied to the
simple and expert drivers, the expert driver <TT>PSSYEVX</TT>
<A NAME="5510"> </A><A NAME="5511"> </A><A NAME="5512"> </A><A NAME="5513"> </A>
satisfies slightly weaker error bounds than <TT>PSSYEV</TT>. The bounds
<IMG WIDTH=138 HEIGHT=34 ALIGN=MIDDLE ALT="tex2html_wrap_inline18736" SRC="img661.gif"> and
<IMG WIDTH=134 HEIGHT=26 ALIGN=MIDDLE ALT="tex2html_wrap_inline18738" SRC="img662.gif"> continue to hold,
but the computed eigenvectors <IMG WIDTH=11 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18638" SRC="img642.gif"> are no longer guaranteed to
be orthogonal to one another. The corresponding LAPACK routine
<TT>SSYEVX</TT> tries to guarantee orthogonality by <EM>reorthogonalizing</EM>
computed eigenvectors against one another provided their corresponding
computed eigenvalues are close enough together:
<IMG WIDTH=213 HEIGHT=34 ALIGN=MIDDLE ALT="tex2html_wrap_inline18742" SRC="img663.gif">, where the threshold
<IMG WIDTH=102 HEIGHT=16 ALIGN=BOTTOM ALT="tex2html_wrap_inline18744" SRC="img664.gif">. If <I>m</I> eigenvalues lie in a cluster satisfying
this closeness criterion, the <TT>SSYEVX</TT> requires <IMG WIDTH=58 HEIGHT=30 ALIGN=MIDDLE ALT="tex2html_wrap_inline18748" SRC="img665.gif"> serial time to
execute. When <I>m</I> is a large fraction of <I>n</I>, this serial bottleneck is
expensive and does not always improve orthogonality.
<P>
ScaLAPACK addresses this problem in two ways. First, it lets the
user use more or less time and space to perform reorthogonalization, rather
than have a fixed criterion. In particular, the user can set the
threshold <TT>ORTOL</TT> used above, decreasing it to make reorthogonalization
less frequent or increasing it to reorthogonalize more.
Furthermore, since each processor computes a subset of the eigenvectors,
ScaLAPACK permits reorthogonalization only with the local eigenvectors;
that is, no communication is allowed. Hence, if a cluster of eigenvalues
is small enough for the corresponding eigenvectors to fit on one processor,
the same reorthogonalization will be done as in LAPACK.
The user can supply more or less workspace to limit the
size of a cluster on one processor. Hence, at one extreme, with a large cluster
and lots of workspace, the algorithm will be essentially equivalent
to <TT>SSYEVX</TT>. At the other extreme, with all small clusters or little
workspace supplied, the algorithm will be perfectly load balanced and
perform minimal communication to compute the eigenvectors.
<P>
The second way ScaLAPACK will deal with reorthogonalization is to introduce
a new algorithm [<A HREF="node189.html#parlettdhillon96">103</A>, <A HREF="node189.html#parlett96">102</A>, <A HREF="node189.html#dhillonthesis">44</A>]
that requires nearly <EM>no</EM> reorthogonalization
to compute orthogonal eigenvectors in a fully parallel way. This algorithm
will be introduced in future ScaLAPACK and LAPACK releases.
<P>
In the special case of a real symmetric tridiagonal matrix <I>T</I>, the eigenvalues
can sometimes be computed much more accurately. PxSYEV (and the other
symmetric eigenproblem drivers) computes the eigenvalues and eigenvectors of
a dense symmetric matrix by first reducing it to tridiagonal form<A NAME="5533"> </A> <I>T</I> and then
finding the eigenvalues and eigenvectors of <I>T</I>.
Reduction of a dense matrix to tridiagonal form<A NAME="5534"> </A> <I>T</I> can introduce
additional errors, so the following bounds for the tridiagonal case do not
apply to the dense case.
<P>
<BLOCKQUOTE> The eigenvalues of the symmetric tridiagonal matrix
<I>T</I> may be computed with small componentwise relative
backward error
<A NAME="5536"> </A>
<A NAME="5537"> </A>
(<IMG WIDTH=31 HEIGHT=26 ALIGN=MIDDLE ALT="tex2html_wrap_inline18764" SRC="img666.gif">) by using subroutine PxSTEBZ
(section <A HREF="node61.html#subseccompsep">3.3.4</A>).
<A NAME="5539"> </A><A NAME="5540"> </A>
To compute tighter error bounds
for the computed eigenvalues <IMG WIDTH=13 HEIGHT=33 ALIGN=MIDDLE ALT="tex2html_wrap_inline18698" SRC="img651.gif"> we must make some assumptions
about <I>T</I>. The bounds discussed here are from [<A HREF="node189.html#barlowdemmel">15</A>].
Suppose <I>T</I> is positive definite, and
write <I>T</I>=<I>DHD</I> where <IMG WIDTH=175 HEIGHT=39 ALIGN=MIDDLE ALT="tex2html_wrap_inline18774" SRC="img667.gif">
and <IMG WIDTH=49 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18776" SRC="img668.gif">.
Then the computed eigenvalues
<IMG WIDTH=13 HEIGHT=33 ALIGN=MIDDLE ALT="tex2html_wrap_inline18698" SRC="img651.gif"> can differ from true eigenvalues <IMG WIDTH=13 HEIGHT=25 ALIGN=MIDDLE ALT="tex2html_wrap_inline18626" SRC="img637.gif"> by
<BR><IMG WIDTH=364 HEIGHT=22 ALIGN=BOTTOM ALT="displaymath18603" SRC="img669.gif"><BR>
where <I>p</I>(<I>n</I>) is a modestly growing function of <I>n</I>.
Thus, if <IMG WIDTH=44 HEIGHT=26 ALIGN=MIDDLE ALT="tex2html_wrap_inline18786" SRC="img670.gif"> is moderate, each eigenvalue will be computed
to high relative accuracy, <A NAME="5547"> </A> no matter how tiny it is.
<P>
</BLOCKQUOTE><HR><A NAME="tex2html3974" HREF="node142.html"><IMG WIDTH=37 HEIGHT=24 ALIGN=BOTTOM ALT="next" SRC="http://www.netlib.org/utk/icons/next_motif.gif"></A> <A NAME="tex2html3972" HREF="node132.html"><IMG WIDTH=26 HEIGHT=24 ALIGN=BOTTOM ALT="up" SRC="http://www.netlib.org/utk/icons/up_motif.gif"></A> <A NAME="tex2html3966" HREF="node140.html"><IMG WIDTH=63 HEIGHT=24 ALIGN=BOTTOM ALT="previous" SRC="http://www.netlib.org/utk/icons/previous_motif.gif"></A> <A NAME="tex2html3976" HREF="node1.html"><IMG WIDTH=65 HEIGHT=24 ALIGN=BOTTOM ALT="contents" SRC="http://www.netlib.org/utk/icons/contents_motif.gif"></A> <A NAME="tex2html3977" HREF="node190.html"><IMG WIDTH=43 HEIGHT=24 ALIGN=BOTTOM ALT="index" SRC="http://www.netlib.org/utk/icons/index_motif.gif"></A> <BR>
<B> Next:</B> <A NAME="tex2html3975" HREF="node142.html">Error Bounds for the </A>
<B>Up:</B> <A NAME="tex2html3973" HREF="node132.html">Accuracy and Stability</A>
<B> Previous:</B> <A NAME="tex2html3967" HREF="node140.html">Error Bounds for Linear </A>
<P><ADDRESS>
<I>Susan Blackford <BR>
Tue May 13 09:21:01 EDT 1997</I>
</ADDRESS>
</BODY>
</HTML>
|