{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Vector operators in SymPy\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sympy import *"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "IPython console for SymPy 1.4 (Python 3.6.9-64-bit) (ground types: gmpy)\n",
      "\n",
      "These commands were executed:\n",
      ">>> from __future__ import division\n",
      ">>> from sympy import *\n",
      ">>> x, y, z, t = symbols('x y z t')\n",
      ">>> k, m, n = symbols('k m n', integer=True)\n",
      ">>> f, g, h = symbols('f g h', cls=Function)\n",
      ">>> init_printing()\n",
      "\n",
      "Documentation can be found at https://docs.sympy.org/1.4/\n",
      "\n"
     ]
    }
   ],
   "source": [
    "init_session()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def scale_coeff(r_vec, coords):\n",
    "    \"\"\"\n",
    "    Compyte scale coefficients for the vector\n",
    "    tranform given by r_vec.\n",
    "    \n",
    "    Parameters\n",
    "    =======\n",
    "    r_vec : Matrix (3, 1)\n",
    "        Transform vector (x, y, z) as a functoin of coordinates\n",
    "        u1, u2, u3.\n",
    "    coords : Tupl (3)\n",
    "        Coordinates for the new reference system.\n",
    "        \n",
    "    Returns\n",
    "    ========\n",
    "    h_vec : Tuple (3)\n",
    "        Scale coefficients.\n",
    "    \"\"\"\n",
    "    if type(r_vec) == list:\n",
    "        r_vec = Matrix(r_vec)\n",
    "    u1, u2, u3 = coords\n",
    "    h1 = simplify((r_vec.diff(u1)).norm())\n",
    "    h2 = simplify((r_vec.diff(u2)).norm())\n",
    "    h3 = simplify((r_vec.diff(u3)).norm())\n",
    "    return h1, h2, h3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As an example, the tranformation for Spherical coordinates."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "r, theta, phi = symbols(\"r theta phi\", positive=True)\n",
    "h_vec = scale_coeff([r*sin(theta)*cos(phi),\n",
    "                     r*sin(theta)*sin(phi),\n",
    "                     r*cos(theta)],\n",
    "                    [r, theta, phi])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAJAAAAAVCAYAAACpO2fTAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAFHklEQVRoBe2ZiVFbMRCGbSYFeJIOKMFAB6QDjgoCHcCkggzpgKQDoAPogOAOoAOODpz/E5JG0nt6l2XjYbwz8pNWe0m7Wh0ez+fz0QaaZ2A8Hk9F8aS5emum/Dy9Xce8FQ5ZTNsqJyFuUzczcKzf3XWbizZ/qX+CP1OfWtx2y3gqY67j8wGkTgReapX9aRG86e4wA5rPqcqrylkH8t4kbf5S/76EXqtcqbyofemUWB9fWBkO3fqt4/MBJG6UneakSNmBChMyydFs8JUZYK6+VbBlEFl/yUdsufQfsu2q3KhOJgp998PS6NMLIj4TQBLMtvVPip5CUShUuVYheklpoQEh6aaezIDmcqYyVjlPuhZu5vwVCCZ4zqX7LcBRJSsZsH13ktUrQ1b4hEDgo8o29VxRP4ognuRoPiteY75Q2V+X8cmWrL/URzJ4TW0VDt+dhHi1SQgVWmgE2TGHfFuKQM4+BE6UfcBtYPgMMK8qnIPYTopBB3+R8Tj3eAhsePFIVeRzMtST+n1mCvtz9ZDvi4gOVO5yxKvC20F8lz4Cmn2WWw9t4FlG/36vDv8tqcM6hYln4X1VYTUfy84d9TEGtn36OX8cqowsz19V6cfJjkbV0Z7KfYdxZv1lx4dstrAQXICkWxo0+B77+saA4xvdijlKbRpEZSsTzVK3MMm/QK8Aex5UDmybem2arbOzCSc5g3TAp+K3MNUJlttUl3CPIU5tbL8OcdShU2GcZ65PdRxPZ9tRIusv8RI4yOAbllfwTlf4FZ4t7yHEURdEY67pN3xkIAyPUhvcqwS7Mu+tTuzhAMrKBUjJC9tXWAfZcVcyOQ+Gq9pflTFckLObrLUrXpdhzRFC8uAhWzQ9pTT5i+2SuTMZD2GA5BIRM9Oo/mAjC6IvGD4CiPQbTkJfQSXoeeWd4RAJY4L8c4LwfVNrzp5iOrBJtjKBPGtgH1nhRvg+22zuzNnmzCZ/MXeRDbKPoAJ+vX8qv9gBX18wfFt9uZZBr4l3AWz26oJB481dgo4dCcdZOJR0/yhnpRlI6CwQgMVAul3guUzuZLMY3VuQwxX7EkBDU1gxIwJBpPRZ0F5GdWEdcpZZsQpK3lp2VNh7cBSPdWYRLMNwK7PNX24xOhM4q+SyDzQsgFw2pD8Hho8AgpnGOsCRjCi1ZeXGU0LHVMJxjAcFEecWzm0muHxH+Uqtv6SfwImCR8FMZmTrjra1xCTsjfiS/lzT8BFArHjS8WAgfarMVbh1DALxYtBEhfNEBVaho6K0GfETmxIS2ukCSGkSlt7NJn9dSRrBzcGZ+STI/UEdfA3g+9TmGrIKyvORcqPrpyLWX+PVx77OldBcBfU1V1N9o6u/2lxNKYNeqsXXaAc2WflL1RGO3dWlN7rSWlvB4SCeN/h6GtVxInOG0RTqBFIdHhxvOywcaJnnytU/sCU7T+JDB/7CJr6tvhAN/pw6+e4rnB+Pw4Vfx2cCRQ2c0vj+EDLn6pLBRCwsJycf/Cp0pPqls3EyU/plt2VPKX+RpWqTR9OY1ef5ttQAmKASf/rtafKGHMiMER1/VqGjoykfRlbKXxz8+9wc3YA9nwkgOZ0DIP/dEFmDwJ4Hngcxd2RahY6OpnwoWUF/8bredMCujNP6wPO5DAShj6oKVzfEUV9juomNqFahI1K4xo1F/cW5LHqx7jjWiM8HkJzP1nOqCIuupx2FcpZqen7vKqaRbhU6MgaQWXl/WRtYxF/Wx6dWRm5MlTHX8f0HCIvZbsCYQvEAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$\\displaystyle \\left( 1, \\  r, \\  r \\left|{\\sin{\\left(\\theta \\right)}}\\right|\\right)$"
      ],
      "text/plain": [
       "(1, r, r⋅│sin(θ)│)"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "h_vec"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And, since $\\theta$ is in $[0, \\pi]$, the last component is $r \\sin(\\theta)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def grad(u, coords=(x, y, z), h_vec=(1, 1, 1)):\n",
    "    \"\"\"\n",
    "    Compute the gradient of a scalara function phi.\n",
    "\n",
    "    Parameters\n",
    "    ==========\n",
    "    u : SymPy expression\n",
    "        Scalar function to compute the gradient from.\n",
    "    coords : Tuple (3), optional\n",
    "        Coordinates for the new reference system. This is an optional\n",
    "        parameters, and it takes a cartesian (x, y, z), as default.\n",
    "    h_vec : Tuple (3), optional\n",
    "        Scale coefficients for the new coordinate system. It takes\n",
    "        (1, 1, 1), as default.\n",
    "        \n",
    "    Returns\n",
    "    ========\n",
    "    gradient: Matrix (3, 1)\n",
    "        Column vector with the components of the gradient.\n",
    "    \"\"\"\n",
    "    return Matrix(3, 1, lambda i, j: u.diff(coords[i])/h_vec[j])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's see an example."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAU0AAABMCAYAAADp2+XVAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAWoklEQVR4Ae2dS5LcNhKG1R06gKYdMXvLF5iwW3tHWLqBbJ9A482sregTKOwbyD6Bbd3A7hPocQNp74XkilnORvN/aIAGQYAEn2BVAxEokngkMn8kEwkQZJ19/Pjxzm0LZ2dnP0jme4r3Fd8pPhUOBx1rWBGBivuK4FbSmyFwLkX+XPFjIv62Bie2zYdr0B6iqbYxmM9lJL9TfKTzC8XroXo1fx4Ce8dd/H0/T8Ja+1QQkC68VYzaRGQ8U/xc8bXiT4pvFf3wToblhZ8w91zM0N4P1mDNJTe6PoCo0iO1j4d5x/KD/J+5NNJrWBaBveMu/hjEGUi/XlbySu3YEJAu/Fs8MxP1wwNdPJZ+nN31UvG+3njXa53+LMJfrUU8k66bllPcTcv9tEwytdhIBHyMd4W7dP8P3SyPuGF0jgNRwy1FINb/1pA+BpLzLXFRw8/V3i9iyt0wWzZv2lLbeJR/eA3j+RJe3RyGfyXHfQvicOFICepCI5J1skmlcc/BXDw+VQc8VdnQy9htv6BHyDaVwRxcptJes15puTEaPA36XEpzJzdSXhHvNKuOymIkaOdebp0tyokfpuvf57alssjxe275VDnRYL34fir/1NMl++a452CuMqxtZut1Xz+JDvfIX9DsKzc1T3RvpS6WkFttMjCpqz7O8jTNTS9CueE7FfxJjRbzMkNGNVrxUOgP8fRjmNdzjdzIMjc8EYFVHrTNZWzt+gVxz8GcNfzJnlsEO7zWTyLpSyTdVl0sKvf5lJ6T0jMaX4ysiyL+PrLOasUlA/zg9WYbQFvnleqYh0hzmLODB+tok5/aqu7jOfXn8D+1rsWwCO45mNu+fWf5nCqmqSdabxTPFJn2Lxosf7vRxUWF6yG2F7lHTc8lD1MCDA5PnLOmqSrn2tjF1Fz88KSUJ/hmacFeDy5PqBxTysWm1KKFF/KX42PsUXVZmF5l6jeWl5zy4rU47uJhEHOVYe09S7dz5F6jjPjblS6uIWOMZim51S42TyxNm56zLWPs08VLNXhQveJTc41UGHCm5S91jqeG4cET6PUeVY7BAtR6y1EmN1g88GowJicdJOMucM/EHIOEzs4K6AxyW9ln0fIrQ5fr26aLe5Hb33Lk90v0XExjbRmFx4YvVCFpbESX0R9D5owq14z0zR5RW+ZK6e8VCawTUaZ5Eq4y3JgYINq6UITOtypD+y5c64T01nqiyri2XbnwiHFt2gozLX/gQ2BPF4YY5YYneH2vNmJrp9Bkb2CStvJWCceOu4f5ZwLorY+v8sCe2dCnkb4dwpytd/egH6mrrJug/KS+2fa5V9BH9Njs/7R12HYHf78qujI6NXrz0peDxEgopotgIn62uldD0YvJHTJCxw8+PVcZwGqmgjofMz1nLTM63VE6ysMTRjaOuukyPJFmpsI6ct2ZjigNuoYnHeGv04bSuJkM3TlHaCuyh69DS+lGkVwe5RTh1/HG+WuX7x9t2WieXy52rrooUdMnsTKpNNU7BdzNEovFobXMAS6KEj/aX/RPEnPlDd4TKpOlb7Sj+FvIh9LQiUZ/yVegTzjpXQJSfhFdtPwVu1dLyW37Bp3R6cc7YzzNK1WYuqB9oQZTniYeX/jmEcqDUhIJlHmh9kMa8PNaox+eA21cMhKq3EHnLkzxjF1d/whPH/wE7xwP+Jl3jedDebeMAZ+pFweg6eTU6WbhqHFXP4PvS4sWr8OGukFaynsfwtzpDzqVCpfKyNG3lM7A76V0FT5NQL8lF+d4p053TF7wU0oXS+tMKbkD+PNGVTya1oMSXY/xNBlVO3vflIaxYGTt5EmBjEenPDfqN56oy+No6zuPg3ZINCO4jr0jtk9n6Fy0GGEfxsopvfWAS9djsDHyxei6NNFjOgTNMCIvfIXpXHe8G4/e0eMu+RrMLQYtj1tp6EErzZO/F3PV4+akflTnPDqD+iYa6GKnL2x6x9tVepJvr93NdVF8FdcZ8bC53B7mjad5LkZ6A56bCjxQ5ZSn1FvfZh50jI3ajNYElC8VUGACNFKBm4DA2iXrhrSFoeHF+6U8TZGLB2ET8gY/3CyLBNHnK0xfhFHE8WCfhen22qyhJRg4etwlo8Fc/YtXdk+x8cyU5vQh5WkmYOkkh/0aFpirbykvNGwn+9rh4lVYShd3rTMryu1BeXOaMz1HKXkKiGvuBzrjYNNzFq/9uu7cvbrIdDYV3LSLGyMVmNYY113gNUsISmN0eA6PSp97Aw1N6QxvamvsDXuhik7GlHxLp58S7gwO7If0DRw6e1BaaqAfwtzpWrJfNtC3vj4voYt70JkScnf6YdDTlOKxlvh1GFFKRTbXkhd7Kuw3hvI5RWzSVQ8aKLYbxZo8TqSYTI/Ip1yz9kMeQfncHAQMOsYKI9kE1cX74Mml81abvAknyMDN1gnwqejkMzxZvk1Z5X3v5Yf14Q35NgungrsFDPxC49a3nkm1IcxdP3OTpsLa+pZql/TNdXEnOrO53LFOOI8lZqY5I5FTnAX7qGFUOp5C56MDMjJMrzGYBL6K9I3SUFQ/UOZHdajzIq8ixgk+Xb5fd+w5vPhbl0x9tYeRxGg7A94y7so3N6hVOlMn+IHmEvwFZAcvjxp3TzpupEYXbX/QF33LI0OYD/WZaz5X3xr+XMWZx1K6WFpnSsnd6i7/e5qsmTkj1SrkX0gpMVQYL2ck8OTY6tOsKfnlObeKjBL/I2Y8lI9S/az4QZH1zU8UeTjUeBAqgyIz9T4oErimjDE4tg0MFvWhRzmm/a29nLqeFCx92mstJXi8MzAQwAGPl3I8kOGBVh82lHmiMoPYQ8sPavuxrnnYNeTp+9Wac4/3o8PdCeHJgK68VwR3g79wafTHleeoOr2YK591cHBtDYABDfQ/qW+iwT3Crgr6iMB98kQRvQ3T2XlBOq/0Qhfd5ZsIGKlOEG3KFNFFtc29VeReLSy3WepTn5jtDXQuT+xaT8eVaZ5eL3UU/eSTr6XaWJuOZMAgL/lEnhtl8j5S1eWGjD4hXhuLLemDkWIW7iqHwes8lXb8QkexF3PqK+4aV2RQzMLEyd53hBY0+8rsIa+U3GoXoykIpr1GqbqTAm8/JEfuSRS3r4SX3TxoWqB5PAtu8qkBLxsP5tRDFHd5Hr8pYuBMsF4Qyt3XR72YWxo4EsnZwU1rxX+jmMzgqheXGXSXrroLubfyNGmn9dbGHkausTxIBpYZZo/wosE0J+kRjeXr1MvHcFca3lbzlpYtk/QQlT+IOfQUO/sq94ivlffW6WIJua1eSA029DTVGGt2fArNrfHo8ijDUiMyD4+ia1ZHicr6TMdwJ413xNmdgMfOfta+9d0czPFSn60vziItxDCZQjgHlyl016pTVO7RD4LmoGCnPtdS7M5T6Dl0t64rOVj/4e2gSVM41cebYaE/+qBia3mOpb05uOdgbsvwdyh90/tdwTUHEwTJwWVXAltmtpbb4sTDN/OiK9Nm1oWynp7PBVCN42nyb5CMFjVUBHaBgL0JuSmOfd19F3ieGhO+0by7tXBSyhdiwGxD0nmJ/YlGZPHAgjJrXHiNeHxM7Q461rAiAjvGHe+yLpes2PenQvq8hCAYTrX7QTcQ+802D/bGxavgg8p4FheK15szcssa3DPu6IK6g2WTGioCvQgUMZq6eVgS4MtEpTzN8GHUM/HDF7bxOmtYD4G94/5GOsBDkRoqAkkEihhNccMbBaWnQr6BdNNyPy0JWs2YhYCP8a5wt4P4In+qNguhWnnXCGxuNDWSszXkFymou2E2B0ht84TU93LxfAmvbg7Dv3ilipOnc9SFxnBLp1OiNO45mItH1jafquy9Y0EePUK2KfzmYDKF7tp1SsuMwZj0GqXqUZe1waxXLlUWI0FbzQdkc+uuWU78sEk6uSk6bFtlkaPz1xphuaFr0WAqOHtz8lA7e82X7JvjnoO5yvBXGdl63Yev6HCP8Apxtn710QvzRHe2LorGUelhCZnVJoOS4J+/ud2ALWK5gcX2n9RwMS8zZFQjFk/R2TPZtyk6rIbcS2yZ4gMO0Lp1oSDuOZjzoHKS55boSLxWPkKzRlhCF3MwWYP3qTSLynw+lWspPaPxxcj6KCKvIe4iSAb4wevNNoC2Dt8Rnb0x3Q4evCUFlpOC6vItz8n1JzU6s5LFsAjuOZjbvl1kbVO0+EDymeLiG+YtjrN1UbzhxMzSw5kqkV19LzKPnp5LQqYEGBw2xWdNU1XOtbOLqbn4YbuT+W8hKY1OzfXgl55UjinlYlNq0cILmfxOvuryRHqVqR+4LB3Fa3HcxcMg5irD2nuWbi+NUS498beYLuZgksvXmuVKyax2sXcSbfr0nP2NY18hvFSjB9VjVCsaNFphwJmWv9Q5nhqGB0+g13tUOQYLkOstR5ncYPHAqymyZzWXzyXKScZd4J6JOQYJnZ0V0BnktrLPouVXhi7XS+liJiY+C5uf70Xmu2MlF+NY3CmfM+N986SxEV1GfwyZM6pcM9I3nz6zZa6U/l6RwDoRZZon4SrDjYkBoq0LReh8qzL+++7XNr21nmgVR1nJgHFt2gpLWf7Ah/BAEUOMcsMTvL5XG7G1U2iyBStJW3mrhGPH3cOcjw/zPcgGX+WBPbOhTyN9O4Q5H5jhYyDMjJxOKqkdlJ/UN9s+9wr6iB6bbXa2Dtvu4O9XRVdGp0Zvcv5zK6mL8Cw66OHSmPAmH7S3uk/VVCsUkbnFgb2g07OenqscgDVTQZ2PmZ6zlhmd7igd5eEJY/OXqTqHL9LMVNhed6YjSoeu4UlH+Ou0obRFPq4KbcXmU2S6mZrpq9KNMrk0yinCr+ON8+in4GzZaJ6jlzqqLorU9EmqXCxd9U4Bd/f3zeDQWuYAF0WJ/nc/uXOl0z9JzJU3eF+oTJa+0Y5i53NzSkMnGv2FNwX6hJPeJSDl9+niWpgU1ZcSMgf6osuPd8Z6mleqNHVB+0JCpzxNPL53ot14lbqmg1BKIoEy/MlbSAN+XmsExHOgjUtGQ5U76NyFKZ6xq+sf4Ym/hogFPGDeLHKBUZ7ybhkDPvFeYgGaTs5Y/lppR427+hl8X1pweB021A3SUt77EOZOf9CpVLhURo6+pXQGfi+lq/BpAvotuTjHO3W6Y/KCn6guroxJaX0pIXMA+83llzr8qfgvdVjjOYXnymckbz0o0fUYT5NRtbP3TWkYC0bWTp7jQXlu1G88UZfH0dZ3oyvtkGhGcB17R2yfztC5aCX/skN5GOoGP12PwcbI59cPz0WPKRE0w4i88BWmc93xbhxd5R097sjgyQMGLY9b1+hBK80r34u56nGDUj+qcx6dQX0TDXSx0xc2vePtKj3Jt9duVBdVdxVMoGv5Knafqv1NZXZYc1T4j+KfnJ/r5L+K/1TkPBrw3JTxQBVSnlK0XpB40PVFkMYlozUB5UsFFJgAjVTgJiCwdsm6Fm1haN6K/6U8TZGLB2ET8gY/3CyLBNHnK0x8vq8VRRwP9lmYbq/7XlU9etwlo8Fc/YtXdk+x8cyU5vQh5WmqeFYI+zWsNFffUl5o2E7W9YqY7FZfVpTZx/x/usBOZk/PUUqeAuKe+wHFPNj0nMVrv647f2VPmM6mgpt2cWOkAtMajOsHgdgsISiNtavn8Kj0uTfQ0JTO8Ka2xt6wF6roZDQ0Nvg5JdwZHNgP6Rs4dPagtNRAP4S507Vkv2ygb31qMKSLS2OyB33ZWuYo/knv0i8txWMt8eswopSKbK4lr3lq6df1zlE+p4hNsupBA8V2I1mTx4kUk+kR+ZRr1n7IIyifm4OAQcdYYSSboLp4Hy8Unbfa5E04QQZutk6AT0Unn+HJ8m3KKo+/ZHD5YX14Q77NwqngbgEDv9C49a1nUm0Ic9fP3KipsLa+pdolPamLttKimOxEXzaVOQX+eSojMz1lBGLVWbCPGkalMyp2PjogI8P0GoNJ+ErxG6WhqH6gzI/qVOdFXkWME3y6fL/u2HN48bcumfpqDyOJ0XYGvGXclW9uUKt4pk7wA80l+AvIDl4eNe6edNxMjS7a/qAv+pZHhjAf6jPXfK6+Nfy5ijOPUV30aK6BSWl9KSGzB+nN6d1OSkaClBJD5YzXQ11jMNjq06wpRcggcHTfm+oxtf5U+T/riDKzvvmJIovOdD4PWPjWIXl8geZAmgKKzVqfMzikwwPGFSXlmmk/D4kMHZ3PCdyEsfVRpi7Gm1W7bHNBufiKEWV5IAP/fdgwmDyh3JYBTMTjMePu4GI5Bt1BL98ruqUepxeunH8cwhxdc1NSv55/3qtv4od7hF0VZjDVNfcJ/YzehunsvCDdvdKLMeY5AroUCylddGUXx2QH+rK5zA7M8EjH8nio9WRcADVPgpc6VxvRp19L0d+CjmTAoC/5RJ4bZfI+UtV9rBh9QrwFHlu1AUaKWbirnBmsUrxBB3qpfNIVGOx2jSsyKG6GSR9eW+WVklntsuwnMae/Rqn6kwJvP7SmrpOolK2EN9M8aFqAFTwLbvKpAW8KL/fUQxR3vDdF480DgM6ZYaDgfX3Ui7mlgTPRNztQdvGwGSbFJf2bgV3IvKWnSVuttza2GqGWbEcyME3IGuH72hUNbvDOPr2+Orc5L4a70vC2mre0bJmkh6j8Qcyhp9jZV7lH7K28LV1U2uKY7En2EjKrzcbT3PR/z9UwngDrOr+oE47WO5IMTO9Yb53lNYsOxpePnyyx3ipSpx1iuCuN9UIGYwJrmfQL6+fRkIO5ymB02BGSpBMlXiBxK0wKiJZssoTMahOjiW4V+d9zRvprNc5C+9EG23EPJcekKZztBD5+XA3mCC2Yg3sO5rYMf4fSN70fwfH6RdfGZH0JxrewtcxWL8oYTeARAzy8eCTFdE8Kx6NWa1QEFkbA3oizZxALs1XJ7QAB32iel+BHxpKp+e9ixGzFKMFDbbMiEEEA7zK1xSdSvCbdRgQm7dNcAihrOJcgVWlUBBZBoM58FoHx5IkU8TRPHtUqYEWgInCyCFSjebJdWwWrCFQE1kCgGs01UK00KwIVgZNFoNia5skiWgWbjIAeDLIdjTc+DorslWTfJXt6d79fUnzWcEsQqEbzlnT0kYh5LT75oLJ58cEaUf7KhO1pdT/rkXTiqbNZp+en3sNHIp8MI29c8Dpg86aYzvE4ueYtshoqArtAoBrNXXRDZUIIsD8y9im3l0rnf8OZutdQESiOQDWaxbugMmARSL3ogLdJSOXf5NbfisBGCFSjuRHQtZk0Ap4X+SFdynwDsye7ZlUEtkGgGs1tcK6t9CNwYbOdV+mXdoaUL/nXUBEojkA1msW7oDLgIVDXLT0w6uk+EahGc5/9ctu4ct5kTG7nhfLfPzVUBIojUI1m8S6oDNitRQDhDKQPivM+6z5NH5V6XgyBajSLQV8bDhBgu5EzkH6WM6Sx7Uh+uXpeEdgEgWo0N4G5NpKBABvY+VvdMPCF/zeeNxrm1+uKwKYIVKO5Kdy1sRQCMor8bcgHbT/iq/4m2K1I3+hi8/+EtyzUQ0Wgg0B997wDSU0oiABe5ZWM5X3LwwMdv5JBrR/sKNgptek2AtVotvGoVwURsFPwo/lDs4JQ1aYLIlCn5wXBr01XBCoCx4dANZrH12eV44pARaAgAtVoFgS/Nl0RqAgcHwLVaB5fn1WOKwIVgYII+A+C+EJ2yMoLLc7X/4EOUanXFYGKwMkiIDvIX63cTwmI0eT1tNQTy7rVI4VcTa8IVAROFQH+pyr2dpqR9//3dslsFoKX4gAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}4 \\left(\\cos^{2}{\\left(x \\right)} + \\cos^{2}{\\left(y \\right)}\\right) \\sin{\\left(x \\right)} \\cos{\\left(x \\right)}\\\\4 \\left(\\cos^{2}{\\left(x \\right)} + \\cos^{2}{\\left(y \\right)}\\right) \\sin{\\left(y \\right)} \\cos{\\left(y \\right)}\\\\0\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "⎡  ⎛   2         2   ⎞              ⎤\n",
       "⎢4⋅⎝cos (x) + cos (y)⎠⋅sin(x)⋅cos(x)⎥\n",
       "⎢                                   ⎥\n",
       "⎢  ⎛   2         2   ⎞              ⎥\n",
       "⎢4⋅⎝cos (x) + cos (y)⎠⋅sin(y)⋅cos(y)⎥\n",
       "⎢                                   ⎥\n",
       "⎣                 0                 ⎦"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "grad(-(cos(x)**2 + cos(y)**2)**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "def grad_vec(A, coords=(x, y, z), h_vec=(1, 1, 1)):\n",
    "    \"\"\"\n",
    "    Gradient of a vector function A.\n",
    "    \n",
    "    Parameters\n",
    "    ==========\n",
    "    A : Matrix (3, 1), list\n",
    "        Vector function to compute the gradient from.\n",
    "    coords : Tuple (3), optional\n",
    "        Coordinates for the new reference system. This is an optional parameter\n",
    "        it takes (x, y, z) as default.\n",
    "    h_vec : Tuple (3), optional\n",
    "        Scale coefficients for the new coordinate system. It takes\n",
    "        (1, 1, 1), as default.\n",
    "        \n",
    "    Returns\n",
    "    =======\n",
    "    gradient: Matrix (3, 3)\n",
    "        Matrix with the components of the gradient. The position (i, j) has as components\n",
    "        diff(A[i], coords[j].\n",
    "    \"\"\" \n",
    "    return Matrix(3, 3, lambda i, j: (S(1)/h_vec[j])*A[i].diff(coords[j]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's see a simple example."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAHwAAABLCAYAAACyYeOrAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAGrElEQVR4Ae2dT5LUNhTGpykOQJGq7NO5QGoI+yyGG0DlBrDJejjD5AYkR2BuAEfIcAOyZ0ExlWU2k+/nlpxnj4wlY89Y7acqtfVfT9+nv2297t3Jycmp7JVsylze3Ny8SEV42DoR2O12HyXZPiWduNw9NBF/yE1ia/62HndXgcCFpHzUk/Sp/M8J28nGEf5EPeADgW6OCwGN+pdq0RtG+IPjapq3ZgwBJ3wMoSOLd8KPjNCx5jjhYwgdWbwTfmSEjjXHHssG02qXxzafnd6Psh+12/s9JlbcXm7O8T8o/NqEn8n9Rvax7AvFvVdaTgQ/y80RcBFjZKV8jiOvZZGRur+T/Wzll/9EeVYvq2lXNge0LWUA4kb2VECcpKziLgiX4Sz3xaaR/1xWQf/nlZ8O8i7EEQ/x2Hc23dxulU+9jayULUMn5buF8+DHfWXrlb8KWSVnEQe9NoKDghpMmp4/SLgSMjqeh8SQ1gcMYjtEyv+oVyFldNLY+LncqqP50iGWF/y0rZFHTzrsPsbzjHExjHjZVckaZCriILYntLGI8JY8VfxFthktsUD5Qa0TFuNCZcwgdIq2HBs/p7tfh/wsNdnkKe0qZbXtknsKB/mER0JUEetcO1oMmYQll4OQ520s466fQd7BzmjlqUHWIGMRB4GnlvCsTZsqwvAS5YMKaDdm8tMJrhV26ytZbTKYPn9V3L28fFH9jFbM+8Nj+LMiWYs4SLX4QSpwIIy1rf8y5ZnCbgEaAHzWJ1vhbwfKniWYemXZhGHojOxU2s6ouHMTTzQ7dDpmFbJKzmwOaFvKlBAO2RHMeJQBVNbn1ghARhabpyu5L2RfygI0O2Q2fYsYlY8sdKiGaD3pjK1RPGAxG7UzlMKqkTU0JIuDttEJx8NE2FAQ59k/BRJkfpblPIjpj/BXArWJU1rWDtIjKOH9tAqezfylki5l96qXoyDTH52tOVnIzWjvn/9rkpUm5HJA2qSZ/Ho0AMmXKE+SJXvg4gjkcqB0DLzm9WjWCFcGpkrOrw258jO1U0hn2pTfzUIIzMVB7hrOWmfXXzrA64Wn6IWgq7bYWTjIGuGC6JXsqXoZayPrM2S3u1/53SyPwCwcZBEeRvKSG67l4aq8hrk4yJ3SK4fLxY8IOOERiY08nfCNEB2b6YRHJDbydMI3QnRsphMekdjI0wnfCNGxmU54RGIjTyd8I0THZjrhEYmNPJ3wjRAdm5n1Xbp5HZp9CV55zlQJb9gey7oiQkTcPA2uhD6VHVSaMGmzOTBVdZy8dhu8eaov7UlcdAle6Xlf7ooIwg78UjZglK00UcqBrVN521urcrsiQg+cvTDJvstu85a4VUe20oTSIpMrIvQBFigdZQf5XRHhMEN3RvgvAuaT7E99AK1f8azJRZfgQx5XRBiY1nv4gu1XlSamcEAdMr/JfsLNpu0f2e9lx3bsRZfgtclwRQSBmmOEFfsozNglkyIODkU2n//qE45HSW5Shw/WEVdEsIh8g5sBEXbelMLsOao0oSTZHFBeyoyNapsn6xJ86K1sSFwRwaJn3MIIgouUJkL2LA5MVbecTOm5JvcSfE2X+2uSFZ5yORjk1BURBqFZf4RmCr7YGlUGUTpXRFg/nV0JRdosyiC5azi7SFdE6HJw175ZOMhdw2e5BH/XCB1ZfbNwkEW4DuycD8fOiEeG77qaMxcHuVP6ulrv0kxGwAmfDF2dGZ3wOnmbLLUTPhm6OjM64XXyNllqJ3wydHVmdMLr5G2y1E74ZOjqzOiE18nbZKmd8MnQ1ZnRCa+Tt8lSZ32Xrldz8XfZsi/BK8+ZpHJFhK9QY3AllSsicMuyxAo0Omb25X7KDnlWrzQR26Unl0NH/5XC4qb0nWvKvGel5UO/eb5XXNEleKXv3w+njFVd7o+A1CCrZCzmILYvdOoiwlvyVPGUX+OnQ/k/IvRmq0RHG1SasGkncpBPeOwpqog1mZnAdoCx2YE8rojQIztiap8B26UUEVrCszZtEgZTdAleGxJXRDjgNvoprBg4mLFLJkUcHIrsfpYcy1hHXBGhi99kHwMi7NIpg5nQFREAIdcIvGb5UPoGPD07P+2teDrstabR61imwhhZVShNBJldESGSp6f/I4IBY8jpighDyFQQrhnKFREq4GmSiCLXFREmIVdvJvYb36wMknssm+USfL1Yr0LyWTjIIlw7W1dEuGfO5+Kg5Bx+z0326udAwI5wFPj7ZV6qZ/HtjptKEBCH/AMk3zkkDYRzmEfRPGX8n4tSqKw7rPk5sCER/wNcaFn0hbXDZwAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}y z & x z & x y\\\\y z & x z & x y\\\\y z & x z & x y\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "⎡y⋅z  x⋅z  x⋅y⎤\n",
       "⎢             ⎥\n",
       "⎢y⋅z  x⋅z  x⋅y⎥\n",
       "⎢             ⎥\n",
       "⎣y⋅z  x⋅z  x⋅y⎦"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "grad_vec([x*y*z, x*y*z, x*y*z])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "def div(A, coords=(x, y, z), h_vec=(1, 1, 1)):\n",
    "    \"\"\"\n",
    "    Divergence of the vector function A.\n",
    "    \n",
    "    Parameters\n",
    "    ==========\n",
    "    A : Matrix, list\n",
    "        Scalar function to compute the divergence from.\n",
    "    coords : Tuple (3), optional\n",
    "        Coordinates for the new reference system. This is an optional parameter\n",
    "        it takes (x, y, z) as default.\n",
    "    h_vec : Tuple (3), optional\n",
    "        Scale coefficients for the new coordinate system. It takes\n",
    "        (1, 1, 1), as default.\n",
    "        \n",
    "    Returns\n",
    "    ========\n",
    "    divergence: SymPy expression\n",
    "        Divergence of A.\n",
    "    \"\"\"  \n",
    "    h = h_vec[0]*h_vec[1]*h_vec[2]\n",
    "    aux = simplify((S(1)/h)*sum(diff(A[k]*h/h_vec[k], coords[k])\n",
    "                                for k in range(3)))\n",
    "    return aux "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's see a simple example."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAHEAAAASCAYAAABsHjEkAAAACXBIWXMAAA7EAAAOxAGVKw4bAAADlUlEQVRYCdWYjVEbMRBGbcYFMKSCQAeEdAAdhEkHpAMoIQMdUAN0AC1AB9BBEjpw3hOSojuffZfjzj87s5a0klafdlcr+abz+XwynU4PJ5PJFSxZl67oe36vbvfvruIfEvc+LrrTmYlpX8Nz+DTJtrUE407iHxI3uiZ3Kqw7CdkfuS4fo8063+DLPrqZt1H8fbEPiXsPZafwE1ynRwT7HPnjeseWtXcV/2C4deIrfIizTEu7SLuKfzDcM1LYlyWeSyfQxQJFR1/E5ldKH0M+hBz7Cf6FvhvKtVFX/AX2I8C9lDjpcw9mo8/I39YBvivuhAWMntxb+AA+Z/4jMu1+kh8zCHOdDif4sLlOcuqe1LKtM1/gcJfF+lMa/z8lc3vfiU3roK8Jf8Ae16rc9cguYVT9s0HXetTX6z6vr4GuBdyOgbT9AyxOWWfKD6G/rihO0jlNL9b8AKI/vWCDjLaOOGzS1yaLcwcxRBN+cblG7HPzlWCjrYGCQdqw1vvVCw+CHT0Ldo+Ys91j2/1kvNMoRPZOHFE3aVj+iKJQIFdRTjW0TT+/kZ2V41bVmaPjjbY6GWmmiZy6iwGvrHFetFdWm/CX2Kn76v6Jzpz2kRnu/i/OsvoijBkVO/ob7d6AwxQqFlPquz+olGk0HNVStqyOEjc+VAQOEs3igW9XYE7pKkc24zWKezleNm+VnHkfxo6OlbjT+owTfyVD2jdDGIhIEMwRwnwCkXlsdXLlhCB345J/Q7aCOuL3RD+zn5xRaGuYN2Qb+TrVEbdf1fTPd3AuZKU9OhygU84YkB2oHHJiIJXApj3JjevcvHH6Lot+u9dGrNuKP4IxKCsBSdvrYCPB2BU34/SD/qk4ELkfOiYzKm7Mxj11c62U7ihTzA3ycIyRq+QertyDUYfRXEY4w8anLvgLFDowZBdlcV/urR68do9KXXEzzgDVL9fULX386B8xZ9w+ULwTmtj/U3SFSTra3C2rJORxSv9qXDiuLzPfSOt1vzKvFX/CxVhxuw+NkfC7wV6v6mibXthZsxNuxuU7nrq29mHm3Pxde+F1SufaiQjTEBpy6etwDFCs64vwhHWXffBoXXZT2Etg2+JET8gBxqzfVyXWD9UxtifQQAkOo+2aRrV3Te87MeoZFXvbxmdtA9bRjxHfWEcek9LdktbQqf437O1AFa0Je8LcWG7FSWxENrCQE+MDRkdKfj/1rsmv6yDd0Z+/pA96ugqJwccAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$\\displaystyle 2 x + 2 y + 2 z$"
      ],
      "text/plain": [
       "2⋅x + 2⋅y + 2⋅z"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "div([x**2 + y*z, y**2+x*z, z**2 + x*y])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "def curl(A, coords=(x, y, z), h_vec=(1, 1, 1)):\n",
    "    \"\"\"\n",
    "    Curl of a function vector A.\n",
    "    \n",
    "    Parameters\n",
    "    ==========\n",
    "    A : Matrix, List\n",
    "        Vector function to compute the curl from.\n",
    "    coords : Tuple (3), optional\n",
    "        Coordinates for the new reference system. This is an optional parameter\n",
    "        it takes (x, y, z) as default.\n",
    "    h_vec : Tuple (3), optional\n",
    "        Scale coefficients for the new coordinate system. It takes\n",
    "        (1, 1, 1), as default.\n",
    "\n",
    "    Returns\n",
    "    =======\n",
    "    curl : Matrix (3, 1)\n",
    "        Column vector with the curl of A.\n",
    "    \"\"\"  \n",
    "    perm = lambda i, j, k: (i - j)*(j - k)*(k - i)/S(2)\n",
    "    h = h_vec[0]*h_vec[1]*h_vec[2]\n",
    "    aux = [(S(1)/h)*sum(perm(i, j, k)*h_vec[i]*diff(A[k]*h_vec[k], coords[j]) for j in range(3) for k in range(3))\n",
    "           for i in range(3)]\n",
    "    return Matrix(aux)  "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's see an example."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAADYAAABLCAYAAAA7+XTCAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAEQ0lEQVRoBe1b0W0UMRTcQ4hPdAoSBVw6OKACSAeBDhI6SEQFKOkgUAHKlRAqgKQDKAApUcQnHxwz5r2Vz/FufOd9ziayJWftZ6894/HavpfdpmmaOeKyI54ul8tmjBF4f3RgBtxl8xiFGj4hwcp++OlnRpY+Ap5pgOkV8rtqU8XmY1RmHUwgtI94QzElOuh1MplwVDm614icEduIX9D7Ba5mwZ+KVp18RcMfQWTBDoToOa47sJlN9UdWbIQEp8ZMSdGGNJUjyVPmrYIpMYB+i3gWAf8NtrmoFynON1kTe9MBkaoxdJX/L834a0bMU+OqB9+spyyryIwYUG0JMlXHB6pkn/nGIdOWxBTnVBMlr5bEVJUYH1XzMlY4hM2MmCzrxKgkfLyq4v3cx8CES72S8Ikp2dhW4NfbOG2mmCDiJvwygu4FbBeeqpEqeSZTYgDOXwxXWPrbE7dsA+9g38uD3n93ibMi1fkAQjOBwp8Wr0H6fh+CZbod9o/v8KWmU3F4uOktVmLpYzWOmlWxceiQjqIqlj5W46hZFRuHDukoHqxi5mdFOfRWh2n6ZOuvaToVoVZ1mPaP//qlpooBTpdDVF1yXeXrMwnuMCMmiwa76/NWzQI8g2XNiAGhOmxUHR+0kq0OU39UUtKWiqkqMRyqpqnD9Cl6/oX4N4ZgU5vnWlMSflPqaxzaYfoEnZBLQ8V+Iz6XNC6DhtIO0z9ATy6O2KBMgsaqwzQYkOys+SEYCKvDNFsmrwHL5d7rpnyyEis/5nk9VsXyxq/83VWx8mOe12NVLG/8yt/9YBUrcVZ0csk/1/V/0TPR8NDqn+xFiKk3GCT4/qILsNE7rG+aDv4iS6mp+BkkVt7rAEmqR0ePyZumpYjRf3iOGAb3CxvqzcOC3HwpYvRtzGRK5mJOur/IM4Zpxx+bsaBKtU4dIU+fPwPf4uGU5WLDuvRDXqK9Y1xvDbyBn4AU/aAA/XF6st8jAHWfmSBN75WfJ0G+q3/AOpI+1/rhFeWsD7Or237bUpoYAa98O4O8+6xDAUueKKcCli+b8XV3NxDhFWUtsaSpKNODHwaoPxDJ3sDVbg8dX8Rqob0T2M9Q/j4o54cHvFcDVWU9Z8PVfZSghX3XJGLScNdz0tf+jTKQOqARbYakaPNJsRofE93UmU8OpVZFBwikOJW2fVKwcbXk4rASYNOFZaPNuxgxAcrvWUKl/Jc0d1FvKgw5DaliO51RduCVrwxEmEmaiuFN6+YBhorwhLFAmgsEAwlsIXLROoadRFiHxy4+SzuIbZA2rlE3nK5tHT9RhBg6JGCSc8+XDwBp3cO+I01CnJqsR4L7SHOhcacWkOKruEmhCDEAunXhESXaQ7KgT9qIY0yLPWOxzi1tlZjl6Fq0XRWzGFXLNv1VkT/Tw74WWK3ClSqscyd5YOUhmltINJAY95Gu81i760fvvluj+yXQBeEfI0cVjvDHMm4AAAAASUVORK5CYII=\n",
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}0\\\\0\\\\- 2 x\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "⎡ 0  ⎤\n",
       "⎢    ⎥\n",
       "⎢ 0  ⎥\n",
       "⎢    ⎥\n",
       "⎣-2⋅x⎦"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "curl([0, -x**2, 0])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "def lap(u, coords=(x, y, z), h_vec=(1, 1, 1)):\n",
    "    \"\"\"\n",
    "    Laplacian of the scalar function phi.\n",
    "    \n",
    "    Parameters\n",
    "    ==========\n",
    "    u : SymPy expression\n",
    "        Scalar function to compute the gradient from.\n",
    "    coords : Tuple (3), optional\n",
    "        Coordinates for the new reference system. This is an optional\n",
    "        parameters, and it takes a cartesian (x, y, z), as default.\n",
    "    h_vec : Tuple (3), optional\n",
    "        Scale coefficients for the new coordinate system. It takes\n",
    "        (1, 1, 1), as default.\n",
    "        \n",
    "    Returns\n",
    "    =======\n",
    "    laplacian: Sympy expression\n",
    "        Laplacian of phi.\n",
    "    \"\"\"\n",
    "    h = h_vec[0]*h_vec[1]*h_vec[2]\n",
    "    return sum([1/h*diff(h/h_vec[k]**2*u.diff(coords[k]), coords[k]) for k in range(3)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Hagamos un ejemplo simple"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAABoAAAAOCAYAAAAxDQxDAAAACXBIWXMAAA7EAAAOxAGVKw4bAAABjklEQVQ4EZWUjU3DMBCFnaoDVHSDskFVNsgIRUwAKyAmQB2hYgIoG7QrkA3IBiXdIHzP+MLJdag46eLz87sfO2dXfd8HL1VVbdL8yDhHX+E0nlOy8ZuBy/eEfqLX6K+vEqVkIn6gtcO2wmz+15h818ZhrnhKuBAWk0QjhD3g1ogO23usZOP3gHb5Gph2GAu13dQAOsOYPXe4NMdPRe5yHtgaVdzZhI/kET1BbOPs/x8VWhL9L0k9/RmDiA0/dMl4h6oR9DNV5YFxVFITaP1rlMRJTR1RvBWBtbsorHXoPdi7YYXxKmFWvadY8rkdnRa1mzfPSvOXrJiMMkzVZaMyoVqrpHW2OajdFWBlQGG0qgtLwXZ7tB2pCSxhyWFRAoW54iyop9ouW0ukmz8ajLVL3aiGsaA+kSU/2D3S/4n9ToX+EutlOLuIniMbGbuww8vig+5w2FgQbFXYofmzoshnzxKYnpucK/+lYlapIuYh0F16MuwItO1n1nWsg8BRwAb8dgAxUmc+YeoOSm7Qwf8b4lRkuuLqfiwAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$\\displaystyle 6.0$"
      ],
      "text/plain": [
       "6.00000000000000"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lap(x**2 +  y**2 + z**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "def lap_vec(A, coords=(x, y, z), h_vec=(1, 1, 1)):\n",
    "    \"\"\"\n",
    "    Laplacian of a vector function A.\n",
    "    \n",
    "    Parameters\n",
    "    ==========\n",
    "    A : Matrix, List\n",
    "        Vector function to compute the curl from.\n",
    "    coords : Tuple (3), optional\n",
    "        Coordinates for the new reference system. This is an optional parameter\n",
    "        it takes (x, y, z) as default.\n",
    "    h_vec : Tuple (3), optional\n",
    "        Scale coefficients for the new coordinate system. It takes\n",
    "        (1, 1, 1), as default.\n",
    "        \n",
    "    Devuelve\n",
    "    ========\n",
    "    laplacian : Matrix (3, 1)\n",
    "        Column vector with the components of the Laplacian.\n",
    "    \"\"\"  \n",
    "    return grad(div(A, coords=coords, h_vec=h_vec), coords=coords, h_vec=h_vec) -\\\n",
    "           curl(curl(A, coords=coords, h_vec=h_vec), coords=coords, h_vec=h_vec)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's check with an example."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAABoAAABLCAYAAABumJecAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACT0lEQVRYCe1ZwVHDMBBUGIYnk+FBAaaDMHQQOoAWoANSA+kAWkhKSAvQARTAB4YnD8xuIjmSLFknP/LhNKNEJ513rdVJvrGNMWaG2mbqqm1bI6m4/i2Dgctbc4xBV57RoLNf3n2j0H7E+DTyuYJ94/rcjGaSO6/xAcEdam9Gjjj5P5lMGgws7CDbLAugvO6a5d/ijABBSYL1gk2puLbz3Cwx1s0I7S4YstLBZ4U6jQHR98ka9zsbYx3REQxJmcPpJeG4Qd8UslKVwSIlYvQ1AIyjahDcH/TD2+8P2pDiMujYG24mxW0gndEe2rYwO8rJ6FviRr56DlHHaCLgPKGuQeJCPoIOzVFEmA1JNiC5DeHyVjURSB4IB5L7PGx/pIoIJDy3LnwS9DEamz502CMmAhgj7NonsVDdoRlCh5YovO0d83RYo82jh4V76gyVJ8qSHUNFRAQAklCe7fpEgMU9RH8REe44t2EjzrwpXqM8hGxEiWQ6JbxUuoQosi6VTqZTwkt0BPE6e7C6p2ljsTSBTIgads1hagIZaDJ6wyIKKSej7x8nkKeY/gfqL2qxQLKaBPIEgMQ2XKNv1HPbxl++gKQ2gfwBGrG3RHlkbwQkmkB6emgCGYhRZ4w+gupoKsK7Fjj21xnFiohtlU4sVeyoCWTwahry6BvIOEZCe+CFBh/vLMW3J6M3rCaQ1dJVJpC7FcRvFdGIBLKeSBPITjPbEJ3eAxs2xsvaVcGQRREMKJFApLTLwaQjPY/6g3zw5bPEvbUisV/Enz5xUeqDb4f1B7ViiLvlVvnoAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}2\\\\2\\\\2\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "⎡2⎤\n",
       "⎢ ⎥\n",
       "⎢2⎥\n",
       "⎢ ⎥\n",
       "⎣2⎦"
      ]
     },
     "execution_count": 17,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "lap_vec([x**2, y**2, z**2])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Plate equation in polar coordinates\n",
    "\n",
    "The equation of a plate in polar coordinates is given by\n",
    "$D\\nabla^4 w(\\mathbf{x})  = -q(\\mathbf{x})$, with $\\nabla^4 u = \\nabla^2 \\nabla^2 u$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "w = symbols(\"w\", cls=Function)\n",
    "rho, phi, z = symbols(\"rho phi z\", positive=True)\n",
    "D, q = symbols(\"D q\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "coords = (rho, phi, z)\n",
    "h_vec = [1, rho, 1]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAVsAAAA6CAYAAAAeEPqmAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAN9klEQVR4Ae2dj5XVNhbGZ+ZQAMx2AB2QbAfQQSAVwHRATirI2e2AUEHCdhCoYAkdhK0gQAez308jGVnPz5b851l+vjpHY0vWn6vvk+67vpY9l7e3txenDpeXl4/V54++X85vJMenU8th/RkChkAZArZ2y/CKS1+upGz/UL9PEUTkPdHhtdKPYsHs3BAwBOpDQOvV1u5IWq5G1pta7Tpq4IPOH4rE+1GenRoChkCdCNjaHcnLvZH1BqtJeT5UIVwEX2W1vosrKP1dlMay/aS8r1GenW4AgT6ONyC+iXgEgT5ebe0eAS0jexHLVmS9Ut8PRcx/dPys9FvFA8vV5/2sMs6lkCGvFakEgVyOKxHXxMhEIJdXW7uZgEbFZvfZigSs2bdStI0PVnlYrzwEexb69mS9UfqF8s2qDcBs4JjL8QaGYiJGCOTyams3Aq3gdAnLll0FN4kM+GV/CHmerH8p/RN5SnNuYTsIDHK8naGYpBECg7za2o3QKjyd3bKlfxGCvxZrlgCBxL8UH8mK/aTrnFMmBHy2jSUcMpc+BjnV969L95W2T99gEecr76XS79L8uEwt5wE7L88Bx6eSM8hhHM6DeMDzGK+6XsXaRb4g62a4l6AXc0ZhgL8WN4JrV+coVSxXNvTeD/lrH71cbGOZdfy57YFJV1nlv0W2rmu15Em+KjgGJ0XjcKY5XAuvOfN8i9zPuhtBvzS4CvDNNlaqzrFksXLZlVCTbxal1viQdV5LeCFB3ivGOzZqkQ1roiaOjcOZZkZlvOaManPcX+WMqqAMFiw7ELrCyW/Vu4QgTxOL2/UP/BAcK7NkvvpnZ8bfXX34H6R3KoP1ODqweKa2caTzKjjW2M6ewyP4L5VdBa85g9sq97P5bAUACuSL4lMpjGZfrfK51cPP86AWy1YyIQ9yrqVssfQ/q/+POh4Ej+X/dP3BwcXMDLWBBYo74t+ZVQaLebmq4HgPHA4SMlOBmnjNGdJWuZ/NstWiDi6Czwlg/GL+FF1PLp82KaJQ/vhpV1G0frQo+k5F62UDy+B+8VXWP0QcrsqxcTjvXKiF15xRbZn7TmWrAT1RfK34SpEXEtg7mxPY8vWzyt/38bXS7DSY1bpS28jl5IuFUt5LxT9Dns65lf6iiNUdAhZfY3mHzK6j6o3F4aA5tfXYt4eyBx9eUUY+p/wPKtzJWKNPeRLHGi9jn8IfUBmHHRNmYtYkXnP63j33UoStp/ECDUuUD8OE3QTc8n4J6aGjyqI88KcRZ3+qrjbdU3wdWXAtuZT+QzHeCYHs7IJ4EuT2ZV6G9LGjyk3CIW5XbfGj8wN5Cg6fcF1pfhwOcFIe+P0ZypUeVRd8XpXWyymvdkdzrLqT+EM+BXg+ew5zuJizjDAdzWuOHGp/19ynipYF/lcMnCegpbDi66c897IEpYUCaykjpZGztQiVZmE2ykzn+GtdG8dk1/XZcFBbLcXg247lQak3PxBBJuWhLFtchGs5R19/EWWb039XGcnEYp7EH+0q7ILDLgy3mmfc317cY+ZGwflXozSnLBBCfCvOE3036+8uLftXE+zS98BDpbDb4bnyfgk9S57g6khdBLgxYv/steoE/3Konh6zcUgrxmnJFKz7eCcGL3bEaXANssfVkTlgH+e3ztUHsmLBp4F2r3U9fDc4vg4mvS4K1ZuNX/U1J3+M46w4jIlZ+nxOXnNkNe6/odQoW5HAgmWB/v7tsjsLiqD1QCcCMSlenlTf9BGUAuedHxNXn05JRrIGxUunyM9e3lixku/qcJITorazcBhoE4WWKv+0yvfKaGGbFuhLa7zulee0jMaBZYwFPcpfrnpBQaZNj06rzcX5Q7itcTga0BEVl+BVeA+uX+P+omXZohiweFLlhBJMrcMDmgU4VhigH3xS8aDwYQa+nPhj4tx6Ny9GHBZ3LyOkMlG/pdi8TP9N6vMk/X6SFycn4RA3pHMUKZuvXfDycAscp8GssdD9JQ7XiukPR3T5tKcT+U2FDRjH48vlj7aMwxTRkelSXlUeN118Z0bPJet3t9xfRRxhGabKCkVABKCjQQSwAR8rCkvz6CcVjzZwp1jC5Q+0pTb7FCKKPbUGyYsXL+1hIcfWL3mUQZEdC6Nx6GiQvuJxpG3jd/61Q0aaYjzpDx/5Jw8z8JvKPIU/2jIOU0RHpEfy2mUExetpaP3ulvsrOPKKDRCILvg8rDIUVqrYfClXF2VMGaeofVmUyJum0MCJ6sSvpqKQuizsuBX30CtkSFbqQDiyuKA8zimXBsYS99dcL8GBsoq3is1Ws6ahbyfc4rOlJgT3IR4SqodsjDO+HspxRMbWj1988VTnknMyvx2yTuGP5ozDDlBLsubktXD97pZ79waZgMe/h2LlV4tzAudsATuqaCmkulhu36tcoxh8Htuyivx+vt57tfdMdVMrle6aoLIodALWH6++cmvzRhG3wT8UeZKf3u4Ef17n/zxTm0U4qHxwCXynvjqtUJVxykryoJS5Q2CyDWKrepTnW7+9+KtMZ/BjGe2zDY2qndn4DW1yVLuj+PN1+XE9ew4Z61JhDK++zvOudYWc/vrg+t0t9wIOnJj4rW1U5OdG1cUi5sk7kYVAmobjLU4sWp6cuy9G6fg4bl9prqPwF/8ymPpASTayBTmUV4yD6qCgD9oKbYYjZRRb29LCtfToy47e9kV7Csg1y9YvL08vv77PXo7TcU5JS6az53AKPjl1C3iFe9Yu65Mj68TW7xEdKWxY6wfrN7y4wPvuoxYm9RSbfaI6pyMIYcU3ilPnWGqOIB0p07yQoHMWKQSS7xZszmQZW0Z9MHmaFzdCO8orxkF13Ebt0Maxo+9zUClTnzYVR/ER+ld9cMzqL9TpOiKH4iC/Xu6jHHe1PSVPMp09h1PwGaqby6vKcSfmjAQdw0sJzC1+7Nz65qho69crX2HRuX6V7xQcK7z1SzVEFtcVsJ4ONfidYo2VabpgIadRxrTh0+QRD9rMkaekjPpwvqNQR2kUPX1n46CyjCNLKapcrlKmzdF3GWE8cxwlRxa/9KXQy/Ec8qRtqE/j0C/wFJu+dC6vHZw2c9hjH15QsfX7TdEeXb/CzN32j1JuqgvIDQGBYOVh4TT5vlzz1pbSuBpUPN9VMXdZ9Y9ybT487WUqwkF1stwCyJ5bVuVaCmTucZe0J1my+PXjo+xJOVZ/xuGINZTLqy8X7kZbhkV8rWROzVV2i9xn+2VTkDRYwEeLNN8d8IsuWIjxLQblYpcCtxzNrWna9qnSkgFZsxXm0nIhCzIt3U9O+5Iji1/aisqenGPwUjQOM5VuxFXvuvW84lZz81FHDKSgeHnoW2SY5My50jLIprgZ7kcrW08GK6112600TvTm1lrn3IqqePNhG0eU8pqFGa7ZcT1Lvwt7eFPs5Zd6CsZxprLrwvnUeQW8YhQ5Zaaju1PVEQXHnWsVRsGpsZvS36SPh2sLB5YYb/68UCTgGP4qgZpXSP02j2vlsyWLwLanar5v6ySyP50I5PBLReO4E77ZMoUvP2Yot1GvXqeC5PJKPZVlTf+tyDrHxcW2yl8ky1cdLRQgMEnZ0o/I4JeOWwzCwX+G1XV8eSjX9E0uV8H+1I3AEL9Ibxwvy6HwnVXZes561206IsnA3eosyj5tey/pycq2DygRhN8Pv081/xKnT167Vo6AcVyOWWmNJZRtiQzqH9fftZRt8+JSSX0re4fA1cJA8Ov50W45FkZ53eaN43XxP0XvrTdET9HhOfaxqGV7joDZmAyBUyOwtmV76vGea39LW7bnipuNyxAwBAyBIgT4UIzbu1NUywobAoZAMQJyp/V+mEkWLE/+w8PmuH2efbCjp+vjTHw5rvcTqDSktm2dA8SKwdwIK4JvXRsCOQiYGyEHpfrLmBuhfo5MQkPAEDgDBO6dwRhsCIaAIbAjBLy7BdcKO2FwrWziJSlTtjuapDZUQ2DrCHhFy+dRnf9aaT4P8F6x87+v1DRecyPUxIbJYggYAkMI8DZdHPhnqY+ldLFyqw6mbKumx4QzBAyBDgRixRq+0RDndVRZP8vcCOtzYBIYAkMI8JosW792H+Q+4ENWceBVYgL/1bfqYFu/qqbHhDMEDIE+BOQ+4ENX+HCr/0iOKds+Ju2aIWAIVIuAf1jGd7FvqhUyEsx8thEYdmoIGALbQECKlm9pb0bRgqr5bLcxt1pSaqI9UQavaHIL9U9FPub8UUcLhsDmERia3/76o2DR+vTn2teAuRE2NjU1sXh/vvlF9xON/+f2YGNDMXENgQMEhua3rvNA7I0iW75CwI3wTGsg7EwI+VUdTdlWRUe/MJpo3DrxtkzzRFZ5bHnBwn2qfPu4cz+EdrViBHLmt8rwzwh4e6wVNPd7P/LTKrxSwtwIKwE/slus2ub/u/k2wv7Cgwk4sg+rZgishcDg/JZS3ewdnD0gW2taFfarX3T8tCjU35OqYZ+h+WwTYCy5HQT2ML9N2W5nPvJAjG+Xpn6pH31+17dOtzM6k3TvCJz9/DY3wnamOJZtyycrawCrllj9Rzi2A7NJuhICZz+/zbJdaWaVdCulivsA3yzRBZ/HF49uZO2aC8HjYoftIbCX+W2W7TbmJr/6hBtNzFd3pxfsSGC7iylaD4gdNovALua3bf3awPyUgn0tMfl30uYu2ABfJmIZAnuZ3+ZGKJsXa5V+ro5/W6tz69cQWBiBXcxvU7YLz6KpzetXHz8tPtvWw7Gp7Vp9Q6AGBPY0v03Z1jDj+mVA2bLly3yz/TjZ1W0isJv5/X8dycEPokhl5gAAAABJRU5ErkJggg==\n",
      "text/latex": [
       "$\\displaystyle \\frac{\\rho \\frac{\\partial^{2}}{\\partial \\rho^{2}} w{\\left(\\rho,\\phi \\right)} + \\frac{\\partial}{\\partial \\rho} w{\\left(\\rho,\\phi \\right)}}{\\rho} + \\frac{\\frac{\\partial^{2}}{\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{2}}$"
      ],
      "text/plain": [
       "    2                            2         \n",
       "   ∂             ∂              ∂          \n",
       "ρ⋅───(w(ρ, φ)) + ──(w(ρ, φ))   ───(w(ρ, φ))\n",
       "    2            ∂ρ              2         \n",
       "  ∂ρ                           ∂φ          \n",
       "──────────────────────────── + ────────────\n",
       "             ρ                       2     \n",
       "                                    ρ      "
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "L1 = lap(w(rho, phi), coords=coords, h_vec=h_vec)\n",
    "L1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAABuoAAACWCAYAAAA4wd63AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAgAElEQVR4Ae2di7XUNteGc1gUQEgFH3QApAKgA/hSQaADWFTAIh2QVEBIB5AK/kAH8FUQQgf539dIczQ+vsge2yPbj9by+KbL1rO35MseyRf//vvvd3sKFxcX71Tfx6r31z3Vu5S6iv8dyXJPyw0tt7W8li4+ak3YOQFsoxwDCLr4KUjkNvtU7fRzORIiyTkIcP08B/VtlUnfsi19UpsyCXA/VaZe5paK/nVuwuvLn76gHJ3RPsvRxVSSoNOpSJLP1ghw7dmaRvPqQ5+Yxykn1vWcSFuJI8N5orrYSUQ4H4FXKvqVXvq/lz4eafs3LXfPJw4ldxGQjqyvlws5trGNLmUse85t9KGLlA080Mp/cLBjffNB9fWfB55uvqIDKyguXD8HMltL9KX7+b32LWuxhzXIST/dqyXup3oRbTLCbu/dNqnNaSpFXzANxylyoX1OQbGsPNBpWfpAmnIIcO0pRxdLSkKfOBHtaxPlU3w2eqj3y+a/tHwpXthtC+jRjO9DFX/UOm5vu9YrrJ1fhEnsN9LXUqNPsY1y7ORmIor7zVuyB4+C3UN4rbraMUkIBLh+btcUztDP77lv2a4hLV8z+ulu5txPdfPZ6ln6161qdny96AvGs5s6Je1zaqLnzw+dnl8HSFAmAa49ZeplbqnoEycivLijTi+FXoWXfhNVoT8blecXzLfkcGCKxX5cs8aw00f6uGM7UEF3tP981gLJfBQB6cdOug9LthlsY5SqZkkkXaSjXP0nh8/WzyyFFZZpsHm/BHYb2H0QB66fG7WCM/Xzu+1bNmpGZ6kW/XQ3dl+v1b651+7GtLmz0vvu+1ffs2h5uznljqwQfcFIcDMko33OAPXMWaLTMyuA4oslwLWnWNXMKhh94nR4LwRzutx6ctKNs50zN1Tm5NOKKe9byvuOlq/K/2iUls4907FfLJ62P2l1V/u7eOnsOi8dunSRyqJ4dgB4eGz6YJlGYfsMBKQXt6PDsOWpRcixD2xjaurN+fXpQuftpPlTi/8Vtatv1KnuHlXnaTD/aKa3raNttqDjXD+3peqqNtLrWft5lb/bvmVKc2prt1OWUXJeqv+u+um6LnL0rzjca9fBrXg/U+e77l/3ZvM5NmGT3xuXczTzHF0ozq7b5zn0ckqZ6PQUesPS5rAeliOx5yaQozPF4T50bkUsnH+f3nWe69yJOllsRJ2UZSedR7XN4aRz438QXqh+Vllvg3H4ptQvo35w+V60fVPLi3Bcm4QpCYhrqy5cjs4/i+VJX3ao+h+/drISyiHwm0SZZaRjl31gG8saQJcuLInO+wJrW7ivtrorJ53rr+A28KLa2vhPmy3oONfP7er+LP28cdK3TGNUbe12mtxXk8tu+um6Rrr0r3Pca9eBbWC/S+exeoqz93u378LzpWdG+BC5bHXdZxM6T1+wkPL7dGExFGf37XMhdUxSDDqdBGNWJjmsszIi0mIEunSmc1x7FtPEsgV16d2S6DzXuQlUsoijTsp6JFmfaPl5ApmbsrCjx99R8s25Xyi/0VK9YNX+Ry3P4+I4Ci99/NvmvL+q+yJOKJXjF6olhFZdhEZ7ePEdZPYIyEWcADvUxWB7ECO31S8zto9G+8A2BqtqigSNunDGQR/+Y0PlsNW+txcJKquIPjO2AcnjPx9sPTTaghloOdv1cwz0UuxnjOxLpRGjs/Tzrp/K9s37rvuWCfXc2G4nzL81q1La2c766bo+GvUf2tjZ7rXrQtb3l7Idl6uySnk2qmMYu9+o85hZ0P2m+1fXNUev6ht+VdQvirvY/WvUw8LrVpsI9nC2vkDlW7bZQ449zC7EtwJadeHTQR+bbp8F6WIqlaPT5a6jnaynUmg9n6X6KZe7l/YR+rqzXXvqOm7aX0rvG9S5cba21aD7s1znLNiW9HrNFZozBGX5n9ueOu3rTGX9rnwfqKxHXkIZRw9nlkOLjaZ6SaRtG9isweWpzpM7oZSvnZ71cEvHJ3uhbNm1PNPifyQeRijWC23Yb9VF0P/zmLfSugO/35DH5Idc5lp1MQWMwNy6fKfFa7eDpvCTDs75bYdG+8A2mlTRf2yAXpsya9RFiOh/Ibuf8VTB/2iJ/ao25wuuT2Ht9KVqu/WXPFZoly1ULxesG8Vb7PppoYaGAu1naBVOij+gPzhLPx8qR99ykpaPEne226OYE+4U2M6K76fF7N2EKohZNer/nPdTUbC29Vy24/KU9+zPRm31GnrcHMLi56wmuduybNR5EnlT/avr1cIn95n3sbJ4ojwWuYdN9LDkZqtNnLMvsH0Xdj+frRPLriXnmbmeZ6suQsRNtU8xauq7cttmnV2p++hUL+Wl68neMXYouo91R9Jxp9zW5+inLM2e28c5rz05ljCX3tekczPQMuY6Z8RdbfUs1zkLtTm9qiF5FNpsi5j54dTf+ZmlDOVth5sNwt++q8oI+7OVGcvpWksGOzs81efk9Va+vqhcyVfHX09RZsjnSSwj7L+L+21rl60FXUg3gdks+m/j33ZcsvgF+0EWbds2PzTF13ELf2hLTXHGHlO+xdmHZFpzO83Wa11n6OJwrejsM90WtMzWJup6Ocd+ibYwhoPqsdq2PKa+9TSqf3Z/MKdNl2hPW7SNc3EukaVkKrqflnx2Ekx6HVF+xd1P1fuk+r5knq2PdlkKsz4b1eszdl9y+jnJn2uwzNbjPzl5hbi7esYKjE7Sq7i5/WUxztFDSXH2aBOq80n20Kc/52+uMZ62W5+ZYxyvnUbLrtqn6jurLlK+59hGp5fvGsWi83n5VP2cg7XK5J5EfddY3Z1DZ2NlTdPNqXflvYo+0XJaf5FLYNL4bjjGiWun08K17tt1f9Z+8ZpAzxbk1fTNsf+BUU2fNlNBBuSpLL/W8nfne5agens032fJNPloup4KmbMbXmeQfDlsbieZVA+VyX7bJrq4JJOli8vos27V/0nqf59f+Tag7MJt1VOR1tvSVMIVZR8baKdZem1RHrr4BqaznYa24H78XgvHLRwuyhbGAC29LY+p04g0Wf0B/fwIsuOSdPYtMUvpI+d+LEZP14u321LbWcn9tJjZifhjqriJthfX/ylyn9F2LHZWWzylfrlpxcHfbPHzob/R7fBFi+/Jc0JxOl+DXsX6D8Hd6hSY2MRly8lq57LZvmtu1r3UZbGHLXRxQFFOn3sp0qgtdHqJbar2dZnj8dairNdw7TrGU+TeojqbgsAZ9Z7VfqaoY2YeY69zzr44vW9Vr7M66qRIO41+DQ/RVuwcwV7dw/fmpCjv+98B8SFojjL78nyhCG/6Io05r/r5od/T0V0JgbOHp1cMrkS4POA8WoPyearFHUoMd7WRwxNdBGIDdBEZz71ObSI64tJjLv+hlupbjzMJU5p9rLqdBh2lOmzTa5M60YWoZLZTX1/cNrYaSrOFMZxLb8tj6jQmTU5/QD8/hmxDmiXuxxqKjYfO0W5Lbmel9tOeGizXERN1m7M+h/5z5GqLM5vtuMCJ2mKb7FMeP+LgexAtv2QWUKLOj+qTWY/saBPq1c+0nmY0vUZmy1FwxF3ZxET20PkOJOg6tZPcZyt0EeBlPlsV3KwOoqHT4TrNaV8HwMnG0qzXcu1KEBW3ubTOpgAwm94nuj5NUcfcPKy/GHKvc45fot43qddrUTtTr2Wsfji1IntHeMWybeBa4rfR/F20o+/MxXi19WvtV15hxfeIIJc3yXfPnF+Q5/DvK+3bEeYRZtWDobb/0VK/KN3RTcrBeei4MTiulqF1rOY5djrl4+/93da255Rt4mMH4WTzSKsM1+2/WlLHnXYbA7o4xjKpLo6zzt+TLd7Wkjpao93UnXLWdaPd1kuzXWgZasez2Ifk2GU7HaDXuvq8jy4uqfS1U4+oi23mMtV2tmaxhTF49tqWx7CqpxnQH6yyn3d9R9rHZu7Hajof3W5HcnTxJbMsrp8WZ1835ppZY7T+a3Y0eHek/bTajgVQnmPuKc/2bDQY2mWC6nlR9bXTKH4HvDp2GaV1a1adb1mvuj56VJ3bYs6zbKsC5jwxsg3MZhMj7cGIWtv6yDou2s4H3EvVzQFdHBPpe7Y6jj3z3kjbQ6fHeplTp6NZj+yrWvspV3mkvSzaVx2r5ix7o3V2qrQjde5iW/W+J52fcJ0zw9n0jl6NNwlS1Oi5abvSqgiP+ur9rlmah+LbAeYGpFXl5MuaU95ptNgxWM35n+Z5yrbyq+aZ1doCVd/t0trTlhzqpe1qWshYjvbtnPwU9+vrEH9QHZXG3/mrvhmndZTJD3ZmfPQ9Me3bafm2Xm66r/MH+dPj9W3Fc11dv2f1c237iosuQpsSi15dtHGc87jkst1c0amO2c6uHG+SRfGKaauSJbaJXbbTqB9xaNRrPF9fK/7kbXWLulCdqn6wzm9L+3PYwhg+W7SfMRymSCOWjf2Bjq+ynzeTofah+Ju5H2uyCdVvVB8+lGNgXzRL1am4fjrh7Pv1w/1Jky7HHFOeo/Q/pqw0TVKvQ510rPXZSOc6bcd5Kwy+p1SaSZ+N0jrOsR04uLKHb6hr23+sbH1mrMuhuLPpXHkPupdW/FXpVfL6PcHBZutsz70v2Qa3AcusMItNKN9B9hBk6bSJMXVUmknbufMbomvFb7yXaspDcdFFoe9ApBvaV9BNartuD1oWf8eYyhC3x7YfpRvUVyl+Zz9leRQG24vSTMoycil5rTrP0uf11VnlDtJ50Gmn3vesc9U9+zoXWM6id/R67Je7LiCTB3lDrTw3hqe5mSuNH/I8b38c0fNF+/6HpR1RX7vyCWliuq6o2edUruV/qbUfoixXlMHfmnBHHIP/XZL+Y9bp0v0Yz/+yGFzHkMZTo/x6yEgblkfnXI7l8z8FY7CcfilwclAZno7lF5Xj0Xu+se2d/k1xrAd08Y3+ZLo4WZkhA+nxlTbfB93Ws72pA9HO6+cO+8qj144Vp7rLOSRq2FCchqP5h1QHZUE7NTFxaNWrzqELMZK9jO0zi2vH+a0kL+Yc/XZeyZexaMuXLE7d6uoPlPfq+nnzGGkfm7kf6+vHdb7XbDZ+zSyqn5Y+7BTwv05nC+fot6duh4ajPHvvKesQQ5qzPBvVZRmw7+dah8MzpHT43ky1+I+m6cwX32LWfufSuWVQUUOfeVv7V4sddDTouT6kmUuvv0sst0m3zdzpRhV1/jCGVZRqDpsYaQ8WqdUmxtRxZnuICFvXKr/12aopEbqY531UE+shx8bYXswfnS6n0zGspds9XLuiORa3HqOzUysxUucudlPXp1M5xvRDr3NON4fe0WvUyOV6FkedsreDzjfavQ8dl6JUadIpKe75nPI4chxIib0vnZM8R22qzAstlbNN5T1WJqkjzI6x9OH7hxg3KexI5uS4ufTWMYnvTaexDH7ouaHV394OoanD8YOgX8RVQWk8bafjpeEwfWdy8I3q0fbg4hvVT8rrkeIcWGgfXVwC7NXFZdTzbElffji9IR22OdBtXzmh145VRv9bw56SJK8d/u4HLNdtLf4X8pETWPu7b6d9ep1QFz9JBw7Wi79jefSHhLgveVbdZ1Y1PP456lOPT51vT5xn73+XqF20zw3bzxIYD2X09QeKWGI/39m3uHIn2MfXA5zjjd7r2HH0au9s92OxnTTINOjQCRxdTjEsGypdTD+tNug25nuto2tkg8y+r19FPx7tL9ZJcg+9zrfZjrEs3habdDH3MbHzn4VcTBOLo+e0pexiT3oN/P0M4etN2/Ou9XOOkNUGlrAL24QWnq16npnRxbD3UedoVEmZtK8ERrJ5lnvaKdrOnq5dib5m2ZxCH7MIVst0Ap07x6b7Lx/P6iMcMQkntZ8kn7NtSvet74aXsgv02q7+uRx1/1WR2U46GYIfav2QkqapPwRWtYjKrHZO/AnGaSfUf5RvW8O1Y84NMTrKLOtf3g8hdZz5kF8WOM5RGFLHo4TH/36306JipPz8sjx6tNMkdtIdXg6oXpWTL42gtK2j44Kc/1P8+0rrB5o0OO9DmEIXlkUZeviy/1XqqSJfaPvuoZDjjVXr4rgqy+6Jrdn5W3XRlr3/Rfupjt0GrthuKmmwj6y2mqYbue22mdrGb9rfpG1ID4PaaeSZqdcY/ZS19VCNqA1lut3aedoUttZOj/rUpgqf45j0cbIz3HJLn4frQdCtHeJtuh1V1YFlbM1+RjEbkyjob3X9/IC+xViG2Meq78fabEB69jXY94BD/xCXZjmEo9MVxTKtSNguqZ822+pb0jU5/V3fD7L3wwwZU/XjtXJG7c7YTzfajoVUmb7nHHNPaX37ntVh8LPRt2Rn+fU9t+tbD4fnNp+Y0i7Q6xHqN9pzO7wlxkfMj2ItuDOkDUxhF7YHM1Bem3/uVh1ne7aaShfKJ3226rr/3to1e5FWdo72tRedjm1fadtxXyxDOPV+1rY0pH1wT5K0vlQfyeGzbPr6tJf2cxbASaFiXT2riHfju+Ep7QK9JuAHbF4bEDcrqhThztYPXr4RzA02lMPouZCHj/2cm8HQeME47XBzZ90V/CAY49TltBf6MMIsZOIb/6aHsHpaczqqo2TySLd3WswvBk/VUaXV+qEazUfH07adFk1TUfrc6IcP5R/r63UMdqB537JMHR6rzOig/VGZx+2mcqJsPlfnuTldNAEYc0z2Yluz0+v/tP3Ii7Y9srNuJ7bzH7R0hTr3K3bclXjgOWyjA9gAvXbkkn3qZhLT/ab7qhvJsXRza+3U9YzXgLSeW9keotuxdR5SxtbsZyyzQekG9Ael9fNDbMNMhtjHqu/HmgxAevY12NPk+d7zs/bfdvTFTVnEY0M4Ok3pLIvpp60bLR51Xi1i5/svh+c6dnDSfTtU1O+QtjjEftpsx5XvvaeUfZ/92WhCLfne++A8UN18P/5RdtH17HNq8ej1kmB8Zq/s7vLwWbd628DE0vFs1QFUbfKOTuc8M3fkkn1qrrZpAdr63V5721if21vfbG3lRUSneZz8Rx3rZor7WZfIPUkm98Kj0X4WUNDC1znXCL2O0Ov1EWn6klQXREUa8tBhh9MfMhp/p8DB/+a/qwcXd7qtwTcSOnlHy6B/FSudH+j9bzr/m6w1/3DCzsKniuePLPqmx9ue+rLaVx4+dgiW2Xlq8bQ3qfw5dbRc/meoRyRWD/PKw+X5339mcyes7UzxiLc0fx2qgkccebrLU4IfIl2uRwu6LHPu1IfijtJF4GUdeiqSO9pvcj7qVBX2qItY91PWfyqxbevILhrsx7bseF0hx46P0mMbs7XTXL0e9HGCLtyvxOA+Pv1uZzwe11trp74eHfXzsaKlrMfq1fKrH+jV7ZD8FdcfJT96GZ1TRsJya/aTVG3Wzdz+oLR+vtf+atSy7cPXONnj2u/HatWv7sf8Zwm3XTvqPELkhRY7AIaEbI7OdAUsi+yn3R8KX3TM/GZ9iWV0FhzpS+dG3UenmeTmYbkkxyL9dIftWPSHWvqe/0p5NkpRj9oWi/hNOj9H+hnOdbvflVmuTp0Heu0ieegzzd12d2T/3SlnPZvTBq4IMMQu0sShPfLcnUI53s69lzqkOkEXQ+5/tnbNPvCbeWPp9oVO8xXqe54p7mddYnb74J6kXUFj+7I0x9w8FG+x+1DL16H3nD5iM/ehQVeDrnO5OnXe6DUQnmIlo7XhTrZIJo+k+2dInor/SYu/f5Yth+L7ZbEbuFbViws7Iewc681D8Z7FeNp22VnpYpq+tfPXclSfUM7RsbZ86mljPOcbt9vWiuNpdToZKI6HFXfGGXJe+Y3WRVpOyKdX/jRN37byLFoXffIveV6s/M/eT11l+rzjdcVJzykuttHQ1sSl184Vp7h2Kpl8o/JBi//oMGUfUmw7DfXt7Xun5DEkL8k3VRtr1O3Q/BW/lZXONZYxpL5NcV2mlqN+SfvZfVU9bSzD+cbttrXi9LbltrTnOO66aimun5dMu7QN1XtQPx84VfaubevSy1Ee2vd5L3YG+A9Qk/TVIc8i25lk83Wpt71OxWLqfCT7yf34kDy6WOnc5G3R5Wk5sh0zVKCf7mif4jPILhS/tQ3oHHr9ZnN+T9F5DZy6fXflZ1m0XGkbPWkG2UVbXirX+Ux6D6P8rrR1Hcuuo+I2snC+bfWIxxWnty6Kc3S9jGnHrpXfybpQHpO3TddHoWhdjGU+JJ0YZNtezFdp0GnDdUlcJm1fys92X9mo1o33s4kdO95k97Sx3KjzuNbxbHuxzDFdunbe6X7TtuL0smxKN+cxyTSF3Wfn0cVJ5+gTG9rgnPpvy3uoXaDXaZ65rY9JHt5TxUo5/2jJvglSXDdE300McpYpvv+lengBoW138J5zvbNOTuN4yWJ5vX/Iqy+PnPPK78BA29l1VFzX48pLcB2z3B4e3lq/kLZyXvbEO8jWFS/3nModpQvnr3B0MdO+D16pf64sTfGU36G+2i5KF03ynutYZNNWfjzvdVuc+nHFxTZqbVZM3MZX1U6tVwW3new/RNRtoW9feRfZTiXX5H1SH4sh5yXf6DYWy1Eerbodmr/jx3zTdVcZabyx28q/SPsZW5+50gU9KPvme4l43uu2OPXjinuSDYYyd9e3mKPCwW7rXOv7iut/H9shddBN2Pd3bap7Q+37RUZ1r6i17xuz8495dK3T/LRdzP2UZCm6n+5i6nMKJ7WhoXm4vCaZdLz1WtAUf8gx5X1ki6EsV/5gz235Kc7sz0ZtZZ/zuOo9yC4cv0leM9YySx+rfFelV8nrZ35har4GLnlccmT3oalcSjfILmppd/3cLXZH9pqyGbN9ii5cnsJsbTPkf6hvKGs3fe6Q+qa6V7rR7WvvOhW7g72lTOvbitd7PxtYznZPW5dV+9n9seJu7p5EdTrJ7oO+svNweXW7CHnQJxZwfxJ1M9Qu0Ot095bXBHOyoKGOblhehkwT5o7a8/N/HShI9e02lRm/ueXkfilxCDr3LCyvta7OqRyX5e9EVEuI/NLHDwmn2XCZ7qwchtTxpmRp4ndPx1unE1VZ5u5v2OVM5VF9NLKSbJqfUboIMr+IIgQdfW2pf4w2Zl2yLsbUZ5Y04u426KnC/G+YpjDEjmN6bCOS0HqN7dTiB7n9cqOaXk373p46FNdOVU9fNzzNZ1OfPHX9x+bX28acsepy5XoYjvva0aXb3vzNyTbhRXndDtvuL6qg/b4yYtRT1sXZzymVmSttSf2864htfDfkfswvLHy/6mt1GvwCPgY76eK9ottgPW6MN3ZdXDuTDa2hn+7j3dvPOgPVtbEfD5l35mFOWs7ZT6e2Y5GH3FMu8WwUMBa16tSpJUWvg/XlERPm1vasMzjDExIMaQNpMTl2caWvUJ19L7b35+4h19yUedv2KF04s6CPrvvvtjKHHE/73SH2toU+d0h9U6boNKER7HTqd4w597OWYs572rRtuKwh9rKF9uE6p6HX7h1Z9nDl2pJk0pmH0p77PtSipnrfu84T1bVudurUqdBrK7vTTuiBvnWE1tBzksQPy/6nztG/tYbm0xdf+btRnfSvYqWPLw4try8Wk47ich0UGv9t0Ve/MedVloeE9/4rdUzeXWlU5km6UHo7M32Tavn9sulOV3ljzynfzetiLJs0XdRDemzstvLCNmr9a+C7xnbqFxvu1OIyy7RByr+odip5fG3oHak9to2cmk6y9bYxlxHq0TjKRudadZuTv+L4ZdfbUE717zjt+/p6+FeltlvLOJVBml7lFGU/qWwlbYtTdb2dQibl1WuDiuN2NNj+ppAv5qHyV28bqoPb0eFe1dtaDtPnaNttMW13nuJt8nsq5VkUS8lTdD8dbbBtLfl725DThnq2taPOPJS2iH5acixmO4HZWZ6N2nQ95LhYdeo01A+91u6z+xiL6yLvKvrkGHs+0y66rrk8dw+0mTZdTaAL7o0n0kWbjoYeR6dX3wuLyeTXUeXZeT9rvSnMfk+rMrgn+ca6934j6KTr2tKZR9An7wsK6/O6+kjprFOnSTtFrzPodTInXVCUOzv3rI1z9nYZwpBzyt8vH47K0L4dd9UDbJDl8CLZcbVUBjSkHOJevVjXmYgrupihYdY5L7kvnXo62JOdSdhGf/tZSq/oYpwupmoLqZ6Vp69Hk/yZJUevLlth1PUwJ/8Qp3IeaLuql9Z21B0cCGn92b5qi2I1mU3k8lWZ9PMru3ZLZ35ZEttY9YJe+4drtbb9pyf/4Sm+iD3cE+faxRrjqb6T2PK56i75e++jLZtCaz/el0c4Tz99QpsXw0X76T6dBpuw7aDXAXoVL9+fuEEdpgwe2/aVx6I2keicdyADdD5Wv33ppP/evltxWvvtvvw5f/V+eW4m6HQZ5uLceT9rPSvs8p4218bFZ7Lrj/Lq7cuCTlr7s748wnnuVwq4dg2wsV67QK/z9ZnXBXfKcDNkNvc0YfbuHqaq1HBL79vpWE33o32/wEhl8DQH1ZRtjkeYlAC6mBRnEZm9lBQvtJzaZrCNItRZCYEuBupC1xE/RDRNMzcwp1mjd+rVJZ94PezNX0V8baihH67cjxDKJUA/X65uGiXTPe4vas+eNsaOOE/J+7gW0fe+j3U8vf+tRdnW7kr66T7ovf1sRj/elwf9dJ8Wyjvfp1NLjF4H6k39o6f4d6r4zmJgDmeP3mkXGX3F2SuwIQHQxYaUGaqCThfQacb9rKXY3T3tAujbiui0eyfKuLb05cH9Shv9co/36dSSo9eZ9De1o+5GkLNJYVNWwcNu/S8Cv7RwJ25H3H0tMTzUhm/E/TLjthZPXXZw7MVIrCchgC4mwVhOJuHm6YPaj//1csoLP2yjHLWiiwG6kO37WvZU9u/rR8mhT6+W/ZTrYU7+P6uMF2L2t9b+Pp0dnB7BHr+TpV1CaQTo50vTSJ480pvvZVvvZ3X+lGt2nhCFxFpRP91HLKef7evH+/Kgn+7TQnnn+3RqidHreL3Fdxbjc/FbtzAAACAASURBVDhPyj676OsrziP1NktFF9vTKzpdSKd997MWY0/3tAthbyumz+6dru/a0pcH9ytt9Ms93qdTS45eZ9LfhTrAK1nr4bf6B4NOeP7gH7V4VEHri4GYgdLFb43dnrtjVVmeY/6eFv+r+OiFoM55Gsxd/as46uAca3RxDurzlimd+h8Uftl+95SSsI1T6E2bFl3k8xQrTx/3fI7rmPL2n0zsBP8lX6L2mF16dSqdP+l62Jd/lEzx/MeY31Wvuf+oE4vczFrsJrWJXDAql34+FxbxiiMg+52tn166sqpL6zONZdH53n68L4+QD/30SOWK7+L9dI5O0etwhYqrp8v1+4NTn3EWt4mg79b+QnXr7SuGEyNFG4GuNoou2qiVfRydlq0fpPtGQHY66fWny+5dYk5/1pdHyIf70G8qXMVvjk7R6zyqvDKiTsqws83fvPAotTjM9U9tfu/9nhD/nfalJ97JpyUf/yo+meI0GaCLaTiWlIt0+ll9wc9a7KyrT6+VLSq2kY1q9ojoIg+xbN7/HvKfU1YxMqVPr671KXXJyT+SVVycdBHGCta2C/r5FSgKEa8QWFs/faUCtQM5/azbay3Z0W5OHk6gePTTR+TK3cnVKXodpcP4zmJU4nMm6rMLne/sK84p+9bKRhdb02h1jeQd4/bUSo16CPT1ZU7ed23JySPkw31ojz5KOZ2rU/Q6vcaORtTpwdcebo8iOEz3pWP+x7VH1j3U8aORa3VxFNf/bvUHjqvJ3+vn2YcABNZFQG3aTnsuputSG9KeQGBum1f+k/4D7oSqkrQQAue2ibltvhDMiLEhAtjshpS5kqqcu59eCaZViCld+r3GTT3f5PwJubVO2EQrGk5AAAIQgMCMBLj+zAiXrCFQAIH6iDqPpntek8uOOodZ/nmmTubq3JvfyuMXApsh0Oe8Vjtw+/LI1dx2ZufZz8q3c0raKdqX8tiMHqgIBHII5Nh8Rpv29dTTSNeD2/hNlfFT/YT2PRVT5wjWKdp0Q7kcmohAl11Id8XaRI7NT4SIbCAwCYE+m+1qixZA6Yu975oEEJk0Esiwi2L76cYKcbAi0KfXUzCVfO0+pV6khcASBPraJtfiJbRAGaUSyGgf3JOUqrwT5OrS+1x9osVV3vheTtBbX9IuvfalrZ8/jKiT0vxC8Z2W71WAnQBV0PFn2nAH0fvdOcX1lGEelXeUh/YJEIAABCCQQUD9qL89EZ043n6qPplpbDLYrSGK9MuIupqi9m7z2ETNIDawG2z6nqpiZ5BnqXitfrzzjzUbqDZV2DCBvds0/XT1gsfvCfwN+MN7gjWavHTpEXWexuswg9CYemATY6iRBgIQWDuBcD/Au4ozKpLrT3VPwjuzM9rgUkWH/mZ3z9TpiDr/i9//5q/ffLsT9vGcF8Ux7U2lidtL6ZByIAABCGyBwCv1tw9dEV2Y4h8oTnqZsAUo1GHTBLD5Tat3l5XzH9xs1+/Dw/Rv2r+7SxJUeisEsOmtaHJEPdSP+Y+4flGyheD3FDnvNbZQV+oAAQhAYGoCPLdNTZT8xhDADsdQW1+aXT5/pI46vxA++gZd8F7aU537cuHvoPcb69M/EkMAAhAogoBfIMTwlzZuqS/mW4GRCOstEsDmt6jVfdcpHXXyo1Ac3V/vGw21XykBbHqlijtVbN2D+h2B70e/nJpXQem3VJeCsCIKBCCwAwI8t+1AySuoIna4AiVNIOIunz8qR51fAgugv0UXv0cXv93wVsc87VrudD1xFF3aaCbQDVlAAAIQ2AcB9bfpHyP8cqRppPM+YFDLXRDA5neh5l1VUjb9NfzZzbNS3NF+NUp6VxCo7KYIYNObUmd2ZeI7Aunfo4Oz0xUe0e89GFFXuJIQDwIQKJOArge8qyhTNbuSCjvch7r3+vxROeqkYr8Mdniqm3B/k87BU63Ze5nrpHMa/p1mCgQIQAACJxIIL0deKBte8J7IkuTrIIDNr0NPSJlHINw/f5RdP9DyofZAmZcJsSBQEAFsuiBlLCfKE+n9l+WKm7ekcJ/hQuKfi+ctkNwhAAEIbJQAz20bVezKqoUdrkxhI8Td4/NHdNT5RfBHAfC/y065GY//TrujfJjmZ4QRzp0k6cg8TamnY3oZDH/uoskfAhAIBPraYTjvbxrdV/vkZQKWs3oC2PzqVUgFagS6bFrnnqnvru6ntfZIlDtabmk73ifXcmMXAucngE2fXwfnkKBN7+63JM8PWvv7IA6eMeeF9t+oLxvyR94qcSE/cfag+LmOQsRCjLkJtNn53OWSPwTWSKCvvYTzvKtYo3JXJDN2uCJlnSBql551bpfP1Be60fY0l/+Iqx02pzjpKtUor3+18avyenqCrkg6EwHp54Oy/tkPWNr2w4r/5f39TMWRLQQg0ECgqx3qnKfkqT6aqrVHKb9QG33ekA2HVkgg6PemdLqrF/bYfLux7tUm2oms40ybTQd9/i/eW2nfL7v/jPvrqB1S7pEANt2u9dCuN3ntbtN7nYbi+X3Bf9SXrfYPZKrDI9XBn/bwrEF/1Os4ZH/LNjGEw1ri5tr5WuqDnBCYk0BXewl9H+8q5lRAT957uf5ghz2GsJHTbXoOdv6/+Ayt/d08U3uyeTtrPmm5KwAn/ztO8JzXV+WVzl2sQ4RzE5BuPK3pj9LNY8sSDN8PXd/r2GofulwXAgTWQqCvHYY+NP7j19XyN+o8FTEBAqskgM2vUm0I3UEgw6afKLn7bWYv6ODIqXIIYNPl6GJJSfr0blkUx38g81Tsfo78Vcsr3Zeu8s9GqotfLrset9daB8lOGEggx84HZkl0CGyWQF970Xm/7+VdxWYtoIyKYYdl6GFuKTL0vMtnak996U7WL4JPdtIFJXrKS8MklEfAoxzTkTn3LKJ0j5OuPF0h0XYJdLZDtUecctvV/V5rhs3vVfPbrXefTftlNgECayKATa9JW9PJ2ql3FxOeE/38mD5DTifBsjlVL5dVp1U6GpdFtanSeu18U7WlMhA4jUBne+FdxWlwSZ1NADvMRrXqiH163uUz9TV1tO8n7mw9tWL8992qLWZLwstT7X9D+uEk/XagR9adNO3HlhhRFwjMTYB2ODdh8i+NADZfmkaQ51QC2PSpBElfGgFsujSNLCPPTvX+QHTTZ+FlYFPK2Qjs1M7PxpuC102A9rJu/W1FeuxwK5rsrgd6bufjEXVTh79Chr4Rxgk0Nd3x+Vkf8V+RdqR6flcfY4pSgyFAYBkCtMNlOFNKOQSw+XJ0gSTTEMCmp+FILuUQwKbL0cWSkuxK73r29R9W/cfVd0tCpqyzE9iVnZ+dNgKsnQDtZe0a3Ib82OE29NhXC/TcQuhCo+laTo0/rBthf/fMI/Wqb6GNz4mUUxGQTl4rr5ta/i/k6en1nktHTHsZgLCCwNwEaIdzEyb/0ghg86VpBHlOJYBNn0qQ9KURwKZL08gy8uxN76qvv03nb9TxfbplTKyIUvZm50VAR4jVEqC9rFZ1mxIcO9yUOlsrg55b0Xw3l6POTqEncgJdtBfNmSUJqBH4o692zDHKcUnwlAWBhADtMIHB5i4IYPO7UPOuKolN70rdu6gsNr0LNV+p5N70rvp6JN0tPQvzLegr1rDdA3uz8+1qkpotQYD2sgRlyugjgB32EdrGefTcrse5HHUewuib4Ye6GWYe+Hb+i5xRA/A0Hx7l+L30wQi6RahTCASOCdAOj3mwt30C2Pz2dby3GmLTe9P49uuLTW9fx0013Jvek/r6T6u/NDHh2PYIJHrnHcj21EuNJiZAe5kYKNmNIoAdjsK2ukTouVtl17pPjzsbnHOflZqpL8chnDqV5+T/iJNuaqzkB4FBBGiHg3AReQMEsPkNKJEqHBHApo9wsLMBAtj0BpQ4ogp70/uTwOjXEaxIsl4Ce7Pz9WoKyUsgQHspQQvIgB3uwwbQc4eeZxlR5/LkIfUNsafA5B9MBkKAAAQgAAEIQAACEIAABCAAAQhAYDECei/hmWXe60+r/Il4MeoUBAEIQAACEIAABCAwlMBsjjoLEm6Kf9dN8dOhghEfAhCAAAQgAAEIQAACEIAABCAAAQiMIaD3EY+U7q2W23on4Rl/CBCAAAQgAAEIQAACECiSwNyOunhjzKi6ItWPUBCAAAQgAAEIQAACEIAABCAAge0RkKPug2r1F38c3p5uqREEIAABCEAAAhDYGoFZHXWGxc3x1kyG+kAAAhCAAAQgAAEIQAACEIAABMolEEbT/SYJ/yNH3ddyJUUyCEAAAhCAAAQgAAEIyI+mm9ZZOegG2R8J/KTlrsr6OGthZA4BCEAAAhCAAAQgAAEIQAACEIDAbgnoHcQNVf5/Wh7rHcT73YKg4hCAAAQgAAEIQAACqyFwbW5JdWPsueD94WbPDU+AAAQgAAEIQAACEIAABCAAAQhAAAJzEfC7h99x0s2Fl3whAAEIQAACEIAABIYS8IwPWp61pZvdUeeCdYP8h1bvJQjOujZNcBwCEIAABCAAAQhAAAIQgAAEIACB0QT0zuGJE+sdxNPRmZAQAhCAAAQgAAEIQAACCxNYxFHnOsUb5S6v4cJ1pzgIQAACEIAABCAAAQhAAAIQgAAENkBA7xoeqBpP9e7h4QaqQxUgAAEIQAACEIAABHZEYDFHnZnqhtlTYN4ON9A7wkxVIQABCEAAAhCAAAQgAAEIQAACEJiDgN4x+Lt0HkV3f478yRMCEIAABCAAAQhAAAJzErg+Z+ZNeceRdU3nOAYBCEAAAhCAAAQgAAEIQAACEIAABIYQ0HuGr4rvPwYTIAABCEAAAhCAAAQgsDoCi46oWx0dBIYABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIDATAQu9M+zmbImWwhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhsn4CmZH+lWvrbyfXgqdpvavlcP+FjF/rBU9dAhkMQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQaCOgwXD2s3UGOfAeKcItxf2lKeL1nEyaEnIMAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhAYT4Bv1I1nR0oIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIjCaAo240OhJCAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAYDwBHHXj2ZESAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAqMJ4KgbjY6EEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEBhP4OLff/8dn5qUEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIBAI4GLi4sbOnFT/rjPTRFw1DVR4RgEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEZibA1JczAyZ7CEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCDQRwFHXRIVjEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEJiZAI66mQGTPQQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAASaCFxvOsgxCEAAAhCAAAQgAAEIQAACEDg/AX10/IGkeKzlk5YftbzUB8g/ak2AAAQgAAEIQAACEIAABCAAgQ0QwFG3ASVSBQhAAAIQgAAEIAABCEBgewTkpHulWt2QY+6paxecdn9q83vvEyAAAQhAAAIQgAAEIAABCEBg/QSY+nL9OqQGEIAABCAAAQhAAAIQgMDGCMgp90RVehSddKF6n7W+ERx2G6sx1YEABCAAAQhAAAIQgAAEILBPAjjq9ql3ag0BCEAAAhCAAAQgAAEIlE3Ao+m8pOFW2LmRHmQbAhCAAAQgAAEIQAACEIAABNZLAEfdenWH5BCAAAQgAAEIQAACEIDABgmEEXN2xv1eq96dsM836mpg2IUABCAAAQhAAAIQgAAEILBWAjjq1qo55IYABCAAAQhAAAIQgAAEtkrgsSr2WdNefq1V8Kdw3FNgEiAAAQhAAAIQgAAEIAABCEBgAwSub6AOVAECEIAABCAAAQhAAAIQgMCWCDxQZd6nFdIoO4+m83I3Pc42BCAAAQhAAAIQgAAEIAABCKybACPq1q0/pIcABCAAAQhAAAIQgAAENkRADjlPeelv0XmpQjj2VjtPNcqOaS8DF1YQgAAEIAABCEAAAhCAAAS2QIARdVvQInWAAAQgAAEIQAACEIAABLZCwKPpHJ7KQffs2+Z3t7V+jJMu0GAFAQhAAAIQgAAEIAABCEBgQwQu9LC3oepQFQhAAAIQgAAEIAABCEAAAuslIOfca0l/T89pTHG5XjUiOQQgAAEIQAACEIAABCAAgWwCTH2ZjYqIEIAABCAAAQhAAAIQgAAEZifwX5XwZvZSKAACEIAABCAAAQhAAAIQgAAEiiCAo64INSAEBCAAAQhAAAIQgAAEILB3AhpN5+/S+Rt17/fOgvpDAAIQgAAEIAABCEAAAhDYCwEcdXvRNPWEAAQgAAEIQAACEIAABEonYEfdZ75FV7qakA8CEIAABCAAAQhAAAIQgMB0BFb1jTr9w/Sdqu6PqH+dDgE5zUlAOruj/O9p8T+Db2t5zYsHUdhQoF1uSJlUBQIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABBYlsJoRdXIGPBEZO3wI6yLwSuL6X8G/aG1H62/rEh9puwjQLrvocA4CEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQg0E3gevfpMs7KGfBAkvyl5UsZEiHFAALpCMgflY7vbQyAV3JU2mXJ2kE2CEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQWAOB4h11cgZ4ysRbGpH1XttrYIqMCQFPUyq9efrLn7Tc0f7D5DSbKyVAu1yp4hAbAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEiiJQvKNOtJ6EaROLAocw+QSkv4+K/VHOnQdaPmj/bn5qYhZKgHZZqGIQCwIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABNZD4EJOkyKklQPnlgTxyKuvkqmaHlHH4kisKKO/U/erljfB+ROPsz4zgSb9WSQdfyZd+ft0VdC+De62jn0Oh1gVSqBDp7TLQnWGWBCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEILAuAkU46uQQ8DfoPL3lr8E58Er7P2v/a4pT5/7R/n/qx9M459y27Es5oFSWp5H0SLWzB8nSqL8g2P8k5/fetsxa/Rn3w/mzrfaqrxzgbTqV7o7apPNS3EXb5VJ6UznFtLEcnREHAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAE1kfgWiEiezTdX5ZFjgCPtHqj5YX3HfTC/IYWO+/8vbpXflHv4yUFyxdkn1ws5e2RhPVwS8ftICshNOovOHWem42WZxLUOr1fgsCWacf6ylFBo07ThGK4eLucS2/Kt/Q2lqJnGwIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABDZCoJRv1P0unk/0sjx1wHn0VRWiw0c7z8OholaS++3Mst2uV1hM/lC5r7V8nsvhVC+zY79Vf5LNU5UWFdBXljpadRpTL90uZ9Zb6W0sYmcNAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACGyJwdkedXr7bOWdH1/3w4t8j6D5ovxphVzpryWqH4rmcZXZc/qblcRcnOzjEtjNOV/quc2vT3xr01cV7iXMl6vSMestqY0vohTIgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEBgewSuFVCl15LhZXTSJfLYebeG4OkcPVXnLEEOCk/3+akp88DMU2CmIxGbojqPucLa9LcGfc2lq9x8S9TpbHqbqI3lsiUeBCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQOBAoNVR55fXWqpvi2n9VsthKspD6mk27GT6GLNSOZXTSU6o9/HYUmuV/UCL63xwEloeLR7hF7+V94/2U8fXHcl6kL8uq+NqGcxRaTwVqL/L5xFzt7XtaS6bdGAn4Tm/VXc2/YkH+qob3DT7s+p0ar0pv623sWm0Si4QgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCBQHIFrHRL9qXNv5IT6RWtP/+b9OYJH7zxyxn6Br5WdU/e9f4bwUPV1XR/55X8o37J98bbOfdXqs5Z73lecyqno7Y4wmKPyfef8giyeVtMyebHDNMrlKA6W52G1dZ6fc+oPfc2j87l1OrXett7G5tEyuUIAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgMDZCTQ66uQMeibJ7CCKI8XsqPKolbqT6OQKqAw7At8r7yfOTPuPtdghtmhQ+Xa6vdTazkLXPcrwo/Yrx1kQyCPY7BxzcJq4XR1If5TXYI4hzVeV/2uaV5DHZdVHz1nOyfWSlt21LbnOoj/01aWV087NqdOp9baHNnaaNkkNAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACJRO43iLcUx33CK4YqhFkwVlUHdML8n/jySnXynfK7HrzUp2qArWuHG4q/7ES/ZEktGPMI4xi+CHGDQeiQy+eT9e9HNPIYdtpLEM11aZWf3s7hCbHoJ2oN2MEye9pOx0vDYfpO5ODcbRkdWgqfc6tv63pK9HHWTZz9D6FTmfU2+Jt7CyKolAIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhDYJIErjjq9lPfoLDt60m/E1Z1XHvk2iUdN5b1TXtXUjdqunGLavz0l7YFlWAa//I+OMvP4y/shpI6zaqRhPJGuczmmacK2nW7R+WcHaaUH5Vd9n05s4ijHmNTxD6P6dL5y8sWTXqf1T4+n20o3iT7TPMdup/Jqu88mVq2vsYymSje13qUvj4z19LX/Ud7RjpvEPVlvKiurr2oo/KQ21pAfhyAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACowhca0jlF+jxe2x28thB5GM/+/gMwS/NY7BDzKO//AJ+yjCkDDsX7IBzqLOwEyIdbWcHWX30mtM51NNe4ah6uq7vavX9PabV2t/y+uh42v7N+1rqwecOjrr6yZXuo68VKk52apt3G47tp6sWU7Qz2lgXYc5BAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIBA8QSuN0hoZ9Afeunu76s5eHTb3a7RMcGRZEeUv62WjsRz+s6g+HeTCH7xnn4frjo1JH/FfaI869936y0jkcEOyafK55NlCdue+rLaV94Hp5i2vypeNfLO24qThhyOdkh61Nx/tVQyKx+X/UqL+d8J6x+0fb+hDB3+znXzdJeTBZVp598ofUYhcvNQPPQVoZ15fYrOlNa2fEs26u9N5tRkina22jaWA4g4EIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgsH0CTY46O8ue64V7OnKslYReyju+X9D/6hf9Wuw0+rnFqdSVj1/0v9Dil++HMCJ/59MYlFdjGWlkye2pJaupL9PjHdt24plBnVcvx1DW95LrUZq/jlffB9Rxj2z8JT3XsH1PcYbI25DF5aERvC8Th62BeaCvKwSXPzCBzuxw7bPVQ8UUd4p2tso2doDABgQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCOyewLWUgF7W22ni0VRDRsU5vqe7s1PJo83eaLHDLTuEcj21o0eNHUashQyG5v+1qeCeMpqSZB0LzokjR9kQjoprJ52dFkdBxz2i7crxNFJIa0fhlGEo76ayh+SBvpoILn9stM6Crf6gtUeC+vt0nrr0RTg+SU3q7Ux5Z/dViltaG5uECZlAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIDA+glcr1XBL+s/6qV4o/OkFjfu+ptqT/Qy3GljsJPpKOh8OpXma5VROaHCC3e/3I+jyF7pXLUdMujNX3m4vJ9igdp3fi6jcvpllBGTjl2/VhnpFI5DON6MctYK90i5oyk80/OhTv6G3ZGTMI2TbOfEidF7eTuiym/UZ8ikMw+lRV+Rdjnr0TqTDbotH5zKbgvaf6njQ/qRHBJpOyutjeXITxwIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAwBGBC71MPzowZEcv5P2y3FNdHr6fpmMftP9X6kDSMY/6eqtj/n6VnTR2xlVTXGrf335LnXz+Rp2/i2dnUG/+iuPp754qzWNtV84yrT3axuX1luFypggqsxq1ozLrIwKnyP5KHirPjrJfVd5kzhDl2cvbgihelz4781Ba9HVFm+c90Kf3HJ0Fu3C782jayja1djuftD1IlsXamcqavI2ZEwECEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgEAnUR9TF47lrO2yaRs7YeZeGB3phH0d12ZFzcC7peOWUSyMn2zn5e/RdzNuOghhuxo2eMmK0k9Yqo/6NupPy60us8n7pizPifA5vZ9uqT53rywN9jVBMmiQ4q/xdyKls4GSdWT7J43Zt/XqZJSzZzibkOwsLMoUABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAIH1EzjVUWenWzrlnff9wv7wjbswGicdVWOnWu6L/N78ldfB6eeyQ/DUly/jDutsAr28M/TZlwf6ylbHYhHR2WKoKQgCEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgcEng1Kkvq++UeeRJcODYCfez9g/OGB2308yOgHdaPHruXerI035rUNqc/G8ogxda/tbi/D2Vpr+zd3AWap+QQSCTd6c++/LQefSVoYuuKGLo6R8nG1GHzrpocw4CEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgMB+Bkxx1Fksv+f3NuXtaMMdfiQAACspJREFU/G25K84xnfc36x7rXDqqTofyQl/+MRfFe6Lt31XOwUkYz7HOJ9DHO0effXlYGvSVr5N6TLGb1FEX9NHZjkMc2lhdGexDAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABE4gcOrUl57m0lNfHqa/bJJlrJPOeeXkH8vESRdJjF/n8O7TZ04elhB9jdfT1CnR2dREyQ8CEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQg0E/g5BF1/UUQAwIQmJLAHCPqppSPvCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEMgjcOSokwPg37xkxIIABE4loFFsF115qD36e4APGuL4O383tTRNJ+spaB83pDkcop0fUIza6NPbqExJBAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAArskcOSo2yUBKg2BlRHYw4g61dHfzPspqMbbT+Uga3JMrkx7iAsBCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQuCRw8jfqLrNiCwIQgMBkBF7JMffQuclp51GF77Tc9j4BAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEILAVAte2UhHqAQEIbIqAp/aM4S9t3JLDzlN+EiAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACmyGAo24zqqQiENgOAY2mu5vUxiPq/O29r8kxNiEAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACqyeAo271KqQCENgugTCK7oVqWE2Dud2aUjMIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhDYIwG+UbdyrSeOjL9VlR+1vNTIo48rrxbi74BAn+2G878JxX1G0+3AIKgiBCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAIEdErjQC/AdVns7VZYz44Nq87Odc9q+pe0P2v5+OzWkJnUCwYF1U3r+XD+3pv0u2w11fKX6ePmi5YXq+3xN9UNWCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgEAfARx1fYQKPi9nxjOJ96McGI8tZnBu/KPN73WM73kZCqFIAn22q/OfJLgdzzH4G3W34w5rCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgMAWCPCNunVr8anEf5NU4Z63cdIlRNgslUCn7dopp8V/JIgLTrpSNYlcEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgMJoAI+pGoztvwqbRczr2WlJ5SsRqhN15JaR0CDQTwHabuXAUAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAE9kfg+v6qvJkaP3BN5JSrpriU8+OOdn3sro8TIFAwAWy3YOUgGgQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCCxHgBF1y7GetKQ4ek6Z/l/I2FMDPo+Ou0kLIzMITEgA250QJllBAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIDAqgngqFup+uTs+CTR7Zj7Y6VVQOydEsB2d6p4qg0BCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAJXCOCou4Kk/ANydNyQlP9o+Z4RdOXrCwkvCWC7lyzYggAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhcA8EqCdyS1B9x0q1Sd3sXGtvduwVQfwhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEDgQYETdAQUbEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEFiOACPqlmNNSRCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhA4ELjQ1r+HPTYgAAEIQCCLgKaedf9JgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQiMJsDUl6PRkRACEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAAC4wkw9eV4dqSEAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAwGgC10enJCEEFiZwcXHxSkXe0HJLy2ctzzX94FetCRCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEBgdQSY+nJ1KtunwMFJ91qOOTvovtP+W61uaf/uPolQawhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABNZOgKkv167B/cj/qFbVl9q/I4edR9cRIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgMDqCOCoW53Kdi1w6pSLU16mx3YNh8pDAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCCwLgJ8o25d+tqttJri8nat8nfC/l+14+xCAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAgVUQYETdKtSEkA0EXunYcznw4si6higcggAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAALlEriQo6Nc6ZAMAg0E9F06O+luyHafNpzmEAQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhBYBQFG1K1CTQgZCchJ90TbOOkiENYQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQislgDfqFut6rYnuJxwD1Srx1o+aflRy0uNmvuodRXC+dtxJF3Y/5LGiXFZQwACEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAIHSCTD1Zeka2ol8crodTWcZnHBv5YT73gi0f0er37S89H4InvryseLwnbpIhDUEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAKrIYCjbjWq2q6gcsJ5OsvncrjdjrXUsVva9si6hzr+Xvv/aPtGPB/XOncRt1lDAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAgTURYOrLNWlru7J6NN3zWvXsqHOonHNyyFUj674d4hcCEIAABCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIrJ/AtfVXgRqsmYBGyvm7dHbG/V6rh6e6dDh8o+7bLr8QgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAYBsEcNRtQ49rrsVjCf9ZI+bq35n7KRz/vObKITsEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQaCPA1JdtZDi+FAGPqHufFqZRdh5N5+VuepxtCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEILAlAoyo25I2V1YXOeQ85aW/ReelCuHYW+081Sg7pr0MXFhBAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACEIAABCCwPQKMqNueTtdUI4+mc3gqB92zb5vf3db6MU66QIMVBCAAAQhAAAIQgAAEIAABCEAAAhCAAAQgAAEIQAACmyVwIYfIZitHxcomIOfca0l4TzbIFJdlqwrpIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQmIEAU1/OAJUsswn8VzHfZMcmIgQgAAEIQAACEIAABCAAAQhAAAIQgAAEIAABCEAAAhDYEAEcdRtS5pqqotF0/i6dv1H3fk1yIysEIAABCEAAAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQmIoAjrqpSJLPUAJ21H3mW3RDsREfAhCAAAQgAAEIQAACEIAABCAAAQhAAAIQgAAEIACBrRD4f/ALC8d6Yt7cAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$\\displaystyle \\frac{\\rho \\left(\\frac{\\rho \\frac{\\partial^{4}}{\\partial \\rho^{4}} w{\\left(\\rho,\\phi \\right)} + 3 \\frac{\\partial^{3}}{\\partial \\rho^{3}} w{\\left(\\rho,\\phi \\right)}}{\\rho} - \\frac{2 \\left(\\rho \\frac{\\partial^{3}}{\\partial \\rho^{3}} w{\\left(\\rho,\\phi \\right)} + 2 \\frac{\\partial^{2}}{\\partial \\rho^{2}} w{\\left(\\rho,\\phi \\right)}\\right)}{\\rho^{2}} + \\frac{\\frac{\\partial^{4}}{\\partial \\rho^{2}\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{2}} + \\frac{2 \\left(\\rho \\frac{\\partial^{2}}{\\partial \\rho^{2}} w{\\left(\\rho,\\phi \\right)} + \\frac{\\partial}{\\partial \\rho} w{\\left(\\rho,\\phi \\right)}\\right)}{\\rho^{3}} - \\frac{4 \\frac{\\partial^{3}}{\\partial \\rho\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{3}} + \\frac{6 \\frac{\\partial^{2}}{\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{4}}\\right) + \\frac{\\rho \\frac{\\partial^{3}}{\\partial \\rho^{3}} w{\\left(\\rho,\\phi \\right)} + 2 \\frac{\\partial^{2}}{\\partial \\rho^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho} - \\frac{\\rho \\frac{\\partial^{2}}{\\partial \\rho^{2}} w{\\left(\\rho,\\phi \\right)} + \\frac{\\partial}{\\partial \\rho} w{\\left(\\rho,\\phi \\right)}}{\\rho^{2}} + \\frac{\\frac{\\partial^{3}}{\\partial \\rho\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{2}} - \\frac{2 \\frac{\\partial^{2}}{\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{3}}}{\\rho} + \\frac{\\frac{\\rho \\frac{\\partial^{4}}{\\partial \\rho^{2}\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)} + \\frac{\\partial^{3}}{\\partial \\rho\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho} + \\frac{\\frac{\\partial^{4}}{\\partial \\phi^{4}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{2}}}{\\rho^{2}}$"
      ],
      "text/plain": [
       "  ⎛    4                3              ⎛    3                2         ⎞      \n",
       "  ⎜   ∂                ∂               ⎜   ∂                ∂          ⎟      \n",
       "  ⎜ρ⋅───(w(ρ, φ)) + 3⋅───(w(ρ, φ))   2⋅⎜ρ⋅───(w(ρ, φ)) + 2⋅───(w(ρ, φ))⎟   ───\n",
       "  ⎜    4                3              ⎜    3                2         ⎟     2\n",
       "  ⎜  ∂ρ               ∂ρ               ⎝  ∂ρ               ∂ρ          ⎠   ∂ρ \n",
       "ρ⋅⎜─────────────────────────────── - ─────────────────────────────────── + ───\n",
       "  ⎜               ρ                                    2                      \n",
       "  ⎝                                                   ρ                       \n",
       "──────────────────────────────────────────────────────────────────────────────\n",
       "                                                                              \n",
       "                                                                              \n",
       "\n",
       " 4                ⎛    2                       ⎞        3                  2  \n",
       "∂                 ⎜   ∂             ∂          ⎟       ∂                  ∂   \n",
       "────(w(ρ, φ))   2⋅⎜ρ⋅───(w(ρ, φ)) + ──(w(ρ, φ))⎟   4⋅──────(w(ρ, φ))   6⋅───(w\n",
       "   2              ⎜    2            ∂ρ         ⎟          2                2  \n",
       " ∂φ               ⎝  ∂ρ                        ⎠     ∂ρ ∂φ               ∂φ   \n",
       "───────────── + ──────────────────────────────── - ───────────────── + ───────\n",
       "     2                          3                           3                 \n",
       "    ρ                          ρ                           ρ                 ρ\n",
       "──────────────────────────────────────────────────────────────────────────────\n",
       "                                                     ρ                        \n",
       "                                                                              \n",
       "\n",
       "       ⎞       3                2                2                            \n",
       "       ⎟      ∂                ∂                ∂             ∂               \n",
       "(ρ, φ))⎟   ρ⋅───(w(ρ, φ)) + 2⋅───(w(ρ, φ))   ρ⋅───(w(ρ, φ)) + ──(w(ρ, φ))   ──\n",
       "       ⎟       3                2                2            ∂ρ              \n",
       "       ⎟     ∂ρ               ∂ρ               ∂ρ                           ∂ρ\n",
       "───────⎟ + ─────────────────────────────── - ──────────────────────────── + ──\n",
       "4      ⎟                  ρ                                2                  \n",
       "       ⎠                                                  ρ                   \n",
       "──────────────────────────────────────────────────────────────────────────────\n",
       "                                                                              \n",
       "                                                                              \n",
       "\n",
       " 3                  2                  4                 3                4   \n",
       "∂                  ∂                  ∂                 ∂                ∂    \n",
       "────(w(ρ, φ))   2⋅───(w(ρ, φ))   ρ⋅───────(w(ρ, φ)) + ──────(w(ρ, φ))   ───(w(\n",
       "   2                2                2   2                 2              4   \n",
       " ∂φ               ∂φ               ∂ρ  ∂φ             ∂ρ ∂φ             ∂φ    \n",
       "───────────── - ──────────────   ──────────────────────────────────── + ──────\n",
       "      2                3                          ρ                           \n",
       "     ρ                ρ                                                      ρ\n",
       "────────────────────────────── + ─────────────────────────────────────────────\n",
       "                                                           2                  \n",
       "                                                          ρ                   \n",
       "\n",
       "      \n",
       "      \n",
       "ρ, φ))\n",
       "      \n",
       "      \n",
       "──────\n",
       "2     \n",
       "      \n",
       "──────\n",
       "      \n",
       "      "
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "L2 = lap(L1, coords=coords, h_vec=h_vec)\n",
    "L2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAssAAAB2CAYAAAA+w95ZAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAgAElEQVR4Ae2di/HcNBfFIZMCklAB0MFHqADogEcFhA5gUkEGOgAqCKEDoAIeHQAVENIB3/npLzmyVl7bu37v0YzWtiRLV0fH8vX1tfbN//777w0HI3BUBN58883/qW+PFR8ovqv4rTj/h7YON45A5MZnEQZ48oW48deNw+LuL4yAePiTmvxE3Hu1cNNuzggYgYEI3B9YzsWMwF4R+FqCf60b0c+6KX2s/e8V39trZyz3pAjAi4+oUdz4UBuUFh6oHIzAIgiId0/UEA/zDkbACGwYASvLGx4cizYJArnF5n3V+PMktbqSIyDwKOvEb9p/R8rLA1v4MlS8OxsC8QEN3r2crRFXbASMwCQIWFmeBEZXsiYCuum8o/Z5jf4KC3IuC4qP8snjdfv/kiUxL+P94yLQw438DQOW5b+sKB+XC2v0rIt/Ssct7B3mK+2vIZrb3CkCXZzaaXd2I/a93UhqQY1ABQFNHF8qmZvOj9q+1PELRW5ETVDeH4pfKeFr5f3eZHjn0AgM4QYARL481W5wyTg0KO7cYgj08O+J5qTvFhPGDR0CgR5OHaKPW+3Em7pgtyqb5TICZxHQxIHF+IU43PiZKg0LIR9qfcLJTC7a/4Z9go4h/LtK84dcAZFj/gzhBj1XOR6s8GP/XJzwB1aA4nA1Auf4p8qfKaYPS2kLv2UU5+fioD8+BhGHEwTOcUq8Cfe7k5OcMBkCVpYng9IVLY1AVHQea6JoXC9i2r9K025QhP7W/kNk0zHK9S/peGl53d5yCMSx7+QGksQy4QNQHeI3+lTc4A2EgxG4CoEh/EsNqOy/2n9b3PPDWgLF2xMExnDq5GQnXI2A3TCuhtAVrIVAvLn8pUnkSYxYlR8hj45xzeDm85X2cb/AXYNX7R8oLh6QRxEL0uKBtstGkaWWXpbb63EfN2K/cMlhTP5URGFhtZTFA+PAeCzesBqscQBZaulryLfXNofwTxg/UORhjbcbzFEn1+mW+4+8iubtQoM0hFMLiRLmjVsb+11ZljU4P4kw9itc6oqotMMEqeRkfUuT+1cal+rrQ5VHSWUViqSIYAVOFhRuEoSg1KqOUa4RqW6dl1wukOcLRdp8qPTUjg7XCxEz1ndehbtqnyXS0pg1QCj9hQ4Yu1G4NxVMsJPGUFVNyo9Ur/pmbpwZpy1yQzJtZo45A93ZrL3w72wnzmTGMfKcVmCUxl3Jk85nNJPqXntOu9mxF/Bv7CFG8uFvytJOu5D5aHKCvSI+wg3+OsYywrh8mKeX+7HM72U6xwpYJ9j5uJZfS6Os4p9lntKwFuKG0ci49n6UCUv3KjKpfZTlk7aVznhWx6RWfs40ycH4V2VR+ih+qLy5URnv2vgJq01xQ/JsZo6p4TUkbU/8G9KfWhn1kXnWc1rHdSZsJpvPwF9hM3OaZLnJsb/HKGw96EmGCRTrpMO6CIQPoXIRdCFjsXyliJWyGjR++AoTnt9t2r+qg49biNQ/NKCkswJGLWzmK3P1HUXvN/VxFettvHb+qYEkmRg3lq7CEr9amIEf5saA0dwoN7Y0xwxAsVpkF/yrSj4gUbzxnHYGpxnmM1rbBKdueux1wzyxOC2VJgJw0XGjhgjE6pNqLIPCzCOWLcsrjZmwx6+zZs1FUWZsWMf4hE9KZ4w78zlHAX9jdlhS6aSOPE1lEhda1myl8/qWOjbDEcmCP2yV13mf5tpX2+BaHRfaVADLVS3xan8yfsT+0DFzo/862hw34KLi6nPMpdfjnvh3RR89p525tsSByeazbI7exJzGtal4k/ez1SzLekJBOf5ZZPhGEesky+l8q3Qm8CbEp7RVrHKNEN5JCDAO4V/OUsLAbfDV1ThX/ZpjHWmMm2XguupWPVhECaxgkAc4hQ9uys/zFt8Xd1HeUfxT3xaXQQ1+dA73iBUfSbauu4UFnYwf2dibG/2DuEVubGKO6YeuXmIv/KtL35/qOa0fI5WYbD6jta1w6tbHfhVlWaBjUcby0ygRkRAozeXr/M+U1/W6HS45LISAxuE9RT4KLZXR5GbRjGchEopYs7xbkZcOg2Kpg7RN6V1bPuR7Ki7xRTnxWx3zD2zNmspdJw5Nj/V+Sd2KWAuaoGNWDMB3KwTtf6z4ryKW2hTwM+vrdyir8z5UDO1oyx+rJExTXYO3nKtIfWAJNjzgIF8XtsgYPoQb3Mi0Bafmh7nRMT7iwKa5sbE5pgPF3uTZ+dcrQaWAxp654Jr5jFo9p1WwLZKmns+o/ipOeeyLEbrkUJPT2Vfec+RLThSb6qtfpfO6IbyKZ6sYTP7aooSQt5lX7HNgs7c6NR5MDIxL14dCKH3k8+cgnVwjP5ar1lM7V+VR/uBIw5NauUvTVG+QRVtuEC2+6vgnxeZjR+0nHJrX/7HMELcSLOJ8WR7wiXW12kt5fVvqUQwfSmob8Enn6Lj6YYbSwa/6cV06d66t2p2FH6rX3CiuN2GyK24kzknudG1V54a5OJTav2Q7N/8ulOmq+Yw2FZj3PKcV11YaD2Ezy3wWsb94TpNcHvuOMUtj17e9JxAXDxKKf1gLfxRRaRzr5Hs8CWmLYtxlrayc6qQVEOAG/KPG6WRpsigLNzpCn4U1fcDJZDwowA3F72KclCfiHxPTr1EQXquV9dOvRlbJQP+IeTnqKN0BlPQ6qB0UVZRbLAcpUAdWoIRdSj+7VXnkQelNb2I4P8edfRTzMiBjbhEv8+c8Tn3M5ay1N4ofwsDcyFDcKTdSD1abY5IAY7dz8m+sLJTX+E8xn1GV5zRQ6A6zzGc0dymnPPbdgzUqRwPQae2bO0+C8hSWlIVkVePxlZs+1jwmyTySx3Hvk+3csrv+YGUIY3MOC40VY6ki53lGGcWqNTXmkT97THKqreYNBnIpNpZx7SfrQetDB6U31mHqUeC8xtKc6s63sUyLz5yjSAVjltLjOmp9GKXjliVOx4xXq0yUM/Qnl2upfclzFT90/uycSG0kTHRsbvRczwmra7fCevY5Jo3vEbcJ/yk4S10KntPOcF/4bGY+89if1zkSPkO392H/0kFPOtycv1dkIvxB8ZHSsHg9VyRgFcI6lixk6ckYhWAzH3Ah6K0GjVfw4dU45RbRGhwofmethqqLcSVUrdNq48277Gl+I/8eqzasqXxQiJL7R167jl9xrLLIT7mGi9on7ZXK5FZkJTV/tsJ+b8jq5hrIA9cHoSXTXVLnLz7HZ3FWPn0eU2dnYxNmXMWPmbjxWewf48BbsNY4mxsTjv6ZqnR9LDLHrMGhM90enRXns9U5i+Ce08K94ew8LIw673dTc5ExWWK+uomxF5CLWpYFKjdHnk5PlrNSWngq07ax4iGfAuRKeXwAONjitnT/bqE98Fcsrai8nistrcFiqfTWeOYYKQ9FFD7w74yDuKiytIUMZ622XfXRVjo31tPpr6v8E2tsPL/xV6YdhSBT3qbSWGank6vKO6k71oVv8YkFOK+73Fd5MGws1NpHnvIYQU/kURrX5Kj2yvYvOVabk/NDdV7NjdSXPlyUfzJ+SoNb5sbAazlhXW6F4WpzjNpekkOj2pJszTWdMFNaM3dq/+y1rPyLOUt7Cp7TOrgtbDY3nyWOxLHz2HeMXY5T1/492L9U0BMVEwM3k2cSqGbhIp/QylNZ/FJZ5oiVGD5RzK18d2f4dxEEohWDsSgtytzcysDETag+aasuFOVfFHmTEJbbofC5oHOwNqGUw4GXOmblCOoZE+BQkglf2LRfqwNOtvioY9Ja1kYdY4EseUmZR4pdAXxabasvTLhELMVjAm3lOJR1M1FyHZUy0gb9ecXOwmFSfkzEjXy8fhMe55ZKNDdmIEy8BlaZY5bk0IVt1ZbWXIqzjLbntG7Ob3E+y6W9Zr7y2OvmOciaN0U5oZ3+vKLx+cvrVT7Wsarfal7O+8uNWY61xoaLDcsCLjMpooQxrieWSaVhIVUVp/IqnYmFuloW6lrZlKayKJGtdnRMPS1LXip/bhvrog+NVaZWXvko543lWfu017KE6xi5TqzcSqP+av+UjmKLqaZpP6aBSct6FNMp28hRyqq8gGdK13Hjr6x9HlCrclCePMWmfKpj7q3anIwfqmsybqR+q04eAFt8S3lsFcyNyrWdYzR2X5iuNseo7cU4dGlbOu/sdar82TgbOe85rYPzwn7r89nF85XH/r83sNRqjJcJepKmsT/U5ntli8pDeUAJSX9SUhbx8coIaIyYDLih1ALW4WD1UDmUZ256qWxpzWSssVCgwJVWWyXXQ+TIY53TWGMTb5R2kV+zzkfJ5AZ0wskkhcqgTBKwvv6jmP6a+1ftv6WIQnXyF9uxbvp4Yg1SHjc1cCIvWeXZr2Ki8ijRBNa6Ro6ToDLgjcWfccIyjZLcWafyQtB5lP9c9Q4ei3jqRRu1Nzk/VCecmowbsb5fVCdvIeBqNaicuVFF5rJE4bnaHLMkhy5pK57zqfh4MteAdsyfjbOxDeZLz2mAEYNw3/x8lsl60XzF+ernTY/9YspyvJBRhnkVXL7CZyB46uGp9aHyq8qA8hxuHAHxBCWci5aAEkNEkXw3V2oi354qHeUWV4vG9QeuqWzz5yU65iGudb6OJwmqG9l4pdxSuJTOpIVy16mklwLoHJRqHjZbdVXKBYxUrnpTzcurTspi3T5R6PNye9iPfbmKG/RT9aB4f6/IA8Rsc5HaMTcAfENhSQ6NaItvdrg+uVa59uFn66F6Kc6qXa4P8xYgZg4j+AEfqve6qUW86bHXzWARNwwNGgOKUlL92EvpXIBnXzEtJavbWYYTY3GGO4qNy4X2uXnwgAWvWq49Og4WKtpQoFxw79EWHjauPtrHGtscj5Wpr7zq5kbHja11ndGmYvVaKMumY5UfdH2oHG22PrZMdZRb6lQcJUdZxxaO6YPiVdygHwpBEdEWzrA/CPNLMFDd5kZxXVyC41TnaDwW49DQtlSOt0PBLUvbwEVt4SX3yzDnxWMevmfnbLxGzNuZeauxnISLU10bqR7JdbNj37qBJ0Dm2gpoLvzmhpbaIY28dOxtW7EyHkGJqfriiTcoxS1lV8flRMPNpVGotc8FH5REbeHeycosU2Ku+uF9o7yyr4g8g9tVWfowSKlVuUEKXqyz0w96SgzmrEv9mJIbKCGMTYqdPstT9EntmBszKx5DxmlJDg1tS+XKeay5riNvwso22l+Us+AZ2/ecNgN3he1k89kQ7o8tc6tjv7SyjJLAhY2yws2fY56ITyxvYwdwifJR3tbHV0u0Sxu0XbalNHA8SS/LHeFY/YQ3zc0i9UlpKMut9Fg23Egop8Cree2u8xCituF58zAY5RmlhOmcwbwbWhaZFHfPH/XB3BjI7VvjxtBrfkkODW0rlgsP1NpvPSzneUP7OGU5te85beA1Nxb3OLatexp1KGziXnerY3+iLAuI8DSrLaNDDNZgbYP1Nzu++CarOnhywrK3G2VPsrYmh7EXwLXlwatWh9LDxwW1vKOkqY/JMtxacSKOCRxtXDCysnkaD2QnbzSWxCfKOljhnVs2ybOba+8cFtl4mxsT3byPwo1zvMnzluTQyLZw1Qr3WW154E+KM65jox628/5OtY9sip7TJrruGBeFXdzrbnHsT5TldCEJDEau+opW6dxoyW+sd+m8o27BQvHiB4RrcVHbXcoyF1d1nK5tc0vnq4/wreW2oGMeFFquCTrmQUyi31mRtR9uLNo2ynPK83YdS/vUuDPejHNer47NjQlv4jm2R9xfkkMj2uIhPyij2ob5X1sU1FXvRUcc/y31SePr+WyDc1d1NQx98YiCwQXJX0s3qwbouAkqw4X8qfIfNokH3VFfeThgya6TVTyW6LLaRyFm0uwaC6z0/3Tl12RUnSiVKP/VOmvnrJkWx4A/Lvk8ykGfX0n+1l9kqxy8fKT4ayzHV+T+i/QIxhE35sY2R3VPc8ySHBraFqOqsmFu1y5zH295WaqSlX1eaeswMQJb4OxQfqic73UTj/+56rqUZVwxuEhREKtrr2qgeC3Excs/l/UuUaVyuw3qKz6RJ8t/LdWhiPXLM2OBMv238gc/uKjOXSnLYC2ZsaqkpcF+Vn9PllCLY4VyXK7tTBUOB0XA3NjewO5tjlmSQ0PaykdU5XmDtgvDRi733va3wtkh/FAZ9BLf6xYi2f2OdsJfD+virCrK8ZykqOx+fdYODEJyJC3v9VN/zxWfKw9FvWVBzRtS3ivJ+ZcifpvNH3bkZY6wH8eg88FM/eehAYX6sBgcYRzn6IO5MQeqt1XnkhzqaytHXvMab3rP3Yvz4t4/AAJ9/PC9bvlBvtfRJNa7PoUDpYSQtndHG/uFVIpfKn7LNhdPx08UcTcJQfsfK/6riNKVAhbYPixCWZ33oWJoR9sXikxyFwXOVaQ+8KUP7ygiXxfeyMi/tt1yABv+tMOvKG+ZBfW+mxt1XJw6HIG1OMSfFw26Bw3vikvuHIG1uLhz2C4X/155qpSxpODhYnEupHJrWlzPyZfynmqi4fUV/XmaEuMW5TKXHyULRflxzGeDlb1RqLP01q5ww22Fv8bFLYX28Cf6pVVo4IHq4lz8idMEyQd8/J00rgUo4VwoZeCVTC53mX/4Y+FT/Sv1w3fcHexFwNzohcgFehBYi0Nqt/NtWo/Izj4oAmtx8aBwDurWibKssxqf0J4a+AthQp9SfVdqhd+oVKYPvVB6c8UYiehrI78IiHJKzMuhmL5U7Axqhw8AWRkk/wCQOrAIJzw7z88zVB55UI6Tz21p5Uc+FPMyIGNuES/zfWwEjIARMAJGwAgYASMwEoH7lfJD/JU5DfcEViNI1s9QlZQ9lj1ZNUimN6MAfBSXlM5PlfYsCSY5k2W8Jb/yseDmyvIjpfW91kd5LX2Kk/V3sAIrmdLat7kl4V3Jkx9TX5Jdu01A5tRmk6g6ka2msFPPI+V/1hR+vQMGLZcOlVt9XF+L5725EdD4p2soNKXxhy+8KRnK51cq+7nqaflamkcBzt3/5PzQmHqO2f2IHr8D5uzxx7jWw3zcyb/0XlZTlktL5kn7agyljlAqiHwI17rJ3hVb51eycMMGHPrETT4pziSThrKfK8akh3PYGRKyun8oyieFtqUsFGXKQxTUUnkvy+BqMbhO9e9kjKhQcvOwM3jpuC2NK/I7LItAvJbeu7ZV8+haBLd3vueY7Y2JJTqPgDl7Hp8j5156L2u5YUiBSgpe45pQgqYyKJ1YEli6K7d4lkW3dBx8kyVvrhhjQW8ppuobltnktpHk73NvSHWXSjYW29JKners2qIIN/7RUR58kUOIx4zR85iUbx7pIO9fnud9I2AEjIARMAJGwAgYgQsQKC3LWFsJLSXyLqllvkYJDO4aKW+trRTInwbIghJcWmNJyy3NdIGP80pLLAooimhXALMWXpIJhZY41hJHWzyMpFDWzYd/30nGUm7K059SYSfdwQgYgRUQ0DyAUYHrmWuTa9t/jiMQHKZFwDybFk/XdodA5NU/OmJ5YHQg3OpuVse4dwdL84s1FFeKUrFMrgxYPX9T/lglsGlgyh0NJq4ELK+WK5i1JrCUc8MKQeVRQhl8FNoQlMZ+zaIOFtX+xnapN68bWfirXRTvBkfKKv6n2FiOVaYMKOr5R4L4Kwdrsc5DNh5S8vz8fGRsKe15pveNgBFYDgFdryjK33K9KmJYYL65aHWc5aR2S3tDwDzb24jtQ17xCsPcn5q7vok6x0sdo9fcbLhPzwUMIKDwBeUxHpOVAgogShtLozUKYMrs26o+fJyp461YlptIUAL7zu3KV53U935Xfp7OgKv8u4oQ4JUiT0tvK36vNNZeRi6IUXMrQUnlvFpA6SZ8Eethn6ewE5xU9yuVoc8ozQ84pnAelPaz8j5RpD2UatZaRr70l83nsMeF43NFByNgBNZHgAf5fN7g4+LfdT3zncBVc9/6XbMEG0LAPNvQYBxMlPwP59BH0rdqB+vmsO4EZVmTN363swTdHJKFJdwgdIySy1rBX6McXtEoA8cNqPVHI131qa2aRba338goWXmgqN3ksBixti996/0rUpVDYWdyw8p0oiwrLVn1Ub6DtVrn1BR4ijYhlkUBP6dMN+W9YwSMwCIIcA0nxThd73naIkK4kcMjkHPKPDv8cM/fQekSpb5082+ug7I8F/RS4lBo+Qvmxg9Y+1hYOeaV5MNL2tb5WMDTTeiSKsaeg8JfukhQR2s5uoGVvi8Maj7H5elYrYc+TEDs3IpV1lU7pm6UdgcjYAQmRkDXeG6Vofbw1k7b3yZuasvVeY6ZeXTMs8kBNmcLSKVvYeBE1/mgyLqpw3sz95anEZ56W0EXOBZQ3BEuNet/NlDhbLV76YHawrrL3003fYn7kGioQot1mvK4gAwJjb/yucKxTh5Iei3beT0qX1s2Ly/ifSNgBKZDIDxwc91NV+W2a/Ics8r43BzPpkTZnG2jKf2CN/cYNp8Jm5t+cz2rsixw+bily3qMZRhlelTQ4KFgj7Wijmqjo3BpvUVx5oO7MQT6VOWHKrXNknEd8qRk/M173UlSYW+NgBFYFgHNWSgwLLU59NpfVkC3dggEzLNDDOOmOsGcpYiehgtpbQGETck7pzBvCog56w91C2ReQfIBWvOPekqjYW4gg5eg0zlYZp8wgNomS+2/2n2otNktNmofBRkrbq8fMfLNHSQPDw5guKRLytzdcv1G4DAIxGv0PV2jpQ/gYfrojqyPgHm2/hgcWYKo+2DAY/GCIW6kh4Njbp9llOTvFbEE/6DI3ytjZXmuSKgqeSqDvy7WUgbnfcX0CoD0tKqFdpvwtc75fW4lNiqlm1CU6fnc/W3QveGdM1y8YVTcdRDo40bMx50qKMrxGIPBmLdRBvvGETDPbpwAM3W/i1dKxyj5t+IHlbnqZr9zms2yHAcCF4ETwJWHOR/Fl0X6W68mlYcyzcoO+Q3mhY5P3DlUFksvCvUilmW143BDCIzh4g3B4q4KgT5uKD8ZCp5lgDGnYZmZ/S1Y1qZ3d4yAebbjwduw6AN4hV71keaqYNBUeXyXnyq+favz1yyWZQGLEotCjDJcs6KQT2jl6TzcCj7WOfmX5AwWHwPi/tB8TBfLJl9d1kt+rvybfD0AkA7TIjCGi9O27Nq2jsBAbvyifmChaS3krznKivLWB3gj8plnGxmIg4kxkFfoVixhy4IEbymis+FOdrPz1yzKskDFOkzoclnAlM9qDI3yG0rfndcsMxfTkmLNjacJOpe6u+pvynnHCFyIABwexMUL6/dp+0Wglxuan07ehO23u5Z8JQTMs5WAP3izQ3iFIbNlzDw4Jr3du9db4rIC/PEGf9Zx8hSiJxWUXmJL0VU6bhmk49ucB15nEjxwdzj4d2YEzMWZAd5x9ebGjgdvR6KbZzsarB2Jal5dPliTK8saDBReQtfi+2lt5dyXj/KY/VmK7RUHWfhM+6RXPwbMynnXCEyFgLk4FZLHq8fcON6YbrFH5tkWR2X/MplXF47hXG4YiIODeC3wkQtr95VKMZbllluGFG+sysTR6zHrHAcjcCkC5uKlyB3/PHPj+GO8hR6aZ1sYhePJYF5dOKaTW5ajEozS+34pk5RfPnbBStzyBY3WaHyTiSHENMrzxyZ2wYi4eDMvAubivPjuuXZzY8+jtx/ZzbP9jNWeJDWvrhut+9ed3nk21uOfNDi4XIT1lbVFQWaNUV4DlIGnHQJfX7JECYEVMVhmyYpygMM/CyFgLi4E9A6bMTd2OGg7FNk82+Gg7UBk8+qKQZptnWVkkuLLh35YmHHJ6PynOZXjT0seSzG2u4WAcFgPAXNxPey33rK5sfUROoZ85tkxxnFrvTCvrhuR+9edfv5sKb+sezxk7eNPVa784O985c41AvMgYC7Og+sRajU3jjCK2++Debb9MdqjhObVFaM2uc/yWFn0tIOfMitotD7uG1uPyxuBaxEwF69F8LjnmxvHHdst9cw829JoHEcW8+r6sVxdWVYXUJb56M++ydePp2u4DgFz8Tr8jny2uXHk0d1O38yz7YzFkSQxr64czVl9lq+UzacbASNgBIyAETACRsAIGIFVEZjVZ3nVnhWN6zUE6zU/VsTlg5U2vrU1Wyg4rIJA5CN/uEOAmyyR6D/eCXD4ZywC4tNPOofVg8r168dW5fJGwAgYASNQIHAzyrL6zf+hf62byc+6sbBKx/eKXn1DIDisggBc/IiWxUeW9EHZ4SHOwQiMQkD8YYlODAEORsAIGAEjMAMCt6Qs51YXlrPzB4UzEMpVDkbgUVaSv4Z/R0rPA1sGM1S824tAfNCCPy97C7uAETACRsAIXITAzSjLKCG6sfC6m1ff/0tWvYtQ80lG4EoExL/8rQaWZT5y9Sv0K3G9pdN5uFJ/3xFveFt2S113X42AETACiyJwKGVZNwy++EQhfsUNpERSaay48YfKfaj4u45zhaUs7mMjcDUCfZxUPgrPU8XgknF1g67gUAj08OeJ5rBvDtVhd8YIGAEjsEEE7m1QpotE0k3lS52IlYU/QXmp4xdREQn1xfywrzIo0v9TGsq1gxGYBYHIuXOcRFHGd/4DcdIf980yCvut9Bx/lIdR4C1tvyZqH7eepzF9v5225EbACBiBDSJwiKXj4g3ihRSO5gMppfFqmxUGPtE+Ssnf2n/IGMTyv6Rj0hyMwJQIDORk+OhU7eJv+lR8/GpKGVzXfhHo40/ZM5X/V2lvi0OvyjwfGwEjYASMwHUIHEVZRhl+rBtF43qhmwdp/yotOPPpmC/GUab/UeQDv2fK8x+hCAiH6RGI/OvkpPL/VKv5mw18lpuHveklco17QqCPP6kvsRxuPLxZ+06RVVb8liIB5K0RMAJGYAIEDqEsg4NuGigeWJMJ3CyIKCTv7vHmkfoj2bkBLhpou8RMaTxs/FymLyrYzhpLYxjF3iwnk5zm2rYIlsZl6/xBviSrObQtDlkaI2AEpoqfjvcAABLHSURBVEHg5AO/OOml18HJ8vWVJsGqFVblsWhgqWXtYgLW3fQqEOsugWPqmMXikWRQ/Z/QWOzDF+wr8Ip7VyHKz5+mrPXRF9glDgTsJMt3kgs/8NnGccpBkpz4dOKGk3jQql75s/I21b91TkpOrnFzrcWOYQdpjFV68rkv1b11/oCUOTSMLy5lBIzAjhHQZPxGiuoGyi2+v3kafpX/KX6Yp5f7sQwrTDTnpn3lYZWkjo9T2lRb6lT8s6xPab8r4oZxIs/W06LsWHdXkV3t8yr3pG2lw4/qGNfKr5kmOXmr8FOfDCoDL6t9UvpFvNV5u+EkfVc01yp87+NOyhd+k3JoT/wBAwVz6Ar+JB55e3rPMSbGZCscuMdMlwW+zP88O0ZhwsKIZfhFnp7vy7KAFY/w/G7T/lUduBIQqX/qgDLPChi1sLgLQ02IMWnCEgXtN2E2ixW+Txa1j0KMX/dJkEzwgDVdscpuNkT5WB3gbFC5uXi7C06q/+baWYb0Z87EoV3wB3TMoX6OuIQRMAL7R6BUlvH5xUpQBlwrHmQ3hjI/+Qo3H9iVBXSMsk0d3KAnCaoLxY7XyPxVcBOUThqK0LMmcT87PJxws1wrPFbD58YRTPmgaJMhjj1K/ZCHjcl5uzNOmmvXs3hSDu2MP6BnDl3PIddgBIzAxhEolWUUjPC3uyPlDr61sjxW/ZpjXUl5meyL/2jppPrSLxllE99alKbdhKjoYc1PWK0h+0fnxjFi+pdkTUrCGjKeaxM/5aFvFCbnbcQH+TbNSXPtHIVG5U3Kob3wB4TMoVE8cWEjYAR2jEBLWdZE/Z4iK2SUSmZ6Xd2lxKE4nbNGAhHWXkLa3h1d/8sHXCzGj9Wa+K2OWYZr8X+2iu1/iQyKLVcFHT9RbKz22v9Y8V9FrOMp4Ovah2Moq/P4F8LQjrZ8eJfGKNU1eMu5itTH2IAhD0zI1zVWyBg+phzcyAIFJS9vLRj/oWEu3s7OSfWVcTLXho70fOXm4NDs/AEOc2g+UrhmI2AEjoXA/b7uaELlZoDS9E1FiWbCTUpayxWiUm8q16VwV07pT5JMrNKA8vZpLL3mOqPhjyUkD0ov/tm5wo5ymfedBxIU5dztAStVp2+48kJQ/VjOH6jv3FQZA8boF8XwpyukDQ06F+WSD+F+1D7jzMduyIn1+HfFT+JxXiUfz4W288Q19yUnWIJJjnGnSCqf+Dg5byXDEpw01zpHd5mMuTi0EH8AyRxahipuxQgYgZ0jcG+A/ChTP2oCxzetFlDUCH0WUZaXI/QpJ3elRvyiIHGDiXGQsjSi+kFFdeNE0fw1FkbpLeUAp6bvkhW8iHk56ihf3yvpdVA7WE9ZVSRXVqkDS2Mai9cnnNlTeeRBOU4fSHJ+Po7s1/ynkTG3iOtw9cCNP3846RMoYZX3t3bORbyVLLNx0lyrDdMqabNxaE7+gJQ5tApf3KgRMAI7ReD+Obk1oaIo80cUuWJWnjLEZ49zsLa+Ul0nyonaCesPlRXv4Vj9Cf8QKFlfaj8pnVi5nyX51b9kxSz7jkKVK8uPdA4W53Mh+GMXBVCyCYMVWMmE0s2SYbl/L3/gkh9TX5Kd+lNA5tRmSlttq77ArecjBdgdbzU25trIQZ65+FUcWnrey/gDLJ6vZiaHqzcCRuA4CHQqy5rIg8+tJthzijJIYF0plcAWQlExI61qnS4m8da5ax9IdpTFz6Ic7PMBWa7ghiylBSVX5cEDJTMpzuSTxoNCeV44hwJDQlb3D0X5pND+UaSfO8Qt5Oy4KR8XkTF1nmtvljxhAtbvC9sc7yFtbY635tq2uVYh1VUcEmfTw0+l6uuTxCfeHOFGdTLPpLRsTsmvH89X18PvGoyAETgQAvdrfdEEiqUOK2OjKCstWBKV1ih8SktKWuNeUNanMigzWEOxUOdWy7Jocxzbou6qJbopOGBnTF0q+6QiIz7QwYKkfG4i9PXcih7BNznHSeU5v6WYRrmS24ayQ+hzb0h1lzc/lPnSSh2r7NygCDf+0VEefJFDiMeMQWMhj1lsHik2PMjS19hlTPggselLFCLwJ6b/qvFoXDSURh5hNd5KBnPtbgySS8AeuBYlDm4Mk3NInGCODbwVX1vzRdNwttPBoVCCPO1wjfeFNKfk17Pnqz7UnG8EjMBNIXCv7K0mWSZrlg9rFOVYBgW6DCgqhOrErrpQlPnwDEUuKJwUPhd0DhZt3AOwdLzUMSs9UM/ocEFdNSX4Udbwb9rvW1qPG15pjSUtvxlRJRbq3JpDGmXy9kjLA3i3sFYfGS8iN70xgbZyXMu6ccHBD7yUkTboT6mwk754QD5FrGetKEGQjz93Ib1RlKOAW+CtufaaLbvg2mtxw96kHLpgrkKIGofSB7/MVTx89wXPV30IOd8IGIGbR+BejoAmbCZOLHSvtP91jCxPRlrN/SC4J0gZKZXDNGGzVBoKy3va9ga1g9JHO0EhjPVyI/2+9+SiwFR1FbJzg0TxP6coYq0ExxAkB+egANO3EKJsNasmOFax0jkottSb101aGptmDCir+J9is1SdypUBl5j8gYg3CUGh13nIRj/z/Px8ZGwp7XnmRvbBpitskrfmWnW4tsq1yTik623KeS/ME+JSMx9UUX2d6PnqNRbeMwJGwAhUEWBN5SYjKleNUtdk3O2gPAVLhsqhoKG0pbKl9ZEJG8Xr2xGTNgo25z3WOY0iFtP+Vdoo/76xdcXyn6qdqqtIzMdKjqUyKJXarwaVRcEnvFL8R5E6v1f8VfEtxT9r7eg8FGswO7EYKQ/LPriTl6z87FcxVvnkUsHa2chxElQm3KSVgVKNZZobZ2edygtB51H+c9U79IYcz5x/I9lw+aFfyfIHN1kaj+XcNsFbyQHPzbUdcm0ODkU+jJr3ujik9C/F9fAmRfvMAZ3Xv/JCUDnPVwkMb42AETACFQRaynIlf/EkTdwo4UnRQSklMuk3lk/tJ8X6qXZRRt9XfFYqb0PqUhl8+1AQaZe2UGRaCqjKkIayi4JYVTyVN0lQW/QVN5iWQq50bmjcUKuW51rjOgel+o+yrrKsygXMVa76oJCXj2VRPsHMISKQMIyHVd6qjLlmrlWvmSH84cRzHFIeD4nB4h0bgW9c0891vc7yYKs2PV9FsL0xAkbgwAhoEuXvlTcRBTP+yi+SPNpHicNSiPmbP5xo5NQx1s3/kaZAOazPeX5vXToHSyofWmnzBh/ysUUx5gYQ2ovHKKq0QV4ol7c15b7q5waHst70Jcr1r9KxGrXSzx0PlTW2iZ94b93UOVaOIfXuuQx4KJ7lrfLNtbvrC36ba9m1NoQ/XB9DOJRfRyrPnNGaN/P8KfZVv+erbCynwNR19N+HjJExWpoD9zTZbSLIQoEVFH/l5kM17WOhw8rMqhiNRVdlUU5wC0nWkpc6xk8XZRbrS29dsQ7qbVlTYzupXapDKeeGgALNzSe5QGh3+hDl4SNClPMQ4j59a9xTYlbnRudQHqv7kNCy2nedEOv8UDKWH8x1nXL4dGFirplrF/N8CH+oXOWY84bMV+Gtm8rzUMscwLcnzVyi40mD56tJ4XRlRsAIbBSBzSjLwofJvfR9TrC1FFol8uHZ85Sp7WP2o6LL7pC6qOMZhXUzKRVLbi7BDUJ1okjirpLiEu4HyJb8CBExyCMZ0sMBaX0Bn9ihSi0PAkMCPr/Nw8yQE26gjLl2539trl1G9iH8oeZB8xUFdd2jVH+lyJxV+zCbYlMGz1dToum6jIAR2BwCm/BZjsoqVlt8dRvrabSIoMg9VPor0MvK5mkolo9U5pMs/2xdKkd7fPzyl/axXvOPVn9oH78/XqkvoRSrqXqQHCjIWHHLB4X6CTOnSh6s66yVHR4iZm5uF9ULEx6yenmrcubaiBG9Fa4N5Q/QmUMjCHSHl+ercZC5tBEwAmcQ2IRlWQpYUIQlJ+4UecDqgoUk5ZMXPv5LabqJpJUPPiczpWu3r64fUl3aolijKKOg8iHfoDWhVW62IHlwM9mEokwnkQWZZuvwDisWHomX5tqE43crXBvBH9D1fDWCY7fCoRGQuKgRMAJXILAJyzLyS1HFEoCSGpRebVGUw+tEbZugcsGKrASWYSNgAW4p1CPqog38emmXD7BY1o1VNZISpEOHORHQWOHzywdfQ1/jzynO6LrNtdGQrXbCFrk2lD+AprKer1Zjz13DW+TQypC4eSNwEwhsRlkGbU1EwfUgIl995a8yuGWgHHf5N4fTh9QV26HdZm3SlObtMggI+10ry6Bkri3DlWtb2SrXxvAn8s3z1bVkuPD8rXLowu74NCNgBAYicH9guUWKSQHmNX+n64EmqgfKR6Fu/Jq7BOurK52nOnHjGPPhXDrVWyMQEDDXTIRrEBjKH9rwfHUN0j7XCBgBI3AZApvwWR4hOooyf7IxpZtE65+zRsjiokZgLALm2ljEXL5EwBwqEdnZsR54WBqUJf2eKH6pyLKnbMP3ODvrjsU1AjeBwKYsy32IS0nGAjz4H+z66iNfdXZasoec7zJGYCgC5tpQpFyuCwFzqAuZfaRHhZglOMNKTEitNI5xR3vIsYMRMALbQ2BvluXtIWiJjIARMAJGwAj0ICClGDdCFGO+uSlXFpr6jWmPNM42AkZgDALNB366kMP/qY452WWNQEJAk/+bab+2Fb/4kr/2mpEbyCPF8uZBNSyfd/ZPUMxbYLqtYK7d1njP0ds1OKS5ipWcnpRtK5112Fma86s5+uo6jYARuB6BRlm+virXYATGI6Abxe5Xwxjfa5+xBgLm2hqoH6vNazgUleLfpBQ36/grje9wWOGp9Sdax0LNvTEC+0fAbhj7H0P3wAgYASNgBDaMgJRi3qARy5WXwts2KdC9KzxtuHsWzQgcHgEryxcMsSY+vmTm6+VvFV8oMgk6GIHVERAX+XMdByMwKQLm1WRwpj/TShXiZlYq0CnPWyNgBDaCgJXlkQOhmwZ+Z3/KEvCN4hfaf6nIRxsORmBVBMRNXFo+1NYPb6uOxLEaN6+uH0/dK16pFpRi3C5CEK78ay2WZVuV7yDxrxHYLAK7WjpuQyjyF9sp/K4dJj0HI7AaAlFBfn81AdzwIREwryYd1g9U2/fC9MtYa/qo2W+DJoXZlRmB6RGwsjwS02hNzs9i3WdbBnJEvL8GAjywPVNMN+I1ZHCbx0PAvJpoTKN1uVndJ1qWWevf94+JMHY1RmAuBKwsX4FstLp8qiqwGDhchgA3CpaOc7gQAfGQv2xPVqoLa7mJ08y1EcNsXlXBmpJDrIphf+UqzE40AttCwD7LF46HbiRY8H5RfCbLgCe8C3HE2qJoRe9C/OJpnwnDH6+r4vhnm2ujx9i8KiCbmEP2Vy7w9aER2CoCVpYvHBlNmnzghwvGu1Kc7XN2IY4+7ToExD1ek/PRqYMRmAwB82oyKE8qEraspsS3LnyI+7H2ff2eoOQEI7AtBOyGcf148M90f2rC+1jKs6171+PpGgYiIM5xs30g3tkyPxAzF+tHwLzqx+iaErpe/U991wDoc43ACghYWa6ArpsFr8f4EIN/VmKFgeBqEW8if+v4A014peuF/W4FjMP0CHTxUS3BU95slJapYLkSR7+bXhrXeBQEzKujjKT7YQSMwNwI+O+uC4R1A8FSjLWONZTfiDeUFzp+GI/TX5MGa57y8V1+qvi2yrCWpoMRmAyBPj7mDaksa7jCz4fmYo6M90sEzKsSER8bASNgBLoRsGU5w0Y3EPw/cafI11FGKX6gvA+VzpfQWJy/0PE/2r6liILynpUToeAwKQID+RjajGXTslSs5fpcnLRb0KQjcozKzKtjjKN7YQSMwHIIWFluY41VufQnQxkm4B/Kmpi4X5QuGGQ5GIGpEejlY2pQvMTlwm4XCRBvzyFgXp1Dx3lGwAgYgQIBr4YRAZG1Bf9PFOIfCoxYw5ZgBfkOB/8ugID5uADIN9iEeXWDg+4uGwEjcDUCVpZfQ8gr7L8q7hSfxXSvOPAaK+/Nj4D5OD/Gt9iCeXWLo+4+GwEjcBUCdsN4DR+WZXySmyArDFZlIuspOxiBJREwH5dE+3baMq9uZ6zdUyNgBCZCwJZlASmlGPcLfJOJIcS0Fzr4IvopxxxvjMC8CJiP8+J7q7WbV7c68u63ETAC1yJgy/IdglhbCKxywVJwBFbE+MSKcsDCP8siYD4ui/ettGZe3cpIu59GwAhMioDXWRacUpD5U4fHUoztbjEpvVzZJQiYj5eg5nP6EDCv+hByvhEwAkagjoDdMO5w+VSb53WInGoEFkfAfFwc8pto0Ly6iWF2J42AEZgagZtXlmVtwU8Zn+XWx31TA+36jMAQBMzHISi5zFgEzKuxiLm8ETACRuA1AjevLAsKlGWWjPM6yq954b31EDAf18P+yC2bV0ceXffNCBiBWRH4P7ty98JYqM6tAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$\\displaystyle D \\frac{\\partial^{4}}{\\partial \\rho^{4}} w{\\left(\\rho,\\phi \\right)} + \\frac{2 D \\frac{\\partial^{3}}{\\partial \\rho^{3}} w{\\left(\\rho,\\phi \\right)}}{\\rho} - \\frac{D \\frac{\\partial^{2}}{\\partial \\rho^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{2}} + \\frac{2 D \\frac{\\partial^{4}}{\\partial \\rho^{2}\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{2}} + \\frac{D \\frac{\\partial}{\\partial \\rho} w{\\left(\\rho,\\phi \\right)}}{\\rho^{3}} - \\frac{2 D \\frac{\\partial^{3}}{\\partial \\rho\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{3}} + \\frac{4 D \\frac{\\partial^{2}}{\\partial \\phi^{2}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{4}} + \\frac{D \\frac{\\partial^{4}}{\\partial \\phi^{4}} w{\\left(\\rho,\\phi \\right)}}{\\rho^{4}} + q$"
      ],
      "text/plain": [
       "                       3                2                    4                \n",
       "                      ∂                ∂                    ∂                 \n",
       "                 2⋅D⋅───(w(ρ, φ))   D⋅───(w(ρ, φ))   2⋅D⋅───────(w(ρ, φ))     \n",
       "    4                  3                2                  2   2            D⋅\n",
       "   ∂                 ∂ρ               ∂ρ                 ∂ρ  ∂φ               \n",
       "D⋅───(w(ρ, φ)) + ──────────────── - ────────────── + ──────────────────── + ──\n",
       "    4                   ρ                  2                   2              \n",
       "  ∂ρ                                      ρ                   ρ               \n",
       "\n",
       "                     3                    2                4             \n",
       "                    ∂                    ∂                ∂              \n",
       "∂             2⋅D⋅──────(w(ρ, φ))   4⋅D⋅───(w(ρ, φ))   D⋅───(w(ρ, φ))    \n",
       "──(w(ρ, φ))            2                  2                4             \n",
       "∂ρ                ∂ρ ∂φ                 ∂φ               ∂φ              \n",
       "─────────── - ─────────────────── + ──────────────── + ────────────── + q\n",
       "     3                  3                   4                 4          \n",
       "    ρ                  ρ                   ρ                 ρ           "
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "display(expand(D*L2 + q))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "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.6.9"
  },
  "varInspector": {
   "cols": {
    "lenName": 16,
    "lenType": 16,
    "lenVar": 40
   },
   "kernels_config": {
    "python": {
     "delete_cmd_postfix": "",
     "delete_cmd_prefix": "del ",
     "library": "var_list.py",
     "varRefreshCmd": "print(var_dic_list())"
    },
    "r": {
     "delete_cmd_postfix": ") ",
     "delete_cmd_prefix": "rm(",
     "library": "var_list.r",
     "varRefreshCmd": "cat(var_dic_list()) "
    }
   },
   "types_to_exclude": [
    "module",
    "function",
    "builtin_function_or_method",
    "instance",
    "_Feature"
   ],
   "window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}