Skip to main content

Section 2.1 Solving and Graphing Systems of Linear Equations

Systems of linear equations can be solved using multiple methods. In this section we will learn how to solve systems of linear equations using a process called row reduction.
 1 
Row reduction is sometimes called Gaussian elimination, after the mathematician Gauss, whose notation was adopted by . The method was actually popularized in Europe by Cambridge University who published notes of Newton. Soon algebra textbooks across Europe included the method. Long before Gauss or Newton, Chinese scholars demonstrated this method of solving systems of linear equations in chapter 8 of The Nine Chapters on the Mathematical Art.

Subsection Row Reduction

Note that row reduction involves symbolic manipulation. Since Numpy is focused on number crunching, we will need a different Python library called Sympy to do the symbolic manipulation of row reduction to reduced echelon form and row reduced echelon form.
In order to row reduce a system of linear equations we convert our system into an augmented matrix.

Example 2.1. Augmented Matrix.

For example, the system
\begin{align*} 2x \amp -3y \amp= 4 \\ x \amp +4y \amp =2 \end{align*}
becomes
\begin{equation*} \left[ \begin{matrix} 2 \amp -3 \amp | \amp 4\\ 1 \amp 4 \amp | \amp 2 \end{matrix} \right]. \end{equation*}
In Sympy, matrices are entered using the function Matrix.

You Try 2.2.

Edit the above code to print the augmented matrix for the system of linear equations in Example 2.1.
Note that sympy will automatically convert decimals to floating point numbers. Often it is helpful to keep fractions exact using sympy’s Rational(a,b) command.

You Try 2.3.

Edit the above code to print a matrix with decimals like \(0.9\) or \(0.1\text{.}\) Also print a matrix with fractions like \(\frac{9}{10}\) and \(\frac{1}{10}\) using sympy’s Rational(a,b) command.
Hint.
Note that Rational(a,b) is not a Python command, but rather a command specific to Sympy.
The Sympy functions echelon_form() and rref() can be used to row reduce a matrix to reduced echelon form and row reduced echelon form, respectively. If you need a refresher on these two forms, read Section 1.2
 2 
Interactive Linear Algebra 1.2 Row Reduction
of Interactive Linear Algebra.

You Try 2.4.

Edit the code below to print the reduced echelon form and the row reduced echelon form for the augmented matrix of the system of linear equations from Example 2.1.
Note that rref() returns not only the row reduced echelon form of a matrix, but also a tuple that tells you which columns have leading values of 1 (where the left-most column is column 0). These columns are called pivot columns, and give you a quick way to see any free variables (which are non-pivot columns).

You Try 2.5.

Edit the code above for various augmented matrices to see how the tuple returned by rref() relates to the pivot columns.
Recall that a system of linear equations either has
  1. no solution,
  2. a unique solution, or
  3. infinitely many solutions.
If you need a refresher on how to interpret the reduced echelon form of an augmented matrix to determine if there is no solution, write down the unique solution, or give the infinitely many solutions in parametric form, see MTH 264 class notes Lesson 3.

Subsection Graphing solutions to systems of equations

In order to graph solutions to systems of linear equations with two or three variables you can graph the corresponding lines, planes, or coordinates.

Example 2.6.

The system
\begin{align*} x \amp +2y \amp= 6 \\ \frac{1}{2}x \amp -y \amp =1 \end{align*}
has solutions \(x=4\) and \(y=1\text{.}\)
Linear equations with two variables can be graphed as lines in the xy-plane, either by plotting two points on the line or solving for \(y\text{.}\) The solution is one point in a scatter plot.

You Try 2.7.

Edit the code in Example 2.6 to clearly indicate the solution \(x=4\) and \(y=1\) on the graph.
Systems of equations with three variables must be graphed in 3D. As mentioned in Section 1.4, Matplotlib paints 2D projections of each surface on top of the other which doesn’t adequately graph their relationships. We can partially fix this by making the surfaces transparent and explicitly plotting their intersections. The intersection of any two non-parallel planes is a line which can be found by solving the system of the two plane equations and writing the solution in parametric form.

Example 2.8.

Graph solutions to the system
\begin{align*} x +y +z \amp= 1 \\ -y +z \amp =1\\ -x +y -z\amp =0 \end{align*}
by finding the equations of the lines of intersection of the planes and plotting the planes, the lines, and the point of intersection.
Solution 1 walks through the process, followed by a Python implementation of the rref().
Solution 1.
First row reducing the augmented matrix for the system of all three equations we obtain a unique solution \(x=-1\text{,}\) \(y=1/2\text{,}\) \(z=3/2\text{.}\)
Row reducing the augmented matrix for the first two equations we obtain
\begin{align*} x +2z \amp= 2 \\ y -z \amp =-1 \end{align*}
which can be written in parametric form,
\begin{equation*} \left[ \begin{matrix} x\\ y\\z \end{matrix} \right] = \left[\begin{matrix} 2-2z\\ -1+z\\z \end{matrix} \right]= \left[\begin{matrix} 2\\ -1\\0 \end{matrix} \right] +z\left[\begin{matrix} -2\\ 1\\1 \end{matrix} \right]. \end{equation*}
Row reducing the augmented matrix for the last two equations we obtain
\begin{align*} x \amp= -1 \\ y -z \amp =-1 \end{align*}
which can be written in parametric form,
\begin{equation*} \left[ \begin{matrix} x\\ y\\z \end{matrix} \right] = \left[\begin{matrix} -1\\ -1+z\\z \end{matrix} \right]= \left[\begin{matrix} -1\\ -1\\0 \end{matrix} \right] +z\left[\begin{matrix} 0\\ 1\\1 \end{matrix} \right]. \end{equation*}
Row reducing the augmented matrix for the first and last equations we obtain
\begin{align*} x +z\amp= 1/2 \\ y \amp =1/2 \end{align*}
which can be written in parametric form,
\begin{equation*} \left[ \begin{matrix} x\\ y\\z \end{matrix} \right] = \left[\begin{matrix} 1/2-z\\ 1/2\\z \end{matrix} \right]= \left[\begin{matrix} 1/2\\ 1/2\\0 \end{matrix} \right] +z\left[\begin{matrix} -1\\ 0\\1 \end{matrix} \right]. \end{equation*}
Next we will use Matplotlib to graph these results. First note that if we are careless with our input domains for the planes and lines, our graph will be inconsistent.
Note: you can change the viewing angle to try and get a better graph.
The real problem in the code above is that we graph our surfaces with domain \(-2\leq x \leq 2\) and \(-2\leq y\leq 2\text{,}\) but use \(-2\leq z\leq 2\) for the parametric equation of our lines. There are several fixes for this.
One fix is to adjust the parametric input values for the lines to better match the planes. For example, you can make \(-1\leq z \leq 3\) to get a more consistent graph. Try several adjustments for the \(z\) range in the above code. To get a fully consistent graph, each line would need a different parametric domain that matches the domain of the planes.
Another fix is to solve each pair of equations for two points on the line in the domain of the planes. You can then plot the line between the two points.
Solution 2 walks through the process of finding the points, followed by an alternative Python implementation to graph the solutions.
Solution 2.
Recall the first two equations intersect in the line
\begin{equation*} \left[ \begin{matrix} x\\ y\\z \end{matrix} \right] = \left[\begin{matrix} 2-2z\\ -1+z\\z \end{matrix} \right]= \left[\begin{matrix} 2\\ -1\\0 \end{matrix} \right] +z\left[\begin{matrix} -2\\ 1\\1 \end{matrix} \right]. \end{equation*}
Solving for when \(x=2\text{,}\) we see \(t=0\text{.}\) Plugging \(t=0\) into the parametric equation for the line we obtain the point \((2,-1,0)\text{.}\) Similarly when \(x=-2\text{,}\) we find \(t=2\text{,}\) giving the point \((-2,1,2)\text{.}\)
The last two equations intersect in the line
\begin{equation*} \left[ \begin{matrix} x\\ y\\z \end{matrix} \right] = \left[\begin{matrix} -1\\ -1+z\\z \end{matrix} \right]= \left[\begin{matrix} -1\\ -1\\0 \end{matrix} \right] +z\left[\begin{matrix} 0\\ 1\\1 \end{matrix} \right]. \end{equation*}
Here \(x\) is always \(-1\text{.}\) Instead we solve for when \(y=-2\) and \(2\text{.}\) We find \(t=-1\) and \(3\text{,}\) giving us the two points \((-1,2,-1)\) and \((-1,2,3)\text{.}\)
The first and last equations intersect in the line
\begin{equation*} \left[ \begin{matrix} x\\ y\\z \end{matrix} \right] = \left[\begin{matrix} 1/2-z\\ 1/2\\z \end{matrix} \right]= \left[\begin{matrix} 1/2\\ 1/2\\0 \end{matrix} \right] +z\left[\begin{matrix} -1\\ 0\\1 \end{matrix} \right]. \end{equation*}
Here \(x=-2,2\) when \(t=5/2,-3/2\) giving us the points \((-2,1/2,5/2)\) and \((2,1/2,-3/2)\text{.}\)
Recall to plot the lines between two points we create arrays of the \(x\)-values, \(y\)-values, and \(z\)-values, then plot the line between them.
Can you find other fixes? If so, share them so we can add them here.

Subsection Graphing vector equations and their solutions

Recall that solving a vector equation is equivalent to solving a system of equations.

Example 2.9.

For example the vector equation
\begin{equation*} x\left[ \begin{matrix} 1\\ 1/2 \end{matrix} \right] + y\left[\begin{matrix} 2\\-1 \end{matrix} \right]= \left[\begin{matrix} 6\\1 \end{matrix} \right] \end{equation*}
corresponds to the linear system found in Example 2.6 and still has solution \(x=4\) and \(y=1\text{.}\)
The solution to the vector equation can be graphed using linear combinations of vectors. That is, \(4\vec{u}+\vec{v}=\vec{b}.\)

Subsection Converting between Numpy arrays and Sympy Matrices

Since Numpy and Sympy are different libraries, we will need to either enter our augmented matrix as a Sympy matrix, or convert a numpy array into a sympy matrix in order to call the row reducing function in Sympy.
Note that Numpy matrices can be coded as 2D Numpy arrays.

Subsection Round-off error in Sympy’s rref()

Whenever you use a computer to perform arithmetic you have to consider round-off error. Numbers in Python are either stored as an integer data type, a floating point number data type, or a complex number data type. Although \(0.1\) has a finite decimal representation, in binary it requires an infinite repeating floating point representation which cannot be stored without being truncated or rounded.

Example 2.10.

Consider the system of equations
\begin{align*} 0.8x -0.8y -0.4z \amp= 0 \\ -0.3x +0.9y -0.4z \amp =0\\ -0.4x -0.1y +0.8z\amp=0 \end{align*}
If you multiply both sides of each equation by 10, you don’t change the solutions to the equations.
\begin{align*} 8x -8y -4z \amp= 0 \\ -3x +9y -4z \amp =0\\ -4x -y +8z\amp=0 \end{align*}
If you use Sympy to row reduce the system, however, you get a different answer.
A possible fix for the floating point round-off error is to tell Sympy explicitly when to round to zero by making a user defined iszerofunc.

You Try 2.11.

Sympy’s rational number data type stores fractions as two integers to avoid floating point representation. Fix the row reduction round-off error for the original system in Example 2.10 by using Sympy’s Rational() function instead.

Summary.

  • TBD
  • TBD

Exercises Exercises

1.

The system
\begin{align*} x +2y +z \amp= 6 \\ x -y \amp =3\\ -x +y -z\amp =0 \end{align*}
has a unique solution.
Use Python code to find the solution. Graph the system of equations and the solution using Matplot3D. Lable your axes. Note: part of the challenge is making sure you have a good viewing window and perspective so the graph helps illuminate the solution.

2.

Consider the system
\begin{align*} x +2y +z \amp= 6 \\ x -y \amp =3\\ -x +y -z\amp =0 \end{align*}
Convert the system of linear equations to a vector equation.
Find the unique solution and graph your answer as a vector equation using Matplotlib3D. (That is, use vectors instead of planes to visualize the solution to this system.)