Developer Reference for Intel® oneAPI Math Kernel Library for Fortran

ID 766686
Date 7/13/2023
Public

A newer version of this document is available. Customers should click here to go to the newest version.

Document Table of Contents

?trsyl

Solves Sylvester equation for real quasi-triangular or complex triangular matrices.

Syntax

call strsyl(trana, tranb, isgn, m, n, a, lda, b, ldb, c, ldc, scale, info)

call dtrsyl(trana, tranb, isgn, m, n, a, lda, b, ldb, c, ldc, scale, info)

call ctrsyl(trana, tranb, isgn, m, n, a, lda, b, ldb, c, ldc, scale, info)

call ztrsyl(trana, tranb, isgn, m, n, a, lda, b, ldb, c, ldc, scale, info)

call trsyl(a, b, c, scale [, trana] [,tranb] [,isgn] [,info])

Include Files

  • mkl.fi, lapack.f90

Description

The routine solves the Sylvester matrix equation op(A)*X±X*op(B) = α*C, where op(A) = A or AH, and the matrices A and B are upper triangular (or, for real flavors, upper quasi-triangular in canonical Schur form); α 1 is a scale factor determined by the routine to avoid overflow in X; A is m-by-m, B is n-by-n, and C and X are both m-by-n. The matrix X is obtained by a straightforward process of back substitution.

The equation has a unique solution if and only if αi±βi 0, where {αi} and {βi} are the eigenvalues of A and B, respectively, and the sign (+ or -) is the same as that used in the equation to be solved.

Input Parameters

trana

CHARACTER*1. Must be 'N' or 'T' or 'C'.

If trana = 'N', then op(A) = A.

If trana = 'T', then op(A) = AT (real flavors only).

If trana = 'C' then op(A) = AH.

tranb

CHARACTER*1. Must be 'N' or 'T' or 'C'.

If tranb = 'N', then op(B) = B.

If tranb = 'T', then op(B) = BT (real flavors only).

If tranb = 'C', then op(B) = BH.

isgn

INTEGER. Indicates the form of the Sylvester equation.

If isgn = +1, op(A)*X + X*op(B) = alpha*C.

If isgn = -1, op(A)*X - X*op(B) = alpha*C.

m

INTEGER. The order of A, and the number of rows in X and C (m 0).

n

INTEGER. The order of B, and the number of columns in X and C (n 0).

a, b, c

REAL for strsyl

DOUBLE PRECISION for dtrsyl

COMPLEX for ctrsyl

DOUBLE COMPLEX for ztrsyl.

Arrays:

a(lda,*) contains the matrix A.

The second dimension of a must be at least max(1, m).

b(ldb,*) contains the matrix B.

The second dimension of b must be at least max(1, n).

c(ldc,*) contains the matrix C.

The second dimension of c must be at least max(1, n).

lda

INTEGER. The leading dimension of a; at least max(1, m).

ldb

INTEGER. The leading dimension of b; at least max(1, n).

ldc

INTEGER. The leading dimension of c; at least max(1, m) .

Output Parameters

c

Overwritten by the solution matrix X.

scale

REAL for single-precision flavors

DOUBLE PRECISION for double-precision flavors.

The value of the scale factor α.

info

INTEGER.

If info = 0, the execution is successful.

If info = -i, the i-th parameter had an illegal value.

If info = 1, A and B have common or close eigenvalues; perturbed values were used to solve the equation.

LAPACK 95 Interface Notes

Routines in Fortran 95 interface have fewer arguments in the calling sequence than their FORTRAN 77 counterparts. For general conventions applied to skip redundant or restorable arguments, see LAPACK 95 Interface Conventions.

Specific details for the routine trsyl interface are the following:

a

Holds the matrix A of size (m,m).

b

Holds the matrix B of size (n,n).

c

Holds the matrix C of size (m,n).

trana

Must be 'N', 'C', or 'T'. The default value is 'N'.

tranb

Must be 'N', 'C', or 'T'. The default value is 'N'.

isgn

Must be +1 or -1. The default value is +1.

Application Notes

Let X be the exact, Y the corresponding computed solution, and R the residual matrix: R = C - (AY±YB). Then the residual is always small:

||R||F = O(ε)*(||A||F +||B||F)*||Y||F.

However, Y is not necessarily the exact solution of a slightly perturbed equation; in other words, the solution is not backwards stable.

For the forward error, the following bound holds:

||Y - X||F||R||F/sep(A,B)

but this may be a considerable overestimate. See [Golub96] for a definition of sep(A, B).

The approximate number of floating-point operations for real flavors is m*n*(m + n). For complex flavors it is 4 times greater.