linSolPosdefComplex

Solves a complex Hermitian positive definite system of linear equations \(Ax = b\). Using optional arguments, any of several related computations can be performed. These extra tasks include computing the Cholesky factor, L, of A such that \(A = LL^H\) or computing the solution to \(Ax = b\) given the Cholesky factor, L.

Synopsis

linSolPosdefComplex (a,  b)

Required Arguments

complex a[[]] (Input)

Array of size n × n containing the matrix.

complex b[] (Input)

Array of size n containing the right-hand side.

Return Value

The solution x of the Hermitian positive definite linear system Ax = b. If no solution was computed, then None is returned.

Optional Arguments

factor, complex pFactor (Output)

An array of size n × n containing the \(LL^H\) factorization of A. The lower‑triangular part of this array contains L, and the upper-triangular part contains \(L^H\).

condition (Output)

A scalar containing an estimate of the \(L_1\) norm condition number of the matrix A. Do not use this option with solveOnly.

factorOnly

Compute the Cholesky factorization \(LL^H\) of A. If factorOnly is used, factor is required. The argument b is then ignored, and the returned value of linSolPosdefComplex is None.

solveOnly

Solve \(Ax = b\) given the \(LL^H\) factorization previously computed by linSolPosdefComplex. By default, the solution to \(Ax = b\) is pointed to by linSolPosdefComplex. If solveOnly is used, argument factor is required and argument a is ignored.

Description

The function linSolPosdefComplex solves a system of linear algebraic equations having a Hermitian positive definite coefficient matrix A. The function first computes the \(LL^H\) factorization of A. The solution of the linear system is then found by solving the two simpler systems, \(y = L^{-1}b\) and \(x = L^{-H}y\). When the solution to the linear system is required, an estimate of the \(L_1\) condition number of A is computed using the algorithm in Dongarra et al. (1979). If the estimated condition number is greater than 1∕ɛ (where ɛ is the machine precision), a warning message is issued. This indicates that very small changes in A may produce large changes in the solution x. The function linSolPosdefComplex fails if L, the lower-triangular matrix in the factorization, has a zero diagonal element.

Examples

Example 1

A system of five linear equations with a Hermitian positive definite coefficient matrix is solved in this example. The equations are as follows:

\[2x_1 +(−1 + i)x_2 = 1 +5i\]

\[(−1 − i)x_1 +4x_2 + (1 + 2i)x_3 = 12 − 6i\]

\[(1 − 2i)x_2 +10x_3 + 4ix_4 = 1 − 16i\]

\[−4ix_3 + 6x_4 + (1 + i)x_5 = −3 − 3i\]

\[(1 − i)x_4 + 9x_5 = 25 + 16i\]

from numpy import *
from pyimsl.math.linSolPosdefComplex import linSolPosdefComplex
from pyimsl.math.writeMatrixComplex import writeMatrixComplex

factor = []
a = array([[2.0 + 0.0j, -1.0 + 1.0j, 0.0 + 0.0j, 0.0 + 0.0j, 0.0 + 0.0j],
           [-1.0 - 1.0j, 4.0 + 0.0j, 1.0 + 2.0j, 0.0 + 0.0j, 0.0 + 0.0j],
           [0.0 + 0.0j, 1.0 - 2.0j, 10.0 + 0.0j, 0.0 + 4.0j, 0.0 + 0.0j],
           [0.0 + 0.0j, 0.0 + 0.0j, 0.0 - 4.0j, 6.0 + 0.0j, 1.0 + 1.0j],
           [0.0 + 0.0j, 0.0 + 0.0j, 0.0 + 0.0j, 1.0 - 1.0j, 9.0 + 0.0j]])
b = array([1.0 + 5.0j, 12.0 - 6.0j, 1.0 - 16.0j, -3.0 - 3.0j, 25.0 + 16.0j])

# Solve Ax=b for x
x = linSolPosdefComplex(a, b)
writeMatrixComplex("Solution, x, of Ax = b", x, column=True)

Output

 
   Solution, x, of Ax = b
1  (          2,          1)
2  (          3,          0)
3  (         -1,         -1)
4  (          0,         -2)
5  (          3,          2)

Example 2

This example solves the same system of five linear equations as in the first example. This time, the \(LL^H\) factorization of A and the solution x are returned.

from numpy import *
from pyimsl.math.linSolPosdefComplex import linSolPosdefComplex
from pyimsl.math.writeMatrixComplex import writeMatrixComplex

factor = []
a = array([[2.0 + 0.0j, -1.0 + 1.0j, 0.0 + 0.0j, 0.0 + 0.0j, 0.0 + 0.0j],
           [-1.0 - 1.0j, 4.0 + 0.0j, 1.0 + 2.0j, 0.0 + 0.0j, 0.0 + 0.0j],
           [0.0 + 0.0j, 1.0 - 2.0j, 10.0 + 0.0j, 0.0 + 4.0j, 0.0 + 0.0j],
           [0.0 + 0.0j, 0.0 + 0.0j, 0.0 - 4.0j, 6.0 + 0.0j, 1.0 + 1.0j],
           [0.0 + 0.0j, 0.0 + 0.0j, 0.0 + 0.0j, 1.0 - 1.0j, 9.0 + 0.0j]])
b = array([1.0 + 5.0j, 12.0 + -6.0j, 1.0 - 16.0j, -3.0 + -3.0j, 25.0 + 16.0j])

# Solve Ax=b for x
x = linSolPosdefComplex(a, b, factor=factor)
writeMatrixComplex("Solution, x, of Ax = b", x, column=True)
writeMatrixComplex("Cholesky factor L, and ctrans(L), of A", factor)

Output

 
   Solution, x, of Ax = b
1  (          2,          1)
2  (          3,          0)
3  (         -1,         -1)
4  (          0,         -2)
5  (          3,          2)
 
        Cholesky factor L, and ctrans(L), of A
                           1                          2
1  (      1.414,      0.000)  (     -0.707,      0.707)
2  (     -0.707,     -0.707)  (      1.732,      0.000)
3  (      0.000,      0.000)  (      0.577,     -1.155)
4  (      0.000,      0.000)  (      0.000,      0.000)
5  (      0.000,      0.000)  (      0.000,      0.000)
 
                           3                          4
1  (      0.000,     -0.000)  (      0.000,     -0.000)
2  (      0.577,      1.155)  (      0.000,     -0.000)
3  (      2.887,      0.000)  (      0.000,      1.386)
4  (      0.000,     -1.386)  (      2.020,      0.000)
5  (      0.000,      0.000)  (      0.495,     -0.495)
 
                           5
1  (      0.000,     -0.000)
2  (      0.000,     -0.000)
3  (      0.000,     -0.000)
4  (      0.495,      0.495)
5  (      2.917,      0.000)

Warning Errors

IMSL_HERMITIAN_DIAG_REAL_1	The diagonal of a Hermitian matrix must be real. Its imaginary part is set to zero.
IMSL_HERMITIAN_DIAG_REAL_2	The diagonal of a Hermitian matrix must be real. The imaginary part will be used as zero in the algorithm.
IMSL_ILL_CONDITIONED	The input matrix is too ill-conditioned. An estimate of the reciprocal of its \(L_1\) condition number is “rcond” = #. The solution might not be accurate.

Fatal Errors

IMSL_NONPOSITIVE_MATRIX	The leading # by # minor matrix of the input matrix is not positive definite.
IMSL_HERMITIAN_DIAG_REAL	During the factorization the matrix has a large imaginary component on the diagonal. Thus, it cannot be positive definite.
IMSL_SINGULAR_TRI_MATRIX	The triangular matrix is singular. The index of the first zero diagonal term is #.