aula-09-eliminacao-gauss.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Álgebra linear com Python: Eliminação Gaussiana e Condicionamento"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Solução de sistemas lineares"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Métodos adequados para a resolução de sistemas lineares e realizar operações no escopo da Álgebra Linear são encontrados no submódulo `linalg` do `scipy`. Importamos essas funcionalidades com: \n",
    "\n",
    "```python\n",
    "from scipy import linalg\n",
    "```\n",
    "\n",
    "Vamos calcular a solução do sistema linear ${\\bf A}{\\bf x} = {\\bf b}$ com\n",
    "\n",
    "$${\\bf A} = \\begin{bmatrix}\n",
    "4   & -2  & -3   & 6   &  \\\\\n",
    "-6  & 7   & 6.5  & -6  &  \\\\\n",
    "1   & 7.5 & 6.25 & 5.5 &  \\\\\n",
    "-12 & 22  & 15.5 & -1 \n",
    "\\end{bmatrix},\n",
    "\\quad\n",
    "{\\bf b} = \\begin{bmatrix}\n",
    "12   \\\\\n",
    "-6.5 \\\\\n",
    "16   \\\\\n",
    "17 \n",
    "\\end{bmatrix}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Vamos importar os módulos e escrever a matriz ${\\bf A}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[  4.    -2.    -3.     6.  ]\n",
      " [ -6.     7.     6.5   -6.  ]\n",
      " [  1.     7.5    6.25   5.5 ]\n",
      " [-12.    22.    15.5   -1.  ]]\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from scipy import linalg\n",
    "\n",
    "A = np.array([[4,-2,-3,6],[-6,7,6.5,-6],[1,7.5,6.25,5.5],[-12,22,15.5,-1]])\n",
    "print(A)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Agora, vamos escrever o vetor ${\\bf b}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[12.  -6.5 16.  17. ]\n"
     ]
    }
   ],
   "source": [
    "b = np.array([12,-6.5,16,17])\n",
    "print(b)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Podemos checar as dimensões com"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(4, 4)"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# dimensões de A\n",
    "A.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(4,)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# dimensão de b\n",
    "b.shape"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A solução do sistema pode ser obtida através do método `linalg.solve`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 2.   4.  -3.   0.5]\n"
     ]
    }
   ],
   "source": [
    "x = linalg.solve(A,b)\n",
    "print(x)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Inversão de matrizes\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Matematicamente, a solução do sistema anterior é dada por ${\\bf x} = {\\bf A}^{-1}{\\bf b}$. Podemos até invocar a matriz inversa aqui com `linalg.inv(A).dot(b)` e a solução é a mesma que no caso anterior."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 2.   4.  -3.   0.5]\n"
     ]
    }
   ],
   "source": [
    "x2 = linalg.inv(A).dot(b)\n",
    "print(x2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    " Por outro lado, este método é ineficiente. Computacionalmente, a inversão de matrizes **não** é aconselhável. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Verificação da solução \n",
    "\n",
    "Podemos usar o fato de que ${\\bf A}{\\bf A}^{-1}{\\bf b} - {\\bf b} = {\\bf 0}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 5.32907052e-15 -3.19744231e-14 -4.44089210e-14 -1.20792265e-13]\n"
     ]
    }
   ],
   "source": [
    "x3 = A.dot(linalg.inv(A).dot(b)) - b\n",
    "print(x3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note que o vetor é próximo do vetor nulo, mas não identicamente nulo."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Podemos também computar a **norma do resíduo (erro)**: $||{\\bf r}|| = ||{\\bf b} - {\\bf A}{\\bf x}|| =  \\langle {\\bf b} - {\\bf A}{\\bf x}, {\\bf b} - {\\bf A}{\\bf x} \\rangle^{1/2}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.464821375527116e-14"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "r = b - A.dot(x) \n",
    "np.sqrt(r.dot(r))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Como a norma do resíduo é próxima de zero, a solução do sistema linear é assumida como correta."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Eliminação Gaussiana\n",
    "\n",
    "A Eliminação Gaussiana (EG) é um algoritmo utilizado para resolver sistemas de equações lineares ao reduzir a matriz plena associada do sistema a uma matriz triangular. Este processo também é chamado de _escalonamento_. Abaixo, usaremos uma matriz genérica 3x3 para exemplificação."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Passos"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Escrever o sistema linear na forma de _matriz estendida_ usando os coeficientes das variáveis como elementos da matriz e o vetor independente como sendo a última coluna;\n",
    "\n",
    "$$\\begin{array}{c}a_{11}x_1+a_{12}x_2+a_{13}x_3 = b_1 \\quad(L_1)\\\\\n",
    "a_{21}x_1+a_{22}x_2+a_{23}x_3 = b_2 \\quad(L_2)\\\\\n",
    "a_{31}x_1+a_{32}x_2+a_{33}x_3 = b_3 \\quad(L_3)\\\\\n",
    "\\downarrow\\\\\n",
    "\\left[\n",
    "\\begin{array}{ccc|c}\n",
    "\ta_{11}&a_{12}&a_{13}&b_{1}\\\\\n",
    "\ta_{21}&a_{22}&a_{23}&b_{2}\\\\\n",
    "\ta_{31}&a_{32}&a_{33}&b_{3}\n",
    "\\end{array}\n",
    "\\right]\n",
    "\\end{array}$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Realizar operações elementares de combinação linear e permutação entre linhas;\n",
    "    - Multiplicação por escalar:    \n",
    "$$\n",
    "\\begin{array}{c}\n",
    "L_2 \\leftarrow L_2 .w\\ \\Rightarrow\\ \n",
    "\\begin{bmatrix}\n",
    "a_{11}&a_{12}\\\\\n",
    "a_{21}&a_{22}\\\\\n",
    "\\end{bmatrix} \\rightarrow\n",
    "\\begin{bmatrix}\n",
    "a_{11}&a_{12}\\\\\n",
    "w.a_{21}&w.a_{22}\\\\\n",
    "\\end{bmatrix}\n",
    "\\end{array}\n",
    "$$\n",
    "\n",
    "    - Combinação linear:    \n",
    "$$\n",
    "\\begin{array}{c}\n",
    "L_2 \\leftarrow L_2 - L_1.w\\ \\Rightarrow\\ \n",
    "\\begin{bmatrix}\n",
    "a_{11}&a_{12}\\\\\n",
    "a_{21}&a_{22}\\\\\n",
    "\\end{bmatrix} \\rightarrow\n",
    "\\begin{bmatrix}\n",
    "a_{11}&a_{12}\\\\\n",
    "a_{21}-a_{11}.w&a_{22}-a_{12}.w\\\\\n",
    "\\end{bmatrix}\n",
    "\\end{array}\n",
    "$$\n",
    "\n",
    "    - Permutação:\n",
    "    \n",
    "$$\n",
    "\\begin{array}{c}\n",
    "L_2 \\leftarrow L_1\\ e\\ L_1\\leftarrow L_2\\ \\Rightarrow\\ \n",
    "\\begin{bmatrix}\n",
    "a_{11}&a_{12}\\\\\n",
    "a_{21}&a_{22}\\\\\n",
    "\\end{bmatrix} \\rightarrow\n",
    "\\begin{bmatrix}\n",
    "a_{21}&a_{22}\\\\\n",
    "a_{11}&a_{12}\\\\\n",
    "\\end{bmatrix}\n",
    "\\end{array}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Anular todos os elementos na porção triangular inferior da matriz original, isto é, todas as entradas exatamente abaixo das entradas dispostas na diagonal principal;\n",
    "\n",
    "$$\\left[\n",
    "\\begin{array}{ccc|c}\n",
    "\ta_{11}&a_{12}&a_{13}&b_{1}\\\\\n",
    "\ta_{21}&a_{22}&a_{23}&b_{2}\\\\\n",
    "\ta_{31}&a_{32}&a_{33}&b_{3}\n",
    "\\end{array}\n",
    "\\right]\\ \\rightarrow\\ \n",
    "\\left[\n",
    "\\begin{array}{ccc|c}\n",
    "\ta_{11}^{(k)}&a_{12}^{(k)}&a_{13}^{(k)}&b_{1}^{(k)}\\\\\n",
    "\t0&a_{22}^{(k)}&a_{23}^{(k)}&b_{2}^{(k)}\\\\\n",
    "\t0&0&a_{33}^{(k)}&b_{3}^{(k)}\n",
    "\\end{array}\n",
    "\\right]$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* A partir da forma triangular, realizar a substituição regressiva.\n",
    "\n",
    "$$\n",
    "\\left[\n",
    "\\begin{array}{ccc|c}\n",
    "\ta_{11}^{(k)}&a_{12}^{(k)}&a_{13}^{(k)}&b_{1}^{(k)}\\\\\n",
    "\t0&a_{22}^{(k)}&a_{23}^{(k)}&b_{2}^{(k)}\\\\\n",
    "\t0&0&a_{33}^{(k)}&b_{3}^{(k)}\n",
    "\\end{array}\n",
    "\\right]\\ \\rightarrow\\ \n",
    "\\left\\{\n",
    "\\begin{array}{c}\n",
    "x_3 = \\frac{b_3^{(k)}}{a_{33}^{(k)}},\\ a_{33}^{(k)} \\neq 0\\\\\\\\\n",
    "x_2 = \\frac{b_2^{(k)}-a_{23}^{(k)}.x_3}{a_{22}^{(k)}},\\ a_{22}^{(k)} \\neq 0\\\\\\\\\n",
    "x_1 = \\frac{b_1^{(k)}-a_{12}^{(k)}.x_2-a_{13}^{(k)}.x_3}{a_{11}^{(k)}},\\ a_{11}^{(k)} \\neq 0\\\\\\\\\n",
    "\\end{array}\n",
    "\\right.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Vejamos um exemplo numérico de como funciona a Eliminação Gaussiana."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.   1.5 -2. ]\n",
      " [ 2.   1.  -1. ]\n",
      " [ 3.  -1.   2. ]]\n"
     ]
    }
   ],
   "source": [
    "# matriz\n",
    "M = np.array([[1.0,1.5,-2.0],[2.0,1.0,-1.0],[3.0,-1.0,2.0]])\n",
    "print(M)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.   1.5 -2. ]\n",
      " [ 0.  -2.   3. ]\n",
      " [ 3.  -1.   2. ]]\n"
     ]
    }
   ],
   "source": [
    "# zeramento da segunda linha \n",
    "m1 = M[1,0]/M[0,0]\n",
    "M[1,:] += - m1*M[0,:]\n",
    "print(M)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.   1.5 -2. ]\n",
      " [ 0.  -2.   3. ]\n",
      " [ 0.  -5.5  8. ]]\n"
     ]
    }
   ],
   "source": [
    "# zeramento da terceira linha\n",
    "m2 = M[2,0]/M[0,0]\n",
    "M[2,:] += - m2*M[0,:]\n",
    "print(M)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.    1.5  -2.  ]\n",
      " [ 0.   -2.    3.  ]\n",
      " [ 0.    0.   -0.25]]\n"
     ]
    }
   ],
   "source": [
    "# eliminação\n",
    "M = np.array([[1.0,1.5,-2.0],[2.0,1.0,-1.0],[3.0,-1.0,2.0]])\n",
    "b = np.array([-2,3,1])\n",
    "\n",
    "m,n = M.shape\n",
    "for i in range(m):\n",
    "    for j in range(i+1,n):\n",
    "        pivo = M[j,i]/M[i,i]                        \n",
    "        for k in range(m):\n",
    "            M[j,k] += -pivo*M[i,k]\n",
    "\n",
    "print(M)            "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "# função simples para eliminação\n",
    "def eliminacao(M):\n",
    "    m,n = M.shape\n",
    "    for i in range(m):\n",
    "        for j in range(i+1,n):\n",
    "            pivo = M[j,i]/M[i,i]                        \n",
    "            for k in range(m):\n",
    "                M[j,k] += -pivo*M[i,k]\n",
    "    return M"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 4.42389663e-01  6.77549133e-01  6.62310314e-01  8.92758362e-01\n",
      "   8.67074978e-01]\n",
      " [ 0.00000000e+00 -8.94767973e-01 -1.17650434e+00 -1.23381565e+00\n",
      "  -1.54014118e+00]\n",
      " [-1.11022302e-16  0.00000000e+00  1.14215572e-01 -6.53327797e-01\n",
      "  -5.87564401e-02]\n",
      " [-2.58946634e-16  0.00000000e+00  0.00000000e+00 -1.60708761e+00\n",
      "  -2.20279294e-01]\n",
      " [ 3.32024317e-17  0.00000000e+00  0.00000000e+00  0.00000000e+00\n",
      "  -9.80998229e-01]]\n"
     ]
    }
   ],
   "source": [
    "# matriz randômica 5x5\n",
    "M2 = np.random.rand(5,5)\n",
    "print(eliminacao(M2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Condicionamento"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Vamos ver como pequenas \"perturbações\" em matrizes podem provocar mudanças drásticas \n",
    "nas soluções de sistemas. Isto ocorre quando temos um problema *mal condicionado*."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "matriz\n",
      "[[1.  2. ]\n",
      " [1.1 2. ]]\n",
      "vetor\n",
      "[10.  10.4]\n"
     ]
    }
   ],
   "source": [
    "A1 = np.array([[1,2],[1.1,2]])\n",
    "b1 = np.array([10,10.4])\n",
    "print('matriz')\n",
    "print(A1)\n",
    "print('vetor')\n",
    "print(b1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[4. 3.]\n"
     ]
    }
   ],
   "source": [
    "# solução do sistema A1x1 = b1\n",
    "x1 = linalg.solve(A1,b1)\n",
    "print(x1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "matriz\n",
      "[[1.    2.   ]\n",
      " [1.055 2.   ]]\n",
      "vetor\n",
      "[10.  10.4]\n"
     ]
    }
   ],
   "source": [
    "d = 0.045\n",
    "A2 = np.array([[1,2],[1.1-d,2]])\n",
    "b2 = np.array([10,10.4])\n",
    "print('matriz')\n",
    "print(A2)\n",
    "print('vetor')\n",
    "print(b2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[7.27272727 1.36363636]\n"
     ]
    }
   ],
   "source": [
    "# solução do sistema perturbado A2x1 = b2\n",
    "x2 = linalg.solve(A2,b2)\n",
    "print(x2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A solução muda drasticamente aqui! Isto se deve à quase dependência linear em que a matriz se encontra. Ou seja, $\\det({\\bf A}_2) \\approx 0$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "-0.20000000000000018 -0.11000000000000032\n"
     ]
    }
   ],
   "source": [
    "print(linalg.det(A1),linalg.det(A2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "169.28388045827452"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "linalg.norm(A)*linalg.norm(linalg.inv(A))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Método de Gauss-Jordan\n",
    "\n",
    "O método de Gauss-Jordan é uma variação da eliminação de Gauss, em que não apenas as entradas da porção inferior da matriz plena do sistema são anuladas, mas também as entradas da porção superior, resultando em uma matriz diagonal. \n",
    "\n",
    "Além disso, todas as linhas são normalizadas através da sua divisão pelo respectivo elemento pivô. Por exemplo, $a_{11}$ é o elemento pivô da primeira equação, $a_{22}$ da segunda, e assim por diante). A partir daí, a obtenção dos valores das variáveis é imediata. \n",
    "\n",
    "O método é melhor ilustrado no seguinte exemplo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[  3.    -0.1   -0.2    7.85]\n",
      " [  0.1    7.    -0.3  -19.3 ]\n",
      " [  0.3   -0.2   10.    71.4 ]]\n"
     ]
    }
   ],
   "source": [
    "# Matriz aumentada\n",
    "AB = np.array([[3. , -0.1, -0.2, 7.85], [0.1, 7., -0.3, -19.3], [0.3, -0.2, 10., 71.4]])\n",
    "print(AB)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.00000000e+00 -3.33333333e-02 -6.66666667e-02  2.61666667e+00]\n",
      " [ 1.00000000e-01  7.00000000e+00 -3.00000000e-01 -1.93000000e+01]\n",
      " [ 3.00000000e-01 -2.00000000e-01  1.00000000e+01  7.14000000e+01]]\n"
     ]
    }
   ],
   "source": [
    "# Normalização da 1a. linha\n",
    "AB[0,:] /= AB[0,0] # L1 <- L1/a11\n",
    "print(AB)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.00000000e+00 -3.33333333e-02 -6.66666667e-02  2.61666667e+00]\n",
      " [ 0.00000000e+00  7.00333333e+00 -2.93333333e-01 -1.95616667e+01]\n",
      " [ 0.00000000e+00 -1.90000000e-01  1.00200000e+01  7.06150000e+01]]\n"
     ]
    }
   ],
   "source": [
    "# Eliminação de x1 da 2a. e 3a. linhas\n",
    "m1 = AB[1,0]\n",
    "AB[1,:] -= m1*AB[0,:] # L2 <- L2 - m1*L1\n",
    "m2 = AB[2,0]\n",
    "AB[2,:] -= m2*AB[0,:] # L3 <- L3 - m2*L1\n",
    "print(AB)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.00000000e+00 -3.33333333e-02 -6.66666667e-02  2.61666667e+00]\n",
      " [ 0.00000000e+00  1.00000000e+00 -4.18848168e-02 -2.79319372e+00]\n",
      " [ 0.00000000e+00 -1.90000000e-01  1.00200000e+01  7.06150000e+01]]\n"
     ]
    }
   ],
   "source": [
    "# Normalização da 2a. linha\n",
    "AB[1,:] /= AB[1,1] # L2 <- L2/a22\n",
    "print(AB)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.00000000e+00  0.00000000e+00 -6.80628272e-02  2.52356021e+00]\n",
      " [ 0.00000000e+00  1.00000000e+00 -4.18848168e-02 -2.79319372e+00]\n",
      " [ 0.00000000e+00  0.00000000e+00  1.00120419e+01  7.00842932e+01]]\n"
     ]
    }
   ],
   "source": [
    "# Eliminação de x2 da 1a. e 3a. linhas\n",
    "m3 = AB[0,1]\n",
    "AB[0,:] -= m3*AB[1,:] # L1 <- L1 - m3*L2  \n",
    "m4 = AB[2,1]\n",
    "AB[2,:] -= m4*AB[1,:] # L3 <- L3 - m4*L2 \n",
    "print(AB)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.          0.         -0.06806283  2.52356021]\n",
      " [ 0.          1.         -0.04188482 -2.79319372]\n",
      " [ 0.          0.          1.          7.        ]]\n"
     ]
    }
   ],
   "source": [
    "# Normalização da 3a. linha\n",
    "AB[2,:] /= AB[2,2] # L3 <- L3/a33\n",
    "print(AB)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.   0.   0.   3. ]\n",
      " [ 0.   1.   0.  -2.5]\n",
      " [ 0.   0.   1.   7. ]]\n"
     ]
    }
   ],
   "source": [
    "# Eliminação de x3 da 1a. e 2a. linhas\n",
    "m5 = AB[0,2]\n",
    "AB[0,:] -= m5*AB[2,:] # L1 <- L1 - m5*L3\n",
    "m6 = AB[1,2]\n",
    "AB[1,:] -= m6*AB[2,:] # L2 <- L2 - m5*L3\n",
    "print(AB)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Do último resultado, vemos que a matriz identidade é obtida, apontando para o vetor solução $[3 \\ \\ -2.5 \\ \\ 7]^T$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Notas\n",
    "\n",
    "Para aqueles acostumados com a notação para matrizes do Matlab, o método `np.mat` pode ajudar. No exemplo a seguir, as linhas da matriz são passadas como uma expressão do tipo `str` e separadas com `;`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1, 2],\n",
       "       [3, 4]])"
      ]
     },
     "execution_count": 29,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# cria matriz\n",
    "np.array(np.mat('1 2; 3 4'))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[1, 2, 3]])"
      ]
     },
     "execution_count": 30,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "np.array(np.mat('1 2 3'))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Tarefa\n",
    "\n",
    "Implemente o algoritmo pleno para a eliminação gaussiana. Verifique exceções (erros que devem ser observados para evitar falhas no algoritmo, tal como pivôs nulos e matrizes singulares, por exemplo) e use o seu método para resolver os sistemas lineares da forma ${\\bf A}{\\bf x} = {\\bf b}$ da lista de exercícios."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.8.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}