chol - Cholesky factorization - MATLAB (original) (raw)

Syntax

Description

[R](#mw%5F3e727a04-bd5b-42cd-b5bb-a3f7b411912e) = chol([A](#mw%5F91dc7b0f-768f-4a6b-ab21-391fd8ee5851)) factorizes symmetric positive definite matrix A into an upper triangularR that satisfies A = R'*R. If A is nonsymmetric, then chol treats the matrix as symmetric and uses only the diagonal and upper triangle of A.

example

[R](#mw%5F3e727a04-bd5b-42cd-b5bb-a3f7b411912e) = chol([A](#mw%5F91dc7b0f-768f-4a6b-ab21-391fd8ee5851),[triangle](#mw%5Fffc71c97-36a0-4335-8e05-82ecff3e6312)) specifies which triangular factor of A to use in computing the factorization. For example, if triangle is 'lower', then chol uses only the diagonal and lower triangular portion ofA to produce a lower triangular matrix R that satisfies A = R*R'. The default value of triangle is'upper'.

example

[[R](#mw%5F3e727a04-bd5b-42cd-b5bb-a3f7b411912e),[flag](#mw%5Fda905fd8-e0ca-4b83-884e-1216badb7607)] = chol(___) also returns the output flag indicating whether A issymmetric positive definite. You can use any of the input argument combinations in previous syntaxes. When you specify the flag output, chol does not generate an error if the input matrix is not symmetric positive definite.

example

[[R](#mw%5F3e727a04-bd5b-42cd-b5bb-a3f7b411912e),[flag](#mw%5Fda905fd8-e0ca-4b83-884e-1216badb7607),[P](#mw%5F9048ccd2-eccf-45df-a23b-5fa1ce84d320)] = chol([S](#mw%5Fbd8d494d-1956-4117-a10b-1adaa40cfc91)) additionally returns a permutation matrix P, which is a preordering of sparse matrix S obtained by amd. If flag = 0, then S is symmetric positive definite andR is an upper triangular matrix satisfying R'*R = P'*S*P.

example

[[R](#mw%5F3e727a04-bd5b-42cd-b5bb-a3f7b411912e),[flag](#mw%5Fda905fd8-e0ca-4b83-884e-1216badb7607),[P](#mw%5F9048ccd2-eccf-45df-a23b-5fa1ce84d320)] = chol(___,[outputForm](#mw%5F2bc2adb8-4c9b-4515-9939-7769773a490d)) specifies whether to return the permutation information P as a matrix or vector, using any of the input argument combinations in previous syntaxes. This option is only available for sparse matrix inputs. For example, if outputForm is'vector' and flag = 0, then S(p,p) = R'*R. The default value of outputForm is'matrix' such that R'*R = P'*S*P.

example

Examples

collapse all

Solve Linear System with Symmetric Positive Definite Matrix

Use chol to factorize a symmetric coefficient matrix, and then solve a linear system using the Cholesky factor.

Create a symmetric matrix with positive values on the diagonal.

A = [1 0 1; 0 2 0; 1 0 3]

A = 3×3

 1     0     1
 0     2     0
 1     0     3

Calculate the Cholesky factor of the matrix.

R = 3×3

1.0000         0    1.0000
     0    1.4142         0
     0         0    1.4142

Create a vector for the right-hand side of the equation Ax=b.

Since A=RTR with the Cholesky decomposition, the linear equation becomes RTR x=b. Solve for x using the backslash operator.

x = 3×1

1.0000
1.0000
1.0000

Cholesky Factorization of Matrix

Calculate the upper and lower Cholesky factorizations of a matrix and verify the results.

Create a 6-by-6 symmetric positive definite test matrix using the gallery function.

Calculate the Cholesky factor using the upper triangle of A.

R = 6×6

1.0000    0.5000    0.3333    0.2500    0.2000    0.1667
     0    0.8660    0.5774    0.4330    0.3464    0.2887
     0         0    0.7454    0.5590    0.4472    0.3727
     0         0         0    0.6614    0.5292    0.4410
     0         0         0         0    0.6000    0.5000
     0         0         0         0         0    0.5528

Verify that the upper triangular factor satisfies R'*R - A = 0, within round-off error.

Now, specify the 'lower' option to calculate the Cholesky factor using the lower triangle of A.

L = 6×6

1.0000         0         0         0         0         0
0.5000    0.8660         0         0         0         0
0.3333    0.5774    0.7454         0         0         0
0.2500    0.4330    0.5590    0.6614         0         0
0.2000    0.3464    0.4472    0.5292    0.6000         0
0.1667    0.2887    0.3727    0.4410    0.5000    0.5528

Verify that the lower triangular factor satisfies L*L' - A = 0, within round-off error.

Suppress Errors for Nonsymmetric Positive Definite Matrices

Use chol with two outputs to suppress errors when the input matrix is not symmetric positive definite.

Create a 5-by-5 matrix of binomial coefficients. This matrix is symmetric positive definite, so subtract 1 from the last element to ensure it is no longer positive definite.

A = pascal(5); A(end) = A(end) - 1

A = 5×5

 1     1     1     1     1
 1     2     3     4     5
 1     3     6    10    15
 1     4    10    20    35
 1     5    15    35    69

Calculate the Cholesky factor for A. Specify two outputs to avoid generating an error if A is not symmetric positive definite.

R = 4×4

 1     1     1     1
 0     1     2     3
 0     0     1     3
 0     0     0     1

Since flag is nonzero, it gives the pivot index where the factorization fails. chol is able to calculate q = flag-1 = 4 rows and columns correctly before failing when it encounters the part of the matrix that changed.

Verify that R'*R returns four rows and columns that agree with A(1:q,1:q).

ans = 4×4

 1     1     1     1
 1     2     3     4
 1     3     6    10
 1     4    10    20

ans = 4×4

 1     1     1     1
 1     2     3     4
 1     3     6    10
 1     4    10    20

Cholesky Factor of Sparse Matrix

Calculate the Cholesky factor of a sparse matrix, and use the permutation output to create a Cholesky factor with fewer nonzeros.

Create a sparse positive definite matrix based on the west0479 matrix.

load west0479 A = west0479; S = A'*A;

Calculate the Cholesky factor of the matrix two different ways. First specify two outputs, and then specify three outputs to enable row and column reordering.

[R,flag] = chol(S); [RP,flagP,P] = chol(S);

For each calculation, check that flag = 0 to confirm the calculation is successful.

if ~flag && ~flagP disp('Factorizations successful.') else disp('Factorizations failed.') end

Factorizations successful.

Compare the number of nonzeros in chol(S) vs. the reordered matrix chol(P'*S*P). Best practice is to use the three output syntax of chol with sparse matrices, since reordering the rows and columns can greatly reduce the number of nonzeros in the Cholesky factor.

subplot(1,2,1) spy(R) title('Nonzeros in chol(S)') subplot(1,2,2) spy(RP) title('Nonzeros in chol(P''SP)')

Figure contains 2 axes objects. Axes object 1 with title Nonzeros in chol(S), xlabel nz = 59551 contains a line object which displays its values using only markers. Axes object 2 with title Nonzeros in chol(P'*S*P), xlabel nz = 7637 contains a line object which displays its values using only markers.

Reorder Sparse Matrix with Permutation Vector

Use the 'vector' option of chol to return the permutation information as a vector rather than a matrix.

Create a sparse finite element matrix.

S = gallery('wathen',10,10); spy(S)

Figure contains an axes object. The axes object with xlabel nz = 4861 contains a line object which displays its values using only markers.

Calculate the Cholesky factor for the matrix, and specify the 'vector' option to return a permutation vector p.

[R,flag,p] = chol(S,'vector');

Verify that flag = 0, indicating the calculation is successful.

if ~flag disp('Factorization successful.') else disp('Factorization failed.') end

Factorization successful.

Verify that S(p,p) = R'*R, within round-off error.

norm(S(p,p) - R'*R,'fro')

Input Arguments

collapse all

A — Input matrix

matrix

Input matrix. Argument A can use full or sparse storage, but must be square and symmetric positive definite.

chol assumes that A is symmetric for real matrices or Hermitian for complex matrices. chol uses only the upper or lower triangle of A to perform its computations, depending on the value of triangle.

Data Types: single | double
Complex Number Support: Yes

S — Sparse input matrix

sparse matrix

Sparse input matrix. S must be square and symmetric positive definite.

chol assumes that S is symmetric for real matrices or Hermitian for complex matrices. chol uses only the upper or lower triangle of S to perform its computations, depending on the value of triangle.

Data Types: double
Complex Number Support: Yes

triangle — Triangular factor of input matrix

'upper' (default) | 'lower'

Triangular factor of input matrix, specified as 'upper' or'lower'. Use this option to specify that chol should use the upper or lower triangle of the input matrix to compute the factorization.chol assumes that the input matrix is symmetric for real matrices or Hermitian for complex matrices. chol uses only the upper or lower triangle to perform its computations.

Using the 'lower' option is equivalent to callingchol with the 'upper' option and the transpose of the input matrix, and then transposing the output R.

Example: R = chol(A,'lower')

outputForm — Shape of permutation output

'matrix' (default) | 'vector'

Shape of permutation output, specified as 'matrix' or'vector'. This flag controls whether the permutation outputP is returned as a permutation matrix or permutation vector.

The Cholesky factor of P'*S*P (if P is a matrix) or S(p,p) (if p is a vector) tends to be sparser than the Cholesky factor of S.

Example: [R,flag,p] = chol(S,'vector')

Output Arguments

collapse all

R — Cholesky factor

matrix

Cholesky factor, returned as a matrix.

flag — Symmetric positive definite flag

scalar

Symmetric positive definite flag, returned as a scalar.

P — Permutation for sparse matrices

matrix | vector

Permutation for sparse matrices, returned as a matrix or vector depending on the value of outputForm. See outputForm for a description of the identities that this output satisfies.

This permutation matrix is based on the approximate minimum degree ordering computed by amd. However, this preordering can differ from the one obtained directly by amd since chol slightly changes the ordering for increased performance.

More About

collapse all

Symmetric Positive Definite Matrix

A symmetric positive definite matrix is a symmetric matrix with all positive eigenvalues.

For any real invertible matrix A, you can construct a symmetric positive definite matrix with the product B = A'*A. The Cholesky factorization reverses this formula by saying that any symmetric positive definite matrixB can be factored into the product R'*R.

A symmetric positive semi-definite matrix is defined in a similar manner, except that the eigenvalues must all be positive or zero.

The line between positive definite and positive semi-definite matrices is blurred in the context of numeric computation. It is rare for eigenvalues to be exactly equal to zero, but they can be numerically zero (on the order of machine precision). For this reason,chol might be able to factorize one positive semi-definite matrix, but could fail with another matrix that has very similar eigenvalues.

Tips

References

[1] Anderson, E., ed. LAPACK Users’ Guide. 3rd ed. Software, Environments, Tools. Philadelphia: Society for Industrial and Applied Mathematics, 1999. https://doi.org/10.1137/1.9780898719604.

[2] Chen, Yanqing, Timothy A. Davis, William W. Hager, and Sivasankaran Rajamanickam. “Algorithm 887: CHOLMOD, Supernodal Sparse Cholesky Factorization and Update/Downdate.” ACM Transactions on Mathematical Software 35, no. 3 (October 2008): 1–14. https://doi.org/10.1145/1391989.1391995.

Extended Capabilities

C/C++ Code Generation

Generate C and C++ code using MATLAB® Coder™.

Usage notes and limitations:

GPU Code Generation

Generate CUDA® code for NVIDIA® GPUs using GPU Coder™.

Usage notes and limitations:

Thread-Based Environment

Run code in the background using MATLAB® backgroundPool or accelerate code with Parallel Computing Toolbox™ ThreadPool.

This function fully supports thread-based environments. For more information, see Run MATLAB Functions in Thread-Based Environment.

GPU Arrays

Accelerate code by running on a graphics processing unit (GPU) using Parallel Computing Toolbox™.

The chol function supports GPU array input with these usage notes and limitations:

For more information, see Run MATLAB Functions on a GPU (Parallel Computing Toolbox).

Distributed Arrays

Partition large arrays across the combined memory of your cluster using Parallel Computing Toolbox™.

Usage notes and limitations:

For more information, see Run MATLAB Functions with Distributed Arrays (Parallel Computing Toolbox).

Version History

Introduced before R2006a