/src/igraph/vendor/lapack/dlaexc.c

Source (jump to first uncovered line)
/*  -- translated by f2c (version 20191129).
   You must link the resulting object file with libf2c:
  on Microsoft Windows system, link with libf2c.lib;
  on Linux or Unix systems, link with .../path/to/libf2c.a -lm
  or, if you install libf2c.a in a standard place, with -lf2c -lm
  -- in that order, at the end of the command line, as in
    cc *.o -lf2c -lm
  Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,

    http://www.netlib.org/f2c/libf2c.zip
*/

#include "f2c.h"

/* Table of constant values */

static integer c__1 = 1;
static integer c__4 = 4;
static logical c_false = FALSE_;
static integer c_n1 = -1;
static integer c__2 = 2;
static integer c__3 = 3;

/* > \brief \b DLAEXC swaps adjacent diagonal blocks of a real upper quasi-triangular matrix in Schur canonica
l form, by an orthogonal similarity transformation.   

    =========== DOCUMENTATION ===========   

   Online html documentation available at   
              http://www.netlib.org/lapack/explore-html/   

   > \htmlonly   
   > Download DLAEXC + dependencies   
   > <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaexc.
f">   
   > [TGZ]</a>   
   > <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaexc.
f">   
   > [ZIP]</a>   
   > <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaexc.
f">   
   > [TXT]</a>   
   > \endhtmlonly   

    Definition:   
    ===========   

         SUBROUTINE DLAEXC( WANTQ, N, T, LDT, Q, LDQ, J1, N1, N2, WORK,   
                            INFO )   

         LOGICAL            WANTQ   
         INTEGER            INFO, J1, LDQ, LDT, N, N1, N2   
         DOUBLE PRECISION   Q( LDQ, * ), T( LDT, * ), WORK( * )   


   > \par Purpose:   
    =============   
   >   
   > \verbatim   
   >   
   > DLAEXC swaps adjacent diagonal blocks T11 and T22 of order 1 or 2 in   
   > an upper quasi-triangular matrix T by an orthogonal similarity   
   > transformation.   
   >   
   > T must be in Schur canonical form, that is, block upper triangular   
   > with 1-by-1 and 2-by-2 diagonal blocks; each 2-by-2 diagonal block   
   > has its diagonal elemnts equal and its off-diagonal elements of   
   > opposite sign.   
   > \endverbatim   

    Arguments:   
    ==========   

   > \param[in] WANTQ   
   > \verbatim   
   >          WANTQ is LOGICAL   
   >          = .TRUE. : accumulate the transformation in the matrix Q;   
   >          = .FALSE.: do not accumulate the transformation.   
   > \endverbatim   
   >   
   > \param[in] N   
   > \verbatim   
   >          N is INTEGER   
   >          The order of the matrix T. N >= 0.   
   > \endverbatim   
   >   
   > \param[in,out] T   
   > \verbatim   
   >          T is DOUBLE PRECISION array, dimension (LDT,N)   
   >          On entry, the upper quasi-triangular matrix T, in Schur   
   >          canonical form.   
   >          On exit, the updated matrix T, again in Schur canonical form.   
   > \endverbatim   
   >   
   > \param[in] LDT   
   > \verbatim   
   >          LDT is INTEGER   
   >          The leading dimension of the array T. LDT >= max(1,N).   
   > \endverbatim   
   >   
   > \param[in,out] Q   
   > \verbatim   
   >          Q is DOUBLE PRECISION array, dimension (LDQ,N)   
   >          On entry, if WANTQ is .TRUE., the orthogonal matrix Q.   
   >          On exit, if WANTQ is .TRUE., the updated matrix Q.   
   >          If WANTQ is .FALSE., Q is not referenced.   
   > \endverbatim   
   >   
   > \param[in] LDQ   
   > \verbatim   
   >          LDQ is INTEGER   
   >          The leading dimension of the array Q.   
   >          LDQ >= 1; and if WANTQ is .TRUE., LDQ >= N.   
   > \endverbatim   
   >   
   > \param[in] J1   
   > \verbatim   
   >          J1 is INTEGER   
   >          The index of the first row of the first block T11.   
   > \endverbatim   
   >   
   > \param[in] N1   
   > \verbatim   
   >          N1 is INTEGER   
   >          The order of the first block T11. N1 = 0, 1 or 2.   
   > \endverbatim   
   >   
   > \param[in] N2   
   > \verbatim   
   >          N2 is INTEGER   
   >          The order of the second block T22. N2 = 0, 1 or 2.   
   > \endverbatim   
   >   
   > \param[out] WORK   
   > \verbatim   
   >          WORK is DOUBLE PRECISION array, dimension (N)   
   > \endverbatim   
   >   
   > \param[out] INFO   
   > \verbatim   
   >          INFO is INTEGER   
   >          = 0: successful exit   
   >          = 1: the transformed matrix T would be too far from Schur   
   >               form; the blocks are not swapped and T and Q are   
   >               unchanged.   
   > \endverbatim   

    Authors:   
    ========   

   > \author Univ. of Tennessee   
   > \author Univ. of California Berkeley   
   > \author Univ. of Colorado Denver   
   > \author NAG Ltd.   

   > \date September 2012   

   > \ingroup doubleOTHERauxiliary   

    =====================================================================   
   Subroutine */ int igraphdlaexc_(logical *wantq, integer *n, doublereal *t, 
  integer *ldt, doublereal *q, integer *ldq, integer *j1, integer *n1, 
  integer *n2, doublereal *work, integer *info)
{
    /* System generated locals */
    integer q_dim1, q_offset, t_dim1, t_offset, i__1;
    doublereal d__1, d__2, d__3;

    /* Local variables */
    doublereal d__[16]  /* was [4][4] */;
    integer k;
    doublereal u[3], x[4] /* was [2][2] */;
    integer j2, j3, j4;
    doublereal u1[3], u2[3];
    integer nd;
    doublereal cs, t11, t22, t33, sn, wi1, wi2, wr1, wr2, eps, tau, tau1, 
      tau2;
    integer ierr;
    doublereal temp;
    extern /* Subroutine */ int igraphdrot_(integer *, doublereal *, integer *, 
      doublereal *, integer *, doublereal *, doublereal *);
    doublereal scale, dnorm, xnorm;
    extern /* Subroutine */ int igraphdlanv2_(doublereal *, doublereal *, 
      doublereal *, doublereal *, doublereal *, doublereal *, 
      doublereal *, doublereal *, doublereal *, doublereal *), igraphdlasy2_(
      logical *, logical *, integer *, integer *, integer *, doublereal 
      *, integer *, doublereal *, integer *, doublereal *, integer *, 
      doublereal *, doublereal *, integer *, doublereal *, integer *);
    extern doublereal igraphdlamch_(char *), igraphdlange_(char *, integer *, 
      integer *, doublereal *, integer *, doublereal *);
    extern /* Subroutine */ int igraphdlarfg_(integer *, doublereal *, doublereal *,
       integer *, doublereal *), igraphdlacpy_(char *, integer *, integer *, 
      doublereal *, integer *, doublereal *, integer *), 
      igraphdlartg_(doublereal *, doublereal *, doublereal *, doublereal *, 
      doublereal *), igraphdlarfx_(char *, integer *, integer *, doublereal *,
       doublereal *, doublereal *, integer *, doublereal *);
    doublereal thresh, smlnum;


/*  -- LAPACK auxiliary routine (version 3.4.2) --   
    -- LAPACK is a software package provided by Univ. of Tennessee,    --   
    -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--   
       September 2012   


    =====================================================================   


       Parameter adjustments */
    t_dim1 = *ldt;
    t_offset = 1 + t_dim1;
    t -= t_offset;
    q_dim1 = *ldq;
    q_offset = 1 + q_dim1;
    q -= q_offset;
    --work;

    /* Function Body */
    *info = 0;

/*     Quick return if possible */

    if (*n == 0 || *n1 == 0 || *n2 == 0) {
  return 0;
    }
    if (*j1 + *n1 > *n) {
  return 0;
    }

    j2 = *j1 + 1;
    j3 = *j1 + 2;
    j4 = *j1 + 3;

    if (*n1 == 1 && *n2 == 1) {

/*        Swap two 1-by-1 blocks. */

  t11 = t[*j1 + *j1 * t_dim1];
  t22 = t[j2 + j2 * t_dim1];

/*        Determine the transformation to perform the interchange. */

  d__1 = t22 - t11;
  igraphdlartg_(&t[*j1 + j2 * t_dim1], &d__1, &cs, &sn, &temp);

/*        Apply transformation to the matrix T. */

  if (j3 <= *n) {
      i__1 = *n - *j1 - 1;
      igraphdrot_(&i__1, &t[*j1 + j3 * t_dim1], ldt, &t[j2 + j3 * t_dim1], 
        ldt, &cs, &sn);
  }
  i__1 = *j1 - 1;
  igraphdrot_(&i__1, &t[*j1 * t_dim1 + 1], &c__1, &t[j2 * t_dim1 + 1], &c__1, 
    &cs, &sn);

  t[*j1 + *j1 * t_dim1] = t22;
  t[j2 + j2 * t_dim1] = t11;

  if (*wantq) {

/*           Accumulate transformation in the matrix Q. */

      igraphdrot_(n, &q[*j1 * q_dim1 + 1], &c__1, &q[j2 * q_dim1 + 1], &c__1, 
        &cs, &sn);
  }

    } else {

/*        Swapping involves at least one 2-by-2 block.   

          Copy the diagonal block of order N1+N2 to the local array D   
          and compute its norm. */

  nd = *n1 + *n2;
  igraphdlacpy_("Full", &nd, &nd, &t[*j1 + *j1 * t_dim1], ldt, d__, &c__4);
  dnorm = igraphdlange_("Max", &nd, &nd, d__, &c__4, &work[1]);

/*        Compute machine-dependent threshold for test for accepting   
          swap. */

  eps = igraphdlamch_("P");
  smlnum = igraphdlamch_("S") / eps;
/* Computing MAX */
  d__1 = eps * 10. * dnorm;
  thresh = max(d__1,smlnum);

/*        Solve T11*X - X*T22 = scale*T12 for X. */

  igraphdlasy2_(&c_false, &c_false, &c_n1, n1, n2, d__, &c__4, &d__[*n1 + 1 + 
    (*n1 + 1 << 2) - 5], &c__4, &d__[(*n1 + 1 << 2) - 4], &c__4, &
    scale, x, &c__2, &xnorm, &ierr);

/*        Swap the adjacent diagonal blocks. */

  k = *n1 + *n1 + *n2 - 3;
  switch (k) {
      case 1:  goto L10;
      case 2:  goto L20;
      case 3:  goto L30;
  }

L10:

/*        N1 = 1, N2 = 2: generate elementary reflector H so that:   

          ( scale, X11, X12 ) H = ( 0, 0, * ) */

  u[0] = scale;
  u[1] = x[0];
  u[2] = x[2];
  igraphdlarfg_(&c__3, &u[2], u, &c__1, &tau);
  u[2] = 1.;
  t11 = t[*j1 + *j1 * t_dim1];

/*        Perform swap provisionally on diagonal block in D. */

  igraphdlarfx_("L", &c__3, &c__3, u, &tau, d__, &c__4, &work[1]);
  igraphdlarfx_("R", &c__3, &c__3, u, &tau, d__, &c__4, &work[1]);

/*        Test whether to reject swap.   

   Computing MAX */
  d__2 = abs(d__[2]), d__3 = abs(d__[6]), d__2 = max(d__2,d__3), d__3 = 
    (d__1 = d__[10] - t11, abs(d__1));
  if (max(d__2,d__3) > thresh) {
      goto L50;
  }

/*        Accept swap: apply transformation to the entire matrix T. */

  i__1 = *n - *j1 + 1;
  igraphdlarfx_("L", &c__3, &i__1, u, &tau, &t[*j1 + *j1 * t_dim1], ldt, &
    work[1]);
  igraphdlarfx_("R", &j2, &c__3, u, &tau, &t[*j1 * t_dim1 + 1], ldt, &work[1]);

  t[j3 + *j1 * t_dim1] = 0.;
  t[j3 + j2 * t_dim1] = 0.;
  t[j3 + j3 * t_dim1] = t11;

  if (*wantq) {

/*           Accumulate transformation in the matrix Q. */

      igraphdlarfx_("R", n, &c__3, u, &tau, &q[*j1 * q_dim1 + 1], ldq, &work[
        1]);
  }
  goto L40;

L20:

/*        N1 = 2, N2 = 1: generate elementary reflector H so that:   

          H (  -X11 ) = ( * )   
            (  -X21 ) = ( 0 )   
            ( scale ) = ( 0 ) */

  u[0] = -x[0];
  u[1] = -x[1];
  u[2] = scale;
  igraphdlarfg_(&c__3, u, &u[1], &c__1, &tau);
  u[0] = 1.;
  t33 = t[j3 + j3 * t_dim1];

/*        Perform swap provisionally on diagonal block in D. */

  igraphdlarfx_("L", &c__3, &c__3, u, &tau, d__, &c__4, &work[1]);
  igraphdlarfx_("R", &c__3, &c__3, u, &tau, d__, &c__4, &work[1]);

/*        Test whether to reject swap.   

   Computing MAX */
  d__2 = abs(d__[1]), d__3 = abs(d__[2]), d__2 = max(d__2,d__3), d__3 = 
    (d__1 = d__[0] - t33, abs(d__1));
  if (max(d__2,d__3) > thresh) {
      goto L50;
  }

/*        Accept swap: apply transformation to the entire matrix T. */

  igraphdlarfx_("R", &j3, &c__3, u, &tau, &t[*j1 * t_dim1 + 1], ldt, &work[1]);
  i__1 = *n - *j1;
  igraphdlarfx_("L", &c__3, &i__1, u, &tau, &t[*j1 + j2 * t_dim1], ldt, &work[
    1]);

  t[*j1 + *j1 * t_dim1] = t33;
  t[j2 + *j1 * t_dim1] = 0.;
  t[j3 + *j1 * t_dim1] = 0.;

  if (*wantq) {

/*           Accumulate transformation in the matrix Q. */

      igraphdlarfx_("R", n, &c__3, u, &tau, &q[*j1 * q_dim1 + 1], ldq, &work[
        1]);
  }
  goto L40;

L30:

/*        N1 = 2, N2 = 2: generate elementary reflectors H(1) and H(2) so   
          that:   

          H(2) H(1) (  -X11  -X12 ) = (  *  * )   
                    (  -X21  -X22 )   (  0  * )   
                    ( scale    0  )   (  0  0 )   
                    (    0  scale )   (  0  0 ) */

  u1[0] = -x[0];
  u1[1] = -x[1];
  u1[2] = scale;
  igraphdlarfg_(&c__3, u1, &u1[1], &c__1, &tau1);
  u1[0] = 1.;

  temp = -tau1 * (x[2] + u1[1] * x[3]);
  u2[0] = -temp * u1[1] - x[3];
  u2[1] = -temp * u1[2];
  u2[2] = scale;
  igraphdlarfg_(&c__3, u2, &u2[1], &c__1, &tau2);
  u2[0] = 1.;

/*        Perform swap provisionally on diagonal block in D. */

  igraphdlarfx_("L", &c__3, &c__4, u1, &tau1, d__, &c__4, &work[1])
    ;
  igraphdlarfx_("R", &c__4, &c__3, u1, &tau1, d__, &c__4, &work[1])
    ;
  igraphdlarfx_("L", &c__3, &c__4, u2, &tau2, &d__[1], &c__4, &work[1]);
  igraphdlarfx_("R", &c__4, &c__3, u2, &tau2, &d__[4], &c__4, &work[1]);

/*        Test whether to reject swap.   

   Computing MAX */
  d__1 = abs(d__[2]), d__2 = abs(d__[6]), d__1 = max(d__1,d__2), d__2 = 
    abs(d__[3]), d__1 = max(d__1,d__2), d__2 = abs(d__[7]);
  if (max(d__1,d__2) > thresh) {
      goto L50;
  }

/*        Accept swap: apply transformation to the entire matrix T. */

  i__1 = *n - *j1 + 1;
  igraphdlarfx_("L", &c__3, &i__1, u1, &tau1, &t[*j1 + *j1 * t_dim1], ldt, &
    work[1]);
  igraphdlarfx_("R", &j4, &c__3, u1, &tau1, &t[*j1 * t_dim1 + 1], ldt, &work[
    1]);
  i__1 = *n - *j1 + 1;
  igraphdlarfx_("L", &c__3, &i__1, u2, &tau2, &t[j2 + *j1 * t_dim1], ldt, &
    work[1]);
  igraphdlarfx_("R", &j4, &c__3, u2, &tau2, &t[j2 * t_dim1 + 1], ldt, &work[1]
    );

  t[j3 + *j1 * t_dim1] = 0.;
  t[j3 + j2 * t_dim1] = 0.;
  t[j4 + *j1 * t_dim1] = 0.;
  t[j4 + j2 * t_dim1] = 0.;

  if (*wantq) {

/*           Accumulate transformation in the matrix Q. */

      igraphdlarfx_("R", n, &c__3, u1, &tau1, &q[*j1 * q_dim1 + 1], ldq, &
        work[1]);
      igraphdlarfx_("R", n, &c__3, u2, &tau2, &q[j2 * q_dim1 + 1], ldq, &work[
        1]);
  }

L40:

  if (*n2 == 2) {

/*           Standardize new 2-by-2 block T11 */

      igraphdlanv2_(&t[*j1 + *j1 * t_dim1], &t[*j1 + j2 * t_dim1], &t[j2 + *
        j1 * t_dim1], &t[j2 + j2 * t_dim1], &wr1, &wi1, &wr2, &
        wi2, &cs, &sn);
      i__1 = *n - *j1 - 1;
      igraphdrot_(&i__1, &t[*j1 + (*j1 + 2) * t_dim1], ldt, &t[j2 + (*j1 + 2) 
        * t_dim1], ldt, &cs, &sn);
      i__1 = *j1 - 1;
      igraphdrot_(&i__1, &t[*j1 * t_dim1 + 1], &c__1, &t[j2 * t_dim1 + 1], &
        c__1, &cs, &sn);
      if (*wantq) {
    igraphdrot_(n, &q[*j1 * q_dim1 + 1], &c__1, &q[j2 * q_dim1 + 1], &
      c__1, &cs, &sn);
      }
  }

  if (*n1 == 2) {

/*           Standardize new 2-by-2 block T22 */

      j3 = *j1 + *n2;
      j4 = j3 + 1;
      igraphdlanv2_(&t[j3 + j3 * t_dim1], &t[j3 + j4 * t_dim1], &t[j4 + j3 * 
        t_dim1], &t[j4 + j4 * t_dim1], &wr1, &wi1, &wr2, &wi2, &
        cs, &sn);
      if (j3 + 2 <= *n) {
    i__1 = *n - j3 - 1;
    igraphdrot_(&i__1, &t[j3 + (j3 + 2) * t_dim1], ldt, &t[j4 + (j3 + 2)
       * t_dim1], ldt, &cs, &sn);
      }
      i__1 = j3 - 1;
      igraphdrot_(&i__1, &t[j3 * t_dim1 + 1], &c__1, &t[j4 * t_dim1 + 1], &
        c__1, &cs, &sn);
      if (*wantq) {
    igraphdrot_(n, &q[j3 * q_dim1 + 1], &c__1, &q[j4 * q_dim1 + 1], &
      c__1, &cs, &sn);
      }
  }

    }
    return 0;

/*     Exit with INFO = 1 if swap was rejected. */

L50:
    *info = 1;
    return 0;

/*     End of DLAEXC */

} /* igraphdlaexc_ */


Coverage Report

Created: 2023-06-07 06:06

Line	Count	Source (jump to first uncovered line)
1		/* -- translated by f2c (version 20191129).
2		You must link the resulting object file with libf2c:
3		on Microsoft Windows system, link with libf2c.lib;
4		on Linux or Unix systems, link with .../path/to/libf2c.a -lm
5		or, if you install libf2c.a in a standard place, with -lf2c -lm
6		-- in that order, at the end of the command line, as in
7		cc *.o -lf2c -lm
8		Source for libf2c is in /netlib/f2c/libf2c.zip, e.g.,
9
10		http://www.netlib.org/f2c/libf2c.zip
11		*/
12
13		#include "f2c.h"
14
15		/* Table of constant values */
16
17		static integer c__1 = 1;
18		static integer c__4 = 4;
19		static logical c_false = FALSE_;
20		static integer c_n1 = -1;
21		static integer c__2 = 2;
22		static integer c__3 = 3;
23
24		/* > \brief \b DLAEXC swaps adjacent diagonal blocks of a real upper quasi-triangular matrix in Schur canonica
25		l form, by an orthogonal similarity transformation.
26
27		=========== DOCUMENTATION ===========
28
29		Online html documentation available at
30		http://www.netlib.org/lapack/explore-html/
31
32		> \htmlonly
33		> Download DLAEXC + dependencies
34		> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaexc.
35		f">
36		> [TGZ]</a>
37		> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaexc.
38		f">
39		> [ZIP]</a>
40		> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaexc.
41		f">
42		> [TXT]</a>
43		> \endhtmlonly
44
45		Definition:
46		===========
47
48		SUBROUTINE DLAEXC( WANTQ, N, T, LDT, Q, LDQ, J1, N1, N2, WORK,
49		INFO )
50
51		LOGICAL WANTQ
52		INTEGER INFO, J1, LDQ, LDT, N, N1, N2
53		DOUBLE PRECISION Q( LDQ, * ), T( LDT, * ), WORK( * )
54
55
56		> \par Purpose:
57		=============
58		>
59		> \verbatim
60		>
61		> DLAEXC swaps adjacent diagonal blocks T11 and T22 of order 1 or 2 in
62		> an upper quasi-triangular matrix T by an orthogonal similarity
63		> transformation.
64		>
65		> T must be in Schur canonical form, that is, block upper triangular
66		> with 1-by-1 and 2-by-2 diagonal blocks; each 2-by-2 diagonal block
67		> has its diagonal elemnts equal and its off-diagonal elements of
68		> opposite sign.
69		> \endverbatim
70
71		Arguments:
72		==========
73
74		> \param[in] WANTQ
75		> \verbatim
76		> WANTQ is LOGICAL
77		> = .TRUE. : accumulate the transformation in the matrix Q;
78		> = .FALSE.: do not accumulate the transformation.
79		> \endverbatim
80		>
81		> \param[in] N
82		> \verbatim
83		> N is INTEGER
84		> The order of the matrix T. N >= 0.
85		> \endverbatim
86		>
87		> \param[in,out] T
88		> \verbatim
89		> T is DOUBLE PRECISION array, dimension (LDT,N)
90		> On entry, the upper quasi-triangular matrix T, in Schur
91		> canonical form.
92		> On exit, the updated matrix T, again in Schur canonical form.
93		> \endverbatim
94		>
95		> \param[in] LDT
96		> \verbatim
97		> LDT is INTEGER
98		> The leading dimension of the array T. LDT >= max(1,N).
99		> \endverbatim
100		>
101		> \param[in,out] Q
102		> \verbatim
103		> Q is DOUBLE PRECISION array, dimension (LDQ,N)
104		> On entry, if WANTQ is .TRUE., the orthogonal matrix Q.
105		> On exit, if WANTQ is .TRUE., the updated matrix Q.
106		> If WANTQ is .FALSE., Q is not referenced.
107		> \endverbatim
108		>
109		> \param[in] LDQ
110		> \verbatim
111		> LDQ is INTEGER
112		> The leading dimension of the array Q.
113		> LDQ >= 1; and if WANTQ is .TRUE., LDQ >= N.
114		> \endverbatim
115		>
116		> \param[in] J1
117		> \verbatim
118		> J1 is INTEGER
119		> The index of the first row of the first block T11.
120		> \endverbatim
121		>
122		> \param[in] N1
123		> \verbatim
124		> N1 is INTEGER
125		> The order of the first block T11. N1 = 0, 1 or 2.
126		> \endverbatim
127		>
128		> \param[in] N2
129		> \verbatim
130		> N2 is INTEGER
131		> The order of the second block T22. N2 = 0, 1 or 2.
132		> \endverbatim
133		>
134		> \param[out] WORK
135		> \verbatim
136		> WORK is DOUBLE PRECISION array, dimension (N)
137		> \endverbatim
138		>
139		> \param[out] INFO
140		> \verbatim
141		> INFO is INTEGER
142		> = 0: successful exit
143		> = 1: the transformed matrix T would be too far from Schur
144		> form; the blocks are not swapped and T and Q are
145		> unchanged.
146		> \endverbatim
147
148		Authors:
149		========
150
151		> \author Univ. of Tennessee
152		> \author Univ. of California Berkeley
153		> \author Univ. of Colorado Denver
154		> \author NAG Ltd.
155
156		> \date September 2012
157
158		> \ingroup doubleOTHERauxiliary
159
160		=====================================================================
161		Subroutine / int igraphdlaexc_(logical wantq, integer n, doublereal t,
162		integer ldt, doublereal q, integer ldq, integer j1, integer *n1,
163		integer n2, doublereal work, integer *info)
164	0	{
165		/* System generated locals */
166	0	integer q_dim1, q_offset, t_dim1, t_offset, i__1;
167	0	doublereal d__1, d__2, d__3;
168
169		/* Local variables */
170	0	doublereal d__[16] /* was [4][4] */;
171	0	integer k;
172	0	doublereal u[3], x[4] /* was [2][2] */;
173	0	integer j2, j3, j4;
174	0	doublereal u1[3], u2[3];
175	0	integer nd;
176	0	doublereal cs, t11, t22, t33, sn, wi1, wi2, wr1, wr2, eps, tau, tau1,
177	0	tau2;
178	0	integer ierr;
179	0	doublereal temp;
180	0	extern /* Subroutine / int igraphdrot_(integer , doublereal , integer ,
181	0	doublereal , integer , doublereal , doublereal );
182	0	doublereal scale, dnorm, xnorm;
183	0	extern /* Subroutine / int igraphdlanv2_(doublereal , doublereal *,
184	0	doublereal , doublereal , doublereal , doublereal ,
185	0	doublereal , doublereal , doublereal , doublereal ), igraphdlasy2_(
186	0	logical , logical , integer , integer , integer *, doublereal
187	0	, integer , doublereal , integer , doublereal , integer ,
188	0	doublereal , doublereal , integer , doublereal , integer *);
189	0	extern doublereal igraphdlamch_(char ), igraphdlange_(char , integer *,
190	0	integer , doublereal , integer , doublereal );
191	0	extern /* Subroutine / int igraphdlarfg_(integer , doublereal , doublereal ,
192	0	integer , doublereal ), igraphdlacpy_(char , integer , integer *,
193	0	doublereal , integer , doublereal , integer ),
194	0	igraphdlartg_(doublereal , doublereal , doublereal , doublereal ,
195	0	doublereal ), igraphdlarfx_(char , integer , integer , doublereal *,
196	0	doublereal , doublereal , integer , doublereal );
197	0	doublereal thresh, smlnum;
198
199
200		/* -- LAPACK auxiliary routine (version 3.4.2) --
201		-- LAPACK is a software package provided by Univ. of Tennessee, --
202		-- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
203		September 2012
204
205
206		=====================================================================
207
208
209		Parameter adjustments */
210	0	t_dim1 = *ldt;
211	0	t_offset = 1 + t_dim1;
212	0	t -= t_offset;
213	0	q_dim1 = *ldq;
214	0	q_offset = 1 + q_dim1;
215	0	q -= q_offset;
216	0	--work;
217
218		/* Function Body */
219	0	*info = 0;
220
221		/* Quick return if possible */
222
223	0	if (n == 0 \|\| n1 == 0 \|\| *n2 == 0) {
224	0	return 0;
225	0	}
226	0	if (j1 + n1 > *n) {
227	0	return 0;
228	0	}
229
230	0	j2 = *j1 + 1;
231	0	j3 = *j1 + 2;
232	0	j4 = *j1 + 3;
233
234	0	if (n1 == 1 && n2 == 1) {
235
236		/* Swap two 1-by-1 blocks. */
237
238	0	t11 = t[j1 + j1 * t_dim1];
239	0	t22 = t[j2 + j2 * t_dim1];
240
241		/* Determine the transformation to perform the interchange. */
242
243	0	d__1 = t22 - t11;
244	0	igraphdlartg_(&t[j1 + j2 t_dim1], &d__1, &cs, &sn, &temp);
245
246		/* Apply transformation to the matrix T. */
247
248	0	if (j3 <= *n) {
249	0	i__1 = n - j1 - 1;
250	0	igraphdrot_(&i__1, &t[j1 + j3 t_dim1], ldt, &t[j2 + j3 * t_dim1],
251	0	ldt, &cs, &sn);
252	0	}
253	0	i__1 = *j1 - 1;
254	0	igraphdrot_(&i__1, &t[j1 t_dim1 + 1], &c__1, &t[j2 * t_dim1 + 1], &c__1,
255	0	&cs, &sn);
256
257	0	t[j1 + j1 * t_dim1] = t22;
258	0	t[j2 + j2 * t_dim1] = t11;
259
260	0	if (*wantq) {
261
262		/* Accumulate transformation in the matrix Q. */
263
264	0	igraphdrot_(n, &q[j1 q_dim1 + 1], &c__1, &q[j2 * q_dim1 + 1], &c__1,
265	0	&cs, &sn);
266	0	}
267
268	0	} else {
269
270		/* Swapping involves at least one 2-by-2 block.
271
272		Copy the diagonal block of order N1+N2 to the local array D
273		and compute its norm. */
274
275	0	nd = n1 + n2;
276	0	igraphdlacpy_("Full", &nd, &nd, &t[j1 + j1 * t_dim1], ldt, d__, &c__4);
277	0	dnorm = igraphdlange_("Max", &nd, &nd, d__, &c__4, &work[1]);
278
279		/* Compute machine-dependent threshold for test for accepting
280		swap. */
281
282	0	eps = igraphdlamch_("P");
283	0	smlnum = igraphdlamch_("S") / eps;
284		/* Computing MAX */
285	0	d__1 = eps * 10. * dnorm;
286	0	thresh = max(d__1,smlnum);
287
288		/* Solve T11X - XT22 = scaleT12 for X. /
289
290	0	igraphdlasy2_(&c_false, &c_false, &c_n1, n1, n2, d__, &c__4, &d__[*n1 + 1 +
291	0	(n1 + 1 << 2) - 5], &c__4, &d__[(n1 + 1 << 2) - 4], &c__4, &
292	0	scale, x, &c__2, &xnorm, &ierr);
293
294		/* Swap the adjacent diagonal blocks. */
295
296	0	k = n1 + n1 + *n2 - 3;
297	0	switch (k) {
298	0	case 1: goto L10;
299	0	case 2: goto L20;
300	0	case 3: goto L30;
301	0	}
302
303	0	L10:
304
305		/* N1 = 1, N2 = 2: generate elementary reflector H so that:
306
307		( scale, X11, X12 ) H = ( 0, 0, * ) */
308
309	0	u[0] = scale;
310	0	u[1] = x[0];
311	0	u[2] = x[2];
312	0	igraphdlarfg_(&c__3, &u[2], u, &c__1, &tau);
313	0	u[2] = 1.;
314	0	t11 = t[j1 + j1 * t_dim1];
315
316		/* Perform swap provisionally on diagonal block in D. */
317
318	0	igraphdlarfx_("L", &c__3, &c__3, u, &tau, d__, &c__4, &work[1]);
319	0	igraphdlarfx_("R", &c__3, &c__3, u, &tau, d__, &c__4, &work[1]);
320
321		/* Test whether to reject swap.
322
323		Computing MAX */
324	0	d__2 = abs(d__[2]), d__3 = abs(d__[6]), d__2 = max(d__2,d__3), d__3 =
325	0	(d__1 = d__[10] - t11, abs(d__1));
326	0	if (max(d__2,d__3) > thresh) {
327	0	goto L50;
328	0	}
329
330		/* Accept swap: apply transformation to the entire matrix T. */
331
332	0	i__1 = n - j1 + 1;
333	0	igraphdlarfx_("L", &c__3, &i__1, u, &tau, &t[j1 + j1 * t_dim1], ldt, &
334	0	work[1]);
335	0	igraphdlarfx_("R", &j2, &c__3, u, &tau, &t[j1 t_dim1 + 1], ldt, &work[1]);
336
337	0	t[j3 + j1 t_dim1] = 0.;
338	0	t[j3 + j2 * t_dim1] = 0.;
339	0	t[j3 + j3 * t_dim1] = t11;
340
341	0	if (*wantq) {
342
343		/* Accumulate transformation in the matrix Q. */
344
345	0	igraphdlarfx_("R", n, &c__3, u, &tau, &q[j1 q_dim1 + 1], ldq, &work[
346	0	1]);
347	0	}
348	0	goto L40;
349
350	0	L20:
351
352		/* N1 = 2, N2 = 1: generate elementary reflector H so that:
353
354		H ( -X11 ) = ( * )
355		( -X21 ) = ( 0 )
356		( scale ) = ( 0 ) */
357
358	0	u[0] = -x[0];
359	0	u[1] = -x[1];
360	0	u[2] = scale;
361	0	igraphdlarfg_(&c__3, u, &u[1], &c__1, &tau);
362	0	u[0] = 1.;
363	0	t33 = t[j3 + j3 * t_dim1];
364
365		/* Perform swap provisionally on diagonal block in D. */
366
367	0	igraphdlarfx_("L", &c__3, &c__3, u, &tau, d__, &c__4, &work[1]);
368	0	igraphdlarfx_("R", &c__3, &c__3, u, &tau, d__, &c__4, &work[1]);
369
370		/* Test whether to reject swap.
371
372		Computing MAX */
373	0	d__2 = abs(d__[1]), d__3 = abs(d__[2]), d__2 = max(d__2,d__3), d__3 =
374	0	(d__1 = d__[0] - t33, abs(d__1));
375	0	if (max(d__2,d__3) > thresh) {
376	0	goto L50;
377	0	}
378
379		/* Accept swap: apply transformation to the entire matrix T. */
380
381	0	igraphdlarfx_("R", &j3, &c__3, u, &tau, &t[j1 t_dim1 + 1], ldt, &work[1]);
382	0	i__1 = n - j1;
383	0	igraphdlarfx_("L", &c__3, &i__1, u, &tau, &t[j1 + j2 t_dim1], ldt, &work[
384	0	1]);
385
386	0	t[j1 + j1 * t_dim1] = t33;
387	0	t[j2 + j1 t_dim1] = 0.;
388	0	t[j3 + j1 t_dim1] = 0.;
389
390	0	if (*wantq) {
391
392		/* Accumulate transformation in the matrix Q. */
393
394	0	igraphdlarfx_("R", n, &c__3, u, &tau, &q[j1 q_dim1 + 1], ldq, &work[
395	0	1]);
396	0	}
397	0	goto L40;
398
399	0	L30:
400
401		/* N1 = 2, N2 = 2: generate elementary reflectors H(1) and H(2) so
402		that:
403
404		H(2) H(1) ( -X11 -X12 ) = ( * * )
405		( -X21 -X22 ) ( 0 * )
406		( scale 0 ) ( 0 0 )
407		( 0 scale ) ( 0 0 ) */
408
409	0	u1[0] = -x[0];
410	0	u1[1] = -x[1];
411	0	u1[2] = scale;
412	0	igraphdlarfg_(&c__3, u1, &u1[1], &c__1, &tau1);
413	0	u1[0] = 1.;
414
415	0	temp = -tau1 * (x[2] + u1[1] * x[3]);
416	0	u2[0] = -temp * u1[1] - x[3];
417	0	u2[1] = -temp * u1[2];
418	0	u2[2] = scale;
419	0	igraphdlarfg_(&c__3, u2, &u2[1], &c__1, &tau2);
420	0	u2[0] = 1.;
421
422		/* Perform swap provisionally on diagonal block in D. */
423
424	0	igraphdlarfx_("L", &c__3, &c__4, u1, &tau1, d__, &c__4, &work[1])
425	0	;
426	0	igraphdlarfx_("R", &c__4, &c__3, u1, &tau1, d__, &c__4, &work[1])
427	0	;
428	0	igraphdlarfx_("L", &c__3, &c__4, u2, &tau2, &d__[1], &c__4, &work[1]);
429	0	igraphdlarfx_("R", &c__4, &c__3, u2, &tau2, &d__[4], &c__4, &work[1]);
430
431		/* Test whether to reject swap.
432
433		Computing MAX */
434	0	d__1 = abs(d__[2]), d__2 = abs(d__[6]), d__1 = max(d__1,d__2), d__2 =
435	0	abs(d__[3]), d__1 = max(d__1,d__2), d__2 = abs(d__[7]);
436	0	if (max(d__1,d__2) > thresh) {
437	0	goto L50;
438	0	}
439
440		/* Accept swap: apply transformation to the entire matrix T. */
441
442	0	i__1 = n - j1 + 1;
443	0	igraphdlarfx_("L", &c__3, &i__1, u1, &tau1, &t[j1 + j1 * t_dim1], ldt, &
444	0	work[1]);
445	0	igraphdlarfx_("R", &j4, &c__3, u1, &tau1, &t[j1 t_dim1 + 1], ldt, &work[
446	0	1]);
447	0	i__1 = n - j1 + 1;
448	0	igraphdlarfx_("L", &c__3, &i__1, u2, &tau2, &t[j2 + j1 t_dim1], ldt, &
449	0	work[1]);
450	0	igraphdlarfx_("R", &j4, &c__3, u2, &tau2, &t[j2 * t_dim1 + 1], ldt, &work[1]
451	0	);
452
453	0	t[j3 + j1 t_dim1] = 0.;
454	0	t[j3 + j2 * t_dim1] = 0.;
455	0	t[j4 + j1 t_dim1] = 0.;
456	0	t[j4 + j2 * t_dim1] = 0.;
457
458	0	if (*wantq) {
459
460		/* Accumulate transformation in the matrix Q. */
461
462	0	igraphdlarfx_("R", n, &c__3, u1, &tau1, &q[j1 q_dim1 + 1], ldq, &
463	0	work[1]);
464	0	igraphdlarfx_("R", n, &c__3, u2, &tau2, &q[j2 * q_dim1 + 1], ldq, &work[
465	0	1]);
466	0	}
467
468	0	L40:
469
470	0	if (*n2 == 2) {
471
472		/* Standardize new 2-by-2 block T11 */
473
474	0	igraphdlanv2_(&t[j1 + j1 * t_dim1], &t[j1 + j2 t_dim1], &t[j2 + *
475	0	j1 * t_dim1], &t[j2 + j2 * t_dim1], &wr1, &wi1, &wr2, &
476	0	wi2, &cs, &sn);
477	0	i__1 = n - j1 - 1;
478	0	igraphdrot_(&i__1, &t[j1 + (j1 + 2) * t_dim1], ldt, &t[j2 + (*j1 + 2)
479	0	* t_dim1], ldt, &cs, &sn);
480	0	i__1 = *j1 - 1;
481	0	igraphdrot_(&i__1, &t[j1 t_dim1 + 1], &c__1, &t[j2 * t_dim1 + 1], &
482	0	c__1, &cs, &sn);
483	0	if (*wantq) {
484	0	igraphdrot_(n, &q[j1 q_dim1 + 1], &c__1, &q[j2 * q_dim1 + 1], &
485	0	c__1, &cs, &sn);
486	0	}
487	0	}
488
489	0	if (*n1 == 2) {
490
491		/* Standardize new 2-by-2 block T22 */
492
493	0	j3 = j1 + n2;
494	0	j4 = j3 + 1;
495	0	igraphdlanv2_(&t[j3 + j3 * t_dim1], &t[j3 + j4 * t_dim1], &t[j4 + j3 *
496	0	t_dim1], &t[j4 + j4 * t_dim1], &wr1, &wi1, &wr2, &wi2, &
497	0	cs, &sn);
498	0	if (j3 + 2 <= *n) {
499	0	i__1 = *n - j3 - 1;
500	0	igraphdrot_(&i__1, &t[j3 + (j3 + 2) * t_dim1], ldt, &t[j4 + (j3 + 2)
501	0	* t_dim1], ldt, &cs, &sn);
502	0	}
503	0	i__1 = j3 - 1;
504	0	igraphdrot_(&i__1, &t[j3 * t_dim1 + 1], &c__1, &t[j4 * t_dim1 + 1], &
505	0	c__1, &cs, &sn);
506	0	if (*wantq) {
507	0	igraphdrot_(n, &q[j3 * q_dim1 + 1], &c__1, &q[j4 * q_dim1 + 1], &
508	0	c__1, &cs, &sn);
509	0	}
510	0	}
511
512	0	}
513	0	return 0;
514
515		/* Exit with INFO = 1 if swap was rejected. */
516
517	0	L50:
518	0	*info = 1;
519	0	return 0;
520
521		/* End of DLAEXC */
522
523	0	} /* igraphdlaexc_ */
524