{ "cells": [ { "cell_type": "markdown", "id": "5d32decb", "metadata": {}, "source": [ "# Métodos iterativos \n", "#### https://meet.noysi.com/metodosnumericos1" ] }, { "cell_type": "markdown", "id": "9d51f4b8", "metadata": {}, "source": [ "Vamos a aplicar el método de Jacobi para resolver el sistema\n", "\\begin{cases}\n", "3 x - y + z = 3,\\\\\n", "x + 5 y - z = 2,\\\\\n", "x + y + 4 z = 1.\n", "\\end{cases}\n" ] }, { "cell_type": "markdown", "id": "da65535b", "metadata": {}, "source": [ "Definimos las matrices del método:" ] }, { "cell_type": "code", "execution_count": 1, "id": "3f720df4", "metadata": {}, "outputs": [], "source": [ "A = matrix(RDF,[[3,-1,1],[1,5,-1],[1,1,4]])\n", "M = matrix(RDF,[[3,0,0],[0,5,0],[0,0,4]])\n", "N = M - A\n", "b = vector(RDF,[3,2,1])" ] }, { "cell_type": "markdown", "id": "5fd0f394", "metadata": {}, "source": [ "En este caso está justificado calcular la inversa, por ser inmediata y no añadir mucho error." ] }, { "cell_type": "code", "execution_count": 2, "id": "0f5fe495", "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ], "text/latex": [ "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left(\\begin{array}{rrr}\n", "0.3333333333333333 & 0.0 & 0.0 \\\\\n", "0.0 & 0.2 & 0.0 \\\\\n", "0.0 & 0.0 & 0.25\n", "\\end{array}\\right)$$" ], "text/plain": [ "[0.3333333333333333 0.0 0.0]\n", "[ 0.0 0.2 0.0]\n", "[ 0.0 0.0 0.25]" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Se podría calcular directamente\n", "Mi = matrix.diagonal([k^-1 for k in M.diagonal()], sparse=true)\n", "show(Mi)" ] }, { "cell_type": "markdown", "id": "12d43b88", "metadata": {}, "source": [ "Tomamos como aproximación inicial el vector nulo:" ] }, { "cell_type": "code", "execution_count": 3, "id": "01c8ca84", "metadata": {}, "outputs": [], "source": [ "x0 = vector([0,0,0])" ] }, { "cell_type": "markdown", "id": "c4dea0ed", "metadata": {}, "source": [ "Y aplicamos el método" ] }, { "cell_type": "code", "execution_count": 6, "id": "a27970e0", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(1.0, 0.4, 0.25)" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x1 = Mi*(N*x0+b)\n", "x1" ] }, { "cell_type": "code", "execution_count": 7, "id": "494391b7", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(1.0, 0.4, 0.25)" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x1 = M\\(N*x0+b)\n", "x1" ] }, { "cell_type": "code", "execution_count": 8, "id": "b1274220", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(1.0499999999999998, 0.25, -0.09999999999999998)" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x2 = Mi*(N*x1+b)\n", "x2" ] }, { "cell_type": "code", "execution_count": 9, "id": "e3495f57", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(1.1166666666666667, 0.17000000000000004, -0.07499999999999996)" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x3 = Mi*(N*x2+b)\n", "x3" ] }, { "cell_type": "markdown", "id": "55540f05", "metadata": {}, "source": [ "Comparamos con la solución exacta." ] }, { "cell_type": "code", "execution_count": 11, "id": "a9fadeab", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "((1.078125, 0.17187500000000003, -0.062499999999999986),\n", " (1.1166666666666667, 0.17000000000000004, -0.07499999999999996),\n", " 0.040561381811329435)" ] }, "execution_count": 11, "metadata": {}, "output_type": "execute_result" } ], "source": [ "A\\b, x3, (A\\b - x3).norm()" ] }, { "cell_type": "code", "execution_count": null, "id": "042b07c5", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "bdefdaff", "metadata": {}, "source": [ "
\n", " Ejercicio 1. \n", " a) Aplica 5 pasos del método de Jacobi para aproximar una solución del sistema y compara los valores obtenidos con la solución exacta $$\\left(\\begin{array}{rrr} 15.0 & -1.0 & 1.0 \\\\ 2.0 & -5.0 & 1.0 \\\\ 1.0 & 1.0 & -3.0 \\end{array}\\right)x = \\left(\\begin{array}{r} 1.0 \\\\ 2.0 \\\\ 3.0 \\end{array}\\right).$$\n", " \n", " b) Aplica 5 pasos del método de Gauss-Seidel para aproximar una solución del sistema anterior y compara los valores obtenidos con la solución exacta y con el método de Jacobi.\n", " \n", " c) Encuentra una matriz tal que el método de Jacobi converja y el de Gauss-Seidel no. \n", " \n", " d) Encuentra una matriz tal que el método de Gauss-Seidel converja y el de Jacobi no. \n", " \n", " e) En los apartados a) y b), aplica los métodos hasta que el error estimado con el criterio de parada sea menor que $10^{-3}$. Compara el error estimado con el error real.\n", " \n", " f) (Para valientes) Elige una matriz estrictamente diagonal dominante y un término independiente. Representa en una gráfica la estimación del error del criterio de parada para los métodos de Jacobi y Gauss-Seidel y el error exacto, para un número de pasos entre 3 y 20.\n", "
" ] }, { "cell_type": "code", "execution_count": 16, "id": "0de88ba5", "metadata": {}, "outputs": [], "source": [ "A = matrix(RDF,[[15,-1,1],[2,-5,1],[1,1,-3]])\n", "M = matrix(RDF,[[15,0,0],[0,-5,0],[0,0,-3]])\n", "N = M - A\n", "b = vector(RDF,[1,2,3])" ] }, { "cell_type": "code", "execution_count": 22, "id": "4b21ea31", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(0.10459390946502056, -0.5897758024691359, -1.161679012345679)" ] }, "execution_count": 22, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x0 = vector([0,0,0])\n", "Mi = matrix(RDF,[[1/15,0,0],[0,-1/5,0],[0,0,-1/3]])\n", "for _ in range(5):\n", " x0 = Mi*(N*x0+b)\n", "x0" ] }, { "cell_type": "code", "execution_count": 26, "id": "741b0591", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(0.10459390946502059, -0.5897758024691357, -1.161679012345679)" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x0 = vector([0,0,0])\n", "for _ in range(5):\n", " x0 = M.solve_right(N*x0+b)\n", "x0" ] }, { "cell_type": "code", "execution_count": 27, "id": "9cf7d71f", "metadata": {}, "outputs": [], "source": [ "A = matrix(RDF,[[15,-1,1],[2,-5,1],[1,1,-3]])\n", "M = matrix(RDF,[[15,0,0],[2,-5,0],[1,1,-3]])\n", "N = M - A\n", "b = vector(RDF,[1,2,3])" ] }, { "cell_type": "code", "execution_count": 28, "id": "71e7d5d2", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(0.10476883150901793, -0.5904553535027791, -1.1618955073312536)" ] }, "execution_count": 28, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x0 = vector([0,0,0])\n", "for _ in range(5):\n", " x0 = M.solve_right(N*x0+b)\n", "x0" ] }, { "cell_type": "code", "execution_count": 47, "id": "d84e381f", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "(0.29367498919688845, 0.10467975714331931)" ] }, "execution_count": 47, "metadata": {}, "output_type": "execute_result" } ], "source": [ "a = 1\n", "A = matrix(RDF,[[15,-1,a],[2,-5,1],[1,1,-3]])\n", "MJ = matrix(RDF,[[15,0,0],[0,-5,0],[0,0,-3]])\n", "NJ = MJ - A\n", "MG = matrix(RDF,[[15,0,0],[2,-5,0],[1,1,-3]])\n", "NG = MG - A\n", "max([k.abs() for k in (MJ.inverse()*NJ).eigenvalues()]),max([k.abs() for k in (MG.inverse()*NG).eigenvalues()])" ] }, { "cell_type": "code", "execution_count": null, "id": "5da1ac87", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "code", "execution_count": null, "id": "7d914261", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "a9c92afe", "metadata": {}, "source": [ "# Representación de matrices\n", "Un método para visualizar matrices cuando la dimensión es muy grande es mediante su representación como imagen. \n", "\n", "Para ello, la matriz $A\\in\\mathcal{M}_n$, consideramos una cuadrícula de $nxn$ píxeles. Numeramos los píxeles comenzando desde la esquina superior izquierda (ojo - esto depende fuertemente del lenguaje y librerías que estemos utilizando). Así, la posición $(1,1)$ será la esquina superior izquierda y la $(1,n)$ la superior derecha.\n", "\n", "En cada posición dibujaremos un color según el valor de dicha posición en la matriz. \n", "\n", "Sage utiliza el color negro para el mayor valor, el blanco para el más pequeño y tonos de grises para los intermedios. (Esto también depende del lenguaje/librería que usemos, usualmente se emplea el blanco para el valor más grande)." ] }, { "cell_type": "code", "execution_count": 48, "id": "60aab1c2", "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ], "text/latex": [ "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left(\\begin{array}{rrrr}\n", "1 & 2 & 3 & 4 \\\\\n", "1 & 1 & 1 & 1 \\\\\n", "4 & 4 & 4 & 4 \\\\\n", "1 & 1 & 1 & 4\n", "\\end{array}\\right)$$" ], "text/plain": [ "[1 2 3 4]\n", "[1 1 1 1]\n", "[4 4 4 4]\n", "[1 1 1 4]" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "A = matrix([[1,2,3,4],[1,1,1,1],[4,4,4,4],[1,1,1,4]])\n", "show(A)" ] }, { "cell_type": "code", "execution_count": 50, "id": "8782f013", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcwAAAHWCAYAAAD6lrl7AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAT/ElEQVR4nO3dX4hc9fnA4XdMyCq4MxBsAiGjhBZKhhDBZCgrKsW2C3sh5q5XaUrTi+AakNypFy29WUEoLTgGpcU7iRQbFVpDF2o2iggbMVQMCIKQFE3VgjPrQjcYz+/Cbuj+stF3kpkzZ7PPA3Oxw+z5vmcPyYdz5l+tKIoiAIBvdNOoBwCAtUAwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwh+Tpp5+OHTt2xM033xx79uyJ119/fdQjrWunTp2KBx54ILZt2xa1Wi1eeumlUY+07s3MzES73Y7x8fHYsmVL7Nu3L95///1Rj7WuHT16NHbv3h31ej3q9XpMTEzEq6++OuqxKkMwh+CFF16IRx55JB5//PF455134t57742pqak4d+7cqEdbtxYXF+POO++Mp556atSj8F9zc3MxPT0db731VszOzsaXX34Zk5OTsbi4OOrR1q3t27fHE088EadPn47Tp0/H/fffHw8++GC89957ox6tEmo+fH3wfvCDH8Rdd90VR48evXzfzp07Y9++fTEzMzPCyYiIqNVqcfz48di3b9+oR+F/fPrpp7Fly5aYm5uL++67b9Tj8F+bN2+OJ598Mg4ePDjqUUbOGeaAXbx4Md5+++2YnJxccf/k5GS8+eabI5oKqq/b7UbE1/9BM3qXLl2KY8eOxeLiYkxMTIx6nErYOOoBbjSfffZZXLp0KbZu3bri/q1bt8aFCxdGNBVUW1EUceTIkbjnnnti165dox5nXXv33XdjYmIi/vOf/8Stt94ax48fj1arNeqxKkEwh6RWq634uSiKK+4Dvvbwww/HP/7xj3jjjTdGPcq69/3vfz/OnDkTn3/+ebz44otx4MCBmJubE80QzIG77bbbYsOGDVecTX7yySdXnHUCEYcPH45XXnklTp06Fdu3bx/1OOvepk2b4nvf+15EROzduzfm5+fj97//fTzzzDMjnmz0PIc5YJs2bYo9e/bE7OzsivtnZ2fj7rvvHtFUUD1FUcTDDz8cf/7zn+Pvf/977NixY9QjsYqiKGJpaWnUY1SCM8whOHLkSOzfvz/27t0bExMT8eyzz8a5c+fi0KFDox5t3friiy/igw8+uPzzhx9+GGfOnInNmzfH7bffPsLJ1q/p6el4/vnn4+WXX47x8fHLV2UajUbccsstI55ufXrsscdiamoqms1mLCwsxLFjx+LkyZNx4sSJUY9WDQVD0el0ijvuuKPYtGlTcddddxVzc3OjHmlde+2114qIuOJ24MCBUY+2bq12PCKieO6550Y92rr1i1/84vL/W9/5zneKH/3oR8Xf/va3UY9VGd6HCQAJnsMEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBHJKlpaX49a9/7SOlKsZxqR7HpHock9X54IIh6fV60Wg0otvtRr1eH/U4/JfjUj2OSfU4Jqtb02eYnU6n0tsbhqrvc9W3NwxV3+dh/A2rflzW49+w6sck4gbY59F+Mt/12blzZ2W31+12i4gout3uwLZZFNXe57WwvWEcl6rv86C3N+htrsdjMoxtVv2YFEW19zmj9G8r+eqrr+Kjjz6K8fHx6/5C5UuXLkWv1xvQZIPd3vJ2BjlfRLX3eS1sbxjHper7POjtDXqb6/GYDGObVT8mEdXc56IoYmFhIbZt2xY33fTNF11Lfw7zn//8ZzSbzTKXBIBvdP78+W/9AvPSzzDHx8cj4uvhPJkMwCj1er1oNpuX2/RNSg/m8mXYer0umABUQuYpwjX9KlkAKItgAkCCYAJAgmACQELpL/rJuN73ZzJ4hw4dGvUIrOJnP/vZqEfg/5mYmBj1CAyJM0wASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASCgtmJ1OJ1qtVrTb7bKWBICBKS2Y09PTcfbs2Zifny9rSQAYGJdkASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIKG0YHY6nWi1WtFut8taEgAGprRgTk9Px9mzZ2N+fr6sJQFgYFySBYAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSAhGsK5tNPPx07duyIm2++Ofbs2ROvv/76oOcCgErpO5gvvPBCPPLII/H444/HO++8E/fee29MTU3FuXPnhjEfAFRC38H87W9/GwcPHoxf/vKXsXPnzvjd734XzWYzjh49Ooz5AKAS+grmxYsX4+23347JyckV909OTsabb7656u8sLS1Fr9dbcQOAtaavYH722Wdx6dKl2Lp164r7t27dGhcuXFj1d2ZmZqLRaFy+NZvNa58WAEbkml70U6vVVvxcFMUV9y179NFHo9vtXr6dP3/+WpYEgJHa2M+Db7vtttiwYcMVZ5OffPLJFWedy8bGxmJsbOzaJwSACujrDHPTpk2xZ8+emJ2dXXH/7Oxs3H333QMdDACqpK8zzIiII0eOxP79+2Pv3r0xMTERzz77bJw7dy4OHTo0jPkAoBL6DuZPf/rT+Pe//x2/+c1v4uOPP45du3bFX//617jjjjuGMR8AVELfwYyIeOihh+Khhx4a9CwAUFk+SxYAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABJKC2an04lWqxXtdrusJQFgYEoL5vT0dJw9ezbm5+fLWhIABsYlWQBIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgoLZidTidarVa02+2ylgSAgSktmNPT03H27NmYn58va0kAGBiXZAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgoVYURVHmgr1eLxqNRnS73ajX66sPVauVORIAN7irpS7TpGXOMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgobRgdjqdaLVa0W63y1oSAAbGZ8kCcMPzWbIAUBLBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgITSgtnpdKLVakW73S5rSQAYmFpRFEWZC/Z6vWg0GtHtdqNer68+VK1W5kgA3OCulrpMk5a5JAsACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACX0H89SpU/HAAw/Etm3bolarxUsvvTSEsQCgWvoO5uLiYtx5553x1FNPDWMeAKikjf3+wtTUVExNTQ1jFgCoLM9hAkBC32eY/VpaWoqlpaXLP/d6vWEvCQADN/QzzJmZmWg0GpdvzWZz2EsCwMANPZiPPvpodLvdy7fz588Pe0kAGLihX5IdGxuLsbGxYS8DAEPVdzC/+OKL+OCDDy7//OGHH8aZM2di8+bNcfvttw90OACojKJPr732WhERV9wOHDiQ+v1ut1tERNHtdq/6mNW27+bm5ubmdq2362nSsr7PMH/4wx/G100DgPXD+zABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASChtGB2Op1otVrRbrfLWhIABqZWlPzllr1eLxqNRnS73ajX66sPVauVORIAN7irpS7TpGUuyQJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkBCacHsdDrRarWi3W6XtSQADEytKIqizAV7vV40Go3odrtRr9dXH6pWK3MkAG5wV0tdpknLXJIFgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBICEjaMeYDUlf5YCAHwrZ5gAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkFBaMDudTrRarWi322UtCQADUytK/uDWXq8XjUYjut1u1Ov1MpcGgBX6aZJLsgCQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJBQWjA7nU60Wq1ot9tlLQkAA1MriqIoc8FerxeNRiO63W7U6/UylwaAFfppkkuyAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQ0FcwZ2Zmot1ux/j4eGzZsiX27dsX77///rBmA4DK6CuYc3NzMT09HW+99VbMzs7Gl19+GZOTk7G4uDis+QCgEq7rC6Q//fTT2LJlS8zNzcV9992X+h1fIA1AVfTTpI3Xs1C3242IiM2bN1/1MUtLS7G0tLRiOABYa675RT9FUcSRI0finnvuiV27dl31cTMzM9FoNC7fms3mtS4JACNzzZdkp6en4y9/+Uu88cYbsX379qs+brUzzGaz6ZIsACM39Euyhw8fjldeeSVOnTr1jbGMiBgbG4uxsbFrWQYAKqOvYBZFEYcPH47jx4/HyZMnY8eOHcOaCwAqpa9gTk9Px/PPPx8vv/xyjI+Px4ULFyIiotFoxC233DKUAQGgCvp6DrNWq616/3PPPRc///nPU9vwthIAqmJoz2Fex1s2AWBN81myAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJBQWjA7nU60Wq1ot9tlLQkAA9PX92EOgu/DBKAq+mmSS7IAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkFBaMDudTrRarWi322UtCQADUyuKoihzwV6vF41GI7rdbtTr9TKXBoAV+mmSS7IAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkLBx1AOsplarjXoEgGtS8lvbKZEzTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIKC2YnU4nWq1WtNvtspYEgIGpFSV/8GGv14tGoxHdbjfq9frqQ/ksWWCN8lmya0umSctckgWABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBICE0oLZ6XSi1WpFu90ua0kAGJhaURRFmQv2er1oNBrR7XajXq+vPlStVuZIAANT8n+pXKdMk5a5JAsACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAl9BfPo0aOxe/fuqNfrUa/XY2JiIl599dVhzQYAldFXMLdv3x5PPPFEnD59Ok6fPh33339/PPjgg/Hee+8Naz4AqITr/gLpzZs3x5NPPhkHDx5MPd4XSAM3Ml8gvbb08wXSG691kUuXLsWf/vSnWFxcjImJias+bmlpKZaWllYMBwBrTd8v+nn33Xfj1ltvjbGxsTh06FAcP348Wq3WVR8/MzMTjUbj8q3ZbF7XwAAwCn1fkr148WKcO3cuPv/883jxxRfjD3/4Q8zNzV01mqudYTabTZdkgRuSS7JrSz+XZK/7Ocwf//jH8d3vfjeeeeaZgQ0nmMBaJZhrSz/BvO73YRZFseIMEgBuRH296Oexxx6LqampaDabsbCwEMeOHYuTJ0/GiRMnhjUfAFRCX8H817/+Ffv374+PP/44Go1G7N69O06cOBE/+clPhjUfAFRCX8H84x//OKw5AKDSfJYsACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAklBbMTqcTrVYr2u12WUsCwMDUiqIoylyw1+tFo9GIbrcb9Xp99aFqtTJHAhiYkv9L5TplmrTMJVkASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIKC2YnU4nWq1WtNvtspYEgIGpFUVRlLlgr9eLRqMR3W436vX66kPVamWOBDAwJf+XynXKNGmZS7IAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJCwcdQDrMYbfwGoGmeYAJAgmACQIJgAkCCYAJAgmACQIJgAkFD620qW3zLS6/XKXhoAVlhuUebtjKUHc2FhISIims1m2UsDwKoWFhai0Wh842NqRcmfEvDVV1/FRx99FOPj41Gr1a5rW+12O+bn5wc02WC31+v1otlsxvnz57/1W7z7UeV9XgvbG8Zxqfo+D3p7g97mejwmw9hm1Y9JRDX3uSiKWFhYiG3btsVNN33zs5Sln2HedNNNsX379oFsa8OGDQM9mIPeXkREvV6v9IzrbXvLBnlcqr7Pw/gbVv3fynr8G1b9mERUd5+/7cxy2Zp+0c/09HSltzcMVd/nqm9vGKq+z8P4G1b9uKzHv2HVj0nE2t/n0i/Jrhe9Xi8ajUZ0u92hnCFxbRyX6nFMqscxWd2aPsOssrGxsfjVr34VY2Njox6F/+G4VI9jUj2OyeqcYQJAgjNMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABI+D9RYLlwQGPeBgAAAABJRU5ErkJggg==\n", "text/plain": [ "Graphics object consisting of 1 graphics primitive" ] }, "execution_count": 50, "metadata": {}, "output_type": "execute_result" } ], "source": [ "matrix_plot(-A)" ] }, { "cell_type": "markdown", "id": "2f8a521e", "metadata": {}, "source": [ "## Matrices dispersas\n", "Muchos sistemas lineales tienen gran parte de los coeficientes nulos. En ese tipo de sistemas, no merece la pena guardar todos esos ceros, sino la lista de valores no nulos. Para indicarle al ordenador que trabaje de ese modo, debemos indicar que la matriz es dispersa (sparse)." ] }, { "cell_type": "code", "execution_count": 51, "id": "08505bd0", "metadata": {}, "outputs": [ { "data": { "text/html": [ "" ], "text/latex": [ "$$\\newcommand{\\Bold}[1]{\\mathbf{#1}}\\left(\\begin{array}{rrrr}\n", "1 & 0 & 0 & 1 \\\\\n", "0 & 2 & 0 & 0 \\\\\n", "1 & 0 & 3 & 0 \\\\\n", "0 & 0 & 0 & 4\n", "\\end{array}\\right)$$" ], "text/plain": [ "[1 0 0 1]\n", "[0 2 0 0]\n", "[1 0 3 0]\n", "[0 0 0 4]" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "A = matrix(QQ,[[1,0,0,1],[0,2,0,0],[1,0,3,0],[0,0,0,4]],sparse=true)\n", "show(A)" ] }, { "cell_type": "markdown", "id": "3670f067", "metadata": {}, "source": [ "En ese caso, al representarla únicamente marca los valores no nulos, sin discriminar el valor. " ] }, { "cell_type": "code", "execution_count": 52, "id": "60f9c0da", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcwAAAHWCAYAAAD6lrl7AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAWq0lEQVR4nO3dX4hc9fn48eckkhhiZkhIDYRda6jgz6msYHaQtVba2obmQgwEtTchbb2xNUJJL76NXrR3ayqlFDqGSIvtjSaCf28MDWiiRQRXDJGKhX5JyQa1Gkhn4pKumJzvRWvo/syfZ3TmzGz29YJhcTN7Ps/scec9c2ZnT1GWZRkAwAUtGvQAADAfCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCGafPPLII7Fu3bq4/PLLY/369fHKK68MeqQF7eWXX47bb7891q5dG0VRxLPPPjvokRa8ycnJaDabsWLFirjyyitj06ZN8de//nXQYy1ou3btirGxsajValGr1WJiYiJeeOGFQY81NASzD/bu3Rs/+clP4sEHH4w333wzvv71r8fGjRvj6NGjgx5twZqZmYkbbrghfvvb3w56FP7j4MGDcd9998Vrr70W+/fvj08++SQ2bNgQMzMzgx5twRoZGYmHHnoopqamYmpqKr71rW/FHXfcEX/5y18GPdpQKPzx9d676aab4sYbb4xdu3ad/dx1110XmzZtisnJyQFORkREURTxzDPPxKZNmwY9Cv/lww8/jCuvvDIOHjwYt95666DH4T9WrVoVDz/8cNxzzz2DHmXgPMPssY8//jjeeOON2LBhw5zPb9iwIV599dUBTQXDr91uR8S/76AZvNOnT8eePXtiZmYmJiYmBj3OULhs0ANcao4fPx6nT5+ONWvWzPn8mjVr4v333x/QVDDcyrKM7du3xy233BLXX3/9oMdZ0N56662YmJiIf/3rX3HFFVfEM888E41GY9BjDQXB7JOiKOb8d1mWn/kc8G/btm2Lw4cPx5///OdBj7LgXXvttXHo0KH45z//GU899VRs3bo1Dh48KJohmD23evXqWLx48WeeTX7wwQefedYJRNx///3x/PPPx8svvxwjIyODHmfBW7JkSVxzzTURETE+Ph6vv/56/OY3v4ndu3cPeLLB8xpmjy1ZsiTWr18f+/fvn/P5/fv3x8033zygqWD4lGUZ27Zti6effjpefPHFWLdu3aBH4hzKsozZ2dlBjzEUPMPsg+3bt8eWLVtifHw8JiYm4tFHH42jR4/GvffeO+jRFqyPPvoo/va3v5397yNHjsShQ4di1apVcdVVVw1wsoXrvvvui8cffzyee+65WLFixdmjMvV6PZYtWzbg6RamBx54IDZu3Bijo6Nx8uTJ2LNnTxw4cCD27ds36NGGQ0lftFqt8stf/nK5ZMmS8sYbbywPHjw46JEWtJdeeqmMiM9ctm7dOujRFqxz7Y+IKB977LFBj7Zg/fCHPzx7v/WlL32pvO2228o//elPgx5raHgfJgAkeA0TABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEs09mZ2fjF7/4hT8pNWTsl+Fjnwwf++Tc/OGCPul0OlGv16PdbketVhv0OPyH/TJ87JPhY5+c27x+htlqtYZ6e/0w7Ld52LfXD8N+m/vxPRz2/bIQv4fDvk8iLoHbPNi/zPfFXHfddUO7vXa7XUZE2W63e7bNshzu2zwftteP/TLst7nX2+v1NhfiPunHNod9n5TlcN/mjMrPVnLmzJl49913Y8WKFV/4hMqnT5+OTqfTo8l6u71Pt9PL+SKG+zbPh+31Y78M+23u9fZ6vc2FuE/6sc1h3ycRw3mby7KMkydPxtq1a2PRogsfdK38Ncxjx47F6OholUsCwAVNT09f9ATmlT/DXLFiRUT8ezgvJgMwSJ1OJ0ZHR8+26UIqD+anh2FrtZpgAjAUMi8RzuvfkgWAqggmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkVBbMVqsVjUYjms1mVUsCQM8UZVmWVS7Y6XSiXq9Hu912AmkABqqbJjkkCwAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmVBbPVakWj0Yhms1nVkgDQM0VZlmWVC3Y6najX69Fut6NWq1W5NADM0U2THJIFgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBICEzxXMRx55JNatWxeXX355rF+/Pl555ZVezzWvHTk+Ezv3vRP3P/Fm7Nz3Thw5PjPokQBS3H+dX9fnw9y7d29s2bIlHnnkkfja174Wu3fvjt/97nfx9ttvx1VXXXXRr7/Uz4f55NR0/Oypw1EURZRlefbjzs1jcef46KDHAzivhXj/1U2Tug7mTTfdFDfeeGPs2rXr7Oeuu+662LRpU0xOTvZ0uPnmyPGZuO1XB+LMOb6ji4qIF3/6jbh69fLqBwO4iIV6/9W3E0h//PHH8cYbb8SGDRvmfH7Dhg3x6quvnvNrZmdno9PpzLlcqp6cmo6iKM75b0VRxN6p6YonAshx/3VxXQXz+PHjcfr06VizZs2cz69Zsybef//9c37N5ORk1Ov1s5fR0UvzaX1ExLETp+J8T9jLsoxjJ05VPBFAjvuvi/tcv/Tz/z8K+fRY97ns2LEj2u322cv09KX7KGVk5bILPkIbWbms4okActx/XVxXwVy9enUsXrz4M88mP/jgg8886/zU0qVLo1arzblcqu4aH73gI7S7L9EXzYH5z/3XxXUVzCVLlsT69etj//79cz6/f//+uPnmm3s62Hy0bvXy2Ll5LBYVEYsXFXM+7tw8dkm+YA5cGtx/Xdxl3X7B9u3bY8uWLTE+Ph4TExPx6KOPxtGjR+Pee+/tx3zzzp3jo9G8elXsnZqOYydOxcjKZXH3+Kj/2YCh5/7rwrp+W0nEv/9wwS9/+ct477334vrrr49f//rXceutt6a+9lJ+WwkA80tf34f5RQkmAMOib+/DBICFSjABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIKGyYLZarWg0GtFsNqtaEgB6xvkwAViwnA8TAHpMMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIKGyYLZarWg0GtFsNqtaEgB6pijLsqxywU6nE/V6PdrtdtRqtSqXBoA5ummSQ7IAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQUFkwW61WNBqNaDabVS0JAD1TlGVZVrlgp9OJer0e7XY7arValUsDwBzdNMkhWQBIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIqCyYrVYrGo1GNJvNqpYEgJ4pyrIsq1yw0+lEvV6PdrsdtVqtyqUBYI5umuSQLAAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkdB3Ml19+OW6//fZYu3ZtFEURzz77bB/GAoDhclm3XzAzMxM33HBD/OAHP4jNmzf3YyboiyPHZ+LJqek4duJUjKxcFneNj8a61csHPRYwT3QdzI0bN8bGjRv7MQv0zZNT0/Gzpw5HURRRlmUURRG7D/5v7Nw8FneOjw56PGAe8Boml7wjx2fiZ08djjNlxOkz5ZyP//PU4fj78ZlBjwjMA30P5uzsbHQ6nTkXqNKTU9NRFMU5/60oitg7NV3xRMB81PdgTk5ORr1eP3sZHXX4i2odO3EqyrI857+VZRnHTpyqeCJgPup7MHfs2BHtdvvsZXrao3mqNbJy2QWfYY6sXFbxRMB81PdgLl26NGq12pwLVOmu8dELPsO82y/9AAldB/Ojjz6KQ4cOxaFDhyIi4siRI3Ho0KE4evRor2eDnli3enns3DwWi4qIxYuKOR93bh6Lq721BEgoyvM99D6PAwcOxDe/+c3PfH7r1q3xhz/84aJf3+l0ol6vR7vd9myTSv39+Ezs/a/3Yd49PiqWsMB106Sug/lFCSYAw6KbJnkfJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJFQWzFarFY1GI5rNZlVLAkDPOB8mAAuW82ECQI8JJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJFQWzFarFY1GI5rNZlVLAkDPFGVZllUu2Ol0ol6vR7vdjlqtVuXSADBHN01ySBYAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEioLZqvVikajEc1ms6olAaBnirIsyyoX7HQ6Ua/Xo91uR61Wq3JpAJijmyY5JAsACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJlQWz1WpFo9GIZrNZ1ZIA0DNFWZZllQt2Op2o1+vRbrejVqtVuTQAzNFNkxySBYAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSAhMu6ufLk5GQ8/fTT8c4778SyZcvi5ptvjp07d8a1117br/nmpSPHZ+LJqek4duJUjKxcFneNj8a61csHPRYMHT8rzCddnQ/zu9/9bnzve9+LZrMZn3zySTz44IPx1ltvxdtvvx3Ll+f+J7/Uz4f55NR0/Oypw1EURZRlefbjzs1jcef46KDHg6HhZ4Vh0E2TvtAJpD/88MO48sor4+DBg3Hrrbf2fLj55sjxmbjtVwfizDm+o4uKiBd/+o242qNn8LPC0KjsBNLtdjsiIlatWnXe68zOzkan05lzuVQ9OTUdRVGc89+Kooi9U9MVTwTDyc8K89HnDmZZlrF9+/a45ZZb4vrrrz/v9SYnJ6Ner5+9jI5euodajp04Fed7wl6WZRw7cariiWA4+VlhPvrcwdy2bVscPnw4nnjiiQteb8eOHdFut89epqcv3UeOIyuXXfBR88jKZRVPBMPJzwrz0ecK5v333x/PP/98vPTSSzEyMnLB6y5dujRqtdqcy6XqrvHRCz5qvtsvMkBE+FlhfuoqmGVZxrZt2+Lpp5+OF198MdatW9evuealdauXx87NY7GoiFi8qJjzcefmMb/EAP/hZ4X5qKvfkv3xj38cjz/+eDz33HNz3ntZr9dj2bLcIZRL+bdkP/X34zOx97/eW3b3+Kg7ADgHPysMWt/eVnK+1xwee+yx+P73v9/z4QCgn7ppUld/6ecLvGUTAOY1f0sWABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIqC2ar1YpGoxHNZrOqJQGgZ7o6H2YvOB8mAMOimyY5JAsACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACZUFs9VqRaPRiGazWdWSANAzRVmWZZULdjqdqNfr0W63o1arVbk0AMzRTZMckgWABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBICEyoLZarWi0WhEs9msakkA6JmiLMuyygU7nU7U6/Vot9tRq9WqXBoA5uimSQ7JAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQEJlwWy1WtFoNKLZbFa1JAD0TFGWZVnlgp1OJ+r1erTb7ajValUuDQBzdNMkh2QBIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASChq2Du2rUrxsbGolarRa1Wi4mJiXjhhRf6NRsAFTtyfCZ27nsn7n/izdi57504cnxm0CMNjcu6ufLIyEg89NBDcc0110RExB//+Me444474s0334yvfvWrfRkQgGo8OTUdP3vqcBRFEWVZRlEUsfvg/8bOzWNx5/jooMcbuC98AulVq1bFww8/HPfcc0/q+k4gDTB8jhyfidt+dSDOnKMIi4qIF3/6jbh69fLqB+uzSk4gffr06dizZ0/MzMzExMTEea83OzsbnU5nzgWA4fLk1HQURXHOfyuKIvZOTVc80fDpOphvvfVWXHHFFbF06dK4995745lnnolGo3He609OTka9Xj97GR31tB5g2Bw7cSrOd8CxLMs4duJUxRMNn66Dee2118ahQ4fitddeix/96EexdevWePvtt897/R07dkS73T57mZ72KAVg2IysXHbBZ5gjK5dVPNHw+cKvYX7729+Or3zlK7F79+7U9b2GCTB8vIbZx9cwP1WWZczOzn7RzQAwQOtWL4+dm8diURGxeFEx5+POzWOXZCy71dXbSh544IHYuHFjjI6OxsmTJ2PPnj1x4MCB2LdvX7/mA6Aid46PRvPqVbF3ajqOnTgVIyuXxd3jo2L5H10F8x//+Eds2bIl3nvvvajX6zE2Nhb79u2L73znO/2aD4AKXb16efzPd//foMcYSl0F8/e//32/5gCAoeZvyQJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQEJlwWy1WtFoNKLZbFa1JAD0TFGWZVnlgp1OJ+r1erTb7ajValUuDQBzdNMkh2QBIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgobJgtlqtaDQa0Ww2q1oSAHqmKMuyrHLBTqcT9Xo92u121Gq1KpcGgDm6aZJDsgCQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJBwWdULfnr6zU6nU/XSADDHpy3KnBq68mCePHkyIiJGR0erXhoAzunkyZNRr9cveJ2izGS1h86cORPvvvturFixIoqi+ELbajab8frrr/dost5ur9PpxOjoaExPT1/0LN7dGObbPB+214/9Muy3udfb6/U2F+I+6cc2h32fRAznbS7LMk6ePBlr166NRYsu/Cpl5c8wFy1aFCMjIz3Z1uLFi3u6M3u9vYiIWq021DMutO19qpf7Zdhvcz++h8P+s7IQv4fDvk8ihvc2X+yZ5afm9S/93HfffUO9vX4Y9ts87Nvrh2G/zf34Hg77flmI38Nh3ycR8/82V35IdqHodDpRr9ej3W735RkSn4/9Mnzsk+Fjn5zbvH6GOcyWLl0aP//5z2Pp0qWDHoX/Yr8MH/tk+Ngn5+YZJgAkeIYJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJ/wdIj1eQHvlwNQAAAABJRU5ErkJggg==\n", "text/plain": [ "Graphics object consisting of 1 graphics primitive" ] }, "execution_count": 52, "metadata": {}, "output_type": "execute_result" } ], "source": [ "matrix_plot(A)" ] }, { "cell_type": "markdown", "id": "d09ac686", "metadata": {}, "source": [ "
\n", " Ejercicio 2. \n", "Consideremos la matriz\n", "$$A=\\left(\\begin{array}{rrrrr} \n", " 4 & -1 & 0 & 0 & -1 \\\\\n", " 0 & 4 & -1 & -1 & 0 \\\\\n", "-1 & 0 & 4 & 0 & 0 \\\\\n", "-1 & 0 & 0 & 4 & 0 \\\\\n", "-1 & 0 & 0 & 0 & 4 \\\\\n", "\\end{array}\\right)$$\n", "\n", "a) Calcula la factorización LU de $A$. Cuenta cuántos valores no nulos hay en la matriz $A$ y cuántos hay en las matrices $L$ y $U$. \n", " \n", "b) Aplica tres pasos del método de Jacobi para resolver el sistema $Ax = b$ con $b=(1,1,1,1,1)$,\n", " usando matrices dispersas. \n", " \n", "c) Idem para el método de Gauss-Seidel. \n", "
" ] }, { "cell_type": "code", "execution_count": null, "id": "8e42302f", "metadata": {}, "outputs": [], "source": [ "A = matrix([[4,-1,0,0,-1],\n", " [0,4,-1,-1,0],\n", " [-1,0,4,0,0],\n", " [-1,0,0,4,0],\n", " [-1,0,0,0,4]])\n", "show(A)" ] }, { "cell_type": "code", "execution_count": null, "id": "a3eac119", "metadata": {}, "outputs": [], "source": [] }, { "cell_type": "markdown", "id": "fd1a3c00", "metadata": {}, "source": [ "## Difuminando y perfilando una imagen (opcional) \n", "\n", "Vamos a ver un ejemplo de cómo se puede aplicar el método iterativo de Gauss-Seidel. Concretamente, veremos cómo funciona (de modo básico) un modelo de difuminado y perfilado. Estos modelos se aplican por ejemplo a imágenes de satélite para eliminar los efectos de difusión de la atmósfera, al análisis de imágenes por ordenador, para retocar fotografías, etc.\n", "\n", "Cada imagen viene cargada en una matriz (por simplicidad asumiremos que está en blanco y negro). Vamos a necesitar algunas funciones que nos pasen los datos de la matriz a un vector y viceversa. Para ello, basta ejecutar la siguiente celda." ] }, { "cell_type": "code", "execution_count": 29, "id": "e81a2b19", "metadata": {}, "outputs": [], "source": [ "def matrix_to_vector(A):\n", " n,m=A.dimensions();\n", " v=vector(RDF,n*m)\n", " for i in range(0,n):\n", " for j in range(0,m):\n", " v[m*i+j]=A[i,j];\n", " return v;\n", "def vector_to_matrix(v):\n", " n=sqrt(v.degree());\n", " A=matrix(RDF,n);\n", " for i in range(0,n):\n", " for j in range(0,n):\n", " A[i,j]=v[n*i+j];\n", " return A;" ] }, { "cell_type": "markdown", "id": "298b5008", "metadata": {}, "source": [ "Para pasar la matriz a un vector, numeramos las celdas comenzando por la primera fila:\n", "\n", "$\\left(\\begin{array}{rrr}\n", "1 & 2 & 3 \\\\\n", "4 & 5 & 6 \\\\\n", "7 & 8 & 9\n", "\\end{array}\\right)\\to(1,2,3,4,5,6,7,8,9)$\n", "\n", "Definimos la matriz en la que tenemos la imagen. Cuando creamos la matriz estará a cero (color negro). Vamos a poner algunas posiciones a uno (blanco)." ] }, { "cell_type": "code", "execution_count": 30, "id": "37d9ab7f", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWMAAAFtCAYAAADSyAuRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAR2UlEQVR4nO3dX4hU9/nA4XdUHEOyO8SkCuIkSFuSLGIgui0bkv6JdmEJIfaqF0HsnxvLKoo3bZKLlt5soFcFNxLbkl6U1FBaTS8ayUKrpgRhNVkiKQQCAbck1qbQ3XUhE7KZ38WvLmxN053VM+edmeeBuZjxzJ73LPLheObM10qz2WwGAKVaVfYAAIgxQApiDJCAGAMkIMYACYgxQAJiDJCAGAMkIMYACYgxQAI9H+Nnn302tmzZEuvWrYvt27fHq6++WvZIhTl79mw89thjsWnTpqhUKnHy5MmyRyrU2NhYDA4ORl9fX2zYsCF2794db7/9dtljFero0aOxbdu26O/vj/7+/hgaGoqXX3657LHaamxsLCqVShw6dKjsUVrS0zF+8cUX49ChQ/H000/HG2+8EQ8//HCMjIzEpUuXyh6tEPPz83H//ffHkSNHyh6lLc6cOROjo6Nx7ty5mJiYiI8//jiGh4djfn6+7NEKs3nz5njmmWfi/Pnzcf78+XjkkUfi8ccfj7feeqvs0dpicnIyjh07Ftu2bSt7lNY1e9iXvvSl5r59+5a8du+99zZ/+MMfljRR+0RE88SJE2WP0VZXrlxpRkTzzJkzZY/SVrfffnvzF7/4RdljFG5ubq75xS9+sTkxMdH86le/2jx48GDZI7WkZ8+MP/roo7hw4UIMDw8veX14eDhee+21kqaiSDMzMxERsX79+pInaY+FhYU4fvx4zM/Px9DQUNnjFG50dDQeffTR2LVrV9mjrMiasgcoywcffBALCwuxcePGJa9v3LgxLl++XNJUFKXZbMbhw4fjoYceiq1bt5Y9TqEuXrwYQ0ND8eGHH8Ztt90WJ06ciIGBgbLHKtTx48fj9ddfj8nJybJHWbGejfE1lUplyfNms3nda3S+/fv3x5tvvhl/+ctfyh6lcPfcc09MTU3Fv/71r/jd734Xe/fujTNnznRtkKenp+PgwYPxyiuvxLp168oeZ8V6NsZ33nlnrF69+rqz4CtXrlx3tkxnO3DgQPzhD3+Is2fPxubNm8sep3Br166NL3zhCxERsWPHjpicnIyf/exn8dxzz5U8WTEuXLgQV65cie3bty++trCwEGfPno0jR45Eo9GI1atXlzjh8vTsNeO1a9fG9u3bY2JiYsnrExMT8eCDD5Y0FTdTs9mM/fv3x+9///v405/+FFu2bCl7pFI0m81oNBplj1GYnTt3xsWLF2NqamrxsWPHjnjiiSdiamqqI0Ic0cNnxhERhw8fjj179sSOHTtiaGgojh07FpcuXYp9+/aVPVohrl69Gu+8887i83fffTempqZi/fr1cdddd5U4WTFGR0fjhRdeiJdeein6+voW/xVUq9XilltuKXm6Yjz11FMxMjIS9Xo95ubm4vjx43H69Ok4depU2aMVpq+v77rPAW699da44447OuvzgXJv5ijf+Ph48+67726uXbu2+cADD3T1bU9//vOfmxFx3WPv3r1lj1aITzvWiGg+//zzZY9WmO9+97uLf58/97nPNXfu3Nl85ZVXyh6r7Trx1rZKs+k/JAUoW89eMwbIRIwBEhBjgATEGCABMQZIQIwBEhBjgAR6PsaNRiN+/OMfd/XXRf+TY+4Njrmz9PyXPmZnZ6NWq8XMzEz09/eXPU5bOGbH3K06+Zg74sx4fHy80O1b1Y55eu2Y2/E7KvrnZ9t+JbIdQzcc87KV+23s5bnvvvsK235mZqYZEc2ZmZkU87RjHxmPuejfUTccc6vbO+YcMy1X21dt++STT+K9996Lvr6+ZS/ivrCwELOzs8veRyvbX9uuqJ+/ku2L3kfGYy76d9QNx9zq9o65/JmazWbMzc3Fpk2bYtWqz74Q0fZrxn/729+iXq+3c5cApZqenv6f/7FB28+M+/r6IuL/h+u0C+wArZidnY16vb7Yvc/S9hhfuzTR398vxkBPWM4l2Y64mwKg24kxQAJiDJCAGAMkIMYACYgxQAJiDJCAGAMkIMYACbQtxuPj4zEwMBCDg4Pt2iVAx2j7QkGdvPgzQCta6Z3LFAAJiDFAAmIMkIAYAyQgxgAJiDFAAmIMkIAYAyQgxgAJiDFAAtamAEjA2hQABbE2BUCHEWOABMQYIAExBkhAjAESEGOABFYU42effTa2bNkS69ati+3bt8err756s+cC6Cktx/jFF1+MQ4cOxdNPPx1vvPFGPPzwwzEyMhKXLl0qYj6AntDylz6+/OUvxwMPPBBHjx5dfO2+++6L3bt3x9jY2P98vy99AL2isC99fPTRR3HhwoUYHh5e8vrw8HC89tprn/qeRqMRs7OzSx4ALNVSjD/44INYWFiIjRs3Lnl948aNcfny5U99z9jYWNRqtcVHvV5f+bQAXWpFH+BVKpUlz5vN5nWvXfPkk0/GzMzM4mN6enoluwToamta2fjOO++M1atXX3cWfOXKlevOlq+pVqtRrVZXPiFAD2jpzHjt2rWxffv2mJiYWPL6xMREPPjggzd1MIBe0tKZcUTE4cOHY8+ePbFjx44YGhqKY8eOxaVLl2Lfvn1FzAfQE1qO8be+9a345z//GT/5yU/i/fffj61bt8Yf//jHuPvuu4uYD6AnWFweoCAWlwfoMGIMkIAYAyQgxgAJiDFAAm2L8fj4eAwMDMTg4GC7dgnQMdzaBlAQt7YBdBgxBkhAjAESEGOABMQYIAExBkhAjAESEGOABMQYIAExBkhAjAESsFAQQAIWCgIoiIWCADqMGAMkIMYACawpe4BeV6lUyh4BCtHmj6M6njNjgATEGCABMQZIQIwBEhBjgATEGCABa1MAJGBtipK5z5hu5T5ja1MAdBwxBkhAjAESEGOABMQYIAExBkhAjAESEGOABMQYIAExBkjA2hQACVibomTWpqBbWZvC2hQAHUeMARIQY4AExBggATEGSECMARJoOcZnz56Nxx57LDZt2hSVSiVOnjxZwFgAvaXlGM/Pz8f9998fR44cKWIegJ60ptU3jIyMxMjISBGzAPSslmPcqkajEY1GY/H57Oxs0bsE6DiFf4A3NjYWtVpt8VGv14veJUDHKTzGTz75ZMzMzCw+pqeni94lQMcp/DJFtVqNarVa9G4AOpr7jAESaPnM+OrVq/HOO+8sPn/33Xdjamoq1q9fH3fddddNHQ6gV7Qc4/Pnz8fXv/71xeeHDx+OiIi9e/fGr371q5s2GEAvaTnGX/va1ywaDXCTuWYMkIAYAyQgxgAJiDFAAmIMkEDbYjw+Ph4DAwMxODjYrl0CdIxKs833qc3OzkatVouZmZno7+9v565TqlQqZY8AhXALbGu9c5kCIAExBkhAjAESEGOABMQYIAExBkhAjAESEGOABMQYIAExBkjA2hQACVibomTWpqBbWZvC2hQAHUeMARIQY4AExBgggTVlD0Dv8cEOXM+ZMUACYgyQgBgDJCDGAAmIMUACYgyQgIWCABKwUFDJenGhIPcZ0yssFATQYcQYIAExBkhAjAESEGOABMQYIAExBkhAjAESEGOABMQYIAFrUwAkYG2KklmbArqXtSkAOowYAyQgxgAJiDFAAmIMkEBLMR4bG4vBwcHo6+uLDRs2xO7du+Ptt98uajaAntFSjM+cOROjo6Nx7ty5mJiYiI8//jiGh4djfn6+qPkAesIN3Wf8j3/8IzZs2BBnzpyJr3zlK8t6j/uMl3KfMXSvtt1nPDMzExER69evv5EfA9Dz1qz0jc1mMw4fPhwPPfRQbN269b9u12g0otFoLD6fnZ1d6S4ButaKz4z3798fb775ZvzmN7/5zO3GxsaiVqstPur1+kp3CdC1VnTN+MCBA3Hy5Mk4e/ZsbNmy5TO3/bQz43q97prxv7lmDN2rlWvGLV2maDabceDAgThx4kScPn36f4Y4IqJarUa1Wm1lNwA9p6UYj46OxgsvvBAvvfRS9PX1xeXLlyMiolarxS233FLIgAC9oKXLFP/tn9TPP/98fPvb317Wz3Br21IuU0D3KvQyBQA3n7UpABIQY4AExBggATEGSECMARIQY4AE2hbj8fHxGBgYiMHBwXbtEqBj3NB6xivhSx9L+dIHdK+2rWcMwM0hxgAJiDFAAmIMkIAYAyQgxgAJiDFAAmIMkIAYAyQgxgAJWJsCIAFrU5TM2hTQvaxNAdBhxBggATEGSECMARJYU/YA9J5e/NCyF/mgtjXOjAESEGOABMQYIAExBkhAjAESsDYFQALWpiiZ27zoVm5tszYFQMcRY4AExBggATEGSECMARIQY4AExBggATEGSECMARIQY4AExBggAQsFASRgoaCSWSiIbmWhIAsFAXQcMQZIQIwBEhBjgATEGCCBlmJ89OjR2LZtW/T390d/f38MDQ3Fyy+/XNRsAD2jpRhv3rw5nnnmmTh//nycP38+HnnkkXj88cfjrbfeKmo+gJ5ww/cZr1+/Pn7605/G9773vWVt7z7jpdxnTLdyn3FrvVuz0p0sLCzEb3/725ifn4+hoaGV/hgAYgUxvnjxYgwNDcWHH34Yt912W5w4cSIGBgb+6/aNRiMajcbi89nZ2ZVNCtDFWr6b4p577ompqak4d+5cfP/734+9e/fGX//61/+6/djYWNRqtcVHvV6/oYEButENXzPetWtXfP7zn4/nnnvuU//8086M6/W6a8b/5pox3co14zZdM76m2Wwuie1/qlarUa1Wb3Q3AF2tpRg/9dRTMTIyEvV6Pebm5uL48eNx+vTpOHXqVFHzAfSElmL897//Pfbs2RPvv/9+1Gq12LZtW5w6dSq+8Y1vFDUfQE9oKca//OUvi5oDoKdZmwIgATEGSECMARIQY4AExBggATEGSKBtMR4fH4+BgYEYHBxs1y4BOsYNr03RKusZL2VtCrqVtSla653LFAAJiDFAAmIMkIAYAyQgxgAJiDFAAmIMkIAYAyQgxgAJiDFAAtamAEjA2hQlszYF3craFNamAOg4YgyQgBgDJCDGAAmsKXuAXudDDiDCmTFACmIMkIAYAyQgxgAJiDFAAtamAEjA2hQABbE2BUCHEWOABMQYIAExBkhAjAESEGOABMQYIAExBkhAjAESEGOABKxNAZCAtSkACmJtCoAOI8YACYgxQAJiDJCAGAMkIMYACdxQjMfGxqJSqcShQ4du0jgAvWnFMZ6cnIxjx47Ftm3bbuY8AD1pRTG+evVqPPHEE/Hzn/88br/99ps9E0DPWVGMR0dH49FHH41du3bd7HkAetKaVt9w/PjxeP3112NycnJZ2zcajWg0GovPZ2dnW90lQNdr6cx4eno6Dh48GL/+9a9j3bp1y3rP2NhY1Gq1xUe9Xl/RoADdrKWFgk6ePBnf/OY3Y/Xq1YuvLSwsRKVSiVWrVkWj0VjyZxGffmZcr9ctFAR0vVYWCmrpMsXOnTvj4sWLS177zne+E/fee2/84Ac/uC7EERHVajWq1WoruwHoOS3FuK+vL7Zu3brktVtvvTXuuOOO614HYPl8Aw8ggZbvpvhPp0+fvgljAPQ2Z8YACYgxQAJiDJCAGAMkIMYACYgxQAJti/H4+HgMDAzE4OBgu3YJ0DFaWpviZmjlu9oAnayV3rlMAZCAGAMkIMYACYgxQAJiDJCAGAMkIMYACYgxQAJiDJCAGAMkYG0KgASsTQFQEGtTAHQYMQZIQIwBEhBjgATEGCABMQZIQIwBEhBjgATEGCCBNe3e4bUv/M3OzrZ71wBtda1zy/mic9tjPDc3FxER9Xq93bsGKMXc3FzUarXP3Kbta1N88skn8d5770VfX19UKpVlvWdwcDAmJyeXvY9Wtp+dnY16vR7T09PLXiujyHnasY+Mx1z076gbjrnV7R1z+cfcbDZjbm4uNm3aFKtWffZV4bafGa9atSo2b97c0ntWr17d0qJCrW4fEdHf37/s97Rjnl475nb8jiI6+5hXcrwRjrnsmf7XGfE1HfEB3ujoaKHbt6od8/TaMbfjd1T0z8+2/UpkO4ZuOOblavtlimx6cUlPx+yYu1UnH3NHnBkXqVqtxo9+9KOoVqtlj9I2jrk3OObO0vNnxgAZ9PyZMUAGYgyQgBgDJCDGAAmIMUACYgyQgBgDJCDGAAn8H1L1wA2ggz6eAAAAAElFTkSuQmCC\n", "text/plain": [ "Graphics object consisting of 1 graphics primitive" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "n=5\n", "B = matrix.zero(RDF,n)\n", "B[1,1:4]=Matrix([1,1,1]) \n", "B[3,1:4]=Matrix([1,1,1]) \n", "B[2,1]=1\n", "show(matrix_plot(B), figsize=5)" ] }, { "cell_type": "markdown", "id": "d93a1ce8", "metadata": {}, "source": [ "Creamos la matriz de difuminado (está libremente basado en el difuminado gaussiano). Para ello cada posición pasa a tomar \"un poco del color de los vecinos\". Estamos considerando cuatro vecinos (arriba, abajo, derecha e izquierda). Si A(i,j) es el color del punto (i,j), después de la transformación será: $$A(i,j)\\to (1-4d) A(i,j)+d\\,A(i-1,j)+d\\,A(i+1,j)+d\\,A(i,j-1)+d\\,A(i,j+1)$$\n", "\n", "Esto puede hacerse escribiendo la matriz $A$ en forma de vector y multiplicando el vector por cierta matriz (dada por la expresión anterior). Guardaremos en $N$ la matriz por la que hay que multiplicar (matriz de difuminado)." ] }, { "cell_type": "code", "execution_count": 31, "id": "7c457499", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAdoAAAHWCAYAAADQJkjUAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAkFElEQVR4nO3df2xV9f3H8VdFuDLWdkPsr1n67RbYRBiJ0IHMobjZ0GTEKvkOtn/AZEZnIWFogsQsQHRWXHT/WB1zG2qmkyWoM9FoukCLC5p0pARi/RKmZdRI05XhrVR3CXC+fwh3XtrbH/ecTz+fcz7PR3KD3N7Tns/7nvhqL/3cV1EQBIEAAIARl9k+AQAAkoygBQDAIIIWAACDCFoAAAwiaAEAMIigBQDAIIIWAACDCFoAAAwiaAEAMIigBQDAoMQF7ZNPPqna2lpdccUVWrBggd566y3bpxQbW7duVVFRUc6toqLC9mk5a9++fVqxYoWqqqpUVFSkV155JefjQRBo69atqqqq0tSpU3XTTTfp3XfftXOyDhptfmvXrh1yPS5evNjOyTqoublZdXV1Ki4uVllZmRobG3XkyJGcx3AN5jeW+UV1DSYqaHft2qUNGzbogQceUGdnp773ve+poaFBx48ft31qsXHttdfqxIkT2dvhw4dtn5KzBgcHNX/+fD3xxBPDfvzRRx/V448/rieeeEIdHR2qqKjQLbfcok8++WSCz9RNo81PkpYvX55zPb7++usTeIZua29vV1NTk9555x21trbq7Nmzqq+v1+DgYPYxXIP5jWV+UkTXYJAg3/nOd4K77747575vfetbwf3332/pjOJly5Ytwfz5822fRixJCl5++eXs38+fPx9UVFQEjzzySPa+//znP0FpaWnwm9/8xsIZuu3S+QVBEKxZsya49dZbrZxPHPX19QWSgvb29iAIuAbH69L5BUF012BifqI9c+aMDhw4oPr6+pz76+vrtX//fktnFT9Hjx5VVVWVamtrtXr1an3wwQe2TymWuru71dvbm3M9plIp3XjjjVyP49DW1qaysjLNnj1bd955p/r6+myfkrPS6bQkafr06ZK4Bsfr0vldFMU1mJig7e/v17lz51ReXp5zf3l5uXp7ey2dVbwsWrRIzz33nN588009/fTT6u3t1ZIlS3Ty5EnbpxY7F685rsfCNTQ06Pnnn9eePXv02GOPqaOjQzfffLMymYztU3NOEATauHGjbrjhBs2dO1cS1+B4DDc/Kbpr8PKoT9i2oqKinL8HQTDkPgyvoaEh+9/z5s3T9ddfr2984xt69tlntXHjRotnFl9cj4VbtWpV9r/nzp2rhQsXqqamRq+99ppuv/12i2fmnnXr1unQoUP629/+NuRjXIOjyze/qK7BxPxEO2PGDE2aNGnId2p9fX1DvqPD2EybNk3z5s3T0aNHbZ9K7Fz8bW2ux+hUVlaqpqaG6/ES69ev16uvvqq9e/fq6quvzt7PNTg2+eY3nEKvwcQE7ZQpU7RgwQK1trbm3N/a2qolS5ZYOqt4y2Qyeu+991RZWWn7VGKntrZWFRUVOdfjmTNn1N7ezvVYoJMnT6qnp4fr8YIgCLRu3Tq99NJL2rNnj2pra3M+zjU4stHmN5yCr8HQv07lkBdffDGYPHly8Pvf/z7o6uoKNmzYEEybNi04duyY7VOLhXvvvTdoa2sLPvjgg+Cdd94JfvjDHwbFxcXML49PPvkk6OzsDDo7OwNJweOPPx50dnYG//znP4MgCIJHHnkkKC0tDV566aXg8OHDwY9//OOgsrIyGBgYsHzmbhhpfp988klw7733Bvv37w+6u7uDvXv3Btdff33wta99jfld8LOf/SwoLS0N2traghMnTmRvn376afYxXIP5jTa/KK/BRAVtEARBS0tLUFNTE0yZMiW47rrrcn5VGyNbtWpVUFlZGUyePDmoqqoKbr/99uDdd9+1fVrO2rt3byBpyG3NmjVBEHy+vWLLli1BRUVFkEqlgqVLlwaHDx+2e9IOGWl+n376aVBfXx9cddVVweTJk4OZM2cGa9asCY4fP277tJ0x3OwkBTt37sw+hmswv9HmF+U1WHThCwIAAAMS82+0AAC4iKAFAMAgghYAAIMIWgAADCJoAQAwiKAFAMAgghYAAIMSGbSZTEZbt26l5aNAzC8c5hceMwyH+YUT9fwS+YYVAwMDKi0tVTqdVklJie3TiR3mFw7zC48ZhsP8wol6fon8iTYKLS0tXh8flu3zt318WLbP3/bxUbC9BtvHh2X7/G0fH6kI3zrSGel0OpAUpNPpgj/HNddcE+oc4nw882N+to9nhlyDcZ/fFzlX/H7+/Hl99NFHKi4uLriceGBgIOfPQpw7d87b45kf87N9PDPkGrR5fDqdlvR5HkXBuX+j/fDDD1VdXW37NAAAnnv//ff19a9/PfTnce4n2uLiYklST08P/4gPAJhwAwMDqq6u1pVXXhnJ53MuaC++XFxSUkLQAgCsKfSfLy/Fbx0DAGAQQQsAgEEELQAABhG0AAAYRNACAGAQQQsAgEEELQAABhG0AAAYZCxon3zySdXW1uqKK67QggUL9NZbb5n6UgAAOMvIO0Pt2rVLGzZs0JNPPqnvfve72rFjhxoaGtTV1aWZM2ea+JJZ3f2D+vPfe/Thqc909Ven6kcLq1U7Y9qEHR93vq8fAKJmpFRg0aJFuu666/TUU09l77vmmmvU2Nio5ubmEY8NU7j757/36P7dh1RUVKQgCLJ/bl/5bf3vwtGLCsIeH3e+rx8ApBgUv585c0YHDhxQfX19zv319fXav39/1F8uq7t/UPfvPqTzgXTufJDz56bdh3Ssf9Do8XHn+/oBwJTIg7a/v1/nzp1TeXl5zv3l5eXq7e0d8vhMJqOBgYGcWyH+/PeevG8AXVRUpF1/7zF6fNz5vn4AMMXYL0Nd+j/tiy9FXqq5uVmlpaXZW6FdtB+e+kz5XgUPgkAfnvrM6PFx5/v6AcCUyIN2xowZmjRp0pCfXvv6+ob8lCtJmzdvVjqdzt56egr7yenqr04d8Seyq7861ejxcef7+gHAlMiDdsqUKVqwYIFaW1tz7m9tbdWSJUuGPD6VSmW7Z8N00P5oYfWIP5GtGuWXecIeH3e+rx8ATDHy0vHGjRv1u9/9Tn/4wx/03nvv6ec//7mOHz+uu+++28SXkyTVzpim7Su/rcuKpEmXFeX8uX3lt/U/o2xRCXt83Pm+fgAwxcj2HunzN6x49NFHdeLECc2dO1e//vWvtXTp0lGPC/tr1cf6B7XrC/tAVy2sHldIhD0+7nxfPwBEvb3HWNAWKuoFAgAwHs7vowUAAP9F0AIAYBBBCwCAQQQtAAAGEbQAABhE0AIAYJCRPlqbbPfR+t7n6vv6AeBSidpHa7uP1vc+V9/XDyAZ2Eebh+0+Wt/7XH1fPwDkk5igtd1H63ufq+/rB4B8nAnalpYWzZkzR3V1dQUdb7uP1vc+V9/XDwD5OBO0TU1N6urqUkdHR0HH2+6j9b3P1ff1A0A+zgRtWLb7aH3vc/V9/QCQT2KC1nYfre99rr6vHwDySdT2Hsl+H63vfa6+rx9A/NFHCwCAQeyjBQAgRghaAAAMImgBADCIoAUAwCCCFgAAgwhaAAAMSlwfbVi2+2h973P1ff0Akod9tF9gu4/W9z5X39cPwA3sozXEdh+t732uvq8fQHIRtBfY7qP1vc/V9/UDSC5ngjZsH21Ytvtofe9z9X39AJLLmaAN20cblu0+Wt/7XH1fP4DkciZobbPdR+t7n6vv6weQXATtBbb7aH3vc/V9/QCSi+09l7DdR+t7n6vv6wdgH320AAAYxD5aAABihKAFAMAgghYAAIMIWgAADCJoAQAwiKAFAMAg+mgjZruP1vc+V9/XD8A97KONkO0+Wt/7XH1fP4BosI/WUbb7aH3vc/V9/QDcRdBGxHYfre99rr6vH4C7nAla2320Ydnuo/W9z9X39QNwlzNBa7uPNizbfbS+97n6vn4A7nImaOPOdh+t732uvq8fgLsI2ojY7qP1vc/V9/UDcBfbeyJmu4/W9z5X39cPIDz6aAEAMIh9tAAAxAhBCwCAQQQtAAAGEbQAABhE0AIAYBBBCwCAQfTROoY+WruYH4CosY/WIfTR2sX8AEjso00s+mjtYn4ATCFoHUEfrV3MD4ApzgRt3Ptow6KP1i7mB8AUZ4I27n20YdFHaxfzA2CKM0HrO/po7WJ+AEwhaB1BH61dzA+AKWzvcQx9tHYxPwD00QIAYBD7aAEAiBGCFgAAgwhaAAAMImgBADCIoAUAwCCCFgAAgyLvo926dau2bduWc195ebl6e3uj/lLDok80HPpww/F9/QCGMlL8fu211+qvf/1r9u+TJk0y8WWGGK5PdEf7+/SJjlHY+fk+f9/XD2B4Rl46vvzyy1VRUZG9XXXVVSa+TA76RMOhDzcc39cPID8jQXv06FFVVVWptrZWq1ev1gcffJD3sZlMRgMDAzm3QtAnGg59uOH4vn4A+UUetIsWLdJzzz2nN998U08//bR6e3u1ZMkSnTx5ctjHNzc3q7S0NHurri7sJTb6RMOhDzcc39cPIL/Ig7ahoUErV67UvHnz9IMf/ECvvfaaJOnZZ58d9vGbN29WOp3O3np6CvvOnz7RcOjDDcf39QPIz/j2nmnTpmnevHk6evTosB9PpVIqKSnJuRWCPtFw6MMNx/f1A8jPeNBmMhm99957qqysNPp16BMNhz7ccHxfP4D8Iq/Ju++++7RixQrNnDlTfX19euihh9Te3q7Dhw+rpqZm1OPD1hPRJxoOfbjh+L5+IAmc76NdvXq19u3bp/7+fl111VVavHixHnzwQc2ZM2dMx9NHCwCwKeocivwNK1588cWoPyUAALHFex0DAGAQQQsAgEEELQAABhG0AAAYRNACAGCQkZq8OPO9T9R2Hy3z93v9QBJFvo82LJv7aIfrEw2CwJs+0bDrt3183Pm+fsAVUecQLx1f4HufqO0+Wubv9/qBJCNoL/C9T9R2Hy3z93v9QJI5E7QtLS2aM2eO6urqrHx93/tEbffRMn+/1w8kmTNB29TUpK6uLnV0dFj5+r73idruo2X+fq8fSDJngtY23/tEbffRMn+/1w8kGUF7ge99orb7aJm/3+sHkoztPZfwvU/Udh8t8/d7/YALnO+jDct20AIA/MY+WgAAYoSgBQDAIIIWAACDCFoAAAwiaAEAMIigBQDAIPpoI+Z7n6jtPlrm7/f6ARexjzZCvveJ2u6jZf5+rx+ICvtoHeV7n6jtPlrm7/f6AZcRtBHxvU/Udh8t8/d7/YDLnAla2320YfneJ2q7j5b5+71+wGXOBK3tPtqwfO8Ttd1Hy/z9Xj/gMmeCNu587xO13UfL/P1eP+AygjYivveJ2u6jZf5+rx9wGdt7IuZ7n6jtPlrm7/f6gSjQRwsAgEHsowUAIEYIWgAADCJoAQAwiKAFAMAgghYAAIMIWgAADKKP1jG+94na7qNl/n6vHzCBfbQO8b1P1HYfLfP3e/3AReyjTSjf+0Rt99Eyf7/XD5hE0DrC9z5R2320zN/v9QMmORO0ce+jDcv3PlHbfbTM3+/1AyY5E7Rx76MNy/c+Udt9tMzf7/UDJjkTtL7zvU/Udh8t8/d7/YBJBK0jfO8Ttd1Hy/z9Xj9gEtt7HON7n6jtPlrm7/f6AYk+WgAAjGIfLQAAMULQAgBgEEELAIBBBC0AAAYRtAAAGETQAgBgEH20CeN7nyh9tHYxP2Ao9tEmiO99ovTR2sX8kBTso8WwfO8TpY/WLuYH5EfQJoTvfaL00drF/ID8nAla3/tow/K9T5Q+WruYH5CfM0Hrex9tWL73idJHaxfzA/JzJmgRju99ovTR2sX8gPwI2oTwvU+UPlq7mB+QH9t7Esb3PlH6aO1ifkgC+mgBADCIfbQAAMQIQQsAgEEELQAABhG0AAAYRNACAGAQQQsAgEHj7qPdt2+ffvWrX+nAgQM6ceKEXn75ZTU2NmY/HgSBtm3bpt/+9rc6deqUFi1apJaWFl177bVRnjcMoU80HPpww2MGSJpxB+3g4KDmz5+vO+64QytXrhzy8UcffVSPP/64nnnmGc2ePVsPPfSQbrnlFh05ckTFxcWRnDTMGK5PdEf7+/SJjlHY+TF/ZoBkCvWGFUVFRTk/0QZBoKqqKm3YsEGbNm2SJGUyGZWXl2v79u266667Rv2cvGGFHd39g/r+Y206P8zVcFmRtOfem3iHnxGEnR/zZwZwh9NvWNHd3a3e3l7V19dn70ulUrrxxhu1f//+YY/JZDIaGBjIuWHi0ScaDn244TEDJFWkQdvb2ytJKi8vz7m/vLw8+7FLNTc3q7S0NHurrublIRvoEw2HPtzwmAGSyshvHV/6XenFf2sZzubNm5VOp7O3nh6+a7WBPtFw6MMNjxkgqSIN2oqKCkka8tNrX1/fkJ9yL0qlUiopKcm5YeLRJxoOfbjhMQMkVaRBW1tbq4qKCrW2tmbvO3PmjNrb27VkyZIovxQiRp9oOPThhscMkFTj3t5z+vRp/eMf/8j+vbu7WwcPHtT06dM1c+ZMbdiwQQ8//LBmzZqlWbNm6eGHH9aXvvQl/eQnP4n0xBG9/11Yrbr/mU6faIHCzo/5MwMk07i397S1tWnZsmVD7l+zZo2eeeaZ7BtW7NixI+cNK+bOnTumz8/2HgCATRS/AwBgkNP7aAEAQC6CFgAAgwhaAAAMImgBADCIoAUAwCCCFgAAg8b9hhXASHwv7bZd/O77/CVmAPewjxaRGa60OwgCb0q7w67f9vFJwAwQBfbRwknd/YO6f/chnQ+kc+eDnD837T6kY/2Dtk/RqLDrt318EjADuIqgRSR8L+22Xfzu+/wlZgB3ORO0LS0tmjNnjurq6myfCgrge2m37eJ33+cvMQO4y5mgbWpqUldXlzo6OmyfCgrge2m37eJ33+cvMQO4y5mgRbz5Xtptu/jd9/lLzADuImgRCd9Lu20Xv/s+f4kZwF1s70GkjvUPel3aHXb9to9PAmaAsOijBQDAIPbRAgAQIwQtAAAGEbQAABhE0AIAYBBBCwCAQQQtAAAG0UcLp/jeJWq7j9b3+UvMANFjHy2c4XuXqO0+Wt/nLzEDfI59tEgk37tEbffR+j5/iRnAHIIWTvC9S9R2H63v85eYAcxxJmjpo/Wb712itvtofZ+/xAxgjjNBSx+t33zvErXdR+v7/CVmAHOcCVr4zfcuUdt9tL7PX2IGMIeghRN87xK13Ufr+/wlZgBz2N4Dp/jeJWq7j9b3+UvMAPTRAgBgFPtoAQCIEYIWAACDCFoAAAwiaAEAMIigBQDAIIIWAACD6KNFovjeJWq7j9b3+UvMAEOxjxaJ4XuXqO0+Wt/nLzGDpGAfLTAM37tEbffR+j5/iRkgP4IWieB7l6jtPlrf5y8xA+TnTNDSR4swfO8Std1H6/v8JWaA/JwJWvpoEYbvXaK2+2h9n7/EDJCfM0ELhOF7l6jtPlrf5y8xA+RH0CIRfO8Std1H6/v8JWaA/Njeg0TxvUvUdh+t7/OXmEES0EcLAIBB7KMFACBGCFoAAAwiaAEAMIigBQDAIIIWAACDCFoAAAyijxb4At+7ROmjtY8ZJg/7aIELfO8SpY/WPmboBvbRAgb43iVKH619zDC5CFpAdInSR2sfM0wuZ4KWPlrY5HuXKH209jHD5HImaOmjhU2+d4nSR2sfM0wuZ4IWsMn3LlH6aO1jhslF0AKiS5Q+WvuYYXKxvQf4At+7ROmjtY8Z2kcfLQAABrGPFgCAGCFoAQAwiKAFAMAgghYAAIMIWgAADBp30O7bt08rVqxQVVWVioqK9Morr+R8fO3atSoqKsq5LV68OKrzBQAgVsbdRzs4OKj58+frjjvu0MqVK4d9zPLly7Vz587s36dMmVL4GQIxQpdoeHTihuP7+l007qBtaGhQQ0PDiI9JpVKqqKgo+KSAOBquS3RH+/t0iY5D2Bn6/hz4vn5XGfk32ra2NpWVlWn27Nm688471dfXZ+LLAM6gSzQ8OnHD8X39Los8aBsaGvT8889rz549euyxx9TR0aGbb75ZmUxm2MdnMhkNDAzk3IC4oUs0PDpxw/F9/S4b90vHo1m1alX2v+fOnauFCxeqpqZGr732mm6//fYhj29ubta2bduiPg1gQtElGh6duOH4vn6XGd/eU1lZqZqaGh09enTYj2/evFnpdDp76+nhuy7ED12i4dGJG47v63eZ8aA9efKkenp6VFlZOezHU6mUSkpKcm5A3NAlGh6duOH4vn6XjTtoT58+rYMHD+rgwYOSpO7ubh08eFDHjx/X6dOndd999+ntt9/WsWPH1NbWphUrVmjGjBm67bbboj53wBl0iYZHJ244vq/fZeOuyWtra9OyZcuG3L9mzRo99dRTamxsVGdnpz7++GNVVlZq2bJlevDBB1VdPbbvpqjJQ5zRJRoenbjh+L7+KNBHCwCAQfTRAgAQIwQtAAAGEbQAABhE0AIAYBBBCwCAQQQtAAAGRf5exwAKR5eo/T5a358D39dvAvtoAUcM1yUaBIFXXaJhZ2D7+Ljzff0XsY8WSCC6RO330fr+HPi+fpMIWsABdIna76P1/Tnwff0mORO0LS0tmjNnjurq6myfCjDh6BK130fr+3Pg+/pNciZom5qa1NXVpY6ODtunAkw4ukTt99H6/hz4vn6TnAlawGd0idrvo/X9OfB9/SYRtIAD6BK130fr+3Pg+/pNYnsP4BC6RO330fr+HPi+fok+WgAAjGIfLQAAMULQAgBgEEELAIBBBC0AAAYRtAAAGETQAgBgEH20QILQJWq/j9b358D39Q+HfbRAQtAlar+P1vfnICnrZx8tgCHoErXfR+v7c+D7+kdC0AIJQJeo/T5a358D39c/EmeClj5aoHB0idrvo/X9OfB9/SNxJmjpowUKR5eo/T5a358D39c/EmeCFkDh6BK130fr+3Pg+/pHQtACCUCXqP0+Wt+fA9/XPxK29wAJQpeo/T5a35+DJKyfPloAAAxiHy0AADFC0AIAYBBBCwCAQQQtAAAGEbQAABhE0AIAYBB9tACy6BK130fr+3OQxPWzjxaApOR0iYZhu4/W9+fAlfWzjxZA5OgStd9H6/tzkOT1E7QA6BKV/T5a35+DJK/fmaCljxawhy5R+320vj8HSV6/M0FLHy1gD12i9vtofX8Okrx+Z4IWgD10idrvo/X9OUjy+glaAHSJyn4fre/PQZLXz/YeAFlJ6BINy3Yfre/PgQvrp48WAACD2EcLAECMELQAABhE0AIAYBBBCwCAQQQtAAAGEbQAABhEHy2AyCSxS3S86KO1y8X5sY8WQCRc6RK1iT5au6KaH/toATgnyV2iY0UfrV0uz4+gBRBakrtEx4o+Wrtcnp8zQUsfLRBfSe4SHSv6aO1yeX7OBC19tEB8JblLdKzoo7XL5fk5E7QA4ivJXaJjRR+tXS7Pj6AFEFqSu0THij5au1yeH9t7AETGhS5R2+ijtSuK+dFHCwCAQeyjBQAgRghaAAAMImgBADCIoAUAwCCCFgAAgwhaAAAMGlfQNjc3q66uTsXFxSorK1NjY6OOHDmS85ggCLR161ZVVVVp6tSpuummm/Tuu+9GetIAkqu7f1Db3/g/rf9Tp7a/8X/qprVmXMLOj/lHb1z7aJcvX67Vq1errq5OZ8+e1QMPPKDDhw+rq6tL06Z9viF4+/bt+uUvf6lnnnlGs2fP1kMPPaR9+/bpyJEjKi4uHvVrsI8W8Bd9rOHQhxsNp96w4l//+pfKysrU3t6upUuXKggCVVVVacOGDdq0aZMkKZPJqLy8XNu3b9ddd9016uckaAE/dfcP6vuPten8MP9HuqxI2nPvTbxD0gjCzo/5/5dTb1iRTqclSdOnT5ckdXd3q7e3V/X19dnHpFIp3Xjjjdq/f/+wnyOTyWhgYCDnBsA/LveJxgF9uO4qOGiDINDGjRt1ww03aO7cuZKk3t5eSVJ5eXnOY8vLy7Mfu1Rzc7NKS0uzt+pqf16eAPBfLveJxgF9uO4qOGjXrVunQ4cO6U9/+tOQj136XdHF1/qHs3nzZqXT6eytp4fvmgAfudwnGgf04bqroKBdv369Xn31Ve3du1dXX3119v6KigpJGvLTa19f35Cfci9KpVIqKSnJuQHwj8t9onFAH667xhW0QRBo3bp1eumll7Rnzx7V1tbmfLy2tlYVFRVqbW3N3nfmzBm1t7dryZIl0ZwxgERyuU80DujDdde4fuv4nnvu0QsvvKC//OUv+uY3v5m9v7S0VFOnfv6ywvbt29Xc3KydO3dq1qxZevjhh9XW1sb2HgBjQh9rOPThhmd1e0++1+937typtWvXSvr8p95t27Zpx44dOnXqlBYtWqSWlpbsL0yNhqAFANjk1D5aEwhaAIBNTu2jBQAAIyNoAQAwiKAFAMAgghYAAIMIWgAADLrc9gkAQJS6+wf15y/sA/3RwmrVerQPNOz6bR+fRGzvAZAYvvep2u6jTcr82d4DAMPo7h/U/bsP6XwgnTsf5Py5afchHesftH2KRoVdv+3jk4ygBZAIvvep2u6j9X3+I3EmaFtaWjRnzhzV1dXZPhUAMeR7n6rtPlrf5z8SZ4K2qalJXV1d6ujosH0qAGLI9z5V2320vs9/JM4ELQCE4Xufqu0+Wt/nPxKCFkAi+N6naruP1vf5j4TtPQASxfc+Vdt9tEmYPzV5AAAYxD5aAABihKAFAMAgghYAAIMIWgAADCJoAQAwiKAFAMAg+mgB4At871O13UebxPmzjxYALkhKn2qhbPfRujJ/9tECgAG+96na7qNN8vwJWgAQfaq2+2iTPH9ngpY+WgA2+d6naruPNsnzdyZo6aMFYJPvfaq2+2iTPH9nghYAbPK9T9V2H22S50/QAoDoU7XdR5vk+bO9BwC+IAl9qmHY7qN1Yf700QIAYBD7aAEAiBGCFgAAgwhaAAAMImgBADCIoAUAwCCCFgAAg+ijBYAIJbFPdTxs99G6OH/20QJARFzpU7XFdh9tVPNnHy0AOCjJfapjYbuP1uX5E7QAEIEk96mOhe0+Wpfn70zQ0kcLIM6S3Kc6Frb7aF2evzNBSx8tgDhLcp/qWNjuo3V5/s4ELQDEWZL7VMfCdh+ty/MnaAEgAknuUx0L2320Ls+f7T0AECEX+lRtst1HG8X86aMFAMAg9tECABAjBC0AAAYRtAAAGETQAgBgEEELAIBBBC0AAAbRRwsADnGxT3Ui0Uc7AdhHC8BX9NnSRwsAMMTlPtWJQB8tAMAol/tUJwJ9tBOAPloAPnO5T3Ui0Ec7AeijBeAzl/tUJwJ9tAAAo1zuU50I9NECAIxyuU91ItBHO4HY3gPAZ/TZ0kdrHEELALCJfbQAAMQIQQsAgEEELQAABhG0AAAYRNACAGAQQQsAgEHjCtrm5mbV1dWpuLhYZWVlamxs1JEjR3Ies3btWhUVFeXcFi9eHOlJAwAQF+MK2vb2djU1Nemdd95Ra2urzp49q/r6eg0O5tYPLV++XCdOnMjeXn/99UhPGgCAuLh8PA9+4403cv6+c+dOlZWV6cCBA1q6dGn2/lQqpYqKimjOEACAGAv1b7TpdFqSNH369Jz729raVFZWptmzZ+vOO+9UX19fmC8DAEBsFfwWjEEQ6NZbb9WpU6f01ltvZe/ftWuXvvzlL6umpkbd3d36xS9+obNnz+rAgQNKpVJDPk8mk1Emk8n+PZ1Oa+bMmerp6eEtGAEAE25gYEDV1dX6+OOPVVpaGv4TBgW65557gpqamqCnp2fEx3300UfB5MmTg927dw/78S1btgSSuHHjxo0bN6du77//fqERmaOgn2jXr1+vV155Rfv27VNtbe2oj581a5Z++tOfatOmTUM+dulPtOfPn9e///1vXXnllXlLfEdz8buRMD8V19XVhSqhj/PxzI/52T6eGXIN2jz+4iurp06d0le+8pWCz+Gicf0yVBAEWr9+vV5++WW1tbWNKWRPnjypnp4eVVZWDvvxVCo15CXlKBYmSSUlJQVfZJMmTQr10nXcj5eYH/OzOz+JGXIN2p3fZZdF81YT4/osTU1N+uMf/6gXXnhBxcXF6u3tVW9vrz777DNJ0unTp3Xffffp7bff1rFjx9TW1qYVK1ZoxowZuu222yI54YnS1NTk9fFh2T5/28eHZfv8bR8fBdtrsH18WLbP3/bxURrXS8f5XsrduXOn1q5dq88++0yNjY3q7OzUxx9/rMrKSi1btkwPPvigqqurIzvp0dBpGw7zC4f5hccMw2F+4UQ9v3G/dDySqVOn6s033wx1QlFIpVLasmXLsL/ljNExv3CYX3jMMBzmF07U8yt4ew8AABgdpQIAABhE0AIAYBBBCwCAQQQtAAAGEbQAABhE0AIAYBBBCwCAQQQtAAAGEbQAABhE0AIAYBBBCwCAQf8PZz9j30w+mxkAAAAASUVORK5CYII=\n", "text/plain": [ "Graphics object consisting of 1 graphics primitive" ] }, "execution_count": 31, "metadata": {}, "output_type": "execute_result" } ], "source": [ "A=matrix(RDF,n*n,sparse=true);\n", "d=1/10;\n", "for i in range(0,n-1):\n", " for j in range(0,n):\n", " A[i+j*n,i+1+j*n]=d; # vecino a la derecha (i+1,j)\n", " A[i+1+j*n,i+j*n]=d; # vecino a la izquierda (i-1,j)\n", " A[i*n+j,(i+1)*n+j]=d; # vecino superior (i,j+1)\n", " A[(i+1)*n+j,i*n+j]=d; # vecino inferior (i,j-1)\n", "for i in range(0,n*n):\n", " A[i,i]=1-4*d; # posición original (i,j)\n", "matrix_plot(A)" ] }, { "cell_type": "markdown", "id": "34a680a3", "metadata": {}, "source": [ "Para difuminar la imagen, escribimos la matriz como vector y multiplicamos ese vector por la matriz de difuminado" ] }, { "cell_type": "code", "execution_count": 32, "id": "37e9e05c", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAWMAAAFtCAYAAADSyAuRAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAShklEQVR4nO3dX2hT9//H8VesmIomqX+mTBpFNlE7qbCaQUU3p65QRNyudiEi+3PhiGLpzaZebAxGvBobWItuw10MVxmbfy6m2DG1Dimm1aJYEAShHerEsSa1YMR6vhe/n2WdzvXEnnPeyXk+IBeJJz3vU7Ynx0/TjxHHcRwBAAI1IegBAADEGABMIMYAYAAxBgADiDEAGECMAcAAYgwABhBjADCAGAOAAcQYAAwIfYz37t2r+fPnq7KyUnV1dTp79mzQI3mmo6ND69ev15w5cxSJRHTkyJGgR/JUJpNRKpVSLBbTrFmz9Oabb+rq1atBj+Wp1tZW1dbWKh6PKx6Pq76+XsePHw96LF9lMhlFIhE1NTUFPYoroY7xoUOH1NTUpF27dunixYtauXKlGhsb1dfXF/RonhgaGtLSpUu1Z8+eoEfxxZkzZ5ROp9XZ2an29nY9ePBADQ0NGhoaCno0z1RXV2v37t3q6upSV1eXVq9erQ0bNujKlStBj+aLbDar/fv3q7a2NuhR3HNC7JVXXnG2bNky6rVFixY5H330UUAT+UeSc/jw4aDH8NXt27cdSc6ZM2eCHsVX06ZNc77++uugx/Dc4OCgs2DBAqe9vd157bXXnO3btwc9kiuhvTO+f/++uru71dDQMOr1hoYGnTt3LqCp4KVcLidJmj59esCT+GN4eFhtbW0aGhpSfX190ON4Lp1Oa926dVq7dm3QoxRlYtADBOXOnTsaHh7W7NmzR70+e/Zs3bp1K6Cp4BXHcdTc3KwVK1ZoyZIlQY/jqcuXL6u+vl737t3T1KlTdfjwYdXU1AQ9lqfa2tp04cIFZbPZoEcpWmhj/EgkEhn13HGcx15D6du6dasuXbqk3377LehRPLdw4UL19PRoYGBAP/74ozZv3qwzZ86UbZD7+/u1fft2nTx5UpWVlUGPU7TQxnjmzJmqqKh47C749u3bj90to7Rt27ZNx44dU0dHh6qrq4Mex3OTJk3Siy++KElatmyZstmsvvzyS+3bty/gybzR3d2t27dvq66ubuS14eFhdXR0aM+ePSoUCqqoqAhwwrEJ7ZrxpEmTVFdXp/b29lGvt7e3a/ny5QFNhfHkOI62bt2qn376Sb/++qvmz58f9EiBcBxHhUIh6DE8s2bNGl2+fFk9PT0jj2XLlmnjxo3q6ekpiRBLIb4zlqTm5mZt2rRJy5YtU319vfbv36++vj5t2bIl6NE8cffuXV27dm3k+fXr19XT06Pp06dr7ty5AU7mjXQ6rYMHD+ro0aOKxWIjfwtKJBKaPHlywNN5Y+fOnWpsbFQymdTg4KDa2tp0+vRpnThxIujRPBOLxR77OcCUKVM0Y8aM0vr5QLAf5gheS0uLM2/ePGfSpEnOyy+/XNYfezp16pQj6bHH5s2bgx7NE0+6VknOgQMHgh7NM+++++7If8/PPfecs2bNGufkyZNBj+W7UvxoW8Rx+AdJASBooV0zBgBLiDEAGECMAcAAYgwABhBjADCAGAOAAcQYAAwIfYwLhYI++eSTsv510X/imsOBay4tof+lj3w+r0QioVwup3g8HvQ4vuCaueZyVcrXXBJ3xi0tLZ4e75Yf84Ttmv34Hnn99a0dXwxr11AO1zxmwf429tgsXrzYs+NzuZwjycnlcibm8eMcFq/Z6+9ROVyz2+O5ZhszjZXvu7Y9fPhQN27cUCwWG/Mm7sPDw8rn82M+h5vjHx3n1dcv5nivz2Hxmr3+HpXDNbs9nmsOfibHcTQ4OKg5c+ZowoSnL0T4vmb8+++/K5lM+nlKAAhUf3//f/7DBr7fGcdiMUn/N1ypLbADgBv5fF7JZHKke0/je4wfLU3E43FiDCAUxrIkWxKfpgCAckeMAcCAUP8beBYMDAwEPQLgiaqqqqBHKCncGQOAAcQYAAwgxgBgADEGAAOIMQAYQIwBwADfYtzS0qKamhqlUim/TgkAJcP3jYJKefNnL/A5Y5QrPmfsrncsUwCAAcQYAAwgxgBgADEGAAOIMQAYQIwBwABiDAAGEGMAMIAYA4ABxBgADGBvCgAwgL0pAsbeFChX7E3B3hQAUHKIMQAYQIwBwABiDAAGEGMAMIAYA4ABRcV47969mj9/viorK1VXV6ezZ8+O91wAECquY3zo0CE1NTVp165dunjxolauXKnGxkb19fV5MR8AhILrGH/++ed677339P7772vx4sX64osvlEwm1dra6sV8ABAKrmJ8//59dXd3q6GhYdTrDQ0NOnfu3BPfUygUlM/nRz0AAKO5ivGdO3c0PDys2bNnj3p99uzZunXr1hPfk8lklEgkRh7JZLL4aQGgTBX1A7xIJDLqueM4j732yI4dO5TL5UYe/f39xZwSAMraRDcHz5w5UxUVFY/dBd++ffuxu+VHotGootFo8RMCQAi4ujOeNGmS6urq1N7ePur19vZ2LV++fFwHA4AwcXVnLEnNzc3atGmTli1bpvr6eu3fv199fX3asmWLF/MBQCi4jvHbb7+tP//8U59++qlu3rypJUuW6Oeff9a8efO8mA8AQoHN5QPG5vIoV2wuz+byAFByiDEAGECMAcAAYgwABhBjADDAtxi3tLSopqZGqVTKr1MCQMngo20B46NtKFd8tI2PtgFAySHGAGAAMQYAA4gxABhAjAHAAGIMAAYQYwAwgBgDgAHEGAAMIMYAYAAxBgAD2CgIAAxgo6CAsVEQyhUbBbFREACUHGIMAAYQYwAwYGLQA/xdGNdPp02bFvQI8MGCBQuCHsF358+fD3oE3z3LOjl3xgBgADEGAAOIMQAYQIwBwABiDAAGEGMAMIC9KQDAAN9inE6n1dvbq2w269cpAaBksEwBAAYQYwAwgBgDgAHEGAAMIMYAYAAxBgADiDEAGECMAcAAYgwABhBjADCAvSkAwAD2pgAAA1imAAADiDEAGECMAcAAYgwABhBjADCAGAOAAa5j3NHRofXr12vOnDmKRCI6cuSIB2MBQLi4jvHQ0JCWLl2qPXv2eDEPAITSRLdvaGxsVGNjoxezAEBouY6xW4VCQYVCYeR5Pp/3+pQAUHI8/wFeJpNRIpEYeSSTSa9PCQAlx/MY79ixQ7lcbuTR39/v9SkBoOR4vkwRjUYVjUa9Pg0AlDQ+ZwwABri+M757966uXbs28vz69evq6enR9OnTNXfu3HEdDgDCwnWMu7q69Prrr488b25uliRt3rxZ33777bgNBgBh4jrGq1atkuM4XswCAKHFmjEAGECMAcAAYgwABhBjADCAGAOAAb7FuKWlRTU1NUqlUn6dEgBKhm8xTqfT6u3tVTab9euUAFAyWKYAAAOIMQAYQIwBwABiDAAGEGMAMIAYA4ABxBgADCDGAGAAMQYAA4gxABjA3hQAYAB7UwCAASxTAIABxBgADCDGAGAAMQYAAyYGPQDC55dffgl6BN9VV1cHPYLvqqqqgh6hpHBnDAAGEGMAMIAYA4ABxBgADCDGAGAAMQYAA9goCAAMYKMgADCAZQoAMIAYA4ABxBgADCDGAGAAMQYAA4gxABhAjAHAAGIMAAYQYwAwgBgDgAHsTQEABrA3BQAYwDIFABhAjAHAAGIMAAYQYwAwgBgDgAGuYpzJZJRKpRSLxTRr1iy9+eabunr1qlezAUBouIrxmTNnlE6n1dnZqfb2dj148EANDQ0aGhryaj4ACIWJbg4+ceLEqOcHDhzQrFmz1N3drVdffXVcBwOAMHmmNeNcLidJmj59+rgMAwBh5erO+O8cx1Fzc7NWrFihJUuW/OtxhUJBhUJh5Hk+ny/2lABQtoq+M966dasuXbqk77///qnHZTIZJRKJkUcymSz2lABQtoqK8bZt23Ts2DGdOnVK1dXVTz12x44dyuVyI4/+/v6iBgWAcuZqmcJxHG3btk2HDx/W6dOnNX/+/P98TzQaVTQaLXpAAAgDVzFOp9M6ePCgjh49qlgsplu3bkmSEomEJk+e7MmAABAGrpYpWltblcvltGrVKj3//PMjj0OHDnk1HwCEgutlCgDA+GNvCgAwgBgDgAHEGAAMIMYAYAAxBgADiDEAGOBbjFtaWlRTU6NUKuXXKQGgZPgW43Q6rd7eXmWzWb9OCQAlg2UKADCAGAOAAcQYAAwgxgBgADEGAAOIMQAYQIwBwABiDAAGEGMAMIAYA4AB7E0BAAawNwUAGMAyBQAYQIwBwABiDAAGEGMAMGBi0AMgfNauXRv0CL5bsGBB0CP47vz580GP4Luqqqqi38udMQAYQIwBwABiDAAGEGMAMIAYA4AB7E0BAAawNwUAGMAyBQAYQIwBwABiDAAGEGMAMIAYA4ABxBgADCDGAGAAMQYAA4gxABhAjAHAAGIMAAawURAAGMBGQQBgAMsUAGAAMQYAA4gxABhAjAHAAGIMAAa4inFra6tqa2sVj8cVj8dVX1+v48ePezUbAISGqxhXV1dr9+7d6urqUldXl1avXq0NGzboypUrXs0HAKEw0c3B69evH/X8s88+U2trqzo7O/XSSy+N62AAECauYvx3w8PD+uGHHzQ0NKT6+vrxnAkAQsd1jC9fvqz6+nrdu3dPU6dO1eHDh1VTU/OvxxcKBRUKhZHn+Xy+uEkBoIy5/jTFwoUL1dPTo87OTn3wwQfavHmzent7//X4TCajRCIx8kgmk880MACUo4jjOM6zfIG1a9fqhRde0L59+57450+6M04mk8rlcorH46OOHRgYeJZRStK0adOCHgE+WLBgQdAj+O78+fNBj+C7qqqqUc/z+bwSicQTe/dPRa8ZP+I4zqjY/lM0GlU0Gn3W0wBAWXMV4507d6qxsVHJZFKDg4Nqa2vT6dOndeLECa/mA4BQcBXjP/74Q5s2bdLNmzeVSCRUW1urEydO6I033vBqPgAIBVcx/uabb7yaAwBCjb0pAMAAYgwABhBjADCAGAOAAcQYAAwgxgBggG8xbmlpUU1NjVKplF+nBICS4VuM0+m0ent7lc1m/TolAJQMlikAwABiDAAGEGMAMIAYA4ABxBgADCDGAGAAMQYAA4gxABhAjAHAAGIMAAawNwUAGMDeFABgAMsUAGAAMQYAA4gxABhAjAHAgIlBD/B3VVVVQY/gu7/++ivoEQBPhPH/52fBnTEAGECMAcAAYgwABhBjADCAGAOAAexNAQAGRBzHcfw8YT6fVyKRUC6XUzwe9/PUJg0MDAQ9AuAJPtrmrncsUwCAAcQYAAwgxgBgADEGAAOIMQAYQIwBwABiDAAGEGMAMIAYA4ABxBgADGBvCgAwgL0pAsbeFChX7E3B3hQAUHKIMQAYQIwBwABiDAAGEGMAMIAYA4ABzxTjTCajSCSipqamcRoHAMKp6Bhns1nt379ftbW14zkPAIRSUTG+e/euNm7cqK+++krTpk0b75kAIHSKinE6nda6deu0du3a8Z4HAEJpots3tLW16cKFC8pms2M6vlAoqFAojDzP5/NuTwkAZc/VnXF/f7+2b9+u7777TpWVlWN6TyaTUSKRGHkkk8miBgWAcuZqo6AjR47orbfeUkVFxchrw8PDikQimjBhggqFwqg/k558Z5xMJtko6P+xURDKFRsFudsoyNUyxZo1a3T58uVRr73zzjtatGiRPvzww8dCLEnRaFTRaNTNaQAgdFzFOBaLacmSJaNemzJlimbMmPHY6wCAseM38ADAANefpvin06dPj8MYABBu3BkDgAHEGAAMIMYAYAAxBgADiDEAGECMAcAA32Lc0tKimpoapVIpv04JACXD1d4U48HN72qHAXtToFyxN4W73rFMAQAGEGMAMIAYA4ABxBgADCDGAGAAMQYAA4gxABhAjAHAAGIMAAYQYwAwgL0pAMAA9qYIGHtToFyxNwV7UwBAySHGAGAAMQYAA4gxABgwMegBwo4fcgCQuDMGABOIMQAYQIwBwABiDAAG+P4DvEe/8JfP5/0+NQD46lHnxvKLzr7HeHBwUJKUTCb9PjUABGJwcFCJROKpx/i+N8XDhw9148YNxWIxRSKRMb0nlUopm82O+Rxujs/n80omk+rv7x/zXhlezuPHOSxes9ffo3K4ZrfHc83BX7PjOBocHNScOXM0YcLTV4V9vzOeMGGCqqurXb2noqLC1aZCbo+XpHg8Pub3+DFP2K7Zj++RVNrXXMz1Slxz0DP91x3xIyXxA7x0Ou3p8W75MU/YrtmP75HXX9/a8cWwdg3lcM1j5fsyhTVh3NKTa+aay1UpX3NJ3Bl7KRqN6uOPP1Y0Gg16FN9wzeHANZeW0N8ZA4AFob8zBgALiDEAGECMAcAAYgwABhBjADCAGAOAAcQYAAwgxgBgwP8ACxoFb9E7fEoAAAAASUVORK5CYII=\n", "text/plain": [ "Graphics object consisting of 1 graphics primitive" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "Bv=matrix_to_vector(B) # Necesitamos considerar la imagen como un vector\n", "b=A*Bv # difuminamos\n", "show(matrix_plot(vector_to_matrix(b)), figsize=5) # Para dibujar necesitamos volver a pasar a la forma de matriz" ] }, { "cell_type": "markdown", "id": "fc2d6ab1", "metadata": {}, "source": [ "Y podemos hacer el proceso inverso resolviendo el sistema $Nx=b$, para lo cual Gauss-Seidel es un método apropiado (a cada paso se perfila un poco más la imagen). Véamoslo en el ejemplo.\n", "\n", "Construimos las matrices de Gauss-Seidel" ] }, { "cell_type": "code", "execution_count": 33, "id": "5925f1db", "metadata": {}, "outputs": [], "source": [ "N = matrix(RDF,n*n);\n", "for i in range(0,n*n):\n", " for j in range(i+1,n*n):\n", " N[i,j] =- A[i,j];\n", "M = A + N;" ] }, { "cell_type": "markdown", "id": "3d0ea2dd", "metadata": {}, "source": [ "Y aplicamos el método" ] }, { "cell_type": "code", "execution_count": 34, "id": "f7453d2b", "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAcwAAAHWCAYAAAD6lrl7AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMywgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/P9b71AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAVu0lEQVR4nO3dX4iVdf7A8c+o6xg1c1gLBfEUsgu7DmKLOsRERX9MGCJyr7oIsaUbYZJkbpbqYpdYmKCbXXCU3F3qYgljabVgS3ZgUYsI1JKVLoIgcJY0K2hmHNjjOp7fxY/mh7+0Pmeec85zzvh6wVx4mmee73w7+PZ5zjnfb0+9Xq8HAPC9lpQ9AADoBoIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYLZoH379sW6detixYoVsXnz5nj33XfLHlLXOH78eDz66KOxZs2a6OnpicOHD5c9pK4zNjYWg4OD0dfXF6tWrYrt27fHJ598Uvawusr+/ftj48aN0d/fH/39/TE0NBTvvPNO2cPqamNjY9HT0xN79uwpeygtJZgNeP3112PPnj3x/PPPx0cffRT33ntvDA8Px9mzZ8seWleYnZ2NO++8M/bu3Vv2ULrWsWPHYmRkJD744IOYmJiIy5cvx7Zt22J2drbsoXWNtWvXxosvvhgnT56MkydPxoMPPhiPPfZYfPzxx2UPrSudOHEiDhw4EBs3bix7KC3XY/H1vLvuuis2bdoU+/fvn39s/fr1sX379hgbGytxZN2np6cnDh06FNu3by97KF3tyy+/jFWrVsWxY8fivvvuK3s4XWvlypXx0ksvxVNPPVX2ULrKxYsXY9OmTbFv37743e9+F7/4xS/i97//fdnDahlXmEmXLl2KU6dOxbZt2656fNu2bfH++++XNCpudFNTUxHxv3/h07i5ubk4ePBgzM7OxtDQUNnD6TojIyPxyCOPxNatW8seSlssK3sA3eKrr76Kubm5WL169VWPr169Os6fP1/SqLiR1ev1GB0djXvuuSc2bNhQ9nC6ypkzZ2JoaCj+85//xC233BKHDh2KgYGBsofVVQ4ePBgffvhhnDhxouyhtI1gNqinp+eqP9fr9e88Bu3w9NNPx7/+9a947733yh5K1/nZz34Wp0+fjm+++SbeeOON2LlzZxw7dkw0kyYnJ+OZZ56Jf/zjH7FixYqyh9M2gpl02223xdKlS79zNXnhwoXvXHVCq+3evTveeuutOH78eKxdu7bs4XSd5cuXx09/+tOIiNiyZUucOHEi/vCHP8TLL79c8si6w6lTp+LChQuxefPm+cfm5ubi+PHjsXfv3qjVarF06dISR9gaXsNMWr58eWzevDkmJiauenxiYiLuvvvukkbFjaZer8fTTz8df/vb3+Kf//xnrFu3ruwhLQr1ej1qtVrZw+gaDz30UJw5cyZOnz49/7Vly5Z44okn4vTp04sylhGuMBsyOjoaO3bsiC1btsTQ0FAcOHAgzp49G7t27Sp7aF3h4sWL8emnn87/+bPPPovTp0/HypUr4/bbby9xZN1jZGQkXnvttXjzzTejr69v/o5HpVKJm266qeTRdYfnnnsuhoeHo1qtxszMTBw8eDCOHj0aR44cKXtoXaOvr+87r5vffPPNceutty7q19MFswGPP/54fP311/HCCy/EuXPnYsOGDfH222/HHXfcUfbQusLJkyfjgQcemP/z6OhoRETs3LkzXn311ZJG1V2+/UjT/ffff9Xjr7zySjz55JPtH1AX+uKLL2LHjh1x7ty5qFQqsXHjxjhy5Eg8/PDDZQ+NDudzmACQ4DVMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQzAbVarX47W9/axmtAsxhceawOHNY3I02hxYuaND09HRUKpWYmpqK/v7+sofTlcxhceawOHNY3I02h115hTk+Pl7q8UU14/xlz4E5LP/4Zij7dyj7+GYo+3co+/hm6Jrfod6F1q9fX9rxU1NT9YioT01NlXL+Zv0Mc2gOi46h2483h8WPXyxzmNX2xdevXLkSn3/+efT19S144+W5ubmYnp5e8BiKHP/tcWWdv1k/wxyaw6Jj6PbjzWHx4xfDHNbr9ZiZmYk1a9bEkiXff9O17a9h/vvf/45qtdrOUwLA95qcnPzBzdjbfoXZ19cXEf87uBvhRWIAOtf09HRUq9X5Nn2ftgfz29uw/f39gglAR8i8RNiV75IFgHYTTABIEEwASBBMAEho+5t+Mi5fvlz2ELreD32eiB+20M8J83/a/Km1RctzsbhmzKG/VQEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIKFtwRwfH4+BgYEYHBxs1ykBoGnavoH09PR0VCqVmJqauu72Xlb6Kc5KP8VZXaU4K/00h+dicdebw0yTvuVvVQBIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASLAfJgAk2A9zkbIfZnH2ICzOfpjN4blYnP0wAaBNBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASFhTMffv2xbp162LFihWxefPmePfdd5s9LgDoKA0H8/XXX489e/bE888/Hx999FHce++9MTw8HGfPnm3F+ACgIzS8+Ppdd90VmzZtiv37988/tn79+ti+fXuMjY394PEWX28Pi68XZ8Hr4iy+3hyei8W1ffH1S5cuxalTp2Lbtm1XPb5t27Z4//33r3lMrVaL6enpq74AoNs0FMyvvvoq5ubmYvXq1Vc9vnr16jh//vw1jxkbG4tKpTL/Va1WFz5aACjJgu7b/f9L23q9ft3L3WeffTampqbmvyYnJxdySgAo1bJGvvm2226LpUuXfudq8sKFC9+56vxWb29v9Pb2LnyEANABGrrCXL58eWzevDkmJiauenxiYiLuvvvupg4MADpJQ1eYERGjo6OxY8eO2LJlSwwNDcWBAwfi7NmzsWvXrlaMDwA6QsPBfPzxx+Prr7+OF154Ic6dOxcbNmyIt99+O+64445WjA8AOkLDn8Msyucw28PnMIvz2bfifA6zOTwXi2v75zAB4EYlmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJDQtmCOj4/HwMBADA4OtuuUANA01pJdpKwlW5z1O4uzlmxzeC4WZy1ZAGgTwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABPthAkCC/TAXKfthFmcPwuLsh9kcnovF2Q8TANpEMAEgQTABIEEwASBBMAEgYVnZA7gW7/Asbtmyjvxf21W8w7M4z8PmqNVqZQ+h6zXjncbKBAAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACfbDBICEjtwP88qVK+0c0qJkhZXirPRTnOdhc1jpp7jrrSBnP0wAaDLBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAE+2ECQIL9MBcp+xAWZz/M4jwPm8N+mMXZDxMA2kQwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIKHhYB4/fjweffTRWLNmTfT09MThw4dbMCwA6CwNB3N2djbuvPPO2Lt3byvGAwAdqeGVkYeHh2N4eLgVYwGAjtXyrQRqtdpVK+1PT0+3+pQA0HQtf9PP2NhYVCqV+a9qtdrqUwJA07U8mM8++2xMTU3Nf01OTrb6lADQdC2/Jdvb2xu9vb2tPg0AtJTPYQJAQsNXmBcvXoxPP/10/s+fffZZnD59OlauXBm33357UwcHAJ2i4WCePHkyHnjggfk/j46ORkTEzp0749VXX23awACgkzQczPvvvz/q9XorxgIAHctrmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQ0LZgjo+Px8DAQAwODrbrlADQND31Ni8MOz09HZVKJaampqK/v/+a33PlypV2DmlRWras5VudLnrWTC7O87A5arVa2UPoekuWXPv6MNOk+Z/RioEBwGIjmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQYD9MAEiwH+YiZR/C4uyHWZznYXPYD7M4+2ECQJsIJgAkCCYAJAgmACQIJgAkeAsbXId3yRbnHe/N0dPTU/YQCFeYAJAimACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQYD9MAEiwH+YiZR/C4jwPizOHzWGln+KuN4f2wwSAJhNMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEiwHyYAJNgPc5GyH2ZxnofFmcPmsB9mcfbDBIA2EUwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIaCiYY2NjMTg4GH19fbFq1arYvn17fPLJJ60aGwB0jIaCeezYsRgZGYkPPvggJiYm4vLly7Ft27aYnZ1t1fgAoCMUWnz9yy+/jFWrVsWxY8fivvvuSx1j8fX2sPh6cZ6HxZnD5rD4enHNWHy90N+qU1NTERGxcuXK635PrVaLWq121eAAoNss+E0/9Xo9RkdH45577okNGzZc9/vGxsaiUqnMf1Wr1YWeEgBKs+BbsiMjI/H3v/893nvvvVi7du11v+9aV5jVatUt2RZzS7Y4z8PizGFzuCVbXGm3ZHfv3h1vvfVWHD9+/HtjGRHR29sbvb29CzkNAHSMhoJZr9dj9+7dcejQoTh69GisW7euVeMCgI7SUDBHRkbitddeizfffDP6+vri/PnzERFRqVTipptuaskAAaATNPQa5vXuAb/yyivx5JNPpn6Gj5W0h9cwi/M8LM4cNofXMItr+2uYBT6yCQBdzVqyAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJDQtmCOj4/HwMBADA4OtuuUANA0C95AeqEsvt4eFl8vzvOwOHPYHBZfL64Zi6+7JQsACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAmCCQAJggkACYIJAAn2wwSABPthLlL2wyzO87A4c9gc9sMszn6YANAmggkACYIJAAmCCQAJggkACR35VkrvCCtuyRL/FirKHBZnDpvjv//9b9lDIFxhAkCKYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgv0wASChI/fDbPOQFqUf/ehHZQ+h69nLsTgr/TSHlX6Ksx8mALSJYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgv0wASDBfpiLlP0wi7MfZnH2w2wO+2EWZz9MAGgTwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBICEhoK5f//+2LhxY/T390d/f38MDQ3FO++806qxAUDHaCiYa9eujRdffDFOnjwZJ0+ejAcffDAee+yx+Pjjj1s1PgDoCIUXX1+5cmW89NJL8dRTT6W+3+Lr7WHx9eIsvl6cxdebw+LrxTVj8fVlCz353Nxc/PWvf43Z2dkYGhq67vfVarWo1WpXDQ4Auk3D//w7c+ZM3HLLLdHb2xu7du2KQ4cOxcDAwHW/f2xsLCqVyvxXtVotNGAAKEPDt2QvXboUZ8+ejW+++SbeeOON+NOf/hTHjh27bjSvdYVZrVbdkm0xt2SLc0u2OLdkm8Mt2eKacUu28GuYW7dujZ/85Cfx8ssvp77fa5jtIZjFCWZxgtkcgllcR2wgXa/Xr7qCBIDFqKE3/Tz33HMxPDwc1Wo1ZmZm4uDBg3H06NE4cuRIq8YHAB2hoWB+8cUXsWPHjjh37lxUKpXYuHFjHDlyJB5++OFWjQ8AOkJDwfzzn//cqnEAQEfzijwAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJLQtmOPj4zEwMBCDg4PtOiUANE3h/TAbZT/M9rAfZnH2wyzOfpjNYT/M4jpiP0wAuBEIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAk2A8TABLsh7lI2Q+zOPthFmc/zOawH2Zx9sMEgDYRTABIEEwASBBMAEgQTABIEEwASFhW9gBoDW9DB2guV5gAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkGA/TABIsB8mAIue/TABoE0EEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwAS7IcJAAn2wwRg0bMfJgC0iWACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAQqFgjo2NRU9PT+zZs6dJwwGAzrTgYJ44cSIOHDgQGzdubOZ4AKAjLSiYFy9ejCeeeCL++Mc/xo9//ONmjwkAOs6CgjkyMhKPPPJIbN26tdnjAYCOtKzRAw4ePBgffvhhnDhxIvX9tVotarXa/J+np6cbPSUAlK6hK8zJycl45pln4i9/+UusWLEidczY2FhUKpX5r2q1uqCBAkCZGtpA+vDhw/HLX/4yli5dOv/Y3Nxc9PT0xJIlS6JWq1313yKufYVZrVZtIA1A2zRjA+mGbsk+9NBDcebMmase+9WvfhU///nP49e//vV3YhkR0dvbG729vY2cBgA6TkPB7Ovriw0bNlz12M033xy33nrrdx4HgMXESj8AkNDwu2T/v6NHjzZhGADQ2VxhAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkCCYAJAgmACQIJgAkBC24I5Pj4eAwMDMTg42K5TAkDTNLQfZjNk9h6zHyYAzdSM/TDdkgWABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgATBBIAEwQSABMEEgAT7YQJAgv0wAVj07IcJAG0imACQIJgAkCCYAJAgmACQIJgAkLCs7AFcy/Xe/gsAZXGFCQAJggkACYIJAAmCCQAJggkACYIJAAlt/1jJtzuRTE9Pt/vUAHCVb1uU2SWr7cGcmZmJiIhqtdruUwPANc3MzESlUvne72n7fphXrlyJzz//PPr6+ha8QMHg4GCcOHFiwWMocvz09HRUq9WYnJz8wb3TWnH+Zv0Mc2gOi46h2483h8WPXwxzWK/XY2ZmJtasWRNLlnz/q5Rtv8JcsmRJrF27ttDPWLp0aaH/OUWPj4jo7+9f8M9oxvnLngNzWP7xEcXmsBlj6PbjI8yhOYwfvLL8Vle+6WdkZKTU44tqxvnLngNzWP7xzVD271D28c1Q9u9Q9vHN0C2/Q9tvyXa76enpqFQqMTU1VfhfZTcqc1icOSzOHBZ3o81hV15hlqm3tzd+85vfRG9vb9lD6VrmsDhzWJw5LO5Gm0NXmACQ4AoTABIEEwASBBMAEgQTABIEEwASBBMAEgQTABIEEwAS/gfGj7LypkOEZQAAAABJRU5ErkJggg==\n", "text/plain": [ "Graphics object consisting of 1 graphics primitive" ] }, "execution_count": 34, "metadata": {}, "output_type": "execute_result" } ], "source": [ "x0 = vector(RDF,n*n)\n", "m = 3 # Número de iteraciones del método\n", "for i in range(1,m+1):\n", " x0 = M.solve_right(N*x0+b)\n", "matrix_plot(vector_to_matrix(vector(x0)))" ] }, { "cell_type": "code", "execution_count": null, "id": "45c41e22", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "SageMath 9.5", "language": "sage", "name": "sagemath" }, "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.11.2" } }, "nbformat": 4, "nbformat_minor": 5 }