Pasted from IPython

file1.py

{
 "metadata": {
  "name": ""
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# Exemplos do M\u00e9todo da Elimina\u00e7\u00e3o de Gauss\n",
      "\n",
      "Veremos como funciona o algoritmo nas situa\u00e7\u00f5es mais simples, usando a pivota\u00e7\u00e3o para estabilizar os c\u00e1lculos.\n",
      "\n",
      "## Um sistema simples $Ax=b$.\n",
      "Vamos considerar o exemplo simples:\n",
      "$$ \\begin{pmatrix}3 & 2& 1 \\\\\n",
      "1 & 4 & 2 \\\\\n",
      "6 & 1& 0 \n",
      "\\end{pmatrix}\\begin{pmatrix} x \\\\ y \\\\ z \\end{pmatrix}=\\begin{pmatrix} 1 \\\\ 2 \\\\ 1\\end{pmatrix}\n",
      "$$\n",
      "Para executar as opera\u00e7\u00f5es elementares, consideremos as matrizes $A$ e $b$ agregadas da forma $[A|b]$. Aplica-se o algoritmo \u00e0 matriz\n",
      "$$ \\begin{bmatrix}3 & 2& 1 & 1\\\\\n",
      "1 & 4 & 2& 2 \\\\\n",
      "6 & 1& 0 &1 \n",
      "\\end{bmatrix}$$\n",
      "\n",
      "1. Pivota\u00e7\u00e3o na coluna 1 transforma o sistema em: $$\\begin{bmatrix}6 & 1& 0 & 1\\\\ 1 & 4 & 2& 2 \\\\ 3 & 2& 1 &1  \\end{bmatrix}$$\n",
      "\n",
      "2. A seguir aplicamos o primeiro passo da elimina\u00e7\u00e3o de Gauss. Obtemos a matriz $$\\begin{bmatrix}6 & 1& 0 & 1\\\\ 0 & 3.833 & 2& 1.833 \\\\ 0 & 1.5& 1 &0.5  \\end{bmatrix}$$\n",
      "\n",
      "3. Agora novamente a pivota\u00e7\u00e3o e o segundo passo do m\u00e9todo da elimina\u00e7\u00e3o temos a matriz $$\\begin{bmatrix}6 & 1& 0 & 1\\\\ 0 & 3.833 & 2& 1.833 \\\\ 0 & 0& 0.271 & -0.271 \\end{bmatrix}$$\n",
      "\n",
      "O sistema est\u00e1 escalonado e resolvemos com o algoritmo *BackStep* e a solu\u00e7\u00e3o \u00e9 $(0,1,-1)$.\n",
      "\n"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import numpy as np\n",
      "# A is type array from numpy\n",
      "# remember: A[0.0] is the first element\n",
      "# Elementary row operations:\n",
      "\n",
      "def  TrocaLinha(A,i,j):\n",
      "    \"\"\"Troca as linhas i e j da matriz A\"\"\"\n",
      "    B=A.copy()\n",
      "    buffer = B[i].copy()\n",
      "    B[i] = B[j]\n",
      "    B[j] = buffer\n",
      "    return B\n",
      "\n",
      "\n",
      "def SubsLinha(A,i,k,alfa, j=0):\n",
      "    '''soma a linha i com um multiplo da linha k a partir da coluna j'''\n",
      "    B=A.copy()\n",
      "    L,C=shape(B)\n",
      "    B[i,j:C] = B[i,j:C] + alfa*B[k,j:C]\n",
      "    return B\n",
      "\n",
      "def BackStep(A,b):\n",
      "    \"\"\" Se a matriz for triangular superior resolve Ax=b \n",
      "    A \u00e9 uma matriz triangular.\n",
      "    b \u00e9 uma lista de n\u00fameros reais de mesma dimens\u00e3o de A. \"\"\"\n",
      "    #i Vamos testar se A \u00e9 TS e se a a diagonal n\u00e3o tem zero\n",
      "    # Esta parte consome muito tempo e n\u00e3o deve ser usada depois\n",
      "    #n = len(A)\n",
      "    #for i in range(n):\n",
      "    #    teste1 = A[i,i]==0\n",
      "    #    if i==n:\n",
      "    #        teste2 = False\n",
      "    #    else: teste2 = False in (A[i+1:n,i]==0)\n",
      "    #    if teste1 or teste2:\n",
      "    #        raise ValueError(\" Matriz A n\u00e3o pode ser parametro do BackStep \")\n",
      "    # Esta foi s\u00f3 a parte de teste.\n",
      "    n=shape(A)[0]\n",
      "    x=np.zeros(n)\n",
      "    k=n-1\n",
      "    while k>=0:\n",
      "        x[k] = (1/A[k,k])*(b[k]-sum([A[k,j]*x[j] for j in range(k+1,n)]))\n",
      "        k-=1\n",
      "    return x\n",
      "\n",
      "# step k of gauss elimination:\n",
      "def GaussStep(A,k):\n",
      "    '''O k-esimo passo na eliminacao de gauss'''\n",
      "    L=shape(A)[0] # numero de linhas da matriz A\n",
      "    B=A.copy()\n",
      "    for j in range(k+1,L):\n",
      "        m=-B[j,k]/B[k,k]\n",
      "        B=SubsLinha(B,j,k,m,k)\n",
      "        B[j,k]=-m # Aproveitamos a matriz A pra guardar o multiplicador(ta errado!_)\n",
      "    return B\n",
      "\n",
      "# triangle form, without care and pivoting\n",
      "def TriangleForm(A):\n",
      "    '''Retorna a matriz escalonada,\n",
      "    com os multiplicadores na parte triangular inferior'''\n",
      "    B=A.copy()\n",
      "    L=shape(B)[0] # numero de linhas = numero de incognitas\n",
      "    for k in range(L-1):\n",
      "        p=pivo(A,k)\n",
      "        B=TrocaLinha(A,k,p)\n",
      "        B = GaussStep(A,k)\n",
      "     \n",
      "    return B\n",
      "   \n",
      "def pivo(A,n):\n",
      "    ''' retorna o indice de pivotacao da matriz A na coluna n,\n",
      "    a partir da linha n'''\n",
      "    tamanho = shape(A)\n",
      "    l=tamanho[0] # numero de linhas de A\n",
      "    c= tamanho[1] # numero de colunas\n",
      "    if (l<=n) or (c<=n):\n",
      "        raise ValueError(\"n deve ser menor que dimensao de A\")\n",
      "    p=n # inicio o pivo como n, e vou aumentando\n",
      "    for k in range(n,l):\n",
      "        if abs(A[k,n]) > abs(A[p,n]):\n",
      "            p=k\n",
      "    return p"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 43
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = np.array([[3.,2.,1.,1],[1,4,2,2],[6,1,0,1]])\n",
      "A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 44,
       "text": [
        "array([[ 3.,  2.,  1.,  1.],\n",
        "       [ 1.,  4.,  2.,  2.],\n",
        "       [ 6.,  1.,  0.,  1.]])"
       ]
      }
     ],
     "prompt_number": 44
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A1 = TrocaLinha(A,0,pivo(A,0))\n",
      "A1\n",
      "        "
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 45,
       "text": [
        "array([[ 6.,  1.,  0.,  1.],\n",
        "       [ 1.,  4.,  2.,  2.],\n",
        "       [ 3.,  2.,  1.,  1.]])"
       ]
      }
     ],
     "prompt_number": 45
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A2 = GaussStep(A1,0)\n",
      "A2"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 46,
       "text": [
        "array([[ 6.        ,  1.        ,  0.        ,  1.        ],\n",
        "       [ 0.16666667,  3.83333333,  2.        ,  1.83333333],\n",
        "       [ 0.5       ,  1.5       ,  1.        ,  0.5       ]])"
       ]
      }
     ],
     "prompt_number": 46
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A3 = GaussStep(A2,1)\n",
      "A3"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 47,
       "text": [
        "array([[ 6.        ,  1.        ,  0.        ,  1.        ],\n",
        "       [ 0.16666667,  3.83333333,  2.        ,  1.83333333],\n",
        "       [ 0.5       ,  0.39130435,  0.2173913 , -0.2173913 ]])"
       ]
      }
     ],
     "prompt_number": 47
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A4= A3[:,0:3]\n",
      "A4"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 48,
       "text": [
        "array([[ 6.        ,  1.        ,  0.        ],\n",
        "       [ 0.16666667,  3.83333333,  2.        ],\n",
        "       [ 0.5       ,  0.39130435,  0.2173913 ]])"
       ]
      }
     ],
     "prompt_number": 48
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "B1= A3[:,3]\n",
      "B1"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 49,
       "text": [
        "array([ 1.        ,  1.83333333, -0.2173913 ])"
       ]
      }
     ],
     "prompt_number": 49
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x=BackStep(A4,B1)\n",
      "x"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 50,
       "text": [
        "array([  9.25185854e-17,   1.00000000e+00,  -1.00000000e+00])"
       ]
      }
     ],
     "prompt_number": 50
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}