Tutorial 08 (AY24/25 Sem 1)¶

In [1]:

# Required imports
import sympy as sym
from ma1522 import Matrix

# Required imports import sympy as sym from ma1522 import Matrix

Question 1¶

Apply Gram-Schmidt Process to convert

(a)¶

$\left\{ \begin{pmatrix} 1 \\ 1 \\ 1 \\ 1 \end{pmatrix}, \begin{pmatrix} 1 \\ -1 \\ 1 \\ 0 \end{pmatrix}, \begin{pmatrix} 1 \\ 1 \\ -1 \\ -1 \end{pmatrix}, \begin{pmatrix} 1 \\ 2 \\ 0 \\ 1 \end{pmatrix} \right\}$ into an orthonormal basis for $\mathbb{R}^4$.

In [2]:

U = Matrix.from_str("1 1 1 1; 1 -1 1 0; 1 1 -1 -1; 1 2 0 1").T
U

U = Matrix.from_str("1 1 1 1; 1 -1 1 0; 1 1 -1 -1; 1 2 0 1").T U

Out[2]:

$\displaystyle \left[\begin{array}{cccc}1 & 1 & 1 & 1\\1 & -1 & 1 & 2\\1 & 1 & -1 & 0\\1 & 0 & -1 & 1\end{array}\right]$

In [3]:

U.gram_schmidt(factor=True, verbosity=1)

U.gram_schmidt(factor=True, verbosity=1)

$\displaystyle v_1 = \left[\begin{array}{c}1\\1\\1\\1\end{array}\right]$

$\displaystyle v_2 = \left[\begin{array}{c}1\\-1\\1\\0\end{array}\right]- \left(\frac{1}{4}\right) \left[\begin{array}{c}1\\1\\1\\1\end{array}\right] = \left(\frac{1}{4}\right) \left[\begin{matrix}3\\-5\\3\\-1\end{matrix}\right]$

$\displaystyle v_3 = \left[\begin{array}{c}1\\1\\-1\\-1\end{array}\right]- \left(\frac{0}{4}\right) \left[\begin{array}{c}1\\1\\1\\1\end{array}\right]- \left(\frac{-1}{\frac{11}{4}}\right) \left[\begin{array}{c}\frac{3}{4}\\- \frac{5}{4}\\\frac{3}{4}\\- \frac{1}{4}\end{array}\right] = \left(\frac{2}{11}\right) \left[\begin{matrix}7\\3\\-4\\-6\end{matrix}\right]$

$\displaystyle v_4 = \left[\begin{array}{c}1\\2\\0\\1\end{array}\right]- \left(\frac{4}{4}\right) \left[\begin{array}{c}1\\1\\1\\1\end{array}\right]- \left(\frac{-2}{\frac{11}{4}}\right) \left[\begin{array}{c}\frac{3}{4}\\- \frac{5}{4}\\\frac{3}{4}\\- \frac{1}{4}\end{array}\right]- \left(\frac{\frac{14}{11}}{\frac{40}{11}}\right) \left[\begin{array}{c}\frac{14}{11}\\\frac{6}{11}\\- \frac{8}{11}\\- \frac{12}{11}\end{array}\right] = \left(\frac{1}{10}\right) \left[\begin{matrix}1\\-1\\-2\\2\end{matrix}\right]$

Out[3]:

$\left[\begin{array}{cccc}1 & 3 & 7 & 1\\1 & -5 & 3 & -1\\1 & 3 & -4 & -2\\1 & -1 & -6 & 2\end{array}\right]\left[\begin{array}{cccc}\frac{1}{2} & 0 & 0 & 0\\0 & \frac{\sqrt{11}}{22} & 0 & 0\\0 & 0 & \frac{1}{\sqrt{110}} & 0\\0 & 0 & 0 & \frac{1}{\sqrt{10}}\end{array}\right]$

(b)¶

$\left\{ \begin{pmatrix} 1 \\ 2 \\ 2 \\ 1 \end{pmatrix}, \begin{pmatrix} 1 \\ 2 \\ 1 \\ 0 \end{pmatrix}, \begin{pmatrix} 1 \\ 0 \\ 1 \\ 0 \end{pmatrix}, \begin{pmatrix} 1 \\ 0 \\ 2 \\ 1 \end{pmatrix} \right\}$ into an orthonormal set. Is the set obtained an orthonormal basis? Why?

In [4]:

U = Matrix.from_str("1 2 2 1; 1 2 1 0; 1 0 1 0; 1 0 2 1").T
U

U = Matrix.from_str("1 2 2 1; 1 2 1 0; 1 0 1 0; 1 0 2 1").T U

Out[4]:

$\displaystyle \left[\begin{array}{cccc}1 & 1 & 1 & 1\\2 & 2 & 0 & 0\\2 & 1 & 1 & 2\\1 & 0 & 0 & 1\end{array}\right]$

In [5]:

U.gram_schmidt(factor=True, verbosity=1) # UserWarning raised because a zero vector is found

U.gram_schmidt(factor=True, verbosity=1) # UserWarning raised because a zero vector is found

$\displaystyle v_1 = \left[\begin{array}{c}1\\2\\2\\1\end{array}\right]$

$\displaystyle v_2 = \left[\begin{array}{c}1\\2\\1\\0\end{array}\right]- \left(\frac{7}{10}\right) \left[\begin{array}{c}1\\2\\2\\1\end{array}\right] = \left(\frac{1}{10}\right) \left[\begin{matrix}3\\6\\-4\\-7\end{matrix}\right]$

$\displaystyle v_3 = \left[\begin{array}{c}1\\0\\1\\0\end{array}\right]- \left(\frac{3}{10}\right) \left[\begin{array}{c}1\\2\\2\\1\end{array}\right]- \left(\frac{- \frac{1}{10}}{\frac{11}{10}}\right) \left[\begin{array}{c}\frac{3}{10}\\\frac{3}{5}\\- \frac{2}{5}\\- \frac{7}{10}\end{array}\right] = \left(\frac{2}{11}\right) \left[\begin{matrix}4\\-3\\2\\-2\end{matrix}\right]$

$\displaystyle v_4 = \left[\begin{array}{c}1\\0\\2\\1\end{array}\right]- \left(\frac{6}{10}\right) \left[\begin{array}{c}1\\2\\2\\1\end{array}\right]- \left(\frac{- \frac{6}{5}}{\frac{11}{10}}\right) \left[\begin{array}{c}\frac{3}{10}\\\frac{3}{5}\\- \frac{2}{5}\\- \frac{7}{10}\end{array}\right]- \left(\frac{\frac{12}{11}}{\frac{12}{11}}\right) \left[\begin{array}{c}\frac{8}{11}\\- \frac{6}{11}\\\frac{4}{11}\\- \frac{4}{11}\end{array}\right] = 0 \left[\begin{matrix}0\\0\\0\\0\end{matrix}\right]$

/tmp/ipykernel_2427/1175356544.py:1: UserWarning: Vectors are linearly dependent. Note that there is no QR factorisation
  U.gram_schmidt(factor=True, verbosity=1) # UserWarning raised because a zero vector is found

Out[5]:

$\left[\begin{array}{cccc}1 & 3 & 4 & 0\\2 & 6 & -3 & 0\\2 & -4 & 2 & 0\\1 & -7 & -2 & 0\end{array}\right]\left[\begin{array}{cccc}\frac{1}{\sqrt{10}} & 0 & 0 & 0\\0 & \frac{1}{\sqrt{110}} & 0 & 0\\0 & 0 & \frac{1}{\sqrt{33}} & 0\\0 & 0 & 0 & 0\end{array}\right]$

Question 2¶

Let $\mathbf{A} = \begin{pmatrix} 0 & 1 & 1 & 0 \\ 1 & -1 & 1 & -1 \\ 1 & 0 & 1 & 0 \\ 1 & 1 & 1 & 1 \end{pmatrix}$ and $\mathbf{b} = \begin{pmatrix} 6 \\ 3 \\ -1 \\ 1 \end{pmatrix}$.

(a)¶

Is the linear system $\mathbf{A}\mathbf{x} = \mathbf{b}$ inconsistent?

In [6]:

A = Matrix.from_str("0 1 1 0; 1 -1 1 -1; 1 0 1 0; 1 1 1 1")
b = Matrix.from_str("6; 3; -1; 1")
A, b

A = Matrix.from_str("0 1 1 0; 1 -1 1 -1; 1 0 1 0; 1 1 1 1") b = Matrix.from_str("6; 3; -1; 1") A, b

Out[6]:

$\displaystyle \left( \left[\begin{array}{cccc}0 & 1 & 1 & 0\\1 & -1 & 1 & -1\\1 & 0 & 1 & 0\\1 & 1 & 1 & 1\end{array}\right], \ \left[\begin{array}{c}6\\3\\-1\\1\end{array}\right]\right)$

In [7]:

# Check if the system is inconsistent
b.is_subspace_of(A, verbosity=2)

# Check if the system is inconsistent b.is_subspace_of(A, verbosity=2)

Check if span(self) is subspace of span(other)

Before RREF: [other | self]

$\displaystyle \left[\begin{array}{cccc|c}0 & 1 & 1 & 0 & 6\\1 & -1 & 1 & -1 & 3\\1 & 0 & 1 & 0 & -1\\1 & 1 & 1 & 1 & 1\end{array}\right]$

After RREF:

$\text{RREF}\left\{\text{rref} = \left[\begin{array}{cccc|c}1 & 0 & 0 & 1 & 0\\0 & 1 & 0 & 1 & 0\\0 & 0 & 1 & -1 & 0\\0 & 0 & 0 & 0 & 1\end{array}\right], \quad\text{pivots} = (0, 1, 2, 4)\right\}$

Span(self) is not a subspace of span(other).

Out[7]:

False

(b)¶

Find a least squares solution to the system. Is the solution unique?

In [8]:

sol = A.solve_least_squares(b, verbosity=2)
sol

sol = A.solve_least_squares(b, verbosity=2) sol

Before RREF: [self.T @ self | self.T @ rhs]

$\displaystyle \left[\begin{array}{cccc|c}3 & 0 & 3 & 0 & 3\\0 & 3 & 1 & 2 & 4\\3 & 1 & 4 & 0 & 9\\0 & 2 & 0 & 2 & -2\end{array}\right]$

After RREF

$\text{RREF}\left\{\text{rref} = \left[\begin{array}{cccc|c}1 & 0 & 0 & 1 & -6\\0 & 1 & 0 & 1 & -1\\0 & 0 & 1 & -1 & 7\\0 & 0 & 0 & 0 & 0\end{array}\right], \quad\text{pivots} = (0, 1, 2)\right\}$

Out[8]:

$\displaystyle \left[\begin{array}{c}- z - 6\\- z - 1\\z + 7\\z\end{array}\right]$

In [9]:

# separate the solution into particular and general solution
sol.sep_part_gen() 

# separate the solution into particular and general solution sol.sep_part_gen()

Out[9]:

$\left(\left[\begin{array}{c}-6\\-1\\7\\0\end{array}\right] + \left[\begin{array}{c}- z\\- z\\z\\z\end{array}\right]\right)$

(c)¶

Use your answer in (b), compute the projection of $\mathbf{b}$ onto the column space of $\mathbf{A}$. Is the solution unique?

In [10]:

proj = A @ sol
proj

proj = A @ sol proj

Out[10]:

$\displaystyle \left[\begin{array}{c}6\\2\\1\\0\end{array}\right]$

Question 3 (Application)¶

A line $p(x) = a_1 x + a_0$ is said to be the least squares approximating line for a given set of data points $(x_1, y_1), (x_2, y_2), \ldots, (x_m, y_m)$ if the sum $$S = [y_1 - p(x_1)]^2 + [y_2 - p(x_2)]^2 + \cdots + [y_m - p(x_m)]^2$$ is minimized. Writing $$\mathbf{x} = \begin{pmatrix} x_1 \\ x_2 \\ \vdots \\ x_m \end{pmatrix}, \quad \mathbf{y} = \begin{pmatrix} y_1 \\ y_2 \\ \vdots \\ y_m \end{pmatrix}, \quad \text{and} \quad p(\mathbf{x}) = \begin{pmatrix} p(x_1) \\ p(x_2) \\ \vdots \\ p(x_m) \end{pmatrix} = \begin{pmatrix} a_1 x_1 + a_0 \\ a_1 x_2 + a_0 \\ \vdots \\ a_1 x_m + a_0 \end{pmatrix}$$

the problem is now rephrased as finding $a_0, a_1$ such that $$S = ||\mathbf{y} - p(\mathbf{x})||^2$$ is minimized. Observe that if we let $$\mathbf{N} = \begin{pmatrix} 1 & x_1 \\ 1 & x_2 \\ \vdots & \vdots \\ 1 & x_m \end{pmatrix} \quad \text{and} \quad \mathbf{a} = \begin{pmatrix} a_0 \\ a_1 \end{pmatrix}$$

then $\mathbf{N}\mathbf{a} = p(\mathbf{x})$. And so our aim is to find $\mathbf{a}$ that minimizes $||\mathbf{y} - \mathbf{N}\mathbf{a}||^2$.

It is known the equation representing the dependency of the resistance of a cylindrically shaped conductor (a wire) at $20°C$ is given by $$R = \rho \frac{L}{A}$$ where $R$ is the resistance measured in Ohms $\Omega$, $L$ is the length of the material in meters $m$, $A$ is the cross-sectional area of the material in meter squared $m^2$, and $\rho$ is the resistivity of the material in Ohm meters $\Omega m$.

A student wants to measure the resistivity of a certain material. Keeping the cross-sectional area constant at $0.002 m^2$, he connected the power sources along the material at various lengths and measured the resistance and obtained the following data.

$L$	$0.01$	$0.012$	$0.015$	$0.02$
$R$	$2.75 \times 10^{-4}$	$3.31 \times 10^{-4}$	$3.92 \times 10^{-4}$	$4.95 \times 10^{-4}$

It is known that the Ohm meter might not be calibrated. Taking that into account, the student wants to find a linear graph $R = \frac{\rho}{0.002} L + R_0$ from the data obtained to compute the resistivity of the material.

(a)¶

Relabeling, we let $R = y$, $\frac{\rho}{0.002} = a_1$ and $R_0 = a_0$. Is it possible to find a graph $y = a_1 x + a_0$ satisfying the points?

In [11]:

A = Matrix.from_str("1 0.01; 1 0.012; 1 0.015; 1 0.02", aug_pos=1)
# scientific notation supported using E
b = Matrix.from_str("2.75E-4; 3.31E-4; 3.92E-4; 4.95E-4")
A, b

A = Matrix.from_str("1 0.01; 1 0.012; 1 0.015; 1 0.02", aug_pos=1) # scientific notation supported using E b = Matrix.from_str("2.75E-4; 3.31E-4; 3.92E-4; 4.95E-4") A, b

Out[11]:

$\displaystyle \left( \left[\begin{array}{cc}1 & 0.01\\1 & 0.012\\1 & 0.015\\1 & 0.02\end{array}\right], \ \left[\begin{array}{c}0.000275\\0.000331\\0.000392\\0.000495\end{array}\right]\right)$

In [12]:

# Observe that numbers are represented as decimals
# This is bad, because it means that SymPy may not be using 
# full precision arithmetic. Use the `simplify` method to convert
# the numbers to symbolic fractions.
A.simplify(tolerance=1e-10)
b.simplify(tolerance=1e-10)

A, b

# Observe that numbers are represented as decimals # This is bad, because it means that SymPy may not be using # full precision arithmetic. Use the `simplify` method to convert # the numbers to symbolic fractions. A.simplify(tolerance=1e-10) b.simplify(tolerance=1e-10) A, b

Out[12]:

$\displaystyle \left( \left[\begin{array}{cc}1 & \frac{1}{100}\\1 & \frac{3}{250}\\1 & \frac{3}{200}\\1 & \frac{1}{50}\end{array}\right], \ \left[\begin{array}{c}\frac{11}{40000}\\\frac{331}{1000000}\\\frac{49}{125000}\\\frac{99}{200000}\end{array}\right]\right)$

In [13]:

b.is_subspace_of(A, verbosity=2)

b.is_subspace_of(A, verbosity=2)

Check if span(self) is subspace of span(other)

Before RREF: [other | self]

$\displaystyle \left[\begin{array}{cc|c}1 & \frac{1}{100} & \frac{11}{40000}\\1 & \frac{3}{250} & \frac{331}{1000000}\\1 & \frac{3}{200} & \frac{49}{125000}\\1 & \frac{1}{50} & \frac{99}{200000}\end{array}\right]$

After RREF:

$\text{RREF}\left\{\text{rref} = \left[\begin{array}{cc|c}1 & 0 & 0\\0 & 1 & 0\\0 & 0 & 1\\0 & 0 & 0\end{array}\right], \quad\text{pivots} = (0, 1, 2)\right\}$

Span(self) is not a subspace of span(other).

Out[13]:

False

(b)¶

Find the least square approximating line for the data points and hence find the resistivity of the material. Would this material make a good wire?

In [14]:

sol = A.solve_least_squares(b, verbosity=2)
sol

sol = A.solve_least_squares(b, verbosity=2) sol

self.T @ self is invertible. The lest squares solution is unique.

$\displaystyle \mathbf{x} = \left(\mathbf{A}^\top \mathbf{A}\right)^{-1} \mathbf{A}^\top \mathbf{b}$

$\displaystyle \left[\begin{array}{c}\frac{14803}{227000000}\\\frac{4907}{227000}\end{array}\right]$

Out[14]:

$\displaystyle \left[\begin{array}{c}\frac{14803}{227000000}\\\frac{4907}{227000}\end{array}\right]$

In [15]:

# Unlike MATLAB, we can get the exact solution with rational numbers
# However, as it is unwieldy and MA1522 does not require exact solution,
# we can use `evalf(n)` to evaluate the solution to `n` significant values.
sol.evalf(4) 

# Unlike MATLAB, we can get the exact solution with rational numbers # However, as it is unwieldy and MA1522 does not require exact solution, # we can use `evalf(n)` to evaluate the solution to `n` significant values. sol.evalf(4)

Out[15]:

$\displaystyle \left[\begin{matrix}6.521 \cdot 10^{-5}\\0.02162\end{matrix}\right]$

Question 4 (Application)¶

Suppose the equation governing the relation between data pairs is not known. We may want to then find a polynomial $$p(x) = a_0 + a_1 x + a_2 x^2 + \cdots + a_n x^n$$ of degree $n$, $n \leq m - 1$, that best approximates the data pairs $(x_1, y_1), (x_2, y_2), \ldots, (x_m, y_m)$. A least square approximating polynomial of degree $n$ is such that $$||\mathbf{y} - p(\mathbf{x})||^2$$ is minimized. If we write $$\mathbf{x} = \begin{pmatrix} x_1 \\ x_2 \\ \vdots \\ x_m \end{pmatrix}, \quad \mathbf{y} = \begin{pmatrix} y_1 \\ y_2 \\ \vdots \\ y_m \end{pmatrix}, \quad \mathbf{N} = \begin{pmatrix} 1 & x_1 & x_1^2 & \cdots & x_1^n \\ 1 & x_2 & x_2^2 & \cdots & x_2^n \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 1 & x_m & x_m^2 & \cdots & x_m^n \end{pmatrix} \quad \text{and} \quad \mathbf{a} = \begin{pmatrix} a_0 \\ a_1 \\ \vdots \\ a_n \end{pmatrix}$$

then $p(\mathbf{x}) = \mathbf{N}\mathbf{a}$, and the task is to find $\mathbf{a}$ such that $||\mathbf{y} - \mathbf{N}\mathbf{a}||^2$ is minimized. Observe that $\mathbf{N}$ is a matrix minor of the Vandermonde matrix. If at least $n + 1$ of the $x$-values $x_1, x_2, \ldots, x_m$ are distinct, the columns of $\mathbf{N}$ are linearly independent, and thus $\mathbf{a}$ is uniquely determined by $$\mathbf{a} = (\mathbf{N}^T \mathbf{N})^{-1} \mathbf{N}^T \mathbf{y}$$

We shall now find a quartic polynomial $$p(x) = a_0 + a_1 x + a_2 x^2 + a_3 x^3 + a_4 x^4$$ that is a least square approximating polynomial for the following data points

$x$	$4$	$4.5$	$5$	$5.5$	$6$	$6.5$	$7$	$8$	$8.5$
$y$	$0.8651$	$0.4828$	$2.590$	$-4.389$	$-7.858$	$3.103$	$7.456$	$0.0965$	$4.326$

In [16]:

x = Matrix.from_str("4 4.5 5 5.5 6 6.5 7 8 8.5").T
y = Matrix.from_str("0.8651 0.4828 2.590 -4.389 -7.858 3.103 7.456 0.0965 4.326").T

x, y

x = Matrix.from_str("4 4.5 5 5.5 6 6.5 7 8 8.5").T y = Matrix.from_str("0.8651 0.4828 2.590 -4.389 -7.858 3.103 7.456 0.0965 4.326").T x, y

Out[16]:

$\displaystyle \left( \left[\begin{array}{c}4\\4.5\\5\\5.5\\6\\6.5\\7\\8\\8.5\end{array}\right], \ \left[\begin{array}{c}0.8651\\0.4828\\2.59\\-4.389\\-7.858\\3.103\\7.456\\0.0965\\4.326\end{array}\right]\right)$

In [17]:

# Similar to the previous example, we can use `simplify` to 
# convert the numbers to symbolic fractions.
x.simplify(tolerance=1e-10)
y.simplify(tolerance=1e-10)

x, y

# Similar to the previous example, we can use `simplify` to # convert the numbers to symbolic fractions. x.simplify(tolerance=1e-10) y.simplify(tolerance=1e-10) x, y

Out[17]:

$\displaystyle \left( \left[\begin{array}{c}4\\\frac{9}{2}\\5\\\frac{11}{2}\\6\\\frac{13}{2}\\7\\8\\\frac{17}{2}\end{array}\right], \ \left[\begin{array}{c}\frac{8651}{10000}\\\frac{1207}{2500}\\\frac{259}{100}\\- \frac{4389}{1000}\\- \frac{3929}{500}\\\frac{3103}{1000}\\\frac{932}{125}\\\frac{193}{2000}\\\frac{2163}{500}\end{array}\right]\right)$

In [18]:

# Create the Vandermonde matrix
N = Matrix.create_vander(num_rows=x.rows, num_cols=5)
N

# Create the Vandermonde matrix N = Matrix.create_vander(num_rows=x.rows, num_cols=5) N

Out[18]:

$\displaystyle \left[\begin{array}{ccccc}1 & x_{1} & x_{1}^{2} & x_{1}^{3} & x_{1}^{4}\\1 & x_{2} & x_{2}^{2} & x_{2}^{3} & x_{2}^{4}\\1 & x_{3} & x_{3}^{2} & x_{3}^{3} & x_{3}^{4}\\1 & x_{4} & x_{4}^{2} & x_{4}^{3} & x_{4}^{4}\\1 & x_{5} & x_{5}^{2} & x_{5}^{3} & x_{5}^{4}\\1 & x_{6} & x_{6}^{2} & x_{6}^{3} & x_{6}^{4}\\1 & x_{7} & x_{7}^{2} & x_{7}^{3} & x_{7}^{4}\\1 & x_{8} & x_{8}^{2} & x_{8}^{3} & x_{8}^{4}\\1 & x_{9} & x_{9}^{2} & x_{9}^{3} & x_{9}^{4}\end{array}\right]$

In [19]:

# To substitute the values of x into the Vandermonde matrix
N = N.apply_vander(x)
N

# To substitute the values of x into the Vandermonde matrix N = N.apply_vander(x) N

Out[19]:

$\displaystyle \left[\begin{array}{ccccc}1 & 4 & 16 & 64 & 256\\1 & \frac{9}{2} & \frac{81}{4} & \frac{729}{8} & \frac{6561}{16}\\1 & 5 & 25 & 125 & 625\\1 & \frac{11}{2} & \frac{121}{4} & \frac{1331}{8} & \frac{14641}{16}\\1 & 6 & 36 & 216 & 1296\\1 & \frac{13}{2} & \frac{169}{4} & \frac{2197}{8} & \frac{28561}{16}\\1 & 7 & 49 & 343 & 2401\\1 & 8 & 64 & 512 & 4096\\1 & \frac{17}{2} & \frac{289}{4} & \frac{4913}{8} & \frac{83521}{16}\end{array}\right]$

In [20]:

sol = N.solve_least_squares(y, verbosity=2)
sol

sol = N.solve_least_squares(y, verbosity=2) sol

self.T @ self is invertible. The lest squares solution is unique.

$\displaystyle \mathbf{x} = \left(\mathbf{A}^\top \mathbf{A}\right)^{-1} \mathbf{A}^\top \mathbf{b}$

$\displaystyle \left[\begin{array}{c}- \frac{43029219289}{210853500}\\\frac{535177649129}{3162802500}\\- \frac{8061301133}{162195000}\\\frac{972997028}{158140125}\\- \frac{47784941}{175711250}\end{array}\right]$

Out[20]:

$\displaystyle \left[\begin{array}{c}- \frac{43029219289}{210853500}\\\frac{535177649129}{3162802500}\\- \frac{8061301133}{162195000}\\\frac{972997028}{158140125}\\- \frac{47784941}{175711250}\end{array}\right]$

In [21]:

sol.evalf(8)  # Evaluate the solution to 8 significant figures

sol.evalf(8) # Evaluate the solution to 8 significant figures

Out[21]:

$\displaystyle \left[\begin{matrix}-204.07164\\169.20995\\-49.701292\\6.1527524\\-0.27195152\end{matrix}\right]$

Question 5¶

Let $\mathbf{A} = \begin{pmatrix} 1 & 1 & 0 \\ 1 & 1 & 0 \\ 1 & 1 & 1 \\ 0 & 1 & 1 \end{pmatrix}$.

(a)¶

Find a $QR$ factorization of $\mathbf{A}$.

In [22]:

A = Matrix.from_str("1 1 0; 1 1 0; 1 1 1; 0 1 1")
A

A = Matrix.from_str("1 1 0; 1 1 0; 1 1 1; 0 1 1") A

Out[22]:

$\displaystyle \left[\begin{array}{ccc}1 & 1 & 0\\1 & 1 & 0\\1 & 1 & 1\\0 & 1 & 1\end{array}\right]$

In [23]:

A.gram_schmidt(factor=True, verbosity=1)

A.gram_schmidt(factor=True, verbosity=1)

$\displaystyle v_1 = \left[\begin{array}{c}1\\1\\1\\0\end{array}\right]$

$\displaystyle v_2 = \left[\begin{array}{c}1\\1\\1\\1\end{array}\right]- \left(\frac{3}{3}\right) \left[\begin{array}{c}1\\1\\1\\0\end{array}\right] = 1 \left[\begin{matrix}0\\0\\0\\1\end{matrix}\right]$

$\displaystyle v_3 = \left[\begin{array}{c}0\\0\\1\\1\end{array}\right]- \left(\frac{1}{3}\right) \left[\begin{array}{c}1\\1\\1\\0\end{array}\right]- \left(\frac{1}{1}\right) \left[\begin{array}{c}0\\0\\0\\1\end{array}\right] = \left(\frac{1}{3}\right) \left[\begin{matrix}-1\\-1\\2\\0\end{matrix}\right]$

Out[23]:

$\left[\begin{array}{ccc}\frac{1}{3} & 0 & - \frac{1}{6}\\\frac{1}{3} & 0 & - \frac{1}{6}\\\frac{1}{3} & 0 & \frac{1}{3}\\0 & 1 & 0\end{array}\right]\left[\begin{array}{ccc}\sqrt{3} & 0 & 0\\0 & 1 & 0\\0 & 0 & \sqrt{6}\end{array}\right]$

In [24]:

QR = A.QRdecomposition(verbosity=1)
QR

QR = A.QRdecomposition(verbosity=1) QR

Finding orthogonal basis via Gram-Schmidt process:

$\displaystyle v_1 = \left[\begin{array}{c}1\\1\\1\\0\end{array}\right]$

$\displaystyle v_2 = \left[\begin{array}{c}1\\1\\1\\1\end{array}\right]- \left(\frac{3}{3}\right) \left[\begin{array}{c}1\\1\\1\\0\end{array}\right] = \left[\begin{array}{c}0\\0\\0\\1\end{array}\right]$

$\displaystyle v_3 = \left[\begin{array}{c}0\\0\\1\\1\end{array}\right]- \left(\frac{1}{3}\right) \left[\begin{array}{c}1\\1\\1\\0\end{array}\right]- \left(\frac{1}{1}\right) \left[\begin{array}{c}0\\0\\0\\1\end{array}\right] = \left[\begin{array}{c}- \frac{1}{3}\\- \frac{1}{3}\\\frac{2}{3}\\0\end{array}\right]$

Q matrix:

$\displaystyle \left[\begin{array}{ccc}\frac{\sqrt{3}}{3} & 0 & - \frac{\sqrt{6}}{6}\\\frac{\sqrt{3}}{3} & 0 & - \frac{\sqrt{6}}{6}\\\frac{\sqrt{3}}{3} & 0 & \frac{\sqrt{6}}{3}\\0 & 1 & 0\end{array}\right]$

R matrix: Q.T @ self

$\displaystyle \left[\begin{array}{ccc}\sqrt{3} & \sqrt{3} & \frac{\sqrt{3}}{3}\\0 & 1 & 1\\0 & 0 & \frac{\sqrt{6}}{3}\end{array}\right]$

Out[24]:

$\left[\begin{array}{ccc}\frac{\sqrt{3}}{3} & 0 & - \frac{\sqrt{6}}{6}\\\frac{\sqrt{3}}{3} & 0 & - \frac{\sqrt{6}}{6}\\\frac{\sqrt{3}}{3} & 0 & \frac{\sqrt{6}}{3}\\0 & 1 & 0\end{array}\right]\left[\begin{array}{ccc}\sqrt{3} & \sqrt{3} & \frac{\sqrt{3}}{3}\\0 & 1 & 1\\0 & 0 & \frac{\sqrt{6}}{3}\end{array}\right]$

(b)¶

Use your answer in (a) to find the least square solution to $\mathbf{A}\mathbf{x} = \mathbf{b}$, where $\mathbf{b} = \begin{pmatrix} 1 \\ 1 \\ 0 \\ 0 \end{pmatrix}$.

In [25]:

b = Matrix.from_str("1; 1; 0; 0")
b

b = Matrix.from_str("1; 1; 0; 0") b

Out[25]:

$\displaystyle \left[\begin{array}{c}1\\1\\0\\0\end{array}\right]$

In [26]:

Q, R = QR.Q, QR.R

sol = R.inverse() @ Q.T @ b
sol

Q, R = QR.Q, QR.R sol = R.inverse() @ Q.T @ b sol

Out[26]:

$\displaystyle \left[\begin{array}{c}0\\1\\-1\end{array}\right]$

In [27]:

# Or brute force
sol = A.solve_least_squares(b, verbosity=2)
sol

# Or brute force sol = A.solve_least_squares(b, verbosity=2) sol

self.T @ self is invertible. The lest squares solution is unique.

$\displaystyle \mathbf{x} = \left(\mathbf{A}^\top \mathbf{A}\right)^{-1} \mathbf{A}^\top \mathbf{b}$

$\displaystyle \left[\begin{array}{c}0\\1\\-1\end{array}\right]$

Out[27]:

$\displaystyle \left[\begin{array}{c}0\\1\\-1\end{array}\right]$