Loading [MathJax]/jax/output/HTML-CSS/jax.js
Skip to main content
Library homepage
 

Text Color

Text Size

 

Margin Size

 

Font Type

Enable Dyslexic Font
Mathematics LibreTexts

41.2: Solve Linear Systems

( \newcommand{\kernel}{\mathrm{null}\,}\)

Login with LibreOne to run this code cell interactively.

If you have already signed in, please refresh the page.

import numpy as np
import sympy as sym
import numpy as np
import sympy as sym

When we have the same number of equations as unknowns, we have the following system Ax=b with a square matrix A. In this section, we assume that the matrix A has full rank. That is the system has an unique solution. We talked about many ways to solve this system of equations. Some examples are:

  • Jacobian iteration/ Gauss-Seidel iteration
  • x=A1b
  • Gaussian elimination
  • LU decomposition

In this assignment, we will show that we can solve it by QR decomposion.

Consider the following system of equations:

[522434467][x1x2x3]=[123]

Login with LibreOne to run this code cell interactively.

If you have already signed in, please refresh the page.

A = np.matrix([[5, -2, 2], [4, -3, 4], [4,-6,7]])
b = np.matrix([[1],[2],[3]])
display(sym.Matrix(A))
display(sym.Matrix(b))
A = np.matrix([[5, -2, 2], [4, -3, 4], [4,-6,7]])
b = np.matrix([[1],[2],[3]])
display(sym.Matrix(A))
display(sym.Matrix(b))

Back substitution

Let’s first implement the back substitution in Python. The back substitution function back_subst solves the system Rx=b where R is an upper-triangular matrix.

When we solve for x, we start with xn: xn=bnRn,n. Then we solve for xn1 as xn1=bn1Rn1,nxnRn1,n1. Then we can find xn2,,x1: xn2=bn2Rn2,n1xn1Rn2,nxnRn2,n2. We can solve for all components of x in the reversed order. So this is called back substitution.

Do This

Complete the following code for back substitution.

Login with LibreOne to run this code cell interactively.

If you have already signed in, please refresh the page.

def back_subst(R,b):
    n = R.shape[0]; x = np.zeros(n);
    for i in reversed(range(n)):
        x[i] = b[i]
        for j in range(i+1,n):            
## Your code starts ##            
            x[i] =   # Complet this line of code.
## Your code ends ##            
        x[i] = x[i]/R[i,i]
    return np.matrix(x).T
def back_subst(R,b):
    n = R.shape[0]; x = np.zeros(n);
    for i in reversed(range(n)):
        x[i] = b[i]
        for j in range(i+1,n):            
## Your code starts ##            
            x[i] =   # Complet this line of code.
## Your code ends ##            
        x[i] = x[i]/R[i,i]
    return np.matrix(x).T
Do This

When we solve for Ax=b with QR decomposition. We have the following steps:

  • Find the QR decomposition of A
  • From QRx=b, we obtain Rx=Qb
  • Solve for x using back substitution.

Use these steps to solve Ax=b with the given A and b. Compare the result with np.linalg.solve.

Login with LibreOne to run this code cell interactively.

If you have already signed in, please refresh the page.

## Your code starts ##
x =  
## Your code ends ##
print(type(x))   # x is a column vector in the np.matrix type
np.allclose(x, np.linalg.solve(A,b))
## Your code starts ##
x =  
## Your code ends ##
print(type(x))   # x is a column vector in the np.matrix type
np.allclose(x, np.linalg.solve(A,b))

This page titled 41.2: Solve Linear Systems is shared under a CC BY-NC 4.0 license and was authored, remixed, and/or curated by Dirk Colbry via source content that was edited to the style and standards of the LibreTexts platform.

Support Center

How can we help?