File: node141.html

package info (click to toggle)
scalapack-doc 1.5-11
  • links: PTS
  • area: main
  • in suites: bullseye, buster, stretch
  • size: 10,336 kB
  • ctags: 4,931
  • sloc: makefile: 47; sh: 18
file content (230 lines) | stat: -rw-r--r-- 16,234 bytes parent folder | download | duplicates (4)
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">&#160;</A>
<P>
The eigendecomposition<A NAME="5404">&#160;</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">&#160;</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">&#160;</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&nbsp;<A HREF="node46.html#subsecdriveeig">3.2.3</A>).
<A NAME="5410">&#160;</A><A NAME="5411">&#160;</A>
<A NAME="5412">&#160;</A><A NAME="5413">&#160;</A>
<A NAME="5414">&#160;</A><A NAME="5415">&#160;</A>
<P>
The approximate error 
bounds
<A NAME="5417">&#160;</A>  
<A NAME="5418">&#160;</A>  
<A NAME="5419">&#160;</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">&#160;</A>
<A NAME="5426">&#160;</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">&#160;</A>
<P>
The usual error analysis of the 
symmetric<A NAME="5459">&#160;</A> 
eigenproblem 
using ScaLAPACK
driver <TT>PSSYEV</TT><A NAME="5461">&#160;</A><A NAME="5462">&#160;</A>
(see subsection&nbsp;<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">&#160;</A> while small ones may not be.
<A NAME="5479">&#160;</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">&#160;</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">&#160;</A><A NAME="5487">&#160;</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">&#160;</A><A NAME="5490">&#160;</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">&#160;</A>
<A NAME="5493">&#160;</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">&#160;</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">&#160;</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">&#160;</A> for reasons
discussed in section&nbsp;<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">&#160;</A><A NAME="5511">&#160;</A><A NAME="5512">&#160;</A><A NAME="5513">&#160;</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">&#160;</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">&#160;</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">&#160;</A>
<A NAME="5537">&#160;</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">&#160;</A><A NAME="5540">&#160;</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">&#160;</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>