LAPACK  3.8.0
LAPACK: Linear Algebra PACKage
ssytrd_2stage.f
Go to the documentation of this file.
1 *> \brief \b SSYTRD_2STAGE
2 *
3 * @generated from zhetrd_2stage.f, fortran z -> s, Sun Nov 6 19:34:06 2016
4 *
5 * =========== DOCUMENTATION ===========
6 *
7 * Online html documentation available at
8 * http://www.netlib.org/lapack/explore-html/
9 *
10 *> \htmlonly
11 *> Download SSYTRD_2STAGE + dependencies
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ssytrd_2stage.f">
13 *> [TGZ]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ssytrd_2stage.f">
15 *> [ZIP]</a>
16 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ssytrd_2stage.f">
17 *> [TXT]</a>
18 *> \endhtmlonly
19 *
20 * Definition:
21 * ===========
22 *
23 * SUBROUTINE SSYTRD_2STAGE( VECT, UPLO, N, A, LDA, D, E, TAU,
24 * HOUS2, LHOUS2, WORK, LWORK, INFO )
25 *
26 * IMPLICIT NONE
27 *
28 * .. Scalar Arguments ..
29 * CHARACTER VECT, UPLO
30 * INTEGER N, LDA, LWORK, LHOUS2, INFO
31 * ..
32 * .. Array Arguments ..
33 * REAL D( * ), E( * )
34 * REAL A( LDA, * ), TAU( * ),
35 * HOUS2( * ), WORK( * )
36 * ..
37 *
38 *
39 *> \par Purpose:
40 * =============
41 *>
42 *> \verbatim
43 *>
44 *> SSYTRD_2STAGE reduces a real symmetric matrix A to real symmetric
45 *> tridiagonal form T by a orthogonal similarity transformation:
46 *> Q1**T Q2**T* A * Q2 * Q1 = T.
47 *> \endverbatim
48 *
49 * Arguments:
50 * ==========
51 *
52 *> \param[in] VECT
53 *> \verbatim
54 *> VECT is CHARACTER*1
55 *> = 'N': No need for the Housholder representation,
56 *> in particular for the second stage (Band to
57 *> tridiagonal) and thus LHOUS2 is of size max(1, 4*N);
58 *> = 'V': the Householder representation is needed to
59 *> either generate Q1 Q2 or to apply Q1 Q2,
60 *> then LHOUS2 is to be queried and computed.
61 *> (NOT AVAILABLE IN THIS RELEASE).
62 *> \endverbatim
63 *>
64 *> \param[in] UPLO
65 *> \verbatim
66 *> UPLO is CHARACTER*1
67 *> = 'U': Upper triangle of A is stored;
68 *> = 'L': Lower triangle of A is stored.
69 *> \endverbatim
70 *>
71 *> \param[in] N
72 *> \verbatim
73 *> N is INTEGER
74 *> The order of the matrix A. N >= 0.
75 *> \endverbatim
76 *>
77 *> \param[in,out] A
78 *> \verbatim
79 *> A is REAL array, dimension (LDA,N)
80 *> On entry, the symmetric matrix A. If UPLO = 'U', the leading
81 *> N-by-N upper triangular part of A contains the upper
82 *> triangular part of the matrix A, and the strictly lower
83 *> triangular part of A is not referenced. If UPLO = 'L', the
84 *> leading N-by-N lower triangular part of A contains the lower
85 *> triangular part of the matrix A, and the strictly upper
86 *> triangular part of A is not referenced.
87 *> On exit, if UPLO = 'U', the band superdiagonal
88 *> of A are overwritten by the corresponding elements of the
89 *> internal band-diagonal matrix AB, and the elements above
90 *> the KD superdiagonal, with the array TAU, represent the orthogonal
91 *> matrix Q1 as a product of elementary reflectors; if UPLO
92 *> = 'L', the diagonal and band subdiagonal of A are over-
93 *> written by the corresponding elements of the internal band-diagonal
94 *> matrix AB, and the elements below the KD subdiagonal, with
95 *> the array TAU, represent the orthogonal matrix Q1 as a product
96 *> of elementary reflectors. See Further Details.
97 *> \endverbatim
98 *>
99 *> \param[in] LDA
100 *> \verbatim
101 *> LDA is INTEGER
102 *> The leading dimension of the array A. LDA >= max(1,N).
103 *> \endverbatim
104 *>
105 *> \param[out] D
106 *> \verbatim
107 *> D is REAL array, dimension (N)
108 *> The diagonal elements of the tridiagonal matrix T.
109 *> \endverbatim
110 *>
111 *> \param[out] E
112 *> \verbatim
113 *> E is REAL array, dimension (N-1)
114 *> The off-diagonal elements of the tridiagonal matrix T.
115 *> \endverbatim
116 *>
117 *> \param[out] TAU
118 *> \verbatim
119 *> TAU is REAL array, dimension (N-KD)
120 *> The scalar factors of the elementary reflectors of
121 *> the first stage (see Further Details).
122 *> \endverbatim
123 *>
124 *> \param[out] HOUS2
125 *> \verbatim
126 *> HOUS2 is REAL array, dimension LHOUS2, that
127 *> store the Householder representation of the stage2
128 *> band to tridiagonal.
129 *> \endverbatim
130 *>
131 *> \param[in] LHOUS2
132 *> \verbatim
133 *> LHOUS2 is INTEGER
134 *> The dimension of the array HOUS2. LHOUS2 = MAX(1, dimension)
135 *> If LWORK = -1, or LHOUS2=-1,
136 *> then a query is assumed; the routine
137 *> only calculates the optimal size of the HOUS2 array, returns
138 *> this value as the first entry of the HOUS2 array, and no error
139 *> message related to LHOUS2 is issued by XERBLA.
140 *> LHOUS2 = MAX(1, dimension) where
141 *> dimension = 4*N if VECT='N'
142 *> not available now if VECT='H'
143 *> \endverbatim
144 *>
145 *> \param[out] WORK
146 *> \verbatim
147 *> WORK is REAL array, dimension (LWORK)
148 *> \endverbatim
149 *>
150 *> \param[in] LWORK
151 *> \verbatim
152 *> LWORK is INTEGER
153 *> The dimension of the array WORK. LWORK = MAX(1, dimension)
154 *> If LWORK = -1, or LHOUS2=-1,
155 *> then a workspace query is assumed; the routine
156 *> only calculates the optimal size of the WORK array, returns
157 *> this value as the first entry of the WORK array, and no error
158 *> message related to LWORK is issued by XERBLA.
159 *> LWORK = MAX(1, dimension) where
160 *> dimension = max(stage1,stage2) + (KD+1)*N
161 *> = N*KD + N*max(KD+1,FACTOPTNB)
162 *> + max(2*KD*KD, KD*NTHREADS)
163 *> + (KD+1)*N
164 *> where KD is the blocking size of the reduction,
165 *> FACTOPTNB is the blocking used by the QR or LQ
166 *> algorithm, usually FACTOPTNB=128 is a good choice
167 *> NTHREADS is the number of threads used when
168 *> openMP compilation is enabled, otherwise =1.
169 *> \endverbatim
170 *>
171 *> \param[out] INFO
172 *> \verbatim
173 *> INFO is INTEGER
174 *> = 0: successful exit
175 *> < 0: if INFO = -i, the i-th argument had an illegal value
176 *> \endverbatim
177 *
178 * Authors:
179 * ========
180 *
181 *> \author Univ. of Tennessee
182 *> \author Univ. of California Berkeley
183 *> \author Univ. of Colorado Denver
184 *> \author NAG Ltd.
185 *
186 *> \date November 2017
187 *
188 *> \ingroup realSYcomputational
189 *
190 *> \par Further Details:
191 * =====================
192 *>
193 *> \verbatim
194 *>
195 *> Implemented by Azzam Haidar.
196 *>
197 *> All details are available on technical report, SC11, SC13 papers.
198 *>
199 *> Azzam Haidar, Hatem Ltaief, and Jack Dongarra.
200 *> Parallel reduction to condensed forms for symmetric eigenvalue problems
201 *> using aggregated fine-grained and memory-aware kernels. In Proceedings
202 *> of 2011 International Conference for High Performance Computing,
203 *> Networking, Storage and Analysis (SC '11), New York, NY, USA,
204 *> Article 8 , 11 pages.
205 *> http://doi.acm.org/10.1145/2063384.2063394
206 *>
207 *> A. Haidar, J. Kurzak, P. Luszczek, 2013.
208 *> An improved parallel singular value algorithm and its implementation
209 *> for multicore hardware, In Proceedings of 2013 International Conference
210 *> for High Performance Computing, Networking, Storage and Analysis (SC '13).
211 *> Denver, Colorado, USA, 2013.
212 *> Article 90, 12 pages.
213 *> http://doi.acm.org/10.1145/2503210.2503292
214 *>
215 *> A. Haidar, R. Solca, S. Tomov, T. Schulthess and J. Dongarra.
216 *> A novel hybrid CPU-GPU generalized eigensolver for electronic structure
217 *> calculations based on fine-grained memory aware tasks.
218 *> International Journal of High Performance Computing Applications.
219 *> Volume 28 Issue 2, Pages 196-209, May 2014.
220 *> http://hpc.sagepub.com/content/28/2/196
221 *>
222 *> \endverbatim
223 *>
224 * =====================================================================
225  SUBROUTINE ssytrd_2stage( VECT, UPLO, N, A, LDA, D, E, TAU,
226  $ HOUS2, LHOUS2, WORK, LWORK, INFO )
227 *
228  IMPLICIT NONE
229 *
230 * -- LAPACK computational routine (version 3.8.0) --
231 * -- LAPACK is a software package provided by Univ. of Tennessee, --
232 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
233 * November 2017
234 *
235 * .. Scalar Arguments ..
236  CHARACTER VECT, UPLO
237  INTEGER N, LDA, LWORK, LHOUS2, INFO
238 * ..
239 * .. Array Arguments ..
240  REAL D( * ), E( * )
241  REAL A( lda, * ), TAU( * ),
242  $ hous2( * ), work( * )
243 * ..
244 *
245 * =====================================================================
246 * ..
247 * .. Local Scalars ..
248  LOGICAL LQUERY, UPPER, WANTQ
249  INTEGER KD, IB, LWMIN, LHMIN, LWRK, LDAB, WPOS, ABPOS
250 * ..
251 * .. External Subroutines ..
252  EXTERNAL xerbla, ssytrd_sy2sb, ssytrd_sb2st
253 * ..
254 * .. External Functions ..
255  LOGICAL LSAME
256  INTEGER ILAENV2STAGE
257  EXTERNAL lsame, ilaenv2stage
258 * ..
259 * .. Executable Statements ..
260 *
261 * Test the input parameters
262 *
263  info = 0
264  wantq = lsame( vect, 'V' )
265  upper = lsame( uplo, 'U' )
266  lquery = ( lwork.EQ.-1 ) .OR. ( lhous2.EQ.-1 )
267 *
268 * Determine the block size, the workspace size and the hous size.
269 *
270  kd = ilaenv2stage( 1, 'SSYTRD_2STAGE', vect, n, -1, -1, -1 )
271  ib = ilaenv2stage( 2, 'SSYTRD_2STAGE', vect, n, kd, -1, -1 )
272  lhmin = ilaenv2stage( 3, 'SSYTRD_2STAGE', vect, n, kd, ib, -1 )
273  lwmin = ilaenv2stage( 4, 'SSYTRD_2STAGE', vect, n, kd, ib, -1 )
274 * WRITE(*,*),'SSYTRD_2STAGE N KD UPLO LHMIN LWMIN ',N, KD, UPLO,
275 * $ LHMIN, LWMIN
276 *
277  IF( .NOT.lsame( vect, 'N' ) ) THEN
278  info = -1
279  ELSE IF( .NOT.upper .AND. .NOT.lsame( uplo, 'L' ) ) THEN
280  info = -2
281  ELSE IF( n.LT.0 ) THEN
282  info = -3
283  ELSE IF( lda.LT.max( 1, n ) ) THEN
284  info = -5
285  ELSE IF( lhous2.LT.lhmin .AND. .NOT.lquery ) THEN
286  info = -10
287  ELSE IF( lwork.LT.lwmin .AND. .NOT.lquery ) THEN
288  info = -12
289  END IF
290 *
291  IF( info.EQ.0 ) THEN
292  hous2( 1 ) = lhmin
293  work( 1 ) = lwmin
294  END IF
295 *
296  IF( info.NE.0 ) THEN
297  CALL xerbla( 'SSYTRD_2STAGE', -info )
298  RETURN
299  ELSE IF( lquery ) THEN
300  RETURN
301  END IF
302 *
303 * Quick return if possible
304 *
305  IF( n.EQ.0 ) THEN
306  work( 1 ) = 1
307  RETURN
308  END IF
309 *
310 * Determine pointer position
311 *
312  ldab = kd+1
313  lwrk = lwork-ldab*n
314  abpos = 1
315  wpos = abpos + ldab*n
316  CALL ssytrd_sy2sb( uplo, n, kd, a, lda, work( abpos ), ldab,
317  $ tau, work( wpos ), lwrk, info )
318  IF( info.NE.0 ) THEN
319  CALL xerbla( 'SSYTRD_SY2SB', -info )
320  RETURN
321  END IF
322  CALL ssytrd_sb2st( 'Y', vect, uplo, n, kd,
323  $ work( abpos ), ldab, d, e,
324  $ hous2, lhous2, work( wpos ), lwrk, info )
325  IF( info.NE.0 ) THEN
326  CALL xerbla( 'SSYTRD_SB2ST', -info )
327  RETURN
328  END IF
329 *
330 *
331  hous2( 1 ) = lhmin
332  work( 1 ) = lwmin
333  RETURN
334 *
335 * End of SSYTRD_2STAGE
336 *
337  END
subroutine ssytrd_2stage(VECT, UPLO, N, A, LDA, D, E, TAU, HOUS2, LHOUS2, WORK, LWORK, INFO)
SSYTRD_2STAGE
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine ssytrd_sy2sb(UPLO, N, KD, A, LDA, AB, LDAB, TAU, WORK, LWORK, INFO)
SSYTRD_SY2SB
Definition: ssytrd_sy2sb.f:245