libflame  revision_anchor
Functions
dsytrd.c File Reference

(r)

Functions

int dsytrd_fla (char *uplo, integer *n, doublereal *a, integer *lda, doublereal *d__, doublereal *e, doublereal *tau, doublereal *work, integer *lwork, integer *info)
 

Function Documentation

◆ dsytrd_fla()

int dsytrd_fla ( char *  uplo,
integer n,
doublereal a,
integer lda,
doublereal d__,
doublereal e,
doublereal tau,
doublereal work,
integer lwork,
integer info 
)
194 {
195  /* System generated locals */
196  integer a_dim1, a_offset, i__1, i__2, i__3;
197  /* Local variables */
198  integer i__, j, nb, kk, nx, iws;
199  extern logical lsame_(char *, char *);
200  integer nbmin, iinfo;
201  logical upper;
202  extern /* Subroutine */
203  int dsytd2_fla(char *, integer *, doublereal *, integer *, doublereal *, doublereal *, doublereal *, integer *), dsyr2k_(char *, char *, integer *, integer *, doublereal *, doublereal *, integer *, doublereal *, integer *, doublereal *, doublereal *, integer *), dlatrd_(char *, integer *, integer *, doublereal *, integer *, doublereal *, doublereal *, doublereal *, integer *), xerbla_(char *, integer *);
204  extern integer ilaenv_(integer *, char *, char *, integer *, integer *, integer *, integer *);
205  integer ldwork, lwkopt;
206  logical lquery;
207  /* -- LAPACK computational routine (version 3.4.0) -- */
208  /* -- LAPACK is a software package provided by Univ. of Tennessee, -- */
209  /* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- */
210  /* November 2011 */
211  /* .. Scalar Arguments .. */
212  /* .. */
213  /* .. Array Arguments .. */
214  /* .. */
215  /* ===================================================================== */
216  /* .. Parameters .. */
217  /* .. */
218  /* .. Local Scalars .. */
219  /* .. */
220  /* .. External Subroutines .. */
221  /* .. */
222  /* .. Intrinsic Functions .. */
223  /* .. */
224  /* .. External Functions .. */
225  /* .. */
226  /* .. Executable Statements .. */
227  /* Test the input parameters */
228  /* Parameter adjustments */
229  a_dim1 = *lda;
230  a_offset = 1 + a_dim1;
231  a -= a_offset;
232  --d__;
233  --e;
234  --tau;
235  --work;
236  /* Function Body */
237  *info = 0;
238  upper = lsame_(uplo, "U");
239  lquery = *lwork == -1;
240  if (! upper && ! lsame_(uplo, "L"))
241  {
242  *info = -1;
243  }
244  else if (*n < 0)
245  {
246  *info = -2;
247  }
248  else if (*lda < max(1,*n))
249  {
250  *info = -4;
251  }
252  else if (*lwork < 1 && ! lquery)
253  {
254  *info = -9;
255  }
256  if (*info == 0)
257  {
258  /* Determine the block size. */
259  nb = ilaenv_(&c__1, "DSYTRD", uplo, n, &c_n1, &c_n1, &c_n1);
260  lwkopt = *n * nb;
261  work[1] = (doublereal) lwkopt;
262  }
263  if (*info != 0)
264  {
265  i__1 = -(*info);
266  xerbla_("DSYTRD", &i__1);
267  return 0;
268  }
269  else if (lquery)
270  {
271  return 0;
272  }
273  /* Quick return if possible */
274  if (*n == 0)
275  {
276  work[1] = 1.;
277  return 0;
278  }
279  nx = *n;
280  iws = 1;
281  if (nb > 1 && nb < *n)
282  {
283  /* Determine when to cross over from blocked to unblocked code */
284  /* (last block is always handled by unblocked code). */
285  /* Computing MAX */
286  i__1 = nb;
287  i__2 = ilaenv_(&c__3, "DSYTRD", uplo, n, &c_n1, &c_n1, & c_n1); // , expr subst
288  nx = max(i__1,i__2);
289  if (nx < *n)
290  {
291  /* Determine if workspace is large enough for blocked code. */
292  ldwork = *n;
293  iws = ldwork * nb;
294  if (*lwork < iws)
295  {
296  /* Not enough workspace to use optimal NB: determine the */
297  /* minimum value of NB, and reduce NB or force use of */
298  /* unblocked code by setting NX = N. */
299  /* Computing MAX */
300  i__1 = *lwork / ldwork;
301  nb = max(i__1,1);
302  nbmin = ilaenv_(&c__2, "DSYTRD", uplo, n, &c_n1, &c_n1, &c_n1);
303  if (nb < nbmin)
304  {
305  nx = *n;
306  }
307  }
308  }
309  else
310  {
311  nx = *n;
312  }
313  }
314  else
315  {
316  nb = 1;
317  }
318  if (upper)
319  {
320  /* Reduce the upper triangle of A. */
321  /* Columns 1:kk are handled by the unblocked method. */
322  kk = *n - (*n - nx + nb - 1) / nb * nb;
323  i__1 = kk + 1;
324  i__2 = -nb;
325  for (i__ = *n - nb + 1;
326  i__2 < 0 ? i__ >= i__1 : i__ <= i__1;
327  i__ += i__2)
328  {
329  /* Reduce columns i:i+nb-1 to tridiagonal form and form the */
330  /* matrix W which is needed to update the unreduced part of */
331  /* the matrix */
332  i__3 = i__ + nb - 1;
333  dlatrd_(uplo, &i__3, &nb, &a[a_offset], lda, &e[1], &tau[1], & work[1], &ldwork);
334  /* Update the unreduced submatrix A(1:i-1,1:i-1), using an */
335  /* update of the form: A := A - V*W**T - W*V**T */
336  i__3 = i__ - 1;
337  dsyr2k_(uplo, "No transpose", &i__3, &nb, &c_b22, &a[i__ * a_dim1 + 1], lda, &work[1], &ldwork, &c_b23, &a[a_offset], lda);
338  /* Copy superdiagonal elements back into A, and diagonal */
339  /* elements into D */
340  i__3 = i__ + nb - 1;
341  for (j = i__;
342  j <= i__3;
343  ++j)
344  {
345  a[j - 1 + j * a_dim1] = e[j - 1];
346  d__[j] = a[j + j * a_dim1];
347  /* L10: */
348  }
349  /* L20: */
350  }
351  /* Use unblocked code to reduce the last or only block */
352  dsytd2_fla(uplo, &kk, &a[a_offset], lda, &d__[1], &e[1], &tau[1], &iinfo);
353  }
354  else
355  {
356  /* Reduce the lower triangle of A */
357  i__2 = *n - nx;
358  i__1 = nb;
359  for (i__ = 1;
360  i__1 < 0 ? i__ >= i__2 : i__ <= i__2;
361  i__ += i__1)
362  {
363  /* Reduce columns i:i+nb-1 to tridiagonal form and form the */
364  /* matrix W which is needed to update the unreduced part of */
365  /* the matrix */
366  i__3 = *n - i__ + 1;
367  dlatrd_(uplo, &i__3, &nb, &a[i__ + i__ * a_dim1], lda, &e[i__], & tau[i__], &work[1], &ldwork);
368  /* Update the unreduced submatrix A(i+ib:n,i+ib:n), using */
369  /* an update of the form: A := A - V*W**T - W*V**T */
370  i__3 = *n - i__ - nb + 1;
371  dsyr2k_(uplo, "No transpose", &i__3, &nb, &c_b22, &a[i__ + nb + i__ * a_dim1], lda, &work[nb + 1], &ldwork, &c_b23, &a[ i__ + nb + (i__ + nb) * a_dim1], lda);
372  /* Copy subdiagonal elements back into A, and diagonal */
373  /* elements into D */
374  i__3 = i__ + nb - 1;
375  for (j = i__;
376  j <= i__3;
377  ++j)
378  {
379  a[j + 1 + j * a_dim1] = e[j];
380  d__[j] = a[j + j * a_dim1];
381  /* L30: */
382  }
383  /* L40: */
384  }
385  /* Use unblocked code to reduce the last or only block */
386  i__1 = *n - i__ + 1;
387  dsytd2_fla(uplo, &i__1, &a[i__ + i__ * a_dim1], lda, &d__[i__], &e[i__], &tau[i__], &iinfo);
388  }
389  work[1] = (doublereal) lwkopt;
390  return 0;
391  /* End of DSYTRD */
392 }
double doublereal
Definition: FLA_f2c.h:31
int integer
Definition: FLA_f2c.h:25
int logical
Definition: FLA_f2c.h:36
int dsytd2_fla(char *uplo, integer *n, doublereal *a, integer *lda, doublereal *d__, doublereal *e, doublereal *tau, integer *info)
Definition: dsytd2.c:169

References dsytd2_fla().