{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"using SymPy"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"E (generic function with 2 methods)"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Create an n x n (default 3), i,j elimination matrix\n",
"function E(i,j,n=3)\n",
" E = sympy\"eye\"(n) # start with the identity\n",
" E[i,j] = - symbols(\"m$i$j\") # insert the negative multiplier\n",
" E\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}A_{11}&A_{12}&A_{13}\\\\A_{21}&A_{22}&A_{23}\\\\A_{31}&A_{32}&A_{33}\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" A11 A12 A13\n",
" A21 A22 A23\n",
" A31 A32 A33"
]
},
"execution_count": 3,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"# Create a symbolic 3x3 A matrix\n",
" A = [symbols(\"A$i$j\") for i=1:3, j=1:3]"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0\\\\- m_{21}&1&0\\\\0&0&1\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" 1 0 0\n",
" -m21 1 0\n",
" 0 0 1"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(2,1) "
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"\n",
" \n",
" \n",
"
"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"using Interact"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "ca73636d-e2de-40e3-bf33-af0aefec0cb4",
"version_major": 2,
"version_minor": 0
}
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [],
"text/plain": [
"Interact.Options{:SelectionSlider,Any}(1: \"input\" = 2 Any , \"i\", 2, \"2\", 2, Interact.OptionDict(DataStructures.OrderedDict{Any,Any}(\"1\"=>1,\"2\"=>2,\"3\"=>3), Dict{Any,Any}(Pair{Any,Any}(2, \"2\"),Pair{Any,Any}(3, \"3\"),Pair{Any,Any}(1, \"1\"))), Any[], Any[], true, \"horizontal\")"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"application/vnd.jupyter.widget-view+json": {
"model_id": "c223c5da-a42c-4d7f-be1d-acea2be9152c",
"version_major": 2,
"version_minor": 0
}
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [],
"text/plain": [
"Interact.Options{:SelectionSlider,Any}(3: \"input-2\" = 2 Any , \"j\", 2, \"2\", 2, Interact.OptionDict(DataStructures.OrderedDict{Any,Any}(\"1\"=>1,\"2\"=>2,\"3\"=>3), Dict{Any,Any}(Pair{Any,Any}(2, \"2\"),Pair{Any,Any}(3, \"3\"),Pair{Any,Any}(1, \"1\"))), Any[], Any[], true, \"horizontal\")"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0\\\\0&- m_{22}&0\\\\0&0&1\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" 1 0 0\n",
" 0 -m22 0\n",
" 0 0 1"
]
},
"execution_count": 7,
"metadata": {
"comm_id": "e2514bf5-28a5-4eed-90ac-ba5ddfda6c73",
"reactive": true
},
"output_type": "execute_result"
}
],
"source": [
"@manipulate for i=1:3, j=1:3\n",
" E(i,j)\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0\\\\- m_{21}&1&0\\\\0&0&1\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" 1 0 0\n",
" -m21 1 0\n",
" 0 0 1"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(2,1)"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0\\\\m_{21}&1&0\\\\0&0&1\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" 1 0 0\n",
" m21 1 0\n",
" 0 0 1"
]
},
"execution_count": 9,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"inv(E(2,1))"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"subtract m21 times the first row from the second row"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}A_{11}&A_{12}&A_{13}\\\\- A_{11} m_{21} + A_{21}&- A_{12} m_{21} + A_{22}&- A_{13} m_{21} + A_{23}\\\\A_{31}&A_{32}&A_{33}\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" A11 A12 A13\n",
" -A11*m21 + A21 -A12*m21 + A22 -A13*m21 + A23\n",
" A31 A32 A33"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(2,1) * A "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"do the row operation twice"
]
},
{
"cell_type": "code",
"execution_count": 11,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}A_{11}&A_{12}&A_{13}\\\\- 2 A_{11} m_{21} + A_{21}&- 2 A_{12} m_{21} + A_{22}&- 2 A_{13} m_{21} + A_{23}\\\\A_{31}&A_{32}&A_{33}\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" A11 A12 A13\n",
" -2*A11*m21 + A21 -2*A12*m21 + A22 -2*A13*m21 + A23\n",
" A31 A32 A33"
]
},
"execution_count": 11,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(2,1)^2 * A "
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
" subtract m32 times the second row from the third row"
]
},
{
"cell_type": "code",
"execution_count": 12,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}A_{11}&A_{12}&A_{13}\\\\A_{21}&A_{22}&A_{23}\\\\- A_{21} m_{32} + A_{31}&- A_{22} m_{32} + A_{32}&- A_{23} m_{32} + A_{33}\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" A11 A12 A13\n",
" A21 A22 A23\n",
" -A21*m32 + A31 -A22*m32 + A32 -A23*m32 + A33"
]
},
"execution_count": 12,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(3,2) * A"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Note that after computing E(3,2) * A, the first row is untouched so we can apply E(2,1) without any interference from row 2"
]
},
{
"cell_type": "code",
"execution_count": 13,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}A_{11}&A_{12}&A_{13}\\\\- A_{11} m_{21} + A_{21}&- A_{12} m_{21} + A_{22}&- A_{13} m_{21} + A_{23}\\\\- A_{21} m_{32} + A_{31}&- A_{22} m_{32} + A_{32}&- A_{23} m_{32} + A_{33}\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" A11 A12 A13\n",
" -A11*m21 + A21 -A12*m21 + A22 -A13*m21 + A23\n",
" -A21*m32 + A31 -A22*m32 + A32 -A23*m32 + A33"
]
},
"execution_count": 13,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(2,1) * E(3,2) * A"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"However, in the below, row 2 has changed"
]
},
{
"cell_type": "code",
"execution_count": 14,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}A_{11}&A_{12}&A_{13}\\\\- A_{11} m_{21} + A_{21}&- A_{12} m_{21} + A_{22}&- A_{13} m_{21} + A_{23}\\\\A_{31}&A_{32}&A_{33}\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" A11 A12 A13\n",
" -A11*m21 + A21 -A12*m21 + A22 -A13*m21 + A23\n",
" A31 A32 A33"
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(2,1) * A"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"so if apply E(3,2)
\n",
"(meaning: subtract m32 times the second row from the third row)
\n",
"it will happen with the UPDATED row 2"
]
},
{
"cell_type": "code",
"execution_count": 15,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}A_{11}&A_{12}&A_{13}\\\\- A_{11} m_{21} + A_{21}&- A_{12} m_{21} + A_{22}&- A_{13} m_{21} + A_{23}\\\\A_{11} m_{21} m_{32} - A_{21} m_{32} + A_{31}&A_{12} m_{21} m_{32} - A_{22} m_{32} + A_{32}&A_{13} m_{21} m_{32} - A_{23} m_{32} + A_{33}\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" A11 … A13\n",
" -A11*m21 + A21 -A13*m21 + A23\n",
" A11*m21*m32 - A21*m32 + A31 A13*m21*m32 - A23*m32 + A33"
]
},
"execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(3,2) * E(2,1) * A"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"The row interpretation of matrix muliply correctly acounts for m32 times the second row and m21 x m32 times the first row -- the combined effect"
]
},
{
"cell_type": "code",
"execution_count": 16,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0\\\\- m_{21}&1&0\\\\m_{21} m_{32}&- m_{32}&1\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" 1 0 0\n",
" -m21 1 0\n",
" m21*m32 -m32 1"
]
},
"execution_count": 16,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(3,2) * E(2,1)"
]
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0\\\\- m_{21}&1&0\\\\0&- m_{32}&1\\end{bmatrix}"
],
"text/plain": [
"3×3 Array{SymPy.Sym,2}:\n",
" 1 0 0\n",
" -m21 1 0\n",
" 0 -m32 1"
]
},
"execution_count": 17,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(2,1) * E(3,2)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's move to 5x5 matrices"
]
},
{
"cell_type": "code",
"execution_count": 18,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\0&0&1&0&0\\\\0&0&0&1&0\\\\0&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" -m21 1 0 0 0\n",
" 0 0 1 0 0\n",
" 0 0 0 1 0\n",
" 0 0 0 0 1"
]
},
"execution_count": 18,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(2,1,5)"
]
},
{
"cell_type": "code",
"execution_count": 19,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\- m_{31}&0&1&0&0\\\\- m_{41}&0&0&1&0\\\\- m_{51}&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" -m21 1 0 0 0\n",
" -m31 0 1 0 0\n",
" -m41 0 0 1 0\n",
" -m51 0 0 0 1"
]
},
"execution_count": 19,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E1 = E(2,1,5) * E(3,1,5) * E(4,1,5)* E(5,1,5)"
]
},
{
"cell_type": "code",
"execution_count": 20,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\- m_{31}&0&1&0&0\\\\- m_{41}&0&0&1&0\\\\- m_{51}&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" -m21 1 0 0 0\n",
" -m31 0 1 0 0\n",
" -m41 0 0 1 0\n",
" -m51 0 0 0 1"
]
},
"execution_count": 20,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(3,1,5) * E(2,1,5) * E(5,1,5)* E(4,1,5) # Why does the order not matter"
]
},
{
"cell_type": "code",
"execution_count": 23,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\- m_{31}&0&1&0&0\\\\- m_{41}&0&0&1&0\\\\- m_{51}&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" -m21 1 0 0 0\n",
" -m31 0 1 0 0\n",
" -m41 0 0 1 0\n",
" -m51 0 0 0 1"
]
},
"execution_count": 23,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E1 # subtracts multiples of the first row"
]
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\m_{21}&1&0&0&0\\\\m_{31}&0&1&0&0\\\\m_{41}&0&0&1&0\\\\m_{51}&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" m21 1 0 0 0\n",
" m31 0 1 0 0\n",
" m41 0 0 1 0\n",
" m51 0 0 0 1"
]
},
"execution_count": 24,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"inv(E1) # INVERSE means add back the same multiples of the first row"
]
},
{
"cell_type": "code",
"execution_count": 25,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\0&1&0&0&0\\\\0&- m_{32}&1&0&0\\\\0&0&0&1&0\\\\0&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" 0 1 0 0 0\n",
" 0 -m32 1 0 0\n",
" 0 0 0 1 0\n",
" 0 0 0 0 1"
]
},
"execution_count": 25,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E(3,2,5)"
]
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\0&1&0&0&0\\\\0&m_{32}&1&0&0\\\\0&0&0&1&0\\\\0&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" 0 1 0 0 0\n",
" 0 m32 1 0 0\n",
" 0 0 0 1 0\n",
" 0 0 0 0 1"
]
},
"execution_count": 26,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"inv(E(3,2,5))"
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\- m_{31}&- m_{32}&1&0&0\\\\- m_{41}&0&0&1&0\\\\- m_{51}&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" -m21 1 0 0 0\n",
" -m31 -m32 1 0 0\n",
" -m41 0 0 1 0\n",
" -m51 0 0 0 1"
]
},
"execution_count": 27,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E1 * E(3,2,5)"
]
},
{
"cell_type": "code",
"execution_count": 28,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\m_{21}&1&0&0&0\\\\m_{31}&m_{32}&1&0&0\\\\m_{41}&0&0&1&0\\\\m_{51}&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" m21 1 0 0 0\n",
" m31 m32 1 0 0\n",
" m41 0 0 1 0\n",
" m51 0 0 0 1"
]
},
"execution_count": 28,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"inv(E1) * inv(E(3,2,5))"
]
},
{
"cell_type": "code",
"execution_count": 21,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\- m_{31}&- m_{32}&1&0&0\\\\- m_{41}&0&0&1&0\\\\- m_{51}&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" -m21 1 0 0 0\n",
" -m31 -m32 1 0 0\n",
" -m41 0 0 1 0\n",
" -m51 0 0 0 1"
]
},
"execution_count": 21,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E1 * E(3,2,5) # Why does this have a simple looking answer?"
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\m_{21}&1&0&0&0\\\\m_{31}&m_{32}&1&0&0\\\\m_{41}&0&0&1&0\\\\m_{51}&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" m21 1 0 0 0\n",
" m31 m32 1 0 0\n",
" m41 0 0 1 0\n",
" m51 0 0 0 1"
]
},
"execution_count": 22,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"inv(E1) * inv(E(3,2,5)) # Why does this have an even slightly simpler looking answer?"
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\m_{21} m_{32} - m_{31}&- m_{32}&1&0&0\\\\- m_{41}&0&0&1&0\\\\- m_{51}&0&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" -m21 1 0 0 0\n",
" m21*m32 - m31 -m32 1 0 0\n",
" -m41 0 0 1 0\n",
" -m51 0 0 0 1"
]
},
"execution_count": 29,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
" E(3,2,5) * E1 # Why doesn't this have a simple looking answer?"
]
},
{
"cell_type": "code",
"execution_count": 30,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\0&1&0&0&0\\\\0&- m_{32}&1&0&0\\\\0&- m_{42}&0&1&0\\\\0&- m_{52}&0&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" 0 1 0 0 0\n",
" 0 -m32 1 0 0\n",
" 0 -m42 0 1 0\n",
" 0 -m52 0 0 1"
]
},
"execution_count": 30,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E2 = prod( E(i,2,5) for i=3:5)"
]
},
{
"cell_type": "code",
"execution_count": 31,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\0&1&0&0&0\\\\0&0&1&0&0\\\\0&0&- m_{43}&1&0\\\\0&0&- m_{53}&0&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" 0 1 0 0 0\n",
" 0 0 1 0 0\n",
" 0 0 -m43 1 0\n",
" 0 0 -m53 0 1"
]
},
"execution_count": 31,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E3 = prod( E(i,3,5) for i=4:5)"
]
},
{
"cell_type": "code",
"execution_count": 32,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\0&1&0&0&0\\\\0&0&1&0&0\\\\0&0&0&1&0\\\\0&0&0&- m_{54}&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" 0 1 0 0 0\n",
" 0 0 1 0 0\n",
" 0 0 0 1 0\n",
" 0 0 0 -m54 1"
]
},
"execution_count": 32,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E4 = E(5,4,5)"
]
},
{
"cell_type": "code",
"execution_count": 33,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\- m_{31}&- m_{32}&1&0&0\\\\- m_{41}&- m_{42}&- m_{43}&1&0\\\\- m_{51}&- m_{52}&- m_{53}&- m_{54}&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" -m21 1 0 0 0\n",
" -m31 -m32 1 0 0\n",
" -m41 -m42 -m43 1 0\n",
" -m51 -m52 -m53 -m54 1"
]
},
"execution_count": 33,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E1 * E2 * E3 * E4 # Why is this simple?"
]
},
{
"cell_type": "code",
"execution_count": 34,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\m_{21}&1&0&0&0\\\\m_{31}&m_{32}&1&0&0\\\\m_{41}&m_{42}&m_{43}&1&0\\\\m_{51}&m_{52}&m_{53}&m_{54}&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" m21 1 0 0 0\n",
" m31 m32 1 0 0\n",
" m41 m42 m43 1 0\n",
" m51 m52 m53 m54 1"
]
},
"execution_count": 34,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"L = inv(E1) * inv(E2) * inv(E3) * inv(E4) # Why is this simple? This is the L Matrix!!"
]
},
{
"cell_type": "code",
"execution_count": 35,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\m_{21} m_{32} - m_{31}&- m_{32}&1&0&0\\\\- m_{21} \\left(m_{32} m_{43} - m_{42}\\right) + m_{31} m_{43} - m_{41}&m_{32} m_{43} - m_{42}&- m_{43}&1&0\\\\- m_{21} \\left(- m_{32} \\left(m_{43} m_{54} - m_{53}\\right) + m_{42} m_{54} - m_{52}\\right) - m_{31} \\left(m_{43} m_{54} - m_{53}\\right) + m_{41} m_{54} - m_{51}&- m_{32} \\left(m_{43} m_{54} - m_{53}\\right) + m_{42} m_{54} - m_{52}&m_{43} m_{54} - m_{53}&- m_{54}&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 … 0 0 0\n",
" -m21 0 0 0\n",
" m21*m32 - m31 1 0 0\n",
" -m21*(m32*m43 - m42) + m31*m43 - m41 -m43 1 0\n",
" -m21*(-m32*(m43*m54 - m53) + m42*m54 - m52) - m31*(m43*m54 - m53) + m41*m54 - m51 m43*m54 - m53 -m54 1"
]
},
"execution_count": 35,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"E4 * E3 * E2 * E1 # Why is this a mess? Good thing we don't need it in Gaussian Elimination"
]
},
{
"cell_type": "code",
"execution_count": 36,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\- m_{21}&1&0&0&0\\\\m_{21} m_{32} - m_{31}&- m_{32}&1&0&0\\\\m_{21} m_{42} - m_{41} - m_{43} \\left(m_{21} m_{32} - m_{31}\\right)&m_{32} m_{43} - m_{42}&- m_{43}&1&0\\\\m_{21} m_{52} - m_{51} - m_{53} \\left(m_{21} m_{32} - m_{31}\\right) - m_{54} \\left(m_{21} m_{42} - m_{41} - m_{43} \\left(m_{21} m_{32} - m_{31}\\right)\\right)&m_{32} m_{53} - m_{52} - m_{54} \\left(m_{32} m_{43} - m_{42}\\right)&m_{43} m_{54} - m_{53}&- m_{54}&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 … 0 0 0\n",
" -m21 0 0 0\n",
" m21*m32 - m31 1 0 0\n",
" m21*m42 - m41 - m43*(m21*m32 - m31) -m43 1 0\n",
" m21*m52 - m51 - m53*(m21*m32 - m31) - m54*(m21*m42 - m41 - m43*(m21*m32 - m31)) m43*m54 - m53 -m54 1"
]
},
"execution_count": 36,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"inv(L) # Right this is the inv(L) , not how the inverse of a triangular matrix isn't so pretty in general"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## How do we solve Ly=b for the unknown y given b?"
]
},
{
"cell_type": "code",
"execution_count": 37,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}1&0&0&0&0\\\\m_{21}&1&0&0&0\\\\m_{31}&m_{32}&1&0&0\\\\m_{41}&m_{42}&m_{43}&1&0\\\\m_{51}&m_{52}&m_{53}&m_{54}&1\\end{bmatrix}"
],
"text/plain": [
"5×5 Array{SymPy.Sym,2}:\n",
" 1 0 0 0 0\n",
" m21 1 0 0 0\n",
" m31 m32 1 0 0\n",
" m41 m42 m43 1 0\n",
" m51 m52 m53 m54 1"
]
},
"execution_count": 37,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"L"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Notice that the exact answer gets messy:\n"
]
},
{
"cell_type": "code",
"execution_count": 38,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"b (generic function with 1 method)"
]
},
"execution_count": 38,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"b(i) = symbols(\"b$i\")"
]
},
{
"cell_type": "code",
"execution_count": 39,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}b_{1}\\\\b_{2}\\\\b_{3}\\\\b_{4}\\\\b_{5}\\end{bmatrix}"
],
"text/plain": [
"5-element Array{SymPy.Sym,1}:\n",
" b1\n",
" b2\n",
" b3\n",
" b4\n",
" b5"
]
},
"execution_count": 39,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"[b(i) for i=1:5]"
]
},
{
"cell_type": "code",
"execution_count": 40,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}b_{1}\\\\- b_{1} m_{21} + b_{2}\\\\- b_{1} m_{31} + b_{3} - m_{32} \\left(- b_{1} m_{21} + b_{2}\\right)\\\\- b_{1} m_{41} + b_{4} - m_{42} \\left(- b_{1} m_{21} + b_{2}\\right) - m_{43} \\left(- b_{1} m_{31} + b_{3} - m_{32} \\left(- b_{1} m_{21} + b_{2}\\right)\\right)\\\\- b_{1} m_{51} + b_{5} - m_{52} \\left(- b_{1} m_{21} + b_{2}\\right) - m_{53} \\left(- b_{1} m_{31} + b_{3} - m_{32} \\left(- b_{1} m_{21} + b_{2}\\right)\\right) - m_{54} \\left(- b_{1} m_{41} + b_{4} - m_{42} \\left(- b_{1} m_{21} + b_{2}\\right) - m_{43} \\left(- b_{1} m_{31} + b_{3} - m_{32} \\left(- b_{1} m_{21} + b_{2}\\right)\\right)\\right)\\end{bmatrix}"
],
"text/plain": [
"5-element Array{SymPy.Sym,1}:\n",
" b1\n",
" -b1*m21 + b2\n",
" -b1*m31 + b3 - m32*(-b1*m21 + b2)\n",
" -b1*m41 + b4 - m42*(-b1*m21 + b2) - m43*(-b1*m31 + b3 - m32*(-b1*m21 + b2))\n",
" -b1*m51 + b5 - m52*(-b1*m21 + b2) - m53*(-b1*m31 + b3 - m32*(-b1*m21 + b2)) - m54*(-b1*m41 + b4 - m42*(-b1*m21 + b2) - m43*(-b1*m31 + b3 - m32*(-b1*m21 + b2)))"
]
},
"execution_count": 40,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"y = L\\[b(i) for i=1:5]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Forward Substitution \n",
"but if you compute the answer the obvious way, it is really straightforward."
]
},
{
"cell_type": "code",
"execution_count": 41,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"$$b_{1}$$"
],
"text/plain": [
"b1"
]
},
"execution_count": 41,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"y[1] = b(1)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"L[2,1]*y[1] + y[2] = b(2) implies"
]
},
{
"cell_type": "code",
"execution_count": 42,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"$$- b_{1} m_{21} + b_{2}$$"
],
"text/plain": [
"-b1*m21 + b2"
]
},
"execution_count": 42,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"y[2] = b(2) - L[2,1]*y[1]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"L[3,1]*y[1] + L[3,2]*y[2] + y[3] = b(3) implies"
]
},
{
"cell_type": "code",
"execution_count": 43,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"$$- b_{1} m_{31} + b_{3} - m_{32} \\left(- b_{1} m_{21} + b_{2}\\right)$$"
],
"text/plain": [
"-b1*m31 + b3 - m32*(-b1*m21 + b2)"
]
},
"execution_count": 43,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"y[3] = b(3) - L[3,1]*y[1] - L[3,2]*y[2]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Putting this all together"
]
},
{
"cell_type": "code",
"execution_count": 44,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"y = [b(i) for i=1:5]\n",
"for i = 1:5, j=1:i-1 \n",
" y[i] -= L[i,j]*y[j]\n",
"end"
]
},
{
"cell_type": "code",
"execution_count": 45,
"metadata": {},
"outputs": [
{
"data": {
"text/latex": [
"\\begin{bmatrix}b_{1}\\\\- b_{1} m_{21} + b_{2}\\\\- b_{1} m_{31} + b_{3} - m_{32} \\left(- b_{1} m_{21} + b_{2}\\right)\\\\- b_{1} m_{41} + b_{4} - m_{42} \\left(- b_{1} m_{21} + b_{2}\\right) - m_{43} \\left(- b_{1} m_{31} + b_{3} - m_{32} \\left(- b_{1} m_{21} + b_{2}\\right)\\right)\\\\- b_{1} m_{51} + b_{5} - m_{52} \\left(- b_{1} m_{21} + b_{2}\\right) - m_{53} \\left(- b_{1} m_{31} + b_{3} - m_{32} \\left(- b_{1} m_{21} + b_{2}\\right)\\right) - m_{54} \\left(- b_{1} m_{41} + b_{4} - m_{42} \\left(- b_{1} m_{21} + b_{2}\\right) - m_{43} \\left(- b_{1} m_{31} + b_{3} - m_{32} \\left(- b_{1} m_{21} + b_{2}\\right)\\right)\\right)\\end{bmatrix}"
],
"text/plain": [
"5-element Array{SymPy.Sym,1}:\n",
" b1\n",
" -b1*m21 + b2\n",
" -b1*m31 + b3 - m32*(-b1*m21 + b2)\n",
" -b1*m41 + b4 - m42*(-b1*m21 + b2) - m43*(-b1*m31 + b3 - m32*(-b1*m21 + b2))\n",
" -b1*m51 + b5 - m52*(-b1*m21 + b2) - m53*(-b1*m31 + b3 - m32*(-b1*m21 + b2)) - m54*(-b1*m41 + b4 - m42*(-b1*m21 + b2) - m43*(-b1*m31 + b3 - m32*(-b1*m21 + b2)))"
]
},
"execution_count": 45,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"y"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Back Substitution"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"For upper triangular matrices with arbitrary diagonal there is an analagous algorithm known as \"back substitution.\"\n",
"Since the diagonal is not 1, there is one more divide by the diagonal in the obvious alagorithm.\n",
"In the solution to Ax=b, the pivots are on the diagonal and we divide by those."
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 0.6.0",
"language": "julia",
"name": "julia-0.6"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "0.6.0"
}
},
"nbformat": 4,
"nbformat_minor": 2
}