LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ ztrsyl()

subroutine ztrsyl ( character  TRANA,
character  TRANB,
integer  ISGN,
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,
double precision  SCALE,
integer  INFO 
)

ZTRSYL

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

Purpose:
 ZTRSYL solves the complex Sylvester matrix equation:

    op(A)*X + X*op(B) = scale*C or
    op(A)*X - X*op(B) = scale*C,

 where op(A) = A or A**H, and A and B are both upper triangular. A is
 M-by-M and B is N-by-N; the right hand side C and the solution X are
 M-by-N; and scale is an output scale factor, set <= 1 to avoid
 overflow in X.
Parameters
[in]TRANA
          TRANA is CHARACTER*1
          Specifies the option op(A):
          = 'N': op(A) = A    (No transpose)
          = 'C': op(A) = A**H (Conjugate transpose)
[in]TRANB
          TRANB is CHARACTER*1
          Specifies the option op(B):
          = 'N': op(B) = B    (No transpose)
          = 'C': op(B) = B**H (Conjugate transpose)
[in]ISGN
          ISGN is INTEGER
          Specifies the sign in the equation:
          = +1: solve op(A)*X + X*op(B) = scale*C
          = -1: solve op(A)*X - X*op(B) = scale*C
[in]M
          M is INTEGER
          The order of the matrix A, and the number of rows in the
          matrices X and C. M >= 0.
[in]N
          N is INTEGER
          The order of the matrix B, and the number of columns in the
          matrices X and C. N >= 0.
[in]A
          A is COMPLEX*16 array, dimension (LDA,M)
          The upper triangular matrix A.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A. LDA >= max(1,M).
[in]B
          B is COMPLEX*16 array, dimension (LDB,N)
          The upper triangular matrix B.
[in]LDB
          LDB is INTEGER
          The leading dimension of the array B. LDB >= max(1,N).
[in,out]C
          C is COMPLEX*16 array, dimension (LDC,N)
          On entry, the M-by-N right hand side matrix C.
          On exit, C is overwritten by the solution matrix X.
[in]LDC
          LDC is INTEGER
          The leading dimension of the array C. LDC >= max(1,M)
[out]SCALE
          SCALE is DOUBLE PRECISION
          The scale factor, scale, set <= 1 to avoid overflow in X.
[out]INFO
          INFO is INTEGER
          = 0: successful exit
          < 0: if INFO = -i, the i-th argument had an illegal value
          = 1: A and B have common or very close eigenvalues; perturbed
               values were used to solve the equation (but the matrices
               A and B are unchanged).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 155 of file ztrsyl.f.

157 *
158 * -- LAPACK computational routine --
159 * -- LAPACK is a software package provided by Univ. of Tennessee, --
160 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
161 *
162 * .. Scalar Arguments ..
163  CHARACTER TRANA, TRANB
164  INTEGER INFO, ISGN, LDA, LDB, LDC, M, N
165  DOUBLE PRECISION SCALE
166 * ..
167 * .. Array Arguments ..
168  COMPLEX*16 A( LDA, * ), B( LDB, * ), C( LDC, * )
169 * ..
170 *
171 * =====================================================================
172 *
173 * .. Parameters ..
174  DOUBLE PRECISION ONE
175  parameter( one = 1.0d+0 )
176 * ..
177 * .. Local Scalars ..
178  LOGICAL NOTRNA, NOTRNB
179  INTEGER J, K, L
180  DOUBLE PRECISION BIGNUM, DA11, DB, EPS, SCALOC, SGN, SMIN,
181  $ SMLNUM
182  COMPLEX*16 A11, SUML, SUMR, VEC, X11
183 * ..
184 * .. Local Arrays ..
185  DOUBLE PRECISION DUM( 1 )
186 * ..
187 * .. External Functions ..
188  LOGICAL LSAME
189  DOUBLE PRECISION DLAMCH, ZLANGE
190  COMPLEX*16 ZDOTC, ZDOTU, ZLADIV
191  EXTERNAL lsame, dlamch, zlange, zdotc, zdotu, zladiv
192 * ..
193 * .. External Subroutines ..
194  EXTERNAL dlabad, xerbla, zdscal
195 * ..
196 * .. Intrinsic Functions ..
197  INTRINSIC abs, dble, dcmplx, dconjg, dimag, max, min
198 * ..
199 * .. Executable Statements ..
200 *
201 * Decode and Test input parameters
202 *
203  notrna = lsame( trana, 'N' )
204  notrnb = lsame( tranb, 'N' )
205 *
206  info = 0
207  IF( .NOT.notrna .AND. .NOT.lsame( trana, 'C' ) ) THEN
208  info = -1
209  ELSE IF( .NOT.notrnb .AND. .NOT.lsame( tranb, 'C' ) ) THEN
210  info = -2
211  ELSE IF( isgn.NE.1 .AND. isgn.NE.-1 ) THEN
212  info = -3
213  ELSE IF( m.LT.0 ) THEN
214  info = -4
215  ELSE IF( n.LT.0 ) THEN
216  info = -5
217  ELSE IF( lda.LT.max( 1, m ) ) THEN
218  info = -7
219  ELSE IF( ldb.LT.max( 1, n ) ) THEN
220  info = -9
221  ELSE IF( ldc.LT.max( 1, m ) ) THEN
222  info = -11
223  END IF
224  IF( info.NE.0 ) THEN
225  CALL xerbla( 'ZTRSYL', -info )
226  RETURN
227  END IF
228 *
229 * Quick return if possible
230 *
231  scale = one
232  IF( m.EQ.0 .OR. n.EQ.0 )
233  $ RETURN
234 *
235 * Set constants to control overflow
236 *
237  eps = dlamch( 'P' )
238  smlnum = dlamch( 'S' )
239  bignum = one / smlnum
240  CALL dlabad( smlnum, bignum )
241  smlnum = smlnum*dble( m*n ) / eps
242  bignum = one / smlnum
243  smin = max( smlnum, eps*zlange( 'M', m, m, a, lda, dum ),
244  $ eps*zlange( 'M', n, n, b, ldb, dum ) )
245  sgn = isgn
246 *
247  IF( notrna .AND. notrnb ) THEN
248 *
249 * Solve A*X + ISGN*X*B = scale*C.
250 *
251 * The (K,L)th block of X is determined starting from
252 * bottom-left corner column by column by
253 *
254 * A(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
255 *
256 * Where
257 * M L-1
258 * R(K,L) = SUM [A(K,I)*X(I,L)] +ISGN*SUM [X(K,J)*B(J,L)].
259 * I=K+1 J=1
260 *
261  DO 30 l = 1, n
262  DO 20 k = m, 1, -1
263 *
264  suml = zdotu( m-k, a( k, min( k+1, m ) ), lda,
265  $ c( min( k+1, m ), l ), 1 )
266  sumr = zdotu( l-1, c( k, 1 ), ldc, b( 1, l ), 1 )
267  vec = c( k, l ) - ( suml+sgn*sumr )
268 *
269  scaloc = one
270  a11 = a( k, k ) + sgn*b( l, l )
271  da11 = abs( dble( a11 ) ) + abs( dimag( a11 ) )
272  IF( da11.LE.smin ) THEN
273  a11 = smin
274  da11 = smin
275  info = 1
276  END IF
277  db = abs( dble( vec ) ) + abs( dimag( vec ) )
278  IF( da11.LT.one .AND. db.GT.one ) THEN
279  IF( db.GT.bignum*da11 )
280  $ scaloc = one / db
281  END IF
282  x11 = zladiv( vec*dcmplx( scaloc ), a11 )
283 *
284  IF( scaloc.NE.one ) THEN
285  DO 10 j = 1, n
286  CALL zdscal( m, scaloc, c( 1, j ), 1 )
287  10 CONTINUE
288  scale = scale*scaloc
289  END IF
290  c( k, l ) = x11
291 *
292  20 CONTINUE
293  30 CONTINUE
294 *
295  ELSE IF( .NOT.notrna .AND. notrnb ) THEN
296 *
297 * Solve A**H *X + ISGN*X*B = scale*C.
298 *
299 * The (K,L)th block of X is determined starting from
300 * upper-left corner column by column by
301 *
302 * A**H(K,K)*X(K,L) + ISGN*X(K,L)*B(L,L) = C(K,L) - R(K,L)
303 *
304 * Where
305 * K-1 L-1
306 * R(K,L) = SUM [A**H(I,K)*X(I,L)] + ISGN*SUM [X(K,J)*B(J,L)]
307 * I=1 J=1
308 *
309  DO 60 l = 1, n
310  DO 50 k = 1, m
311 *
312  suml = zdotc( k-1, a( 1, k ), 1, c( 1, l ), 1 )
313  sumr = zdotu( l-1, c( k, 1 ), ldc, b( 1, l ), 1 )
314  vec = c( k, l ) - ( suml+sgn*sumr )
315 *
316  scaloc = one
317  a11 = dconjg( a( k, k ) ) + sgn*b( l, l )
318  da11 = abs( dble( a11 ) ) + abs( dimag( a11 ) )
319  IF( da11.LE.smin ) THEN
320  a11 = smin
321  da11 = smin
322  info = 1
323  END IF
324  db = abs( dble( vec ) ) + abs( dimag( vec ) )
325  IF( da11.LT.one .AND. db.GT.one ) THEN
326  IF( db.GT.bignum*da11 )
327  $ scaloc = one / db
328  END IF
329 *
330  x11 = zladiv( vec*dcmplx( scaloc ), a11 )
331 *
332  IF( scaloc.NE.one ) THEN
333  DO 40 j = 1, n
334  CALL zdscal( m, scaloc, c( 1, j ), 1 )
335  40 CONTINUE
336  scale = scale*scaloc
337  END IF
338  c( k, l ) = x11
339 *
340  50 CONTINUE
341  60 CONTINUE
342 *
343  ELSE IF( .NOT.notrna .AND. .NOT.notrnb ) THEN
344 *
345 * Solve A**H*X + ISGN*X*B**H = C.
346 *
347 * The (K,L)th block of X is determined starting from
348 * upper-right corner column by column by
349 *
350 * A**H(K,K)*X(K,L) + ISGN*X(K,L)*B**H(L,L) = C(K,L) - R(K,L)
351 *
352 * Where
353 * K-1
354 * R(K,L) = SUM [A**H(I,K)*X(I,L)] +
355 * I=1
356 * N
357 * ISGN*SUM [X(K,J)*B**H(L,J)].
358 * J=L+1
359 *
360  DO 90 l = n, 1, -1
361  DO 80 k = 1, m
362 *
363  suml = zdotc( k-1, a( 1, k ), 1, c( 1, l ), 1 )
364  sumr = zdotc( n-l, c( k, min( l+1, n ) ), ldc,
365  $ b( l, min( l+1, n ) ), ldb )
366  vec = c( k, l ) - ( suml+sgn*dconjg( sumr ) )
367 *
368  scaloc = one
369  a11 = dconjg( a( k, k )+sgn*b( l, l ) )
370  da11 = abs( dble( a11 ) ) + abs( dimag( a11 ) )
371  IF( da11.LE.smin ) THEN
372  a11 = smin
373  da11 = smin
374  info = 1
375  END IF
376  db = abs( dble( vec ) ) + abs( dimag( vec ) )
377  IF( da11.LT.one .AND. db.GT.one ) THEN
378  IF( db.GT.bignum*da11 )
379  $ scaloc = one / db
380  END IF
381 *
382  x11 = zladiv( vec*dcmplx( scaloc ), a11 )
383 *
384  IF( scaloc.NE.one ) THEN
385  DO 70 j = 1, n
386  CALL zdscal( m, scaloc, c( 1, j ), 1 )
387  70 CONTINUE
388  scale = scale*scaloc
389  END IF
390  c( k, l ) = x11
391 *
392  80 CONTINUE
393  90 CONTINUE
394 *
395  ELSE IF( notrna .AND. .NOT.notrnb ) THEN
396 *
397 * Solve A*X + ISGN*X*B**H = C.
398 *
399 * The (K,L)th block of X is determined starting from
400 * bottom-left corner column by column by
401 *
402 * A(K,K)*X(K,L) + ISGN*X(K,L)*B**H(L,L) = C(K,L) - R(K,L)
403 *
404 * Where
405 * M N
406 * R(K,L) = SUM [A(K,I)*X(I,L)] + ISGN*SUM [X(K,J)*B**H(L,J)]
407 * I=K+1 J=L+1
408 *
409  DO 120 l = n, 1, -1
410  DO 110 k = m, 1, -1
411 *
412  suml = zdotu( m-k, a( k, min( k+1, m ) ), lda,
413  $ c( min( k+1, m ), l ), 1 )
414  sumr = zdotc( n-l, c( k, min( l+1, n ) ), ldc,
415  $ b( l, min( l+1, n ) ), ldb )
416  vec = c( k, l ) - ( suml+sgn*dconjg( sumr ) )
417 *
418  scaloc = one
419  a11 = a( k, k ) + sgn*dconjg( b( l, l ) )
420  da11 = abs( dble( a11 ) ) + abs( dimag( a11 ) )
421  IF( da11.LE.smin ) THEN
422  a11 = smin
423  da11 = smin
424  info = 1
425  END IF
426  db = abs( dble( vec ) ) + abs( dimag( vec ) )
427  IF( da11.LT.one .AND. db.GT.one ) THEN
428  IF( db.GT.bignum*da11 )
429  $ scaloc = one / db
430  END IF
431 *
432  x11 = zladiv( vec*dcmplx( scaloc ), a11 )
433 *
434  IF( scaloc.NE.one ) THEN
435  DO 100 j = 1, n
436  CALL zdscal( m, scaloc, c( 1, j ), 1 )
437  100 CONTINUE
438  scale = scale*scaloc
439  END IF
440  c( k, l ) = x11
441 *
442  110 CONTINUE
443  120 CONTINUE
444 *
445  END IF
446 *
447  RETURN
448 *
449 * End of ZTRSYL
450 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
complex *16 function zdotu(N, ZX, INCX, ZY, INCY)
ZDOTU
Definition: zdotu.f:83
subroutine zdscal(N, DA, ZX, INCX)
ZDSCAL
Definition: zdscal.f:78
complex *16 function zdotc(N, ZX, INCX, ZY, INCY)
ZDOTC
Definition: zdotc.f:83
double precision function zlange(NORM, M, N, A, LDA, WORK)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: zlange.f:115
complex *16 function zladiv(X, Y)
ZLADIV performs complex division in real arithmetic, avoiding unnecessary overflow.
Definition: zladiv.f:64
Here is the call graph for this function:
Here is the caller graph for this function: