![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DLAED1 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download DLAED1 + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaed1.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaed1.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaed1.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE DLAED1( N, D, Q, LDQ, INDXQ, RHO, CUTPNT, WORK, IWORK, 00022 * INFO ) 00023 * 00024 * .. Scalar Arguments .. 00025 * INTEGER CUTPNT, INFO, LDQ, N 00026 * DOUBLE PRECISION RHO 00027 * .. 00028 * .. Array Arguments .. 00029 * INTEGER INDXQ( * ), IWORK( * ) 00030 * DOUBLE PRECISION D( * ), Q( LDQ, * ), WORK( * ) 00031 * .. 00032 * 00033 * 00034 *> \par Purpose: 00035 * ============= 00036 *> 00037 *> \verbatim 00038 *> 00039 *> DLAED1 computes the updated eigensystem of a diagonal 00040 *> matrix after modification by a rank-one symmetric matrix. This 00041 *> routine is used only for the eigenproblem which requires all 00042 *> eigenvalues and eigenvectors of a tridiagonal matrix. DLAED7 handles 00043 *> the case in which eigenvalues only or eigenvalues and eigenvectors 00044 *> of a full symmetric matrix (which was reduced to tridiagonal form) 00045 *> are desired. 00046 *> 00047 *> T = Q(in) ( D(in) + RHO * Z*Z**T ) Q**T(in) = Q(out) * D(out) * Q**T(out) 00048 *> 00049 *> where Z = Q**T*u, u is a vector of length N with ones in the 00050 *> CUTPNT and CUTPNT + 1 th elements and zeros elsewhere. 00051 *> 00052 *> The eigenvectors of the original matrix are stored in Q, and the 00053 *> eigenvalues are in D. The algorithm consists of three stages: 00054 *> 00055 *> The first stage consists of deflating the size of the problem 00056 *> when there are multiple eigenvalues or if there is a zero in 00057 *> the Z vector. For each such occurence the dimension of the 00058 *> secular equation problem is reduced by one. This stage is 00059 *> performed by the routine DLAED2. 00060 *> 00061 *> The second stage consists of calculating the updated 00062 *> eigenvalues. This is done by finding the roots of the secular 00063 *> equation via the routine DLAED4 (as called by DLAED3). 00064 *> This routine also calculates the eigenvectors of the current 00065 *> problem. 00066 *> 00067 *> The final stage consists of computing the updated eigenvectors 00068 *> directly using the updated eigenvalues. The eigenvectors for 00069 *> the current problem are multiplied with the eigenvectors from 00070 *> the overall problem. 00071 *> \endverbatim 00072 * 00073 * Arguments: 00074 * ========== 00075 * 00076 *> \param[in] N 00077 *> \verbatim 00078 *> N is INTEGER 00079 *> The dimension of the symmetric tridiagonal matrix. N >= 0. 00080 *> \endverbatim 00081 *> 00082 *> \param[in,out] D 00083 *> \verbatim 00084 *> D is DOUBLE PRECISION array, dimension (N) 00085 *> On entry, the eigenvalues of the rank-1-perturbed matrix. 00086 *> On exit, the eigenvalues of the repaired matrix. 00087 *> \endverbatim 00088 *> 00089 *> \param[in,out] Q 00090 *> \verbatim 00091 *> Q is DOUBLE PRECISION array, dimension (LDQ,N) 00092 *> On entry, the eigenvectors of the rank-1-perturbed matrix. 00093 *> On exit, the eigenvectors of the repaired tridiagonal matrix. 00094 *> \endverbatim 00095 *> 00096 *> \param[in] LDQ 00097 *> \verbatim 00098 *> LDQ is INTEGER 00099 *> The leading dimension of the array Q. LDQ >= max(1,N). 00100 *> \endverbatim 00101 *> 00102 *> \param[in,out] INDXQ 00103 *> \verbatim 00104 *> INDXQ is INTEGER array, dimension (N) 00105 *> On entry, the permutation which separately sorts the two 00106 *> subproblems in D into ascending order. 00107 *> On exit, the permutation which will reintegrate the 00108 *> subproblems back into sorted order, 00109 *> i.e. D( INDXQ( I = 1, N ) ) will be in ascending order. 00110 *> \endverbatim 00111 *> 00112 *> \param[in] RHO 00113 *> \verbatim 00114 *> RHO is DOUBLE PRECISION 00115 *> The subdiagonal entry used to create the rank-1 modification. 00116 *> \endverbatim 00117 *> 00118 *> \param[in] CUTPNT 00119 *> \verbatim 00120 *> CUTPNT is INTEGER 00121 *> The location of the last eigenvalue in the leading sub-matrix. 00122 *> min(1,N) <= CUTPNT <= N/2. 00123 *> \endverbatim 00124 *> 00125 *> \param[out] WORK 00126 *> \verbatim 00127 *> WORK is DOUBLE PRECISION array, dimension (4*N + N**2) 00128 *> \endverbatim 00129 *> 00130 *> \param[out] IWORK 00131 *> \verbatim 00132 *> IWORK is INTEGER array, dimension (4*N) 00133 *> \endverbatim 00134 *> 00135 *> \param[out] INFO 00136 *> \verbatim 00137 *> INFO is INTEGER 00138 *> = 0: successful exit. 00139 *> < 0: if INFO = -i, the i-th argument had an illegal value. 00140 *> > 0: if INFO = 1, an eigenvalue did not converge 00141 *> \endverbatim 00142 * 00143 * Authors: 00144 * ======== 00145 * 00146 *> \author Univ. of Tennessee 00147 *> \author Univ. of California Berkeley 00148 *> \author Univ. of Colorado Denver 00149 *> \author NAG Ltd. 00150 * 00151 *> \date November 2011 00152 * 00153 *> \ingroup auxOTHERcomputational 00154 * 00155 *> \par Contributors: 00156 * ================== 00157 *> 00158 *> Jeff Rutter, Computer Science Division, University of California 00159 *> at Berkeley, USA \n 00160 *> Modified by Francoise Tisseur, University of Tennessee 00161 *> 00162 * ===================================================================== 00163 SUBROUTINE DLAED1( N, D, Q, LDQ, INDXQ, RHO, CUTPNT, WORK, IWORK, 00164 $ INFO ) 00165 * 00166 * -- LAPACK computational routine (version 3.4.0) -- 00167 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00168 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00169 * November 2011 00170 * 00171 * .. Scalar Arguments .. 00172 INTEGER CUTPNT, INFO, LDQ, N 00173 DOUBLE PRECISION RHO 00174 * .. 00175 * .. Array Arguments .. 00176 INTEGER INDXQ( * ), IWORK( * ) 00177 DOUBLE PRECISION D( * ), Q( LDQ, * ), WORK( * ) 00178 * .. 00179 * 00180 * ===================================================================== 00181 * 00182 * .. Local Scalars .. 00183 INTEGER COLTYP, I, IDLMDA, INDX, INDXC, INDXP, IQ2, IS, 00184 $ IW, IZ, K, N1, N2, ZPP1 00185 * .. 00186 * .. External Subroutines .. 00187 EXTERNAL DCOPY, DLAED2, DLAED3, DLAMRG, XERBLA 00188 * .. 00189 * .. Intrinsic Functions .. 00190 INTRINSIC MAX, MIN 00191 * .. 00192 * .. Executable Statements .. 00193 * 00194 * Test the input parameters. 00195 * 00196 INFO = 0 00197 * 00198 IF( N.LT.0 ) THEN 00199 INFO = -1 00200 ELSE IF( LDQ.LT.MAX( 1, N ) ) THEN 00201 INFO = -4 00202 ELSE IF( MIN( 1, N / 2 ).GT.CUTPNT .OR. ( N / 2 ).LT.CUTPNT ) THEN 00203 INFO = -7 00204 END IF 00205 IF( INFO.NE.0 ) THEN 00206 CALL XERBLA( 'DLAED1', -INFO ) 00207 RETURN 00208 END IF 00209 * 00210 * Quick return if possible 00211 * 00212 IF( N.EQ.0 ) 00213 $ RETURN 00214 * 00215 * The following values are integer pointers which indicate 00216 * the portion of the workspace 00217 * used by a particular array in DLAED2 and DLAED3. 00218 * 00219 IZ = 1 00220 IDLMDA = IZ + N 00221 IW = IDLMDA + N 00222 IQ2 = IW + N 00223 * 00224 INDX = 1 00225 INDXC = INDX + N 00226 COLTYP = INDXC + N 00227 INDXP = COLTYP + N 00228 * 00229 * 00230 * Form the z-vector which consists of the last row of Q_1 and the 00231 * first row of Q_2. 00232 * 00233 CALL DCOPY( CUTPNT, Q( CUTPNT, 1 ), LDQ, WORK( IZ ), 1 ) 00234 ZPP1 = CUTPNT + 1 00235 CALL DCOPY( N-CUTPNT, Q( ZPP1, ZPP1 ), LDQ, WORK( IZ+CUTPNT ), 1 ) 00236 * 00237 * Deflate eigenvalues. 00238 * 00239 CALL DLAED2( K, N, CUTPNT, D, Q, LDQ, INDXQ, RHO, WORK( IZ ), 00240 $ WORK( IDLMDA ), WORK( IW ), WORK( IQ2 ), 00241 $ IWORK( INDX ), IWORK( INDXC ), IWORK( INDXP ), 00242 $ IWORK( COLTYP ), INFO ) 00243 * 00244 IF( INFO.NE.0 ) 00245 $ GO TO 20 00246 * 00247 * Solve Secular Equation. 00248 * 00249 IF( K.NE.0 ) THEN 00250 IS = ( IWORK( COLTYP )+IWORK( COLTYP+1 ) )*CUTPNT + 00251 $ ( IWORK( COLTYP+1 )+IWORK( COLTYP+2 ) )*( N-CUTPNT ) + IQ2 00252 CALL DLAED3( K, N, CUTPNT, D, Q, LDQ, RHO, WORK( IDLMDA ), 00253 $ WORK( IQ2 ), IWORK( INDXC ), IWORK( COLTYP ), 00254 $ WORK( IW ), WORK( IS ), INFO ) 00255 IF( INFO.NE.0 ) 00256 $ GO TO 20 00257 * 00258 * Prepare the INDXQ sorting permutation. 00259 * 00260 N1 = K 00261 N2 = N - K 00262 CALL DLAMRG( N1, N2, D, 1, -1, INDXQ ) 00263 ELSE 00264 DO 10 I = 1, N 00265 INDXQ( I ) = I 00266 10 CONTINUE 00267 END IF 00268 * 00269 20 CONTINUE 00270 RETURN 00271 * 00272 * End of DLAED1 00273 * 00274 END