ChEn-3170: Computational Methods in Chemical Engineering Fall 2018 UMass Lowell; Prof. V. F. de Almeida 01Oct2018
$ \newcommand{\Amtrx}{\boldsymbol{\mathsf{A}}} \newcommand{\Bmtrx}{\boldsymbol{\mathsf{B}}} \newcommand{\Cmtrx}{\boldsymbol{\mathsf{C}}} \newcommand{\Mmtrx}{\boldsymbol{\mathsf{M}}} \newcommand{\Imtrx}{\boldsymbol{\mathsf{I}}} \newcommand{\Pmtrx}{\boldsymbol{\mathsf{P}}} \newcommand{\Lmtrx}{\boldsymbol{\mathsf{L}}} \newcommand{\Umtrx}{\boldsymbol{\mathsf{U}}} \newcommand{\xvec}{\boldsymbol{\mathsf{x}}} \newcommand{\yvec}{\boldsymbol{\mathsf{y}}} \newcommand{\zvec}{\boldsymbol{\mathsf{z}}} \newcommand{\avec}{\boldsymbol{\mathsf{a}}} \newcommand{\bvec}{\boldsymbol{\mathsf{b}}} \newcommand{\cvec}{\boldsymbol{\mathsf{c}}} \newcommand{\rvec}{\boldsymbol{\mathsf{r}}} \newcommand{\norm}[1]{\bigl\lVert{#1}\bigr\rVert} \DeclareMathOperator{\rank}{rank} \DeclareMathOperator{\abs}{abs} $
your name
¶Context | Points |
---|---|
Precision of the answer | 80% |
Answer Markdown readability | 10% |
Code readability | 10% |
Program a solution algorithm, as a Python function, for a system of equations with a lower triangular matrix of coefficients $\Lmtrx\,\xvec=\bvec$ (suggestion: use a test matrix and right side vector of variable size to test your code). The algorithm for $\Lmtrx\,\xvec=\bvec$ is as follows: $x_i = \Bigl(b_i - \sum\limits_{j=1}^{i-1} L_{i,j}\,x_j \Bigr)\,L^{-1}_{i,i} \ \ \forall \ \ i=1,\ldots,m$ for $i$ and $j$ with offset 1. Recall that numpy
and Python have offset 0 for their sequence data types. Provide an output example for $\Lmtrx$, $\bvec$, and $\xvec$.
from chen_3170.help import get_triangular_matrix
def forward_solve(l_mtrx, b_vec):
'''
Performs a forward solve with a lower triangular matrix and right side vector
Parameters
----------
l_mtrx: numpy.ndarray
Lower triangular matrix is required.
b_vec: numpy.ndarray
Vector is required.
Returns
-------
x_vec: numpy.narray
Solution vector returned.
'''
l_mtrx = [[0.901 0. 0. 0. 0. 0. 0. ] [0.688 0.073 0. 0. 0. 0. 0. ] [0.704 0.523 0.01 0. 0. 0. 0. ] [0.092 0.039 0.12 0.073 0. 0. 0. ] [0.899 0.648 0.666 0.696 0.476 0. 0. ] [0.933 0.484 0.739 0.1 0.817 0.799 0. ] [0.378 0.836 0.374 0.427 0.683 0.421 0.412]] b_vec = [0.369 0.479 0.332 0.92 0.359 0.009 0.959] x_vec = [ 0.41 2.691 -131.39 227.19 -151.995 246.354 -119.523]
y_vec = [ 0.41 2.691 -131.39 227.19 -151.995 246.354 -119.523] ||x_vec - y_vec|| = 3.859383e-11
'''Read the image into a matrix'''
from matplotlib import pyplot as plt # import the pyplot function of the matplotlib package
plt.rcParams['figure.figsize'] = [20, 4] # extend the figure size on screen output
img_mtrx = plt.imread('https://raw.githubusercontent.com/dpploy/chen-3170/master/notebooks/images/cermet.png',format='png')
img_mtrx = img_mtrx.astype('float64')
plt.figure(1) # create a figure placeholder
plt.imshow( img_mtrx,cmap='gray')
plt.axis('off')
print('')
'''Set image matrix to lower triangular'''
'''Solve L y = b for b equal to a vector of ones'''
||x_vec - y_vec|| = 1.923128e-13
Program a solution algorithm, as a Python function, for a system of equations with an upper triangular matrix of coefficients $\Umtrx\,\xvec=\bvec$ (suggestion: use a test matrix and right side vector of variable size to test your code). The algorithm for $\Umtrx\,\xvec=\bvec$ is as follows: $x_i = \Bigl(b_i - \sum\limits_{j=i+1}^{m} U_{i,j}\,x_j \Bigr)\,U^{-1}_{i,i} \ \ \forall \ \ i=m,\ldots,1$ for $i$ and $j$ with offset 1. Recall that numpy
and Python have offset 0 for their sequence data types. Provide an output example for $\Umtrx$, $\bvec$, and $\xvec$.
def backward_solve(u_mtrx, b_vec):
'''
Performs a backward solve with an upper triangular matrix and right side vector
Parameters
----------
u_mtrx: numpy.ndarray
Upper triangular matrix is required.
b_vec: numpy.ndarray
Vector is required.
Returns
-------
x_vec: numpy.narray
Solution vector returned.
'''
return x_vec
u_mtrx = [[0.59 0.444 0.905 0.874 0.261 0.019 0.628] [0. 0.443 0.338 0.494 0.573 0.444 0.155] [0. 0. 0.989 0.749 0.566 0.436 0.567] [0. 0. 0. 0.07 0.042 0.784 0.987] [0. 0. 0. 0. 0.146 0.908 0.516] [0. 0. 0. 0. 0. 0.593 0.776] [0. 0. 0. 0. 0. 0. 0.823]] b_vec = [0.56 0.984 0.799 0.86 0.13 0.729 0.264] x_vec = [-6.191 4.605 1.872 1.869 -5.271 0.809 0.32 ]
y_vec = [-6.191 4.605 1.872 1.869 -5.271 0.809 0.32 ] ||x_vec - y_vec|| = 3.330669e-15
'''Set image matrix to upper triangular'''
'''Solve U y = b for b equal to a vector of ones'''
||x_vec - y_vec|| = 2.304705e-13
Program (in a function) a $\Lmtrx\,\Umtrx$ factorization algorithm (without using pivoting) for the given matrix $\Amtrx$ below and compute the $\Lmtrx\,\Umtrx$ factors. Compute the $\max\bigl(\abs(\Amtrx-\Lmtrx\,\Umtrx)\bigr)$. The factorization is obtained by elimination steps $k = 1,\ldots,m-1$ so that $A^{(k+1)}_{i,j} = A^{(k)}_{i,j} - A^{(k)}_{k,j}\, m_{i,k} \ \forall\ i=k+1,\ldots,m \ \text{and}\ j=k+1,\ldots,m$ where the multipliers $m_{i,k}$ are given by $m_{i,k} = \frac{A^{(k)}_{i,k}}{A^{(k)}_{k,k}}$. When $k = m-1$, $A^{(m)}_{i,j}$, is upper triangular, that is, $U_{i,j} = A^{(m-1)}_{i,j}$ . The lower triangular matrix is obtained using the multipliers $m_{i,k}$, that is $L_{i,j} = m_{i,j} \ \forall \ i>j$, $L_{i,i}=1$, and $L_{i,j}=0 \ \forall \ i<j$.
2. Print the smallest diagonal entry of the $\Umtrx$ factor, compute the $\rank(\Amtrx)$ and the $\det(\Amtrx)$ using the $\Lmtrx\,\Umtrx$ factors. 3. Using your codes (Assignment 1 and 2) solve the linear system of algebraic equations $\Amtrx\,\xvec=\bvec$ for $\bvec$ with all ones as elements (*i.e.* $b_i = 1\ \forall \ i=1,\ldots,m$). Solve the same system with `numpy.linalg.solve` and compute the norm of the solutions difference.def lu_factorization( mtrx ):
return (l_mtrx, u_mtrx)
a_mtrx = np.copy(img_mtrx[:120,:120]) # given matrix A
max(A - LU) = 2.91989e-14
min(diag(U)) = -8.51004e-01 rank(A) = 120 det(A) = -1.09249e-148
''' Solve A x = b for b equal to a vector of ones '''
||x_vec - x_vec_gold|| = 3.26639e-10