{"metadata": {"signature": "sha256:a3986859d2cd9957992a15ea7048e4fa966ee90aea26dca3478c72acf3b5f10c", "name": ""}, "worksheets": [{"metadata": {}, "cells": [{"metadata": {}, "source": ["# Gaussian Elimination with Partial Pivoting"], "cell_type": "markdown"}, {"metadata": {}, "outputs": [], "collapsed": false, "input": ["import numpy as np"], "cell_type": "code", "language": "python", "prompt_number": 1}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 0., -7., 5.],\n", " [-1., -2., -5.],\n", " [ 5., -6., 3.]])"], "output_type": "pyout", "prompt_number": 38}], "collapsed": false, "input": ["n = 3\n", "\n", "np.random.seed(22)\n", "A = np.round(5*np.random.randn(n, n))\n", "A[0,0] = 0\n", "A"], "cell_type": "code", "language": "python", "prompt_number": 38}, {"metadata": {}, "source": ["----------------\n", "Now define a function `row_swap_mat(i, j)` that returns a permutation matrix that swaps row i and j:"], "cell_type": "markdown"}, {"metadata": {}, "outputs": [], "collapsed": false, "input": ["def row_swap_mat(i, j):\n", " P = np.eye(n)\n", " P[i] = 0\n", " P[j] = 0\n", " P[i, j] = 1\n", " P[j, i] = 1\n", " return P"], "cell_type": "code", "language": "python", "prompt_number": 39}, {"metadata": {}, "source": ["What do these matrices look like?"], "cell_type": "markdown"}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 0., 1., 0.],\n", " [ 1., 0., 0.],\n", " [ 0., 0., 1.]])"], "output_type": "pyout", "prompt_number": 40}], "collapsed": false, "input": ["row_swap_mat(0,1)"], "cell_type": "code", "language": "python", "prompt_number": 40}, {"metadata": {}, "source": ["Do they work?"], "cell_type": "markdown"}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[-1., -2., -5.],\n", " [ 0., -7., 5.],\n", " [ 5., -6., 3.]])"], "output_type": "pyout", "prompt_number": 41}], "collapsed": false, "input": ["row_swap_mat(0,1).dot(A)"], "cell_type": "code", "language": "python", "prompt_number": 41}, {"metadata": {}, "source": ["--------------"], "cell_type": "markdown"}, {"metadata": {}, "source": ["`U` is the copy of `A` that we'll modify:"], "cell_type": "markdown"}, {"metadata": {}, "outputs": [], "collapsed": false, "input": ["U = A.copy()"], "cell_type": "code", "language": "python", "prompt_number": 63}, {"metadata": {}, "source": ["Create P1 to swap up the right row: "], "cell_type": "markdown"}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 5., -6., 3.],\n", " [-1., -2., -5.],\n", " [ 0., -7., 5.]])"], "output_type": "pyout", "prompt_number": 64}], "collapsed": false, "input": ["P1 = row_swap_mat(0, 2)\n", "U = P1.dot(U)\n", "U"], "cell_type": "code", "language": "python", "prompt_number": 64}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 1. , 0. , 0. ],\n", " [ 0.2, 1. , 0. ],\n", " [ 0. , 0. , 1. ]])"], "output_type": "pyout", "prompt_number": 65}], "collapsed": false, "input": ["M1 = np.eye(n)\n", "M1[1,0] = -U[1,0]/U[0,0]\n", "M1"], "cell_type": "code", "language": "python", "prompt_number": 65}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 5. , -6. , 3. ],\n", " [ 0. , -3.2, -4.4],\n", " [ 0. , -7. , 5. ]])"], "output_type": "pyout", "prompt_number": 66}], "collapsed": false, "input": ["U = M1.dot(U)\n", "U"], "cell_type": "code", "language": "python", "prompt_number": 66}, {"metadata": {}, "source": ["Create `P2` to swap up the right row:"], "cell_type": "markdown"}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 5. , -6. , 3. ],\n", " [ 0. , -7. , 5. ],\n", " [ 0. , -3.2, -4.4]])"], "output_type": "pyout", "prompt_number": 67}], "collapsed": false, "input": ["P2 = row_swap_mat(1,2)\n", "U = P2.dot(U)\n", "U"], "cell_type": "code", "language": "python", "prompt_number": 67}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 1. , 0. , 0. ],\n", " [ 0. , 1. , 0. ],\n", " [ 0. , -0.45714286, 1. ]])"], "output_type": "pyout", "prompt_number": 68}], "collapsed": false, "input": ["M2 = np.eye(n)\n", "M2[2,1] = -U[2,1]/U[1,1]\n", "M2"], "cell_type": "code", "language": "python", "prompt_number": 68}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 5. , -6. , 3. ],\n", " [ 0. , -7. , 5. ],\n", " [ 0. , 0. , -6.68571429]])"], "output_type": "pyout", "prompt_number": 69}], "collapsed": false, "input": ["U = M2.dot(U)\n", "U"], "cell_type": "code", "language": "python", "prompt_number": 69}, {"metadata": {}, "source": ["So we've built $M_2P_2M_1P_1A=U$."], "cell_type": "markdown"}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 5. , -6. , 3. ],\n", " [ 0. , -7. , 5. ],\n", " [ 0. , 0. , -6.68571429]])"], "output_type": "pyout", "prompt_number": 75}], "collapsed": false, "input": ["M2.dot(P2).dot(M1).dot(P1).dot(A)"], "cell_type": "code", "language": "python", "prompt_number": 75}, {"metadata": {}, "source": ["------------------\n", "Can't simply sort the permutation matrices out of the way:"], "cell_type": "markdown"}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 0. , 0. , 1. ],\n", " [ 1. , 0. , 0. ],\n", " [-0.45714286, 1. , 0.2 ]])"], "output_type": "pyout", "prompt_number": 77}], "collapsed": false, "input": ["M2.dot(P2).dot(M1).dot(P1)"], "cell_type": "code", "language": "python", "prompt_number": 77}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 0. , 0. , 1. ],\n", " [ 1. , 0. , 0.2 ],\n", " [-0.45714286, 1. , -0.09142857]])"], "output_type": "pyout", "prompt_number": 80}], "collapsed": false, "input": ["M2.dot(M1).dot(P2).dot(P1)"], "cell_type": "code", "language": "python", "prompt_number": 80}, {"metadata": {}, "source": ["-----------------\n", "But: We now know the right reordering!"], "cell_type": "markdown"}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 0., 0., 1.],\n", " [ 1., 0., 0.],\n", " [ 0., 1., 0.]])"], "output_type": "pyout", "prompt_number": 87}], "collapsed": false, "input": ["P = P2.dot(P1)\n", "P"], "cell_type": "code", "language": "python", "prompt_number": 87}, {"metadata": {}, "outputs": [{"metadata": {}, "text": ["array([[ 5., -6., 3.],\n", " [ 0., -7., 5.],\n", " [-1., -2., -5.]])"], "output_type": "pyout", "prompt_number": 90}], "collapsed": false, "input": ["PA = P.dot(A)\n", "PA"], "cell_type": "code", "language": "python", "prompt_number": 90}]}], "nbformat": 3, "nbformat_minor": 0}