Skip to content
Snippets Groups Projects
test_fd_derivation.ipynb 7.47 KiB
Newer Older
{
 "cells": [
  {
   "cell_type": "code",
   "metadata": {},
   "outputs": [],
   "source": [
Martin Bauer's avatar
Martin Bauer committed
    "from pystencils.session import *\n",
    "from pystencils.fd.derivation import *"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 2D standard stencils\n"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "outputs": [],
   "source": [
    "stencil = [(-1, 0), (1, 0), (0, -1), (0, 1), (0, 0)]\n",
    "standard_2d_00 = FiniteDifferenceStencilDerivation((0,0), stencil)\n",
    "f = ps.fields(\"f: [2D]\")\n",
    "standard_2d_00_res = standard_2d_00.get_stencil()\n",
    "res = standard_2d_00_res.apply(f.center)\n",
    "expected = f[-1, 0] - 2 * f[0, 0] + f[1, 0]\n",
    "assert res == expected"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Finite difference stencil of accuracy 2, isotropic error: False"
      ]
     },
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "assert standard_2d_00_res.accuracy == 2\n",
    "assert not standard_2d_00_res.is_isotropic\n",
    "standard_2d_00_res"
   ]
  },
  {
   "cell_type": "code",
   "metadata": {},
   "outputs": [
    {
     "ename": "AttributeError",
     "evalue": "MutableDenseMatrix has no attribute tomatrix.",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mAttributeError\u001b[0m                            Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-4-ea41cd8e50a0>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mstandard_2d_00\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mget_stencil\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mas_matrix\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m",
      "\u001b[0;32m~/pystencils/pystencils/pystencils/fd/derivation.py\u001b[0m in \u001b[0;36mas_matrix\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m    185\u001b[0m                 \u001b[0;32mfor\u001b[0m \u001b[0mdirection\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mweight\u001b[0m \u001b[0;32min\u001b[0m \u001b[0mzip\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mstencil\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mweights\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    186\u001b[0m                     \u001b[0mresult\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0mmax_offset\u001b[0m \u001b[0;34m-\u001b[0m \u001b[0mdirection\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m1\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmax_offset\u001b[0m \u001b[0;34m+\u001b[0m \u001b[0mdirection\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m0\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mweight\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 187\u001b[0;31m                 \u001b[0mresult\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtomatrix\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mtomatrix\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m    188\u001b[0m                 \u001b[0mprint\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m\"test\"\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m    189\u001b[0m             \u001b[0;32mif\u001b[0m \u001b[0mdim\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;36m3\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;32m~/anaconda3/envs/pystencils/lib/python3.7/site-packages/sympy/matrices/matrices.py\u001b[0m in \u001b[0;36m__getattr__\u001b[0;34m(self, attr)\u001b[0m\n\u001b[1;32m   2139\u001b[0m         \u001b[0;32melse\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   2140\u001b[0m             raise AttributeError(\n\u001b[0;32m-> 2141\u001b[0;31m                 \"%s has no attribute %s.\" % (self.__class__.__name__, attr))\n\u001b[0m\u001b[1;32m   2142\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m   2143\u001b[0m     \u001b[0;32mdef\u001b[0m \u001b[0m__len__\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mAttributeError\u001b[0m: MutableDenseMatrix has no attribute tomatrix."
     ]
    }
   ],
   "source": [
    "standard_2d_00.get_stencil().as_matrix()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 2D isotropic stencils\n",
    "\n",
    "## second x-derivative"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
Martin Bauer's avatar
Martin Bauer committed
    "stencil = [(i, j) for i in (-1, 0, 1) for j in (-1, 0, 1)]\n",
    "isotropic_2d_00 = FiniteDifferenceStencilDerivation((0,0), stencil)\n",
    "isotropic_2d_00_res = isotropic_2d_00.get_stencil(isotropic=True)\n",
    "assert isotropic_2d_00_res.is_isotropic\n",
    "assert isotropic_2d_00_res.accuracy == 2\n",
    "isotropic_2d_00_res"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "source": [
    "isotropic_2d_00_res.as_matrix()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "source": [
    "plt.figure(figsize=(2,2))\n",
    "isotropic_2d_00_res.visualize()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expected_result = sp.Matrix([[1, -2, 1], [10, -20, 10], [1, -2, 1]]) / 12\n",
    "assert expected_result == isotropic_2d_00_res.as_matrix()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(isotropic_2d_00_res.as_matrix())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(expected_result)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Isotropic laplacian"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "source": [
    "isotropic_2d_11 = FiniteDifferenceStencilDerivation((1,1), stencil)\n",
    "isotropic_2d_11_res = isotropic_2d_11.get_stencil(isotropic=True)\n",
    "iso_laplacian = isotropic_2d_00_res.as_matrix() + isotropic_2d_11_res.as_matrix()\n",
    "iso_laplacian"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expected_result = sp.Matrix([[1, 4, 1], [4, -20, 4], [1, 4, 1]]) / 6\n",
    "assert iso_laplacian == expected_result"
   ]
  }
 ],
 "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",
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}