{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Populating the interactive namespace from numpy and matplotlib\n" ] } ], "source": [ "%pylab inline" ] }, { "cell_type": "code", "execution_count": 112, "metadata": {}, "outputs": [], "source": [ "import pandas\n", "import sympy as Ω\n", "import matplotlib.pyplot as plt\n", "import math\n", "import pylab as p\n", "import numpy as np\n", "import decimal as d\n", "from scipy.misc import derivative\n", "import pickle\n", "from tqdm import tqdm_notebook as tqdm\n", "\n", "from dataclasses import dataclass\n", "\n", "import scipy.sparse\n", "import scipy.sparse.linalg\n", "\n", "import functools" ] }, { "cell_type": "code", "execution_count": 3, "metadata": {}, "outputs": [], "source": [ "Ω.init_printing()" ] }, { "cell_type": "code", "execution_count": 4, "metadata": {}, "outputs": [], "source": [ "l_float = 20\n", "T_float = 100\n", "\n", "D_float = 0.6\n", "H_float = 4" ] }, { "cell_type": "code", "execution_count": 113, "metadata": {}, "outputs": [], "source": [ "with open('/Users/abra/projects/study/7semestr/degtyarev/roots128.pickle', 'rb') as file:\n", " roots128 = pickle.load(file)\n", "\n", "roots = roots128\n", "\n", "def u_elem128(x, t, alpha_n):\n", " x = np.float128(x)\n", " t = np.float128(t)\n", " alpha_n = np.float128(alpha_n)\n", " return -4*np.pi**2*(alpha_n*np.cos(alpha_n*x/10) + 40*np.sin(alpha_n*x/10))*(alpha_n*np.sin(2*alpha_n) - 40*np.cos(2*alpha_n) + 40)*np.exp(-3*alpha_n**2*t/500)/(5*(alpha_n**2 - np.pi**2)*(alpha_n**2*(4*alpha_n + np.sin(4*alpha_n)) - 80*alpha_n*(np.cos(4*alpha_n) - 1) + 6400*alpha_n - 1600*np.sin(4*alpha_n)))\n", "\n", "@functools.lru_cache(maxsize=100000)\n", "def u128(x, t, root_n=len(roots)):\n", " return sum(u_elem128(x, t, roots[:root_n]))" ] }, { "cell_type": "code", "execution_count": 6, "metadata": {}, "outputs": [], "source": [ "@dataclass\n", "class Grid:\n", " I: int\n", " K: int\n", " \n", " xs: ndarray\n", " ts: ndarray\n", " \n", " h_x: float\n", " h_t: float\n", " \n", " def __init__(self, I, K):\n", " self.I = I\n", " self.K = K\n", "\n", " self.xs = linspace(0, l_float, I + 1)\n", " self.ts = linspace(0, T_float, K + 1)\n", "\n", " self.h_x = l_float / I\n", " self.h_t = T_float / K\n", " \n", " self.gamma = D_float * self.h_t / (self.h_x ** 2)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": {}, "outputs": [], "source": [ "def ug_mse(ug1, ug2):\n", " return np.mean((ug1 - ug2) ** 2)\n", "\n", "def ug_mae(ug1, ug2):\n", " return np.max(np.abs(ug1 - ug2))" ] }, { "cell_type": "code", "execution_count": 64, "metadata": {}, "outputs": [], "source": [ "n = 5000\n", "\n", "np.random.seed(0)\n", "\n", "a = np.random.uniform(0, 0.1, size=n-1)\n", "b = np.ones(n)\n", "c = np.random.uniform(0, 0.1, size=n-1)\n", "\n", "z = np.random.uniform(0, 0.1, size=n)\n", "\n", "A = np.zeros((n, n))\n", "\n", "for i in range(n):\n", " A[i, i] = b[i]\n", " \n", " if i > 0:\n", " A[i, i - 1] = a[i - 1]\n", " \n", " if i < n - 1:\n", " A[i, i + 1] = c[i]" ] }, { "cell_type": "code", "execution_count": 65, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "1.52 s ± 78.3 ms per loop (mean ± std. dev. of 5 runs, 5 loops each)\n" ] } ], "source": [ "%%timeit -n 5 -r 5\n", "np.linalg.solve(A, z)" ] }, { "cell_type": "code", "execution_count": 66, "metadata": {}, "outputs": [], "source": [ "A_sp = scipy.sparse.diags([a, b, c], [-1, 0, 1], format='csc')" ] }, { "cell_type": "code", "execution_count": 67, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "3.1 ms ± 430 µs per loop (mean ± std. dev. of 5 runs, 5 loops each)\n" ] } ], "source": [ "%%timeit -n 5 -r 5\n", "scipy.sparse.linalg.spsolve(A_sp, z)" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## теория" ] }, { "cell_type": "code", "execution_count": 8, "metadata": {}, "outputs": [], "source": [ "def calc_theoretical(grid):\n", " ug = zeros((grid.K + 1, grid.I + 1), dtype=float)\n", " \n", " for k, t in enumerate(grid.ts):\n", " for i, x in enumerate(grid.xs):\n", " ug[k, i] = u128(x, t, 100)\n", " \n", " return ug" ] }, { "cell_type": "code", "execution_count": 9, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "58.89331758385291\n" ] } ], "source": [ "grid = Grid(10, 50)\n", "print(sum(calc_theoretical(grid)))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## простейшая явная схема" ] }, { "cell_type": "code", "execution_count": 10, "metadata": {}, "outputs": [], "source": [ "def calc_explicit(grid):\n", " if grid.gamma > 1/2:\n", " raise Exception('Bad gamma!')\n", " \n", " ug = zeros((grid.K + 1, grid.I + 1), dtype=float)\n", "\n", " ug[0, :] = sin(np.pi * grid.xs / l_float) ** 2 * 4 / 10\n", "\n", " for k in range(1, grid.K + 1):\n", " # middle\n", " ug[k, 1 : grid.I] = (\n", " grid.gamma * ug[k - 1, 0 : grid.I - 1] +\n", " (-2 * grid.gamma + 1) * ug[k - 1, 1 : grid.I] + \n", " grid.gamma * ug[k - 1, 2 : grid.I + 1]\n", " )\n", "\n", " # left\n", " ug[k, 0] = ug[k, 1] / (H_float * grid.h_x + 1)\n", "\n", " # right\n", " ug[k, grid.I] = ug[k, grid.I - 1] / (H_float * grid.h_x + 1)\n", " \n", " return ug" ] }, { "cell_type": "code", "execution_count": 11, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "58.72948448630509\n", "4.3733803061715146e-07\n" ] } ], "source": [ "grid = Grid(10, 50)\n", "ug = calc_explicit(grid)\n", "print(sum(ug))\n", "print(ug_mse(ug, calc_theoretical(grid)))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## простейшая неявная схема" ] }, { "cell_type": "code", "execution_count": 76, "metadata": {}, "outputs": [], "source": [ "def calc_implicit(grid):\n", " ug = zeros((grid.K + 1, grid.I + 1), dtype=float)\n", "\n", " ug[0, :] = sin(np.pi * grid.xs / l_float) ** 2 * 4 / 10\n", "\n", " for k in range(1, grid.K + 1):\n", " diags = {\n", " -1: np.zeros(grid.I + 1, dtype=float),\n", " 0: np.zeros(grid.I + 1, dtype=float),\n", " +1: np.zeros(grid.I + 1, dtype=float),\n", " }\n", " b = np.zeros(grid.I + 1, dtype=float)\n", "\n", " # left\n", " diags[ 0][0] = H_float + 1 / grid.h_x\n", " diags[+1][0] = -1 / grid.h_x\n", "\n", " # middle\n", " for i in range(1, grid.I):\n", " diags[-1][i] = grid.gamma\n", " diags[ 0][i] = - 2 * grid.gamma - 1\n", " diags[+1][i] = grid.gamma\n", " b[i] = -ug[k - 1, i]\n", "\n", " # right\n", " diags[ 0][grid.I] = H_float + 1 / grid.h_x\n", " diags[-1][grid.I] = -1 / grid.h_x\n", "\n", " A = scipy.sparse.diags([\n", " diags[-1][1:], \n", " diags[ 0], \n", " diags[+1][:-1],\n", " ], [-1, 0, 1], format='csc')\n", " \n", " ug[k, :] = scipy.sparse.linalg.spsolve(A, b)\n", " \n", " return ug" ] }, { "cell_type": "code", "execution_count": 69, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "59.49556361037571\n", "4.006492537092588e-06\n" ] } ], "source": [ "grid = Grid(10, 50)\n", "ug = calc_implicit(grid)\n", "print(sum(ug))\n", "print(ug_mse(ug, calc_theoretical(grid)))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## схема Кранка-Николсон" ] }, { "cell_type": "code", "execution_count": 79, "metadata": {}, "outputs": [], "source": [ "def calc_kn(grid):\n", " ug = zeros((grid.K + 1, grid.I + 1), dtype=float)\n", "\n", " ug[0, :] = sin(np.pi * grid.xs / l_float) ** 2 * 4 / 10\n", "\n", " for k in range(1, grid.K + 1):\n", " diags = {\n", " -1: np.zeros(grid.I + 1, dtype=float),\n", " 0: np.zeros(grid.I + 1, dtype=float),\n", " +1: np.zeros(grid.I + 1, dtype=float),\n", " }\n", " b = np.zeros(grid.I + 1, dtype=float)\n", "\n", " # left\n", " diags[ 0][0] = H_float + 1 / grid.h_x\n", " diags[+1][0] = -1 / grid.h_x\n", "\n", " # middle\n", " for i in range(1, grid.I):\n", " diags[-1][i] = grid.gamma\n", " diags[ 0][i] = - 2 * grid.gamma - 2\n", " diags[+1][i] = grid.gamma\n", " b[i] = -(\n", " ug[k - 1][i - 1] * grid.gamma +\n", " ug[k - 1][i] * (-2 * grid.gamma + 2) +\n", " ug[k - 1][i + 1] * grid.gamma\n", " )\n", "\n", " # right\n", " diags[ 0][grid.I] = H_float + 1 / grid.h_x\n", " diags[-1][grid.I] = -1 / grid.h_x\n", "\n", " A = scipy.sparse.diags([\n", " diags[-1][1:], \n", " diags[ 0], \n", " diags[+1][:-1],\n", " ], [-1, 0, 1], format='csc')\n", " \n", " ug[k, :] = scipy.sparse.linalg.spsolve(A, b)\n", " \n", " return ug" ] }, { "cell_type": "code", "execution_count": 15, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "0.048762671501371764\n", "5.115545125775158e-07\n" ] } ], "source": [ "grid = Grid(10, 50)\n", "ug = calc_kn(grid)\n", "print(mean(ug[-1]))\n", "print(ug_mse(ug, calc_theoretical(grid)))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## простейшая неявная схема повышенного порядка аппроксимации" ] }, { "cell_type": "code", "execution_count": 174, "metadata": {}, "outputs": [], "source": [ "def calc_implicit_improved(grid):\n", " ug = zeros((grid.K + 1, grid.I + 1), dtype=float)\n", "\n", " ug[0, :] = sin(np.pi * grid.xs / l_float) ** 2 * 4 / 10\n", "\n", " for k in range(1, grid.K + 1):\n", " diags = {\n", " -1: np.zeros(grid.I + 1, dtype=float),\n", " 0: np.zeros(grid.I + 1, dtype=float),\n", " +1: np.zeros(grid.I + 1, dtype=float),\n", " }\n", " b = np.zeros(grid.I + 1, dtype=float)\n", "\n", " # left\n", " diags[ 0][0] = (\n", " 2 * H_float * grid.gamma * grid.h_x + \n", " 2 * grid.gamma + \n", " 1\n", " )\n", " diags[+1][0] = - 2 * grid.gamma\n", " b[0] = ug[k - 1, 0]\n", "\n", " # middle\n", " for i in range(1, grid.I):\n", " diags[-1][i] = grid.gamma\n", " diags[ 0][i] = - 2 * grid.gamma - 1\n", " diags[+1][i] = grid.gamma\n", " b[i] = -ug[k - 1, i]\n", "\n", " # right\n", " diags[ 0][grid.I] = (\n", " 2 * H_float * grid.gamma * grid.h_x + \n", " 2 * grid.gamma + \n", " 1\n", " )\n", " diags[-1][grid.I] = - 2 * grid.gamma\n", " b[grid.I] = ug[k - 1, grid.I]\n", "\n", " A = scipy.sparse.diags([\n", " diags[-1][1:], \n", " diags[ 0], \n", " diags[+1][:-1],\n", " ], [-1, 0, 1], format='csc')\n", " \n", " ug[k, :] = scipy.sparse.linalg.spsolve(A, b)\n", " \n", " return ug" ] }, { "cell_type": "code", "execution_count": 177, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "59.486724613775095\n", "4.014426575945732e-06\n" ] } ], "source": [ "grid = Grid(10, 50)\n", "ug = calc_implicit_improved(grid)\n", "print(sum(ug))\n", "print(ug_mse(ug, calc_theoretical(grid)))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## Кранк-Николсон повышенного порядка аппроксимации" ] }, { "cell_type": "code", "execution_count": 174, "metadata": {}, "outputs": [], "source": [ "def calc_kn_improved(grid):\n", " ug = zeros((grid.K + 1, grid.I + 1), dtype=float)\n", "\n", " ug[0, :] = sin(np.pi * grid.xs / l_float) ** 2 * 4 / 10\n", "\n", " for k in range(1, grid.K + 1):\n", " diags = {\n", " -1: np.zeros(grid.I + 1, dtype=float),\n", " 0: np.zeros(grid.I + 1, dtype=float),\n", " +1: np.zeros(grid.I + 1, dtype=float),\n", " }\n", " b = np.zeros(grid.I + 1, dtype=float)\n", "\n", " # left\n", " #diags[ 0][0] = H_float + 1 / grid.h_x\n", " #diags[+1][0] = -1 / grid.h_x\n", " \n", " diags[ 0][0] = (\n", " 2 * H_float * grid.gamma * grid.h_x + \n", " 2 * grid.gamma + \n", " 1\n", " )\n", " diags[+1][0] = - 2 * grid.gamma\n", " b[0] = ug[k - 1, 0]\n", "\n", " # middle\n", " for i in range(1, grid.I):\n", " diags[-1][i] = grid.gamma\n", " diags[ 0][i] = - 2 * grid.gamma - 2\n", " diags[+1][i] = grid.gamma\n", " b[i] = -(\n", " ug[k - 1][i - 1] * grid.gamma +\n", " ug[k - 1][i] * (-2 * grid.gamma + 2) +\n", " ug[k - 1][i + 1] * grid.gamma\n", " )\n", "\n", " # right\n", " diags[ 0][grid.I] = H_float + 1 / grid.h_x\n", " diags[-1][grid.I] = -1 / grid.h_x\n", "\n", " A = scipy.sparse.diags([\n", " diags[-1][1:], \n", " diags[ 0], \n", " diags[+1][:-1],\n", " ], [-1, 0, 1], format='csc')\n", " \n", " ug[k, :] = scipy.sparse.linalg.spsolve(A, b)\n", " \n", " return ug" ] }, { "cell_type": "code", "execution_count": 177, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "59.486724613775095\n", "4.014426575945732e-06\n" ] } ], "source": [ "grid = Grid(10, 50)\n", "ug = calc_kn_improved(grid)\n", "print(sum(ug))\n", "print(ug_mse(ug, calc_theoretical(grid)))" ] }, { "cell_type": "markdown", "metadata": {}, "source": [ "## сравнение" ] }, { "cell_type": "code", "execution_count": 128, "metadata": {}, "outputs": [ { "data": { "image/png": "iVBORw0KGgoAAAANSUhEUgAAAO4AAAAPBAMAAAAR2O42AAAAMFBMVEX///8AAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMpndu3bvImbNiRBUq0Qb3U6NAAAACXBIWXMAAA7EAAAOxAGVKw4bAAADlklEQVQ4Eb2TT4hbVRTGf8nkz5u8JBN0J8jEwT906rSDFUrFYlDRjdqhIGqlTIRS0Rlq6KLKbPo2yuCmEV0IVYnakTEgBksLFiRv5a4kDnQYhdKg6LbS/zNJTc8596UdU9deeCffO/d835dz3r0A98gDC/UXiO2erNcNuXBq3zbwJz8skK6fCOGH5S+lsr6voAS4d6VerxkafT32qAENqpdZ0srhZVK3WY/rdrpIK/D6/f5NQxYokgsZDf2dHCZ7Az7jZIUfO/Eq/NyAY1LeJjsHqX6/JiKaBNHzv2CLwqFlUo4FyQO6m+wxNh2vQNuQhWxAusoH8C2zIZfxymSqPE2i5K9cEItHIMdPU10Y/e4VMbOk6SWKPDvkqa8mZSx5ee8TTeUP0axmZEAFQ+61hD/Ni9Bq/NXx10lIc938dS1nVnxrINNMim/OcpY0vYvSwH8sk6oZS9jmK1WtQMKf8jjUCvJXGpmAm3AxROc81ibfS5S1wlmQn7nbV/Wes5q7QiSlLLxK5Os/oXUleQxpuHBjGf+S+ErhyRrNIvmrYx8tyWGLfBOCrN/FvYUoaXrr9+3twJmXZaab1kBKWZzC+aa/3ipv0h+GLGT6x0lvwHiNhccaHK0R22h+Sioc+C7jfBON9D8Cdfiq568HzJPfyYIk76xICmXJRxrM+ZsONK1KERLO/nGtk5Z+xZfMIY4W1bdLXE6wWuD/LUH7lfWMPO6ji2+/wcONI6H/gG0NQiRlrGzhtm9K5F6zGkWk5rw2s3OD4XCw4eY8jSe9mW+yKnWR7/aKJZ3eZblkndX3X7U/9/EOWbsK8jfdJzPWGZxvukKiZ3fPkIVUhVhPz1UrfBf2hHKusr2RKp5cZfMdKUW+v8L3HUs6vV3iGwhxaJkUxjp7/vy1c7I91lVfX3YMWdChr+lFPNaQwe0JEyW8bq58p99mWyq037dgu/Qmf8bprWq/G0OuOCmMJXsHdD9VYuQ6MZmEIRcq8Lxe9i1a81TFKxOvxuT7lqN+x4vCVF/5lRbdELR2XL4vb8ueJDctk8JYkr0i5w8voDVDVnwNWcgfx5thJPTfYJHsVXiS0wV+43QYWbSKQlffw3hvCrThix65mj8vZ8X/9z1yUhgLTvTX2A2/T7wkV2he2IYs3D+5TU7D0v4C+YkHO3ohf5GDPSFh6p3PAzgSQnx1fY305IqM2SVVj6+mCiQfWhS1zcukjLU5+3/iW7ckZCmz5Ei9AAAAAElFTkSuQmCC\n", "text/latex": [ "$$4.3733803061715146e-07$$" ], "text/plain": [ "4.3733803061715146e-07" ] }, "execution_count": 128, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ug_mse(calc_theoretical(grid), calc_explicit(grid))" ] }, { "cell_type": "code", "execution_count": 169, "metadata": {}, "outputs": [], "source": [ "def print_comparison(Is, Ks, E, i_p, k_p, i_d, k_d):\n", " if (\n", " not 0 <= i_p < len(Is) or\n", " not 0 <= k_p < len(Ks) or\n", " not 0 <= i_p + i_d < len(Is) or\n", " not 0 <= k_p + k_d < len(Ks) or\n", " E[k_p, i_p] < 0 or\n", " E[k_p + k_d, i_p + i_d] < 0\n", " ):\n", " return\n", " \n", " error_ratio = E[k_p, i_p] / E[k_p + k_d, i_p + i_d]\n", " \n", " I1 = Is[i_p]\n", " I2 = Is[i_p + i_d]\n", " K1 = Ks[k_p]\n", " K2 = Ks[k_p + k_d]\n", " \n", " print(f'h_x x{I2//I1} h_t x{K2//K1} @ ({I1:5},{K1:5}): {error_ratio}')" ] }, { "cell_type": "code", "execution_count": 170, "metadata": {}, "outputs": [], "source": [ "def compare_errors(Is, Ks, calc_function):\n", " E = np.zeros((len(Ks), len(Is)))\n", "\n", " for k_p in range(len(Ks)):\n", " for i_p in range(len(Is)):\n", " grid = Grid(Is[i_p], Ks[k_p])\n", "\n", " try:\n", " error = ug_mae(calc_function(grid), calc_theoretical(grid))\n", " except:\n", " error = -1\n", "\n", " E[k_p, i_p] = error\n", " \n", " for k_p in range(len(Ks)):\n", " for i_p in range(len(Is)):\n", " print_comparison(Is, Ks, E, i_p, k_p, +1, +1)\n", " \n", " for k_p in range(len(Ks)):\n", " for i_p in range(len(Is)):\n", " print_comparison(Is, Ks, E, i_p, k_p, +2, +1)\n", " \n", " for k_p in range(len(Ks)):\n", " for i_p in range(len(Is)):\n", " print_comparison(Is, Ks, E, i_p, k_p, +1, +2)" ] }, { "cell_type": "code", "execution_count": 171, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "h_x x2 h_t x4 @ ( 16, 256): 1.919077767270527\n" ] } ], "source": [ "compare_errors([16, 32, 64, 128, 256], [256, 1024], calc_explicit)" ] }, { "cell_type": "code", "execution_count": 172, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "h_x x2 h_t x2 @ ( 16, 16): 1.8945839194805507\n", "h_x x2 h_t x2 @ ( 32, 16): 1.8480652483857702\n", "h_x x2 h_t x2 @ ( 64, 16): 1.825718428756318\n", "h_x x2 h_t x2 @ ( 128, 16): 1.8200879583854386\n", "h_x x2 h_t x2 @ ( 16, 32): 2.011003900865954\n", "h_x x2 h_t x2 @ ( 32, 32): 1.8852164156444766\n", "h_x x2 h_t x2 @ ( 64, 32): 1.8445734893796497\n", "h_x x2 h_t x2 @ ( 128, 32): 1.8438755692851074\n", "h_x x2 h_t x2 @ ( 16, 64): 2.17619657916294\n", "h_x x2 h_t x2 @ ( 32, 64): 1.9562383707853712\n", "h_x x2 h_t x2 @ ( 64, 64): 1.912181287104185\n", "h_x x2 h_t x2 @ ( 128, 64): 1.9058730801460957\n", "h_x x2 h_t x2 @ ( 16, 128): 2.3986401861806734\n", "h_x x2 h_t x2 @ ( 32, 128): 2.066009495516749\n", "h_x x2 h_t x2 @ ( 64, 128): 1.9741429057474853\n", "h_x x2 h_t x2 @ ( 128, 128): 1.955899186659026\n", "h_x x4 h_t x2 @ ( 16, 16): 1.9431575629889533\n", "h_x x4 h_t x2 @ ( 32, 16): 1.8496006164556726\n", "h_x x4 h_t x2 @ ( 64, 16): 1.8261475842007147\n", "h_x x4 h_t x2 @ ( 16, 32): 2.0776417503976705\n", "h_x x4 h_t x2 @ ( 32, 32): 1.8918649575362794\n", "h_x x4 h_t x2 @ ( 64, 32): 1.845407456579902\n", "h_x x4 h_t x2 @ ( 16, 64): 2.326489502477712\n", "h_x x4 h_t x2 @ ( 32, 64): 1.9755444903443895\n", "h_x x4 h_t x2 @ ( 64, 64): 1.9125944713395193\n", "h_x x4 h_t x2 @ ( 16, 128): 2.8326076010600714\n", "h_x x4 h_t x2 @ ( 32, 128): 2.1104815578650395\n", "h_x x4 h_t x2 @ ( 64, 128): 1.9752019588093679\n", "h_x x2 h_t x4 @ ( 16, 16): 3.45714271998841\n", "h_x x2 h_t x4 @ ( 32, 16): 3.396912364589268\n", "h_x x2 h_t x4 @ ( 64, 16): 3.3648762818042726\n", "h_x x2 h_t x4 @ ( 128, 16): 3.3552270369444783\n", "h_x x2 h_t x4 @ ( 16, 32): 3.6798635241236157\n", "h_x x2 h_t x4 @ ( 32, 32): 3.569646703206942\n", "h_x x2 h_t x4 @ ( 64, 32): 3.514763487476464\n", "h_x x2 h_t x4 @ ( 128, 32): 3.512604694405925\n", "h_x x2 h_t x4 @ ( 16, 64): 3.8072301037478553\n", "h_x x2 h_t x4 @ ( 32, 64): 3.780516278348026\n", "h_x x2 h_t x4 @ ( 64, 64): 3.7380285131624063\n", "h_x x2 h_t x4 @ ( 128, 64): 3.726890300676519\n" ] } ], "source": [ "compare_errors([16, 32, 64, 128, 256], [16, 32, 64, 128, 256], calc_implicit)" ] }, { "cell_type": "code", "execution_count": 173, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "h_x x2 h_t x2 @ ( 16, 16): 1.7867407245342732\n", "h_x x2 h_t x2 @ ( 32, 16): 2.0659339289789247\n", "h_x x2 h_t x2 @ ( 64, 16): 2.206577164218172\n", "h_x x2 h_t x2 @ ( 128, 16): 2.205398898130605\n", "h_x x2 h_t x2 @ ( 16, 32): 1.0260466758467883\n", "h_x x2 h_t x2 @ ( 32, 32): 1.985994321847054\n", "h_x x2 h_t x2 @ ( 64, 32): 2.0974055723617107\n", "h_x x2 h_t x2 @ ( 128, 32): 2.1318052621473274\n", "h_x x2 h_t x2 @ ( 16, 64): 1.8854175495471748\n", "h_x x2 h_t x2 @ ( 32, 64): 1.5501374544674071\n", "h_x x2 h_t x2 @ ( 64, 64): 2.0605293935478826\n", "h_x x2 h_t x2 @ ( 128, 64): 2.1241455078974423\n", "h_x x2 h_t x2 @ ( 16, 128): 2.2073516657031678\n", "h_x x2 h_t x2 @ ( 32, 128): 1.8202146022973351\n", "h_x x2 h_t x2 @ ( 64, 128): 1.8615590435162546\n", "h_x x2 h_t x2 @ ( 128, 128): 2.093178235555532\n", "h_x x4 h_t x2 @ ( 16, 16): 1.5416112730228488\n", "h_x x4 h_t x2 @ ( 32, 16): 1.9411758005230835\n", "h_x x4 h_t x2 @ ( 64, 16): 2.132375626592411\n", "h_x x4 h_t x2 @ ( 16, 32): 0.672923482269376\n", "h_x x4 h_t x2 @ ( 32, 32): 1.8096548816821525\n", "h_x x4 h_t x2 @ ( 64, 32): 2.0030692793517577\n", "h_x x4 h_t x2 @ ( 16, 64): 2.09685316770308\n", "h_x x4 h_t x2 @ ( 32, 64): 1.3513796667001667\n", "h_x x4 h_t x2 @ ( 64, 64): 1.9355394149338634\n", "h_x x4 h_t x2 @ ( 16, 128): 5.089679627450866\n", "h_x x4 h_t x2 @ ( 32, 128): 2.07031910687416\n", "h_x x4 h_t x2 @ ( 64, 128): 1.8247920519286784\n", "h_x x2 h_t x4 @ ( 16, 16): 5.410544492760973\n", "h_x x2 h_t x4 @ ( 32, 16): 4.755334696107873\n", "h_x x2 h_t x4 @ ( 64, 16): 4.925531450048806\n", "h_x x2 h_t x4 @ ( 128, 16): 4.865081194253896\n", "h_x x2 h_t x4 @ ( 16, 32): 1.4301343984858494\n", "h_x x2 h_t x4 @ ( 32, 32): 4.6940709147810455\n", "h_x x2 h_t x4 @ ( 64, 32): 4.742894620222057\n", "h_x x2 h_t x4 @ ( 128, 32): 4.7415271368204355\n", "h_x x2 h_t x4 @ ( 16, 64): 1.4883710757235502\n", "h_x x2 h_t x4 @ ( 32, 64): 2.537069293873301\n", "h_x x2 h_t x4 @ ( 64, 64): 4.399957273891053\n", "h_x x2 h_t x4 @ ( 128, 64): 4.733335280216436\n" ] } ], "source": [ "compare_errors([16, 32, 64, 128, 256], [16, 32, 64, 128, 256], calc_kn)" ] }, { "cell_type": "code", "execution_count": 181, "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "h_x x2 h_t x2 @ ( 32, 128): 2.038390823677269\n", "h_x x2 h_t x2 @ ( 64, 128): 1.9736276882133699\n", "h_x x2 h_t x2 @ ( 32, 256): 2.1614872398503246\n", "h_x x2 h_t x2 @ ( 64, 256): 2.025328231523804\n", "h_x x2 h_t x2 @ ( 32, 512): 2.373745419372439\n", "h_x x2 h_t x2 @ ( 64, 512): 2.0878950429120877\n", "h_x x2 h_t x2 @ ( 32, 1024): 2.713172301556474\n", "h_x x2 h_t x2 @ ( 64, 1024): 2.185936854335897\n", "h_x x4 h_t x2 @ ( 32, 128): 2.1195416132190292\n", "h_x x4 h_t x2 @ ( 32, 256): 2.330981592342115\n", "h_x x4 h_t x2 @ ( 32, 512): 2.7365187229137757\n", "h_x x4 h_t x2 @ ( 32, 1024): 3.543382713251337\n", "h_x x2 h_t x4 @ ( 32, 128): 3.8282183811233788\n", "h_x x2 h_t x4 @ ( 64, 128): 3.8442015881699336\n", "h_x x2 h_t x4 @ ( 32, 256): 3.9146871029344097\n", "h_x x2 h_t x4 @ ( 64, 256): 3.921190229218736\n", "h_x x2 h_t x4 @ ( 32, 512): 3.973114349424527\n", "h_x x2 h_t x4 @ ( 64, 512): 3.958967998390918\n" ] } ], "source": [ "compare_errors([32, 64, 128], [128, 256, 512, 1024, 2048], calc_implicit_improved)" ] }, { "cell_type": "code", "execution_count": null, "metadata": {}, "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.7.2" } }, "nbformat": 4, "nbformat_minor": 2 }