---
title: "LinearSolve"
language: "en"
type: "Symbol"
summary: "LinearSolve[m, b] finds an x that solves the matrix equation m . x == b. LinearSolve[m] generates a LinearSolveFunction[...] that can be applied repeatedly to different b. LinearSolve[a, b] finds an x that solves the array equation a . x == b."
keywords: 
- Cholesky
- cofactor expansion
- direct solver methods
- division-free row reduction
- division of matrices
- Gaussian elimination
- Krylov
- Krylov methods
- LINPACK
- matrix equations
- multifrontal
- multifrontal methods
- one-step row reduction
- solution of linear systems
- sparse matrices
- UMFPACK
- GMRES
- PARDISO
- BiCGSTAB
- CHOLSOL
- CRAMER
- GS_ITER
- LINBCG
- LU_COMPLEX
- LUSOL
- SVSOL
- TRISOL
- LinearSolve
- linsolve
- lsolve
- \
- /
- bicg
- bicgstab
- cgs
- gmres
- linsolve
- lsqr
- luinc
- minres
- mldivide
- mrdivide
- pcg
- qmr
- symmlq
canonical_url: "https://reference.wolfram.com/language/ref/LinearSolve.html"
source: "Wolfram Language Documentation"
related_guides: 
  - 
    title: "Linear Systems"
    link: "https://reference.wolfram.com/language/guide/LinearSystems.en.md"
  - 
    title: "Equation Solving"
    link: "https://reference.wolfram.com/language/guide/EquationSolving.en.md"
  - 
    title: "Matrix Operations"
    link: "https://reference.wolfram.com/language/guide/MatrixOperations.en.md"
  - 
    title: "Finite Mathematics"
    link: "https://reference.wolfram.com/language/guide/FiniteMathematics.en.md"
  - 
    title: "GPU Computing"
    link: "https://reference.wolfram.com/language/guide/GPUComputing.en.md"
  - 
    title: "Matrices and Linear Algebra"
    link: "https://reference.wolfram.com/language/guide/MatricesAndLinearAlgebra.en.md"
  - 
    title: "Systems Modeling"
    link: "https://reference.wolfram.com/language/guide/SystemsModeling.en.md"
  - 
    title: "Finite Fields"
    link: "https://reference.wolfram.com/language/guide/FiniteFields.en.md"
  - 
    title: "GPU Programming"
    link: "https://reference.wolfram.com/language/guide/GPUProgramming.en.md"
  - 
    title: "Symbolic Vectors, Matrices and Arrays"
    link: "https://reference.wolfram.com/language/guide/SymbolicArrays.en.md"
  - 
    title: "Structured Arrays"
    link: "https://reference.wolfram.com/language/guide/StructuredArrays.en.md"
related_functions: 
  - 
    title: "Inverse"
    link: "https://reference.wolfram.com/language/ref/Inverse.en.md"
  - 
    title: "Solve"
    link: "https://reference.wolfram.com/language/ref/Solve.en.md"
  - 
    title: "NullSpace"
    link: "https://reference.wolfram.com/language/ref/NullSpace.en.md"
  - 
    title: "CoefficientArrays"
    link: "https://reference.wolfram.com/language/ref/CoefficientArrays.en.md"
  - 
    title: "CholeskyDecomposition"
    link: "https://reference.wolfram.com/language/ref/CholeskyDecomposition.en.md"
  - 
    title: "PseudoInverse"
    link: "https://reference.wolfram.com/language/ref/PseudoInverse.en.md"
  - 
    title: "LeastSquares"
    link: "https://reference.wolfram.com/language/ref/LeastSquares.en.md"
  - 
    title: "RowReduce"
    link: "https://reference.wolfram.com/language/ref/RowReduce.en.md"
  - 
    title: "LinearSolveFunction"
    link: "https://reference.wolfram.com/language/ref/LinearSolveFunction.en.md"
  - 
    title: "MatrixPower"
    link: "https://reference.wolfram.com/language/ref/MatrixPower.en.md"
  - 
    title: "Adjugate"
    link: "https://reference.wolfram.com/language/ref/Adjugate.en.md"
related_tutorials: 
  - 
    title: "Solving Linear Systems"
    link: "https://reference.wolfram.com/language/tutorial/LinearAlgebra.en.md#22511"
  - 
    title: "Implementation notes: Numerical and Related Functions"
    link: "https://reference.wolfram.com/language/tutorial/SomeNotesOnInternalImplementation.en.md#5783"
---
# LinearSolve

LinearSolve[m, b] finds an x that solves the matrix equation m.x == b. 

LinearSolve[m] generates a LinearSolveFunction[…] that can be applied repeatedly to different b. 

LinearSolve[a, b] finds an x that solves the array equation a.x == b.

## Details and Options

* ``LinearSolve`` works on both numerical and symbolic matrices, as well as ``SparseArray`` objects.

* The argument ``b`` can be either a vector or a matrix. »

* The matrix ``m`` can be square or rectangular. »

* ``LinearSolve[m]`` and ``LinearSolveFunction[…]`` provide an efficient way to solve the same approximate numerical linear system many times.

* ``LinearSolve[m, b]`` is equivalent to ``LinearSolve[m][b]``.

* For underdetermined systems, ``LinearSolve`` will return one of the possible solutions; ``Solve`` will return a general solution. »

* For an ``n1``×…×``nk``×``m`` array ``a`` and an ``n1``×…×``nk``×``d1``×…×``dl`` array ``b``, ``LinearSolve[a, b]`` gives an ``m``×``d1``×…×``dl`` array ``x``, such that ``a.x == b``.

* ``LinearSolve`` has the following options and settings:

|          |           |                                             |
| -------- | --------- | ------------------------------------------- |
| Method   | Automatic | method to use                               |
| Modulus  | 0         | prime modulus to use                        |
| ZeroTest | Automatic | test to determine when expressions are zero |

* ``LinearSolve[m, …Modulus -> p]`` solves rational systems modulo the prime ``p``. If ``p`` is zero, ordinary arithmetic is used. »

* The ``ZeroTest`` option only applies to exact and symbolic matrices.

* With ``Method -> Automatic``, the method is automatically selected depending upon input.

* Explicit ``Method`` settings for exact and symbolic matrices include:

|                            |                                               |
| -------------------------- | --------------------------------------------- |
| "CofactorExpansion"        | Laplace cofactor expansion                    |
| "DivisionFreeRowReduction" | Bareiss method of division-free row reduction |
| "OneStepRowReduction"      | standard row reduction                        |

* Explicit ``Method`` settings for approximate numeric matrices include:

|                |                                                          |
| -------------- | -------------------------------------------------------- |
| "Banded"       | banded matrix solver                                     |
| "Cholesky"     | Cholesky method for positive definite Hermitian matrices |
| "Krylov"       | iterative Krylov sparse solver                           |
| "Multifrontal" | direct sparse LU decomposition                           |
| "Mumps"        | parallel direct sparse solver                            |
| "Pardiso"      | parallel direct sparse solver                            |

---

## Examples (53)

### Basic Examples (3)

Solve the matrix-vector equation $m.x=b$ with $m=\left(
\begin{array}{cc}
 r & s \\
 t & u \\
\end{array}
\right)$ and $b=\left(
\begin{array}{c}
 y \\
 z \\
\end{array}
\right)$ :

```wl
In[1]:= LinearSolve[{{r, s}, {t, u}}, {y, z}]

Out[1]= {(u y - s z/-s t + r u), (t y - r z/s t - r u)}
```

Verify the solution:

```wl
In[2]:= {{r, s}, {t, u}}.% == {y, z}//Simplify

Out[2]= True
```

---

Solve the matrix equation $m.x=b$ with $m=\left(
\begin{array}{cc}
 1 & 2 \\
 3 & 4 \\
\end{array}
\right)$ and $b=\left(
\begin{array}{cc}
 5 & 6 \\
 7 & 8 \\
\end{array}
\right)$ :

```wl
In[1]:=
m = {{1, 2}, {3, 4}};
b = {{5, 6}, {7, 8}};
LinearSolve[m, b]//MatrixForm

Out[1]//MatrixForm=
(⁠|    |    |
| -- | -- |
| -3 | -4 |
| 4  | 5  |⁠)
```

Verify the solution:

```wl
In[2]:= m.% == b

Out[2]= True
```

---

Solve a rectangular matrix equation:

```wl
In[1]:=
m = (⁠|   |   |
| - | - |
| 1 | 5 |
| 2 | 6 |
| 3 | 7 |
| 4 | 8 |⁠);b = (⁠|    |
| -- |
| 9  |
| 10 |
| 11 |
| 12 |⁠);

In[2]:= LinearSolve[m, b]

Out[2]= {-1, 2}
```

Verify the solution:

```wl
In[3]:= m.% == b

Out[3]= True
```

### Scope (16)

#### Basic Uses (9)

Solve $m.x=b$ at machine precision:

```wl
In[1]:=
m = (⁠|     |      |      |
| --- | ---- | ---- |
| 1.2 | 3.2  | 3.2  |
| 7.9 | -1.4 | 5.1  |
| 1.1 | 2.5  | -1.5 |⁠);

In[2]:= LinearSolve[m, {.5, 1.3, -2.2}]

Out[2]= {-0.309701, -0.362687, 0.635075}
```

Solve a case where $b$ is a matrix:

```wl
In[3]:= LinearSolve[m, {{1.5, 4.75, -3.2}, {-1.7, 6.7, -9.3}, {4.9, -8.65, 15.4}}]//MatrixForm

Out[3]//MatrixForm=
(⁠|           |          |          |
| --------- | -------- | -------- |
| 0.577099  | -1.30145 | 2.11617  |
| 1.16092   | -1.06494 | 2.59547  |
| -0.908587 | 3.03736  | -4.38903 |⁠)
```

---

Solve $m.x=b$ for a complex matrix:

```wl
In[1]:= LinearSolve[{{1 + I, 2, 3 - 2 I}, {0, 4, 5I}, {1 + I, 6, 3 + 3I}}, {1, 2, 3}]

Out[1]= {0, (1/2), 0}
```

---

Find a solution for an exact, rectangular matrix:

```wl
In[1]:= (m = {{1, 2, 3, 4}, {5, 6, 7, 8}, {9, 10, 11, 12}})//MatrixForm

Out[1]//MatrixForm=
(⁠|   |    |    |    |
| - | -- | -- | -- |
| 1 | 2  | 3  | 4  |
| 5 | 6  | 7  | 8  |
| 9 | 10 | 11 | 12 |⁠)

In[2]:= LinearSolve[m, {13, 14, 15}]

Out[2]= {-(25/2), (51/4), 0, 0}
```

---

Compute a solution at arbitrary precision:

```wl
In[1]:= (m = RandomReal[4, {2, 3}, WorkingPrecision -> 20])//MatrixForm

Out[1]//MatrixForm=
(⁠|                       |                        |                        |
| --------------------- | ---------------------- | ---------------------- |
| 3.6861532821998711058 | 2.3679012393477645585  | 1.3808843904771543609  |
| 1.4893310201074335201 | 0.50491076107417303962 | 0.14079547574411953659 |⁠)

In[2]:= LinearSolve[m, {Pi, E}]

Out[2]= {2.9124348192863460180, -3.2070968115239222589, 0``20.}
```

---

Solve $m.x=b$ for a symbolic matrix:

```wl
In[1]:= LinearSolve[{{a, b, c}, {d, e, f}, {g, h, i}}, {3, 2, 1}]

Out[1]= {(c e - b f - 2 c h + 3 f h + 2 b i - 3 e i/c e g - b f g - c d h + a f h + b d i - a e i), (c d - a f - 2 c g + 3 f g + 2 a i - 3 d i/-c e g + b f g + c d h - a f h - b d i + a e i), (b d - a e - 2 b g + 3 e g + 2 a h - 3 d h/c e g - b f g - c d h + a f h + b d i - a e i)}
```

Solve the system when $b$ is a matrix:

```wl
In[2]:= LinearSolve[{{a, b, c}, {d, e, f}, {g, h, i}}, DiagonalMatrix[{1, 0, -1}]]//MatrixForm

Out[2]//MatrixForm=
(⁠|                                                            |   |                                                             |
| ---------------------------------------------------------- | - | -------------------------------------------------- ... /-c e g + b f g + c d h - a f h - b d i + a e i) | 0 | (-c d + a f/-c e g + b f g + c d h - a f h - b d i + a e i) |
| (e g - d h/c e g - b f g - c d h + a f h + b d i - a e i)  | 0 | (-b d + a e/c e g - b f g - c d h + a f h + b d i - a e i)  |⁠)
```

---

Solve $m.x=b$ over a finite field:

```wl
In[1]:=
ℱ = FiniteField[29, 4];
m = {{ℱ[12], ℱ[23], ℱ[34]}, {ℱ[45], ℱ[56], ℱ[67]}, {ℱ[78], ℱ[89], ℱ[90]}};
b = {{ℱ[123], ℱ[234]}, {ℱ[345], ℱ[456]}, {ℱ[567], ℱ[678]}};
LinearSolve[m, b]//MatrixForm

Out[1]//MatrixForm= [image]

In[2]:= m.% === b

Out[2]= True
```

---

Solve $m.x=b$ for ``CenteredInterval`` matrices:

```wl
In[1]:=
m = Map[CenteredInterval, RandomReal[{-10, 10}, {3, 3}, WorkingPrecision -> 10], {2}];
b = Map[CenteredInterval, RandomReal[{-10, 10}, {3, 2}, WorkingPrecision -> 10], {2}];
(sol = LinearSolve[m, b])//MatrixForm

Out[1]//MatrixForm=
(⁠|                                                                |                                                                 |
| -------------------------------------------------------------- | ---------------------------------------------- ... {-5849645801641, -45, 542762282, -59}, 44}]  | CenteredInterval[{{-9804892802255, -43, 779862720, -59}, 44}]   |
| CenteredInterval[{{-9430695778397, -43, 935092794, -60}, 44}]  | CenteredInterval[{{-13964024333331, -44, 1014801528, -60}, 44}] |⁠)
```

Find random representatives ``mrep`` and ``brep`` of ``m`` and ``b`` :

```wl
In[2]:=
ranrep[e_CenteredInterval] := e["Center"] + RandomInteger[{-1000, 1000}] / 1000 e["Radius"]
(mrep = Map[ranrep, m, {2}])//MatrixForm

Out[2]//MatrixForm=
(⁠|                                                 |                                               |                                                 |
| ----------------------------------------------- | -------------------------------------------- ... 293422105772630194481/36028797018963968000) | (2859375052400080104329/288230376151711744000)  |
| (3002767587093399590663/576460752303423488000)  | (527630747089909040631/57646075230342348800)  | (5428595781884011830961/2305843009213693952000) |⁠)

In[3]:= (brep = Map[ranrep, b, {2}])//MatrixForm

Out[3]//MatrixForm=
(⁠|                                                 |                                                 |
| ----------------------------------------------- | ----------------------------------------------- |
| (334414015954834426303/72057594037927936 ... 8386492610281/57646075230342348800)     |
| -(3895351429657908445429/576460752303423488000) | -(178328624310696070181/36028797018963968000)   |
| -(1733131702738263520209/288230376151711744000) | -(4195039226539635706087/576460752303423488000) |⁠)
```

Verify that ``sol`` contains ``LinearSolve[mrep, brep]`` :

```wl
In[4]:= MapThread[IntervalMemberQ, {sol, LinearSolve[mrep, brep]}, 2]//MatrixForm

Out[4]//MatrixForm=
(⁠|      |      |
| ---- | ---- |
| True | True |
| True | True |
| True | True |⁠)
```

---

Solve $m.x=b$ for $x$ when $b$ is a matrix of different dimensions:

```wl
In[1]:=
m = {{1, 1, 1}, {1, 2, 3}, {1, 4, 9}};
b = {{1, 2}, {3, 4}, {5, 6}};
LinearSolve[m, b]

Out[1]= {{-2, -1}, {4, 4}, {-1, -1}}
```

---

When no right‐hand side $b$ for $m.x=b$ is given, a ``LinearSolveFunction`` is returned:

```wl
In[1]:=
m = {{1, 1, 1}, {1, 2, 3}, {1, 4, 9}};
ls = LinearSolve[m]

Out[1]=
LinearSolveFunction[{3, 3}, {2, True, {{{1, 1, 1}, {1, 1, 2}, {1, 3, 2}}, {1, 2, 3}, 0}, 
  {0, Automatic, Automatic}, 0}]
```

This contains data to solve the problem quickly for a few values of $b$ :

```wl
In[2]:= Map[ls, {{1, 2, 3}, {0, 1, 3}, {1, 0, 2}}]

Out[2]= {{-(1/2), 2, -(1/2)}, {-1, 1, 0}, {4, -5, 2}}
```

#### Special Matrices (6)

Solve $m.x=b$  with sparse matrices:

```wl
In[1]:= m = SparseArray[{{1, 2} -> 1, {3, 1} -> 5, {2, 3} -> 2}, {3, 3}]

Out[1]= SparseArray[Automatic, {3, 3}, 0, {1, {{0, 1, 2, 3}, {{2}, {3}, {1}}}, {1, 2, 5}}]

In[2]:= b = SparseArray[{1, 0, 2}]

Out[2]= SparseArray[Automatic, {3}, 0, {1, {{0, 2}, {{1}, {3}}}, {1, 2}}]
```

As the result is typically not sparse, the result is returned as an ordinary list:

```wl
In[3]:= LinearSolve[m, b]

Out[3]= {(2/5), 1, 0}
```

---

Sparse methods are used to efficiently solve sparse matrices:

```wl
In[1]:= s = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {999, 999}]

Out[1]=
SparseArray[Automatic, {999, 999}, 0, 
 {1, {CompressedData["«2067»"], 
   CompressedData["«3953»"]}, CompressedData["«882»"]}]

In[2]:= b = SparseArray[{499 -> 1}, {999}]

Out[2]= SparseArray[Automatic, {999}, 0, {1, {{0, 1}, {{499}}}, {1}}]

In[3]:= AbsoluteTiming[x = LinearSolve[s, b];]

Out[3]= {0.000895, Null}
```

Visualize the result:

```wl
In[4]:= ListPlot[x]

Out[4]= [image]
```

---

Solve a system with structured matrices:

```wl
In[1]:= m = SymmetrizedArray[{{1, 2} -> 3, {2, 2} -> 5, {2, 3} -> 7, {3, 2} -> 13}, {3, 3}, Symmetric[All]]

Out[1]=
SymmetrizedArray[StructuredArray`StructuredData[{3, 3}, {{{1, 2} -> 3, {2, 2} -> 5, {2, 3} -> 10}, 
   Symmetric[{1, 2}]}]]

In[2]:= b = SymmetrizedArray[{{2} -> 1}, {3}]

Out[2]= SymmetrizedArray[StructuredArray`StructuredData[{3}, {{{2} -> 1}, {}}]]

In[3]:= LinearSolve[m, b]

Out[3]= {(1/3), 0, 0}
```

Use a different type of matrix structure:

```wl
In[4]:= m = QuantityArray[{{1, 2}, {3, 4}}, {"Meters", "Seconds"}]

Out[4]=
QuantityArray[StructuredArray`StructuredData[{2, 2}, {{{1, 2}, {3, 4}}, {"Meters", "Seconds"}, 
   {{1}, {2}}}]]

In[5]:= b = QuantityArray[{5, 6}, "Meters"]

Out[5]= QuantityArray[StructuredArray`StructuredData[{2}, {{5, 6}, "Meters", {{1}}}]]

In[6]:= LinearSolve[m, b]

Out[6]= {-4, Quantity[9/2, "Meters"/"Seconds"]}
```

---

An identity matrix always produces a trivial solution:

```wl
In[1]:= LinearSolve[IdentityMatrix[3], {a, b, c}]

Out[1]= {a, b, c}
```

---

Solve a linear system whose coefficient matrix is a Hilbert matrix:

```wl
In[1]:= LinearSolve[HilbertMatrix[4], {a, b, c, d}]

Out[1]= {4 (4 a - 30 b + 60 c - 35 d), -60 (2 a - 20 b + 45 c - 28 d), 60 (4 a - 45 b + 108 c - 70 d), -140 (a - 12 b + 30 c - 20 d)}
```

---

Solve a $10\times 10$ system whose coefficients are univariate polynomials of degree $100$ :

```wl
In[1]:=
rpoly[n_] := RandomInteger[{-2 ^ 10, 2 ^ 10}, {n + 1}].x ^ Range[0, n]
SeedRandom[1234];
m = Table[rpoly[100], {10}, {10}];
b = Table[rpoly[100], {10}];

In[2]:= LinearSolve[m, b]//Short[#, 3]&//AbsoluteTiming

Out[2]= {0.71148, {(5625322287210941851316952361428 + «1497» + 2578979863996329764873509807266 x^1000) / (3034638622673506403581743718005 + «1485» + 998816316987602849360464925610 x^1000), («1»/«1»), «6», («1»/«1»), («1»/«1»)}}
```

#### Arrays (1)

Solve $a.x=b$ with a 2×3×6 array $a$ and a 2×3×4×5 array $b$:

```wl
In[1]:= a = RandomInteger[{-10, 10}, {2, 3, 6}];

In[2]:= b = RandomInteger[{-10, 10}, {2, 3, 4, 5}];
```

The result is a 6×4×5 array:

```wl
In[3]:= (x = LinearSolve[a, b])//Dimensions

Out[3]= {6, 4, 5}
```

Verify the solution:

```wl
In[4]:= a.x == b

Out[4]= True
```

### Options (7)

#### Method (6)

##### "Banded" (1)

Solve $m.x=b$ using a banded matrix method:

```wl
In[1]:=
m = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {19999, 19999}];
b = Table[0., {19999}];b[[9999]] = 1.;

In[3]:= Timing[x = LinearSolve[m, b, Method -> "Banded"];]

Out[3]= {0., Null}
```

Check a relative error of the computed solution:

```wl
In[4]:= Norm[m.x - b] / Norm[x]

Out[4]= 1.3669827509714134`*^-16
```

##### "Cholesky" (1)

---

Solve $m.x=b$ using the Cholesky decomposition:

```wl
In[1]:=
m = SparseArray[{{i_, i_} -> 2., {i_, j_} /; Abs[i - j] == 1 -> -1.}, {19999, 19999}];
b = Table[0., {19999}];b[[9999]] = 1.;

In[2]:= Timing[x = LinearSolve[m, b, Method -> "Cholesky"];]

Out[2]= {0.062500, Null}
```

Check a relative error of the computed solution:

```wl
In[3]:= Norm[m.x - b] / Norm[x]

Out[3]= 2.898012744861011`*^-16
```

##### "Krylov" (2)

The following suboptions can be specified for the method ``"Krylov"`` :

|                        |                                                                      |
| ---------------------- | -------------------------------------------------------------------- |
| "BasisSize"            | the size of the Krylov basis (GMRES only)                            |
| "MaxIterations"        | the maximum number of iterations                                     |
| "Method"               | methods to be used                                                   |
| "Preconditioner"       | which preconditioner to apply                                        |
| "PreconditionerSide"   | how to apply a preconditioner ("Left" or "Right")                    |
| "ResidualNormFunction" | A norm function that computes a norm of the residual of the solution |
| "StartingVector"       | the initial vector to start iterations                               |
| "Tolerance"            | the tolerance used to terminate iterations                           |

Possible settings for ``"Method"`` include:

|                     |                                                           |
| ------------------- | --------------------------------------------------------- |
| "BiCGSTAB"          | iterative method for arbitrary square matrices            |
| "ConjugateGradient" | iterative method for Hermitian positive definite matrices |
| "GMRES"             | iterative method for arbitrary square matrices            |

Possible settings for ``"Preconditioner"`` include:

|         |                                                                                                 |
| ------- | ----------------------------------------------------------------------------------------------- |
| "ILU0"  | a preconditioner based on an incomplete LU factorization of the original matrix without fill-in |
| "ILUT"  | a variant of ILU0 with fill-in                                                                  |
| "ILUTP" | a variant of ILUT with column permutation                                                       |

Possible suboptions for ``"Preconditioner"`` include:

|     |     |
| --- | --- |
| "FillIn" | upper bound on the number of additional nonzero elements in a row introduced by the ILUT preconditioner |
| "PermutationTolerance" | when to permute columns |
| "Tolerance" | drop tolerance (any element of magnitude smaller than this tolerance is treated as zero) |

---

Solve $m.x=b$ using a Krylov method:

```wl
In[1]:=
m = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {19999, 19999}];
b = Table[0., {19999}];b[[9999]] = 1.;

In[2]:= AbsoluteTiming[x = LinearSolve[m, b, Method -> "Krylov"];]

Out[2]= {0.016112, Null}
```

Check a relative error of the computed solution:

```wl
In[3]:= Norm[m.x - b] / Norm[x]

Out[3]= 1.366982750971415`*^-16
```

##### "Multifrontal" (1)

---

Solve $m.x=b$ using a direct multifrontal method:

```wl
In[1]:=
m = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {19999, 19999}];
b = Table[0., {19999}];b[[9999]] = 1.;

In[2]:= AbsoluteTiming[x = LinearSolve[m, b, Method -> "Multifrontal"];]

Out[2]= {0.059621, Null}
```

Check a relative error of the computed solution:

```wl
In[3]:= Norm[m.x - b] / Norm[x]

Out[3]= 1.3685767137209222`*^-16
```

##### "Pardiso" (1)

---

Solve $m.x=b$ using Pardiso:

```wl
In[1]:=
m = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {19999, 19999}];
b = Table[0., {19999}];b[[9999]] = 1.;

In[2]:= AbsoluteTiming[x = LinearSolve[m, b, Method -> "Pardiso"];]

Out[2]= {0.039811, Null}
```

Check a relative error of the computed solution:

```wl
In[3]:= Norm[m.x - b] / Norm[x]

Out[3]= 1.632657682200382`*^-16
```

#### Modulus (1)

Find the solution ``x`` to ``m.x == b`` modulo 47:

```wl
In[1]:=
m = {{1, 1, 1}, {1, 2, 3}, {1, 4, 9}};
b = {1, 2, 3};
LinearSolve[m, b, Modulus -> 47]

Out[1]= {23, 2, 23}
```

Verify the solution:

```wl
In[2]:= Mod[m.%, 47]

Out[2]= {1, 2, 3}
```

### Applications (11)

#### Spans and Linear Independence (3)

The following three vectors are not linearly independent:

```wl
In[1]:=
v1 = {1, 2, 3};v2 = {4, 5, 6};v3 = {7, 8, 9};
Solve[a v1 + b v2 == v3, {a, b}]

Out[1]= {{a -> -1, b -> 2}}
```

The equation with a generic right-hand side does not have a solution:

```wl
In[2]:= LinearSolve[{v1, v2, v3}, {x, y, z}]
```

LinearSolve::nosol: Linear equation encountered that has no solution.

```wl
Out[2]= LinearSolve[{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}}, {x, y, z}]
```

Equivalently, the equation with the identity matrix on the right-hand side has no solution:

```wl
In[3]:= LinearSolve[{v1, v2, v3}, IdentityMatrix[3]]
```

LinearSolve::nosol: Linear equation encountered that has no solution.

```wl
Out[3]= LinearSolve[{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}}, {{1, 0, 0}, {0, 1, 0}, {0, 0, 1}}]
```

---

The following three vectors are linearly independent:

```wl
In[1]:=
v1 = {1, 2, 3};v2 = {4, 5, 6};v3 = {7, 8, 10};
Solve[a v1 + b v2 == v3, {a, b}]

Out[1]= {}
```

The equation with a generic right-hand side has a solution:

```wl
In[2]:= LinearSolve[{v1, v2, v3}, {x, y, z}]

Out[2]= {(1/3) (-2 x - 4 y + 3 z), (1/3) (-2 x + 11 y - 6 z), x - 2 y + z}
```

Equivalently, the equation with the identity matrix on the right-hand side has a solution:

```wl
In[3]:= LinearSolve[{v1, v2, v3}, IdentityMatrix[3]]

Out[3]= {{-(2/3), -(4/3), 1}, {-(2/3), (11/3), -2}, {1, -2, 1}}
```

The solution is the inverse of $\left\{v_1,v_2,v_3\right\}$ :

```wl
In[4]:= % == Inverse[{v1, v2, v3}]

Out[4]= True
```

---

Determine if the following vectors are linearly independent or not:

```wl
In[1]:=
Subscript[v, 1] = {108, 90, 252, 186};
Subscript[v, 2] = {240, 260, 520, 420};
Subscript[v, 3] = {264, 340, 536, 468};
Subscript[v, 4] = {522, 705, 1038, 929};
```

As $\left\{v_1,v_2,v_3,v_4\right\}.x$ does not have a solution for an arbitrary $b$, they are not linearly independent:

```wl
In[2]:= LinearSolve[{Subscript[v, 1], Subscript[v, 2], Subscript[v, 3], Subscript[v, 4]}, {w, x, y, z}]
```

LinearSolve::nosol: Linear equation encountered that has no solution.

```wl
Out[2]= LinearSolve[{{108, 90, 252, 186}, {240, 260, 520, 420}, {264, 340, 536, 468}, {522, 705, 1038, 929}}, {w, x, y, z}]
```

#### Equation Solving and Invertibility (6)

Solve the following system of equations:

```wl
In[1]:= system = {2x + 3y + 5z == -2, -3x + 4y - z == 7, 4x + y - 6z == 3};
```

Rewrite the system in matrix form:

```wl
In[2]:=
a = {{2, 3, 5}, {-3, 4, -1}, {4, 1, -6}};
b = {-2, 7, 3};
system  == Thread[ a.{x, y, z} == b]

Out[2]= True
```

Use ``LinearSolve`` to find a solution:

```wl
In[3]:= LinearSolve[a, b]

Out[3]= {-(2/3), (73/69), -(53/69)}
```

Show that the solution is unique using ``NullSpace`` :

```wl
In[4]:= NullSpace[a] == {}

Out[4]= True
```

Verify the result using ``SolveValues`` :

```wl
In[5]:= SolveValues[system, {x, y, z}]

Out[5]= {{-(2/3), (73/69), -(53/69)}}
```

---

Find all solutions of the following system of equations:

```wl
In[1]:= system = {10 x + 12 y + 14 z == 1, 22 x + 27 y + 32 z == 0, 34 x + 42 y + 50 z == -1}

Out[1]= {10 x + 12 y + 14 z == 1, 22 x + 27 y + 32 z == 0, 34 x + 42 y + 50 z == -1}
```

First, write the coefficient matrix $m$, variable vector $v$ and constant vector $b$ :

```wl
In[2]:=
m = {{10, 12, 14}, {22, 27, 32}, {34, 42, 50}};
b = {1, 0, -1};
v = {x, y, z};
```

Verify the rewrite:

```wl
In[3]:= system == Thread[m.v == b]

Out[3]= True
```

``LinearSolve`` gives a particular solution:

```wl
In[4]:= x1 = LinearSolve[m, b]

Out[4]= {(9/2), -(11/3), 0}
```

``NullSpace`` gives a basis for solutions to the homogeneous equation $m.b=0$ :

```wl
In[5]:= ns = NullSpace[m]

Out[5]= {{1, -2, 1}}
```

Define $x_0$ to be an arbitrary linear combination of the elements of $*$$\text{ns}$$*$ :

```wl
In[6]:= x0 = Sum[C[i]ns[[i]], {i, Length[ns]}]

Out[6]= {C[1], -2 C[1], C[1]}
```

The general solution is the sum of $x_1$ and $x_0$ :

```wl
In[7]:= x1 + x0

Out[7]= {(9/2) + C[1], -(11/3) - 2 C[1], C[1]}

In[8]:= m.% == b//Simplify

Out[8]= True
```

---

Determine if the following matrix has an inverse:

```wl
In[1]:=
a = (⁠|     |     |      |     |
| --- | --- | ---- | --- |
| 108 | 90  | 252  | 186 |
| 240 | 260 | 520  | 420 |
| 264 | 340 | 536  | 468 |
| 522 | 705 | 1038 | 929 |⁠);
```

Since the system $a.x=I$ has no solution, $a$ does not have an inverse:

```wl
In[2]:= LinearSolve[a, IdentityMatrix[4]]
```

LinearSolve::nosol: Linear equation encountered that has no solution.

```wl
Out[2]= LinearSolve[{{108, 90, 252, 186}, {240, 260, 520, 420}, {264, 340, 536, 468}, {522, 705, 1038, 929}}, {{1, 0, 0, 0}, {0, 1, 0, 0}, {0, 0, 1, 0}, {0, 0, 0, 1}}]
```

Verify the result using ``Inverse`` :

```wl
In[3]:= Inverse[a]
```

Inverse::sing: Matrix {{108,90,252,186},{240,260,520,420},{264,340,536,468},{522,705,1038,929}} is singular.

```wl
Out[3]= Inverse[{{108, 90, 252, 186}, {240, 260, 520, 420}, {264, 340, 536, 468}, {522, 705, 1038, 929}}]
```

---

Determine if the following matrix has a nonzero determinant:

```wl
In[1]:=
a = (|   |    |     |    |
| - | -- | --- | -- |
| 1 | 4  | 2   | -9 |
| 4 | 12 | 2   | 5  |
| 6 | 7  | -11 | 9  |
| 5 | 15 | 10  | 12 |);
```

Since the system $a.x=I$ has a solution, $a$'s determinant must be nonzero:

```wl
In[2]:= LinearSolve[a, IdentityMatrix[4]]

Out[2]= {{(1387/3395), -(2901/3395), (1052/3395), (292/679)}, {-(367/3395), (1311/3395), -(342/3395), -(113/679)}, {(253/3395), -(654/3395), (23/3395), (89/679)}, {-(66/679), (23/679), -(6/679), (2/679)}}
```

Confirm the result using ``Det`` :

```wl
In[3]:= Det[a]

Out[3]= -3395
```

---

Find the inverse of the following matrix:

```wl
In[1]:= m = {{a, b, c}, {d, e, f}, {g, h, i}};
```

To find the inverse, first solve the system $m.x=I$ :

```wl
In[2]:= LinearSolve[m, IdentityMatrix[3]]//MatrixForm

Out[2]//MatrixForm=
(⁠|                                                            |                                                            |                                                            |
| ----------------------------------------------------------  ... (c d - a f/-c e g + b f g + c d h - a f h - b d i + a e i) |
| (e g - d h/c e g - b f g - c d h + a f h + b d i - a e i)  | (b g - a h/-c e g + b f g + c d h - a f h - b d i + a e i) | (b d - a e/c e g - b f g - c d h + a f h + b d i - a e i)  |⁠)
```

Verify the result using ``Inverse`` :

```wl
In[3]:= Simplify[% == Inverse[m]]

Out[3]= True
```

---

Solve the system $m.x=b$, with several different $b$ by means of computing a ``LinearSolveFunction`` :

```wl
In[1]:= m = RandomReal[1, {500, 500}];

In[2]:=
AbsoluteTiming[
	lsf = LinearSolve[m];
	x1 = lsf[UnitVector[500, 1]];
	x250 = lsf[UnitVector[500, 250]];
	x500 = lsf[UnitVector[500, 500]];]

Out[2]= {0.005914, Null}
```

Perform the computation by inverting the matrix and multiplying by the inverse:

```wl
In[3]:=
AbsoluteTiming[inv = Inverse[m];
	y1 = inv.UnitVector[500, 1];
	y250 = inv.UnitVector[500, 250];
	y500 = inv.UnitVector[500, 500];]

Out[3]= {0.026855, Null}
```

The results are practically identical, even though ``LinearSolveFunction`` is multiple times faster:

```wl
In[4]:= Norm[{x1 - y1, x250 - y250, x500 - y500}]

Out[4]= 1.4285339682116405`*^-13
```

#### Calculus (2)

Newton's method for finding a root of a multivariate function:

```wl
In[1]:=
f[{x_, y_, z_}] = {x - Sin[y], z - Cos[y], Cos[x]Sin[z]};
J[{x_, y_, z_}] = Grad[f[{x, y, z}], {x, y, z}];
FixedPoint[(# - LinearSolve[J[#], f[#]])&, {.5, .5, .5}]

Out[1]= {1., 1.5708, 0.}
```

Compare with the answer found by ``FindRoot`` :

```wl
In[2]:= x /. FindRoot[f[x], {x, {.5, .5, .5}}]

Out[2]= {1., 1.5708, -6.468659422812297`*^-29}
```

---

Approximately solve the boundary value problem $u_{x x} + u = f; u(0) = u(1) = 0$ using discrete differences:

```wl
In[1]:=
f[x_] := Exp[-100(x - 1 / 2) ^ 2];
n = 100;h = 1. / n;grid = N[Range[n - 1]] / n;
s = SparseArray[{{i_, i_} -> 1. - 2 / h ^ 2, {i_, j_} /; Abs[i - j] == 1 -> 1 / h ^ 2}, {n - 1, n - 1}];
b = f[grid];
ux = LinearSolve[s, b];
ListPlot[Transpose[{grid, ux}]]

Out[1]= [image]
```

Show the error compared with the exact solution:

```wl
In[2]:=
ex = DSolveValue[{u''[x] + u[x] == f[x], u[0] == u[1] == 0}, u, x][grid];
ListPlot[Transpose[{grid, Abs[ux - ex]}]]

Out[2]= [image]
```

### Properties & Relations (9)

For an invertible matrix $m$, ``LinearSolve[m, b]`` gives the same result as ``SolveValues`` for the corresponding system of equations:

```wl
In[1]:=
m = {{1, 1, 1}, {1, 0, 1}, {0, 1, 1}};
b = {6, 4, 5};

In[2]:= sol = LinearSolve[m, b]

Out[2]= {1, 2, 3}
```

Create the corresponding system of linear equations:

```wl
In[3]:= eqns = Thread[m.{x, y, z} == {6, 4, 5}]

Out[3]= {x + y + z == 6, x + z == 4, y + z == 5}
```

Confirm that ``SolveValues`` gives the same result:

```wl
In[4]:= sol == First[SolveValues[eqns, {x, y, z}]]

Out[4]= True
```

---

``LinearSolve`` always returns the trivial solution $x=0$ to the homogenous equation $m.x=0$ :

```wl
In[1]:= m = {{1, 2, 3}, {4, 5, 6}, {7, 8, 9}};

In[2]:= LinearSolve[m, {0, 0, 0}]

Out[2]= {0, 0, 0}
```

Use ``NullSpace`` to get the complete spanning set of solutions if $m$ is singular:

```wl
In[3]:= NullSpace[m]

Out[3]= {{1, -2, 1}}
```

Compare with the result of ``SolveValues`` :

```wl
In[4]:= SolveValues[{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}}.{x, y, z} == {0, 0, 0}, {x, y, z}]//Quiet

Out[4]= {{x, -2 x, x}}
```

---

If $m$ is nonsingular, the solution $x$ of $m.x = b$ is the inverse of $m$ when $b$ is the identity matrix:

```wl
In[1]:=
m = {{1, 1, 1}, {1, 2, 3}, {1, 4, 9}};
LinearSolve[m, IdentityMatrix[3]] == Inverse[m]

Out[1]= True
```

---

In this case there is no solution to $m.x=b$:

```wl
In[1]:= m = {{1, 2, 3}, {4, 5, 6}, {7, 8, 9}};b = {1, -2, 1};

In[2]:= LinearSolve[m, b]
```

LinearSolve::nosol: Linear equation encountered that has no solution.

```wl
Out[2]= LinearSolve[{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}}, {1, -2, 1}]
```

Use ``LeastSquares`` to minimize $\| m.x-b\|$ :

```wl
In[3]:= x = LeastSquares[m, b]

Out[3]= {0, 0, 0}

In[4]:= Norm[m.x - b]

Out[4]= Sqrt[6]
```

Compare to general minimization:

```wl
In[5]:= Minimize[Norm[m.{x1, x2, x3} - b], {x1, x2, x3}]

Out[5]= {Sqrt[6], {x1 -> 0, x2 -> 0, x3 -> 0}}
```

---

If $m.x=b$ can be solved, ``LeastSquares`` is equivalent to ``LinearSolve`` :

```wl
In[1]:=
m = {{1, 1}, {1, 2}, {1, 9}};
b = {7, 7, 7};

In[2]:= LeastSquares[m, b] == LinearSolve[m, b]

Out[2]= True
```

---

For a square matrix, ``LinearSolve[m, b]`` has a solution for a generic ``b`` iff ``Det[m] != 0`` :

```wl
In[1]:= m = RandomReal[{-1, 1}, {3, 3}];

In[2]:= LinearSolve[m, {Subscript[b, x], Subscript[b, y], Subscript[b, z]}]

Out[2]= {0.893904 (1. Subscript[b, x] + 0.00394667 Subscript[b, y] - 0.548404 Subscript[b, z]), -0.144473 (1. Subscript[b, x] - 9.65098 Subscript[b, y] - 10.0947 Subscript[b, z]), -0.390202 (1. Subscript[b, x] + 0.535568 Subscript[b, y] + 3.48509 Subscript[b, z])}

In[3]:= Det[m] != 0

Out[3]= True
```

---

For a square matrix, ``LinearSolve[m, b]`` has a solution for a generic ``b`` iff ``m`` has full rank:

```wl
In[1]:= m = RandomReal[{-1, 1}, {3, 3}];

In[2]:= LinearSolve[m, {Subscript[b, x], Subscript[b, y], Subscript[b, z]}]

Out[2]= {-1.33419 (1. Subscript[b, x] + 1.97016 Subscript[b, y] + 0.97383 Subscript[b, z]), 1.17314 (1. Subscript[b, x] - 0.0569917 Subscript[b, y] - 0.208923 Subscript[b, z]), -0.464811 (1. Subscript[b, x] - 3.08538 Subscript[b, y] - 4.46255 Subscript[b, z])}

In[3]:= MatrixRank[m] == Length[m]

Out[3]= True
```

---

For a square matrix, ``LinearSolve[m, b]`` has a solution for a generic ``b`` iff ``m`` has an inverse:

```wl
In[1]:= m = {{1, 2, 3}, {4, 5, 6}, {7, 8, 9}};

In[2]:= LinearSolve[m, {Subscript[b, x], Subscript[b, y], Subscript[b, z]}]
```

LinearSolve::nosol: Linear equation encountered that has no solution.

```wl
Out[2]= LinearSolve[{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}}, {Subscript[b, x], Subscript[b, y], Subscript[b, z]}]

In[3]:= Inverse[m]
```

Inverse::sing: Matrix {{1,2,3},{4,5,6},{7,8,9}} is singular.

```wl
Out[3]= Inverse[{{1, 2, 3}, {4, 5, 6}, {7, 8, 9}}]
```

---

For a square matrix, ``LinearSolve[m, b]`` has a solution for a generic ``b`` iff ``m`` has a trivial null space:

```wl
In[1]:= m = RandomReal[{-1, 1}, {3, 3}];

In[2]:= LinearSolve[m, {Subscript[b, x], Subscript[b, y], Subscript[b, z]}]

Out[2]= {0.419861 (1. Subscript[b, x] - 1.39845 Subscript[b, y] + 5.69611 Subscript[b, z]), 1.45569 (1. Subscript[b, x] - 0.099699 Subscript[b, y] + 1.34671 Subscript[b, z]), 1.33585 (1. Subscript[b, x] - 0.550053 Subscript[b, y] - 0.388031 Subscript[b, z])}

In[3]:= NullSpace[m] == {}

Out[3]= True
```

### Possible Issues (4)

Solution found for an underdetermined system is not unique:

```wl
In[1]:= LinearSolve[{{1, 2, 3}, {4, 5, 6}}, {6, 15}]

Out[1]= {0, 3, 0}
```

All solutions are found by ``Solve`` :

```wl
In[2]:= {x, y, z} /. Solve[{{1, 2, 3}, {4, 5, 6}}.{x, y, z} == {6, 15}, {x, y}]

Out[2]= {{z, 3 - 2 z, z}}
```

``LinearSolve`` gave the solution corresponding to $z=0$ :

```wl
In[3]:= % /. z -> 0

Out[3]= {{0, 3, 0}}
```

---

With ill-conditioned matrices, numerical solutions may not be sufficiently accurate:

```wl
In[1]:=
m = HilbertMatrix[20];
b = m.Table[1, {20}];
LinearSolve[N[m], N[b]]
```

LinearSolve::luc: Result for LinearSolve of badly conditioned matrix {{1.,0.5,0.333333,0.25,0.2,<<20>>,0.142857,0.125,0.111111,0.1,<<10>>},<<10>>} may contain significant numerical errors.

```wl
Out[1]= {1., 0.999998, 1.00006, 0.999168, 1.00673, 0.965595, 1.10946, 0.836849, 0.7731, 2.74426, -3.06032, 5.75178, -0.248755, -5.66391, 18.5275, -26.3362, 29.634, -18.121, 8.28211, -0.20038}
```

The solution is more accurate if sufficiently high precision is used:

```wl
In[2]:= LinearSolve[N[m, 30], N[b, 30]]

Out[2]= {1.0000000000000000000000000000, 1.00000000000000000000000000000, 0.99999999999999999999999999999, 1.00000000000000000000000000046, 0.99999999999999999999999998962, 1.00000000000000000000000014777, 0.99999999999999999999999857870, 1.000000000000000 ... 999999799632636, 1.00000000000000000000263945141, 0.99999999999999999999734891393, 1.00000000000000000000199227843, 0.99999999999999999999891596852, 1.00000000000000000000040318336, 0.99999999999999999999990832653, 1.00000000000000000000000960912}
```

---

Some of the linear solvers available are not deterministic. Set up a system of equations:

```wl
In[1]:=
n = 1000;m = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {n, n}];
b = Table[0., {n}];b[[99]] = 1.;
```

Create a solver function:

```wl
In[2]:= solver[method_] := LinearSolve[m, b, Method -> method]
```

The ``"Pardiso"`` solver is not deterministic:

```wl
In[3]:= MinMax[solver["Pardiso"] - solver["Pardiso"]]

Out[3]= {2.6645352591003757`*^-15, 1.7905676941154525`*^-12}
```

The ``Automatic`` solver method is deterministic:

```wl
In[4]:= MinMax[solver[Automatic] - solver[Automatic]]

Out[4]= {0., 0.}
```

---

``LinearSolve`` may not give a result if a non-prime ``Modulus`` is provided:

```wl
In[1]:= LinearSolve[{{2, 0}, {0, 2}}, Modulus -> 4]
```

LinearSolve::modp: Value of option Modulus -> 4 should be a prime number or zero.

```wl
Out[1]= LinearSolve[{{2, 0}, {0, 2}}, Modulus -> 4]
```

### Neat Examples (3)

Solve 100,000 equations using a direct method:

```wl
In[1]:= s = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {10 ^ 5, 10 ^ 5}]

Out[1]= SparseArray[<299998>, {100000, 100000}]

In[2]:= Timing[LinearSolve[s, RandomReal[1, 10 ^ 5]];]

Out[2]= {0.312791, Null}
```

---

Solve a million equations using an iterative method:

```wl
In[1]:= s = SparseArray[{{i_, i_} -> -2., {i_, j_} /; Abs[i - j] == 1 -> 1.}, {10 ^ 6, 10 ^ 6}]

Out[1]= SparseArray[<2999998>, {1000000, 1000000}]

In[2]:= Timing[x1 = LinearSolve[s, b1 = RandomReal[1, 10 ^ 6], Method -> "Krylov"];]

Out[2]= {1.89455, Null}
```

Check a relative error of the solution:

```wl
In[3]:= Norm[s.x1 - b1] / Norm[x1]

Out[3]= 1.5342499634813195`*^-16
```

Solve the same system of equations using a banded matrix method:

```wl
In[4]:= Timing[x2 = LinearSolve[s, b2 = RandomReal[1, 10 ^ 6], Method -> "Banded"];]

Out[4]= {1.73436, Null}
```

Check a relative error of the solution:

```wl
In[5]:= Norm[s.x2 - b2] / Norm[x2]

Out[5]= 1.3029978169639978`*^-16
```

---

## See Also

* [`Inverse`](https://reference.wolfram.com/language/ref/Inverse.en.md)
* [`Solve`](https://reference.wolfram.com/language/ref/Solve.en.md)
* [`NullSpace`](https://reference.wolfram.com/language/ref/NullSpace.en.md)
* [`CoefficientArrays`](https://reference.wolfram.com/language/ref/CoefficientArrays.en.md)
* [`CholeskyDecomposition`](https://reference.wolfram.com/language/ref/CholeskyDecomposition.en.md)
* [`PseudoInverse`](https://reference.wolfram.com/language/ref/PseudoInverse.en.md)
* [`LeastSquares`](https://reference.wolfram.com/language/ref/LeastSquares.en.md)
* [`RowReduce`](https://reference.wolfram.com/language/ref/RowReduce.en.md)
* [`LinearSolveFunction`](https://reference.wolfram.com/language/ref/LinearSolveFunction.en.md)
* [`MatrixPower`](https://reference.wolfram.com/language/ref/MatrixPower.en.md)
* [`Adjugate`](https://reference.wolfram.com/language/ref/Adjugate.en.md)

## Tech Notes

* [Solving Linear Systems](https://reference.wolfram.com/language/tutorial/LinearAlgebra.en.md#22511)
* [Linear Algebra in the Wolfram Language](https://reference.wolfram.com/language/tutorial/LinearAlgebraInMathematicaOverview.en.md)
* [Implementation notes: Numerical and Related Functions](https://reference.wolfram.com/language/tutorial/SomeNotesOnInternalImplementation.en.md#5783)

## Related Guides

* [Linear Systems](https://reference.wolfram.com/language/guide/LinearSystems.en.md)
* [Equation Solving](https://reference.wolfram.com/language/guide/EquationSolving.en.md)
* [Matrix Operations](https://reference.wolfram.com/language/guide/MatrixOperations.en.md)
* [Finite Mathematics](https://reference.wolfram.com/language/guide/FiniteMathematics.en.md)
* [GPU Computing](https://reference.wolfram.com/language/guide/GPUComputing.en.md)
* [Matrices and Linear Algebra](https://reference.wolfram.com/language/guide/MatricesAndLinearAlgebra.en.md)
* [Systems Modeling](https://reference.wolfram.com/language/guide/SystemsModeling.en.md)
* [Finite Fields](https://reference.wolfram.com/language/guide/FiniteFields.en.md)
* [GPU Programming](https://reference.wolfram.com/language/guide/GPUProgramming.en.md)
* [Symbolic Vectors, Matrices and Arrays](https://reference.wolfram.com/language/guide/SymbolicArrays.en.md)
* [Structured Arrays](https://reference.wolfram.com/language/guide/StructuredArrays.en.md)

## Related Links

* [NKS\|Online](http://www.wolframscience.com/nks/search/?q=LinearSolve)
* [A New Kind of Science](http://www.wolframscience.com/nks/)

## History

* Introduced in 1988 (1.0) \| Updated in 1996 (3.0) ▪ 2003 (5.0) ▪ [2014 (10.0)](https://reference.wolfram.com/language/guide/SummaryOfNewFeaturesIn100.en.md) ▪ [2022 (13.2)](https://reference.wolfram.com/language/guide/SummaryOfNewFeaturesIn132.en.md) ▪ [2023 (13.3)](https://reference.wolfram.com/language/guide/SummaryOfNewFeaturesIn133.en.md) ▪ [2024 (14.0)](https://reference.wolfram.com/language/guide/SummaryOfNewFeaturesIn140.en.md) ▪ [2024 (14.1)](https://reference.wolfram.com/language/guide/SummaryOfNewFeaturesIn141.en.md)