{
"cells": [
{
"cell_type": "markdown",
"metadata": {},
"source": [
"# Stationary Iterative Methods for Linear Systems\n",
"\n",
"Copyright (C) 2020 Andreas Kloeckner\n",
"\n",
"\n",
"MIT License
\n",
"Permission is hereby granted, free of charge, to any person obtaining a copy\n",
"of this software and associated documentation files (the \"Software\"), to deal\n",
"in the Software without restriction, including without limitation the rights\n",
"to use, copy, modify, merge, publish, distribute, sublicense, and/or sell\n",
"copies of the Software, and to permit persons to whom the Software is\n",
"furnished to do so, subject to the following conditions:\n",
"\n",
"The above copyright notice and this permission notice shall be included in\n",
"all copies or substantial portions of the Software.\n",
"\n",
"THE SOFTWARE IS PROVIDED \"AS IS\", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR\n",
"IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,\n",
"FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE\n",
"AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER\n",
"LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,\n",
"OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN\n",
"THE SOFTWARE.\n",
" "
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [],
"source": [
"import numpy as np\n",
"import scipy.linalg as la\n",
"\n",
"import matplotlib.pyplot as pt"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Let's solve $u''=-30x^2$ with $u(0)=1$ and $u(1)=-1$."
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [],
"source": [
"n = 50\n",
"\n",
"mesh = np.linspace(0, 1, n)\n",
"h = mesh[1] - mesh[0]"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Set up the system matrix `A` to carry out centered finite differences\n",
"\n",
"$$\n",
"u''(x)\\approx \\frac{u(x+h) - 2u(x) + u(x-h)}{h^2}.\n",
"$$\n",
"\n",
"Use `np.eye(n, k=...)`. What needs to be in the first and last row?"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [],
"source": [
"A = (np.eye(n, k=1) + -2*np.eye(n) + np.eye(n, k=-1))/h**2\n",
"A[0] = 0\n",
"A[-1] = 0\n",
"A[0,0] = 1\n",
"A[-1,-1] = 1"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next, fix the right hand side:"
]
},
{
"cell_type": "code",
"execution_count": 9,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [],
"source": [
"b = -30*mesh**2\n",
"b[0] = 1\n",
"b[-1] = -1"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Compute a reference solution `x_true` to the linear system:"
]
},
{
"cell_type": "code",
"execution_count": 10,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"[]"
]
},
"execution_count": 10,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAEACAYAAABVtcpZAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAGtVJREFUeJzt3XucFOWV//HPATKAAnIRRe7ooqgkIrhcXIOzaBREUWOM\nklXAHzEQNWo2iUTXvECz64r7SqImRoxhvcS4iMFVQYgm6kQIGl0IispFwnA3qFwUBLkM5/fH0+MM\nQw/TTNd0dXd9369Xvbp6pqbqUMzU6fM89Txl7o6IiCRTo7gDEBGR+CgJiIgkmJKAiEiCKQmIiCSY\nkoCISIIpCYiIJFgkScDMpprZRjN7q5bvn2lmW81sYWq5NYrjiohIdppEtJ+HgJ8Djx5km1fcfURE\nxxMRkQhEUgm4+zxgSx2bWRTHEhGR6OSyT2Cgmf3VzJ4zs5NyeFwREalFVM1BdVkAdHP3HWY2DHga\nOD5HxxYRkVrkJAm4+/Zq63PM7Jdm1tbdN1ffzsw0kZGISD24e72a3KNsDjJqafc3s6OrrfcHrGYC\nqOTuWtyZOHFi7DHky6JzoXOhc3HwJRuRVAJm9jhQCrQzszXARKAEcHf/FfA1M/s2sAfYCVwWxXFF\nRCQ7kSQBd/9GHd+/D7gvimOJiEh0NGI4T5WWlsYdQt7Quaiic1FF5yIalm17UpTMzPMpHhGRQmBm\neB50DIuISIFREhARSTAlARGRBFMSEBFJMCUBEZEEUxIQEUkwJQERkQRTEhARSTAlARGRBFMSEBFJ\nMCUBEZEEUxIQEUkwJQERkQRTEhARSTAlARGRBFMSEBFJMCUBEZEEUxIQEUkwJQERkQRTEhARSTAl\nARGRBFMSEBFJMCUBEZEEUxIQEUkwJQERkQRrEncAUjh274Zt22D79qrX7dthxw747LP0y969UFFx\n4LJvHzRqBI0bh6X6+he+AM2aQfPm+y/NmkGLFtC6ddXSsiWYxX1mRAqXkkDCVFTAli3w0Udh2bRp\n//WtW8OyZcv+r598Au7hotuiRdVrixZw2GHhAp1uadIESkqqLvDVL/r79lUlhOoJYs+ekGR27jxw\n2b4dPv64Ks6dO6FVq5AQjjwSjj4aOnQ48LVLF+jcOcQjIlXM3eOO4XNm5vkUT6HYtw8+/BDef79q\n2bgxLB98sP/65s1wxBHhgtmuXXitXNq2hTZtqj5lV18/4ohwMc+3T91794YEVZnY/v738G+t/vr+\n+7B2bXjfsSN0777/csIJcOKJ4d8pUojMDHev11+nkkCe27YN1q+Hdev2f12/HjZsCBe4Dz4In4Y7\ndoRjjgnL0UdXLUcdVbXerl1yPw3v3h2SwapVVUt5OSxdGpZWrUIyOOmkqtdTTw0JUCSfKQkUqM8+\nCxelNWvCa7pl797QjNGpU9XSuXO44FcuHTqET+lSf/v2hQT77ruwZElY3n4b3nornON+/aqWvn2V\nGCS/KAnkqa1bYfXq8Ilz9eqwrFlT9bplS7igd+0a2qzTLa1b518TTJJUVIQq4f/+DxYsCMubb4b/\nty9/GQYPDku3bnFHKkmmJBCTjz8OzQmVzQqVF/vKpoaKinBx6N49vHbrFi74lesdOoQOUiksFRWw\neDHMnQuvvBKWZs2qEsKQIXDccXFHKUmiJNBAduyousBXXyq/tncv9OhxYEdj5UW/bVt9ik8Cd1i+\nPCSDP/0JXnwx9C+cd15YBg+Gpk3jjlKKmZJAPVVUhHbg8nJYuXL/1/Ly0FzTrVu40PfoAcceGy7w\nle91kZd09u2DRYtg9myYMyf0LZSWwrBhMGJE6GMQiVLsScDMpgLnAxvd/Uu1bHMvMAz4FBjj7ovS\nbBN5EtiyperiXn0pLw/t8u3bV13gq7/26BH+WNVcI9natAleeAGeey4kht694bLL4JJLQpOgSLby\nIQmcAWwHHk2XBMxsGHCduw83swHAPe4+MM12h5wE9uwJF/OaF/nKZe/e0D5b/SJfud69e2jLFcmV\nXbtCQpg+HWbNgj59qhJC+/ZxRyeFKvYkkAqiGzCzliQwBXjZ3Z9IvV8ClLr7xhrbHZAE3MOn+eoX\n9r/9rWp9w4bwaaryQl9zaddOTTaSn3buhN//PiSEOXPC3Ubf/CYMH57csRxSP4WQBGYC/+nu81Pv\n/wjc5O4La2zn99/vB3yadw8X9HQX+q5ddY+8FL5PPw3J4Ne/Dk2VY8bA2LG6y0gyk00SyNXnjXTB\npc0+Dz44iTZtwpQF48aVcumlpbRpo0/zUtwOPxyuuios774LU6fCoEGh/+Cb34SvfU0fdqRKWVkZ\nZWVlkewrruagpcCZmTQHiSTVrl3w7LMwZUoYsPad78C4ceEDkkh12VQCUd77YqT/xA/wLDAKwMwG\nAltrJgAR2V/TpnDppWHcwezZYSqL444LyeBvf4s7OikWkSQBM3scmA8cb2ZrzOwqMxtnZt8CcPfZ\nQLmZrQAeAK6J4rgiSXHKKfDII2HMQcuWMHBguKPo1VfjjkwKXaIHi4kUqu3b4eGH4Sc/gZ494fbb\nQ2KQZMqX5iARyZEWLeC662DZstBp/PWvh1tLFyyIOzIpNEoCIgWspAS+9S14770wT9GIEXDRRWGm\nU5FMKAmIFIGmTeHaa2HFijBP0dChYSTyqlVxRyb5TklApIg0bw433hiSQe/ecNppcOutoQ9BJB0l\nAZEidPjh8KMfhdlMV62CXr3gN78JM5yKVKe7g0QS4LXX4IYbwvo99+hOomKju4NE5KAGDgxjCq69\nNtxNNGZMmOJaRElAJCEaNYJRo8IUFK1bhz6DJ54IEzRKcqk5SCShXn01TE533HHwy19C585xRyT1\npeYgETlkgwbBwoXQty+cemqYqE4dx8mjSkBEeOed8PyCpk3DMw169ow7IjkUqgREJCsnnwx//nMY\nbXz66WFeIn0eSwZVAiKyn8WLYeTI0HE8ZUroRJb8pkpARCLzxS/CG2/AkUdCnz4wb17cEUlDUiUg\nIrWaOROuvhrGjw/TTzTJ1QNp5ZDkxYPmo6AkIJJ/NmyA0aNhxw6YNg26dIk7IqlJzUEi0mA6doTn\nnw/TVPfvDxE931zyhCoBEcnYH/4AV1wBt9wC118PVq/PnhI1NQeJSM6Ul8PFF8OXvgQPPBCmr5Z4\nqTlIRHKmRw+YPx8qKuCf/glWr447IsmGkoCIHLLDDoPHHoMrr4QBA+Cll+KOSOpLzUEikpWXXoJv\nfANuvz0871hyT30CIhKrFStg2LDwXOMf/1gdxrmmJCAisfvwQzj//PAoywcfhJKSuCNKDnUMi0js\n2reHl1+GrVth+HD4+OO4I5JMKAmISGQOOwyeeipMRT14MKxfH3dEUhclARGJVOPGcN99obN40CB4\n++24I5KDUZ+AiDSYxx+H734XnnkmPOxeGoY6hkUkb82eDWPGwIwZ8OUvxx1NcVLHsIjkrfPOCxXB\nV78KL74YdzRSk5KAiDS4s88OlcDIkTBnTtzRSHVKAiKSE4MHh76B0aPDq+QHPSdIRHJm0KDQR3D+\n+bB7N1x6adwRiZKAiOTUaaeFh9QMHRoSwb/8S9wRJZuSgIjk3CmnhE7is8+GZs3gkkvijii5IukT\nMLOhZrbUzJab2YQ03x9tZh+Y2cLU8v+iOK6IFK6TTgpNQ9dcEyoDiUfW4wTMrBGwHDgL2AC8AVzu\n7kurbTMa6Ofu19exL40TEEmY+fPhoovCdBNnnBF3NIUp7nEC/YH33H21u+8BpgEXptlOk8uKyAFO\nPx1++9swjmDhwrijSZ4okkAnYG219+tSX6vpq2a2yMymm1nnCI4rIkXiK18JzysePhyWLIk7mmSJ\nIgmk+4Rfs03nWaC7u/cBXgQeieC4IlJELr4YJk+Gc8+FVavijiY5org7aB3Qtdr7zoS+gc+5+5Zq\nbx8EJte2s0mTJn2+XlpaSmlpaQQhikghGDUKtm0Ldw3NnQvHHBN3RPmprKyMsrKySPYVRcdwY2AZ\noWP4feB1YKS7L6m2TQd3/3tq/WLgB+5+epp9qWNYRPiP/4AnnwyJoGXLuKPJf7HPImpmQ4F7CM1L\nU939TjO7DXjD3WeZ2R3ACGAPsBn4trsvT7MfJQERwR3Gj4e1a+HZZ6GJRjQdVOxJICpKAiJSae/e\nML3EsceGh9To4fW1i/sWURGRyDVpAtOnw7x58NOfxh1N8VKRJSJ5q1UreO65MPFcjx5hLIFES81B\nIpL3Fi4ME87NmgX9+8cdTf5Rc5CIFLW+feG//ztML1FeHnc0xUVJQEQKwvnnw803h1HFW7bUvb1k\nRs1BIlJQbrgBli8PTUONG8cdTX5Qc5CIJMZPfgKffQbVJheQLCgJiEhBadIEnngCHnkEnn467mgK\nn5qDRKQgvf566Cd45RXo1SvuaOKl5iARSZz+/eGOO8Lso9u2xR1N4VIlICIFbdw4+Ogj+N3vkju1\nhCoBEUmse++F9evDswjk0KkSEJGCt25daB56+GE455y4o8k9VQIikmidO8O0aXDllXoq2aFSEhCR\nojB4MPzgBzByJOzZE3c0hUNJQESKxr/+KxxxhAaSHQr1CYhIUdm4EU49FR57DIYMiTua3FCfgIhI\nytFHhw7iUaPgww/jjib/qRIQkaI0YQK88w7MnFn84wdUCYiI1PDv/x4qgXvvjTuS/KZKQESK1sqV\nMGAAvPBC6CcoVqoERETSOPbYUAlcfjls3x53NPlJlYCIFL2rrgqvDz0UbxwNRZWAiMhB/Pzn8Oc/\nw4wZcUeSf1QJiEgizJ8Pl1wCb70F7dvHHU20sqkElAREJDFuuinMLTR9etyRREvNQSIiGbj9dli8\nuPiSQDZUCYhIorz+OowYAW++GUYXFwM1B4mIHIJbboF334X//d/iGE2s5iARkUMwcSKsWAGPPx53\nJPFTJSAiibRgAQwbBosWQceOcUeTHVUCIiKHqF8/GD8+PKg+yZ89lQREJLFuvRXWrIFHH407kvio\nOUhEEm3RIvjKV8Ktox06xB1N/ejuIBGRLNx8cxhE9j//E3ck9aMkICKShR07oHdvuP9+OPfcuKM5\ndLF3DJvZUDNbambLzWxCmu+XmNk0M3vPzF41s65RHFdEJAqHHQb33QfXXAM7d8YdTW5lnQTMrBHw\nC+Bc4GRgpJn1qrHZWGCzu/cE7gbuyva4IiJRGjYMTjstPJEsSaKoBPoD77n7anffA0wDLqyxzYXA\nI6n13wFnRXBcEZFI3X03/OpX4dnESRFFEugErK32fl3qa2m3cfcKYKuZtY3g2CIikTnmGJg0KYwf\n2Lcv7mhyI4okkK4zombvbs1tLM02IiKxGz8edu0q3qeQ1dQkgn2sA6p39HYGNtTYZi3QBdhgZo2B\nVu6+Jd3OJk2a9Pl6aWkppaWlEYQoIpKZxo1Dk9A558AFF8BRR8Ud0YHKysooKyuLZF9Z3yKauqgv\nI7Tzvw+8Dox09yXVtrkG6O3u15jZ5cBF7n55mn3pFlERyQvf+x58+GFhjCaOfZyAmQ0F7iE0L011\n9zvN7DbgDXefZWZNgd8ApwKbgMvdfVWa/SgJiEhe2L4dTjoJHn4YhgyJO5qDiz0JREVJQETyyTPP\nwIQJ4bnEJSVxR1O72AeLiYgUoxEjoHv3MJCsWKkSEBE5iCVLYPDg8CSy9u3jjiY9NQeJiDSgG2+E\nzz6DKVPijiQ9JQERkQa0ZQv06gXPPw99+sQdzYHUJyAi0oDatAkjiW+8sfieQqYkICKSgauvhs2b\n4amn4o4kWmoOEhHJ0EsvwdixoZO4efO4o6mi5iARkRwYMgT69oWf/jTuSKKjSkBE5BCsXAn9+8Ob\nb0KnmvMlx0R3B4mI5NAtt8C6dfkzr5CSgIhIDm3bFm4ZnTEDBg6MOxr1CYiI5FTLlnDHHWGm0UL/\n3KokICJSD1dcAZ98AjNnxh1JdpQERETqoXFjuPNOuPlmqKiIO5r6UxIQEamn886Ddu3yp4O4PtQx\nLCKShVdfhcsug2XL4htApo5hEZGYDBoE/foV7jMHVAmIiGRpyRI480xYvhxat8798VUJiIjE6MQT\nw1PI7ror7kgOnSoBEZEIrFsHp5wCixdDx465PbZGDIuI5IEJE2DrVnjggdweV0lARCQPbNkCxx8P\n8+bBCSfk7rjqExARyQNt2sD3vw//9m9xR5I5VQIiIhHauTNUAzNmhCmnc0GVgIhInmjeHH70I5g4\nMe5IMqMkICISsTFjwiMoX3st7kjqpiQgIhKxkpLw4Jnbbos7krqpT0BEpAHs3g09e8ITTzT8g2fU\nJyAikmcKpRpQJSAi0kByVQ2oEhARyUOFUA2oEhARaUC5qAZUCYiI5KmSkvAIynytBlQJiIg0sF27\nQjUwfXrDVAOqBERE8ljTpvnbN6BKQEQkBxqyGoitEjCzNmb2gpktM7PnzeyIWrarMLOFZvZXM3s6\nm2OKiBSifK0GsqoEzGwysMnd7zKzCUAbd/9hmu0+cfdWGexPlYCIFK3KauDJJ2HAgOj2G9tDZcxs\nKXCmu280sw5Ambv3SrPdNndvmcH+lAREpKjdey+UlcFTT0W3zziTwGZ3b1vt/SZ3b5dmu93AImAv\nMNndn6llf0oCIlLUPv0UevSAuXOje/pYNkmgSQY7/wNwdPUvAQ7cegjH6erufzezHsBLZvaWu5en\n23DSpEmfr5eWllJaWnoIhxERyW+HHw7XXgv/9V/w61/Xbx9lZWWUlZVFEk+2lcASoLRac9DL7n5i\nHT/zEDDT3Q8ohlQJiEgSbNoU+gYWL4ZOnbLfX5zjBJ4FxqTWRwMHNPOYWWszK0mtHwmcDryb5XFF\nRApWu3Zw5ZVw991xR5J9JdAWmA50AdYAl7r7VjPrB4xz92+Z2SDgAaCCkHR+5u4P17I/VQIikgir\nV0PfvrBiRXhAfTZi6xiOmpKAiCTJqFHQq1cYP5ANJQERkQL09ttw9tlQXh4eUF9fmjtIRKQA9e4N\n//iP8Mgj8cWgSkBEJEbz5sHo0bBsGTSp86b99FQJiIgUqDPOgA4dYMaMeI6vJCAiErMf/hAmT4Y4\nGkKUBEREYjZ8eJhc7o9/zP2xlQRERGLWqBHcdBPceWcMx879IUVEpKaRI2H5cvjrX3N7XCUBEZE8\nUFIC110H99yT2+PqFlERkTyxeTMcdxwsWRLuGMqUbhEVESkCbdvCZZfBlCm5O6YqARGRPPLuuzBk\nSJhgrmnTzH5GlYCISJE46SQ45RSYNi03x1MSEBHJMzfcEDqIc9EwoiQgIpJnhg4NzyKeO7fhj6Uk\nICKSZxo1guuvz83touoYFhHJQ9u3Q7dusGABdO9+8G3VMSwiUmRatICrroJf/KJhj6NKQEQkT61a\nBf36hdtFW7SofTtVAiIiRah7dygtbdgnj6kSEBHJY3PnwtixsHRp6DBOR5WAiEiROuOM0BT0+983\nzP6VBERE8pgZ3Hgj3H13A+0/n5pf1BwkInKgXbuga9fQNHT88Qd+X81BIiJFrGnTcLvoAw9Ev29V\nAiIiBWDlShgwANasgebN9/+eKgERkSJ37LFw2mnw5JPR7ldJQESkQIwfH/0DZ5QEREQKxPDhsHYt\nvPlmdPtUEhARKRBNmsDVV0dbDahjWESkgKxfD1/8YphPqGXL8DV1DIuIJESnTvDP/wy//W00+1MS\nEBEpMJUdxFE0nCgJiIgUmLPOCo+f/Mtfst+XkoCISIFp1AjGjYP7749gX9n8sJl9zczeNrMKM+t7\nkO2GmtlSM1tuZhOyOaaIiMCYMfDMM7B5c3b7ybYSWAxcDPyptg3MrBHwC+Bc4GRgpJn1yvK4Ra+s\nrCzuEPKGzkUVnYsqST8XRx4JF1yQ/QNnskoC7r7M3d8DDnZrUn/gPXdf7e57gGnAhdkcNwmS/gte\nnc5FFZ2LKjoX0YwgzkWfQCdgbbX361JfExGRLJx+ephhNBt1JgEz+4OZvVVtWZx6vSDDY6SrEjQi\nTEQkS2ahGshqH1GM0DWzl4HvufvCNN8bCExy96Gp9z8E3N0np9lWyUFEpB7qO2K4SYQx1BbAG8A/\nmFk34H3gcmBkug3r+48QEZH6yfYW0YvMbC0wEJhlZnNSXz/GzGYBuHsFcB3wAvAOMM3dl2QXtoiI\nRCGvJpATEZHcimXEcF2Dx8ysxMymmdl7ZvaqmXWNI85cyOBcfNfM3jGzRalO+i5xxJkLmQ4qTA1S\n3HewAYqFLpNzYWZfT/1uLDazx3IdY65k8DfSxcxeMrOFqb+TYXHE2dDMbKqZbTSztw6yzb2p6+Yi\nM+uT0Y7dPacLIfGsALoBXwAWAb1qbPNt4Jep9csITUg5jzVPzsWZQLPU+vgkn4vUdi0IgxPnA33j\njjvG34t/ABYArVLvj4w77hjPxQPAuNT6iUB53HE30Lk4A+gDvFXL94cBz6XWBwCvZbLfOCqBTAaP\nXQhUjoP7HXBWDuPLpTrPhbv/yd0/S719jeIdY5HpoMIfA5OBXbkMLscyORdXA/e5+ycA7v5RjmPM\nlUzOxT6gVWq9NbA+h/HljLvPA7YcZJMLgUdT2/4FOMLMjq5rv3EkgUwGj32+jYeO5a1m1jY34eXU\noQ6kGwvMadCI4lPnuUiVt53dfXYuA4tBJr8XxwMnmNk8M5tvZufmLLrcyuRc3AZcmbpJZRbwnRzF\nlm9qnqv1ZPChMcpbRDOVyeCxmttYmm2KQcYD6czsCqAfoXmoGB30XJiZAT8DRtfxM8Ugk9+LJoQm\nocFAV2CumZ1cWRkUkUzOxUjgIXf/WWpc0mOEecqSpl4Dc+OoBNYRfmkrdQY21NhmLdAFwMwaE9o9\nD1YGFapMzgVmdjZwM3BBqiQuRnWdi5aEP+wyMysn3Jb8TJF2Dmfye7EOeMbd97n7KmAZ0DM34eVU\nJudiLDAdwN1fA5qZ2ZG5CS+vrCN13UxJez2pKY4k8PngMTMrIQwee7bGNjOp+sR3KfBSDuPLpTrP\nhZmdCkwBRrj7phhizJWDngt3/8Tdj3L3Y929B6F/5AJPM0q9CGTyN/I0MAQgdcHrCazMaZS5kcm5\nWA2cDWBmJwJNi7iPxKi9An4WGAWfz9Sw1d031rXDnDcHuXuFmVUOHmsETHX3JWZ2G/CGu88CpgK/\nMbP3gE2E//iik+G5uAs4HHgy1SSy2t0vii/qhpHhudjvRyjS5qBMzoW7P29m55jZO8Be4PvFWC1n\n+HvxfeBBM/suoZN4dO17LFxm9jhQCrQzszXARKCEMA3Pr9x9tpmdZ2YrgE+BqzLab+p2IhERSSA9\nXlJEJMGUBEREEkxJQEQkwZQEREQSTElARCTBlARERBJMSUBEJMGUBEREEuz/A3hI5yO7oWfsAAAA\nAElFTkSuQmCC\n",
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"x_true = la.solve(A, b)\n",
"pt.plot(mesh, x_true)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next, we'll try all the stationary iterative methods we have seen."
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Jacobi"
]
},
{
"cell_type": "code",
"execution_count": 44,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [],
"source": [
"x = np.zeros(n)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"Next, apply a Jacobi step:"
]
},
{
"cell_type": "code",
"execution_count": 42,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [],
"source": [
"x_new = np.empty(n)\n",
"\n",
"for i in range(n):\n",
" x_new[i] = b[i]\n",
" for j in range(n):\n",
" if i != j:\n",
" x_new[i] -= A[i,j]*x[j]\n",
" \n",
" x_new[i] = x_new[i] / A[i,i]\n",
"\n",
"x = x_new"
]
},
{
"cell_type": "code",
"execution_count": 43,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [
{
"data": {
"text/plain": [
""
]
},
"execution_count": 43,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAEACAYAAABVtcpZAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xl8VNX9//HXJwlJgCQQ1rCGRVxRligIX0tToSouqHVB\nrCIUFcS9tXX5asHW+tXaRW2xuCJU/SEVFVAQFxpBEVAUkF1lXwVMWMOS5Pz+mCSEmIRkZjJ3lvfz\n8biPuTP3zLkfrvF+5pxzz73mnENERGJTnNcBiIiId5QERERimJKAiEgMUxIQEYlhSgIiIjFMSUBE\nJIYFJQmY2Ytmtt3MllSy/admlmdmXxYvDwZjvyIiEpiEINUzDvgHMKGKMrOdcwOCtD8REQmCoLQE\nnHOfALnHKWbB2JeIiARPKMcEzjazr8zsXTM7NYT7FRGRSgSrO+h4FgKZzrkDZtYfeBs4MUT7FhGR\nSoQkCTjn9pVZn2Fmz5hZI+fcD2XLmZluZCQi4gfnnF9d7sHsDjIq6fc3s+Zl1nsAVj4BlHDOaXGO\nUaNGeR5DuCw6FjoWOhZVL4EISkvAzF4DsoHGZrYBGAUkAs459xxwpZndAhwB8oGBwdiviIgEJihJ\nwDl37XG2jwHGBGNfIiISPJoxHKays7O9DiFs6FgcpWNxlI5FcFig/UnBZGYunOIREYkEZobzc2A4\nVJeIiojUWLt27Vi/fr3XYYSNzMxM1q1bF9Q61RIQkbBV/AvX6zDCRmXHI5CWgMYERERimJKAiEgM\nUxIQEYlhSgIiIjFMSUBExE/t27dn1qxZXocRECUBEZFaUFhY6HUI1aIkICLih8GDB7NhwwYuvvhi\n0tLSeOKJJ4iLi+Oll14iMzOTvn378vHHH9OmTZtjvle29eCc47HHHuOEE06gadOmXHPNNeTl5YX0\n36EkICLihwkTJtC2bVveffdd9uzZw9VXXw3A7NmzWblyJTNnzgR81/BX5qmnnmLq1KnMmTOHLVu2\nkJ6ezsiRI0MSfwnNGBaRiGYPB+fJtW6Uf5PSyk7eMjMefvhh6tatW63vPvfcc4wZM4YWLVoA8Pvf\n/57MzExeeeUV4uJC8xtdSUBEIpq/J+/a0rp162qXXb9+PZdffnnpCd85R506ddi+fXtpYqhtSgIi\nIn6qqKun7Gf169fnwIEDpe8LCwvZsWNH6fu2bdvy0ksv0atXr9oNtAoaExAR8VNGRgZr1qwBqPAp\nXyeeeCIHDx5kxowZFBQU8Mgjj3D48OHS7cOHD+eBBx5gw4YNAOzYsYOpU6eG7h+AkoCIiN/uu+8+\n/vjHP9KoUSMmT578o5ZBWloazzzzDMOGDaN169akpqYe01105513cumll3LeeefRoEEDevfuzYIF\nC0L6b9BdREUkbOkuosfSXURFRCSolARERGKYkoCISAxTEhARiWFKAiIiMUxJQEQkhmnGsIiErczM\nzCpvwBZrMjMzg16n5gmIiEQ4zRMQERG/qDtIqu1w4WH2HtrLvsP72HvY97rv8D4OHDnAwYKDFS4F\nRQUUFhVS6AqPeS1yRcRZHPFx8cRb/DHrdeLrkJyQTN2EutStU7f0NTkhmZTEFBomNyxdUhNT1V0g\nEgB1B8WYwqJCcg/msvPATnYe2MmuA7uOrufvIu9gHnkH88g9mOt7zfe97jm0B4cjNTGVlMQUUpN8\nrymJKdSrU4/khGTfEp98dD0hmYS4hNKTe9nXOIujyBWVJoSySeJI4REOFhwkvyCf/CP5vtfi9X2H\n97H70O7SOPOP5JOWlEbD5IY0qdeE5inNyaif4XtNyaB5fd9rmwZtaJ3WmoQ4/e6R6BNId5CSQBQo\nckXs2L+Drfu2snXvVrbu28r2fdvZvn873+//nu37t7N9n2/9h/wfaJDcgCb1mtC4bmOa1GtSujSq\n24j05PTSX9npdY+uN0hqQGJ8Ytj96i4oKmDPoT3k5vsS27Z929i+f7vvdd92tu3fxta9W9m4ZyPb\n922nZWpL2jVsd8xyUuOTOKXpKTRMbuj1P0fEL0oCUWzvob1s3ruZTXs2sXlP8evezWzeu5kte7ew\nde9Wvt//PWlJabRMbUmL1Ba0SGlB8/rNaZ7SnOb1m9OsfrPS9cb1Gsfsr+HDhYfZuHsj6/LWlS5r\n89aycudKVu5cSVpSGqc0PYVTm5zqe216Kt0yutEguYHXoYtUSUkgQh0sOMjG3RvZsHsDG/dsZOPu\njb7XMusFRQW0TmtNq9RWtEprRavUVrROa03L1JalS0ZKBonxiV7/cyJakSti055NLN+xnBU7VrBi\n5wqWfr+UJduX0DK1JVkts8hq4Vu6t+iuxCBhRUkgTOUdzGN93nrW5a1j/e71rM9bz4Y9G3yvuzeQ\nezCX1mmtadugLW3S2viWBse+NkxuGHZdMLGksKiQlTtX8sWWL1i4dSELty5k8bbFtE5rzU/a/oQ+\nmX3ok9mHzIbBv35bpLqUBDyy++Bu1uat9XUr5K4tPdmXdDUUukIyG2TSrmE7Mhtkktkwk7YN2pau\nZ6RkEGe6SjfSFBYV8vX3XzNn/Rxmb5jN7PWzSU5I9iWEtn04t/25dGzU0eswJYYoCdSSA0cOlJ7g\n1+atLX0t6UsuKCqgfcP2PxpoLDnpN6rbSL/iY4BzjtW7VjN7/Ww+Xv8xH639iLSkNC484UIu7HQh\nfTL7kJSQ5HWYEsWUBPxUWFTIpj2bWJu3ljW5a1ibu5Y1eWtKT/a5+blkNsykfcP2tG/Yng7pHWjX\nsB3t033vdZKXihS5IhZtW8T0b6Yz49sZLP1+Kdntsul/Qn8GnDSAlqktvQ5RooznScDMXgQuBrY7\n586opMzTQH9gPzDEObeogjJBTwK5+bmlJ/myy9q8tWzYvYGm9ZrSPt13gi850bdv2J726e1pmdpS\n3TUSsF0HdvH+d+/z7jfvMv2b6XRu1pmBpw3kilOvICMlw+vwJAqEQxI4B9gHTKgoCZhZf+A259xF\nZtYTeMo5d3YF5WqcBI4UHmHD7g3HnuTzjq4XFBXQMb3jMSf5DukdaJ/u68ZJTkj2818tUnOHCg7x\n/nfvM2n5JN5Z/Q5dM7r6EsIpV9C0flOvw5MI5XkSKA4iE5hWSRIYC/zXOfd68fsVQLZzbnu5cj9K\nAs45cg/mHnOS/+6H70pP9Fv2biEjJaP0RF9+aVy3sbpsJCzlH8nnvW/fY9LyScz4ZgY/yfwJN3a7\nkYtOvChm53KIfyIhCUwD/s85N7f4/YfA75xzX5Yr5/71+b9+1HXjcHRI71Dhib5tg7a6Rl4i3v7D\n+5m0bBIvfPUCa3PXMqTrEIZ1G6arjKRaAkkCofq5UVFwFWafZ594nsb100lPTmf4ucO5avBVpCen\n69e8RLX6ifUZ2m0oQ7sNZfmO5bz45Yv0erEXnZt15sbuN3LlqVfqx46UysnJIScnJyh1edUdtBL4\naUXdQfPmOXr2DEpIIhHtUMEhpq6aytiFY1m5cyW397id4VnDSa+b7nVoEmbC5XkCRsW/+AGmAoMB\nzOxsIK98Aijx+edBjEgkgiUlJHHVaVfx0eCPmH7tdFbsXEHHpzty+/Tb+e6H77wOT6JEUJKAmb0G\nzAVONLMNZjbUzIab2c0AzrnpwFoz+xZ4FhhZWV1ffBGMiESiS5eMLoy/bDxLRy4lNSmVs188mysm\nXcFnGz/zOjSJcGE3Wey00xxLl3odiUh423d4Hy8vepm/fvZXOjXqxB9+9gfObv2jq64lRoTF1UHB\nYGauXj3H9u2QkuJ1NCLh73DhYV5e9DKPzH6E05ufzh+y/0BWyyyvw5IQC5cxgaDo3Bm++srrKEQi\nQ2J8Ijdn3cw3t3/DhSdcyICJA7hs4mUs3rbY69AkQoRdEjjzTI0LiNRUUkISt/a4lW9v/5bsdtlc\n8OoFDHxjIOvy1nkdmoS5sEsCZ52lJCDir7p16nLX2Xfx7e3f0rlpZ8587kwenPUg+w7v8zo0CVNh\nlwTUEhAJXP3E+jz004dYNGIR6/LWcfI/T+bfi/9NkSvyOjQJM2E3MHzkiKNhQ9i0CRrqud8iQTFv\n0zzufO9OAJ664CldSRRlompgOCEBunaFL788flkRqZ6zW5/NZ8M+49azbuXKSVcy5O0h7Dqwy+uw\nJAyEXRIAjQuI1IY4i2Nwl8GsvG0lDZMb0vlfnXl96euEU2+AhF7YdQc553j1VZgyBSZN8joikej1\n2cbPuHHajXRM78gzFz1D67TWXockfoqq7iDwDQ7rHkIitatXm158efOXdG/RnW7PdmPsF2M1cByD\nwrIlUFQE6enw3XfQpInXUYlEv2XfL2PY1GEkJSTxwiUv0KlxJ69DkhqIupZAXBxkZcHChV5HIhIb\nTmt2Gp/+6lMuO+kyer/Um5cXvayxghgRlkkANF9AJNTi4+K5u9fdzBo8i7/M/QuDJg8i72Ce12FJ\nLQvrJKBxAZHQO7356Xx+0+c0qdeErmO78smGT7wOSWpRWI4JAKxZA336+CaNiYg3pq2axk3TbmLE\nmSN4sM+DJMSF6om0UhNRdSvpknic8w0KL10KLVp4HJhIDNuydws3vH0DB44cYOIVE2nToI3XIUk5\nUTcwDGDm6xLS4LCIt1qmtmTmdTMZcOIAerzQg5x1OV6HJEEUtkkANC4gEi7iLI57z7mXCZdNYOAb\nA3lq3lO6eihKhH0S0BVCIuHj5x1/zrxh8xi3aBw3vH0D+UfyvQ5JAhTWSaDkHkL6wSESPtqnt2fu\nsLkUukL+56X/YX3eeq9DkgCEdRJo1co3NqArhETCS7069Xjl8le4/ozr6flCT2atneV1SOKnsE4C\nJYPDGhcQCT9mxt297ua1K17j2snX8tzC57wOSfwQ1kkANC4gEu7ObX8un/zqE56Y+wQPznpQA8YR\nJuyTgJ4tIBL+Tmh0AnN/NZcP1nzAkClDOFx42OuQpJrCdrJYiW3b4NRTYdcuX/eQiISvA0cOMGjy\nIA4cOcAbV71Bg+QGXocUE6JysliJjAyoX993GwkRCW/16tTjzavfpFOjTvR5uQ+b92z2OiQ5jrBP\nAqBxAZFIEh8Xz5gLx3Bt52vp9WIvln6/1OuQpAoRkQQ0LiASWcyMe8+5l8f6PUbfCX2Zt2me1yFJ\nJSIiCegyUZHIdO3p1zLu0nEM+H8DmLN+jtfhSAXCfmAYfIPC7dtDXp7vqWMiElk+XPMhgyYPYuIV\nE+nboa/X4USdqB4YBmjc2Hdb6dWrvY5ERPzRr0M/Jl89mUGTBzHjmxlehyNlREQSAI0LiES6Ppl9\nmHLNFG54+wamrJzidThSLGKSgK4QEol8vdr0YvovpzP8neH8Z9l/vA5HiLAkoMFhkch3ZsszmXnd\nTO547w5eXfKq1+HEvIgYGAbYvdt3V9G8PEjQY05FIt7yHcvpN6Ef/+j/D6449Qqvw4long8Mm9kF\nZrbSzFab2b0VbL/BzL43sy+Ll1/VdB8NGkDr1rBiRTAiFhGvndr0VKb/cjojp49k5rczvQ4nZgWc\nBMwsDvgncD5wGjDIzE6uoOhE51z34uUlf/alcQGR6NI1oytvDXyL69+6nk82fOJ1ODEpGC2BHsA3\nzrn1zrkjwETg0grKBXz7N40LiESf3m168+ovXuUXr/+CL7d+6XU4MScYSaAVsLHM+03Fn5X3CzNb\nZGaTzKy1PztSS0AkOv2848959uJnuei1i1ixQ32+oRSMJFDRL/zyo7tTgXbOua7AR8B4f3bUrRss\nXQqHdatykahz+SmX83i/xzn/lfNZl7fO63BiRjCus9kEtC3zvjWwpWwB51xumbfPA49XVtno0aNL\n17Ozs8nOzi59X78+dOzoSwTduwcUs4iEocFdBrP30F76TejHnKFzaJHawuuQwlJOTg45OTlBqSvg\nS0TNLB5YBfQFtgILgEHOuRVlymQ457YVr18O/NY517uCuiq9RLTE0KHQqxfcfHNAYYtIGPvT7D/x\nn+X/Yc7QOaQmpXodTtjz9BJR51whcBvwPrAM31VAK8zsYTO7uLjYHWa21My+Ki47xN/9aXBYJPo9\n8JMH6NmqJwPfGEhBUYHX4US1iJksVmL+fBgxAr76KkRBiYgnCooKuPi1i+mQ3oExF47B9HzZSnk+\nWSyUunSBVasgP9/rSESkNiXEJTDpqkl8suET/vbZ37wOJ2pFXBJIToaTT4YlS7yORERqW1pSGu9e\n+y5/n/d33lzxptfhRKWISwKgcQGRWNKmQRumDprKiHdGsGDzAq/DiToRmwQ0aUwkdnRv0Z2XLn2J\nyyZextrctV6HE1UiMgnoATMisefiEy/m/nPu56LXLiI3P/f4X5Bqibirg8A3Yzg9HbZvh5SUEAQm\nImHjzhl3svqH1bwz6B3i4+K9DicsxNTVQQCJidC5sy4TFYlFfz3/rxwsOMjonNFehxIVIjIJgMYF\nRGJVQlwCr1/5OuMXj+ftlW97HU7Ei9gkoHEBkdjVrH4z3rj6DW6edjMrd670OpyIFrFJQC0BkdjW\no1UPHu37KJe/fjl7D+31OpyIFZEDwwAFBdCwIWze7Hv0pIjEpuHThrMzfydvXPVGzN5aIuYGhsH3\nsPmuXWHhQq8jEREvPd3/aTbv2czjn1Z6h3qpQsQmAdC4gIhAUkISb1z9Bk/Pf5r3v3vf63AiTkQn\nAY0LiAhA67TWTLxyIte/db2eSlZDSgIiEhX6ZPbht71/y6DJgzhSeMTrcCJGRCeBTp1g1y7YudPr\nSEQkHPy6169pkNRAE8lqIKKTQFwcZGVpcFhEfOIsjvGXjWfconHMWjvL63AiQkQnAVCXkIgcq3lK\nc16+7GUGvzWYHft3eB1O2Iv4JHDBBTBuHBw86HUkIhIuzut4Hr88/ZcMnTKUcJoLFY4iPgmce67v\nkZOPPeZ1JCISTh459xF2HNjB0/Of9jqUsBaxM4bL2rgRunWDzz7zDRaLiACsyV1Dzxd68v5179Ot\nRTevw6k1MTljuKw2beC+++DWWyGMcpqIeKxDegeevuBprpl8DfsO7/M6nLAUFS0BgCNHoHt3eOgh\nuPrqIAcmIhFt6JShAIy7dJzHkdSOmG8JANSpA2PHwq9/DXv2eB2NiISTf/T/B59u+JTJyyd7HUrY\niZqWQIlhwyA1FZ58MkhBiUhUmLtxLldMuoIlI5bQtH5Tr8MJqkBaAlGXBHbuhNNOg/fe8w0Wi4iU\n+N0Hv2Nd3jomXTXJ61CCSt1BZTRpAo8+CrfcAkVFXkcjIuHkDz/7A19//zWTlkVXEghE1CUBgKFD\nIT4enn/e60hEJJwkJyQz/rLx3DHjDrbv2+51OGEh6rqDSixZAv36wdKl0KxZUKoUkSjxwEcPsHzH\nct4a+FZUPI1M3UEVOOMMuP56+N3vvI5ERMLNqJ+O4tsfvuW1r1/zOhTPRW1LAGDvXjj1VHj1VejT\nJ2jVikgUWLhlIf1f7c+iEYtomdrS63ACopZAJUouFb3lFjh82OtoRCScZLXMYsSZIxj+zvCYvslc\nVCcBgF/8AjIzNW9ARH7swT4PsmH3BiYsnuB1KJ6J6u6gEmvWQI8evofPZGYGvXoRiWCLti3i5//+\nOV/f8jUZKRleh+MXdQcdR4cOcOedvkVEpKyuGV25sduN3D3zbq9D8URMJAHwXSW0YgVMm+Z1JCIS\nbh766UPM3zSfmd/O9DqUkAtKEjCzC8xspZmtNrN7K9ieaGYTzewbM/vMzNoGY781kZQEzzwDt98O\n+/eHeu8iEs7q1anHmAvHMHL6SPKP5HsdTkgFnATMLA74J3A+cBowyMxOLldsGPCDc64T8CTw50D3\n64++faF3b/jTn7zYu4iEs/6d+nNmyzN5ZPYjXocSUgEPDJvZ2cAo51z/4vf3Ac4593iZMu8Vl5lv\nZvHANufcj27jV1sDw2Vt2wannw6zZ8Mpp9TqrkQiknO+pbDQd/+toqKj68d7rUnZsq/VWS/7WXXe\nV/ZZVdv2sZUZ7c6gz5oc6h84rdLvVrYUFByt+6GH4IYbQvPfLJCB4YQg7L8VsLHM+01Aj8rKOOcK\nzSzPzBo5534Iwv5rJCMDRo2CkSNh1iwomTFeNveUrFf1WXW2lS9TslT2vvxrVWUrel9+W1Xfq+5S\nVBTc75X/vOz76m4rWT/ea2XbKls/3nZ/31e0VKdMZeXKnnCrs62y8pWVdc73/0VcnO8eXHFxR5f4\n+KOfVfZa1bbKypavt7L1yrZX9lliYsXlKvt+fHwL0veM5rPUEfyx48fUSYir8LtVLQkJsHUrXHUV\nnHcetGhx/POSl4KRBCrKPuV/zpcvYxWUCZlbboGJE31/BJUpSQ5lbytS/rPqbCtfpmSp7H3516rK\nVvS+/LbjfV7VEhdX8+9U9b3yn5d9X9NtJZ+VL1fZ98qeyCqro/wJrPz3ypePj6+47pJtFX23ou+U\n3VbVvsqfjGuyraK6y54MK4onVl1TNIJeL45nd/txDOs+zK86TjrJ92yT3/wGXgvzO1MEIwlsAsoO\n9LYGtpQrsxFoA2wp7g5Kc87lVlTZ6NGjS9ezs7PJzs4OQojHio+HOXNi+w9dRCoWHxfPc5c8x3n/\nPo9LTrqEZvX9uwPlQw/5nm3y4Ye+m1kGU05ODjk5OUGpKxhjAvHAKqAvsBVYAAxyzq0oU2Yk0Nk5\nN9LMrgEuc85dU0FdtT4mICJSHb+Z+Rt2HNjBhMv9n008dSr89re+uxonJQUxuHI8nSzmnCsEbgPe\nB5YBE51zK8zsYTO7uLjYi0ATM/sGuAu4L9D9iojUpod/9jA563KYtXaW33UMGAAnnwxPPBHEwIIs\nJm4bISLijykrp3Dvh/ey5JYlJMYn+lXH+vWQlQULFvjuXlAbdNsIEZFaMOCkAbRr2I4xC8b4XUdm\npq9L6Lbbjr1iMFwoCYiIVMLM+Pv5f+fRTx5lx/4dftdz992+FsGbbwYxuCBRd5CIyHHc9d5dHCw4\nyNiLx/pdx+zZ8MtfwvLlvmedBFMg3UFKAiIix5Gbn8vJY05m5nUz6ZrR1e96hgyBxo3hr38NXmyg\nJCAiUuv+9fm/eH3Z6/z3hv/i78Ppd+zwzR344APo0iV4sWlgWESklt2UdRM/5P/Amyv879hv2hT+\n+Ee49VbfbTrCgZKAiEg1JMQl8OQFT3LPB/cEdLvpG2/0PfN8Qpg80VJJQESkms5tfy7dW3Tnb5/9\nze864uN9zza57z7IrfDmOaGlMQERkRpYk7uGHs/3YPGIxbRKa+V3PSNH+l6feSbwmDQwLCISQg98\n9ACb9mwK6L5Cubm+Z5pMmwZnnRVYPBoYFhEJofvPuZ+P1n7EvE3z/K4jPR0ee8zXIigsDGJwNaQk\nICJSQ6lJqTx67qP85v3fEEjvxeDBvruLPv98EIOrISUBERE/XHfGdew5tIdpq6f5XUdcnG9M4Pe/\n980h8IKSgIiIH+Lj4nms72Pc/9H9FBb5359zxhlw3XVw771BDK4GlARERPx0YacLaVy3MRMWB3bR\n/+jR8P778OmnwYmrJpQERET8ZGY83u9xRuWMCmgCWVoa/OUv3gwSKwmIiASgV5teZLXMYszn/j9z\nAGDgQMjPh2XLghRYNSkJiIgE6NFzH+XPn/6ZvIN5ftdhBj16wBdfBDGwalASEBEJ0ClNT2HASQP4\n86d/DqieM8+EhQuDFFQ1KQmIiATB6OzRPLvwWbbs3eJ3HVlZoW8J6LYRIiJBcu8H95J3MI9nL3nW\nr+/v2wfNm0NeHtSpU/3v6bYRIiJh4L5z7uPNlW+yaucqv76fkuJ7MP3y5UEOrApKAiIiQZJeN517\net3D/876X7/rCHWXkJKAiEgQ3dHzDuZvns+CzQv8+n6oB4eVBEREgqhunbo81OchRuWM8uv7agmI\niES4IV2HsHzHcr9uNd21q2/C2JEjtRBYBZQERESCLDE+kQfOeYCHP364xt9NSYF27UI3c1hJQESk\nFgztNtTv1kAou4SUBEREakEgrYFQDg4rCYiI1BJ/WwOhbAloxrCISC169otneXvV28z45Yxqf2f/\nfmjWzPcw+sTE45fXjGERkTDlT2ugfn1o3z40g8NKAiIitSgxPpH7z7m/xmMDoeoSUhIQEallQ7sO\nZdn3y2rUGgjV4LCSgIhILUtKSOKBn9TsSiG1BEREokhNWwNdu/ruJnr4cO3GFVASMLN0M3vfzFaZ\n2Uwza1BJuUIz+9LMvjKztwPZp4hIJKppa6BePejYEZYurd24Am0J3Ad86Jw7CZgF3F9Juf3Oue7O\nuW7OucsC3KeISEQqaQ3M3zS/WuVD0SUUaBK4FBhfvD4eqOwE79f1qyIi0SQpIYl7et/D458+Xq3y\noRgcDjQJNHPObQdwzm0DmlZSLsnMFpjZXDO7NMB9iohErGHdhvHJhk+q9fSxULQEEo5XwMw+AJqX\n/QhwwIM12E9b59w2M2sPzDKzJc65tRUVHD16dOl6dnY22dnZNdiNiEh4q59Yn1vPupUn5j7BCwNe\nqLJsly6wYgUcOgRJSUc/z8nJIScnJyjxBHTbCDNbAWQ757abWQbwX+fcKcf5zjhgmnPuzQq26bYR\nIhL1dh3YRad/dOLrW76mVVqrKsuecQaMG+drFVTGy9tGTAWGFK/fAEwpX8DMGppZYvF6E6A3EMLH\nKIuIhJfG9Rpz/RnX8+S8J49btra7hAJNAo8DPzezVUA/4DEAM8sys+eKy5wCfGFmXwEfAf/nnFsZ\n4H5FRCLar3v9mpcWvURufm6V5Wp7cFh3ERUR8cjgtwZzcpOTeeAnD1RaZt48GDkSvvyy8noC6Q5S\nEhAR8cjS75fSb0I/1t65lrp16lZYJj8fGjf23Va67OBwWbqVtIhIBOrcrDNntTqL8YvHV1qmbl3o\n1Am+/rp2YlASEBHx0L3/cy9PzH2CgqKCSsvU5uCwkoCIiIfOaXsOGSkZTF4+udIytTk4rCQgIuKx\n+/7nPh7/9HEqGxNVS0BEJIpddOJFHCo8xIdrPqxwe5cusGqVb+ZwsCkJiIh4LM7i+F3v3/HYp49V\nuD05ufYGh5UERETCwKDTB7F612q+2vpVhdu7d696roC/lARERMJAYnwit511G0/Nf6rC7d27187g\nsJKAiEg3MpkSAAAHOElEQVSYuCnrJqasmsK2fdt+tC0rSy0BEZGo1qhuIwaeNpCxX4z90bYuXXzP\nHD5yJLj7VBIQEQkjd/S8g7FfjOVQwbGXAtWvD+3awbJlwd2fkoCISBg5tempdMnowsSlE3+0rTYG\nh5UERETCzJ097+Sp+U/9aPJYVlbwB4eVBEREwswFJ1zA/iP7mbNhzjGfqyUgIhID4iyOO3rc8aPL\nRbt29U0YK6j8XnM131fwqhIRkWC5oesN5KzLYV3eutLP0tKgVStYGcRnMyoJiIiEoZTEFIZ2Hco/\nF/zzmM+D3SWkJCAiEqZu63Eb4xaNY9/hfaWfBXtwWElARCRMtWvYjux22YxfdPTJY2oJiIjEkLt6\n3sVT85+iyBUB0K0bLFoERUXBqV9JQEQkjJ3T9hxSElN479v3AEhPh2bNYPXq4NSvJCAiEsbMjLvO\nvosn5z1Z+lkwu4SUBEREwtzA0wayePtiVu/y/fwP5uCwkoCISJhLSkhiaNehPPvFs0BwWwJW2YON\nvWBmLpziEREJF2ty19DzhZ5suGsD+3fX5YQT4IcfIC7O12XknDN/6lVLQEQkAnRI78CZLc/kP8v/\nQ5Mm0KABrFkTeL1KAiIiEWJE1ojSB84Eq0tISUBEJEJcdOJFbNyzkcXbFgdtcFhJQEQkQiTEJXBT\n95sY+8XYoLUENDAsIhJBNu/ZzOn/Op0Fg9bTs1sqO3dCXJwGhkVEYkKrtFb8rP3P+PD7V0lOhvXr\nA6tPSUBEJMKUDBB36+4C7hJSEhARiTB9O/Rl/5H9ZGTND3hwWElARCTCxFkcw7OGs6bxv7xtCZjZ\nlWa21MwKzax7FeUuMLOVZrbazO4NZJ8iIgJDug5h4b4pfL70h4DqCbQl8DVwOfBxZQXMLA74J3A+\ncBowyMxODnC/US8nJ8frEMKGjsVROhZHxfqxaFKvCQNOuoRDJ48/fuEqBJQEnHOrnHPfAFVdmtQD\n+MY5t945dwSYCFwayH5jQaz/gZelY3GUjsVROhYw4swRFGWNDaiOUIwJtAI2lnm/qfgzEREJQO82\nvUlJTgqojuMmATP7wMyWlFm+Ln69pJr7qKiVoBlhIiIBMjMubTUisDqCMUPXzP4L/MY596NxajM7\nGxjtnLug+P19gHPOPV5BWSUHERE/+DtjOCGIMVQWwOfACWaWCWwFrgEGVVTQ33+EiIj4J9BLRC8z\ns43A2cA7Zjaj+PMWZvYOgHOuELgNeB9YBkx0zq0ILGwREQmGsLqBnIiIhJYnM4aPN3nMzBLNbKKZ\nfWNmn5lZWy/iDIVqHIu7zWyZmS0qHqRv40WcoVDdSYXFkxSLqpqgGOmqcyzM7Oriv42vzeyVUMcY\nKtX4f6SNmc0ysy+L/z/p70Wctc3MXjSz7Wa2pIoyTxefNxeZWddqVeycC+mCL/F8C2QCdYBFwMnl\nytwCPFO8PhBfF1LIYw2TY/FTILl4fUQsH4vicin4JifOBbp7HbeHfxcnAAuBtOL3TbyO28Nj8Sww\nvHj9FGCt13HX0rE4B+gKLKlke3/g3eL1nsC86tTrRUugOpPHLgVKpsG9AfQNYXyhdNxj4Zz72Dl3\nsPjtPKJ3jkV1JxX+EXgcOBTK4EKsOsfiJmCMc24PgHNuZ4hjDJXqHIsiIK14vSGwOYTxhYxz7hMg\nt4oilwITisvOBxqYWfPj1etFEqjO5LHSMs43sJxnZo1CE15I1XQi3TBgRq1G5J3jHovi5m1r59z0\nUAbmger8XZwInGRmn5jZXDM7P2TRhVZ1jsXDwPXFF6m8A9weotjCTfljtZlq/GgM5iWi1VWdyWPl\ny1gFZaJBtSfSmdl1QBa+7qFoVOWxMDMD/g7ccJzvRIPq/F0k4OsS6gO0BeaY2WklLYMoUp1jMQgY\n55z7e/G8pFfw3acs1vg1MdeLlsAmfH+0JVoDW8qV2Qi0ATCzeHz9nlU1gyJVdY4FZtYPuB+4pLhJ\nHI2OdyxS8f2PnWNma/FdljwlSgeHq/N3sQmY4pwrcs6tA1YBnUITXkhV51gMAyYBOOfmAclm1iQ0\n4YWVTRSfN4tVeD4pz4skUDp5zMwS8U0em1quzDSO/uK7CpgVwvhC6bjHwsy6AWOBAc65XR7EGCpV\nHgvn3B7nXDPnXAfnXHt84yOXuApmqUeB6vw/8jZwLkDxCa8TsCakUYZGdY7FeqAfgJmdAiRF8RiJ\nUXkLeCowGErv1JDnnNt+vApD3h3knCs0s5LJY3HAi865FWb2MPC5c+4d4EXg32b2DbAL33/4qFPN\nY/FnoD7wn+IukfXOucu8i7p2VPNYHPMVorQ7qDrHwjk308zOM7NlQAFwTzS2lqv5d3EP8LyZ3Y1v\nkPiGymuMXGb2GpANNDazDcAoIBHfbXiec85NN7MLzexbYD8wtFr1Fl9OJCIiMUiPlxQRiWFKAiIi\nMUxJQEQkhikJiIjEMCUBEZEYpiQgIhLDlARERGKYkoCISAz7/4iJKvWcGNXOAAAAAElFTkSuQmCC\n",
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"pt.plot(mesh, x)\n",
"pt.plot(mesh, x_true, label=\"true\")\n",
"pt.legend()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"* Ideas to accelerate this?\n",
"* Multigrid"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"## Gauss-Seidel"
]
},
{
"cell_type": "code",
"execution_count": 45,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [],
"source": [
"x = np.zeros(n)"
]
},
{
"cell_type": "code",
"execution_count": 64,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [
{
"data": {
"text/plain": [
""
]
},
"execution_count": 64,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAEACAYAAABVtcpZAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xd4FWXax/Hvk4QkENLoJRCCooDoUlwFVIzgqiAo2ABX\nKWJBFsvaZX1XsC2sr66iKCIIovIKigpIXcUINkAR6YISOoSWQAKEkOR5/5g0YoCQnJw5Oef3ua7n\nmjlhMnNnSOaep80Yay0iIhKYgtwOQERE3KMkICISwJQEREQCmJKAiEgAUxIQEQlgSgIiIgHMI0nA\nGDPBGJNijFl5kn+/3BiTZoxZnlee8sRxRUSkfEI8tJ+JwGvA5FNss8hae52HjiciIh7gkZqAtfYb\nIPU0mxlPHEtERDzHm30C7Y0xPxtjZhtjWnrxuCIichKeag46nZ+AeGvtEWNMV+Az4BwvHVtERE7C\nK0nAWptRZH2uMeYNY0wNa+2BotsZY/QgIxGRMrDWlqnJ3ZPNQYaTtPsbY+oWWb8IMMUTQD5rrYq1\nPP30067H4CtF50LnQufi1KU8PFITMMZMARKBmsaYrcDTQChgrbXjgJuMMfcCx4GjQG9PHFdERMrH\nI0nAWnvraf59DDDGE8cSERHP0YxhH5WYmOh2CD5D56KQzkUhnQvPMOVtT/IkY4z1pXhERCoDYwy2\njB3D3hoiKiJyxpo0acKWLVvcDsNnxMfHs3nzZo/uUzUBEfFZeXe4bofhM052PspTE1CfgIhIAFMS\nEBEJYEoCIiIBTElARCSAKQmIiJRRQkICCxcudDuMclESEBGpADk5OW6HUCpKAiIiZdCvXz+2bt1K\n9+7diYqK4sUXXyQoKIh33nmH+Ph4unTpwtdff02jRo1O+L6itQdrLSNHjuTss8+mdu3a9OnTh7S0\nNK/+HEoCIiJlMHnyZBo3bszs2bM5dOgQt9xyCwCLFi1i/fr1zJ8/H3DG8J/Mq6++ysyZM1m8eDE7\nd+4kNjaWIUOGeCX+fJoxLCKVmhnhmTfX2qfLNimt6OQtYwwjRoygatWqpfrecePGMWbMGOrXrw/A\nP//5T+Lj43n//fcJCvLOPbqSgIhUamW9eFeUuLi4Um+7ZcsWevXqVXDBt9ZSpUoVUlJSChJDRVMS\nEBEpo5Kaeop+LSIigiNHjhR8zsnJYe/evQWfGzduzDvvvEOHDh0qNtBTUJ+AiEgZ1atXj02bNgGU\n+Javc845h8zMTObOnUt2djbPPfccWVlZBf9+zz33MGzYMLZu3QrA3r17mTlzpvd+AJQERETK7Ikn\nnuDZZ5+lRo0aTJ8+/Q81g6ioKN544w0GDRpEXFwckZGRJzQXPfDAA1x//fVcddVVREdH07FjR5Yu\nXerVn0FPERURn6WniJ5ITxEVERGPUhIQEQlgSgIiIgFMSUBEJIApCYiIBDAlARGRAKYZwyLis+Lj\n40/5ALZAEx8f7/F9ap6AiEglp3kCIiJSJmoOklLLyski/Vg6GVkZpGc5y4ysDI4cP0JmdmaJJTs3\nm5zcHHJszgnLXJtLkAkiOCiYYBN8wnqV4CqEh4RTNaQqVatULViGh4RTPbQ6MeExBSUyNFLNBSLl\noOagAJOTm0NqZir7juxj35F97D+yv3D96H7SMtNIy0wjNTPVWR51loeOHcJiiQyNpHpodSLDnGX1\n0OpUq1KN8JBwpwSHF66HhBMSFFJwcS+6DDJB5NrcgoRQNEkczzlOZnYmR7OPcvT4UWeZt56RlcHB\nYwcL4jx6/ChRYVHEhMdQq1ot6lavS72Ies6yej3qRjjLRtGNiIuKIyRI9z3if8rTHKQk4AdybS57\nD+9lV8YudqXvYlfGLlIyUkg5nMKew3tIOZxCSoazfuDoAaLDo6lVrRY1q9akVrVaBaVG1RrEhscW\n3GXHVi1cjw6LJjQ41OfuurNzszl07BCpR53EtjtjNymHU5xlRgq7D+9mV/outh3aRkpGCg0iG9Ak\npskJ5dya59KidgtiwmPc/nFEykRJwI+lH0tnR/oOth/azo5Decv0HexI38HO9J3sSt/FnsN7iAqL\nokFkA+pH1qd+9frUjahL3ep1qRtRlzoRdQrWa1arGbB3w1k5WWw7uI3NaZsLSnJaMuv3rWf9vvVE\nhUXRonYLWtZq6Sxrt6RNvTZEh0e7HbrIKSkJVFKZ2ZlsO7iNrQe3su3QNrYd3OYsi6xn52YTFxVH\nw8iGNIxqSMPIhsRFxdEgskFBqVe9HqHBoW7/OJVars1l+6HtrN27lnV717Fu3zpW71nNypSVNIhs\nQLsG7WhX3ylt67dVYhCfoiTgo9Iy09iStoXNaZvZcnALW9K2sPXQVmd5cCupmanERcXROLoxjaIa\nOSX6xGVMeIzPNcEEkpzcHNbvW8+PO3/kp10/8dOun/hl9y/ERcVxWePL6BTfiU7xnYiP8fz4bZHS\nUhJwycHMgySnJTvNCqnJBRf7/KaGHJtDfHQ8TWKaEB8dT3xMPI2jGxes16tejyCjUbqVTU5uDqv2\nrGLxlsUs2rqIRVsWER4S7iSExp3onNCZs2qc5XaYEkCUBCrIkeNHCi7wyWnJBcv8tuTs3GwSYhL+\n0NGYf9GvUbWG7uIDgLWWDfs3sGjLIr7e8jVfJn9JVFgU3c7uRrdm3egU34mwkDC3wxQ/piRQRjm5\nOWw/tJ3ktGQ2pW4iOTWZTWmbCi72qUdTiY+JJyEmgYSYBJrGNqVJTBMSYp3PushLSXJtLit2r2DO\nxjnM/W0uq/esJrFJIl3P7sp1515Hg8gGbocofsb1JGCMmQB0B1KstRecZJvRQFfgMDDAWruihG08\nngRSj6YWXOSLluS0ZLYe3ErtarVJiHUu8PkX+oSYBBJiE2gQ2UDNNVJu+4/sZ8HvC5i9cTZzNs6h\nVZ1W9D6vNze2vJF61eu5HZ74AV9IApcCGcDkkpKAMaYrMNRae60x5mLgVWtt+xK2O+MkcDznOFsP\nbj3xIp9WuJ6dm81ZsWedcJFvGtuUhFinGSc8JLyMP7XImTuWfYwFvy9g2tppfL7hc1rXa+0khBY3\nUjuittvhSSXlehLICyIemHWSJDAW+MpaOzXv8zog0VqbUmy7PyQBay2pmaknXOR/P/B7wYV+Z/pO\n6lWvV3ChL15qVq2pJhvxSUePH2Xeb/OYtnYaczfO5bL4y7izzZ1ce861ATuXQ8qmMiSBWcC/rLXf\n5X3+AnjMWru82Hb2zWVv/qHpxmJpGtu0xAt94+jGGiMvld7hrMNMWzON8T+PJzk1mQGtBzCozSCN\nMpJSKU8S8NbtRknBlZh9Xn3ubRrWiCU2PJZ7Ot/Dzf1uJjY8Vnfz4tciQiMY2GYgA9sMZO3etUxY\nPoEOEzrQqk4r7mx7Jze1vEk3O1IgKSmJpKQkj+zLreag9cDlJTUHXXqpZdEi0DVfAt2x7GPM/HUm\nY38ay/p967nvovu4p909xFaNdTs08TG+8j4BQ8l3/AAzgX4Axpj2QFrxBJAvIwOmTvVgVCKVVFhI\nGDefdzNf9vuSObfOYd2+dZw1+izum3Mfvx/43e3wxE94anTQFCARqAmkAE8DoYC11o7L2+Z14Bqc\nIaIDi/cH5G1jFy2y/PWvsG4dRESUOzQRv7IzfSevL32dt5e/Taf4TjzS4RE6NOrgdljiMp/oGPaE\n/NFBffrAuefCiBFuRyTimzKyMpi0YhIvff8SzWo045krnqF93B9GXUuA8LsksHUrtGkDy5dDBbxX\nWcRvZOVkMWnFJJ5b9Bzn1z2fZxKfoV2Ddm6HJV7md0kAnFrA6tXw0UcuByVSCRzLPsb45eN54ZsX\n+HODPzMicQR/qvcnt8MSL/HLJHDkCLRsCRMnwhVXuByYSCVx9PhR3vrpLUZ9O4pO8Z0YdeUomsQ0\ncTssqWC+MjrIo6pVgxdfhAcfhOxst6MRqRyqVqnKg+0f5Lf7fqNV7VZcOO5Cnlr4FBlZGW6HJj7K\nZ5MAwE03QWwsvP2225GIVC4RoRH8z+X/w4rBK9ictpnmrzfnvV/eI9fmuh2a+BifbQ7K98svcNVV\nzpDRGjVcCkykkvth+w88MO8BAF695lWNJPIzftknUNS990JICLz2mgtBifiJXJvL+yvfZ9iXw7iy\n6ZW8dNVL1KxW0+2wxAP8sk+gqGefhWnTnCGjIlI2QSaIfn/qx/qh64kJj6HVm62YunoqvnQjKN5X\nKWoCABMmwLhx8N13EBzs5cBE/ND3277nzll3clbsWbxx7RvERcW5HZKUkd/XBAAGDnSahNRJLOIZ\nHRp1YPndy2lbvy1t3mrD2B/HquM4AFWamgDAqlXQpYuzrFvXi4GJ+Lk1e9YwaOYgwkLCGN9jPM1q\nNnM7JDkDAVETADj/fOjfHx591O1IRPzLeXXO49s7vqXnuT3p+E5HJq2YpL6CAFGpagLgPGq6ZUuY\nPBkSE70Tl0ggWZWyir7T+9KqTivGdh9LTHiM2yHJaQRMTQCgenV49VVn2GhWltvRiPif8+uez7K7\nllGrWi1aj23NN1u/cTskqUCVriYAYC1cdx106ADDhnkhMJEANevXWdw16y4GXziYpzo9RUiQt95I\nK2fC7yeLlWTzZrjwQli2DBISKjYukUC2M30n/T/rz5HjR/jwxg9pFN3I7ZCkmIBqDsrXpAk8/DDc\nd59TMxCRitEgsgHzb5vPdedcx0XjLyJpc5LbIYkHVdqaADh9Aq1bw/PPQ69eFRiYiADw39//y22f\n3sawS4dx/8X3Y0yZbj7FwwKyOSjfokVw663O3IHY2AoKTEQKJKcm02tqLy6oewFvdX+LqlWquh1S\nwAvI5qB8nTo5tYC//93tSEQCQ0JsAt8N+o4cm8Ml71zClrQtbock5VDpkwDAyJGweDHMnu12JCKB\noVqVarzf631uv+B2Lh5/MQuTF7odkpRRpW8OypeUBLfdpmYhEW9bmLyQW6ffyjNXPMPd7e52O5yA\nFNB9AkUNHQqHDzvvJRYR7/ntwG90/aArvc/rzbNXPKsOYy9TEsiTkQEXXACvvw7dunkwMBE5rb2H\n99L9/7rTvFZz3u7xNqHBoW6HFDACumO4qOrVYfx4uOceSEtzOxqRwFI7ojZf9f+KtMw0rp1yLQcz\nD7odkpSCXyUBgM6doUcPeOghtyMRCTzVqlTjk1s+oVmNZnSa1Ikdh3a4HZKcht8lAYBRo+Crr2Du\nXLcjEQk8wUHBjOk2hltb3UqHCR1YvWe12yHJKfhVn0BRX34JAwbA6tUQHe2RXYrIGZqyagp/n/93\nZvSZQfu49m6H47fUMXwS994L6enw3nugwQoi7pizcQ4DPhvA9Fumc1n8ZW6H45fUMXwSL70Ey5c7\nL6AREXd0a9aNKTdO4YZpN/Dlpi/dDkeK8euaADjNQVdc4cwobt7co7sWkTOwaMsibpp2E+/2fJeu\nzbq6HY5fUU3gFFq1cp4y2rs3ZGa6HY1I4OoU34kZfWbQ/7P+zFg/w+1wJI/f1wTAed9A795Qp44z\nkUxE3PPjzh/pPqU7r3V9jZvPu9ntcPyCagKnYQyMG+c8YO6TT9yORiSwXdjgQubfNp/7593PBys/\ncDucgBcQNYF8S5Y4E8mWLYP4+Ao7jIiUwtq9a7ly8pW81vU1bmx5o9vhVGqu1wSMMdcYY9YbYzYY\nYx4v4d/7G2P2GGOW55U7PHHcM3XxxfDoo9C3Lxw/7kYEIpKvZe2WzPnrHIbMGcL83+a7HU7AKndN\nwBgTBGwAugA7gWVAH2vt+iLb9AfaWWvvP82+KrQmAJCb6zxcrm1beOGFCj2UiJTCd9u+o+eHPfmk\n9ydc2vhSt8OplNyuCVwEbLTWbrHWHgc+BK4vYTufmK4VFOTMG3j3XT1WQsQXdGzUkQ9u+IAbpt7A\n8l3L3Q4n4HgiCTQEthX5vD3va8XdYIxZYYyZZoyJ88Bxy6xOHZg2Dfr3h3Xr3IxERAD+ctZfeKv7\nW1w75VrW7dUfpTd5IgmUdIdfvE1nJtDEWtsa+BJ41wPHLZdLLoF//9vpKN6/3+1oRKRXi16MunIU\nV79/NZvTNrsdTsAI8cA+tgONi3yOw+kbKGCtTS3y8W1g1Ml2Nnz48IL1xMREEhMTPRBiyQYMgDVr\n4KabYP58CNU7MERc1e9P/Ug/ls6Vk69k8cDF1I+s73ZIPikpKYmkpCSP7MsTHcPBwK84HcO7gKVA\nX2vtuiLb1LPW7s5b7wU8aq3tWMK+KrxjuLicHOjZExo0gLFj9aA5EV/w/KLn+WjtRyweuJjIsEi3\nw/F5rnYMW2tzgKHAAmAN8KG1dp0xZoQxpnveZvcbY1YbY37O23ZAeY/rKcHB8MEH8O23mk0s4iuG\nXTaMixteTO+Pe5Odm+12OH4toCaLnUpyMnTs6IwauuoqV0IQkSKyc7PpPqU7TWObMqbbGL28/hTc\nHiLqFxISnBFDt90G69effnsRqVghQSFMu3ka32z9hpe/f9ntcPyWkkARl10GI0c6I4b27nU7GhGJ\nCoti9q2z+c8P/+GTdXrwV0VQEijmjjugTx+nSSg19fTbi0jFahTdiJl9ZzL488Es3bHU7XD8jpJA\nCZ55BhIT4Zpr4NAht6MRkbb12/LO9e/Q88OeJKcmux2OX1HH8ElY67yjeO1amDcPqlVzOyIReW3J\na7z545t8e8e3xFaNdTscn6EXzVeQ3FxnQtnu3TBzJoSHux2RiDww9wE2HNjA530/Jzgo2O1wfIJG\nB1WQoCB45x2IiYFbbtHjp0V8wUtXv0RmdibDk4a7HYpfUBI4jZAQeP99p3nottucGcYi4p6QoBCm\n3jSVd395l8/Wf+Z2OJWekkAphIbCRx/BgQPO6KFsTWAUcVWdiDp8fMvH3D3rbtbv08Se8lASKKXw\ncPjsM9i1C264AY4ccTsikcB2UcOLeKHLC/Sa2ov0Y+luh1NpKQmcgYgI+PxziI2Fzp01oUzEbXe2\nvZNOjTsxYMYAfGlQSWWiJHCGQkNh0iTo0sV5J8Hvv7sdkUhgG911NDsO7WDUtyd9Qr2cgpJAGRgD\nzz8PDz3kPGpi2TK3IxIJXGEhYXx8y8eMXjKaBb8vcDucSkdJoBwGD3beQdCtG8ye7XY0IoErLiqO\nD2/6kNs/vV1vJTtDSgLldN11MGsW3HknvPWWM5RURLyvU3wnHu34KH2n9+V4jib1lJZmDHvIxo3O\nG8ratoU33oBIvQxJxOtybS7dPuhGu/rteL7L826H4zWaMewDmjVz+gbCwuDCC2HFCrcjEgk8QSaI\nd3u+y8QVE1mYvNDtcCoFJQEPqlYNxo+Hp5+Gv/wFxoxR85CIt9WtXpdJPSfR79N+7D2scdyno+ag\nCrJxI/Tu7byxbMIE5/lDIuI9j//3cdbsXcOsvrP8/tWUag7yQc2awfffQ8OG0KaN8yJ7EfGe5zo/\nx94jexm9ZLTbofg01QS84LPPYMgQZyjpyJFQq5bbEYkEhk2pm7h4/MUsuG0Bbeq3cTucCqOagI/r\n2RPWrXMeO3HeeU6/QW6u21GJ+L+msU0Zfc1o+kzvQ0ZWhtvh+CTVBLzs55+dWgE4Q0nb+O/NiYjP\nGDhjIAATr5/ociQVQzWBSiS/f+DOO513GN93n15oL1LRXuv6Gt9u/Zbpa6e7HYrPURJwQVAQDBrk\nvL84K8vpRH7qKdi/3+3IRPxT9dDqTOo5iaFzh2rYaDFKAi6qWdN51MTSpbBnj5MMHnsMUlLcjkzE\n/3Rs1JHbL7idv835m9uh+BQlAR/QtCmMG+fMMj5yBFq0gAcfhJ073Y5MxL88c8UzrNqzimlrprkd\nis9QEvAhjRvD66/D6tXO46pbtYLbb4ekJM08FvGE8JBw3u35LvfPvZ+UDFW5QaODfNq+fc5L7idM\ngKNHnfcb9+/vTEATkbIb9uUw1u5dy6e9P/WL2cQaHeSnatVymoVWroQpU2DLFjj/fOjeHT75xEkM\nIpWRtYXFDU9f/jS/HfiNKaumuBOAD1FNoJI5fBimT3decfnjj5CY6CSF7t2hQQO3oxN/Y63zO7d/\nv1MOHHCGNOcvi5b0dGfbjIwTy5EjJ7/gBwc77+yuVcsZKFG0NGoErVs7JSrK8z/bTzt/ousHXVkx\neAUNIiv3H095agJKApVYairMm+e81GbePOdhdT16OI+naNsWQkLcjlB8SW4uHDxYeEHft+/Uy/wS\nElJ4Ya5Rw7lox8aeuB4T41yoq1cvLBERzrJaNWdYNDh9XUVbX44fdxJK0ePll+RkZ3LlqlVQv74z\nxya/XHKJZxLDP7/6Jz/v/pmZfWZW6mYhJQEhO9uZhPb55zB3rtN0dOGF0LGjUzp0cP5opfKzFg4d\nKrwjL1ry79ZLurCmpTkX5vwLetG775Ot16wJVau6+/NmZ8OGDU5C+PlnWL4cfvoJrroKbr0VunaF\n8PCy7TsrJ4s/v/1nHmr/EP1b9/ds4F6kJCB/kJoKS5bAd985ZelSp0P54oud5xfll0aNCu/SxHvy\n78rT0gpL0aaV4p+LN8FUq1Z4F170Lj2/5H8uejGvUcN/aocHDjj9YlOmOEOre/VyEkJiotPEdCZW\n7F7BX977C6vuXUW96vUqJN6KpiQgp5WT4ww9XbYM1qxxytq1zoWoRQto2RLOPdcZptq4McTHO30M\n/nLR8KTsbKf9++BB5468aMn/Wlqas55fin5OS3O+v3r1wqaU6OjCi3pJJb/ppUYNZ/sqVdw+C75j\nxw6YOtVJCDt3wsMPw9Chzlv+SuvJL55k88HN/N+N/1dxgVYgJQEps7Q05wmna9Y4L8LZurWwpKRA\nvXpOUqhf32kmKKlERzsXtMhI5w7VF5pWc3IgM9MpR4865fBhp5Myv+R/Lt6ZWfRzerpTDh0qXM/K\ncn7W6GinXToq6o/rJysxMYXt52d6xyqnt3o1DBvmLEeOhJtvLt3v45HjR2j1RivevPZNrj776ooP\n1MNcTwLGmGuAV3CGnE6w1o4q9u+hwGSgHbAP6G2t3VrCfpQEfMjx485d1pYtTkLYt+/Esnevszx0\nqPCCmZnptDtHRjrL8HDnrjU09MRSpUphJ2HxAk5zSU5OYcn/nJ3tXITzy/Hjhev5F/3MTGe78HCn\nPTsszImlWrXCUvRzSZ2Z+SUy8o+lalXfSHRycl99BY884vyuvfSS0y92OnM3zmXo3KGsvnc1Vau4\n3BFyhlxNAsaYIGAD0AXYCSwD+lhr1xfZ5l7gfGvtEGNMb6CXtbZPCftSEqjkcnIK76Tz75pLumBn\nZZ04Vrz4uPHgYKevIji4sAQFOc1TYWF/TCxVqjgX/fwLf36SkcCVmwsffAD/+IfTFzZyJJx11qm/\np/fHvTk79mye7/K8d4L0ELeTQHvgaWtt17zPTwC2aG3AGDMvb5slxphgYLe1tnYJ+1ISEBGPOnIE\nXnkFXn4ZnnsO7rnn5DcIu9J3ccHYC0jqn8R5dc7zbqDl4PaM4YbAtiKft+d9rcRtrLU5QJoxRgMW\nRaTCVavm9BN8/z289hrcdZfTZFiS+pH1GX75cAbPHkyuDYzX/3kiCZSUfYrfzhffxpSwjYhIhWnW\nzBk2ffAgXH45bN9e8naDLxzMsexjTPzZP99CVpwnBgBuBxoX+RyH0zdQ1DagEbAzrzkoylpb4vu0\nhg8fXrCemJhIYmKiB0IUEXE6+6dNg3//Gy66CD78EDp1OnGb4KBgxvUYx1XvXUWPc3tQJ6KOO8Ge\nQlJSEklJSR7Zlyf6BIKBX3E6hncBS4G+1tp1RbYZArTK6xjuA/RUx7CIuGnBAudR7f/4h/Oa1+L9\nBA/Pf5i9R/YyuddkdwI8A74yRPRVCoeIjjTGjACWWWs/N8aEAe8BbYD9OKOHNpewHyUBEfGaTZvg\nhhuch9SNH3/i5MiMrAxajmnJpJ6T6JzQ2b0gS8H1JOApSgIi4m1HjjiJoG5dmDjxxMeozFg/g8e/\neJyV964kNDjUvSBPw+3RQSIilVa1as5ziJKT4f77T3zk9XXnXkeTmCaMWTrGvQArmGoCIiI4o4Y6\nd4ZrroHni8wVW7d3HZ0mdWLtkLXUjvjD9CafoOYgEREP2LfPGT56++3wxBOFX39w3oNkZmcytvtY\n94I7BSUBEREP2bkTLrvMeRrpkCHO11KPptJ8THPm3zaf1vVauxtgCdQnICLiIQ0awBdfwL/+Be+9\n53wttmoswy8fzoPzHsTfblSVBEREiklIcOYRPPYYfPqp87W72t3FgaMH+GTdJ+4G52FqDhIROYnl\ny+Hqq+Hrr50XLy1MXsigmYNYO2StTz1uWs1BIiIVoG1b5xETN9/sPCK9c0Jn2tZvy8vfv+x2aB6j\nmoCIyClYCwMGOJPIJk6ETambuOjti/hl8C80jCr+wGR3qCYgIlJBjIE33nCeQDppEjSNbcrd7e7m\nyS+fdDs0j1BNQESkFNasgcRESEqCxmen03xMc6bfMp32ce3dDk01ARGRinbeefDii07/QFB2JC90\nfoGHFzxc6YeMKgmIiJTSgAHO+4qHDIG/nn8bh44dYtaGWW6HVS5qDhIROQOHDzsvpHnkEahzyWwe\n++IxVg5eSXBQsGsxqTlIRMRLIiLgo4+ciWSNMrtRs2pNJv/i+y+eORnVBEREymDiRBg9GkZ/+j1/\n/bQ3vw791bUJZKoJiIh42YABULMmLPukA+0atGPMssr5zgHVBEREymjjRujQAaZ+uY6+Cy5nw30b\niAmP8XocepS0iIhLnnvOmUhW5647qRtRhxe6vOD1GNQcJCLiksceg99/h/ZHh/PWT2+xM32n2yGd\nEdUERETKafFi6NsXbhr7OEdtGm/1eMurx1dzkIiIy+66C6iaymeNzuGbgd9wbq1zvXZsJQEREZcd\nOOA8WuLGl0exO3gZH9/ysdeOrT4BERGX1agB//u/sOjF+1myfQlLdyx1O6RSUU1ARMRDrHXeRBbd\neRwZjT5l7l/neuW4qgmIiPgAY+DNN2HhywNYuWstP2z/we2QTktJQETEg846Cx4YGkr934Yx4usR\nbodzWmoOEhHxsMOHoVnzLOzQZnz616kV/uIZNQeJiPiQiAgY/j+hRK0cxogk364NKAmIiFSAO+4A\nfh7I8m0RyBtyAAAJ0ElEQVS+3TegJCAiUgFCQmDk86FUWTKM4T5cG1ASEBGpID17QqN9A/lxi+/W\nBpQEREQqiDHw4shQ7KIneXqhb9YGlARERCrQpZdCh6oDWbp5jU/WBpQEREQq2L9fCCP7q2E89YXv\n1QaUBEREKljLlnDz2QNZmux7tYFyJQFjTKwxZoEx5ldjzHxjTPRJtssxxiw3xvxsjPmsPMcUEamM\nnh0eRk7SMJ6c51u1gfLWBJ4AvrDWngssBJ48yXaHrbVtrbVtrLU9y3lMEZFKp2FD+NslA1m2ZQ1L\nti9xO5wC5XpshDFmPXC5tTbFGFMPSLLWNi9hu3RrbWQp9qfHRoiI3zp4EOJuHE37Pkn8985PPLZf\nNx8bUcdamwJgrd0N1D7JdmHGmKXGmO+MMdeX85giIpVSdDQ8duUgFm/5hl/3/ep2OACEnG4DY8x/\ngbpFvwRY4KkzOE5ja+1uY0wCsNAYs9Jam1zShsOHDy9YT0xMJDEx8QwOIyLi2x4YEsHIPn/jqXkv\n8tFt48u0j6SkJJKSkjwST3mbg9YBiUWag76y1rY4zfdMBGZZa/9QF1JzkIgEgsdH7OeVnGZsemQV\nDaMalnt/bjYHzQQG5K33B2YU38AYE2OMCc1brwV0BNaW87giIpXW4/fVxKy6nWcWvOJ2KOVOAqOA\nvxhjfgWuBEYCGGPaGWPG5W3TAvjRGPMz8CXwL2vt+nIeV0Sk0qpRA+5o/hDvrnqH1KOprsail8qI\niLhgzx5odH8/HurXnH91G1aufemlMiIilUydOnBzg8cYvWQ0R48fdS0OJQEREZeMfKgVxzf/mTHf\nvutaDEoCIiIuiYuDrlGP83zSi2TnZrsSg5KAiIiLXnnoUjJ21WPyj9NdOb6SgIiIixISoFPQEzw1\ndxRuDIxREhARcdnrD1zLngPHmLXmC68fW0lARMRlLZoH0S7zMR7+dKTXj60kICLiA964ty+bDm3g\n+80/e/W4SgIiIj6gXetQmqUO5eGPXvXqcZUERER8xL9uuoslaTPYlb7ba8dUEhAR8RE9r65BzPbe\nPDp1rNeOqSQgIuIjjIGHLrmfj5LHciz7mFeOqSQgIuJDHunfElL+xEvzP/TK8ZQERER8SFgY3NL4\nAV7+7lWvTB5TEhAR8TH/O/gaUjMOM/OXxRV+LCUBEREfU7dOEO2D7ufJGRU/XFRJQETEB70ysD/r\nM5P4NWVzhR5HSUBExAf9+YLqxB8YyAMfvF6hx1ESEBHxUSOuHcoX+yaSfiyjwo6hJCAi4qNu79GE\nansTGTat4t48piQgIuKjjIEhbR/knTWvkmtzK+QYSgIiIj7s6QGXkpVRnXFfzauQ/SsJiIj4sKpV\nDd1rP8hz/32lQvavJCAi4uNG392bnTm/sPT3DR7ft5KAiIiPa1Q/jJZZA3ls6lse37eSgIhIJTC8\nx90sTp/M4WNHPbpfJQERkUrgxs5NiTh4IcM//sij+1USEBGpBIyB25sPZsIKz75wRklARKSSeGHg\ntRy025i34heP7VNJQESkkoiODOHPwXcx7BPP1QaUBEREKpF/9xnEiuNT2Xco3SP7UxIQEalEOrVp\nSK2MK3j8gw88sj8lARGRSubeCwcz9fexHnn9pJKAiEgl849bu5CZe5gPvl5S7n0pCYiIVDKhVYLo\nEn0Pz859s9z7KlcSMMbcZIxZbYzJMca0PcV21xhj1htjNhhjHi/PMUVEBF7qN4ANwTPYnHKgXPsp\nb01gFdAL+PpkGxhjgoDXgauB84C+xpjm5Tyu30tKSnI7BJ+hc1FI56JQoJ+LVgm1iD/ag4cml++F\nM+VKAtbaX621GwFzis0uAjZaa7dYa48DHwLXl+e4gSDQf8GL0rkopHNRSOcCHrliMLN3l2/OgDf6\nBBoC24p83p73NRERKYch3TticsPKtY/TJgFjzH+NMSuLlFV5yx6lPEZJtYTyj2sSEQlwQUGG6xsO\nLtc+jCfGmRpjvgIettYuL+Hf2gPDrbXX5H1+ArDW2lElbKvkICJSBtbaUzXLn1SIB2M4WQDLgLON\nMfHALqAP0LekDcv6Q4iISNmUd4hoT2PMNqA98LkxZm7e1+sbYz4HsNbmAEOBBcAa4ENr7bryhS0i\nIp7gkeYgERGpnFyZMXy6yWPGmFBjzIfGmI3GmO+NMY3diNMbSnEu/m6MWWOMWZHXSd/IjTi9obST\nCvMmKeaeaoJiZVeac2GMuSXvd2OVMeZ9b8foLaX4G2lkjFlojFme93fS1Y04K5oxZoIxJsUYs/IU\n24zOu26uMMa0LtWOrbVeLTiJ5zcgHqgCrACaF9vmXuCNvPXeOE1IXo/VR87F5UB43vrgQD4XedtV\nx5mc+B3Q1u24Xfy9OBv4CYjK+1zL7bhdPBdvAffkrbcAkt2Ou4LOxaVAa2DlSf69KzA7b/1i4IfS\n7NeNmkBpJo9dD+RPg/sY6OLF+LzptOfCWvu1tTYz7+MP+O8ci9JOKnwWGAUc82ZwXlaac3EXMMZa\newjAWrvPyzF6S2nORS4QlbceA+zwYnxeY639Bkg9xSbXA5Pztl0CRBtj6p5uv24kgdJMHivYxjod\ny2nGmBreCc+rznQi3SBgboVG5J7Tnou86m2ctXaONwNzQWl+L84BzjXGfGOM+c4Yc7XXovOu0pyL\nEcDteYNUPgfu81Jsvqb4udpBKW4aPTlEtLRKM3ms+DamhG38Qakn0hljbgPa4TQP+aNTngtjjAH+\nA/Q/zff4g9L8XoTgNAl1AhoDi40x5+XXDPxIac5FX2CitfY/efOS3sd5TlmgKdPEXDdqAttxfmnz\nxQE7i22zDWgEYIwJxmn3PFU1qLIqzbnAGHMl8CTQI69K7I9Ody4icf6wk4wxyTjDkmf4aedwaX4v\ntgMzrLW51trNwK9AM++E51WlOReDgGkA1tofgHBjTC3vhOdTtpN33cxT4vWkODeSQMHkMWNMKM7k\nsZnFtplF4R3fzcBCL8bnTac9F8aYNsBY4Dpr7X4XYvSWU54La+0ha20da21Ta20CTv9ID1vCLHU/\nUJq/kc+AzgB5F7xmwCavRukdpTkXW4ArAYwxLYAwP+4jMZy8BjwT6AcFT2pIs9amnG6HXm8Ostbm\nGGPyJ48FAROsteuMMSOAZdbaz4EJwHvGmI3Afpz/eL9TynPxbyAC+CivSWSLtbane1FXjFKeixO+\nBT9tDirNubDWzjfGXGWMWQNkA4/4Y225lL8XjwBvG2P+jtNJ3P/ke6y8jDFTgESgpjFmK/A0EIrz\nGJ5x1to5xphuxpjfgMPAwFLtN284kYiIBCC9XlJEJIApCYiIBDAlARGRAKYkICISwJQEREQCmJKA\niEgAUxIQEQlgSgIiIgHs/wGthhk4tXDXRAAAAABJRU5ErkJggg==\n",
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"x_new = np.empty(n)\n",
"\n",
"for i in range(n):\n",
" x_new[i] = b[i]\n",
" for j in range(i):\n",
" x_new[i] -= A[i,j]*x_new[j]\n",
" for j in range(i+1, n):\n",
" x_new[i] -= A[i,j]*x[j]\n",
" \n",
" x_new[i] = x_new[i] / A[i,i]\n",
"\n",
"x = x_new\n",
"pt.plot(mesh, x)\n",
"pt.plot(mesh, x_true, label=\"true\")\n",
"pt.legend()"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"### And now Successive Over-Relaxation (\"SOR\")"
]
},
{
"cell_type": "code",
"execution_count": 92,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [],
"source": [
"x = np.zeros(n)"
]
},
{
"cell_type": "code",
"execution_count": 106,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [
{
"data": {
"text/plain": [
"(-1.3, 1.3)"
]
},
"execution_count": 106,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAYEAAAD7CAYAAACMlyg3AAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAIABJREFUeJzt3Xd4FNX6wPHvSQ8hBRJDTxEQuCiCCAIiBlBKlCKKwFVB\nREURuBau+NOrAlbUa7kXFBFEFBERRFHwKoihiIBSLBCKAiH0UJKQQCDl/P446SRkSTY7W97P88yz\ns7uzM+9ONvPOKTNHaa0RQgjhmbysDkAIIYR1JAkIIYQHkyQghBAeTJKAEEJ4MEkCQgjhwSQJCCGE\nB/OxOoDilFLSX1UIISpBa60q8zmnKwlorWXSmmeffdbyGJxlkn0h+0L2xYWnqnC6JCCEEMJxJAkI\nIYQHkyTgpOLi4qwOwWnIvigi+6KI7Av7UFWtT7InpZR2pniEEMIVKKXQlWwYdqreQUIIUVxMTAxJ\nSUlWh+E0oqOj2bt3r13XKSUBIYTTyj/DtToMp1He/qhKSUDaBIQQwoNJEhBCCA8mSUAIITyYJAEh\nhPBgkgSEEKKSYmNjWbFihdVhVIl0ERWFtNZknMsgNSuV1KxU0s6mFc6nn00n41wGGecyyDyXaeaz\nzfMz2WfIyskiKyeLs7lnzWOOeczJyyFX55Kn88jNyyVX55KbZ557KS+8lBfeXt7mUZlHX29fAnwC\nCPAJINAnkEDfwML5IL8gQv1DCQsII9Q/lNCAovnwGuHUCapDZFAkIf4hKFWpzhJC2EVubi7e3t5W\nh1Eh6SLqpgoO6Eczj3I08ygpp1MK54+fPs7xM/nT6aLH1KxUAnwCCAsIMwfWYgfYUP9QavrVPG8K\n8guihm8NAnwC8Pf2Lzx4B/gE4O/jj4+XT+EB3tvLu/BA76W80OjChFA8UWTnZRcmlTPZZziTc6Zw\nPuNcBmln00jLSit8TD2bSlpWGsdOH+No5lGOZB4hJy+HyKBIIoMiqRNUh4YhDYkOjSY6LJro0Gii\nQqOoH1wfby/n/yf1ZM7cRXTo0KF8/PHH+Pv74+Pjw9NPP8348eOZMWMGEydOJDY2lokTJ3LnnXeS\nnJxc+LnY2FhmzpxJt27d0FozefJkZsyYQVpaGt27d2fatGmEhYWVuc3q6CIqJQEXk52bzeGMwxw8\ndZDDGYc5lHGIwxmHz5s/mnkUL+VVeCCMDIokskYkETUiqFuzLi0jWxIeGE54jfDCx1oBtfD19nXs\nF6qmY/Dp7NMcyTjC0cyjHM44zP70/SSlJfHrjl9JSksiKTWJ42eO0yC4AU3Dm9I8vDnNIprRPKI5\nzcKbUT+4vpQkxAV9+OGHrF69mvfff5+uXbuSlJTE+PHjWbVqFdu3b8fLy4t169Zd8Hf01ltvsXjx\nYlavXk1ERARjx45l1KhRzJ0712HfQ5KAk9Bak342nf3p+zlw6gD70/eb+fQDHMw4aB5PHeT4meNE\nBkVSr2Y96gfXp27NutStWZcr61xJz8Y9C59HBkUS5Bdk9deyTA3fGsTWiiW2Vmy5y5zNOUtyejI7\nj+9k+7Ht/Hr4Vz7d+inbj23nTPYZmkc0p3Xd1rSp24Y29drQqk4ravjWcOC3ELZQE+2TrPWzlStx\nFD8zV0oxceJEAgMDbfrs9OnTmTp1KvXq1QPgmWeeITo6mjlz5uDl5ZgmW0kCDnI6+zTJacnsS9tH\ncnr+Y1oyyelm2p++H4CGIQ0LpwbBDWhTrw03Bd9E/eD6NAhuQGRQpFRh2Im/jz9NajehSe0mxDeN\nL/HeyTMn2ZayjS2Ht7Dx0EZmbJ5BYkoisbViaVO3DVfXv5prG11L67qtHV96EiVU9uBdXRo2bGjz\nsklJSdxyyy2FB3ytNb6+vhw5cqQwMVQ3SQJ2oLXmaOZRktKS2Je2j6RU87gvvWg+MzuThiENiQqN\nolFII6JCo+jQsAMDQwfSKKQRjUIbEeIfYvVXEflqBdbi2qhruTbq2sLXzuWeY1vKNjYf2syGAxt4\nf/P77EndQ7v67egc1ZnOUZ3p0LCD/B09SFlVPcVfCwoK4vTp04XPc3NzSUlJKXweFRXF+++/T8eO\nHas30AuQJGCD3LxcDmUcYm/qXpJSk8xjWhJJaWZ+X9o+gnyDSjQ6xoTF0CW6C1GhUUSHRXNJjUuk\njtnF+Xn70bpua1rXbc3wNsMBSM1K5afkn1izbw0vrH6BjQc30jyiOT0b96Rnk550bNhRSgpurG7d\nuuzevbuwkbd0o+1ll11GVlYW33zzDTfeeCMvvPAC586dK3x/5MiRPPnkk8yePZuoqChSUlL46aef\n6Nu3r8O+gyQBSh7k95zcw97UvWZKMwf95PRkwgPDCw/yMWExtKnbhv7N+xMTFkNUaBQ1/Wpa/TWE\nBcICwujdtDe9m/YGTGlh3f51fPvntzz67aP8eeJP4mLiCpPCpbUutThiYU9PPPEEY8aM4fHHH+ep\np54670QvJCSEt99+mxEjRpCXl8fjjz9eorroH//4BwA9evTg0KFDREZGMmjQIIcmAbt0EVVKzQRu\nBo5orVuVs8x/gN5AJnC31npLGctUSxdRrTWHMw6bg3zqnqKDfZo52CenJVM7sDYxYTGFU2xYLNFh\n0YUH+QCfALvHJdxfSmYKy3Yv49u/vuXbP78lvEY4t7a4lVtb3EqrOq2kdFgBZ+4iaoXq6CJqryTQ\nGcgAPiwrCSilegOjtdY3KaWuAd7SWncoY7lKJQGtNcfPHC88i9+TuqfwIL/n5B6S0pII9gsmtlZs\n4QG++GNUaBSBvra15gtRWXk6jw0HNrBw20IWJC7AW3mbhPC3W2lXv50khDJIEijJaZNAfhDRwFfl\nJIFpwA9a60/znycCcVrrI6WWKzcJpJ9NZ8/JPSXO5PekFj338fIhNsx0CYwJjTnvgO/J3SWF89Fa\ns/nw5sKEkJWTxZDLh3B367tpHtHc6vCchiSBklw5CXwFvKS1Xpv/fDnwuNZ6U6nl9JKdS0qe0ecf\n5LNyss47yBc+D4shLKDsK+yEcHZaa/44+gcf/vohc36fQ0xYDMNbD2dQy0GEBoRaHZ6lJAmU5MpJ\n4GvgxVJJ4J9a682lltMNb7qUyOAwagXWotN1nbi5x83EhMVI7xrhEXLycvj2z2+ZtWUWy3cvJ75p\nPMNbD6f7pd3xUp53v0dJAiUV7I+EhAQSEhIKX584caLTJ4HS1UHbgevLqg5q106zdi34SL8l4eGO\nnz7OJ398wszNMzmdfZox7ccw7MphBPsHWx2aw0gSKMnZh5dU+VNZFgNDAZRSHYDU0gmgQHAwvP66\nHaMSwkWF1whndPvRbLp/EzP7zmRl0kpi3orhkf89wl8n/rI6POEm7NU7aC4QB4QDR4BnAT9Aa62n\n5y8zBeiF6SI6vHR7QP4yevduTbt28OOP0KxZlUMTwq3sS9vH2z+/zczNM+nQsAMPX/Mw3WK7uW1V\naUxMDElJSVaH4TSio6PZu3fvea87RZuAPRT0DpoyBT75BFatAhe4HbcQDnc6+zRzf5/L6z+9TrB/\nMM90eYb4pvFumwzEhbldEsjLg7g4GDAAHn7Y6qiEcF55Oo+F2xby/Orn8fHy4ekuT9O3WV+PbET2\nZG6XBAB27YKOHWH9emjc2OLAhHByeTqPxTsW89yq58jOzebpLk9z699ulWTgIdwyCYBpIF68GFas\nAAfdWlsIl6a1ZumupTy36jkyszOZfMNkejfpLdVEbs5tk0BuLnTuDHfdBaNGWRiYEC5Ga81XO7/i\n8WWP0zCkIa/e+Cpt6rWxOixRTdw2CQAkJsJ118Evv0BMjDVxCeGqsnOzmbFpBpNWTaJH4x483/V5\nGoU2sjosYWfOcp1AtWjRAsaNg/vug7w8q6MRwrX4evvyYLsH2TF6B41CGtH63dY8+f2TpJ9Ntzo0\n4SScPgmASQKnTsHbb1sdiRCuKcQ/hOe7Pc+vD/zKwVMHaTG1BfO3zpercYXzVwcV2LkTOnWSi8iE\nsIe1yWsZ+fVIGoY0ZGr8VBnsxsW5dXVQgcsug0mTTCNxTo7V0Qjh2jo16sSm+zfRNaYr7d9rz4ur\nX+Rc7rmKPyjcjsuUBAC0ht69TYngmWccGJgQbmxv6l5GLx3N7pO7mXbzNLpEd7E6JHGR3Lp3UGkH\nDsBVV8GSJXD11Q4KTAg3p7Vm0fZF/ON//yC+STz/7vlvGTfbhXhEdVCBBg3grbdMtdCZM1ZHI4R7\nUEoxoMUAto7ayrm8c7Se1pq1yWutDks4gMuVBAoMGQJ16sCbb1ZzUEJ4oEWJi3hwyYPc0+YeJsRN\nwM/bz+qQxAV4VHVQgRMnoFUrmD0bunev5sCE8EBHMo5w31f3sS9tH3MGzOHyyMutDkmUw6OqgwrU\nrg0zZ8Lw4ZCaanU0QrifOjXr8OXgLxl7zVi6zu7Ka2tfI0/LFZvuxmVLAgVGjYK0NJgzB+QeWUJU\njz0n9zD0i6EE+QYxZ8AcImpEWB2SKMYjSwIFXnsNfv0VPvjA6kiEcF+xtWL5YdgPtKrTirbT27J+\n/3qrQxJ24vIlAYCtW80gNKtWmXsNCSGqzxfbv+D+r+7n6S5PM7r9aLlNtRPwyIbh0mbMgP/8xwxC\nExho58CEECX8deIvbvvsNppHNGf6zdMJ9g+2OiSP5tHVQQVGjIC//Q0efdTqSIRwf41rN2btPWup\n6VuT9jPas/XoVqtDEpXkNklAKZg+HZYtg88+szoaIdxfoG8g7/V9j/HXjidudhyLEhdZHZKoBLep\nDirwyy8QH2+qhWJj7RSYEOKCfjn4C/3n9WfsNWP5Z6d/SjuBg0mbQClvvAHz5sGaNeDra4fAhBAV\n2p++nz6f9OGqulfxzs3vyFXGDiRJoBStoU8f00bwyit2CEwIYZOMcxnc8fkdpJ9NZ+HtC6kdWNvq\nkDyCNAyXopS5buCTT2DpUqujEcJz1PSryee3f87V9a6mw4wO7Dq+y+qQRAXcMgkARESYKqG774Zd\n8jsUwmG8vbx5tcerPH7t41w36zoS9iZYHZK4ALesDipu2jT4739h3ToIlq7MQjjUij0rGLJwCFPj\np3Lb326zOhy3JW0CF6A1jBwJKSmwcCF4uW3ZRwjn9OvhX4mfG8+kuEmMuGqE1eG4JUkCFTh7Frp2\nNUNTPv203VcvhKjAruO76DGnBw+1e4hxncZZHY7bkSRgg0OHoF07eOcd03NICOFY+9P30+OjHvRv\n3p8Xur0g1xLYkSQBG61bB337mhvNNW9ebZsRQpTj2Olj9JrTi/YN2jMlfgpeSupn7UG6iNqoQwd4\n8UXo39+MQSCEcKyIGhGsGLaCbSnbuPPzO8nOzbY6JI/nUSWBAg89BPv2wZdfSkOxEFY4k32GQQsG\nodEsGLgAfx9/q0NyaVIddJHOnYMbbjBtBP/+d7VvTghRhuzcbIYsHEJOXg6fDfwMX2+5x0tlWV4d\npJTqpZTarpTaqZQaX8b7w5RSR5VSm/Kne+yx3cry84MvvjBXE7/xhpWRCOG5fL19mXvrXDSaIQuH\nSNWQRapcElBKeQE7ge7AQeBnYLDWenuxZYYBbbXWYytYl0NKAgWSkuDaa80QlYMHO2yzQohizuac\nZcD8AQT7BTNnwBx8vHysDsnlWF0SaA/s0lonaa2zgXlAvzKWc7r+YNHRpjQwdiz88IPV0Qjhmfx9\n/Fl4+0JOZp3k7i/uJjcv1+qQPIo9kkADILnY8/35r5U2QCm1RSk1XynV0A7btYtWreDTT2HQIPjt\nN6ujEcIzBfgE8MWgLziUcYh7v7qXPJ1ndUgewx7lrrLO8EvX6SwG5mqts5VSI4HZmOqj80yYMKFw\nPi4ujri4ODuEeGFdu5rxiW+6CX78EaKiqn2TQohSAn0DWTx4MfFz4xn51Uje7fOuXEdQjoSEBBIS\nEuyyLnu0CXQAJmite+U/fwLQWuvJ5SzvBZzQWoeV8Z5D2wRKe+MNeO89MxhNbbkNuhCWOHX2FL0+\n7kXbem15q9dbcmWxDaxuE/gZaKKUilZK+QGDMWf+xQOsW+xpP2CbHbZrd488Yoam7NsXTp+2Ohoh\nPFOwfzBL/76UlUkreXnNy1aH4/aqnAS01rnAaOA7YCswT2udqJSaqJS6OX+xsUqpP5RSm/OXvbuq\n260ur7wCjRubqqHMTKujEcIzhQaE8s0d3/Duxnf5YMsHVofj1jzyYrGK5ObCvffC7t2wZAnUrGl1\nREJ4pu3HthP3QRyz+s2id9PeVofjtKyuDnI73t4wcyY0aWJuP33qlNURCeGZmkc0Z9GgRQz9Yigb\nDmywOhy3JEmgHF5eppG4RQuTCNLTrY5ICM/UsVFHZvadSb95/WTM4mogSeACvLzM8JSXXw69esmd\nR4WwSt9mfZkUN4mec3pyOOOw1eG4FUkCFfDygrffhtatoWdPSQRCWOW+tvcx7MphxH8cz6mzUkdr\nL9IwbCOtze0l1q+Hr7+GyEirIxLC82itGfn1SA6eOsiXg7/E28vb6pCcgjQMO4BS5qrinj2hY0dI\nTLQ6IiE8j1KKqfFTyczO5Mnvn7Q6HLcgSeAiKAXPPWcGq4+Lk5vOCWEFX29fPhv4GZ9t+4w5v82x\nOhyXJ0mgEu6+Gz75xNx++oMPrI5GCM8TUSOCxUMW88i3j7B+/3qrw3Fp0iZQBYmJ5sriO+6ASZNM\nSUEI4TiLdyxm1JJRrL93PQ1Cyrp5sWeQ4SUtdPQo9OsHsbHw/vsQEGB1REJ4lpdWv8Tn2z9n1d2r\nCPQNtDocS0jDsIUiI2HFCsjJMbek3rfP6oiE8CxPdH6CprWbMmLxCFztJNIZSBKwg8BAmDcPBgww\ng9cvWmR1REJ4DqUUM/vOZNeJXUz+scw72IsLkOogO1u/HoYMMbea+Pe/pXpICEc5kH6Aa2Zcw/Q+\n04lvGm91OA4l1UFO5JprYNMmSEkx83I9gRCO0SCkAfNum8fwL4ezL03qZW0lSaAahIWZcYtHj4Yu\nXWDWLHPFsRCienWO6sy4juMYtGAQ53LPWR2OS5DqoGq2dasZxL5FC3jrLahf3+qIhHBtv/5qTrSi\no8t+P0/n0X9ef5rUbsLrPV93bHAWkeogJ9ayJfz8MzRrBldeaW49kZtrdVRCuK5x4+Ddd8t/30t5\n8UH/D1i0fRGLEqWXRkWkJOBAiYkwapQZm2DaNNOTSAhhu4wMqF3b3Lblu+8uvOyGAxu4ee7N/DTi\nJxrXbuyQ+KwiJQEX0aKFuabg4YfNYPajR0NqqtVRCeE6vv/elKp/+aXidrb2Ddrzry7/4vYFt5OV\nk+WYAF2QJAEHUwruusu0FeTkwN/+ZoayzM62OjIhnN/SpTB8OAQFwZ49FS8/pv0YYsNiefTbR6s/\nOBclScAitWubKqFFi2DuXFNK+PBDkxiEEOfT2iSB3r3h6qtNaaAiBReSfffXd3zy+yfVH6QLkiRg\nsWuuMUXcGTPM1LKlSQrSeCzc3caNMPkiLvD94w/w8YHmzaFtW9uSAEBoQCifDfyMsf8by58n/qxc\nsG5MkoCTiIuDlSvNUJZTpsAVV8D8+ZIMhHvaswf69IHnn4cjR2z7zNKlEB9vqlRtLQkUaFOvDf+6\n7l/c+fmd5ORJcbs4SQJORCno3h1+/BHeeANefx2aNIGXXzZ3KxXCHZw4Yap0nnoKbr0V5tg4LkxB\nEgBTEti4EfLybN/u1XljCPEP4cXVL1580G5MkoATUsoMY7lunSkN7NwJl11mxi1Ys0auPhau6+xZ\nuOUWUwp46CHTyPv++xX/plNTze1YunY1zy+5BGrVgj9trN05dQrirvdiXJNZTP15KhsObKjaF3Ej\nkgScXLt25p9kzx4zP2IEtGplqo2OHbM6OiFsl5dnRuWrU6eoLaBLF5MYfv75wp9dtgw6d4YaNYpe\nu/pqUxqwxcqVptPFn5sbMKX3FO78/E4yz2VW6nu4G0kCLqJWLXN9wfbt5vYTK1dC48Zwww2ml9Hh\nw1ZHKDyJ1qZaJzHRjLW9aFHFJyVPPgnJyaYXnFf+kUcpUxqYNevCny1eFVTgYtoFli+Hpk1h7VoY\n2HIgHRp2YNx342z7sJuTK4Zd2OnT8O23sGCB+Se54gq47TZT3G7UyOrohKvJzTUH9iNHzHT0aNnz\nhw+b54GBULeuObOvUcMcYLt0gTvvNNU9xc/ap00z7Vxr10J4eMntJiebW6ocOGDWWVpenrnn1po1\npo2swPLl8Nxz5oSoIi1bwvjxMGEC7N4NaVlpXDntSqbET+Hmy26u1P5yJjK8pCAry/xTLFgAX31l\n6ky7djVTXJwZAU14Fq1NXfrRo+bW5ikp588Xn06cgJAQc1AvPkVGlnxet655rfRYGadOmRLBxx/D\nhg3Qv79JCJmZMHKkOYg3LufuDb16wdCh8Pe/n//exo3m9R07Sr5+8qS5idzJk+DtXf5+OHjQnCAd\nOWLi3roV6tWDVUmrGLxgMFse2EJkkGv/g0gSECXk5sJvv5li+g8/wOrV0LChSQhdusBVV8Gll5qi\nuHAdp0/D8eOm2uXYsZLzKSnnPx4/bs7GL7mk5BQZaR4LDvAFU0SE6YdvD4cOmdH2Pv7YVGF+/725\nJqY88+fD9OnmRKa0554zCeqNN85/r0kTc9LTokX56/7oI/jyS3OCdNNNpl1twADz3hPLnyDxWCJf\nDPoC5cL/EJIExAXl5MDmzSYhrFlj5k+dgtatoU2boql5c/D1tTpa93f2rDl7PXHi/On48aKp9PO8\nPHPwDg83B+zwcDMVHNwjIko+hoeDv7/V39bcEqWi31VWljlR+eUXiIkp+V6nTjBxItx44/mfGzzY\nHNjvuqv8dQ8datbxwAPwwgtm37/2mnnvXO45rplxDQ+1e4h7r7r3or6XM5EkIC7asWMmGRSf9u41\nbQlNmpw/RUWVXV/riXJyzJ1g09JMdUvpx9RUc6ApeCyYCp6fO2duG1K7tmnwL3gsOKiHh5vXis9H\nRJizehc+Wa3QmDHmez77bNFrx46ZUmtKStkJ7bXXTJvCW2+VvU6toUEDUxpu3NicCD31lGmbKLD1\n6Fau/+B6Nt6/keiwcgYpcHKSBIRdnD1rEsGff5acdu2C/ftNHXC9eudPERFmkI/SU3BwUS8Qq+Tk\nwJkzpirl9GlTP136MSPDlIxOnSqaL3hMTy864BfMZ2WZ7xYWBqGhZT/WqlVyKnitdm1z8zN3PphX\n1ubNplPD7t1Fv5u5c80ofV9+WfZnEhLMQf3HH8t+f9s2U1LYvdvs88xMU/V1/HjJNo2XVr/EyqSV\nfHPHNy5ZLVSVJGCnGkDhDvz9zW16mzU7/z2tzVnsoUMlp/37TftD6TPg1FTzDxcQYEoQpaeAANOY\n5+NjHgsmHx9zAMjLM9vUumg+L8+0d2RnF03nzpWcP3PGHKQLJq3NtoKCzJl0WY81a5qDes2a5iDe\nsGHRayEhZgoNLZqXg3j1aNPGJMuEBOjWzbxWVtfQ4q66yow0lpNTdnvGsmWmG3XB3ysoyLQfbNpk\nqogKjOs0jvnb5jPntzncdeUF6pbckF1KAkqpXsCbmOsOZmqtJ5d63w/4EGgLHAMGaa3PGwlaSgLu\nJTfXHJTLmrKyzPs5Oeax9OTlZf5xSz96e5v65YLJz69o3t/fHPALEk9AgP0aOoVj/Oc/pmfRnDnm\nd1CnjjlgR0WV/5lmzUyj7xVXnP9enz6mveD224teGzvWrG9cqcsENh3aRO+Pe/PbA79Rp2Yd+3wh\nB7G0Okgp5QXsBLoDB4GfgcFa6+3FlnkQuEJrPUopNQi4RWs9uIx1SRIQwoMdP27q7vfuNb2K7rsP\nfv/9wp+54w5ztj98eMnXs7NNVeVff5nHAvPmmd5In39+/rrGLxvP3rS9fHrbp1X+Lo5k9chi7YFd\nWuskrXU2MA/oV2qZfsDs/PkFmIQhhBAlhIebXkCfflpxVVCB8q4c3rDBJJTiCQBMNdCPP5Z9v6IJ\ncRPYdGgTX24vpxHCDdkjCTQAkos935//WpnLaK1zgVSlVG07bFsI4WbuucfcL6uqSWDZsrK7lTZq\nZKoPd+8+/71A30Bm9JnBQ0sfIjXLM8Z+tUeNaVlFkNI5tvQyqoxlAJgwYULhfFxcHHFxcVUITQjh\nanr0MNVAGRklG2/L06aNqTI6d860ERVYvrxkd9MCSpn1rl1b9hXM18dcz82X3czjyx5nep/plf8i\n1SghIYGEhAS7rMsebQIdgAla6175z58AdPHGYaXUN/nLrFdKeQOHtNbnXactbQJCCICnnzZdk+fN\ns235li1NY3KbNuZ5erq5PqDgHkelvfmmuQ3FO++Uvb60rDQuf+dyZvefTbfYbpX7Eg5kdZvAz0AT\npVR0fi+gwcDiUst8BQzLnx8IrLDDdoUQburZZ81wq7YqXSW0cqW5TUV5FzgWlATKExoQytvxb3P/\nV/dzOvu07YG4oCongfw6/tHAd8BWYJ7WOlEpNVEpVXB7vplAhFJqF/Aw8ERVtyuEcF8+PuZaDVuV\nTgLLl5seQ+Vp3dr0GkpPL3+ZPs360K5BO579oYw6JTciVwwLIVzeTz/B6NFFg8y0bAmzZ5vkUJ4u\nXUy1U1mNxwWOZh6l5dstWTF0BVfUKeNCBCdhdXWQEEJY6sorzQA3Z8+acQkOHy5qHyjPtddeuEoI\nIDIokgnXT2D0N6Nx1xNUSQJCCJdXo4YZOez3381tq7t1u/AYA1BxuwDAH3/AmdUPcOrsKeb+Ptd+\nATsRuaheCOEWCtoF1q69cHtAgY4dzS0lcnPLThiZmTBwICQlefP5hqmM+O42+jTrQ4h/iP2Dt5CU\nBIQQbqFtWzNg/fLlF67nLxARYe5NtG1b2e8/8ohJLP36wV8rO9KzcU8mJEywa8zOQJKAEMItXH21\nuR9QQIAZg8AWBbeQKG3BAlixAqZOhWHDTCPzyze8zEe/fcTvRyq4mZGLkSQghHALrVqZq4xtqQoq\nUFa7wL59MGqUGcsgJMSUKvbvh+P7IpkYN9HtGoklCQgh3EJAgOkl1KOH7Z8pnQRycsxdSR97DNq3\nN695e8Odd5rSwMi2I92ukViuExBCuI1Dh8zIYRX1DCqQl2fuXLp9u2kfmDjRDEX53XclR8X74w/o\n1QuSkmBkfi14AAAN4klEQVTDwZ+47bPbSHwo0WkaieU6ASGEwAx3amsCAHOg79DBXGy2erW5l9CH\nH54/LOrll5sk8f330LFRR3o17uU2jcSSBIQQHq1TJ1iyxFT5zJgB9euXvdywYSZBgHs1Ekt1kBDC\no33/vWlMHjPGDG9ZnpQUc0Havn2mwXjqhqksTFzI90O/t3xweqkOEkKISurY0SSAV1658HKXXAJx\ncab7KMDIq0dyOOMwS3YtqfYYq5OUBIQQwkaLFpmxCFauNM+X7FzCuGXj+O2B3/D19rUsLikJCCGE\nA9x0k7nCeM8e8zy+aTwNghvw3qb3rA2sCiQJCCGEjfz8YPDgogZipRSv9XiNSSsnkZaVZm1wlSRJ\nQAghLkJBL6GCmuvWdVvTu2lvXlrzkrWBVZIkASGEuAht24K/P6xZU/Ta812f571N75GUmmRdYJUk\nSUAIIS6CUiWvGQBoENKA0e1G8+SKJ60LrJKkd5AQQlykAwfgiivMY8Fg9hnnMmg2pRmLBi2ifYP2\nDo1HegcJIYQDNWgA7dqZLqMFavrVZFLcJB777jGXusuoJAEhhKiExx4zA8/8739Fr93d+m7Sz6az\naPui8j/oZCQJCCFEJfToAQsXwj33mKuNtQZvL29eu/E1xi8fz7ncc1aHaBNJAkIIUUmdO8P69TB/\nvhmH4PRpuLHxjTSt3ZRpv0yzOjybSMOwEEJU0ZkzcP/9ZtyBL76Ak/5b6P1xb/4c8ydBfkHVvn1p\nGBZCCAsFBpouo3fdZcYnSNvRmi7RXfjvhv9aHVqFpCQghBB2tGyZGZtgwpTtPLP3OnaN2UVYQFi1\nbrMqJQFJAkIIYWfffQdjx0KHl4YTFdaISV0nVev2pDpICCGcyI03mrGLW6c9y9Sfp5KSmWJ1SOWS\nkoAQQlSD5cvhoYeg+78fooZfIK/1eK3atiUlASGEcDLdu5vRyFoce4pZW2Zx8NRBq0Mqk5QEhBCi\nmnz/PTz4IPT57z85k5PJ2ze9XS3bkZKAEEI4oW7doE4daHxoPJ9u/ZQ9J/dYHdJ5pCQghBDVaMUK\nGDkSBk97luRTSXzQ/wO7b0NKAkII4aS6doV69aBB8qMs3bWUxJREq0MqoUolAaVULeBTIBrYC9yu\ntT5voE2lVC7wK6CAJK11/3LWJyUBIYTb+eEHc1uJe2ZOZvORjcwfON+u67eyJPAEsFxr3QxYAfxf\nOctlaq2v0lq3KS8BCCGEu4qLg/r1Ifyv0azZt4Yth7dYHVKhqpYEtgPXa62PKKXqAgla6+ZlLHdK\nax1sw/qkJCCEcEsJCXDvvTBy9utsOPgTnw38zG7rtrIkEKm1PgKgtT4MXFLOcv5KqQ1KqbVKqX5V\n3KYQQricuDho1AhCdoxkVdIqtqVsszokAHwqWkAptQyoU/wlQAP/uojtRGmtDyulYoEVSqnftNZl\n9pWaMGFC4XxcXBxxcXEXsRkhhHBeEybAPfcEMeb9h3lpzUt8dMtHlVpPQkICCQkJdompqtVBiUBc\nseqgH7TWLSr4zCzgK63152W8J9VBQgi31r079LktnedPXcq6e9fRpHaTKq/TyuqgxcDd+fPDgC9L\nL6CUClNK+eXPRwCdAOcoBwkhhIO98gpMnhTCvVc+xMtrXrY6nCqXBGoD84FGwD5goNY6VSnVFhip\ntb5fKdUReBfIxSSdN7TWH5SzPikJCCHc3l13QWTMcWbVbMqWB7YQFRpVpfXJeAJCCOFC9u2DNm1g\n0Hvj8QrIZEr8lCqtT5KAEEK4mCefhL8OH2FZsxZsHbWVesH1Kr0uSQJCCOFi0tPhssug2yv/oH5d\n3yqNNyBJQAghXNA778CcxftJjGvFzjE7iagRUan1yA3khBDCBd13H5zY25Brggfy5ro3LYlBkoAQ\nQljExwdefRV2zRrPtF+mkZqV6vAYJAkIIYSFbroJooIvpam+mSkbqtZLqDKkTUAIISy2eTP0GLID\n7rmOPQ/vpqZfzYv6vLQJCCGEC2vTBuKvaUZEZhfe3/y+Q7ctJQEhhHACyclwea91hI4Ywu6Hd+Hj\nVeH9PQtJSUAIIVxco0Zw9w0dyEttwKLERQ7briQBIYRwEuPGQdo343hp1as4qlZEkoAQQjiJRo1g\n4JV9SDqSypp9axyyTUkCQgjhRP5vvDdnEx7l5dWVv43ExZAkIIQQTqRpU+hdfyir/lrHjmM7qn17\nkgSEEMLJPP1EDfjlQV5d83q1b0uSgBBCOJlWraCT7yjm/jafo5lHq3VbkgSEEMIJTXw8Eu/EQfxn\n3dRq3Y5cLCaEEE6qw0072NbhOg4/sZcavjXKXU4uFhNCCDf0wqPN0Ps6MWvz7GrbhpQEhBDCSWkN\nLXuv4eT1w9k/fjveXt5lLiclASGEcENKwUsPXkv64XAW7/iqerbhTGfeUhIQQoiS8vIgJn4BwTe+\nwdbHfixzGSkJCCGEm/LyghfvuoW/Ug6w8eBG+6/f7msUQghhV4MHeVNj6yieWWL/kcckCQghhJPz\n8YF/dh/BsuQvSMlMseu6JQkIIYQLGHNvOF7bBzB5+Qy7rleSgBBCuICaNeHvTUfz7sZ3yMnLsdt6\nJQkIIYSLeGF0G84cieLjXxbbbZ2SBIQQwkXUqwedfUcz8X/2ayCWJCCEEC7kzZEDSMrYzsbkP+yy\nPkkCQgjhQlpf4celJ0by2Kf2KQ1IEhBCCBfz4q0jWX3yU06cTq3yuiQJCCGEi7mtV11CjvRm/LxZ\nVV6XJAEhhHAxSsE/Oozh451TydN5VVpXlZKAUuo2pdQfSqlcpdRVF1iul1Jqu1Jqp1JqfFW2KYQQ\nAp4c2oGczFD++83/qrSeqpYEfgduAVaWt4BSyguYAvQEWgJDlFLNq7hdIYTwaH5+ilvqj2HyD1Vr\nIK5SEtBa79Ba7wIudAvT9sAurXWS1jobmAf0q8p2hRBCwFv3D+Kw1y9VWocj2gQaAMnFnu/Pf00I\nIUQV1A0PpJ33iCqtw6eiBZRSy4A6xV8CNPCU1tqWoW7KKiWUO3LMhAkTCufj4uKIi4uzYRNCCOE5\nEhISSEhIAOCKtCw2VGFddhlZTCn1A/CY1npTGe91ACZorXvlP38C0FrryWUsKyOLCSHERXKWkcXK\nC+BnoIlSKlop5QcMBux39yMhhBCVVtUuov2VUslAB+BrpdQ3+a/XU0p9DaC1zgVGA98BW4F5WuvE\nqoUthBDCHmSgeSGEcHHOUh0khBDCxUgSEEIIDyZJQAghPJgkASGE8GCSBIQQwoNJEhBCCA8mSUAI\nITyYJAEhhPBgkgSEEMKDSRIQQggPJklACCE8mCQBIYTwYJIEhBDCg0kSEEIIDyZJQAghPJgkASGE\n8GCSBIQQwoNJEnBSCQkJVofgNGRfFJF9UUT2hX1IEnBS8gMvIvuiiOyLIrIv7EOSgBBCeDBJAkII\n4cGU1trqGAoppZwnGCGEcCFaa1WZzzlVEhBCCOFYUh0khBAeTJKAEEJ4MEuSgFKql1Jqu1Jqp1Jq\nfBnv+yml5imldimlflJKRVkRpyPYsC8eUUptVUptUUotU0o1siJOR6hoXxRb7jalVJ5S6ipHxudI\ntuwLpdTt+b+N35VScxwdo6PY8D/SSCm1Qim1Kf//pLcVcVY3pdRMpdQRpdRvF1jmP/nHzS1KqdY2\nrVhr7dAJk3j+BKIBX2AL0LzUMg8Cb+fPDwLmOTpOJ9oX1wMB+fMPePK+yF+uJrASWAtcZXXcFv4u\nmgAbgZD85xFWx23hvngXGJk/3wLYY3Xc1bQvOgOtgd/Keb83sCR//hpgnS3rtaIk0B7YpbVO0lpn\nA/OAfqWW6QfMzp9fAHR3YHyOVOG+0Fqv1Fpn5T9dBzRwcIyOYsvvAuA5YDJw1pHBOZgt++I+YKrW\nOh1Aa33MwTE6ii37Ig8IyZ8PAw44MD6H0VqvAU5eYJF+wIf5y64HQpVSdSparxVJoAGQXOz5fs4/\nsBUuo7XOBVKVUrUdE55D2bIvihsBfFOtEVmnwn2RX7xtqLVe6sjALGDL7+IyoJlSao1Saq1SqqfD\nonMsW/bFROAupVQy8DUwxkGxOZvS++oANpw0+lRbOOUrqy9r6X6qpZdRZSzjDmzZF2ZBpe4E2mKq\nh9zRBfeFUkoBbwDDKviMO7Dld+GDqRLqAkQBq5VSLQtKBm7Eln0xBJiltX5DKdUBmAO0rPbInI/N\nx5PirCgJ7Mf8aAs0BA6WWiYZaASglPLG1HteqBjkqmzZFyilbgD+D+iTXyR2RxXti2DMP3aCUmoP\n0AH40k0bh235XewHvtRa52mt9wI7gKaOCc+hbNkXI4D5AFrrdUCAUirCMeE5lf3kHzfzlXk8Kc2K\nJPAz0EQpFa2U8gMGA4tLLfMVRWd8A4EVDozPkSrcF0qpNsA0oK/W+rgFMTrKBfeF1jpdax2ptb5U\nax2LaR/po7XeZFG81cmW/5EvgG4A+Qe8psBuh0bpGLbsiyTgBgClVAvA343bSBTll4AXA0MB8ktE\nqVrrIxWt0OHVQVrrXKXUaOA7TBKaqbVOVEpNBH7WWn8NzAQ+UkrtAo5j/vBux8Z98QoQBHyWXyWS\npLXub13U1cPGfVHiI7hpdZAt+0Jr/a1SqodSaiuQA4xzx9Kyjb+LccB7SqlHMI3Ew8pfo+tSSs0F\n4oBwpdQ+4FnAD9Ba6+la66VKqXil1J9AJjDcpvXmdycSQgjhgeSKYSGE8GCSBIQQwoNJEhBCCA8m\nSUAIITyYJAEhhPBgkgSEEMKDSRIQQggPJklACCE82P8DeNVn3NiwUnIAAAAASUVORK5CYII=\n",
"text/plain": [
""
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"x_new = np.empty(n)\n",
"\n",
"for i in range(n):\n",
" x_new[i] = b[i]\n",
" for j in range(i):\n",
" x_new[i] -= A[i,j]*x_new[j]\n",
" for j in range(i+1, n):\n",
" x_new[i] -= A[i,j]*x[j]\n",
" \n",
" x_new[i] = x_new[i] / A[i,i]\n",
"\n",
"direction = x_new - x\n",
"omega = 1.5\n",
"x = x + omega*direction\n",
"\n",
"pt.plot(mesh, x)\n",
"pt.plot(mesh, x_true, label=\"true\")\n",
"pt.legend()\n",
"pt.ylim([-1.3, 1.3])"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": false,
"jupyter": {
"outputs_hidden": false
}
},
"outputs": [],
"source": []
}
],
"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.3"
}
},
"nbformat": 4,
"nbformat_minor": 4
}