sgltrsylv - Man Page

Name

sgltrsylv — Single Precision

— Single Precision routines for triangular standard Sylvester equations.

subroutine sla_trsylv2_dag (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-3 Bartels-Stewart Algorithm for the discrete time Sylvester equation with DAG based parallelization.
subroutine sla_trsylv_dag (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-3 Bartels-Stewart Algorithm for the generalized Sylvester equation with DAG parallelization.
subroutine sla_trsylv2_kernel_44nn (sgn, m, n, a, lda, b, ldb, x, ldx, scale, info)
Solver for a 4x4 discrete time Sylvester equation (TRANSA = N, TRANSB = N)
subroutine sla_trsylv2_kernel_44nt (sgn, m, n, a, lda, b, ldb, x, ldx, scale, info)
Solver for a 4x4 discrete time Sylvester equation (TRANSA = N, TRANSB = T)
subroutine sla_trsylv2_kernel_44tn (sgn, m, n, a, lda, b, ldb, x, ldx, scale, info)
Solver for a 4x4 discrete time Sylvester equation (TRANSA = T, TRANSB = N)
subroutine sla_trsylv2_kernel_44tt (sgn, m, n, a, lda, b, ldb, x, ldx, scale, info)
Solver for a 4x4 discrete time Sylvester equation (TRANSA = T, TRANSB = T)
subroutine sla_trsylv_kernel_44nn (sgn, m, n, a, lda, b, ldb, x, ldx, scale, info)
Solver for a 4x4 Sylvester equation (TRANSA = N, TRANSB = N)
subroutine sla_trsylv_kernel_44nt (sgn, m, n, a, lda, b, ldb, x, ldx, scale, info)
Solver for a 4x4 Sylvester equation (TRANSA = N, TRANSB = T)
subroutine sla_trsylv_kernel_44tn (sgn, m, n, a, lda, b, ldb, x, ldx, scale, info)
Solver for a 4x4 Sylvester equation (TRANSA = T, TRANSB = N)
subroutine sla_trsylv_kernel_44tt (sgn, m, n, a, lda, b, ldb, x, ldx, scale, info)
Solver for a 4x4 Sylvester equation (TRANSA = T, TRANSB = T)
subroutine sla_trsylv2_l2 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)
subroutine sla_trsylv2_l2_local_copy (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)
subroutine sla_trsylv2_l2_local_copy_128 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)
subroutine sla_trsylv2_l2_local_copy_32 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)
subroutine sla_trsylv2_l2_local_copy_64 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)
subroutine sla_trsylv2_l2_local_copy_96 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)
subroutine sla_trsylv2_l2_reorder (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)
subroutine sla_trsylv2_l2_unopt (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (unoptimized)
subroutine sla_trsylv_l2 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the Sylvester equation (optimized)
subroutine sla_trsylv_l2_local_copy (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the Sylvester equation (optimized)
subroutine sla_trsylv_l2_local_copy_128 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the generalized Sylvester equation (Optimized)
subroutine sla_trsylv_l2_local_copy_32 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the generalized Sylvester equation (Optimized)
subroutine sla_trsylv_l2_local_copy_64 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the generalized Sylvester equation (Optimized)
subroutine sla_trsylv_l2_local_copy_96 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the generalized Sylvester equation (Optimized)
subroutine sla_trsylv_l2_reorder (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the Sylvester equation (optimized)
subroutine sla_trsylv_l2_unopt (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-2 Bartels-Stewart Algorithm for the Sylvester equation (unoptimized)
subroutine sla_trsylv2_l3_2s (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Two level blocked Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)
subroutine sla_trsylv2_l3 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-3 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)
subroutine sla_trsylv2_l3_unopt (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-3 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Not Optimized)
subroutine sla_trsylv_l3_2s (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-3 Bartels-Stewart Algorithm for the generalized Sylvester equation.
subroutine sla_trsylv_l3 (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-3 Bartels-Stewart Algorithm for the generalized Sylvester equation.
subroutine sla_trsylv_l3_unopt (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-3 Bartels-Stewart Algorithm for the generalized Sylvester equation (unoptimized)
recursive subroutine sla_trsylv2_recursive (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-3 Recursive Blocked Algorithm for the discrete time Sylvester equation.
recursive subroutine sla_trsylv_recursive (transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)
Level-3 Recursive Blocked Solver for the Sylvester equation.

Detailed Description

Single Precision routines for triangular standard Sylvester equations.

This section contains the solvers for the standard Sylvester equation with (quasi) triangular coefficient matrices. The coefficient matrices are normally generated with the help of the Schur decomposition from LAPACK. All codes use single precision arithmetic.

Function Documentation

subroutine sla_trsylv2_dag (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, m) a, integer lda, real, dimension(ldb, n) b, integer ldb, real, dimension(ldx, n) x, integer ldx, real scale, real, dimension(*) work, integer info)

Level-3 Bartels-Stewart Algorithm for the discrete time Sylvester equation with DAG based parallelization.

Purpose:

 SLA_TRSYLV2_DAG solves a discrete time Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Remarks

The algorithm is implemented using BLAS level 3 operations and OpenMP 4.0 DAG scheduling.

Attention

Due to the parallel nature of the algorithm the scaling is not applied to the right hand. If the problem is ill-posed and scaling appears you have to solve the equation again with a solver with complete scaling support like the SLA_TRSYLV2_L3 routine.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDX >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.
          If SCALE .NE. 1 the problem is no solved correctly in this case
          one have to use an other solver.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 181 of file sla_trsylv2_dag.f90.

subroutine sla_trsylv2_kernel_44nn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

Solver for a 4x4 discrete time Sylvester equation (TRANSA = N, TRANSB = N)

Purpose:

 SLA_TRSYLV2_KERNEL_44NN solves a discrete time Sylvester equation of the following form

    A * X * B + SGN * X  = SCALE * Y                              (1)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N quasi upper
 triangular matrix.  The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.
 The algorithm is implemented using BLAS level 2
 operations. Thereby the order of M and N is at most 4. Furthermore, for fast execution
 the function does not check the input arguments.

Parameters

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.

          M is INTEGER
          The order of the matrix A.  4 >= M >= 0.

          N is INTEGER
          The order of the matrix B.  4 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1)
          as selected by SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCAL <= 1 holds true.

INFO

          INFO is INTEGER
          On input:
          = 1 :  Skip the input data checks
          <> 1:  Check input data like normal LAPACK like routines.
          On output:
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 135 of file sla_trsylv2_kernel_44_nn.f90.

subroutine sla_trsylv2_kernel_44nt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

Solver for a 4x4 discrete time Sylvester equation (TRANSA = N, TRANSB = T)

Purpose:

 SLA_TRSYLV2_KERNEL_44NT solves a discrete time Sylvester equation of the following form

    A * X * B**T + SGN * X  = SCALE * Y                              (1)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N quasi upper
 triangular matrix.  The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.
 The algorithm is implemented using BLAS level 2
 operations. Thereby the order of M and N is at most 4. Furthermore, for fast execution
 the function does not check the input arguments.

Parameters

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.

          M is INTEGER
          The order of the matrix A.  4 >= M >= 0.

          N is INTEGER
          The order of the matrix B.  4 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1)
          as selected by SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCAL <= 1 holds true.

INFO

          INFO is INTEGER
          On input:
          = 1 :  Skip the input data checks
          <> 1:  Check input data like normal LAPACK like routines.
          On output:
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 134 of file sla_trsylv2_kernel_44_nt.f90.

subroutine sla_trsylv2_kernel_44tn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

Solver for a 4x4 discrete time Sylvester equation (TRANSA = T, TRANSB = N)

Purpose:

 SLA_TRSYLV2_KERNEL_44TN solves a discrete time Sylvester equation of the following form

    A**T * X * B + SGN * X  = SCALE * Y                              (1)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N quasi upper
 triangular matrix.  The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.
 The algorithm is implemented using BLAS level 2
 operations. Thereby the order of M and N is at most 4. Furthermore, for fast execution
 the function does not check the input arguments.

Parameters

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.

          M is INTEGER
          The order of the matrix A.  4 >= M >= 0.

          N is INTEGER
          The order of the matrix B.  4 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1)
          as selected by SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCAL <= 1 holds true.

INFO

          INFO is INTEGER
          On input:
          = 1 :  Skip the input data checks
          <> 1:  Check input data like normal LAPACK like routines.
          On output:
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 135 of file sla_trsylv2_kernel_44_tn.f90.

subroutine sla_trsylv2_kernel_44tt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

Solver for a 4x4 discrete time Sylvester equation (TRANSA = T, TRANSB = T)

Purpose:

 SLA_TRSYLV2_KERNEL_44TT solves a discrete time Sylvester equation of the following form

    A**T * X * B**T + SGN * X  = SCALE * Y                              (1)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N quasi upper
 triangular matrix.  The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.
 The algorithm is implemented using BLAS level 2
 operations. Thereby the order of M and N is at most 4. Furthermore, for fast execution
 the function does not check the input arguments.

Parameters

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.

          M is INTEGER
          The order of the matrix A.  4 >= M >= 0.

          N is INTEGER
          The order of the matrix B.  4 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1)
          as selected by SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCAL <= 1 holds true.

INFO

          INFO is INTEGER
          On input:
          = 1 :  Skip the input data checks
          <> 1:  Check input data like normal LAPACK like routines.
          On output:
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 135 of file sla_trsylv2_kernel_44_tt.f90.

subroutine sla_trsylv2_l2 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV2_L2 solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV and SGER by Fortran intrinsic,

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 182 of file sla_trsylv2_l2.f90.

subroutine sla_trsylv2_l2_local_copy (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV2_L2_LOCAL_COPY solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV and SGER by Fortran intrinsic,
Use local copies of A, B, C, D, and X (M, N <=128) .

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 184 of file sla_trsylv2_l2_opt_local_copy.f90.

subroutine sla_trsylv2_l2_local_copy_128 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV2_L2_LOCAL_COPY_128 solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Attention

The size of the Problem is limited by M,N <= 128

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV by Fortran intrinsic,
Use local copies of A, B, C, D, and X.
Align local copies and fix problem size to <= 128

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 184 of file sla_trsylv2_l2_opt_local_copy_128.f90.

subroutine sla_trsylv2_l2_local_copy_32 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV2_L2_LOCAL_COPY_32 solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Attention

The size of the Problem is limited by M,N <= 32

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV by Fortran intrinsic,
Use local copies of A, B, C, D, and X.
Align local copies and fix problem size to <= 32

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 184 of file sla_trsylv2_l2_opt_local_copy_32.f90.

subroutine sla_trsylv2_l2_local_copy_64 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV2_L2_LOCAL_COPY_64 solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Attention

The size of the Problem is limited by M,N <= 64

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV by Fortran intrinsic,
Use local copies of A, B, C, D, and X.
Align local copies and fix problem size to <= 64

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 184 of file sla_trsylv2_l2_opt_local_copy_64.f90.

subroutine sla_trsylv2_l2_local_copy_96 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV2_L2_LOCAL_COPY_96 solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Attention

The size of the Problem is limited by M,N <= 96

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV by Fortran intrinsic,
Use local copies of A, B, C, D, and X.
Align local copies and fix problem size to <= 96

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 184 of file sla_trsylv2_l2_opt_local_copy_96.f90.

subroutine sla_trsylv2_l2_reorder (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV2_L2_UNOPT solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV and SGER by Fortran intrinsic,

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 183 of file sla_trsylv2_l2_opt_reorder.f90.

subroutine sla_trsylv2_l2_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-2 Bartels-Stewart Algorithm for the discrete time Sylvester equation (unoptimized)

Purpose:

 SLA_TRSYLV2_L2_UNOPT solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Attention

The algorithm is implemented using BLAS level 2 operations without further optimizations. For a faster implementation see SLA_TRSYLV2_L2.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Nothing.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 181 of file sla_trsylv2_l2_unopt.f90.

subroutine sla_trsylv2_l3 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-3 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV2_L3 solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 173 of file sla_trsylv2_l3_opt.f90.

subroutine sla_trsylv2_l3_2s (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Two level blocked Bartels-Stewart Algorithm for the discrete time Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV2_L3_2S solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Remarks

The algorithm used level-3 BLAS operations and a DAG scheduled inner solver.

Attention

Due to the parallel nature of the inner solvers the scaling is turned off and SCALE is set to ONE.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 174 of file sla_trsylv2_l3_2stage.f90.

subroutine sla_trsylv2_l3_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-3 Bartels-Stewart Algorithm for the discrete time Sylvester equation (Not Optimized)

Purpose:

 SLA_TRSYLV2_L3_UNOPT solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Attention

This function iterates column first over the result and from this fact it will be a bit slower than SLA_TRSYLV2_L3.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 175 of file sla_trsylv2_l3_unopt.f90.

recursive subroutine sla_trsylv2_recursive (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-3 Recursive Blocked Algorithm for the discrete time Sylvester equation.

Purpose:

 SLA_TRSYLV2_RECURSIVE solves a generalized Sylvester equation of the following forms

    op1(A) * X * op2(B) +  X  = SCALE * Y                              (1)

 or

    op1(A) * X * op2(B) -  X  = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix and B is a N-by-N upper
 quasi triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are created by SGEES from LAPACK.

Remarks

The algorithm uses recursive blocking instead of the Bartels-Stewart approach.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension M*N
          Workspace for the algorithm

INFO

          INFO is INTEGER
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 163 of file sla_trsylv2_recursive.f90.

subroutine sla_trsylv_dag (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, n) x, integer ldx, real scale, real, dimension(*) work, integer info)

Level-3 Bartels-Stewart Algorithm for the generalized Sylvester equation with DAG parallelization.

Purpose:

 SLA_TRSYLV_DAG solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Remarks

The algorithm is implemented using BLAS level 3 operations and OpenMP 4.0 DAG Scheduling.

Attention

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 177 of file sla_trsylv_dag.f90.

subroutine sla_trsylv_kernel_44nn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

Solver for a 4x4 Sylvester equation (TRANSA = N, TRANSB = N)

Purpose:

 SLA_TRSYLV_KERNEL_44NN solves a Sylvester equation of the following form

    A * X  + SGN * X * B = SCALE * Y                              (1)

 where A is a M-by-M quasi upper triangular matrix and B is N-by-N quasi
 upper triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are create by SGEES from LAPACK.
 The algorithm is implemented without using BLAS level 2
 operations. Thereby the order of M and N is at most 4. Furthermore, for fast execution
 the function does not check the input arguments.

Parameters

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.

          M is INTEGER
          The order of the matrices A and C.  4 >= M >= 0.

          N is INTEGER
          The order of the matrices B and D.  4 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCAL <= 1 holds true.

INFO

          INFO is INTEGER
          On output:
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 134 of file sla_trsylv_kernel_44_nn.f90.

subroutine sla_trsylv_kernel_44nt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

Solver for a 4x4 Sylvester equation (TRANSA = N, TRANSB = T)

Purpose:

 SLA_TRSYLV_KERNEL_44NT solves a Sylvester equation of the following form

    A * X  + SGN * X * B**T = SCALE * Y                              (1)

 where A is a M-by-M quasi upper triangular matrix and B is N-by-N quasi
 upper triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are create by SGEES from LAPACK.
 The algorithm is implemented using without BLAS level 2
 operations. Thereby the order of M and N is at most 4. Furthermore, for fast execution
 the function does not check the input arguments.

Parameters

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.

          M is INTEGER
          The order of the matrices A and C.  4 >= M >= 0.

          N is INTEGER
          The order of the matrices B and D.  4 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCAL <= 1 holds true.

INFO

          INFO is INTEGER
          On output:
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 133 of file sla_trsylv_kernel_44_nt.f90.

subroutine sla_trsylv_kernel_44tn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

Solver for a 4x4 Sylvester equation (TRANSA = T, TRANSB = N)

Purpose:

 SLA_TRSYLV_KERNEL_44TN solves a Sylvester equation of the following form

    A**T * X  + SGN * X * B = SCALE * Y                              (1)

 where A is a M-by-M quasi upper triangular matrix and B is N-by-N quasi
 upper triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are create by SGEES from LAPACK.
 The algorithm is implemented without using BLAS level 2
 operations. Thereby the order of M and N is at most 4. Furthermore, for fast execution
 the function does not check the input arguments.

Parameters

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.

          M is INTEGER
          The order of the matrices A and C.  4 >= M >= 0.

          N is INTEGER
          The order of the matrices B and D.  4 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCAL <= 1 holds true.

INFO

          INFO is INTEGER
          On output:
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 133 of file sla_trsylv_kernel_44_tn.f90.

subroutine sla_trsylv_kernel_44tt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

Solver for a 4x4 Sylvester equation (TRANSA = T, TRANSB = T)

Purpose:

 SLA_TRSYLV_KERNEL_44TT solves a Sylvester equation of the following form

    A**T * X  + SGN * X * B**T = SCALE * Y                              (1)

 where A is a M-by-M quasi upper triangular matrix and B is N-by-N quasi
 upper triangular matrices. The right hand side Y and the solution X
 M-by-N matrices. Typically the matrices A and B are create by SGEES from LAPACK.
 The algorithm is implemented without using BLAS level 2
 operations. Thereby the order of M and N is at most 4. Furthermore, for fast execution
 the function does not check the input arguments.

Parameters

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.

          M is INTEGER
          The order of the matrices A and C.  4 >= M >= 0.

          N is INTEGER
          The order of the matrices B and D.  4 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular. If the matrix D is already
          quasi-upper triangular the matrix B must be upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCAL <= 1 holds true.

INFO

          INFO is INTEGER
          On output:
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 133 of file sla_trsylv_kernel_44_tt.f90.

subroutine sla_trsylv_l2 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-2 Bartels-Stewart Algorithm for the Sylvester equation (optimized)

Purpose:

 SLA_TRSYLV_L2 solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV and SGER by Fortran intrinsic,

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 180 of file sla_trsylv_l2.f90.

subroutine sla_trsylv_l2_local_copy (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-2 Bartels-Stewart Algorithm for the Sylvester equation (optimized)

Purpose:

 SLA_TRSYLV_L2_LOCAL_COPY solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV and SGER by Fortran intrinsic,
Use local copies of A, B, C, D, and X (M, N <=128) .

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 181 of file sla_trsylv_l2_opt_local_copy.f90.

subroutine sla_trsylv_l2_local_copy_128 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

Level-2 Bartels-Stewart Algorithm for the generalized Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV_L2_LOCAL_COPY_128 solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  128 >= M >= 0.

          N is INTEGER
          The order of the matrices B and D.  128 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV by Fortran intrinsic,
Use local copies of A, B, C, D, and X.
Align local copies and fix problem size to <= 64

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 183 of file sla_trsylv_l2_opt_local_copy_128.f90.

subroutine sla_trsylv_l2_local_copy_32 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

Level-2 Bartels-Stewart Algorithm for the generalized Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV_L2_LOCAL_COPY_32 solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  32 >= M >= 0.

          N is INTEGER
          The order of the matrices B and D.  32 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV by Fortran intrinsic,
Use local copies of A, B, C, D, and X.
Align local copies and fix problem size to <= 64

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 183 of file sla_trsylv_l2_opt_local_copy_32.f90.

subroutine sla_trsylv_l2_local_copy_64 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

Level-2 Bartels-Stewart Algorithm for the generalized Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV_L2_LOCAL_COPY_64 solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  64 >= M >= 0.

          N is INTEGER
          The order of the matrices B and D.  64 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV by Fortran intrinsic,
Use local copies of A, B, C, D, and X.
Align local copies and fix problem size to <= 64

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 182 of file sla_trsylv_l2_opt_local_copy_64.f90.

subroutine sla_trsylv_l2_local_copy_96 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

Level-2 Bartels-Stewart Algorithm for the generalized Sylvester equation (Optimized)

Purpose:

 SLA_TRSYLV_L2_LOCAL_COPY_96 solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  96 >= M >= 0.

          N is INTEGER
          The order of the matrices B and D.  96 >= N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV by Fortran intrinsic,
Use local copies of A, B, C, D, and X.
Align local copies and fix problem size to <= 64

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 181 of file sla_trsylv_l2_opt_local_copy_96.f90.

subroutine sla_trsylv_l2_reorder (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-2 Bartels-Stewart Algorithm for the Sylvester equation (optimized)

Purpose:

 SLA_TRSYLV_L2_REORDER solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Attention

The algorithm is implemented using BLAS level 2 operations without further optimizations. For a faster implementation see SLA_TRSYLV_L2 .

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Replaced level-3 by level-2 and level-1 calls,
Reorder the solution order to column first style,
Replaced SAXPY operation by Fortran intrinsic
Replaced all BLAS calls except of STRMV and SGER by Fortran intrinsic,

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 182 of file sla_trsylv_l2_opt_reorder.f90.

subroutine sla_trsylv_l2_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-2 Bartels-Stewart Algorithm for the Sylvester equation (unoptimized)

Purpose:

 SLA_TRSYLV_L2_UNOPT solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Attention

The algorithm is implemented using BLAS level 2 operations without further optimizations. For a faster implementation see SLA_TRSYLV_L2 .

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A and C:
          == 'N':  op1(A) = A, op1(C) = C (No transpose for A and C)
          == 'T':  op1(A) = A**T, op1(C) = C **T (Transpose A and C)

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B and D:
          == 'N':  op2(B) = B, op2(D) = D (No transpose for B and D)
          == 'T':  op2(B) = B**T, op2(D) = D **T (Transpose B and D)

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Optimizations:

Nothing.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 180 of file sla_trsylv_l2_unopt.f90.

subroutine sla_trsylv_l3 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-3 Bartels-Stewart Algorithm for the generalized Sylvester equation.

Purpose:

 SLA_TRSYLV_L3 solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Attention

The algorithm is implemented using BLAS level 3 operations.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 172 of file sla_trsylv_l3_opt.f90.

subroutine sla_trsylv_l3_2s (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-3 Bartels-Stewart Algorithm for the generalized Sylvester equation.

Purpose:

 SLA_TRSYLV_L3_2S solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Remarks

The algorithm is implemented using BLAS level 3 operations and a two level inner solver consisting of a massive parallel DAG scheduled solver from SLA_TRSYLV_DAG and optimized kernel solvers.

Attention

Due to the parallel nature of the inner solvers the scaling is turned off and SCALE is set to ONE.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 174 of file sla_trsylv_l3_2stage.f90.

subroutine sla_trsylv_l3_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-3 Bartels-Stewart Algorithm for the generalized Sylvester equation (unoptimized)

Purpose:

 SLA_TRSYLV_L3_UNOPT solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Attention

The algorithm is implemented using BLAS level 3 operations.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension LWORK
          Workspace for the algorithm.
          The workspace needs to queried before the running the computation.
          The query is performed by calling the subroutine with INFO == -1 on input.
          The required workspace is then returned in INFO.

INFO

          INFO is INTEGER

          On input:
            == -1 : Perform a workspace query
            <> -1: normal operation

          On exit, workspace query:
            < 0 :  if INFO = -i, the i-th argument had an illegal value
            >= 0:  The value of INFO is the required number of elements in the workspace.

          On exit, normal operation:
            == 0:  successful exit
            < 0:  if INFO = -i, the i-th argument had an illegal value
            > 0:  The equation is not solved correctly. One of the arising inner
                  system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 172 of file sla_trsylv_l3_unopt.f90.

recursive subroutine sla_trsylv_recursive (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Level-3 Recursive Blocked Solver for the Sylvester equation.

Purpose:

 SLA_TRSYLV_RECURSIVE solves a Sylvester equation of the following forms

    op1(A) * X  +  X * op2(B) = SCALE * Y                              (1)

 or

    op1(A) * X  -  X * op2(B) = SCALE * Y                              (2)

 where A is a M-by-M quasi upper triangular matrix, B is a N-by-N quasi upper triangular
 matrix. The right hand side Y and the solution X are M-by-N matrices.  Typically the matrices
 A and B are generated via SGEES form LAPACK.

Attention

The algorithm is implemented using BLAS level 3 operations and recursive blocking.

Parameters

TRANSA

          TRANSA is CHARACTER(1)
          Specifies the form of the system of equations with respect to A:
          == 'N':  op1(A) = A
          == 'T':  op1(A) = A**T

TRANSB

          TRANSB is CHARACTER(1)
          Specifies the form of the system of equations with respect to B:
          == 'N':  op2(B) = B,
          == 'T':  op2(B) = B**T

SGN

          SGN is REAL, allowed values: +/-1
          Specifies the sign between the two parts of the Sylvester equation.
          = 1 :  Solve Equation (1)
          == -1:  Solve Equation (2)

          M is INTEGER
          The order of the matrices A and C.  M >= 0.

          N is INTEGER
          The order of the matrices B and D.  N >= 0.

          A is REAL array, dimension (LDA,M)
          The matrix A must be (quasi-) upper triangular.

LDA

          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).

          B is REAL array, dimension (LDB,N)
          The matrix B must be (quasi-) upper triangular.

LDB

          LDB is INTEGER
          The leading dimension of the array B.  LDB >= max(1,N).

          X is REAL array, dimension (LDX,N)
          On input, the matrix X contains the right hand side Y.
          On output, the matrix X contains the solution of Equation (1) or (2)
          as selected by TRANSA, TRANSB, and SGN.
          Right hand side Y and the solution X are M-by-N matrices.

LDX

          LDX is INTEGER
          The leading dimension of the array X.  LDB >= max(1,M).

SCALE

          SCALE is REAL
          SCALE is a scaling factor to prevent the overflow in the result.
          If INFO == 0 then SCALE is 1.0 otherwise if one of the inner systems
          could not be solved correctly, 0 < SCALE <= 1 holds true.

WORK

          WORK is REAL array, dimension 1
          Workspace for the algorithm

INFO

          INFO is INTEGER
          == 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value
          > 0:  The equation is not solved correctly. One of the arising inner
                system got singular.

Author

Martin Koehler, MPI Magdeburg

Date

January 2024

Definition at line 160 of file sla_trsylv_recursive.f90.

Author

Generated automatically by Doxygen for MEPACK from the source code.

Info

Fri Feb 2 2024 00:00:00 Version 1.1.1 MEPACK

sgltrsylv - Man Page

Name

Synopsis

Functions

Detailed Description

Function Documentation

subroutine sla_trsylv2_dag (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, m) a, integer lda, real, dimension(ldb, n) b, integer ldb, real, dimension(ldx, n) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv2_kernel_44nn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_kernel_44nt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_kernel_44tn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_kernel_44tt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_l2 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l2_local_copy (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l2_local_copy_128 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv2_l2_local_copy_32 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv2_l2_local_copy_64 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv2_l2_local_copy_96 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv2_l2_reorder (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l2_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l3 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l3_2s (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l3_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_dag (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, n) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv_kernel_44nn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_kernel_44nt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_kernel_44tn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_kernel_44tt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_l2 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l2_local_copy (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l2_local_copy_128 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv_l2_local_copy_32 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv_l2_local_copy_64 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv_l2_local_copy_96 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv_l2_reorder (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l2_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l3 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l3_2s (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l3_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

Author

Info

Name

Synopsis

Functions

Detailed Description

Function Documentation

subroutine sla_trsylv2_dag (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, m) a, integer lda, real, dimension(ldb, n) b, integer ldb, real, dimension(ldx, n) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv2_kernel_44nn (real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_kernel_44nt (real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_kernel_44tn (real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_kernel_44tt (real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_l2 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv2_l2_local_copy_128 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv2_l2_local_copy_32 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv2_l2_local_copy_64 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv2_l2_local_copy_96 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv2_l2_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv2_l3 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv2_l3_2s (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv2_l3_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv_dag (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, n) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv_kernel_44nn (real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_kernel_44nt (real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_kernel_44tn (real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_kernel_44tt (real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_l2 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv_l2_local_copy_128 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv_l2_local_copy_32 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv_l2_local_copy_64 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv_l2_local_copy_96 ( transa, transb, sgn, m, n, a, lda, b, ldb, x, ldx, scale, work, info)

subroutine sla_trsylv_l2_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv_l3 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv_l3_2s (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv_l3_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, *) a, integer lda, real, dimension(ldb, *) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, real, dimension(*) work, integer info)

Author

Info

subroutine sla_trsylv2_kernel_44nn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_kernel_44nt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_kernel_44tn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_kernel_44tt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv2_l2 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l2_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l3 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l3_2s (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv2_l3_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_dag (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, n) x, integer ldx, real scale, real, dimension(*) work, integer info)

subroutine sla_trsylv_kernel_44nn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_kernel_44nt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_kernel_44tn (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_kernel_44tt (real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, *) x, integer ldx, real scale, integer info)

subroutine sla_trsylv_l2 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l2_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l3 (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l3_2s (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)

subroutine sla_trsylv_l3_unopt (character, dimension(1) transa, character, dimension(1) transb, real sgn, integer m, integer n, real, dimension(lda, ) a, integer lda, real, dimension(ldb, ) b, integer ldb, real, dimension(ldx, ) x, integer ldx, real scale, real, dimension() work, integer info)