System of linear equations¶
Linear systems · Chapter 6 · Lecture notes (PDF) · Exercises: Systems of linear equations · Slides (PDF) · App: Linear systems · App: Rank
1. Existence and uniqueness of solutions¶
Consider a system of \(m\) linear equations with \(n\) variables:
where:
-
\({\boldsymbol A} \in \R^{m \times n}\) is the coefficient matrix (\(m\) equations, \(n\) variables)
-
\({\boldsymbol x} \in \R^{n \times 1}\) is the vector of variables
-
\({\boldsymbol b} \in \R^{m \times 1}\) is the right-hand side vector
-
\(({\boldsymbol A} | {\boldsymbol b}) \in \R^{m \times (n+1)}\) is the augmented matrix
Rouché–Capelli theorem. The existence and uniqueness of solutions to the system \({\boldsymbol A} \; {\boldsymbol x} = {\boldsymbol b}\) depend on the rank of the coefficient matrix \({\boldsymbol A}\) and the rank of the augmented matrix \(({\boldsymbol A} | {\boldsymbol b})\).
Case 1 - No solution (inconsistent system):
The system is inconsistent and has no solution.
Case 2 - Unique solution:
The system has a unique solution.
Case 3 - Infinite solutions:
The system has infinitely many solutions with \(n - \text{rank}({\boldsymbol A})\) degrees of freedom (free variables). If \(r=\text{rank}({\boldsymbol A})\), the solutions can be described by means of \(n-r\) free parameters.
-
Since \(({\boldsymbol A} | {\boldsymbol b})\) is obtained from \({\boldsymbol A}\) by adding one column, we always have
\[ \text{rank}({\boldsymbol A}) \le \text{rank}({\boldsymbol A} | {\boldsymbol b}) \le \text{rank}({\boldsymbol A}) + 1, \]so the three cases above cover all the possibilities. In particular, the system has at least one solution (it is consistent) if and only if \(\text{rank}({\boldsymbol A}) = \text{rank}({\boldsymbol A} | {\boldsymbol b})\).
-
Since \(\text{rank}({\boldsymbol A}) \le \min\{m,n\}\), a unique solution is possible only if \(m \ge n\).
-
If \({\boldsymbol A} \in \R^{n \times n}\) is square, recall that \(\text{rank}({\boldsymbol A})=n\) if and only if \(\det({\boldsymbol A}) \neq 0\). Hence, if \(\det({\boldsymbol A}) \neq 0\), then also \(\text{rank}({\boldsymbol A} | {\boldsymbol b})=n\) and the system has a unique solution, for every right-hand side \({\boldsymbol b}\).
2. Gaussian Elimination Method¶
The Gaussian elimination method is an algorithm for solving systems of linear equations by transforming the augmented matrix \(({\boldsymbol A} | {\boldsymbol b})\) into row echelon form (upper triangular form) through elementary row operations.
-
Given the system of linear equations:
\[ {\boldsymbol A} \; {\boldsymbol x} = {\boldsymbol b} \]we form the augmented matrix:
\[ ({\boldsymbol A} | {\boldsymbol b}) \] -
The method consists of two phases:
Phase 1 - Forward elimination: Transform the augmented matrix into upper triangular form using elementary row operations:
-
Swap two rows
-
Multiply a row by a non-zero scalar
-
Add a multiple of one row to another row
The goal is to create zeros below the diagonal, obtaining:
\[ \left(\begin{array}{cccc|c} \tilde{a}_{11} & \tilde{a}_{12} & \cdots & \tilde{a}_{1n} & \tilde{b}_1 \\ 0 & \tilde{a}_{22} & \cdots & \tilde{a}_{2n} & \tilde{b}_2 \\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 & 0 & \cdots & \tilde{a}_{nn} & \tilde{b}_n \end{array}\right) \]Phase 2 - Backward substitution: Solve the upper triangular system from bottom to top (assuming \(\tilde{a}_{ii} \neq 0\) for all \(i\)):
\[\begin{align*} x_n &= \frac{\tilde{b}_n}{\tilde{a}_{nn}}\\ x_{n-1} &= \frac{\tilde{b}_{n-1} - \tilde{a}_{n-1,n} \; x_n}{\tilde{a}_{n-1,n-1}}\\ &\vdots\\ x_i &= \frac{\tilde{b}_i - \sum_{j=i+1}^{n} \tilde{a}_{ij} \; x_j}{\tilde{a}_{ii}} \end{align*}\] -
Elementary row operations:
-
\(R_i \leftrightarrow R_j\) : swap row \(i\) with row \(j\)
-
\(R_i \leftarrow k \; R_i\) : multiply row \(i\) by scalar \(k \neq 0\)
-
\(R_i \leftarrow R_i + k \; R_j\) : add \(k\) times row \(j\) to row \(i\)
These operations do not change the solution set of the system.
Example 1: Solving a system using Gaussian elimination - Part 1
Consider the system of linear equations:
Augmented matrix:
Step 1: Eliminate \(x_1\) from rows 2 and 3
\(R_2 \leftarrow R_2 - 2 \; R_1\) (subtract \(2\) times row 1 from row 2):
\(R_3 \leftarrow R_3 - 4 \; R_1\) (subtract \(4\) times row 1 from row 3):
Step 2: Eliminate \(x_2\) from row 3
\(R_3 \leftarrow R_3 - 3 \; R_2\) (subtract \(3\) times row 2 from row 3):
Upper triangular form achieved! The corresponding system is:
Example 2: Solving a system using Gaussian elimination - Part 2
Backward substitution:
From the third equation:
From the second equation:
From the first equation:
Solution:
Verification:
3. Gaussian elimination and the Rouché–Capelli theorem¶
-
Gaussian elimination can be applied to any system \({\boldsymbol A} \; {\boldsymbol x} = {\boldsymbol b}\), with \({\boldsymbol A} \in \R^{m \times n}\). In general, forward elimination transforms \(({\boldsymbol A} | {\boldsymbol b})\) into row echelon form: in each nonzero row, the first nonzero entry (called pivot) is strictly to the right of the pivot of the row above, and the rows of zeros (if any) are at the bottom.
-
Elementary row operations do not change the solution set of the system, and they do not change the rank of \({\boldsymbol A}\) nor the rank of \(({\boldsymbol A} | {\boldsymbol b})\).
-
The rank of a matrix in row echelon form is equal to the number of its nonzero rows. Therefore, at the end of forward elimination we can read \(\text{rank}({\boldsymbol A})\) and \(\text{rank}({\boldsymbol A} | {\boldsymbol b})\) and apply the Rouché–Capelli theorem:
-
a row of the form \(\left(\begin{array}{ccc|c} 0 & \cdots & 0 & c \end{array}\right)\) with \(c \neq 0\) corresponds to the impossible equation \(0=c\): then \(\text{rank}({\boldsymbol A}) < \text{rank}({\boldsymbol A} | {\boldsymbol b})\) and the system has no solution;
-
otherwise, the variables whose columns contain a pivot are computed by backward substitution, and the remaining \(n-r\) variables are the free parameters.
-
When the system has infinitely many solutions with one free parameter \(t \in \R\), the solutions can be written in parametric form
where \({\boldsymbol x}_0\) is a particular solution (obtained for \(t=0\)) and \({\boldsymbol v}\) is a nonzero vector such that \({\boldsymbol A} \; {\boldsymbol v} = {\boldsymbol 0}\).
Example 3: Rouché–Capelli: unique solution
Consider the system of linear equations:
The entry in position \((1,1)\) is zero, so we first swap rows 1 and 2.
\(R_1 \leftrightarrow R_2\):
\(R_3 \leftarrow R_3 - 2 \; R_1\):
\(R_3 \leftarrow R_3 + R_2\):
There are \(3\) nonzero rows both in \({\boldsymbol A}\) and in \(({\boldsymbol A} | {\boldsymbol b})\): \(\text{rank}({\boldsymbol A}) = \text{rank}({\boldsymbol A} | {\boldsymbol b}) = 3 = n\), hence the solution is unique. By backward substitution:
that is,
Example 4: Rouché–Capelli: infinitely many solutions
Consider the system of linear equations:
\(R_2 \leftarrow R_2 - 2 \; R_1\):
\(R_3 \leftarrow R_3 - 3 \; R_1\):
\(R_3 \leftarrow R_3 - R_2\):
There are \(2\) nonzero rows both in \({\boldsymbol A}\) and in \(({\boldsymbol A} | {\boldsymbol b})\): \(\text{rank}({\boldsymbol A}) = \text{rank}({\boldsymbol A} | {\boldsymbol b}) = 2 < 3 = n\), hence there are infinitely many solutions with \(n - r = 3 - 2 = 1\) free parameter. The pivots are in the columns of \(x_1\) and \(x_2\), so \(x_3\) is free: we set \(x_3 = t\), \(t \in \R\). By backward substitution:
In parametric form:
Verification: \({\boldsymbol A} \; {\boldsymbol x}_0 = (-3+4, \; -6+10, \; -9+14)' = (1, 4, 5)' = {\boldsymbol b}\) and \({\boldsymbol A} \; {\boldsymbol v} = (7-6-1, \; 14-15+1, \; 21-21)' = {\boldsymbol 0}\), hence \({\boldsymbol A}({\boldsymbol x}_0 + t \; {\boldsymbol v}) = {\boldsymbol b}\) for every \(t \in \R\) \(\checkmark\)
Example 5: Rouché–Capelli: no solution
We change only the last right-hand side of the previous system:
With the same row operations \(R_2 \leftarrow R_2 - 2 \; R_1\), \(R_3 \leftarrow R_3 - 3 \; R_1\):
\(R_3 \leftarrow R_3 - R_2\):
The last row corresponds to the equation \(0 \; x_1 + 0 \; x_2 + 0 \; x_3 = 1\), which is impossible. Indeed \({\boldsymbol A}\) has \(2\) nonzero rows while \(({\boldsymbol A} | {\boldsymbol b})\) has \(3\):
hence the system is inconsistent and has no solution.
Try it in the lab
the same computation step by step: change the matrix and watch the steps change.
4. Homogeneous systems¶
Definition 1: homogeneous system
A system of linear equations \({\boldsymbol A} \; {\boldsymbol x} = {\boldsymbol b}\) with \({\boldsymbol A} \in \R^{m \times n}\) is called homogeneous if \({\boldsymbol b} = {\boldsymbol 0}\), i.e., if it has the form
Observation 1: solutions of a homogeneous system
Let \({\boldsymbol A} \in \R^{m \times n}\). Then:
-
the homogeneous system \({\boldsymbol A} \; {\boldsymbol x} = {\boldsymbol 0}\) always has the solution \({\boldsymbol x} = {\boldsymbol 0}\), called the trivial solution;
-
it has nontrivial solutions (\({\boldsymbol x} \neq {\boldsymbol 0}\)) if and only if \(\text{rank}({\boldsymbol A}) < n\); in this case it has infinitely many solutions, with \(n - \text{rank}({\boldsymbol A})\) free parameters;
-
if \({\boldsymbol A} \in \R^{n \times n}\) is square, it has nontrivial solutions if and only if \(\det({\boldsymbol A}) = 0\).
Proof
Adding a column of zeros does not change the rank, hence \(\text{rank}({\boldsymbol A} | {\boldsymbol 0}) = \text{rank}({\boldsymbol A})\) and the system is always consistent (as shown by \({\boldsymbol x} = {\boldsymbol 0}\)). By the Rouché–Capelli theorem, the solution \({\boldsymbol x} = {\boldsymbol 0}\) is the unique one if and only if \(\text{rank}({\boldsymbol A}) = n\); otherwise there are infinitely many solutions with \(n - \text{rank}({\boldsymbol A})\) free parameters. For a square matrix, \(\text{rank}({\boldsymbol A}) < n\) if and only if \(\det({\boldsymbol A}) = 0\). □
-
In particular, if \(m < n\) (fewer equations than variables), then \(\text{rank}({\boldsymbol A}) \le m < n\) and the homogeneous system always has nontrivial solutions.
-
If \({\boldsymbol x}_0\) is a solution of \({\boldsymbol A} \; {\boldsymbol x} = {\boldsymbol b}\) and \({\boldsymbol v}\) is a solution of \({\boldsymbol A} \; {\boldsymbol x} = {\boldsymbol 0}\), then \({\boldsymbol A}({\boldsymbol x}_0 + {\boldsymbol v}) = {\boldsymbol b} + {\boldsymbol 0} = {\boldsymbol b}\). This explains the parametric form \({\boldsymbol x} = {\boldsymbol x}_0 + t \; {\boldsymbol v}\) of Example Example 4.
Example 6: a homogeneous system
Consider the homogeneous system with the coefficient matrix of Example Example 4:
By the Sarrus rule,
so there are nontrivial solutions. With the same row operations as in Example Example 4 (the last column stays equal to zero):
Setting \(x_3 = t\): \(x_2 = -3t\) and \(x_1 = -2 \; x_2 + x_3 = 6t + t = 7t\). Hence the solutions are
and the trivial solution is obtained for \(t = 0\).
5. LU Factorization Method¶
The LU factorization (or LU decomposition) is a method for solving systems of linear equations by decomposing the coefficient matrix \({\boldsymbol A}\) into the product of two triangular matrices:
where \({\boldsymbol L}\) is a lower triangular matrix (with 1's on the diagonal) and \({\boldsymbol U}\) is an upper triangular matrix.
-
Given the system of linear equations:
\[ {\boldsymbol A} \; {\boldsymbol x} = {\boldsymbol b} \]if we have the LU factorization \({\boldsymbol A} = {\boldsymbol L} \; {\boldsymbol U}\), we can substitute:
\[ {\boldsymbol L} \; {\boldsymbol U} \; {\boldsymbol x} = {\boldsymbol b} \] -
We introduce the auxiliary variable \({\boldsymbol y}\) defined as:
\[ {\boldsymbol y} \;:=\; {\boldsymbol U} \; {\boldsymbol x} \]so that \({\boldsymbol L} \; {\boldsymbol U} \; {\boldsymbol x} = {\boldsymbol b}\) becomes \({\boldsymbol L} \; {\boldsymbol y} = {\boldsymbol b}\).
Why does this work? We are splitting the original system \({\boldsymbol A}\,{\boldsymbol x}={\boldsymbol b}\) into two simpler triangular systems:
-
Find \({\boldsymbol y}\) such that \({\boldsymbol L}\,{\boldsymbol y} = {\boldsymbol b}\).
-
Find \({\boldsymbol x}\) such that \({\boldsymbol U}\,{\boldsymbol x} = {\boldsymbol y}\).
If both steps succeed, then \({\boldsymbol A}\,{\boldsymbol x} = {\boldsymbol L}\,{\boldsymbol U}\,{\boldsymbol x} = {\boldsymbol L}\,{\boldsymbol y} = {\boldsymbol b}\), so \({\boldsymbol x}\) is indeed a solution of the original system. Since both \({\boldsymbol L}\) and \({\boldsymbol U}\) are triangular, each of the two systems can be solved in \(O(n^2)\) operations by simple substitution.
-
-
We then solve the system in two steps:
Step 1 - Forward substitution: Solve \({\boldsymbol L} \; {\boldsymbol y} = {\boldsymbol b}\) for \({\boldsymbol y}\)
Since \({\boldsymbol L}\) is lower triangular, we can solve this system easily from top to bottom:
\[\begin{align*} y_1 &= b_1\\ y_2 &= b_2 - \ell_{21} \; y_1\\ y_3 &= b_3 - \ell_{31} \; y_1 - \ell_{32} \; y_2\\ &\vdots\\ y_i &= b_i - \sum_{j=1}^{i-1} \ell_{ij} \; y_j \end{align*}\]Step 2 - Backward substitution: Solve \({\boldsymbol U} \; {\boldsymbol x} = {\boldsymbol y}\) for \({\boldsymbol x}\)
Since \({\boldsymbol U}\) is upper triangular, we can solve this system easily from bottom to top:
\[\begin{align*} x_n &= \frac{y_n}{u_{nn}}\\ x_{n-1} &= \frac{y_{n-1} - u_{n-1,n} \; x_n}{u_{n-1,n-1}}\\ &\vdots\\ x_i &= \frac{y_i - \sum_{j=i+1}^{n} u_{ij} \; x_j}{u_{ii}} \end{align*}\]
Advantages of LU factorization:
-
Once the factorization \({\boldsymbol A} = {\boldsymbol L} \; {\boldsymbol U}\) is computed, we can solve the system for different right-hand sides \({\boldsymbol b}\) efficiently
-
Both forward and backward substitution require only \(O(n^2)\) operations
-
The factorization itself requires \(O(n^3)\) operations, but it needs to be done only once
Example 7: Solving a system using LU factorization - Part 1
Consider the system of linear equations:
In matrix form: \({\boldsymbol A} \; {\boldsymbol x} = {\boldsymbol b}\) where
Given the LU factorization of \({\boldsymbol A}\):
Verification:
Example 8: Solving a system using LU factorization - Part 2
Step 1 - Forward substitution: Solve \({\boldsymbol L} \; {\boldsymbol y} = {\boldsymbol b}\)
From the first equation:
From the second equation:
From the third equation:
Thus:
Step 2 - Backward substitution: Solve \({\boldsymbol U} \; {\boldsymbol x} = {\boldsymbol y}\)
From the third equation:
From the second equation:
From the first equation:
Solution:
Verification:
6. Cramer's rule¶
Cramer's rule is an explicit formula for the solution of a system of \(m\) linear equations with \(m\) variables (valid whenever the system has a unique solution).
-
Given a column vector \({\boldsymbol b} \in \R^{m \times 1}\) of \(m\) rows, a matrix \({\boldsymbol A } \in \R^{m \times m}\) of \(m\) rows and \(m\) columns and a column vector \({\boldsymbol x} \in \R^{m \times 1}\) of \(m\) rows containing the \(m\) variables:
\[ {\boldsymbol A}= \begin{pmatrix} a_{11} & a_{12} & \dots & a_{1m}\\ a_{21} & a_{22} & \dots & a_{2m}\\ \vdots & \vdots & \dots & \vdots \\ a_{m1} & a_{m2} & \dots & a_{mm}\\ \end{pmatrix} \qquad {\boldsymbol b}= \begin{pmatrix} b_{1} \\ b_{2} \\ \vdots \\ b_{m} \\ \end{pmatrix} \qquad {\boldsymbol x}= \begin{pmatrix} x_{1} \\ x_{2} \\ \vdots \\ x_{m} \\ \end{pmatrix} \]if \(\det({\boldsymbol A}) \neq 0\), the solution \({\tilde{\boldsymbol x}}\) of the system of \(m\) linear equations:
\[ {\boldsymbol A} \; {\boldsymbol x}= {\boldsymbol b} \]is given by the formula:
\[\begin{equation} \label{CCCC} \tilde{x}_j = \frac{\det({\boldsymbol A}_j) }{\det({\boldsymbol A})},~~~~~ \forall j \in \{1,2,\dots,m\} \end{equation}\]where \({\boldsymbol A}_j\) is the matrix formed by replacing the \(j\)-th column of \({\boldsymbol A}\) by the column vector \({\boldsymbol b}\).
6.1 Systems of two equations and two variables¶
With \(m=2\), we have:
and
If \(\det({\boldsymbol A}) \neq 0\), we then have:
Example 9: solution of a system of two equations and two variables
-
Given
\[ {\boldsymbol A}= \begin{pmatrix} -1 & 1 \\[1ex] 8 & 2 \end{pmatrix} \qquad {\boldsymbol b}= \begin{pmatrix} 2 \\[1ex] 19 \end{pmatrix} \]we have
\[ {\boldsymbol A_1}= \begin{pmatrix} 2 & 1 \\[1ex] 19 & 2 \end{pmatrix} ~~~ {\boldsymbol A_2}= \begin{pmatrix} -1 & 2 \\[1ex] 8 & 19 \end{pmatrix} \]and
\[\begin{align*} \det({\boldsymbol A}) &= (-1) \cdot 2 - 1 \cdot 8 =-10\\[1ex] \det({\boldsymbol A}_1) &= 2 \cdot 2 - 1 \cdot 19 = -15\\[1ex] \det({\boldsymbol A}_2) &= (-1) \cdot 19 - 2 \cdot 8 = -35 \end{align*}\]Since \(\det({\boldsymbol A}) \neq 0\), we then have:
\[\begin{equation*} \begin{cases} \begin{array}{rrrrrrrrrrrrr} -x_1 & + & x_2 & = &2\\[2ex] 8\;x_1 & + & 2\;x_2 & = &19 \end{array} \end{cases} ~~\Longrightarrow~~ (\tilde{x}_1, \tilde{x}_2) = \left(\frac{-15}{-10}, \frac{-35}{-10}\right) = \left(\frac{3}{2}, \frac{7}{2}\right) \end{equation*}\] -
Given
\[ {\boldsymbol A}= \begin{pmatrix} -1 & -1 \\[1ex] 1 & -1 \end{pmatrix} \qquad {\boldsymbol b}= \begin{pmatrix} -2 \\[1ex] 0 \end{pmatrix} \]we have
\[ {\boldsymbol A_1}= \begin{pmatrix} -2 & -1 \\[1ex] 0 & -1 \end{pmatrix} ~~~ {\boldsymbol A_2}= \begin{pmatrix} -1 & -2 \\[1ex] 1 & 0 \end{pmatrix} \]and
\[\begin{align*} \det({\boldsymbol A}) &= (-1) \cdot (-1) - (-1) \cdot 1 =2\\[1ex] \det({\boldsymbol A}_1) &= (-2) \cdot (-1) - (-1) \cdot 0 = 2\\[1ex] \det({\boldsymbol A}_2) &= (-1) \cdot 0 - (-2) \cdot 1 = 2 \end{align*}\]Since \(\det({\boldsymbol A}) \neq 0\), we then have:
\[\begin{equation*} \begin{cases} \begin{array}{rrrrrrrrrrrrr} -x_1 & - & x_2 & = &-2\\[2ex] x_1 & - & x_2 & = &0 \end{array} \end{cases} ~~\Longrightarrow~~ (\tilde{x}_1, \tilde{x}_2) = \left(\frac{2}{2}, \frac{2}{2}\right) = \left(1, 1\right) \end{equation*}\]
6.2 Systems of three equations and three variables¶
With \(m=3\), we have:
and (Sarrus rule)
If \(\det({\boldsymbol A}) \neq 0\), we then have:
Example 10: solution of a system of three equations and three variables
We solve with Cramer's rule the system of Example Example 1:
We have
and
Since \(\det({\boldsymbol A}) \neq 0\), we then have:
which is the solution found with Gaussian elimination.
6.3 Computational cost¶
-
Cramer's rule gives an explicit formula, which is very useful for small systems (\(m=2\) or \(m=3\)) and for theoretical purposes. However, it is not a practical method for large systems.
-
To apply Cramer's rule we need \(m+1\) determinants of order \(m\): \(\det({\boldsymbol A})\), \(\det({\boldsymbol A}_1)\), \(\dots\), \(\det({\boldsymbol A}_m)\). If each determinant is computed with the Laplace expansion, the number of arithmetic operations grows roughly like \(m!\) for each determinant.
-
Gaussian elimination, instead, solves the system with a number of operations which grows roughly like \(m^3\) (about \(\frac{2}{3} m^3\) operations).
-
For example, with \(m=10\) we have \(10! = 3\,628\,800\), while \(\frac{2}{3} \cdot 10^3 \approx 667\); with \(m=20\), \(20!\) is larger than \(2 \cdot 10^{18}\), while \(\frac{2}{3}\cdot 20^3 \approx 5\,333\).
-
Even if the determinants are computed in a more efficient way (for instance with Gaussian elimination itself), Cramer's rule still needs \(m+1\) of them, and it remains more expensive than solving the system directly with Gaussian elimination.
Try it in the lab
the same computation step by step: change the matrix and watch the steps change.
Exercises and lab¶
- Exercises · the exercise sheet of this chapter: 11 exercises with worked solutions
- Lab · Linear systems — Gauss and back substitution, Rouché–Capelli, Cramer, or the LU factorization.
- Lab · Rank — the rank of \(\boldsymbol A\) and of the augmented matrix \((\boldsymbol A\mid\boldsymbol b)\) for Rouché–Capelli