#include "f2c.h" #include "blaswrap.h" /* Table of constant values */ static integer c__4 = 4; static integer c__8 = 8; /* Subroutine */ int zlarot_(logical *lrows, logical *lleft, logical *lright, integer *nl, doublecomplex *c__, doublecomplex *s, doublecomplex *a, integer *lda, doublecomplex *xleft, doublecomplex *xright) { /* System generated locals */ integer i__1, i__2, i__3, i__4; doublecomplex z__1, z__2, z__3, z__4, z__5, z__6; /* Builtin functions */ void d_cnjg(doublecomplex *, doublecomplex *); /* Local variables */ integer j, ix, iy, nt; doublecomplex xt[2], yt[2]; integer iyt, iinc, inext; doublecomplex tempx; extern /* Subroutine */ int xerbla_(char *, integer *); /* -- LAPACK auxiliary test routine (version 3.1) -- */ /* Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. */ /* November 2006 */ /* .. Scalar Arguments .. */ /* .. */ /* .. Array Arguments .. */ /* .. */ /* Purpose */ /* ======= */ /* ZLAROT applies a (Givens) rotation to two adjacent rows or */ /* columns, where one element of the first and/or last column/row */ /* November 2006 */ /* for use on matrices stored in some format other than GE, so */ /* that elements of the matrix may be used or modified for which */ /* no array element is provided. */ /* One example is a symmetric matrix in SB format (bandwidth=4), for */ /* which UPLO='L': Two adjacent rows will have the format: */ /* row j: * * * * * . . . . */ /* row j+1: * * * * * . . . . */ /* '*' indicates elements for which storage is provided, */ /* '.' indicates elements for which no storage is provided, but */ /* are not necessarily zero; their values are determined by */ /* symmetry. ' ' indicates elements which are necessarily zero, */ /* and have no storage provided. */ /* Those columns which have two '*'s can be handled by DROT. */ /* Those columns which have no '*'s can be ignored, since as long */ /* as the Givens rotations are carefully applied to preserve */ /* symmetry, their values are determined. */ /* Those columns which have one '*' have to be handled separately, */ /* by using separate variables "p" and "q": */ /* row j: * * * * * p . . . */ /* row j+1: q * * * * * . . . . */ /* The element p would have to be set correctly, then that column */ /* is rotated, setting p to its new value. The next call to */ /* ZLAROT would rotate columns j and j+1, using p, and restore */ /* symmetry. The element q would start out being zero, and be */ /* made non-zero by the rotation. Later, rotations would presumably */ /* be chosen to zero q out. */ /* Typical Calling Sequences: rotating the i-th and (i+1)-st rows. */ /* ------- ------- --------- */ /* General dense matrix: */ /* CALL ZLAROT(.TRUE.,.FALSE.,.FALSE., N, C,S, */ /* A(i,1),LDA, DUMMY, DUMMY) */ /* General banded matrix in GB format: */ /* j = MAX(1, i-KL ) */ /* NL = MIN( N, i+KU+1 ) + 1-j */ /* CALL ZLAROT( .TRUE., i-KL.GE.1, i+KU.LT.N, NL, C,S, */ /* A(KU+i+1-j,j),LDA-1, XLEFT, XRIGHT ) */ /* [ note that i+1-j is just MIN(i,KL+1) ] */ /* Symmetric banded matrix in SY format, bandwidth K, */ /* lower triangle only: */ /* j = MAX(1, i-K ) */ /* NL = MIN( K+1, i ) + 1 */ /* CALL ZLAROT( .TRUE., i-K.GE.1, .TRUE., NL, C,S, */ /* A(i,j), LDA, XLEFT, XRIGHT ) */ /* Same, but upper triangle only: */ /* NL = MIN( K+1, N-i ) + 1 */ /* CALL ZLAROT( .TRUE., .TRUE., i+K.LT.N, NL, C,S, */ /* A(i,i), LDA, XLEFT, XRIGHT ) */ /* Symmetric banded matrix in SB format, bandwidth K, */ /* lower triangle only: */ /* [ same as for SY, except:] */ /* . . . . */ /* A(i+1-j,j), LDA-1, XLEFT, XRIGHT ) */ /* [ note that i+1-j is just MIN(i,K+1) ] */ /* Same, but upper triangle only: */ /* . . . */ /* A(K+1,i), LDA-1, XLEFT, XRIGHT ) */ /* Rotating columns is just the transpose of rotating rows, except */ /* for GB and SB: (rotating columns i and i+1) */ /* GB: */ /* j = MAX(1, i-KU ) */ /* NL = MIN( N, i+KL+1 ) + 1-j */ /* CALL ZLAROT( .TRUE., i-KU.GE.1, i+KL.LT.N, NL, C,S, */ /* A(KU+j+1-i,i),LDA-1, XTOP, XBOTTM ) */ /* [note that KU+j+1-i is just MAX(1,KU+2-i)] */ /* SB: (upper triangle) */ /* . . . . . . */ /* A(K+j+1-i,i),LDA-1, XTOP, XBOTTM ) */ /* SB: (lower triangle) */ /* . . . . . . */ /* A(1,i),LDA-1, XTOP, XBOTTM ) */ /* Arguments */ /* ========= */ /* LROWS - LOGICAL */ /* If .TRUE., then ZLAROT will rotate two rows. If .FALSE., */ /* then it will rotate two columns. */ /* Not modified. */ /* LLEFT - LOGICAL */ /* If .TRUE., then XLEFT will be used instead of the */ /* corresponding element of A for the first element in the */ /* second row (if LROWS=.FALSE.) or column (if LROWS=.TRUE.) */ /* If .FALSE., then the corresponding element of A will be */ /* used. */ /* Not modified. */ /* LRIGHT - LOGICAL */ /* If .TRUE., then XRIGHT will be used instead of the */ /* corresponding element of A for the last element in the */ /* first row (if LROWS=.FALSE.) or column (if LROWS=.TRUE.) If */ /* .FALSE., then the corresponding element of A will be used. */ /* Not modified. */ /* NL - INTEGER */ /* The length of the rows (if LROWS=.TRUE.) or columns (if */ /* LROWS=.FALSE.) to be rotated. If XLEFT and/or XRIGHT are */ /* used, the columns/rows they are in should be included in */ /* NL, e.g., if LLEFT = LRIGHT = .TRUE., then NL must be at */ /* least 2. The number of rows/columns to be rotated */ /* exclusive of those involving XLEFT and/or XRIGHT may */ /* not be negative, i.e., NL minus how many of LLEFT and */ /* LRIGHT are .TRUE. must be at least zero; if not, XERBLA */ /* will be called. */ /* Not modified. */ /* C, S - COMPLEX*16 */ /* Specify the Givens rotation to be applied. If LROWS is */ /* true, then the matrix ( c s ) */ /* ( _ _ ) */ /* (-s c ) is applied from the left; */ /* if false, then the transpose (not conjugated) thereof is */ /* applied from the right. Note that in contrast to the */ /* output of ZROTG or to most versions of ZROT, both C and S */ /* are complex. For a Givens rotation, |C|**2 + |S|**2 should */ /* be 1, but this is not checked. */ /* Not modified. */ /* A - COMPLEX*16 array. */ /* The array containing the rows/columns to be rotated. The */ /* first element of A should be the upper left element to */ /* be rotated. */ /* Read and modified. */ /* LDA - INTEGER */ /* The "effective" leading dimension of A. If A contains */ /* a matrix stored in GE, HE, or SY format, then this is just */ /* the leading dimension of A as dimensioned in the calling */ /* routine. If A contains a matrix stored in band (GB, HB, or */ /* SB) format, then this should be *one less* than the leading */ /* dimension used in the calling routine. Thus, if A were */ /* dimensioned A(LDA,*) in ZLAROT, then A(1,j) would be the */ /* j-th element in the first of the two rows to be rotated, */ /* and A(2,j) would be the j-th in the second, regardless of */ /* how the array may be stored in the calling routine. [A */ /* cannot, however, actually be dimensioned thus, since for */ /* band format, the row number may exceed LDA, which is not */ /* legal FORTRAN.] */ /* If LROWS=.TRUE., then LDA must be at least 1, otherwise */ /* it must be at least NL minus the number of .TRUE. values */ /* in XLEFT and XRIGHT. */ /* Not modified. */ /* XLEFT - COMPLEX*16 */ /* If LLEFT is .TRUE., then XLEFT will be used and modified */ /* instead of A(2,1) (if LROWS=.TRUE.) or A(1,2) */ /* (if LROWS=.FALSE.). */ /* Read and modified. */ /* XRIGHT - COMPLEX*16 */ /* If LRIGHT is .TRUE., then XRIGHT will be used and modified */ /* instead of A(1,NL) (if LROWS=.TRUE.) or A(NL,1) */ /* (if LROWS=.FALSE.). */ /* Read and modified. */ /* ===================================================================== */ /* .. Local Scalars .. */ /* .. */ /* .. Local Arrays .. */ /* .. */ /* .. External Subroutines .. */ /* .. */ /* .. Intrinsic Functions .. */ /* .. */ /* .. Executable Statements .. */ /* Set up indices, arrays for ends */ /* Parameter adjustments */ --a; /* Function Body */ if (*lrows) { iinc = *lda; inext = 1; } else { iinc = 1; inext = *lda; } if (*lleft) { nt = 1; ix = iinc + 1; iy = *lda + 2; xt[0].r = a[1].r, xt[0].i = a[1].i; yt[0].r = xleft->r, yt[0].i = xleft->i; } else { nt = 0; ix = 1; iy = inext + 1; } if (*lright) { iyt = inext + 1 + (*nl - 1) * iinc; ++nt; i__1 = nt - 1; xt[i__1].r = xright->r, xt[i__1].i = xright->i; i__1 = nt - 1; i__2 = iyt; yt[i__1].r = a[i__2].r, yt[i__1].i = a[i__2].i; } /* Check for errors */ if (*nl < nt) { xerbla_("ZLAROT", &c__4); return 0; } if (*lda <= 0 || ! (*lrows) && *lda < *nl - nt) { xerbla_("ZLAROT", &c__8); return 0; } /* Rotate */ /* ZROT( NL-NT, A(IX),IINC, A(IY),IINC, C, S ) with complex C, S */ i__1 = *nl - nt - 1; for (j = 0; j <= i__1; ++j) { i__2 = ix + j * iinc; z__2.r = c__->r * a[i__2].r - c__->i * a[i__2].i, z__2.i = c__->r * a[ i__2].i + c__->i * a[i__2].r; i__3 = iy + j * iinc; z__3.r = s->r * a[i__3].r - s->i * a[i__3].i, z__3.i = s->r * a[i__3] .i + s->i * a[i__3].r; z__1.r = z__2.r + z__3.r, z__1.i = z__2.i + z__3.i; tempx.r = z__1.r, tempx.i = z__1.i; i__2 = iy + j * iinc; d_cnjg(&z__4, s); z__3.r = -z__4.r, z__3.i = -z__4.i; i__3 = ix + j * iinc; z__2.r = z__3.r * a[i__3].r - z__3.i * a[i__3].i, z__2.i = z__3.r * a[ i__3].i + z__3.i * a[i__3].r; d_cnjg(&z__6, c__); i__4 = iy + j * iinc; z__5.r = z__6.r * a[i__4].r - z__6.i * a[i__4].i, z__5.i = z__6.r * a[ i__4].i + z__6.i * a[i__4].r; z__1.r = z__2.r + z__5.r, z__1.i = z__2.i + z__5.i; a[i__2].r = z__1.r, a[i__2].i = z__1.i; i__2 = ix + j * iinc; a[i__2].r = tempx.r, a[i__2].i = tempx.i; /* L10: */ } /* ZROT( NT, XT,1, YT,1, C, S ) with complex C, S */ i__1 = nt; for (j = 1; j <= i__1; ++j) { i__2 = j - 1; z__2.r = c__->r * xt[i__2].r - c__->i * xt[i__2].i, z__2.i = c__->r * xt[i__2].i + c__->i * xt[i__2].r; i__3 = j - 1; z__3.r = s->r * yt[i__3].r - s->i * yt[i__3].i, z__3.i = s->r * yt[ i__3].i + s->i * yt[i__3].r; z__1.r = z__2.r + z__3.r, z__1.i = z__2.i + z__3.i; tempx.r = z__1.r, tempx.i = z__1.i; i__2 = j - 1; d_cnjg(&z__4, s); z__3.r = -z__4.r, z__3.i = -z__4.i; i__3 = j - 1; z__2.r = z__3.r * xt[i__3].r - z__3.i * xt[i__3].i, z__2.i = z__3.r * xt[i__3].i + z__3.i * xt[i__3].r; d_cnjg(&z__6, c__); i__4 = j - 1; z__5.r = z__6.r * yt[i__4].r - z__6.i * yt[i__4].i, z__5.i = z__6.r * yt[i__4].i + z__6.i * yt[i__4].r; z__1.r = z__2.r + z__5.r, z__1.i = z__2.i + z__5.i; yt[i__2].r = z__1.r, yt[i__2].i = z__1.i; i__2 = j - 1; xt[i__2].r = tempx.r, xt[i__2].i = tempx.i; /* L20: */ } /* Stuff values back into XLEFT, XRIGHT, etc. */ if (*lleft) { a[1].r = xt[0].r, a[1].i = xt[0].i; xleft->r = yt[0].r, xleft->i = yt[0].i; } if (*lright) { i__1 = nt - 1; xright->r = xt[i__1].r, xright->i = xt[i__1].i; i__1 = iyt; i__2 = nt - 1; a[i__1].r = yt[i__2].r, a[i__1].i = yt[i__2].i; } return 0; /* End of ZLAROT */ } /* zlarot_ */