LAPACK  3.6.0
LAPACK: Linear Algebra PACKage
Collaboration diagram for complex16:

Functions

subroutine zlaesy (A, B, C, RT1, RT2, EVSCAL, CS1, SN1)
 ZLAESY computes the eigenvalues and eigenvectors of a 2-by-2 complex symmetric matrix. More...
 
double precision function zlansy (NORM, UPLO, N, A, LDA, WORK)
 ZLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex symmetric matrix. More...
 
subroutine zlaqsy (UPLO, N, A, LDA, S, SCOND, AMAX, EQUED)
 ZLAQSY scales a symmetric/Hermitian matrix, using scaling factors computed by spoequ. More...
 
subroutine zsymv (UPLO, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
 ZSYMV computes a matrix-vector product for a complex symmetric matrix. More...
 
subroutine zsyr (UPLO, N, ALPHA, X, INCX, A, LDA)
 ZSYR performs the symmetric rank-1 update of a complex symmetric matrix. More...
 
subroutine zsyswapr (UPLO, N, A, LDA, I1, I2)
 ZSYSWAPR More...
 
subroutine ztgsy2 (TRANS, IJOB, M, N, A, LDA, B, LDB, C, LDC, D, LDD, E, LDE, F, LDF, SCALE, RDSUM, RDSCAL, INFO)
 ZTGSY2 solves the generalized Sylvester equation (unblocked algorithm). More...
 

Detailed Description

This is the group of complex16 auxiliary functions for SY matrices

Function Documentation

subroutine zlaesy ( complex*16  A,
complex*16  B,
complex*16  C,
complex*16  RT1,
complex*16  RT2,
complex*16  EVSCAL,
complex*16  CS1,
complex*16  SN1 
)

ZLAESY computes the eigenvalues and eigenvectors of a 2-by-2 complex symmetric matrix.

Download ZLAESY + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZLAESY computes the eigendecomposition of a 2-by-2 symmetric matrix
    ( ( A, B );( B, C ) )
 provided the norm of the matrix of eigenvectors is larger than
 some threshold value.

 RT1 is the eigenvalue of larger absolute value, and RT2 of
 smaller absolute value.  If the eigenvectors are computed, then
 on return ( CS1, SN1 ) is the unit eigenvector for RT1, hence

 [  CS1     SN1   ] . [ A  B ] . [ CS1    -SN1   ] = [ RT1  0  ]
 [ -SN1     CS1   ]   [ B  C ]   [ SN1     CS1   ]   [  0  RT2 ]
Parameters
[in]A
          A is COMPLEX*16
          The ( 1, 1 ) element of input matrix.
[in]B
          B is COMPLEX*16
          The ( 1, 2 ) element of input matrix.  The ( 2, 1 ) element
          is also given by B, since the 2-by-2 matrix is symmetric.
[in]C
          C is COMPLEX*16
          The ( 2, 2 ) element of input matrix.
[out]RT1
          RT1 is COMPLEX*16
          The eigenvalue of larger modulus.
[out]RT2
          RT2 is COMPLEX*16
          The eigenvalue of smaller modulus.
[out]EVSCAL
          EVSCAL is COMPLEX*16
          The complex value by which the eigenvector matrix was scaled
          to make it orthonormal.  If EVSCAL is zero, the eigenvectors
          were not computed.  This means one of two things:  the 2-by-2
          matrix could not be diagonalized, or the norm of the matrix
          of eigenvectors before scaling was larger than the threshold
          value THRESH (set below).
[out]CS1
          CS1 is COMPLEX*16
[out]SN1
          SN1 is COMPLEX*16
          If EVSCAL .NE. 0,  ( CS1, SN1 ) is the unit right eigenvector
          for RT1.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 117 of file zlaesy.f.

117 *
118 * -- LAPACK auxiliary routine (version 3.4.2) --
119 * -- LAPACK is a software package provided by Univ. of Tennessee, --
120 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
121 * September 2012
122 *
123 * .. Scalar Arguments ..
124  COMPLEX*16 a, b, c, cs1, evscal, rt1, rt2, sn1
125 * ..
126 *
127 * =====================================================================
128 *
129 * .. Parameters ..
130  DOUBLE PRECISION zero
131  parameter( zero = 0.0d0 )
132  DOUBLE PRECISION one
133  parameter( one = 1.0d0 )
134  COMPLEX*16 cone
135  parameter( cone = ( 1.0d0, 0.0d0 ) )
136  DOUBLE PRECISION half
137  parameter( half = 0.5d0 )
138  DOUBLE PRECISION thresh
139  parameter( thresh = 0.1d0 )
140 * ..
141 * .. Local Scalars ..
142  DOUBLE PRECISION babs, evnorm, tabs, z
143  COMPLEX*16 s, t, tmp
144 * ..
145 * .. Intrinsic Functions ..
146  INTRINSIC abs, max, sqrt
147 * ..
148 * .. Executable Statements ..
149 *
150 *
151 * Special case: The matrix is actually diagonal.
152 * To avoid divide by zero later, we treat this case separately.
153 *
154  IF( abs( b ).EQ.zero ) THEN
155  rt1 = a
156  rt2 = c
157  IF( abs( rt1 ).LT.abs( rt2 ) ) THEN
158  tmp = rt1
159  rt1 = rt2
160  rt2 = tmp
161  cs1 = zero
162  sn1 = one
163  ELSE
164  cs1 = one
165  sn1 = zero
166  END IF
167  ELSE
168 *
169 * Compute the eigenvalues and eigenvectors.
170 * The characteristic equation is
171 * lambda **2 - (A+C) lambda + (A*C - B*B)
172 * and we solve it using the quadratic formula.
173 *
174  s = ( a+c )*half
175  t = ( a-c )*half
176 *
177 * Take the square root carefully to avoid over/under flow.
178 *
179  babs = abs( b )
180  tabs = abs( t )
181  z = max( babs, tabs )
182  IF( z.GT.zero )
183  $ t = z*sqrt( ( t / z )**2+( b / z )**2 )
184 *
185 * Compute the two eigenvalues. RT1 and RT2 are exchanged
186 * if necessary so that RT1 will have the greater magnitude.
187 *
188  rt1 = s + t
189  rt2 = s - t
190  IF( abs( rt1 ).LT.abs( rt2 ) ) THEN
191  tmp = rt1
192  rt1 = rt2
193  rt2 = tmp
194  END IF
195 *
196 * Choose CS1 = 1 and SN1 to satisfy the first equation, then
197 * scale the components of this eigenvector so that the matrix
198 * of eigenvectors X satisfies X * X**T = I . (No scaling is
199 * done if the norm of the eigenvalue matrix is less than THRESH.)
200 *
201  sn1 = ( rt1-a ) / b
202  tabs = abs( sn1 )
203  IF( tabs.GT.one ) THEN
204  t = tabs*sqrt( ( one / tabs )**2+( sn1 / tabs )**2 )
205  ELSE
206  t = sqrt( cone+sn1*sn1 )
207  END IF
208  evnorm = abs( t )
209  IF( evnorm.GE.thresh ) THEN
210  evscal = cone / t
211  cs1 = evscal
212  sn1 = sn1*evscal
213  ELSE
214  evscal = zero
215  END IF
216  END IF
217  RETURN
218 *
219 * End of ZLAESY
220 *
double precision function zlansy ( character  NORM,
character  UPLO,
integer  N,
complex*16, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( * )  WORK 
)

ZLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex symmetric matrix.

Download ZLANSY + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZLANSY  returns the value of the one norm,  or the Frobenius norm, or
 the  infinity norm,  or the  element of  largest absolute value  of a
 complex symmetric matrix A.
Returns
ZLANSY
    ZLANSY = ( max(abs(A(i,j))), NORM = 'M' or 'm'
             (
             ( norm1(A),         NORM = '1', 'O' or 'o'
             (
             ( normI(A),         NORM = 'I' or 'i'
             (
             ( normF(A),         NORM = 'F', 'f', 'E' or 'e'

 where  norm1  denotes the  one norm of a matrix (maximum column sum),
 normI  denotes the  infinity norm  of a matrix  (maximum row sum) and
 normF  denotes the  Frobenius norm of a matrix (square root of sum of
 squares).  Note that  max(abs(A(i,j)))  is not a consistent matrix norm.
Parameters
[in]NORM
          NORM is CHARACTER*1
          Specifies the value to be returned in ZLANSY as described
          above.
[in]UPLO
          UPLO is CHARACTER*1
          Specifies whether the upper or lower triangular part of the
          symmetric matrix A is to be referenced.
          = 'U':  Upper triangular part of A is referenced
          = 'L':  Lower triangular part of A is referenced
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.  When N = 0, ZLANSY is
          set to zero.
[in]A
          A is COMPLEX*16 array, dimension (LDA,N)
          The symmetric matrix A.  If UPLO = 'U', the leading n by n
          upper triangular part of A contains the upper triangular part
          of the matrix A, and the strictly lower triangular part of A
          is not referenced.  If UPLO = 'L', the leading n by n lower
          triangular part of A contains the lower triangular part of
          the matrix A, and the strictly upper triangular part of A is
          not referenced.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(N,1).
[out]WORK
          WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK)),
          where LWORK >= N when NORM = 'I' or '1' or 'O'; otherwise,
          WORK is not referenced.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015

Definition at line 125 of file zlansy.f.

125 *
126 * -- LAPACK auxiliary routine (version 3.6.0) --
127 * -- LAPACK is a software package provided by Univ. of Tennessee, --
128 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
129 * November 2015
130 *
131 * .. Scalar Arguments ..
132  CHARACTER norm, uplo
133  INTEGER lda, n
134 * ..
135 * .. Array Arguments ..
136  DOUBLE PRECISION work( * )
137  COMPLEX*16 a( lda, * )
138 * ..
139 *
140 * =====================================================================
141 *
142 * .. Parameters ..
143  DOUBLE PRECISION one, zero
144  parameter( one = 1.0d+0, zero = 0.0d+0 )
145 * ..
146 * .. Local Scalars ..
147  INTEGER i, j
148  DOUBLE PRECISION absa, scale, sum, value
149 * ..
150 * .. External Functions ..
151  LOGICAL lsame, disnan
152  EXTERNAL lsame, disnan
153 * ..
154 * .. External Subroutines ..
155  EXTERNAL zlassq
156 * ..
157 * .. Intrinsic Functions ..
158  INTRINSIC abs, sqrt
159 * ..
160 * .. Executable Statements ..
161 *
162  IF( n.EQ.0 ) THEN
163  VALUE = zero
164  ELSE IF( lsame( norm, 'M' ) ) THEN
165 *
166 * Find max(abs(A(i,j))).
167 *
168  VALUE = zero
169  IF( lsame( uplo, 'U' ) ) THEN
170  DO 20 j = 1, n
171  DO 10 i = 1, j
172  sum = abs( a( i, j ) )
173  IF( VALUE .LT. sum .OR. disnan( sum ) ) VALUE = sum
174  10 CONTINUE
175  20 CONTINUE
176  ELSE
177  DO 40 j = 1, n
178  DO 30 i = j, n
179  sum = abs( a( i, j ) )
180  IF( VALUE .LT. sum .OR. disnan( sum ) ) VALUE = sum
181  30 CONTINUE
182  40 CONTINUE
183  END IF
184  ELSE IF( ( lsame( norm, 'I' ) ) .OR. ( lsame( norm, 'O' ) ) .OR.
185  $ ( norm.EQ.'1' ) ) THEN
186 *
187 * Find normI(A) ( = norm1(A), since A is symmetric).
188 *
189  VALUE = zero
190  IF( lsame( uplo, 'U' ) ) THEN
191  DO 60 j = 1, n
192  sum = zero
193  DO 50 i = 1, j - 1
194  absa = abs( a( i, j ) )
195  sum = sum + absa
196  work( i ) = work( i ) + absa
197  50 CONTINUE
198  work( j ) = sum + abs( a( j, j ) )
199  60 CONTINUE
200  DO 70 i = 1, n
201  sum = work( i )
202  IF( VALUE .LT. sum .OR. disnan( sum ) ) VALUE = sum
203  70 CONTINUE
204  ELSE
205  DO 80 i = 1, n
206  work( i ) = zero
207  80 CONTINUE
208  DO 100 j = 1, n
209  sum = work( j ) + abs( a( j, j ) )
210  DO 90 i = j + 1, n
211  absa = abs( a( i, j ) )
212  sum = sum + absa
213  work( i ) = work( i ) + absa
214  90 CONTINUE
215  IF( VALUE .LT. sum .OR. disnan( sum ) ) VALUE = sum
216  100 CONTINUE
217  END IF
218  ELSE IF( ( lsame( norm, 'F' ) ) .OR. ( lsame( norm, 'E' ) ) ) THEN
219 *
220 * Find normF(A).
221 *
222  scale = zero
223  sum = one
224  IF( lsame( uplo, 'U' ) ) THEN
225  DO 110 j = 2, n
226  CALL zlassq( j-1, a( 1, j ), 1, scale, sum )
227  110 CONTINUE
228  ELSE
229  DO 120 j = 1, n - 1
230  CALL zlassq( n-j, a( j+1, j ), 1, scale, sum )
231  120 CONTINUE
232  END IF
233  sum = 2*sum
234  CALL zlassq( n, a, lda+1, scale, sum )
235  VALUE = scale*sqrt( sum )
236  END IF
237 *
238  zlansy = VALUE
239  RETURN
240 *
241 * End of ZLANSY
242 *
double precision function zlansy(NORM, UPLO, N, A, LDA, WORK)
ZLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm, or the element of largest absolute value of a complex symmetric matrix.
Definition: zlansy.f:125
logical function disnan(DIN)
DISNAN tests input for NaN.
Definition: disnan.f:61
subroutine zlassq(N, X, INCX, SCALE, SUMSQ)
ZLASSQ updates a sum of squares represented in scaled form.
Definition: zlassq.f:108
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine zlaqsy ( character  UPLO,
integer  N,
complex*16, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( * )  S,
double precision  SCOND,
double precision  AMAX,
character  EQUED 
)

ZLAQSY scales a symmetric/Hermitian matrix, using scaling factors computed by spoequ.

Download ZLAQSY + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZLAQSY equilibrates a symmetric matrix A using the scaling factors
 in the vector S.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          Specifies whether the upper or lower triangular part of the
          symmetric matrix A is stored.
          = 'U':  Upper triangular
          = 'L':  Lower triangular
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is COMPLEX*16 array, dimension (LDA,N)
          On entry, the symmetric matrix A.  If UPLO = 'U', the leading
          n by n upper triangular part of A contains the upper
          triangular part of the matrix A, and the strictly lower
          triangular part of A is not referenced.  If UPLO = 'L', the
          leading n by n lower triangular part of A contains the lower
          triangular part of the matrix A, and the strictly upper
          triangular part of A is not referenced.

          On exit, if EQUED = 'Y', the equilibrated matrix:
          diag(S) * A * diag(S).
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(N,1).
[in]S
          S is DOUBLE PRECISION array, dimension (N)
          The scale factors for A.
[in]SCOND
          SCOND is DOUBLE PRECISION
          Ratio of the smallest S(i) to the largest S(i).
[in]AMAX
          AMAX is DOUBLE PRECISION
          Absolute value of largest matrix entry.
[out]EQUED
          EQUED is CHARACTER*1
          Specifies whether or not equilibration was done.
          = 'N':  No equilibration.
          = 'Y':  Equilibration was done, i.e., A has been replaced by
                  diag(S) * A * diag(S).
Internal Parameters:
  THRESH is a threshold value used to decide if scaling should be done
  based on the ratio of the scaling factors.  If SCOND < THRESH,
  scaling is done.

  LARGE and SMALL are threshold values used to decide if scaling should
  be done based on the absolute size of the largest matrix element.
  If AMAX > LARGE or AMAX < SMALL, scaling is done.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 136 of file zlaqsy.f.

136 *
137 * -- LAPACK auxiliary routine (version 3.4.2) --
138 * -- LAPACK is a software package provided by Univ. of Tennessee, --
139 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
140 * September 2012
141 *
142 * .. Scalar Arguments ..
143  CHARACTER equed, uplo
144  INTEGER lda, n
145  DOUBLE PRECISION amax, scond
146 * ..
147 * .. Array Arguments ..
148  DOUBLE PRECISION s( * )
149  COMPLEX*16 a( lda, * )
150 * ..
151 *
152 * =====================================================================
153 *
154 * .. Parameters ..
155  DOUBLE PRECISION one, thresh
156  parameter( one = 1.0d+0, thresh = 0.1d+0 )
157 * ..
158 * .. Local Scalars ..
159  INTEGER i, j
160  DOUBLE PRECISION cj, large, small
161 * ..
162 * .. External Functions ..
163  LOGICAL lsame
164  DOUBLE PRECISION dlamch
165  EXTERNAL lsame, dlamch
166 * ..
167 * .. Executable Statements ..
168 *
169 * Quick return if possible
170 *
171  IF( n.LE.0 ) THEN
172  equed = 'N'
173  RETURN
174  END IF
175 *
176 * Initialize LARGE and SMALL.
177 *
178  small = dlamch( 'Safe minimum' ) / dlamch( 'Precision' )
179  large = one / small
180 *
181  IF( scond.GE.thresh .AND. amax.GE.small .AND. amax.LE.large ) THEN
182 *
183 * No equilibration
184 *
185  equed = 'N'
186  ELSE
187 *
188 * Replace A by diag(S) * A * diag(S).
189 *
190  IF( lsame( uplo, 'U' ) ) THEN
191 *
192 * Upper triangle of A is stored.
193 *
194  DO 20 j = 1, n
195  cj = s( j )
196  DO 10 i = 1, j
197  a( i, j ) = cj*s( i )*a( i, j )
198  10 CONTINUE
199  20 CONTINUE
200  ELSE
201 *
202 * Lower triangle of A is stored.
203 *
204  DO 40 j = 1, n
205  cj = s( j )
206  DO 30 i = j, n
207  a( i, j ) = cj*s( i )*a( i, j )
208  30 CONTINUE
209  40 CONTINUE
210  END IF
211  equed = 'Y'
212  END IF
213 *
214  RETURN
215 *
216 * End of ZLAQSY
217 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the caller graph for this function:

subroutine zsymv ( character  UPLO,
integer  N,
complex*16  ALPHA,
complex*16, dimension( lda, * )  A,
integer  LDA,
complex*16, dimension( * )  X,
integer  INCX,
complex*16  BETA,
complex*16, dimension( * )  Y,
integer  INCY 
)

ZSYMV computes a matrix-vector product for a complex symmetric matrix.

Download ZSYMV + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZSYMV  performs the matrix-vector  operation

    y := alpha*A*x + beta*y,

 where alpha and beta are scalars, x and y are n element vectors and
 A is an n by n symmetric matrix.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
           On entry, UPLO specifies whether the upper or lower
           triangular part of the array A is to be referenced as
           follows:

              UPLO = 'U' or 'u'   Only the upper triangular part of A
                                  is to be referenced.

              UPLO = 'L' or 'l'   Only the lower triangular part of A
                                  is to be referenced.

           Unchanged on exit.
[in]N
          N is INTEGER
           On entry, N specifies the order of the matrix A.
           N must be at least zero.
           Unchanged on exit.
[in]ALPHA
          ALPHA is COMPLEX*16
           On entry, ALPHA specifies the scalar alpha.
           Unchanged on exit.
[in]A
          A is COMPLEX*16 array, dimension ( LDA, N )
           Before entry, with  UPLO = 'U' or 'u', the leading n by n
           upper triangular part of the array A must contain the upper
           triangular part of the symmetric matrix and the strictly
           lower triangular part of A is not referenced.
           Before entry, with UPLO = 'L' or 'l', the leading n by n
           lower triangular part of the array A must contain the lower
           triangular part of the symmetric matrix and the strictly
           upper triangular part of A is not referenced.
           Unchanged on exit.
[in]LDA
          LDA is INTEGER
           On entry, LDA specifies the first dimension of A as declared
           in the calling (sub) program. LDA must be at least
           max( 1, N ).
           Unchanged on exit.
[in]X
          X is COMPLEX*16 array, dimension at least
           ( 1 + ( N - 1 )*abs( INCX ) ).
           Before entry, the incremented array X must contain the N-
           element vector x.
           Unchanged on exit.
[in]INCX
          INCX is INTEGER
           On entry, INCX specifies the increment for the elements of
           X. INCX must not be zero.
           Unchanged on exit.
[in]BETA
          BETA is COMPLEX*16
           On entry, BETA specifies the scalar beta. When BETA is
           supplied as zero then Y need not be set on input.
           Unchanged on exit.
[in,out]Y
          Y is COMPLEX*16 array, dimension at least
           ( 1 + ( N - 1 )*abs( INCY ) ).
           Before entry, the incremented array Y must contain the n
           element vector y. On exit, Y is overwritten by the updated
           vector y.
[in]INCY
          INCY is INTEGER
           On entry, INCY specifies the increment for the elements of
           Y. INCY must not be zero.
           Unchanged on exit.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 159 of file zsymv.f.

159 *
160 * -- LAPACK auxiliary routine (version 3.4.2) --
161 * -- LAPACK is a software package provided by Univ. of Tennessee, --
162 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
163 * September 2012
164 *
165 * .. Scalar Arguments ..
166  CHARACTER uplo
167  INTEGER incx, incy, lda, n
168  COMPLEX*16 alpha, beta
169 * ..
170 * .. Array Arguments ..
171  COMPLEX*16 a( lda, * ), x( * ), y( * )
172 * ..
173 *
174 * =====================================================================
175 *
176 * .. Parameters ..
177  COMPLEX*16 one
178  parameter( one = ( 1.0d+0, 0.0d+0 ) )
179  COMPLEX*16 zero
180  parameter( zero = ( 0.0d+0, 0.0d+0 ) )
181 * ..
182 * .. Local Scalars ..
183  INTEGER i, info, ix, iy, j, jx, jy, kx, ky
184  COMPLEX*16 temp1, temp2
185 * ..
186 * .. External Functions ..
187  LOGICAL lsame
188  EXTERNAL lsame
189 * ..
190 * .. External Subroutines ..
191  EXTERNAL xerbla
192 * ..
193 * .. Intrinsic Functions ..
194  INTRINSIC max
195 * ..
196 * .. Executable Statements ..
197 *
198 * Test the input parameters.
199 *
200  info = 0
201  IF( .NOT.lsame( uplo, 'U' ) .AND. .NOT.lsame( uplo, 'L' ) ) THEN
202  info = 1
203  ELSE IF( n.LT.0 ) THEN
204  info = 2
205  ELSE IF( lda.LT.max( 1, n ) ) THEN
206  info = 5
207  ELSE IF( incx.EQ.0 ) THEN
208  info = 7
209  ELSE IF( incy.EQ.0 ) THEN
210  info = 10
211  END IF
212  IF( info.NE.0 ) THEN
213  CALL xerbla( 'ZSYMV ', info )
214  RETURN
215  END IF
216 *
217 * Quick return if possible.
218 *
219  IF( ( n.EQ.0 ) .OR. ( ( alpha.EQ.zero ) .AND. ( beta.EQ.one ) ) )
220  $ RETURN
221 *
222 * Set up the start points in X and Y.
223 *
224  IF( incx.GT.0 ) THEN
225  kx = 1
226  ELSE
227  kx = 1 - ( n-1 )*incx
228  END IF
229  IF( incy.GT.0 ) THEN
230  ky = 1
231  ELSE
232  ky = 1 - ( n-1 )*incy
233  END IF
234 *
235 * Start the operations. In this version the elements of A are
236 * accessed sequentially with one pass through the triangular part
237 * of A.
238 *
239 * First form y := beta*y.
240 *
241  IF( beta.NE.one ) THEN
242  IF( incy.EQ.1 ) THEN
243  IF( beta.EQ.zero ) THEN
244  DO 10 i = 1, n
245  y( i ) = zero
246  10 CONTINUE
247  ELSE
248  DO 20 i = 1, n
249  y( i ) = beta*y( i )
250  20 CONTINUE
251  END IF
252  ELSE
253  iy = ky
254  IF( beta.EQ.zero ) THEN
255  DO 30 i = 1, n
256  y( iy ) = zero
257  iy = iy + incy
258  30 CONTINUE
259  ELSE
260  DO 40 i = 1, n
261  y( iy ) = beta*y( iy )
262  iy = iy + incy
263  40 CONTINUE
264  END IF
265  END IF
266  END IF
267  IF( alpha.EQ.zero )
268  $ RETURN
269  IF( lsame( uplo, 'U' ) ) THEN
270 *
271 * Form y when A is stored in upper triangle.
272 *
273  IF( ( incx.EQ.1 ) .AND. ( incy.EQ.1 ) ) THEN
274  DO 60 j = 1, n
275  temp1 = alpha*x( j )
276  temp2 = zero
277  DO 50 i = 1, j - 1
278  y( i ) = y( i ) + temp1*a( i, j )
279  temp2 = temp2 + a( i, j )*x( i )
280  50 CONTINUE
281  y( j ) = y( j ) + temp1*a( j, j ) + alpha*temp2
282  60 CONTINUE
283  ELSE
284  jx = kx
285  jy = ky
286  DO 80 j = 1, n
287  temp1 = alpha*x( jx )
288  temp2 = zero
289  ix = kx
290  iy = ky
291  DO 70 i = 1, j - 1
292  y( iy ) = y( iy ) + temp1*a( i, j )
293  temp2 = temp2 + a( i, j )*x( ix )
294  ix = ix + incx
295  iy = iy + incy
296  70 CONTINUE
297  y( jy ) = y( jy ) + temp1*a( j, j ) + alpha*temp2
298  jx = jx + incx
299  jy = jy + incy
300  80 CONTINUE
301  END IF
302  ELSE
303 *
304 * Form y when A is stored in lower triangle.
305 *
306  IF( ( incx.EQ.1 ) .AND. ( incy.EQ.1 ) ) THEN
307  DO 100 j = 1, n
308  temp1 = alpha*x( j )
309  temp2 = zero
310  y( j ) = y( j ) + temp1*a( j, j )
311  DO 90 i = j + 1, n
312  y( i ) = y( i ) + temp1*a( i, j )
313  temp2 = temp2 + a( i, j )*x( i )
314  90 CONTINUE
315  y( j ) = y( j ) + alpha*temp2
316  100 CONTINUE
317  ELSE
318  jx = kx
319  jy = ky
320  DO 120 j = 1, n
321  temp1 = alpha*x( jx )
322  temp2 = zero
323  y( jy ) = y( jy ) + temp1*a( j, j )
324  ix = jx
325  iy = jy
326  DO 110 i = j + 1, n
327  ix = ix + incx
328  iy = iy + incy
329  y( iy ) = y( iy ) + temp1*a( i, j )
330  temp2 = temp2 + a( i, j )*x( ix )
331  110 CONTINUE
332  y( jy ) = y( jy ) + alpha*temp2
333  jx = jx + incx
334  jy = jy + incy
335  120 CONTINUE
336  END IF
337  END IF
338 *
339  RETURN
340 *
341 * End of ZSYMV
342 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine zsyr ( character  UPLO,
integer  N,
complex*16  ALPHA,
complex*16, dimension( * )  X,
integer  INCX,
complex*16, dimension( lda, * )  A,
integer  LDA 
)

ZSYR performs the symmetric rank-1 update of a complex symmetric matrix.

Download ZSYR + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZSYR   performs the symmetric rank 1 operation

    A := alpha*x*x**H + A,

 where alpha is a complex scalar, x is an n element vector and A is an
 n by n symmetric matrix.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
           On entry, UPLO specifies whether the upper or lower
           triangular part of the array A is to be referenced as
           follows:

              UPLO = 'U' or 'u'   Only the upper triangular part of A
                                  is to be referenced.

              UPLO = 'L' or 'l'   Only the lower triangular part of A
                                  is to be referenced.

           Unchanged on exit.
[in]N
          N is INTEGER
           On entry, N specifies the order of the matrix A.
           N must be at least zero.
           Unchanged on exit.
[in]ALPHA
          ALPHA is COMPLEX*16
           On entry, ALPHA specifies the scalar alpha.
           Unchanged on exit.
[in]X
          X is COMPLEX*16 array, dimension at least
           ( 1 + ( N - 1 )*abs( INCX ) ).
           Before entry, the incremented array X must contain the N-
           element vector x.
           Unchanged on exit.
[in]INCX
          INCX is INTEGER
           On entry, INCX specifies the increment for the elements of
           X. INCX must not be zero.
           Unchanged on exit.
[in,out]A
          A is COMPLEX*16 array, dimension ( LDA, N )
           Before entry, with  UPLO = 'U' or 'u', the leading n by n
           upper triangular part of the array A must contain the upper
           triangular part of the symmetric matrix and the strictly
           lower triangular part of A is not referenced. On exit, the
           upper triangular part of the array A is overwritten by the
           upper triangular part of the updated matrix.
           Before entry, with UPLO = 'L' or 'l', the leading n by n
           lower triangular part of the array A must contain the lower
           triangular part of the symmetric matrix and the strictly
           upper triangular part of A is not referenced. On exit, the
           lower triangular part of the array A is overwritten by the
           lower triangular part of the updated matrix.
[in]LDA
          LDA is INTEGER
           On entry, LDA specifies the first dimension of A as declared
           in the calling (sub) program. LDA must be at least
           max( 1, N ).
           Unchanged on exit.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
September 2012

Definition at line 137 of file zsyr.f.

137 *
138 * -- LAPACK auxiliary routine (version 3.4.2) --
139 * -- LAPACK is a software package provided by Univ. of Tennessee, --
140 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
141 * September 2012
142 *
143 * .. Scalar Arguments ..
144  CHARACTER uplo
145  INTEGER incx, lda, n
146  COMPLEX*16 alpha
147 * ..
148 * .. Array Arguments ..
149  COMPLEX*16 a( lda, * ), x( * )
150 * ..
151 *
152 * =====================================================================
153 *
154 * .. Parameters ..
155  COMPLEX*16 zero
156  parameter( zero = ( 0.0d+0, 0.0d+0 ) )
157 * ..
158 * .. Local Scalars ..
159  INTEGER i, info, ix, j, jx, kx
160  COMPLEX*16 temp
161 * ..
162 * .. External Functions ..
163  LOGICAL lsame
164  EXTERNAL lsame
165 * ..
166 * .. External Subroutines ..
167  EXTERNAL xerbla
168 * ..
169 * .. Intrinsic Functions ..
170  INTRINSIC max
171 * ..
172 * .. Executable Statements ..
173 *
174 * Test the input parameters.
175 *
176  info = 0
177  IF( .NOT.lsame( uplo, 'U' ) .AND. .NOT.lsame( uplo, 'L' ) ) THEN
178  info = 1
179  ELSE IF( n.LT.0 ) THEN
180  info = 2
181  ELSE IF( incx.EQ.0 ) THEN
182  info = 5
183  ELSE IF( lda.LT.max( 1, n ) ) THEN
184  info = 7
185  END IF
186  IF( info.NE.0 ) THEN
187  CALL xerbla( 'ZSYR ', info )
188  RETURN
189  END IF
190 *
191 * Quick return if possible.
192 *
193  IF( ( n.EQ.0 ) .OR. ( alpha.EQ.zero ) )
194  $ RETURN
195 *
196 * Set the start point in X if the increment is not unity.
197 *
198  IF( incx.LE.0 ) THEN
199  kx = 1 - ( n-1 )*incx
200  ELSE IF( incx.NE.1 ) THEN
201  kx = 1
202  END IF
203 *
204 * Start the operations. In this version the elements of A are
205 * accessed sequentially with one pass through the triangular part
206 * of A.
207 *
208  IF( lsame( uplo, 'U' ) ) THEN
209 *
210 * Form A when A is stored in upper triangle.
211 *
212  IF( incx.EQ.1 ) THEN
213  DO 20 j = 1, n
214  IF( x( j ).NE.zero ) THEN
215  temp = alpha*x( j )
216  DO 10 i = 1, j
217  a( i, j ) = a( i, j ) + x( i )*temp
218  10 CONTINUE
219  END IF
220  20 CONTINUE
221  ELSE
222  jx = kx
223  DO 40 j = 1, n
224  IF( x( jx ).NE.zero ) THEN
225  temp = alpha*x( jx )
226  ix = kx
227  DO 30 i = 1, j
228  a( i, j ) = a( i, j ) + x( ix )*temp
229  ix = ix + incx
230  30 CONTINUE
231  END IF
232  jx = jx + incx
233  40 CONTINUE
234  END IF
235  ELSE
236 *
237 * Form A when A is stored in lower triangle.
238 *
239  IF( incx.EQ.1 ) THEN
240  DO 60 j = 1, n
241  IF( x( j ).NE.zero ) THEN
242  temp = alpha*x( j )
243  DO 50 i = j, n
244  a( i, j ) = a( i, j ) + x( i )*temp
245  50 CONTINUE
246  END IF
247  60 CONTINUE
248  ELSE
249  jx = kx
250  DO 80 j = 1, n
251  IF( x( jx ).NE.zero ) THEN
252  temp = alpha*x( jx )
253  ix = jx
254  DO 70 i = j, n
255  a( i, j ) = a( i, j ) + x( ix )*temp
256  ix = ix + incx
257  70 CONTINUE
258  END IF
259  jx = jx + incx
260  80 CONTINUE
261  END IF
262  END IF
263 *
264  RETURN
265 *
266 * End of ZSYR
267 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine zsyswapr ( character  UPLO,
integer  N,
complex*16, dimension( lda, n )  A,
integer  LDA,
integer  I1,
integer  I2 
)

ZSYSWAPR

Download ZSYSWAPR + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZSYSWAPR applies an elementary permutation on the rows and the columns of
 a symmetric matrix.
Parameters
[in]UPLO
          UPLO is CHARACTER*1
          Specifies whether the details of the factorization are stored
          as an upper or lower triangular matrix.
          = 'U':  Upper triangular, form is A = U*D*U**T;
          = 'L':  Lower triangular, form is A = L*D*L**T.
[in]N
          N is INTEGER
          The order of the matrix A.  N >= 0.
[in,out]A
          A is COMPLEX*16 array, dimension (LDA,N)
          On entry, the NB diagonal matrix D and the multipliers
          used to obtain the factor U or L as computed by ZSYTRF.

          On exit, if INFO = 0, the (symmetric) inverse of the original
          matrix.  If UPLO = 'U', the upper triangular part of the
          inverse is formed and the part of A below the diagonal is not
          referenced; if UPLO = 'L' the lower triangular part of the
          inverse is formed and the part of A above the diagonal is
          not referenced.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,N).
[in]I1
          I1 is INTEGER
          Index of the first row to swap
[in]I2
          I2 is INTEGER
          Index of the second row to swap
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011

Definition at line 104 of file zsyswapr.f.

104 *
105 * -- LAPACK auxiliary routine (version 3.4.0) --
106 * -- LAPACK is a software package provided by Univ. of Tennessee, --
107 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
108 * November 2011
109 *
110 * .. Scalar Arguments ..
111  CHARACTER uplo
112  INTEGER i1, i2, lda, n
113 * ..
114 * .. Array Arguments ..
115  COMPLEX*16 a( lda, n )
116 *
117 * =====================================================================
118 *
119 * ..
120 * .. Local Scalars ..
121  LOGICAL upper
122  INTEGER i
123  COMPLEX*16 tmp
124 *
125 * .. External Functions ..
126  LOGICAL lsame
127  EXTERNAL lsame
128 * ..
129 * .. External Subroutines ..
130  EXTERNAL zswap
131 * ..
132 * .. Executable Statements ..
133 *
134  upper = lsame( uplo, 'U' )
135  IF (upper) THEN
136 *
137 * UPPER
138 * first swap
139 * - swap column I1 and I2 from I1 to I1-1
140  CALL zswap( i1-1, a(1,i1), 1, a(1,i2), 1 )
141 *
142 * second swap :
143 * - swap A(I1,I1) and A(I2,I2)
144 * - swap row I1 from I1+1 to I2-1 with col I2 from I1+1 to I2-1
145  tmp=a(i1,i1)
146  a(i1,i1)=a(i2,i2)
147  a(i2,i2)=tmp
148 *
149  DO i=1,i2-i1-1
150  tmp=a(i1,i1+i)
151  a(i1,i1+i)=a(i1+i,i2)
152  a(i1+i,i2)=tmp
153  END DO
154 *
155 * third swap
156 * - swap row I1 and I2 from I2+1 to N
157  DO i=i2+1,n
158  tmp=a(i1,i)
159  a(i1,i)=a(i2,i)
160  a(i2,i)=tmp
161  END DO
162 *
163  ELSE
164 *
165 * LOWER
166 * first swap
167 * - swap row I1 and I2 from I1 to I1-1
168  CALL zswap( i1-1, a(i1,1), lda, a(i2,1), lda )
169 *
170 * second swap :
171 * - swap A(I1,I1) and A(I2,I2)
172 * - swap col I1 from I1+1 to I2-1 with row I2 from I1+1 to I2-1
173  tmp=a(i1,i1)
174  a(i1,i1)=a(i2,i2)
175  a(i2,i2)=tmp
176 *
177  DO i=1,i2-i1-1
178  tmp=a(i1+i,i1)
179  a(i1+i,i1)=a(i2,i1+i)
180  a(i2,i1+i)=tmp
181  END DO
182 *
183 * third swap
184 * - swap col I1 and I2 from I2+1 to N
185  DO i=i2+1,n
186  tmp=a(i,i1)
187  a(i,i1)=a(i,i2)
188  a(i,i2)=tmp
189  END DO
190 *
191  ENDIF
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine zswap(N, ZX, INCX, ZY, INCY)
ZSWAP
Definition: zswap.f:52

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine ztgsy2 ( character  TRANS,
integer  IJOB,
integer  M,
integer  N,
complex*16, dimension( lda, * )  A,
integer  LDA,
complex*16, dimension( ldb, * )  B,
integer  LDB,
complex*16, dimension( ldc, * )  C,
integer  LDC,
complex*16, dimension( ldd, * )  D,
integer  LDD,
complex*16, dimension( lde, * )  E,
integer  LDE,
complex*16, dimension( ldf, * )  F,
integer  LDF,
double precision  SCALE,
double precision  RDSUM,
double precision  RDSCAL,
integer  INFO 
)

ZTGSY2 solves the generalized Sylvester equation (unblocked algorithm).

Download ZTGSY2 + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZTGSY2 solves the generalized Sylvester equation

             A * R - L * B = scale * C               (1)
             D * R - L * E = scale * F

 using Level 1 and 2 BLAS, where R and L are unknown M-by-N matrices,
 (A, D), (B, E) and (C, F) are given matrix pairs of size M-by-M,
 N-by-N and M-by-N, respectively. A, B, D and E are upper triangular
 (i.e., (A,D) and (B,E) in generalized Schur form).

 The solution (R, L) overwrites (C, F). 0 <= SCALE <= 1 is an output
 scaling factor chosen to avoid overflow.

 In matrix notation solving equation (1) corresponds to solve
 Zx = scale * b, where Z is defined as

        Z = [ kron(In, A)  -kron(B**H, Im) ]             (2)
            [ kron(In, D)  -kron(E**H, Im) ],

 Ik is the identity matrix of size k and X**H is the conjuguate transpose of X.
 kron(X, Y) is the Kronecker product between the matrices X and Y.

 If TRANS = 'C', y in the conjugate transposed system Z**H*y = scale*b
 is solved for, which is equivalent to solve for R and L in

             A**H * R  + D**H * L   = scale * C           (3)
             R  * B**H + L  * E**H  = scale * -F

 This case is used to compute an estimate of Dif[(A, D), (B, E)] =
 = sigma_min(Z) using reverse communicaton with ZLACON.

 ZTGSY2 also (IJOB >= 1) contributes to the computation in ZTGSYL
 of an upper bound on the separation between to matrix pairs. Then
 the input (A, D), (B, E) are sub-pencils of two matrix pairs in
 ZTGSYL.
Parameters
[in]TRANS
          TRANS is CHARACTER*1
          = 'N', solve the generalized Sylvester equation (1).
          = 'T': solve the 'transposed' system (3).
[in]IJOB
          IJOB is INTEGER
          Specifies what kind of functionality to be performed.
          =0: solve (1) only.
          =1: A contribution from this subsystem to a Frobenius
              norm-based estimate of the separation between two matrix
              pairs is computed. (look ahead strategy is used).
          =2: A contribution from this subsystem to a Frobenius
              norm-based estimate of the separation between two matrix
              pairs is computed. (DGECON on sub-systems is used.)
          Not referenced if TRANS = 'T'.
[in]M
          M is INTEGER
          On entry, M specifies the order of A and D, and the row
          dimension of C, F, R and L.
[in]N
          N is INTEGER
          On entry, N specifies the order of B and E, and the column
          dimension of C, F, R and L.
[in]A
          A is COMPLEX*16 array, dimension (LDA, M)
          On entry, A contains an upper triangular matrix.
[in]LDA
          LDA is INTEGER
          The leading dimension of the matrix A. LDA >= max(1, M).
[in]B
          B is COMPLEX*16 array, dimension (LDB, N)
          On entry, B contains an upper triangular matrix.
[in]LDB
          LDB is INTEGER
          The leading dimension of the matrix B. LDB >= max(1, N).
[in,out]C
          C is COMPLEX*16 array, dimension (LDC, N)
          On entry, C contains the right-hand-side of the first matrix
          equation in (1).
          On exit, if IJOB = 0, C has been overwritten by the solution
          R.
[in]LDC
          LDC is INTEGER
          The leading dimension of the matrix C. LDC >= max(1, M).
[in]D
          D is COMPLEX*16 array, dimension (LDD, M)
          On entry, D contains an upper triangular matrix.
[in]LDD
          LDD is INTEGER
          The leading dimension of the matrix D. LDD >= max(1, M).
[in]E
          E is COMPLEX*16 array, dimension (LDE, N)
          On entry, E contains an upper triangular matrix.
[in]LDE
          LDE is INTEGER
          The leading dimension of the matrix E. LDE >= max(1, N).
[in,out]F
          F is COMPLEX*16 array, dimension (LDF, N)
          On entry, F contains the right-hand-side of the second matrix
          equation in (1).
          On exit, if IJOB = 0, F has been overwritten by the solution
          L.
[in]LDF
          LDF is INTEGER
          The leading dimension of the matrix F. LDF >= max(1, M).
[out]SCALE
          SCALE is DOUBLE PRECISION
          On exit, 0 <= SCALE <= 1. If 0 < SCALE < 1, the solutions
          R and L (C and F on entry) will hold the solutions to a
          slightly perturbed system but the input matrices A, B, D and
          E have not been changed. If SCALE = 0, R and L will hold the
          solutions to the homogeneous system with C = F = 0.
          Normally, SCALE = 1.
[in,out]RDSUM
          RDSUM is DOUBLE PRECISION
          On entry, the sum of squares of computed contributions to
          the Dif-estimate under computation by ZTGSYL, where the
          scaling factor RDSCAL (see below) has been factored out.
          On exit, the corresponding sum of squares updated with the
          contributions from the current sub-system.
          If TRANS = 'T' RDSUM is not touched.
          NOTE: RDSUM only makes sense when ZTGSY2 is called by
          ZTGSYL.
[in,out]RDSCAL
          RDSCAL is DOUBLE PRECISION
          On entry, scaling factor used to prevent overflow in RDSUM.
          On exit, RDSCAL is updated w.r.t. the current contributions
          in RDSUM.
          If TRANS = 'T', RDSCAL is not touched.
          NOTE: RDSCAL only makes sense when ZTGSY2 is called by
          ZTGSYL.
[out]INFO
          INFO is INTEGER
          On exit, if INFO is set to
            =0: Successful exit
            <0: If INFO = -i, input argument number i is illegal.
            >0: The matrix pairs (A, D) and (B, E) have common or very
                close eigenvalues.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2015
Contributors:
Bo Kagstrom and Peter Poromaa, Department of Computing Science, Umea University, S-901 87 Umea, Sweden.

Definition at line 261 of file ztgsy2.f.

261 *
262 * -- LAPACK auxiliary routine (version 3.6.0) --
263 * -- LAPACK is a software package provided by Univ. of Tennessee, --
264 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
265 * November 2015
266 *
267 * .. Scalar Arguments ..
268  CHARACTER trans
269  INTEGER ijob, info, lda, ldb, ldc, ldd, lde, ldf, m, n
270  DOUBLE PRECISION rdscal, rdsum, scale
271 * ..
272 * .. Array Arguments ..
273  COMPLEX*16 a( lda, * ), b( ldb, * ), c( ldc, * ),
274  $ d( ldd, * ), e( lde, * ), f( ldf, * )
275 * ..
276 *
277 * =====================================================================
278 *
279 * .. Parameters ..
280  DOUBLE PRECISION zero, one
281  INTEGER ldz
282  parameter( zero = 0.0d+0, one = 1.0d+0, ldz = 2 )
283 * ..
284 * .. Local Scalars ..
285  LOGICAL notran
286  INTEGER i, ierr, j, k
287  DOUBLE PRECISION scaloc
288  COMPLEX*16 alpha
289 * ..
290 * .. Local Arrays ..
291  INTEGER ipiv( ldz ), jpiv( ldz )
292  COMPLEX*16 rhs( ldz ), z( ldz, ldz )
293 * ..
294 * .. External Functions ..
295  LOGICAL lsame
296  EXTERNAL lsame
297 * ..
298 * .. External Subroutines ..
299  EXTERNAL xerbla, zaxpy, zgesc2, zgetc2, zlatdf, zscal
300 * ..
301 * .. Intrinsic Functions ..
302  INTRINSIC dcmplx, dconjg, max
303 * ..
304 * .. Executable Statements ..
305 *
306 * Decode and test input parameters
307 *
308  info = 0
309  ierr = 0
310  notran = lsame( trans, 'N' )
311  IF( .NOT.notran .AND. .NOT.lsame( trans, 'C' ) ) THEN
312  info = -1
313  ELSE IF( notran ) THEN
314  IF( ( ijob.LT.0 ) .OR. ( ijob.GT.2 ) ) THEN
315  info = -2
316  END IF
317  END IF
318  IF( info.EQ.0 ) THEN
319  IF( m.LE.0 ) THEN
320  info = -3
321  ELSE IF( n.LE.0 ) THEN
322  info = -4
323  ELSE IF( lda.LT.max( 1, m ) ) THEN
324  info = -6
325  ELSE IF( ldb.LT.max( 1, n ) ) THEN
326  info = -8
327  ELSE IF( ldc.LT.max( 1, m ) ) THEN
328  info = -10
329  ELSE IF( ldd.LT.max( 1, m ) ) THEN
330  info = -12
331  ELSE IF( lde.LT.max( 1, n ) ) THEN
332  info = -14
333  ELSE IF( ldf.LT.max( 1, m ) ) THEN
334  info = -16
335  END IF
336  END IF
337  IF( info.NE.0 ) THEN
338  CALL xerbla( 'ZTGSY2', -info )
339  RETURN
340  END IF
341 *
342  IF( notran ) THEN
343 *
344 * Solve (I, J) - system
345 * A(I, I) * R(I, J) - L(I, J) * B(J, J) = C(I, J)
346 * D(I, I) * R(I, J) - L(I, J) * E(J, J) = F(I, J)
347 * for I = M, M - 1, ..., 1; J = 1, 2, ..., N
348 *
349  scale = one
350  scaloc = one
351  DO 30 j = 1, n
352  DO 20 i = m, 1, -1
353 *
354 * Build 2 by 2 system
355 *
356  z( 1, 1 ) = a( i, i )
357  z( 2, 1 ) = d( i, i )
358  z( 1, 2 ) = -b( j, j )
359  z( 2, 2 ) = -e( j, j )
360 *
361 * Set up right hand side(s)
362 *
363  rhs( 1 ) = c( i, j )
364  rhs( 2 ) = f( i, j )
365 *
366 * Solve Z * x = RHS
367 *
368  CALL zgetc2( ldz, z, ldz, ipiv, jpiv, ierr )
369  IF( ierr.GT.0 )
370  $ info = ierr
371  IF( ijob.EQ.0 ) THEN
372  CALL zgesc2( ldz, z, ldz, rhs, ipiv, jpiv, scaloc )
373  IF( scaloc.NE.one ) THEN
374  DO 10 k = 1, n
375  CALL zscal( m, dcmplx( scaloc, zero ),
376  $ c( 1, k ), 1 )
377  CALL zscal( m, dcmplx( scaloc, zero ),
378  $ f( 1, k ), 1 )
379  10 CONTINUE
380  scale = scale*scaloc
381  END IF
382  ELSE
383  CALL zlatdf( ijob, ldz, z, ldz, rhs, rdsum, rdscal,
384  $ ipiv, jpiv )
385  END IF
386 *
387 * Unpack solution vector(s)
388 *
389  c( i, j ) = rhs( 1 )
390  f( i, j ) = rhs( 2 )
391 *
392 * Substitute R(I, J) and L(I, J) into remaining equation.
393 *
394  IF( i.GT.1 ) THEN
395  alpha = -rhs( 1 )
396  CALL zaxpy( i-1, alpha, a( 1, i ), 1, c( 1, j ), 1 )
397  CALL zaxpy( i-1, alpha, d( 1, i ), 1, f( 1, j ), 1 )
398  END IF
399  IF( j.LT.n ) THEN
400  CALL zaxpy( n-j, rhs( 2 ), b( j, j+1 ), ldb,
401  $ c( i, j+1 ), ldc )
402  CALL zaxpy( n-j, rhs( 2 ), e( j, j+1 ), lde,
403  $ f( i, j+1 ), ldf )
404  END IF
405 *
406  20 CONTINUE
407  30 CONTINUE
408  ELSE
409 *
410 * Solve transposed (I, J) - system:
411 * A(I, I)**H * R(I, J) + D(I, I)**H * L(J, J) = C(I, J)
412 * R(I, I) * B(J, J) + L(I, J) * E(J, J) = -F(I, J)
413 * for I = 1, 2, ..., M, J = N, N - 1, ..., 1
414 *
415  scale = one
416  scaloc = one
417  DO 80 i = 1, m
418  DO 70 j = n, 1, -1
419 *
420 * Build 2 by 2 system Z**H
421 *
422  z( 1, 1 ) = dconjg( a( i, i ) )
423  z( 2, 1 ) = -dconjg( b( j, j ) )
424  z( 1, 2 ) = dconjg( d( i, i ) )
425  z( 2, 2 ) = -dconjg( e( j, j ) )
426 *
427 *
428 * Set up right hand side(s)
429 *
430  rhs( 1 ) = c( i, j )
431  rhs( 2 ) = f( i, j )
432 *
433 * Solve Z**H * x = RHS
434 *
435  CALL zgetc2( ldz, z, ldz, ipiv, jpiv, ierr )
436  IF( ierr.GT.0 )
437  $ info = ierr
438  CALL zgesc2( ldz, z, ldz, rhs, ipiv, jpiv, scaloc )
439  IF( scaloc.NE.one ) THEN
440  DO 40 k = 1, n
441  CALL zscal( m, dcmplx( scaloc, zero ), c( 1, k ),
442  $ 1 )
443  CALL zscal( m, dcmplx( scaloc, zero ), f( 1, k ),
444  $ 1 )
445  40 CONTINUE
446  scale = scale*scaloc
447  END IF
448 *
449 * Unpack solution vector(s)
450 *
451  c( i, j ) = rhs( 1 )
452  f( i, j ) = rhs( 2 )
453 *
454 * Substitute R(I, J) and L(I, J) into remaining equation.
455 *
456  DO 50 k = 1, j - 1
457  f( i, k ) = f( i, k ) + rhs( 1 )*dconjg( b( k, j ) ) +
458  $ rhs( 2 )*dconjg( e( k, j ) )
459  50 CONTINUE
460  DO 60 k = i + 1, m
461  c( k, j ) = c( k, j ) - dconjg( a( i, k ) )*rhs( 1 ) -
462  $ dconjg( d( i, k ) )*rhs( 2 )
463  60 CONTINUE
464 *
465  70 CONTINUE
466  80 CONTINUE
467  END IF
468  RETURN
469 *
470 * End of ZTGSY2
471 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zgetc2(N, A, LDA, IPIV, JPIV, INFO)
ZGETC2 computes the LU factorization with complete pivoting of the general n-by-n matrix...
Definition: zgetc2.f:113
subroutine zscal(N, ZA, ZX, INCX)
ZSCAL
Definition: zscal.f:54
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine zaxpy(N, ZA, ZX, INCX, ZY, INCY)
ZAXPY
Definition: zaxpy.f:53
subroutine zlatdf(IJOB, N, Z, LDZ, RHS, RDSUM, RDSCAL, IPIV, JPIV)
ZLATDF uses the LU factorization of the n-by-n matrix computed by sgetc2 and computes a contribution ...
Definition: zlatdf.f:171
logical function lde(RI, RJ, LR)
Definition: dblat2.f:2945
subroutine zgesc2(N, A, LDA, RHS, IPIV, JPIV, SCALE)
ZGESC2 solves a system of linear equations using the LU factorization with complete pivoting computed...
Definition: zgesc2.f:117

Here is the call graph for this function:

Here is the caller graph for this function: