# decomposition

Matrix decomposition for solving linear systems

## Description

`decomposition` creates reusable matrix decompositions (LU, LDL, Cholesky, QR, and more) that enable you to solve linear systems (Ax = b or xA = b) more efficiently. For example, after computing `dA = decomposition(A)` the call `dA\b` returns the same vector as `A\b`, but is typically much faster. `decomposition` objects are well-suited to solving problems that require repeated solutions, since the decomposition of the coefficient matrix does not need to be performed multiple times.

You can use a `decomposition` object `dA` with many of the same operators you might use on the original coefficient matrix `A`:

• Complex conjugate transpose `dA'`

• Negation `-dA`

• Multiply or divide by a scalar using `c*dA` or `dA/c`.

• Solve a linear system Ax = b using ```x = dA\b```.

• Solve a linear system xA = b using ```x = b/dA```.

## Creation

### Syntax

``dA = decomposition(A)``
``dA = decomposition(A,type)``
``dA = decomposition(A,type,triangularFlag)``
``dA = decomposition(___,Name,Value)``

### Description

example

````dA = decomposition(A)` returns a decomposition of matrix `A` that you can use to solve linear systems more efficiently. The decomposition type is automatically chosen based on the properties of the input matrix.```

example

````dA = decomposition(A,type)` specifies the type of decomposition to perform. `type` can be `'qr'`, `'cod'`, `'lu'`, `'ldl'`, `'chol'`, `'triangular'`, `'permutedTriangular'`, `'banded'`, `'hessenberg'`, or `'diagonal'`.```

example

````dA = decomposition(A,type,triangularFlag)` specifies that only the upper or lower triangular portion of `A` is to be used in the decomposition. `triangularFlag` can be `'upper'` or `'lower'`. With this syntax, the decomposition type must be `'ldl'`, `'chol'`, or `'triangular'`.```

example

````dA = decomposition(___,Name,Value)` specifies additional options using one or more `Name,Value` pair arguments using any of the previous syntaxes. For example, ```dA = decomposition(A,'CheckCondition',false)``` specifies not to throw a warning based on the condition of `A` while solving `dA\b`.```

### Input Arguments

expand all

Coefficient matrix. The coefficient matrix appears in the system of linear equations on the left as Ax = b or on the right as xA = b.

Data Types: `single` | `double`
Complex Number Support: Yes

Decomposition type, specified as one of the options in these tables.

These options work for any coefficient matrix.

Value

Matrix Decomposition of `A`

Notes

`'auto'` (default)

N/A

Automatic selection of the matrix decomposition type based on the properties of the coefficient matrix. For information about how the decomposition is chosen, see the Algorithms section of `mldivide`.

`'qr'`

`$AP=QR$`

`Q` is unitary, `R` is upper triangular, and `P` is a permutation matrix.

QR decompositions give a least-squares solution.

If `type` is `'qr'`, then you cannot solve `A'\B` or `B/A`. Instead, use `'cod'` for problems with those forms.

`'cod'`

`$A=QR{Z}^{*}$`

`R` is upper triangular, and both `Q` and `Z` have orthonormal columns.

Complete orthogonal decompositions give a minimum norm least-squares solution.

For square coefficient matrices, you also can use these options.

Value

Matrix Decomposition of `A`

Notes

`'lu'`

Dense matrices:

`$PA=LU$`

`L` is lower triangular, `U` is upper triangular, and `P` is a permutation matrix.

Sparse matrices:

`$P\left(R\A\right)Q=LU$`

`P` and `Q` are permutation matrices and `R` is a diagonal scaling matrix.

`'ldl'`

Dense matrices:

`${P}^{*}AP=LD{L}^{*}$`

`L` is a lower triangular matrix with 1s on the diagonal, `D` is a diagonal matrix, and `P` is a permutation matrix.

Sparse matrices:

`${P}^{*}SASP=LD{L}^{*}$`

`S` is a scaling matrix.

`A` must be symmetric.

`'chol'`

Dense matrices:

`$A=L{L}^{*}$`

`L` is lower triangular.

Sparse matrices:

`$A=PL{L}^{*}{P}^{*}$`

`P` is a permutation matrix.

`A` must be symmetric positive definite.

`'triangular'`

`$A=T$`

`T` is triangular.

`A` must be triangular.

`'permutedTriangular'`

`$A=PT{Q}^{*}$`

`T` is triangular, and both `P` and `Q` are permutation matrices.

`A` must be a permutation of a triangular matrix.

`'banded'`

`$A={P}^{*}LU$`

`P` is a permutation matrix and both `L` and `U` are banded.

Most effective for matrices with a low bandwidth. See `bandwidth` for more information.

`'hessenberg'`

`$A={P}^{*}LU$`

`P` is a permutation matrix, `L` is banded, and `U` is upper triangular.

`A` must have zeros below the first subdiagonal.

`'diagonal'`

`$A=D$`

`D` is diagonal.

`A` must be diagonal.

Flag to use only upper or lower triangular portion of coefficient matrix, specified as either `'upper'` or `'lower'`. This option supports the `'triangular'`, `'chol'`, and `'ldl'` decomposition types.

• `'triangular'` — If both an upper and lower triangular matrix are stored in the same matrix, then use `triangularFlag` to specify only one of the triangular matrices.

• `'chol'` and `'ldl'` — Use `triangularFlag` to avoid symmetrizing a nearly symmetric coefficient matrix.

#### Name-Value Pair Arguments

Specify optional comma-separated pairs of `Name,Value` arguments. `Name` is the argument name and `Value` is the corresponding value. `Name` must appear inside quotes. You can specify several name and value pair arguments in any order as `Name1,Value1,...,NameN,ValueN`.

Example: `dA = decomposition(A,'qr','CheckCondition',false)` performs a QR decomposition of `A` and turns off warnings about the condition of the coefficient matrix when it is used to solve a linear system.

#### General Parameters

expand all

Toggle to check condition of coefficient matrix, specified as the comma-separated pair consisting of `'CheckCondition'` and either logical `1` (`true`) or logical `0` (`false`). If `CheckCondition` is `true` and the coefficient matrix is badly conditioned or of low rank, then solving linear systems using `mldivide (\)` or `mrdivide (/)` produces warnings.

Data Types: `logical`

Rank tolerance, specified as a nonnegative scalar. Specifying the tolerance can help prevent the solution from being susceptible to random noise in the coefficient matrix.

`decomposition` computes the rank of `A` as the number of diagonal elements in the `R` matrix of the QR decomposition `[Q,R,p] = qr(A,0)` with absolute value larger than `tol`. If the rank of `A` is `k`, then a low-rank approximation of `A` is formed by multiplying the first `k` columns of `Q` by the first `k` rows of `R`. Changing the tolerance affects this low-rank approximation of `A`.

### Note

This option only applies when `'Type'` is `'qr'` or `'cod'`, or when `'Type'` is `'auto'` and `A` is rectangular. Otherwise, this option is ignored.

#### Sparse Matrix Parameters

expand all

Band density threshold, specified as a scalar value in the range `[0 1]`. The value of `'BandDensity'` determines how dense a sparse, banded coefficient matrix must be for the banded solver to be used by `mldivide (\)` or ```mrdivide (/)``` when solving a system of equations. If the band density of the coefficient matrix is larger than the specified band density, then the banded solver is used.

The band density is defined as: (# nonzeros in the band) / (# elements in the band). A value of `1.0` indicates to never use the banded solver.

Pivot tolerance for LDL factorization, specified as a scalar value in the interval `[0 0.5]`. Using smaller values of the pivot tolerance can give faster factorization times and fewer entries, but also can result in a less stable factorization.

This pivot tolerance is the same that `ldl` uses for real sparse matrices.

Pivot tolerance for LU factorization, specified as a scalar or vector. Specify a scalar value to change the first element in the tolerance vector, or specify a two-element vector to change both values. Smaller pivot tolerances tend to lead to sparser LU factors, but the solution can become inaccurate. Larger values can lead to a more accurate solution, but not always, and usually increase the total work and memory usage.

This pivot tolerance is the same that `lu` uses for sparse matrices.

## Properties

expand all

Size of coefficient matrix, returned as a two-element row vector.

Data Types: `double`

Decomposition type, returned as `'qr'`, `'cod'`, `'lu'`, `'ldl'`, `'chol'`, `'triangular'`, `'permutedTriangular'`, `'banded'`, `'hessenberg'`, or `'diagonal'`.

Data Types: `char`

Toggle to check condition of coefficient matrix, specified as either logical `1` (`true`) or logical `0` (`false`). If `CheckCondition` is `true` and the coefficient matrix is badly conditioned or of low rank, then solving linear systems using `mldivide (\)` or ```mrdivide (/)``` produces warnings.

Data Types: `logical`

Data type of coefficient matrix, returned as either `'double'` or `'single'`.

Data Types: `char`

Indicator that coefficient matrix is complex conjugate transposed, returned as either logical `1` (`true`) or logical `0` (`false`). This indicator is `false` by default for any `decomposition` object constructed from the coefficient matrix. However, the value is `true` if you use the `ctranspose` operator on a `decomposition` object in an expression, such as `dA'\b`. In that case, `dA'` is the same `decomposition` object as `dA`, but with a value of `true` for `IsConjugateTransposed`.

Data Types: `logical`

Indicator that coefficient matrix is real, returned as either logical `1` (`true`) or logical `0` (`false`). A value of `false` indicates that the coefficient matrix contains complex numbers.

Data Types: `logical`

Indicator that coefficient matrix is sparse, returned as either logical `1` (`true`) or logical `0` (`false`).

Data Types: `logical`

Multiplicative scale factor for coefficient matrix, returned as a scalar. The default value of `1` indicates that the coefficient matrix is not scaled. However, when you multiply or divide the `decomposition` object by a scalar, the value of `ScaleFactor` changes. For example, `3*dA` is a `decomposition` object equivalent to `dA`, but with a value of `3` for `ScaleFactor`.

Data Types: `double`
Complex Number Support: Yes

## Object Functions

The primary functions and operators that you can use with `decomposition` objects are related to solving linear systems of equations. If the decomposition type is `'qr'`, then you cannot solve `A'\B` or `B/A`. Instead, use `'cod'` for problems with those forms.

 `ctranspose` Complex conjugate transpose `mldivide` Solve systems of linear equations Ax = B for x `mrdivide` Solve systems of linear equations xA = B for x `isIllConditioned` Determine whether matrix is ill conditioned

You also can check the condition number or rank of the underlying matrix of `decomposition` objects. Since different algorithms are employed, the results of using these functions on the `decomposition` object can differ compared to using the same functions directly on the coefficient matrix.

 `rank` Only the 1-input form `rank(dA)` is supported.The decomposition type must be `'qr'` or `'cod'`.The value of the rank depends on the choice of `RankTolerance`, if specified. `rcond` Runs the same condition check that backslash `\` uses to determine whether to issue a warning.Supports all decomposition types except for `'qr'` and `'cod'`.

## Examples

collapse all

Show how using `decomposition` objects can improve the efficiency of solving $\mathrm{Ax}=\mathit{b}$ with many right-hand sides.

The inverse iteration is an iterative eigenvalue algorithm that solves linear systems with many right-hand sides. It is a method to iteratively compute an eigenvalue of a matrix starting from a guess of the corresponding eigenvector. Each iteration computes `x = A\x`, and then scales `x` by its norm.

Create a sparse matrix `A` and random starting vectors `x1` and `x2`.

```n = 1e3; rng default % for reproducibility A = sprandn(n,n,0.2) + speye(n); x1 = randn(n,1); x2 = x1;```

Apply 100 iterations of the inverse iteration algorithm using backslash to calculate an eigenvalue of `A`.

```tic for ii=1:100 x1 = A \ x1; x1 = x1 / norm(x1); end toc```
```Elapsed time is 18.987628 seconds. ```
`lambda = x1'*A*x1`
```lambda = -0.6707 ```

Now use a `decomposition` object to solve the same problem.

```tic dA = decomposition(A); for ii=1:100 x2 = dA \ x2; x2 = x2 / norm(x2); end toc```
```Elapsed time is 1.176988 seconds. ```
`lambda = x2'*A*x2`
```lambda = -0.6707 ```

The performance of the algorithm improves dramatically because the matrix `A` does not need to be factorized during each iteration. Also, even though the backslash algorithm can be improved by performing an LU decomposition of `A` before the `for`-loop, the `decomposition` object gives access to all of the same performance gains without requiring that you write complex code.

Choose a decomposition type to override the automatic default selection based on the input matrix.

Create a coefficient matrix and decompose the matrix using the default selection of decomposition type.

```A = ones(3); dA = decomposition(A)```
```dA = decomposition with properties: MatrixSize: [3 3] Type: 'ldl' Show all properties ```

Solve the linear system using a vector of ones for the right-hand side.

```b = ones(3,1); x = dA\b```
```Warning: Matrix is singular to working precision. ```
```x = 3×1 NaN NaN NaN ```

Specify the decomposition type to use the `'qr'` method instead of the default `'ldl'` method. This forces backslash (\) to find a least-squares solution to the problem instead of returning a vector of `NaN`s.

`dA_qr = decomposition(A,'qr')`
```dA_qr = decomposition with properties: MatrixSize: [3 3] Type: 'qr' Show all properties ```
`x = dA_qr\b`
```Warning: Rank deficient, rank = 1, tol = 1.153778e-15. ```
```x = 3×1 1.0000 0 0 ```

Specify `'upper'` to use only the upper triangular portion of an input matrix in the decomposition.

Create a coefficient matrix. Construct a triangular decomposition for the matrix using only the upper triangular portion. This option can be useful in cases where both an upper triangular and lower triangular matrix are stored in the same matrix.

`A = randi([0 5],10)`
```A = 10×10 4 0 3 4 2 1 4 5 2 0 5 5 0 0 2 4 1 1 4 0 0 5 5 1 4 3 3 4 3 3 5 2 5 0 4 0 4 1 3 4 3 4 4 0 1 0 5 5 5 5 0 0 4 4 2 2 5 2 1 0 1 2 4 4 2 5 3 1 4 3 3 5 2 1 3 2 0 1 4 2 5 4 3 5 4 3 0 3 2 0 5 5 1 0 4 1 1 2 3 2 ```
`dA = decomposition(A,'triangular','upper')`
```dA = decomposition with properties: MatrixSize: [10 10] Type: 'triangular' Show all properties ```

Use the `'CheckCondition'` name-value pair to turn off warnings based on the condition of the coefficient matrix when solving a linear system using `decomposition`.

Create a coefficient matrix that is ill conditioned. In this matrix, averaging together the first two columns produces the third column.

`A = [1 2 1.5; 3 4 3.5; 5 6 5.5]`
```A = 3×3 1.0000 2.0000 1.5000 3.0000 4.0000 3.5000 5.0000 6.0000 5.5000 ```

Solve a linear system $\mathrm{Ax}=\mathit{b}$ using a vector of 1s for the right-hand side. `mldivide` produces a warning about the conditioning of the coefficient matrix.

```b = ones(3,1); x = A\b```
```Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND = 8.326673e-18. ```
```x = 3×1 -0.0556 1.9444 -1.8889 ```

Now create a `decomposition` object for the matrix and solve the same linear system. Specify `'CheckCondition'` as `false` so that `mldivide` does not check the condition of the coefficient matrix. Even though the same solution is returned, `mldivide` does not display the warning message.

```dA = decomposition(A,'CheckCondition',false); x = dA\b```
```x = 3×1 -0.0556 1.9444 -1.8889 ```

Use the `isIllConditioned` function to check whether the `decomposition` object is based on an ill-conditioned matrix.

`tf = isIllConditioned(dA)`
```tf = logical 1 ```