{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "# General Linear Elliptic Equation" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The `GeneralMG2d` class implements support for a general elliptic equation of the form:\n", "\n", "$$\\alpha \\phi + \\nabla \\cdot (\\beta \\nabla \\phi) + \\gamma \\cdot \\nabla \\phi = f$$\n", "\n", "with inhomogeneous boundary conditions.\n", "\n", "It subclasses the `CellCenterMG2d` class, and the basic interface is the same" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "We will solve the above with\n", "\n", "\\begin{align}\n", "\\alpha &= 10 \\\\\n", "\\beta &= xy + 1 \\\\\n", "\\gamma &= \\hat{x} + \\hat{y}\n", "\\end{align}\n", "\n", "and\n", "\n", "\\begin{align}\n", "f = &-\\frac{\\pi}{2}(x + 1)\\sin\\left( \\frac{\\pi y}{2}\\right) \\cos\\left(\\frac{\\pi x}{2}\\right ) \\\\\n", "{} &-\\frac{\\pi}{2}(y + 1)\\sin\\left( \\frac{\\pi x}{2}\\right) \\cos\\left(\\frac{\\pi y}{2}\\right ) \\\\\n", "{} &+ \\left( \\frac{-\\pi^2 (xy+1)}{2} + 10\\right ) \\cos \\left ( \\frac{\\pi x}{2} \\right ) \\cos\\left ( \\frac{\\pi y}{2} \\right )\n", "\\end{align}\n", "\n", "\n", "on $[0, 1] \\times [0,1]$ with boundary conditions:\n", "\n", "\\begin{align}\n", "\\phi(x=0) &= \\cos(\\pi y/2) \\\\\n", "\\phi(x=1) &= 0 \\\\\n", "\\phi(y=0) &= \\cos(\\pi x/2) \\\\\n", "\\phi(y=1) &= 0\n", "\\end{align}\n", "\n", "\n", "This has the exact solution:\n", "\n", "$$\\phi = \\cos(\\pi x/2) \\cos(\\pi y/2)$$" ] }, { "cell_type": "code", "execution_count": 1, "metadata": { "tags": [] }, "outputs": [], "source": [ "import numpy as np\n", "import matplotlib.pyplot as plt" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Setting up the solver" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "tags": [] }, "outputs": [], "source": [ "import pyro.mesh.boundary as bnd\n", "import pyro.mesh.patch as patch\n", "import pyro.multigrid.general_MG as gMG" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "For reference, we'll define a function providing the analytic solution" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "tags": [] }, "outputs": [], "source": [ "def true(x,y): \n", " return np.cos(np.pi*x/2.0)*np.cos(np.pi*y/2.0) " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now the coefficents--note that since $\\gamma$ is a vector, we have a different function for each component" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "tags": [] }, "outputs": [], "source": [ "def alpha(x,y): \n", " return 10.0*np.ones_like(x) \n", " \n", "def beta(x,y): \n", " return x*y + 1.0 \n", " \n", "def gamma_x(x,y): \n", " return np.ones_like(x) \n", " \n", "def gamma_y(x,y): \n", " return np.ones_like(x) " ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "and the righthand side function" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "tags": [] }, "outputs": [], "source": [ "def f(x,y):\n", " return -0.5*np.pi*(x + 1.0)*np.sin(np.pi*y/2.0)*np.cos(np.pi*x/2.0) - \\\n", " 0.5*np.pi*(y + 1.0)*np.sin(np.pi*x/2.0)*np.cos(np.pi*y/2.0) + \\\n", " (-np.pi**2*(x*y+1.0)/2.0 + 10.0)*np.cos(np.pi*x/2.0)*np.cos(np.pi*y/2.0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Our inhomogeneous boundary conditions require a function that can be evaluated on the boundary to give the value" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "tags": [] }, "outputs": [], "source": [ "def xl_func(y):\n", " return np.cos(np.pi*y/2.0)\n", "\n", "def yl_func(x):\n", " return np.cos(np.pi*x/2.0)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we can setup our grid object and the coefficients, which are stored as a `CellCenter2d` object. Note, the coefficients do not need to have the same boundary conditions as $\\phi$ (and for real problems, they may not). The one that matters the most is $\\beta$, since that will need to be averaged to the edges of the domain, so the boundary conditions on the coefficients are important.\n", "\n", "Here we use Neumann boundary conditions" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "tags": [] }, "outputs": [], "source": [ "nx = 128\n", "ny = 128\n", "\n", "g = patch.Grid2d(nx, ny, ng=1)\n", "d = patch.CellCenterData2d(g)\n", "\n", "bc_c = bnd.BC(xlb=\"neumann\", xrb=\"neumann\",\n", " ylb=\"neumann\", yrb=\"neumann\")\n", "\n", "d.register_var(\"alpha\", bc_c)\n", "d.register_var(\"beta\", bc_c)\n", "d.register_var(\"gamma_x\", bc_c)\n", "d.register_var(\"gamma_y\", bc_c)\n", "d.create()\n", "\n", "a = d.get_var(\"alpha\")\n", "a[:,:] = alpha(g.x2d, g.y2d)\n", "\n", "b = d.get_var(\"beta\")\n", "b[:,:] = beta(g.x2d, g.y2d)\n", "\n", "gx = d.get_var(\"gamma_x\")\n", "gx[:,:] = gamma_x(g.x2d, g.y2d)\n", "\n", "gy = d.get_var(\"gamma_y\")\n", "gy[:,:] = gamma_y(g.x2d, g.y2d)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "Now we can setup the multigrid object" ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "cc data: nx = 2, ny = 2, ng = 1\n", " nvars = 7\n", " variables:\n", " v: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " f: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " r: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " alpha: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " beta: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_x: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_y: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", "\n", "cc data: nx = 4, ny = 4, ng = 1\n", " nvars = 7\n", " variables:\n", " v: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " f: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " r: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " alpha: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " beta: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_x: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_y: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", "\n", "cc data: nx = 8, ny = 8, ng = 1\n", " nvars = 7\n", " variables:\n", " v: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " f: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " r: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " alpha: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " beta: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_x: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_y: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", "\n", "cc data: nx = 16, ny = 16, ng = 1\n", " nvars = 7\n", " variables:\n", " v: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " f: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " r: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " alpha: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " beta: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_x: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_y: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", "\n", "cc data: nx = 32, ny = 32, ng = 1\n", " nvars = 7\n", " variables:\n", " v: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " f: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " r: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " alpha: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " beta: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_x: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_y: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", "\n", "cc data: nx = 64, ny = 64, ng = 1\n", " nvars = 7\n", " variables:\n", " v: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " f: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " r: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " alpha: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " beta: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_x: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_y: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", "\n", "cc data: nx = 128, ny = 128, ng = 1\n", " nvars = 7\n", " variables:\n", " v: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " f: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " r: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: dirichlet +x: dirichlet -y: dirichlet +y: dirichlet \n", " alpha: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " beta: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_x: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", " gamma_y: min: 0.0000000000 max: 0.0000000000\n", " BCs: -x: neumann +x: neumann -y: neumann +y: neumann \n", "\n" ] } ], "source": [ "a = gMG.GeneralMG2d(nx, ny,\n", " xl_BC_type=\"dirichlet\", yl_BC_type=\"dirichlet\",\n", " xr_BC_type=\"dirichlet\", yr_BC_type=\"dirichlet\",\n", " xl_BC=xl_func,\n", " yl_BC=yl_func,\n", " coeffs=d,\n", " verbose=1, vis=0, true_function=true)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "just as before, we specify the righthand side and initialize the solution" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Source norm = 1.775181492337501\n" ] } ], "source": [ "a.init_zeros()\n", "a.init_RHS(f(a.x2d, a.y2d))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Solving the system" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "tags": [ "nbval-ignore-output" ] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "source norm = 1.775181492337501\n", "<<< beginning V-cycle (cycle 1) >>>\n", "\n", " level = 6, nx = 128, residual change: 1.77518 → 188.933\n", " level = 5, nx = 64, residual change: 129.938 → 56.2871\n", " level = 4, nx = 32, residual change: 38.8869 → 18.7228\n", " level = 3, nx = 16, residual change: 12.9261 → 6.74186\n", " level = 2, nx = 8, residual change: 4.64648 → 2.06513\n", " level = 1, nx = 4, residual change: 1.37453 → 0.0224452\n", " bottom solve\n", " level = 1, nx = 4, residual change: 0.0312525 → 0.0312525\n", " level = 2, nx = 8, residual change: 2.80598 → 2.80598\n", " level = 3, nx = 16, residual change: 8.7724 → 8.7724\n", " level = 4, nx = 32, residual change: 19.591 → 19.591\n", " level = 5, nx = 64, residual change: 50.4641 → 50.4641\n", " level = 6, nx = 128, residual change: 160.213 → 160.213\n", "cycle 1: relative err = 0.9999999999999981, residual err = 2.323786088373021\n", "\n", "<<< beginning V-cycle (cycle 2) >>>\n", "\n", " level = 6, nx = 128, residual change: 4.12514 → 2.42473\n", " level = 5, nx = 64, residual change: 1.69154 → 1.04862\n", " level = 4, nx = 32, residual change: 0.728342 → 0.455482\n", " level = 3, nx = 16, residual change: 0.316533 → 0.221286\n", " level = 2, nx = 8, residual change: 0.153325 → 0.0747197\n", " level = 1, nx = 4, residual change: 0.0497494 → 0.000813357\n", " bottom solve\n", " level = 1, nx = 4, residual change: 0.00113252 → 0.00113252\n", " level = 2, nx = 8, residual change: 0.101526 → 0.101526\n", " level = 3, nx = 16, residual change: 0.298147 → 0.298147\n", " level = 4, nx = 32, residual change: 0.521885 → 0.521885\n", " level = 5, nx = 64, residual change: 0.991063 → 0.991063\n", " level = 6, nx = 128, residual change: 2.04419 → 2.04419\n", "cycle 2: relative err = 0.036315310129800826, residual err = 0.03283823443993396\n", "\n", "<<< beginning V-cycle (cycle 3) >>>\n", "\n", " level = 6, nx = 128, residual change: 0.0582938 → 0.0417201\n", " level = 5, nx = 64, residual change: 0.0292467 → 0.0233563\n", " level = 4, nx = 32, residual change: 0.0163063 → 0.0129066\n", " level = 3, nx = 16, residual change: 0.00901111 → 0.00731526\n", " level = 2, nx = 8, residual change: 0.0050815 → 0.00256253\n", " level = 1, nx = 4, residual change: 0.00170641 → 2.79124e-05\n", " bottom solve\n", " level = 1, nx = 4, residual change: 3.88653e-05 → 3.88653e-05\n", " level = 2, nx = 8, residual change: 0.00348191 → 0.00348191\n", " level = 3, nx = 16, residual change: 0.010065 → 0.010065\n", " level = 4, nx = 32, residual change: 0.0160323 → 0.0160323\n", " level = 5, nx = 64, residual change: 0.0243037 → 0.0243037\n", " level = 6, nx = 128, residual change: 0.0377753 → 0.0377753\n", "cycle 3: relative err = 0.0012532978372415558, residual err = 0.0006216334987521017\n", "\n", "<<< beginning V-cycle (cycle 4) >>>\n", "\n", " level = 6, nx = 128, residual change: 0.00110351 → 0.000889832\n", " level = 5, nx = 64, residual change: 0.00062574 → 0.00060774\n", " level = 4, nx = 32, residual change: 0.000426042 → 0.000397674\n", " level = 3, nx = 16, residual change: 0.000278462 → 0.000242683\n", " level = 2, nx = 8, residual change: 0.000168818 → 8.63435e-05\n", " level = 1, nx = 4, residual change: 5.75013e-05 → 9.40799e-07\n", " bottom solve\n", " level = 1, nx = 4, residual change: 1.30997e-06 → 1.30997e-06\n", " level = 2, nx = 8, residual change: 0.000117324 → 0.000117324\n", " level = 3, nx = 16, residual change: 0.000338509 → 0.000338509\n", " level = 4, nx = 32, residual change: 0.000524953 → 0.000524953\n", " level = 5, nx = 64, residual change: 0.000708087 → 0.000708087\n", " level = 6, nx = 128, residual change: 0.000918517 → 0.000918517\n", "cycle 4: relative err = 4.257466296364851e-05, residual err = 1.4967652930826935e-05\n", "\n", "<<< beginning V-cycle (cycle 5) >>>\n", "\n", " level = 6, nx = 128, residual change: 2.65703e-05 → 2.30982e-05\n", " level = 5, nx = 64, residual change: 1.62749e-05 → 1.79061e-05\n", " level = 4, nx = 32, residual change: 1.25859e-05 → 1.28807e-05\n", " level = 3, nx = 16, residual change: 9.03506e-06 → 8.103e-06\n", " level = 2, nx = 8, residual change: 5.6415e-06 → 2.90121e-06\n", " level = 1, nx = 4, residual change: 1.93217e-06 → 3.16168e-08\n", " bottom solve\n", " level = 1, nx = 4, residual change: 4.40233e-08 → 4.40233e-08\n", " level = 2, nx = 8, residual change: 3.94227e-06 → 3.94227e-06\n", " level = 3, nx = 16, residual change: 1.14059e-05 → 1.14059e-05\n", " level = 4, nx = 32, residual change: 1.7696e-05 → 1.7696e-05\n", " level = 5, nx = 64, residual change: 2.28172e-05 → 2.28172e-05\n", " level = 6, nx = 128, residual change: 2.72045e-05 → 2.72045e-05\n", "cycle 5: relative err = 1.437223355768636e-06, residual err = 4.2910353907176844e-07\n", "\n", "<<< beginning V-cycle (cycle 6) >>>\n", "\n", " level = 6, nx = 128, residual change: 7.61737e-07 → 6.88796e-07\n", " level = 5, nx = 64, residual change: 4.8583e-07 → 5.69884e-07\n", " level = 4, nx = 32, residual change: 4.01145e-07 → 4.28873e-07\n", " level = 3, nx = 16, residual change: 3.01132e-07 → 2.72291e-07\n", " level = 2, nx = 8, residual change: 1.89676e-07 → 9.77049e-08\n", " level = 1, nx = 4, residual change: 6.50717e-08 → 1.06486e-09\n", " bottom solve\n", " level = 1, nx = 4, residual change: 1.48271e-09 → 1.48271e-09\n", " level = 2, nx = 8, residual change: 1.32767e-07 → 1.32767e-07\n", " level = 3, nx = 16, residual change: 3.85631e-07 → 3.85631e-07\n", " level = 4, nx = 32, residual change: 6.03884e-07 → 6.03884e-07\n", " level = 5, nx = 64, residual change: 7.68242e-07 → 7.68242e-07\n", " level = 6, nx = 128, residual change: 8.86509e-07 → 8.86509e-07\n", "cycle 6: relative err = 4.849259894834445e-08, residual err = 1.3530556515124825e-08\n", "\n", "<<< beginning V-cycle (cycle 7) >>>\n", "\n", " level = 6, nx = 128, residual change: 2.40192e-08 → 2.21253e-08\n", " level = 5, nx = 64, residual change: 1.56138e-08 → 1.88696e-08\n", " level = 4, nx = 32, residual change: 1.32927e-08 → 1.44857e-08\n", " level = 3, nx = 16, residual change: 1.01772e-08 → 9.19808e-09\n", " level = 2, nx = 8, residual change: 6.40947e-09 → 3.30184e-09\n", " level = 1, nx = 4, residual change: 2.19906e-09 → 3.59875e-11\n", " bottom solve\n", " level = 1, nx = 4, residual change: 5.01092e-11 → 5.01092e-11\n", " level = 2, nx = 8, residual change: 4.48679e-09 → 4.48679e-09\n", " level = 3, nx = 16, residual change: 1.30812e-08 → 1.30812e-08\n", " level = 4, nx = 32, residual change: 2.0705e-08 → 2.0705e-08\n", " level = 5, nx = 64, residual change: 2.62808e-08 → 2.62808e-08\n", " level = 6, nx = 128, residual change: 2.99444e-08 → 2.99444e-08\n", "cycle 7: relative err = 1.6392149576904378e-09, residual err = 4.458207725000789e-10\n", "\n", "<<< beginning V-cycle (cycle 8) >>>\n", "\n", " level = 6, nx = 128, residual change: 7.91413e-10 → 7.35586e-10\n", " level = 5, nx = 64, residual change: 5.1922e-10 → 6.36466e-10\n", " level = 4, nx = 32, residual change: 4.4855e-10 → 4.92822e-10\n", " level = 3, nx = 16, residual change: 3.4637e-10 → 3.11941e-10\n", " level = 2, nx = 8, residual change: 2.17418e-10 → 1.11945e-10\n", " level = 1, nx = 4, residual change: 7.45572e-11 → 1.22015e-12\n", " bottom solve\n", " level = 1, nx = 4, residual change: 1.69894e-12 → 1.69894e-12\n", " level = 2, nx = 8, residual change: 1.52121e-10 → 1.52121e-10\n", " level = 3, nx = 16, residual change: 4.44914e-10 → 4.44914e-10\n", " level = 4, nx = 32, residual change: 7.10977e-10 → 7.10977e-10\n", " level = 5, nx = 64, residual change: 9.034e-10 → 9.034e-10\n", " level = 6, nx = 128, residual change: 1.0238e-09 → 1.0238e-09\n", "cycle 8: relative err = 5.555097426033948e-11, residual err = 1.5072807373286882e-11\n", "\n" ] } ], "source": [ "a.solve(rtol=1.e-10)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Checking the result\n", "\n", "We can compare to the true solution" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "tags": [] }, "outputs": [], "source": [ "v = a.get_solution()\n", "b = true(a.x2d, a.y2d)\n", "e = v - b" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "The norm of the error is" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "tags": [] }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ " 1.671934405e-05\n" ] } ], "source": [ "print(f\"{e.norm():20.10g}\")" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "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.11.2" } }, "nbformat": 4, "nbformat_minor": 4 }