{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "6d240197-fd38-4a27-9c6f-87c78d4f0901",
   "metadata": {},
   "source": [
    "# Complex numbers arising from cubics"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "36f80441-2092-4f8c-8c24-deaeb6dadf26",
   "metadata": {},
   "source": [
    "In this notebook, we're going to use a little computer algebra to talk about how imaginary numbers arise in the context of cubics. It might seem that complex numbers arise more naturally in the solutions of quadratics like\n",
    "$$x^2 + 1 = 0$$\n",
    "This is not how complex numbers arose historcially, though, simply because it's easy to assume that such equations are meaningless. As it turns out, though, complex numbers arise in solution of the general cubic - even when the solutions are manifestly real. This played an important role in their acceptance."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ef314ada-9b3a-4fb1-a2de-b07b7f3eba4d",
   "metadata": {},
   "source": [
    "## A surprising cubic"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "55394eed-375e-4cbf-abc0-24365d21a8b5",
   "metadata": {},
   "source": [
    "Let's ask SageMath to solve \n",
    "$$x^3 - 3x - 1 = 0:$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "id": "05207390-cb23-4f7a-8a6e-7c52e9d64e68",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = -\\frac{1}{2} \\, {\\left(i \\, \\sqrt{3} + 1\\right)} {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} - \\frac{-i \\, \\sqrt{3} + 1}{2 \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = -\\frac{1}{2} \\, {\\left(i \\, \\sqrt{3} + 1\\right)} {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} - \\frac{-i \\, \\sqrt{3} + 1}{2 \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}$"
      ],
      "text/plain": [
       "x == -1/2*(I*sqrt(3) + 1)*(1/2*I*sqrt(3) + 1/2)^(1/3) - 1/2*(-I*sqrt(3) + 1)/(1/2*I*sqrt(3) + 1/2)^(1/3)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = -\\frac{1}{2} \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} {\\left(-i \\, \\sqrt{3} + 1\\right)} - \\frac{i \\, \\sqrt{3} + 1}{2 \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = -\\frac{1}{2} \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} {\\left(-i \\, \\sqrt{3} + 1\\right)} - \\frac{i \\, \\sqrt{3} + 1}{2 \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}$"
      ],
      "text/plain": [
       "x == -1/2*(1/2*I*sqrt(3) + 1/2)^(1/3)*(-I*sqrt(3) + 1) - 1/2*(I*sqrt(3) + 1)/(1/2*I*sqrt(3) + 1/2)^(1/3)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} + \\frac{1}{{\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} + \\frac{1}{{\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}$"
      ],
      "text/plain": [
       "x == (1/2*I*sqrt(3) + 1/2)^(1/3) + 1/(1/2*I*sqrt(3) + 1/2)^(1/3)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "surprising_solutions = solve(x^3 - 3*x - 1 == 0, x)\n",
    "for sol in surprising_solutions:\n",
    "    pretty_print(sol)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e51271a6-40ca-4af7-86dd-8401ff9a4af8",
   "metadata": {},
   "source": [
    "Note that all three roots are expressed in terms of the imaginary unit $i$. Of course, every cubic with real coefficients has at least one real root. In fact, a simple look at a plot shows that all three roots of this polynomial are real:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "id": "2cd43856-4698-488c-912d-7bf0a0bfe04f",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAnIAAAHUCAYAAAC+g8X7AAAAOnRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjEwLjUsIGh0dHBzOi8vbWF0cGxvdGxpYi5vcmcvWftoOwAAAAlwSFlzAAAPYQAAD2EBqD+naQAAWdlJREFUeJzt3Qd0lNXWxvGH0BEIKtKRJlVAEaUpl6YIIldsn1hQFAtNBcWCFa5XEVQQBbti74q9XwG5CgpKFamCIoIUkS51vrXnGC5gEpIwM2/7/9Z6VyAwmZMzk8yec87eO18sFosJAAAAgZPm9QAAAACQNwRyAAAAAUUgBwAAEFAEcgAAAAFFIAcAABBQBHIAAAABRSAHAAAQUARyACLBSmauX78+/hEAwoJADkAkbNiwQenp6fGPABAWBHIAAAABRSAHAAAQUARyAAAAAUUgBwAAEFAEcgBS7osvvlDnzp1VoUIF5cuXT2+99dZ+bzNhwgQ1btxYRYoUUfXq1fXII4+kZKwA4GcEcgBSbtOmTTrqqKM0atSoHP3/xYsX65RTTlHLli01bdo03XTTTbrqqqv0xhtvJH2sAOBn+WIUVQLg5S+hfPk0duxYdenSJcv/c8MNN+idd97RDz/8sPtzPXv21IwZMzRp0qQc3Y/VkLPyI+vWrVPJkiUTMnYA8FoBrwcA7GnXLmnJEmnBAnf9/LO0bp29CLt/K1xYKlJEqlhRqlxZOvxw6cgjpfLlmccws2Ctffv2e33u5JNP1pNPPqnt27erYMGCf7vN1q1b49eegRwAeOW666SaNaXLL0/s1yWQg+d++kl6803p88+lL7+U1q51n7fX5ipVpFKlpBIlpPz57cVZ2rxZevddacWK/32NChWkJk2kdu2kjh2lGjU8+3aQBCtWrFDZsmX3+pz9fceOHVq9erXKZxLJDxkyRIMHD+bxAOALr7wiXXBB4r8ugRw8sXKlNGaM9Prr0tSpbqWtZUvp6qulpk2l2rXdapsFb1mxoO6XX6QZM6QpU2zVRurfX7rySqlOHen8890PTdWqqfzOkMwt2D1lnArZ9/MZBg4cqGuuuWavFbnKtowLACm2aZO0dKl7bUs0AjmklB1xGj5ceu45KS1NOuUU6dprpU6d3KpbbljwZytvdp1xhvucdV/6z3/cCt/dd0u33upW6ez1vEMHd58InnLlysVX5fa0cuVKFShQQIceemimtylcuHD8AgCvLVzoPiYjkONlDSlhq2adO0v16knvvy8NGuRW02xFrmvX3AdxWbGvY2fmn33Wbb3aRzsaZYFi/frSU09JO3Yk5r6QOs2bN9enn3661+c++eQTHXvssZmejwMAP5k3z32sVSvxX5tADknfQr34YqlRI2n+fOmZZ1wyw403Soccktz7Ll5c6tZN+vpraeJE906oRw8X0I0da1tzyb1/ZG3jxo2aPn16/MooL2J//tmyW/7aFr3wwgv3ylD96aef4lullrn61FNPxRMdBgwYwDQDCEQgV7p0cl73COSQFJZh+sQT7qzaO+9IVi5s9mzJXpsLFUrtpNsRqhNOcMHbd9+5BArbim3RQvrmm9SOBc7UqVPVqFGj+GUsQLM/33bbbfG/L1++fHdQZ6pVq6YPPvhA48eP19FHH6077rhDDzzwgM4880ymFIDv2UJGMrZVDXXkkHC2pXnRRbb1JXXvLt1zj3sn4ieffeZSwW3Lt08f6c47JUqLhRt15AB4xaoqZBzvSTRW5JBQdv6tYUMXIH30kctM9VsQZ0480WW63nefG6Od3XvvPa9HBQAIm1jMba0ma0WOQA4JsXOnVd+XTj3VvfOYOdMKtvp7cgsUcOVK5syRjj7aJWP07Stt2eL1yAAAYfHbby7pjkAOvmWdF/75T+nee90KlxXrLVNGgWH16mzMo0e7c31Wx86COwAAEnE+zhDIwZfsPHrz5tJXX0kffujqtWVRn9XXbMy9e7vtVkvUaNxYevxxr0cFAAi6efNcDdPq1ZPz9dlaxQEV9z3+eLcVaSU+9mmFGUgNGrhgzpI1rB+eBXfbt3s9KhyI0aNHq169ejruuOOYSACeBHLVqrki9slA1iryxIId62lqLS4//tj1Og0bW5GzjFYLVl97zZ9JG8g5slYBeMGOHtk5cksGTAZW5JBr1tO0bVtXoXrChHAGceayy1y7r++/dwkcnJsDAORlRS4ZHR0yEMghV6ygrq3EWR1X65iU7O4MXmvZ0q0+WpcIKypsQSwAADlhR3N+/DF5iQ6GQA45Zp0Z7BycPSGt5tpBB0Vj8qwTxBdfuGKO7dq5pA4AAPZn8WLX35tADp5buNAV0a1c2RX6jVoXhFKl3FnAk05y5x2ef97rEQEAgrCtathahadWr3bbqenpru3WwQdH8wEpWlR64w2pWzd3Pfyw1yMCAPg9kLOjOck8S14geV8aYfDnn1KXLq7o7+TJ0mGHKdKsG8STT7qg1kqTWM05y2wFACCrRIdk1lclkEO2/eGsnpolOIwbl7xihkFjP5DDh7sCj9bSyxDMAQAy6+qQzG1VQyCHLA0ZIr36qttOtLZV2DuYs5ZkxoI5C+p69WKGAAB7r8i1aaOkIpBDpj74QLrlFum226QzzmCSsgvmMrZXLSHi3HOZKwCA9Mcf0m+/SXXqJHc2COSQaYbqeedJp54q3X47E7S/YO6++6S1a6ULL3SJIB06MGd+a9Fl104rrQ4AKTJ3rvuY7ECOFl3Yy9atUosW1s5ImjrVHerH/lmdoDPPdEWSP/vMzSH8hRZdAFLp6aeliy+WNm2SihVL3v1QEBh7ufFGadYs6ZVXCOJym8368suS9WXv1MnNIQAg2ityVaokN4gzBHLYzbo13H+/dM890jHHMDF5qTP3zjtS1arSySdLP//MHAJAlAO5OkneVjUEcohbtkzq3t2di7vqKiYlr2wr2jpfFCrkOkBs3MhcAkAUzSWQQ6pY1qV1KihcWBozJrmFC6OgbFm3ummNks8/X+KMPQBEy/bt0qJFrMghRR56yBX8ffZZqXRppj0R6td3Z+YsoBs4kDkFgChZtMglwdWtm/z7Yms14mzV6IYbXDHbdu28Hk24nHKKK01iZw5tpRMAEA1zU1R6xFBHLuJbqj16uP6pQ4d6PZpwuvpq6YcfpCuucC3OWrXyekQAgFQEclYkvkyZpN8VgVyUPfKINH689J//SCVKeD2acLLzhqNGuSLL1iFjyhR61gJAVBId8qXgzDlbqxG1eLF0/fVSz55S27ZejybcChaUXnvNvTs7+2zpzz+9HlG0WFeHevXq6Tgr8gcAKWA7ManYVjV0doigWMyd35ozR5o9m9W4VJk2TWreXLroIunRR1N2t/gLnR0ApOo11t6433STO4OebKzIRdDYsa7W2YMPEsSlUqNGbpv1scdchjAAIHxWrHBtLlO1IkcgFzFWoNYO4FvhXytYi9Sy5BIrvGxb2rTxAoDwmZvCjFVDIBcxd9whrV4tPfCA1yOJJjv4Onq0VLOmdOaZ7l0bACBcgVyBAqlLbCOQixA7Ezd8uHTzzVK1al6PJrqsgfLrr7vld1uhs/MUAIDwBHJHHOES3VKBQC4iLFjo29cFcNdd5/VoYCtyTz/tAjpWRwEgPOamqMdqBgK5iHj7bdeGa+RI11MV3rO6ctdcIw0Y4OrLAQCCb26KAznKj0TAtm2u96etxn38sdejwb6PzfHHS3/84cqTFC/O/CQL5UcAJNumTe73+DPPSBdeqJRgRS4iHRysge+993o9EuyrUCHpxRel5culq65ifgAgyObPdx/ZWkXCrF0rDR7sDtU3aMDE+vW8nNWXGzNGeuUVr0cDADjQ0iO1aytlWJELuX//223f/etfXo8E2bFuD+ecI11xhfTTT8xVItGiC0AqW3OVLy+lp6fsLjkjF2Y//uiWd2+7TbrlFq9Hg/2xc3JHHSXVqCF99pmUxtushOKMHIBk+7//c7VaP/9cKcNLRYgNGiSVLu0yI+F/1pvvqadcdvFDD3k9GgCA3zNWDYFciIv/Pv+8K/5rBWgRDO3aSX36uEbLCxZ4PRoAQE7t3OmSHQjkkBC33y4dfrh06aVMaNAMHerOWFhPVvvFAADwv59+krZuJZBDAlg9MusYYGfjKP4bPAcd5Lo+TJokjRjh9WgAALnJWGVFDgfs1ltdSYtUFSNE4p1wgjvbaEkqtk0OAPB/IFesmFSpUmrvlzNyIWOrOO+/72rHFSjg9WhwIO64w3XjsNIk27eHby4feughVatWTUWKFFHjxo01ceLELP/v+PHjlS9fvr9dczPeAgOAx+xNt9WPS3XFAQK5kLHtVCv8azXJEGxFi7o2L7ZVfvfdCpVXXnlF/fr1080336xp06apZcuW6tixo37++edsbzdv3jwtX75891XTlp4BwCeB3JFHpv5+CeRCthpn9ccsmKMGWTg0aSLdeKMr6DxzpkJj+PDh6tGjhy699FLVrVtX999/vypXrqyHH34429uVKVNG5cqV233lz58/ZWMGgKzEYgRySIA775Tq1pXOOIPpDBMLzGvVki6/PBxZrNu2bdO3336r9u3b7/V5+/tXX32V7W0bNWqk8uXLq127dhpnBfeysXXr1ngR4D0vAEgG65e9bp1Ur55SjhW5kLDtNzsbd9NNrMaFTaFC0mOPSV9/HY5CwatXr9bOnTtVtmzZvT5vf1+xYkWmt7Hg7bHHHtMbb7yhN998U7Vr144Hc1988UWW9zNkyBClp6fvvmzFDwCSISMpzYtAjuPwIXHXXVL16lLXrl6PBMlw/PFSr14uUO/SRQpDTGLJCnuKxWJ/+1wGC9zsytC8eXMtXbpU9957r/7xj39kepuBAwfqmj3amtiKHMEcgGQFclbuyxLUUo0VuZA06X3jDXeWikzV8BoyRCpZ0nV+sPMYQVW6dOn42bZ9V99Wrlz5t1W67DRr1kwLsml/UbhwYZUsWXKvCwCSFchZ/Tgvju0SyIWAZTRWrEjduLBLT5dGjZLefdcVfA6qQoUKxcuNfPrpp3t93v7eokWLHH8dy3a1LVcA8EMg58W2qmFrNeCsWsMLL0j33UcXhyg4/XR3XXmldOKJ0sEHK5Bsy7Nbt2469thj49ukdv7NSo/07Nlz97bosmXL9Oyzz8b/blmtVatW1ZFHHhlPlnj++efj5+XsAgAv2Q7J999bwpY3908gF3APPiiVKCH16OH1SJDKx9ze+d1wg0uCCKJzzjlHa9as0b/+9a94Pbj69evrgw8+UJUqVeL/bp/bs6acBW8DBgyIB3dFixaNB3Tvv/++TjnlFA+/CwCQVq6Ufv/duxW5fDE7YYxAsmoKdujdDsGHrWAssmfl1nr3to4HUqtWzFZOWLKDZa+uW7eO83IAEsYqIbVt686rp7rPquGMXIA9+aS0ebPbZkO0XHGFZMfJLIgPY/suAAjS+biCBaUaNby5fwK5gNqxQxo50pUbsUQHRIt17rBVufnz7fyY16MBgOia81ePVQvmvEAgF1B2xvunn6Rrr/V6JPBKw4ZS377S4MHSL7/wOABA1DJWDYFcANmpRstStT35o4/2ejTwkgVxluxCQA8A3iCQQ67997/SlCm8eMPVlrv3XunVV6XPPmNGMjN69GjVq1dPxx13HBMEIKFWr3ZZq16uyJG1GkBWR2zuXFe3xs5KIdpshbZ1a8kaJcycST3BrJC1CiDRJk6UrEvg7NnSkUfKE4QBAfPjj9Lbb1tBVYI4ONaedPRoadEiacQIZgUAUsUWVKwtV82a8gyBXMA8+qjbTjv/fK9HAj+pX1/q10+64w7X7QMAkJrzcRbEFSokzxDIBciff7racRdfLBUr5vVo4De33+6C/Ouu83okABANczzOWDUEcgHy2mvSmjXSX+0ogb1Y9qp1+LDEBzu3AQBILgI55MpDD0knnSTVqsXEIXMXXCBZcqZts+7axSwBQLKsXWt9oVmRQw599500ebLrrwlkxbKYrdODPV+efpp5AoBksd6qxqts1QxsrQZoNa5yZenUU70eCfzOerCee650001WcsPr0QBAeLdV09K83yUjkAvI8u2LL7pG6QUKeD0aBMHQoS6Iu+sur0cCAOEN5GrUkIoU8XYcBHIB8Mwz0o4dUo8eXo8EQWGrt9df7+rKWX05AEDia8h5nbFqCOR8zg6s27bqmWdK5cp5PRoEiQVyZcq4j1FGiy4AYc1YNQRyPjdhgrRggdSrl9cjQdBYrUErR/Lmm9L48YqsPn36aM6cOZpiDYoBIAHs6MovvxDIIQeeespVjW7ZkulC7p13ntSsmStHsnMnMwgAicxYZUUO2frjD+n116VLLnH9NIHcsufNyJHSjBmuKwgAIDHbqvb7tU4deY6tVR97+WVp+3bpoou8HgmCrEkTqVs36ZZbpA0bvB4NAIQjkKta1R/tMgnkfMxWUDp2lMqX93okCLo773RB3LBhXo8EAIJvjk8SHQyBnE/NnClNnUrJESSuHImdk7vvPunXX5lVADjQQM7rjg4ZCOR8nORgpSM6dfJ6JAiLG2902wC33eb1SAAguDZtkpYskerWlS8QyPnQ1q3S889LF14oFSzo9WgQFunpLogbM0aaPdvr0QBAcAsBm/r15QsEcj70zjvSmjUuWxVIpJ49perVpRtuYF4BIC/sjbBlrHJGDtluqzZv7p9lW4RHoULSkCHSBx9In3+uSKCzA4BEB3LWY9UPGasmXywWi3k9CPzP0qVSlSrSY49Jl17KzCDx7Ce+RQtp2zbJmh2kRWRdfv369UpPT9e6detUsmRJr4cDIKDat5cOOkgaO1a+EJFf4cFhZ+OKFJHOOcfrkSCsbEvg3nul776TXnrJ69EAQPBW5Or75HycIZDz2UqJBXJdukglSng9GoTZ8cdLp53mkh9sZQ4AsH92fn35cgI5ZMHaKFltmgsuYIqQfP/+t7R4sfTEE8w2AORERsZ/gwbyDVbkfOSFF6TSpaWTTvJ6JIgC2xqwNw133OHqIgEA9h/IWVmwmjXlGwRyPrFzp/Tii1LXrtSOQ+oMHuy2Ch58kFkHgJwEcnXq+Ot1mkDOJyZMcK2Tzj/f65EgSqpVk664Qho6VFq71uvRAIC/zfZZooMhkPMJS3KwujRNm3o9EkTNLbe4hIdhw7weCQD4OyFx1iwCOWRiyxbpjTfceSUrDQGkUtmyUr9+0siRLhsLAPB3y5ZJ69b5K9HBsCLnA++9Z8VK2VaFd667ztUvtMQHAEDWGatsrSLTbdUmTfyVBYNoKVVKuvFG6fHHpUWLFCq06AKQqEDOOjpY9yU/oUWXxyxjsHx5V2n/qqu8Hg2ibPNm6YgjpLZt3ZuLsKFFF4ADcdFF0ty50tdfy1fYWvXYa69Ju3bRkgveswbQ1unByuDMnOn1aADAfytyDXx2Ps4QyHnMel2eeKI7cA54rUcPqXp16eabvR4JAPir1uucOf47H2cI5DxkGYITJ7IaB/+wIpeW8GAJOJMmeT0aAPCHH3+U/vyTQA77sJIj+fO75uWAX5xzjnTkkdLtt3s9EgDwh9k+zVg1rMh5fD7O+qoecoiXowD2lpYmDRokffqpWzEGgKibNcv1QvfjMSgCOY+3Vc8+26sRAFk74wypYUNW5QBgz9ZcfizaTyDn4bZqgQJSly5ejQDIflVu8GBp3Dh3AUCUzfZhj9UMBHIeefVVl6168MFejQDInp3dbNTIrcpZj0EAiKKtW6X58/0byBXI6w1jsZg2bNiQ2NFEbFv1oYdcay7Ar6zbgyU/vPuu1Lq1AmXr1q3xK0PG7ysrDAwAuTkfZ+VHqlVL/Wt2iRIllG8/+7l57uyQUSUdAAAAibdu3TqVLFkyOYFcblfkLPCrXLmyli5dut9B7eu4447TlClT8jDKvN82mffZoYNUvLj0+uuJmZ8DGa8f5yczQZsfL+4zWT9jlr161lnuXKcdB0jUeJN9u31X5JYvX64mTZpozpw5qlixYkrGeiC35WeM+Qnacy+sr/O33CK9/bZbmUv1/ORkRS7PW6v2hfPygmq3ye3t8ufPn6f7OpDbJus+ly2TJk+WnnrK5iIx87O/+/TT7Q70tkGZH6/uMxk/Y5bB2ry5NHSodPrpf8/aCtJzL+MXY5R/B4XlZ4z58efchvF1ft486eij//ea7cX8BD7ZoU+fPim/bbLuMyNbNdFFgFP9fXrxmByIIH2ffpsfC9z+9S/pm2+kDz7I3W3zep/JuN2BCNrzIEhzxPwkZ368mtu88vP3OXOmdNRRB35/yZLnrdW8nqnLyX5vmLVs6aL699/f+/PMT/aYH2/nyH5LtGolbdokTZ3qz1pK+/PLL7/s3taoVKmS18PxHX7GmB+eQ3+3cqUrAmwF/O2IiR9/xlK2Ile4cGHdfvvt8Y9R9dtv0pdfSmee+fd/Y36yx/x4O0cZq3LffSe9844CKWNeovw7KDv8jDE/PIcyX40zViDdrz9jKVuRg/TEE9IVV0grVkiHHcaMIHjatpV+/12aNi14q3KsOAHIreHDXbKD5XZab3Q/CsQZubB46y3phBMI4hBcVhx4xgxXVw4Awm7mTKlBA/8GcYZALkUsmv/sM1pyIdjsnNw//uG2WVnLBxB2M2YkZls1mQjkUuTjj12bj0RnqwKpdttt0rffSh99xNwDCK/t26U5cwjksMe2qkX11aszJQj+ObkWLaTBg1mVAxBe8+dL27btXXokUityS5YsUY8ePVStWjUVLVpUNWrUiGdrbLNZyYblXgwaNEgVKlSI365169b6/vvvFfSo/r33XDHVPd15551q0aKFihUrplKlSuXoa3Xv3j1ejHnPq1mzZgqjvMxPGJ8/WVm7dq26desWT2e3y/78xx9/JP35Y0kOt94qff21Oy7gd6NHj1a9evXi1dij7qGHHor/Ti5SpIgaN26sidb0OQvjx4//23PFrrlz5yqMvvjiC3Xu3Dn+u8O+z7fs3fd+TJgwIT6PNp/Vq1fXI488orDK7fyE4fkz86+MVTsjtz9DhgyJ/46xguNlypRRly5dNM8qCafgOZS0QM4erF27dunRRx+Nv5COGDEiPsCbbrop29sNGzZMw4cP16hRo+ItMMqVK6eTTjopV+3A/GbCBOuX9vfzcRbUnn322erVq1euvl6HDh3i7YYyrg8yq9IaAnmZnzA+f7Jy3nnnafr06froo4/il/3ZgrlUPH9OPtna1ATjrJwV77S2XHlt/xMWr7zyivr166ebb75Z06ZNU8uWLdWxY0f9/PPP2d7OXoz2fL7UrFlTYbRp0yYdddRR8d8dObF48WKdcsop8Xm0+bTXtquuukpvWNX3EMrt/ITh+TNjhlS5snTwwcpRQGa/ayZPnqxPP/1UO3bsUPv27ePzlvTnUCyFhg0bFqtWrVqW/75r165YuXLlYnfffffuz/3555+x9PT02COPPBILqj59YrEqVez7y/zfx4wZE/8ec+Kiiy6KnXbaabEoyen8hPX5k5k5c+ZY+BSbPHny7s9NmjQp/rm5c+em5Pnz7rsWwsVi48bFAmHdunXx+bGPUdSkSZNYz5499/pcnTp1YjfeeGOm/3/cuHHx+Vq7dm0sauz7Hjt2bLb/5/rrr4/P356uuOKKWLNmzWJhl5P5CcPzp2PHWKxTp7zdduXKlfHvf8KECUl/DqU02cGqGR9yyCHZRqcrVqyIR7EZrHBeq1at9NVXXymI7ClvK9C2Gpeoulu2ZG1Lt7Vq1dJll12mlVZ6GqF8/mRl0qRJ8e3Upk2b7v6cbZHa5/b3vSbq+dOpk9SokVuVg/9Xt7/99tu9fjaM/X1/z5dGjRqpfPnyateuncaNG5fkkQbrZ3Df+Tz55JM1depUbbfzNAj882fmPq25chvvmOxinkQ9h1IWyC1atEgPPvigevbsmeX/sRdhU9b6YezB/p7xb0Fj7YyWLUtc2RHbCnnhhRf0+eef67777otvF7Vt21ZbLSU24sL4/MmKfT8WjO3LPpfd95rI50/GWTn73ZzNUSv4wOrVq7Vz585c/WzYi+9jjz0W3+Z58803Vbt27fiLsZ2VgvsZzGw+bUvN5jvqgv78WbPGvXbnpfSILVpec801OuGEE1S/fv2kP4dyHcjZQfLMDjDueVk0uadff/01fi7Hzjtdeuml+70P+xr7Tsq+n/OrfeenSZM77SmhNm0KZDk/uXHOOeeoU6dO8SeHHTz98MMPNX/+fL2/b/NWn8rL8ye3wvT8yW5+Mvue9ve9Jvr5Y+V07CDwHXfk6eZIsdz8bNgLr63YHnPMMWrevHk8UcKeO/fee2+KRhvM+czs81EU9OfPzANozdW3b1/NnDlTL730UkqeQwXyMsCuXbtm+3+qVq26VxDXpk2b+ANp0Xl27GB6RpRq0XwG2/rZN2r1q33n59RTq6lBgz81ZMjsTOfnQNk8ValSRQsWLFAQ5Pb5kxthfP5kNT/2S+I3a967j1WrVuXqez3Q509amluV+7//kyZPtu3dPH0ZJFnp0qWVP3/+v62+5fZnw7bvn3/++SSMMHjs901m81mgQAEdeuihno3Lz4L0/Jk5047mSLnNzbjyyiv1zjvvxFceK1WqlJLnUIG8/EKwKyeWLVsWD+IstXbMmDFKs9/62bC0ePvGLOPD9tUzznZYNsjQoUMVBHvOz8KFtqUs3XNPYdWpk56U+1uzZo2WLl26V+DiZ7l5/uRW2J4/2bE3RnYG45tvvlGTJk3in/v666/jn7OSLal8/px5plSvnluVC8jCcOQUKlQo/nvYfjZO36MOkv39tFxUKbfMuqD8rkk2+xl8d59edZ988omOPfZYFSxY0LNx+VmQnj8zZ0q2K1ogh1GSraRZEDd27Nj4OWR7PUrZcyiWJMuWLYsdccQRsbZt28Z++eWX2PLly3dfe6pdu3bszTff3P13yzi0LEP73KxZs2LnnnturHz58rH169fHgub++2OxQoVisQ0bMv/3n376KTZt2rTY4MGDY8WLF4//2a4Ne9xgz/mxz1977bWxr776KrZ48eJ4VlDz5s1jFStWDOT87E9u5ydsz5/96dChQ6xhw4bxbFW7GjRoEDv11FP3+j+pev688ILLYJ0yJeZbUc9affnll2MFCxaMPfnkk/Gs5379+sUOOuig2JIlS+L/btmr3bp12/3/R4wYEc9MnD9/fmz27Nnxf7f5e+ONN2JhZD8fGb9j7PscPnx4/M/2eyiz+fnxxx9jxYoVi/Xv3z8+nzavNr+vv/56LIxyOz9Bf/40bhyLXXxxzv9/r1694q8948eP3yve2bx58+7/k6znkJJZMsIetMyuvQYgxf/vniUkbr/99ngZicKFC8f+8Y9/xF+Qg+ikk2Kx9u1j2ZaCyGx+7AU2s/mxJ0T79u1jhx12WPzBPvzww+Nf4+eff46FUW7nJ2zPn/1Zs2ZN7Pzzz4+VKFEiftmf9031T9XzZ8eOWKxWrVjMz5Vxoh7ImdGjR8eqVKkSK1SoUOyYY47ZqzSCPRdatWq1++9Dhw6N1ahRI1akSJHYwQcfHDvhhBNi77//fiysMspl7HvZvGQ2P8ZetBs1ahSfz6pVq8YefvjhWFjldn6C/PzZvj0WK1LEgtGc3yareGfP16dkPYfy/TUAJJjVn7Udsnvuka66iulF+I0ZI11yiTR7tnTkkfJVZwe7LGvTEjts+7lkyZJeDwuAT82dK9WtK/3nP64lod8RyCXJ2LHSGWe4c3I1aiTrXgD/sO57RxwhtWolPfecfGf9+vXxOnsEcgCy88orkuWcrVrlFmT8LqUFgaPEDn3Xrk0Qh+goVEi67jrJMu5//NHr0QBA3hMdKlQIRhBnCOSSwDarrX3lqacm46sD/tWjh1Uyd0cKACCogVzDPNSP8wqBXBJMmyYtX+5aGAFRUqyY1L+/9NRT7mcAAIJmxoy8t+byAoFckrZV7Sz1CSck46sD/ta7t1SkiDR8uNcjAYDc+f13aelSVuQizwI564NLTUhEUXq6daiQHn7Y/VIEgKCYPt19/KumfCCwIpdgK1dK33zDtiqi7eqrpV27pFGjvB4JAOQukCtaVKpVS4FBIJdgH37okh06dkz0VwaCo0wZ6dJLpZEjpY0bvR4NAOQ8kLNEh/z5FRgEcknYVj3uOCkgPdqBpBkwwGq3SY89xiQDCAZLVjz6aAUKgVwCbd8uffwx26qAOfxwqVs36d57pa1bmRMA/vbnn9IPPwTrfJwhkEugL790KxCUHQGcG26QVqyQnnnGuxmx9lz16tXTcbZUDgBZsPaCO3cGb0WOFl0JftGyF6xff5XSCJGBuP/7P+nbb6V586QCBbybFFp0AcjOE09IV1zheqVbTcygINxIINtWtbIjBHHA/wwc6Fp2vfoqswLA3+fjatcOVhBnCOQSxLaPrBr0yScn6isC4WDnTSyLe8gQV5IEAPyasdooYOfjDIFcgnzyift40kmJ+opAuFbl7PyJZXUDgN/s2uUWY4J2Ps4QyCVwW/WYY1z9LAB7a9nStay76y5XZxEA/GThQmnTJgK5SEfytiLHtiqQ/arc5MnSxInMEgB/tuY6mhW56B6QXL2aQA7Ijp2Tq19fGjaMeQLgv9fxihWlww5T4LC1mqBt1eLFpebNE/HVgHDKl0+67jp3Ts7OywGAX0wPaKKDIZBLUCDXtq1UqFAivhoQXueeK1WuLN1zj9cjAYC9A7kgbqsaArkDZIUDv/qKbVUgJwoWlPr3l158UVq6lDkD4I/yYStWEMhF1rhx0o4dBHJATl16qTuKMGJEauaMFl0AcpLowNZqhLdVa9RwF4D9K1FC6tNHeuwxae3a5M9Ynz59NGfOHE2ZMoWHB0CmiQ4lS0pVqyqQ2FpNQCBH2REgd6680q1kP/wwMwfAH+fj0gIaEQV02P6waJG7COSA3ClbVureXRo5UtqyhdkD4O2K3NEBTXQwBHIH4NNPpQIFpNatE/eAAFExYIC0apX07LNejwRAVK1fLy1YENzzcYZA7gB89pnUtKnbWweQO0ccIZ15pnTvvdLOncweAG9W40zjxgosArkDaMtlGavt2iX2AQGi5PrrXY/DsWO9HgmAKPruO6loUaluXQUWgdwBHI78/XcCOSC31q5dq27duik9PV0nnpiusmXnaMiQHYrFsr5N9+7dlS9fvr2uZs2aMfkADsi330pHHeWOSQUVgVwe/ec/UrFiEq8lQO6cd955mj59uj766KP4VaTIA/ruuwKaMCH723Xo0EHLly/ffX3wwQdMPYADDuSCvK1qAhyDeh/ItWxJWy4gN3744Yd48DZ58mQ1tQOmkl56KZ9atJih226roS++KJ7lbQsXLqxy5cox4QAS1plp3jzXAzrIWJHLg61bpYkT2VYFcmvSpEnxLdWMIM40b95MRYuO0sSJxTVzZta3HT9+vMqUKaNatWrpsssu08qVK/fzc7pV69ev3+sCgD2PSNmRjqCvyBHI5cHkydLmzdKJJyb+AQHCbMWKFfFgbF8VK36pUqXWadiwzG/XsWNHvfDCC/r888913333xbs0tG3bNh6sZWXIkCHxoDHjqly5ciK/FQAh2FYtXFiqV0+BRiCXx23VQw91ByQBSIMGDfpbMsK+19SpU+NTZX/+u+064YQpevll6aef/v6v55xzjjp16qT69eurc+fO+vDDDzV//ny9//77WU7/wIEDtW7dut3X0qVLeagA7JWxaq/jBQsq0Dgjl8dArk2b4LbzABKtb9++6tq1a7b/p2rVqpo5c6Z+++23v/3bqlWrdO21y/XVV9KIEdL992d/f+XLl1eVKlW0wCp5ZnOmzi4AyGpFrlUrBR6BXB4OR37zjfTgg8l5QIAgKl26dPzan+bNm8dXx7755hs1adIk/rmvv/46/rk2bZqoTx/pvvukW291q95ZWbNmTXyFzQI6AMitTZukuXOla65R4LGmlEtffOGafVMIGMi9unXrxsuIWLKCZa7aZX8+9dRTVbt2bfXt64pt1649UmP/qhK8ceNGDRgwIJ4osWTJknjSg22vWuB4+umn8zAAyFOig/2uOeYYBR6BXB62Ve3MtLUXApB7lrTQoEEDtW/fPn41bNhQzz33XPzfLA/ikktsxe08rVq1Mf65/Pnza9asWTrttNPiGasXXXRR/KMFdiVKlOAhAJCnbdVChaQjj1Tg5YvFsqunjn01bCgde6z01FPMDZAMixZJtWq54wu9eyfu61r5EctetW3ckjRIBiLtooukOXOkKVMUeKzI5YKVrZo1i21VIJlq1JDOOksaPlzauZO5BpB434ago0MGArlcGDfOfWzbNkmPBoC4AQPcytzbbzMhABJr82brMkMgF0njx0t16ljpA69HAoTbcce5Fnj33uv1SACEzYwZLtGBFbmIBnKtW3s9CiA6q3KTJileWw4AErmtWrBgOBIdDFurObRihas5E4bigUAQnHqqS3qwunIHYvTo0apXr56Os2U+AJH37bdSgwauPVcYEMjlon6cIZADUsM6p1ixTisnt3Bh3r9Onz59NGfOnHh/VgD4NkSJDoZALhfbqrVrcz4OSKULL3QdHvbXsgsAcmLLFld2hEAuooEcq3FAahUtan1cXd3GNWuYfQAHnuiwcyeBXORYj29LVSbRAUg9KwpsZcsffpjZB3BgvvtOKlBAql8/PDPJ1moOcD4O8M5hh7kq7KNGSX/+ySMB4MDOx9WvLxUpEp5ZJJDL4bZqzZpShQrJf0AA/F3//q6zygsvMDsA8m7KFNdmM0wI5HKA+nGAtyzRqHNnV4rECnkCQG5t2iR9/73UpEm45o5Abj9sFcAyXDgfB3hfINjOqn70EY8EgLydj9u1i0AucjgfB/jDCSe4X8C07QKQF9984zLhw9LRIQMrcjnYVj3iCKlixdQ8IAAyly+fdO210rhx7p11TtHZAUDG+bhjjnFZq2GSLxazxH5kxbJbmjeXHn+cOQK8tmOHSzxq0SL3iQ/r169Xenq61q1bp5IlSyZriAB8qnp1qUsXafhwhQorctlYtcodjKQQMOAP9k7aMlhfeUX6+WevRwMgKFavlhYvlsLYcplALhucjwP855JLpBIlpJEjvR4JgKCY8ler5bBlrBoCuf0EcrYUW7ly6h4QANkrXlzq2dMdd1i3jtkCkLNEh0MOca/pYUMgl42JE6WWLVP3YADImSuvdF0eOLsKIKcrcrataklTYUMgl4X1611zXSt5AMBfrMvK+edL998vbdvm9WgA+Fks5lbkwritagjksjB5siscSCAH+NM110jLlkmvvur1SAD42c8/u+TFMCY6GAK5LPz3v9Khh7rWQAD8p0ED6eSTXdsuiigByIqtxhkCuQgGcrYaF8b9dCBMbbumT5c+/9zrkQDwcyB3+OFSuXIKJVbkMrF9u9taJdEB8Ld27aSGDWnbBWD/iQ5hRSCXiWnTpC1bOB8H+J2tmNuq3EcfueLdmaFFFxBdO3dKU6eGN9HBEMhlUXbEGus2apT6BwRA7pxzjstizartTp8+fTRnzhxNyagICiAy5s6VNm1iRS6S5+OaNpUKFfJ6JAD2x35Or7pKev55acUK5gvA3ufjbOW+cWOFFity+7Dst4xEBwDBcPnlUsGC0qhRXo8EgN8Cubp1pZIlFVoEcvuYP9811yWQA4Lj4IOlHj2khx922ygAEIVEB0Mgtw9bjUtLk5o39+YBAZA3/fpJf/whPfMMMwhA8TZ+1qEpzIkOhkAuk0DuqKPCvQwLhFG1atIZZ0gjRrhMNQDRNmOGtGMHK3KRw/k4ILiuvVZauFB6912vRwLAa19/7ZKhrNZkmLEitwfLeLMXAc7HAcHUrJnUooVr2wUg2iZNctmqhQsr1Ajk9lmNM8cf79GjASAhq3L2s5zRXxFANE2e7N7chR2B3B7sl7+ds6lY0bsHBMCBOe00qUYNVuWAqO+wLVkSjcRFArk9cD4OCL78+aX+/aXXX3e/yGnRBURzNc6wIhchGza4HqucjwOCr3t3KT1dGjmSFl1AVAO5ihWlypUVeqzI7ZHdsmuX1LKltw8IgAN30EFSr17SE0+42nIAopfo0CwC5+MMgdxfvvzSVYevXdvbBwRAYvTtK23bJj3+ODMKRMmOHa6jQxTOxxkCuX2id+vqACD4ypeXzjtPeuABaft2r0cDIFVmzZK2bGFFLlJsS9X206MSvQNRcc010i+/SK++6vVIAKRyYaZgQemYY6Ix56w/SZo7V1q3jkAOCJsGDaT27V0pkljM69EASIXJk6Wjj5aKFo3GfBPI/RW958sX/sa6gB/ceeedatGihYoVK6ZSpUrl6DaxWEyDBg1ShQoVVLRoUbVu3Vrff/99jgsEW0b6xIkHOHAAgXlNbx6hHTYCub8e9Pr1pZIlvX44gPDbtm2bzj77bPWytNIcGjZsmIYPH65Ro0ZpypQpKleunE466SRtsLpB+3HSSW5l7sEHD3DgAHxv9WrXajMqGauGQO6vZdgoRe+AlwYPHqz+/furgUVXOVyNu//++3XzzTfrjDPOUP369fXMM89o8+bNevHFF/d7e1ttt7Nyn3ySgMEDCEQh4OYRek2PfCBnZ+PmzInWgw4EyeLFi7VixQq1t8NufylcuLBatWqlr776Ksvbbd26VevXr49fnTqt12GH7UrRiAF4GciVLStVqRKdxyDygZwVArZD0ARygD9ZEGfK2m/nPdjfM/4tM0OGDFF6enr8KlMmXatW3RH//MqVSR4wAE8DuWbN3Ep8VEQ+kLPzcYccItWq5fVDAQSXJSLky5cv22vq1KkHdB/2Nfbdct33c3saOHCg1q1bt/uaObN//PPW7QFA+Ozc6bo0RW1hpoAiLqMQcJSidyDR+vbtq65du2b7f6pWrZqnr22JDcZW38pbld+/rFy58m+rdHuy7Ve7MmRstVinh9tvj05pAiAq5syRNm6MVqKDoh7IZRQCHjDA65EAwVa6dOn4lQzVqlWLB3OffvqpGjVqtDvzdcKECRo6dGiuv97vv0vPPitdcUUSBgvA04WZ/PmlY4+N1oMQ6a1VCgEDqffzzz9r+vTp8Y87d+6M/9mujfZW+i916tTR2LFj43+27dN+/frprrvuin9u9uzZ6t69e7wO3XnWgyuXOneWhg93b+QAhCuQa9hQOuggRUqBqD/o1luVQsBA6tx2223x8iEZMlbZxo0bFy/0a+bNmxc/15bh+uuv15YtW9S7d2+tXbtWTZs21SeffKISJUrk+v6vvNJ1e3j/fRfUAQiHL7+UOnRQ5OSL2YnhiLr0UmnKFGnGDK9HAiDZrAyJZbBagHjyySVlx+fGj2fegTD47Tc7Tyu9/LJ0zjmKlEhvrUatjQeA/7XtmjBB+vZbZgQIg6/+Kil5/PGKnMgGcn/8QSFgIKpOP92SKKT77vN6JAAS4b//dZnplSpFbz4jG8hZrRnDihwQPZbZ1q+f9Oqrlnzh9WgAJOJ83PERXI2LdCBnZUesEHDNml6PBEAyjR49WvXq1dNxxx231+cvuUSyXIkHHmD+gSDbvNkdk4hqIBfZZIdOnVz5gQ8/9HokAFKd7FCyZMn45268UXr4YWnpUumvTwEImAkTJEt4t8RFKz8SNZFckbPQ9ZtvpKZNvR4JAC9ZKRJ7N0/bLiDY26olS0pHHqlIimQgt3ixtHo19eOAqKtYUTr3XGnkSGn7dq9HAyCviQ7Nm7uzr1EUyUDOVuPMPkdmAESQlSKxhIfXX/d6JABya9cuV3rkhBOiO3eRDeSqV5cOO8zrkQDw2lFHSe3auVIk0TwxDATX999L1gQmqokOkQ7kaMsFYM9VOct6++IL5gQI2vm4/Pmj/ZoeuUDOzsHYL+woP+gA9mb9GevVo0AwEMRA7phjpIMOUmRFLpCbPVv6808yVgH8T758blXu3XelefOYGSBIiQ7HR3hbNZKBnG2r2jJso0ZejwSAn5x/vlS2rDRihNcjAZATy5ZJS5YQyKVFsTWXFQwsWtTrkQDwk8KFpb59pWeekVat8no0AHKyrWpYkYsYEh2AaMmqRVdmevVy26zW7QGA/wO56tWl8uUVaZFq0bVhg5SeLj35pHTxxV6PBoDXLboy07u39MYb0k8/SUWKpHSIAHLBkhwaNHCr6FEWqa3VqVNdnSgyVgFkpX9/t7X6/PPMEeBXVjtu+nSpVSuvR+K9tKhtqxYvLtWp4/VIAPhVzZrSP/8pDR/uqsYD8Ge2qi3M/OMfXo/Ee5EL5OyYTFT7sQHIGStF8sMP0kcfMWOAH02YIFWoINWo4fVIvBe5QI5tVQD7Y30b7U2fte0C4D/WhcVW4/Ll83ok3otMIPfrr9IvvxDIAdi/jALBn3/uzuEA8I+NG92Zd87HRSyQs9U407Sp1yMBEARnnilVqcKqHOA3kyZJO3dyPi5ygdyUKa7WTMWKXo8EQBAUKCBdfbX08stuNR+Af87HHXaYVLeu1yPxh8gEcrYMe+yxXo8CQJD06CEVKyY9+KDXIwGQgfNxEQzkLEWZQA5Ablnd4Msvlx591BUUB+CtLVtcq03KjkQskLOmur//zoocEEW5adGVmauukjZtch1hAHh/3n3bNhIdIhfI2WqcadzY65EASLU+ffpozpw5mmIHZfOgcmXpnHOk+++XduxI+PAA5PJ8XKlSUv36TFukArlvv5UqVZLKlvV6JACCyEqRWO/VN9/0eiRAtFkg17Ilhf0jF8hxPg7AgWjUSGrTxpUisTO3AFLPtlSt9Ajn4yIWyJHoACBRq3J2PufLL5lPwKtFGUt2oBBwxAK5RYukdetIdABwYDp2lOrUoUAw4GXZkeLF3Qo5IhTIkegAIBHS0qRrrpHefltasIA5Bbw4H3f88a5YNyIUyFmig7XZKV3a65EACLpu3dzvEstgBZA6ljFuxxo4HxfBQI5EBwCJUqSIlTORxoyR1qxhXoFUvpZbUe62bZnzSAVyu3a5FTlacwFIlN69XRLVI48wp0CqfP65VKIEr+eRC+TsHItF8ARyABLFmnVfeKE0apS0dSvzCqQqkLNtVc7HRSyQI9EBwIG26MpM//7SihXSiy8yv0Cy/fmnOx/Htmrm8sVi4S1vmZFhZiVIAETb+vXrlZ6ernXr1qlkyZIH/PU6d3Z9nGfOlPLlS8gQAWRi/HhXkHvaNOnoo5miyK3Isa0KIFkFgmfPlj75hPkFkr2tesghUsOGzHOkArmdO6XvviOQA5AcVl3+mGMoEAwk27hxbkXOajni70I7LfPmSZs2EcgBSA7bTrVVuU8/lWbMYJaBZLDX8cmTOR8XyUAuI9HB3jEDQDKcfbZUubI0fDjzCyTDf//rigGT6BDBQM7qx9WsKaWnez0SAGFVsKB09dXSSy9Jy5Z5PRognOfjypeXatf2eiT+FdpAjkQHwJ/uvPNOtWjRQsWKFVOpUqVydJvu3bsrX758e13NmjWTH1x2mVS0qPTAA16PBAhnIGercWSGRyyQs0SH6dOlxo29HgmAfW3btk1nn322evXqlavJ6dChg5YvX777+uCDD3wxuVbJpGdP1+lh3TqvRwOEx9q1LmmRbdUIBnLz50ubN3M+DvCjwYMHq3///mrQoEGuble4cGGVK1du93WI1SPwCdtetaKljz7q9UiA8PjiC9dqk0AugoGcRfCGwoFAeIwfP15lypRRrVq1dNlll2nlypXZ/v+tW7fGiwDveSVLhQpSt27S/ffTtgtI5LZqtWpS1arMaeQCOav+bA/8wQd7PRIAidCxY0e98MIL+vzzz3XfffdpypQpatu2bTxYy8qQIUPinRwyrsqWXppE113n2nY9/3xS7waIjP/8x9WPQ0QDOcqOAKkzaNCgvyUj7HtNzagJlAfnnHOOOnXqpPr166tz58768MMPNX/+fL3//vtZ3mbgwIHxdlwZ19KlS5VMllV32mnSPfe47SAAeffrr9L330snncQs7k8BhYx1jrWtVSvUCSA1+vbtq65du2b7f6omcH+kfPnyqlKlihYsWJDtmTq7UumGG6TmzaV33pG6dEnpXQOhYoW2LVP1xBO9Hon/hS6Q++kn6Y8/pEaNvB4JEB2lS5eOX6myZs2a+AqbBXR+YhVRWraUhg51q3OUTADyxnoY285aCn+tBFZaGLdVDVurgD/9/PPPmj59evzjzp0743+2a+PGjbv/T506dTR27Nj4n+3zAwYM0KRJk7RkyZJ40oNtr1rgePrpp8tvrr/etRSyivQAcs+OJtiKXPv2zF4kV+RsW7VsWVcJGoD/3HbbbXrmmWd2/73RX8vn48aNU+vWreN/njdvXvxcm8mfP79mzZqlZ599Vn/88Ud8Fa5NmzZ65ZVXVKJECfnNKadIRx4pDRvmVucA5I71Ll61ikAup/LFYnaqLDxOPdUVBP7wQ69HAsBPrPyIZa9agFjSqvgmkcWp3btLs2ZJ9esn9a6A0LGjCXfcYUco7Kyr16Pxv1BurbKtCsBL554rVaok3XsvjwOQl/NxtjhPEBfBQO6331zKMokOALxUqJDUv7/0wgtSkqueAKGyaZM7X8r5uIgGchmJDgRyALx22WVS8eKu2wOAnLfl2raNQC7SgVx6ulS9utcjAeAXo0ePVr169XTcccel9H4tD6N3b+mxx1zzbwA521a1JixWYBsRTHY4+2yX6TJ+vNcjARDlZIc9j3tUqWKZutJNN6XkLoFAs4xvK6r9xBNejyQ4QrciR6IDAL+wUkiWvTpypPTnn16PBvC3X36R5sxhWzWygZyVnFq0iPNxAPxlwAC3U/D0016PBAhGW6527bweSbCEJpCbPt19JNEBgJ8ccYR01lnSPfdIO3Z4PRrA3+fjGjeWDj3U65EES1qYtlWLFLHWPl6PBAD2NnCg9OOP0quvMjNAZmjLlXehCuQaNpQKhK7pGICgs52Cjh2lIUPcCxaAv7+GWycH6sdFOJCzHqtsqwLwK8tanT1beu89r0cC+HNb9aCDXMYqIhjIWTbYDz8QyAHwrxNOcNddd0nhKfoEJMbHH0tt2riuKIhgIPf999LOndLRR3s9EgDIflXu66+pdQns6Y8/XFuuU05hXiIbyM2Y4VKW69f3eiQAkLUOHdwbTluVA/C/bVVbjCGQi3ggV7Om218HAD+06MqMveG0DNbPPpOmTPF6NIA/fPCB6+hgXVAQ0RZdrVtLZcqQ2g/AXy26MmMrD3Xruh2EN9/0bBiAL1gWd/nyrgPK0KFejyaYAr8iZ2GorcgddZTXIwGA/cufX7rxRmnsWNeOCIiyqVOllSulTp28HklwBT6QW7rUHZQkkAMQFBdcIFWqJN19t9cjAbzfVk1Pp+xIpAM5W40zBHIAgsJKLFgP1hdflJYs8Xo0gHfef186+WSpYEEehUgHcgcf7N7dAkBQXHqp+91lPViBKFqxwm2tsq16YEIRyNlqnGWDAUBQWJZ9v37Sk0+6FzQgaj76yL12W1ke5F1oAjkACJo+fdw26/33ez0SwJttVasKZFUnENFAbtMmaeFCAjkAwVSqlNS7t/TQQ9LatV6PBkid7dtdIWC2VSMeyM2a5cqPsCIHIKj693cvaqNHez0SIHW+/NJqOxLIKeqBnG2rWk2mevW8HgkA5E3ZslKPHtLIkW6XAYjKtqo99xs18nokwRf4QK5OHalIEa9HAsCv/NSiKytWisS2Vh9/3OuRAKkL5Ky3alqgoxB/CHSLruOPl6pWlV54weuRAPA7v7ToyspFF7kerIsW8eYU4bZ4sVS9uvTaa9JZZ3k9muBLC3J/tpkzOR8HIBxuusmVIXnqKa9HAiTX22+7bG0rBIwIB3IW0W/cSCAHIBxq15a6dpWGDJG2bvV6NEDyvPWW1K6dVKIEsxzpQI7WXADC5pZbpGXLpGee8XokQHKsWSNNnCh16cIMJ0qgAzkrIliunNcjAYDEqFtX+r//k+66y5UkAcLmvfdc2bB//tPrkYRHoAM56scBCOOq3E8/Sc895/VIgORsqzZtyiJMIhHIAYCP1K8vnXmmdOedrMohXDZvlj7+mG3VRAtkILdunbRkCStyAMK7Kvfjj9KLL3o9EiBxrLzOli0EcokWyEDOyo4YtlYBhNHRR0unneZW5Xbs8Ho0QOLKjlh2tl2IeCBn5+OsBo11dQCAMLr1VmnBAumVV7weCXDgdu6U3nmH1bhkCGwgZ/1VCxb0eiQA/C4ILboy07ixayj+73+7F0EgyL76Slq9mkAuGQLZoqtJExfIPf201yMBEBR+b9GVmW++cRl+L78snXOO16MBDqyfsLXTtDqJ9FeN+IqcvTOdPZvzcQDCz960dugg3XGHa0sIBJEtF1nZEasdRxCXeIEL5BYudFkvDRt6PRIASL7bbpO+/156801mG8Fkz99Fi9hWTZbABXKzZrmPDRp4PRIASL7mzaUTT2RVDsE1dqxUvLjUtq3XIwmnQAZy1prLLgCIyqqclV2yrD8gaF57TercWSpc2OuRhFPgAjk7H8dqHIAoadlSat1aGjSIs3IIlnnz3ALM2Wd7PZLwClwgZ08Ia2EDIHiWLFmiHj16qFq1aipatKhq1Kih22+/Xdu2bcv2dpZcP2jQIFWoUCF+u9atW+t7O3gTIYMHu9JLtk0FBMXrr7ttVUvaQXKkBa1PmyU7sCIHBNPcuXO1a9cuPfroo/FAbMSIEXrkkUd00003ZXu7YcOGafjw4Ro1apSmTJmicuXK6aSTTtKGDRsUFf/4hzsrd/vtrMohWNuqp54qFS3q9UjCK1B15KZOlaym59dfu7R8AMF3zz336OGHH9aP1lw0E/Yrylbi+vXrpxtuuCH+ua1bt6ps2bIaOnSorrjiitDWkdvX5Mku+eGll6SuXb0eDZA960xSq5ZblTvzTGYrWdKCdj7OWDFgAOFggdUhhxyS5b8vXrxYK1asUPv27Xd/rnDhwmrVqpW+snLxWbBgz4K3Pa+ga9ZMOuUUd1aObg/wOwvgihWTOnb0eiThlha083HVq7v9dgDBt2jRIj344IPq2bNnlv/HgjhjK3B7sr9n/FtmhgwZEl+By7gqV66ssJyVswPkL77o9UiAnG2rWjCH5AlcIMf5OMB/LBEhX7582V5T7WzEHn799Vd16NBBZ599ti699NL93od9jX23XPf93J4GDhwYX+3LuJYuXaowOPZY6bTTXEC3Y4fXowEyZwWAp00jWzUVCihggVwOft8DSLG+ffuq634ObVWtWnWvIK5NmzZq3ry5HnvssWxvZ4kNxlbfypcvv/vzK1eu/Nsq3Z5s+9WuMLIg7uijpWeflS65xOvRAJlvq1qCA9uqyReYQG71avtFTukRwI9Kly4dv3Ji2bJl8SCucePGGjNmjNL203zRSpVYMPfpp5+qUaNG8c9ZuZIJEybEkx2i6KijpLPOct0eLrhAKlTI6xEBf99W7dRJOuggZibZ0oKW6MDWKhBcthJnNeDsvNq9996rVatWxVfa9j3rVqdOHY39q2CabZ9axupdd90V/9zs2bPVvXt3FStWTOedd56iyhIefvpJGjPG65EAe1u8WPr2W7ZVU6VAkLZV7V1nzZpejwRAXn3yySdauHBh/KpUqdJe/7ZnJaR58+bFz7VluP7667Vlyxb17t1ba9euVdOmTeNfq0SJEpF9MI480pUg+fe/pe7daX8E/7DyOJbgYCtySL7A1JG7/HLpm2+k6dO9HgmAIApDHbl9WfaqlWMaMUK66iqvRwPYGzJ3BMrOcL7wAjOSCoHaWqU1FwD8T+3abjXOVuUi1OQCPt89mzNHivCph5RLC0qEb4Ec5+MAYG/WsstqHd9/PzMD79kq3KGHSnvU70aSBSKQswO99m6TQA4A9nb44VLv3tbqzGX3A17Ztcudj/u//5MKFuRxSJW0IGWssrUKILdGjx6tevXq6Thr1BxSAwe6j3ff7fVIEGVffilZ3W22VVMrLSh77unpUkg67ABIoT59+mjOnDmaMmVKaOf9sMOka6+VRo1yL6SAV9uqVapILVow/6kUmEDOVuOy6cYDAJF2zTWSVWP517+8HgmiaNs2VwT43HOl/dT4RoIFYrrpsQoA2bMg7uabpaeecmVJgFT65BPp99/ZVvVCWhCi/LlzOR8HAPvTs6dUsaJ0663MFVK/rWo7ZyQlpp7vA7n586UdO3hyAMD+FCkiDR7strisRRKQChs3Sm+/LZ1/PvPthbQgbKsaMlYBYP+6dbNetdJNNzFbSA1ri7xli2sZh9RLC0LpkQoVpEMO8XokAOB/BQpId97pziyNG+f1aBAFTz8ttW4tVa3q9UiiKRCBHKtxAJBzp58uWdk8qy8XjG7aCCor2P/5565VHLzh+0Du+++lI4/0ehQAEBxWqsmKA3/9tfTGG16PBmH2zDNS8eLSWWd5PZLo8nUgZ3vuP/5IIAcAudW2rdSxo3TjjS77H0hGSy7bVj37bOmgg5hfr/g6kLOyI7YtUK+e1yMBEFRRaNGVlWHDpMWLpUce8XokCKOJE93zi21Vb+WLxfx7guL5510G1h9/uBZdAJBX69evV3p6utatW6eSJUtGZiIvvVR66y1p4UKpVCmvR4Mwufhi6Ysv3HOLzkveSfP7+bhKlQjiACCvrK7c5s3S0KHMIRJbO87qFdpqHEGct3wdyM2Zw/k4ADgQ1unh2mul+++Xli5lLpEYr7/u3iBcdBEz6jXfr8hxPg4ADsx117lerLfcwkwiMSzJwRJqDj+cGfWabwM5i/TJWAWAA2dHAgcNkp57Tpo2jRnFgVmwQJowwZ2Rg/fS/J6xSg05ADhwl10m1a4t9e9PkWAcmMcec92WzjyTmfSDND+fjzNsrQLAgStYUBo+3K2kWBYrkBdbt0pjxrgkhyJFmEM/SPPz+bjKld2WAADgwFmB4A4dpAED3AsykFvWKWTNGunyy5k7v/B1IMdqHAAk1n33uf6YDzzAzCL3Hn1Uat3abdPDH3wdyHE+DgASy94g9+ol/fvf0sqVzC5y7ocfXAHgK65g1vwkza8Zq9b2g0AOwIGKcouurFgGa1qadNttXo8EQUtyOOww6fTTvR4JfN+i67vvpMaNpcmTpaZNvR4NgDCIaouurIwcKV1zjStH0rCh16OB323Z4opLW/YzXUL8Jc2v26qmbl2vRwIA4dS7t1Szpgvm/Pd2Hn5j7bjWrnWBHPzFt4EcGasAkNxyJJb48J//SG+/zUxj/0kOJ54oHXEEM+U3vgzk6LEKAMl3yimuHEm/fu5sMpCZ6dOlr76SevZkfvzIl4EcGasAkHz58rkyJMuXS3ffzYwjcw8+KFWqJJ12GjPkR2l+zVilhhwAJJ+dk7vuOmnYMGnhQmYce1u9WnrhBXemskABZseP0vxYp4YeqwCQOjfdJJUtK119NYkP2NsTT7iPJDn4l+8COXqsAkBqFSsm3X+/9MEH0rvvMvtwduywOozS+edLpUszK36V5sfzcYcfLpUo4fVIACA6unRxiQ9XXUXiA5y33pJ++UW68kpmxM98GcjR0QEAUovEB+xrxAipZUvp6KOZGz/zZSBHogOARKFFV86R+IAMVm7ErgEDmBO/81WLLstYLV7cHa685BKvRwMgTGjRlfPfw9ZVp3596b333Eodosf6qVryoZ1bt7688C9fPTzz5rmMKVbkAMAbJD5g/nzX7ePaawnigsBXgZxF/4YeqwDgfeKDHXLfuJFHImqGD5fKlJG6dfN6JAhcIDd3rlS+vJSe7vVIACC6bDt11Chp1Srpttu8Hg1SaeVK6emnXRBfpAhzHwS+W5GrU8frUQAAatSQBg2SRo6Upk5lPqLC6sZZB4devbweCQIbyLGtCgD+cM01UsOGrqr/9u1ejwapSHSxQK5HD+mQQ5jvoEjzUwVpO2BJIAeE15IlS9SjRw9Vq1ZNRYsWVY0aNXT77bdr27Zt2d6ue/fuypcv315Xs2bNUjbuqLKVmccfl2bOdDXFEG5PPimtXSv17+/1SJAbvmmB++OP7h0fgRwQXnPnztWuXbv06KOP6ogjjtDs2bN12WWXadOmTbr33nuzvW2HDh00ZsyY3X8vVKhQCkaMY491PVhtm/XMMxXfckX4bN0qDR3q2nFVrer1aBDIQI6MVSD8LBizK0P16tU1b948Pfzww/sN5AoXLqxy5crl+L62bt0av/asI4e8+de/pDfflHr2lD75hNpyYWTvkX79Vbr5Zq9HgsBurVogV7Kky1oFEB3r1q3TITk4kDN+/HiVKVNGtWrViq/irbT0umwMGTJE6enpu6/KlSsncNTRYoXaH35Y+uwz6bnnvB4NEs1ONgwZInXtKtWuzfwGjW86O1x0kSsIPHmy1yMBkCqLFi3SMccco/vuu0+XXnpplv/vlVdeUfHixVWlShUtXrxYt956q3bs2KFvv/02vlKX0xU5C+YscCxp7xqRa+edJ330kWulyJvucJ2Nsx+/2bPpdR5EvgnkmjRxT6A9jsAACIhBgwZp8ODB2f6fKVOm6Fg7cPWXX3/9Va1atYpfT1hfvlxYvnx5PKh7+eWXdcYZZ+ToNrToOnBr1rjf08cdJ73zDlusYWBn020VrnFj6bXXvB4NAntGzkJJKwZ81llejwRAXvTt21ddbV8mG1X3OEFtQVybNm3UvHlzPfbYY7m+v/Lly8cDuQULFuRpvMibQw+V7OE67TTpmWcsm5iZDLoXX5QWL5bGjvV6JAh0IGcHLDdsoBgwEFSlS5eOXzmxbNmyeBDXuHHjeBZqWh46cq9Zs0ZLly6NB3RIrX/+U7rwQpfJ2q6dxNHDYK/G/fvfLjA/6iivR4NAJzuQsQpEg63EtW7dOn5WzbJUV61apRUrVsSvPdWpU0dj/1oi2LhxowYMGKBJkybF69BZ0kPnzp3jgePpp5/u0XcSbdbtoUQJd67KH4dzkBd2lGnhQpeVjOAq4JdAzkpCVavm9UgAJNMnn3yihQsXxq9KlSrt9W97Hte1kiSWlGDy58+vWbNm6dlnn9Uff/wRX4WzFT1LgChh0QRSrlQpyY41duzotlqvuIIHIWi2bHEB3Lnnuu4dCC5fJDv07i1NnCjNmuX1SACEFckOiXf55e6Mlf3u5o14sNx3n3TDDW4hpWZNr0eDUGyt0tEBAIIXDNjRyIsvlnbt8no0yCmrjW1146ynKkFc8BHIAQDyxHa27ZzVhAnSgw8yiUExfLidPZVuvdXrkSAUgZw16P3tN1bkACCI2rSRrrrKbdPNnOn1aJCTKhHWDa9vX2mfY6oIKM8DOasfZ9haBZAMo0ePVr169XScVbFFUlizdSsqa6UEN29mkv3sllukIkXcR4SD58kOTz3lUtg3bZKKFvVyJADCjGSH5J91tu4AF1zgMlnhP9OmucfItsH79PF6NAjNipz98FvBd4I4AAgu21V54AHp8cdp9eRHtmRzzTWu8D7lYsLF8zpyZKwCQDhYFuSnn0qXXeb6Z1ep4vWIkMF6444fL73/vlTA81d+hG5FjvNxABB8+fJJjz7qCgafd560Y4fXI4LZtk0aMEBq394VcUa4pHldWdqa9RLIAUA4WBD30kvS119LgwZ5PRpk1Puz11r7aME2wsXTQG7+fLdvTyAHAOHRvLlr/3TXXdLHH3s9mmizAO6OO6R+/aT69b0eDUIXyNm2qiGQA4BwufFGqUMHt8X6009ejya6rr5aOvRQVkfDLM3rGnJly0oHH+zlKAAAiZaWJj3/vOv+cPbZ0tatzHGqvf229O670siRUvHizH9Yeb4ix2ocAITTIYdIb7whzZjhtvaQOlab1TpuWHLD6acz82HmeSBnNW0AAOFkBWhHjZIeeUR68kmvRxMd1rlh5Uo39yQ4hJtngdyuXdKCBQRyAJKLFl3es+49VoS2Vy9p4kSvRxN+Nse2nXrnnVL16l6PBqFt0bVkiVStmvThh+5ALAAkEy26vLV9u6tjNnu2NGWK6+iD5GypHnWUVK6cNGGClD8/sxx2nq3IzZvnPlqjZQBAuBUs6Fp3WfLDP/8pbdjg9YjCaeBA6ddfpTFjCOKiwtNArnBh6fDDvRoBACCVSpd2WZRW26xbN3fEBoljLbgefFC6+26pZk1mNio8DeTsicayLwBEx5FHus4P1vvTDuQjMTZulC6+WPrHP6S+fZnVKCngZSDHtioARM+pp0rDhknXXSdVqeISIXBg+vSRVq2SPvvM1fBDdHgayF14oVf3DgDw0rXXuo4PvXu7wvBduvB45NWzz/7vqlGDeYyaNK+WgH/5hRU5AIgqq212//3SGWdI554rffml1yMKJlsUsWD4oovcuUNEjyeB3Pz57iNbqwAQXXZG+rnnpKZNpc6dpTlzvB5RsPz5p3TOOVKlSq7wL6LJk0CO0iMAAFOkiPTWWy4YOflkV2MU+2cVYC+/3L2evvoqvVSjzLNArkwZqVQpL+4dAOAn9lpgxeGtJFWbNtLPP3s9Iv8bMcKtZj71lNSwodejQSQDObZVAaQCLbqCoWJFadw4d3bOgrmlS70ekX99/LHL+L3hBne+ENHmSYuuY45xjZQffzzV9wwgqmjRFQyWydqqlesEYQVuLcDD/1iP8iZNpBYtXC0+arEi5StyFjZasgMrcgCAfVldOVuZ27pVattWWr6cOcpgdeKsBp+Va3nxRYI4eBTILVvmmvoSyAEAMlOtmgvm7LWidWvOzBnrTXvKKdK6ddL770vp6Tx34FEgR8YqAGB/rLCtba1u2yYdf7z0ww/RnTNbnbR6e7ab9dFHFP2FDwK5AgXcOy4AALJyxBGuULBltbZsKU2ZEr252rnTFfudOFF6+23p6KO9HhH8xpNAzt5p2UFWAACyU6GCNGGCVKuW22YdOzY687Vjh2tl+frr7kycff+ALwI5zscBAHLqkENcM/hOndwW49ChLnEuzLZvly64QHrlFemll9z3DWSGQA4A4HvFikkvvyzdeqt0443SxRe7s2Nhbb3Vtav0xhuua8PZZ3s9IvhZSgO5LVtcjSBW5AAAuZWWJv3rX9ILL7igLoxdINaskU46SfrgA+nNN1mJg88CuYUL3XI4gRwAIK/OO8+dm7NyVnb435IAwmDRIlfod+5cV36lc2evR4QgSGkgR+kRAKlGi65watpUmj7dJQB06SJdfXWwt1qt7VazZm6xY9Ik92fAl4HcwQdLpUun8l4BRFmfPn00Z84cTYli7YqQs9cTO0f24IPSI4+44G7qVAUuM3XgQKlDB+nYY10QZ2VXAN8Gcratak2RAQA4UPZ60revNHmy+7MFc9dcI23c6P+5/eUXd87vnnuku+92HRsOPdTrUSFoPAnkAABIpEaNpG++cQGRrc7Vry+9954/y5TYmMaMkRo2lJYscef9brjBJXMAuZWWyicugRwAIFms0Px110mzZrkCwpYsYCtetl3pFzNnSu3aSZdc4sZn5/ysBRng+0Bu5UrX7JcVOQBAMln3IEseePddae1alwl62mkuiPLK0qXSZZe5lUPbUrWeqc88w1YqAhTIkbEKwPzzn//U4YcfriJFiqh8+fLq1q2bfv3112wnJxaLadCgQapQoYKKFi2q1q1b6/vvv2dCkSU7L3fqqdK0aa7u3OzZ0lFHuRU6S5CwJINUsKep9UqtXt3VhRsxwo3l5JN58BDAQM72/8nGAaKtTZs2evXVVzVv3jy98cYbWrRokc4666xsbzNs2DANHz5co0aNimeflitXTieddJI2bNiQsnEjmOx1x+rOWW02KyJsAZw93apWdV0iZsxI/Dm6zZvdfVmwZmf1Pv/cnsOuIP5VV0mFCiX2/hBt+WL2VjcFBgxwzY6t4CEAZHjnnXfUpUsXbd26VQXtkNM+7FeUrcT169dPN9iJcFm9sK0qW7ashg4dqiuuuCJHk7l+/Xqlp6dr3bp1KlmyJA9AhFnwNnq0a39lR35sK9bOq7Vs6c6rlS2bu69nr6Lz50tffOE6Mti2rnUysq91+eWu3RbBGwIfyNkPyc6d7kkOAOb3339Xr169tGzZMv33v//NdFJ+/PFH1ahRQ999950a2QGjv5x22mkqVaqUnrGDRpmwYM+uPQO5ypUrE8hht23bXAcF22r99FOXQWoqV5Zq1nSXrdylp0slSrhkCgvQ7LJz3/b/Fy92Z+/sLJ6t/jVp4goUW5N7uz2QbAWUIgsWSB07pureAPiZrazZNunmzZvVrFkzvWd1IrKwYsWK+EdbgduT/f0n26vKwpAhQzR48OAEjhphY6tktv2ZcV7NkhDs/YQFZvaaZbXpbNVu/Xq3EJGhQAFX2N6CvCpVpBNPdMV8bQXOAj4glCtyf/7prlKlUnFvAFLJEhH2FzTZ2bZj7dVO0urVq+OrcRaI2e1sy9OCuXyZVAv/6quvdPzxx8cTIiw5IsNll12mpUuX6iNL/8sEK3JIFHuVtNcvO19XtKgL5AC/SNnTsUgRdwEIn759+6qrHQTKRlVbvvhL6dKl41etWrVUt27d+Jbn5MmT1bx587/dzhIbMlbm9gzkVq5c+bdVuj0VLlw4fgEHyt5fWAAH+BHvKwAcsIzALC8yNgX2PM+2p2rVqsWDuU8//XT3Gblt27ZpwoQJ8WQHAIgyGoIASJlvvvkmfjZu+vTp8W3VcePG6bzzzosnM+y5GlenTh2NtTT3+GpIvnjG6l133RX/3OzZs9W9e3cVK1YsflsAiDJW5ACkjBXzffPNN3X77bdr06ZN8a3SDh066OWXX95rG9RqzFmZkAzXX3+9tmzZot69e2vt2rVq2rSpPvnkE5XgZDmAiEtZsgMAeIk6cgDCiK1VAACAgCKQAwAACCgCOQAAgIDijByASLDjwBs2bIgnSGRWeBgAgohADgAAIKDYWgUAAAgoAjkAAICAIpADAAAIKAI5AACAgCKQAwAACCgCOQAAgIAikAMAAFAw/T9Dqid0EQvPkgAAAABJRU5ErkJggg==",
      "text/plain": [
       "Graphics object consisting of 1 graphics primitive"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "plot(x^3 - 3*x - 1, [x,-2,2])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "69d3070c-4769-4df3-a421-8ae1b1a6321a",
   "metadata": {},
   "source": [
    "So, why do these expressions involve $i$? Can the expression be rewritten in terms of only real expressions?\n",
    "\n",
    "If you restrict yourself to just [algebraic expressions](https://en.wikipedia.org/wiki/Algebraic_expression), the answer is *no*! Closed form algebraic expressions for the roots of this polynomial *must* involve imaginary numbers."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c3875ed2-5b70-4550-b00d-75bd3526aad6",
   "metadata": {},
   "source": [
    "### Viete's formula"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "0cab8259-3733-4fae-a9fd-2aabca8e74a3",
   "metadata": {},
   "source": [
    "It's worth mentioning that the roots can be expressed in terms of trig functions without reference to $i$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "id": "73de72cc-444e-4de7-842d-421b9fe5f319",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = \\sqrt{3} \\sin\\left(\\frac{1}{9} \\, \\pi\\right) - \\cos\\left(\\frac{1}{9} \\, \\pi\\right)\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = \\sqrt{3} \\sin\\left(\\frac{1}{9} \\, \\pi\\right) - \\cos\\left(\\frac{1}{9} \\, \\pi\\right)$"
      ],
      "text/plain": [
       "x == sqrt(3)*sin(1/9*pi) - cos(1/9*pi)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = -\\sqrt{3} \\sin\\left(\\frac{1}{9} \\, \\pi\\right) - \\cos\\left(\\frac{1}{9} \\, \\pi\\right)\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = -\\sqrt{3} \\sin\\left(\\frac{1}{9} \\, \\pi\\right) - \\cos\\left(\\frac{1}{9} \\, \\pi\\right)$"
      ],
      "text/plain": [
       "x == -sqrt(3)*sin(1/9*pi) - cos(1/9*pi)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = 2 \\, \\cos\\left(\\frac{1}{9} \\, \\pi\\right)\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = 2 \\, \\cos\\left(\\frac{1}{9} \\, \\pi\\right)$"
      ],
      "text/plain": [
       "x == 2*cos(1/9*pi)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "surprising_solutions = solve(x^3 - 3*x - 1 == 0, x)\n",
    "for sol in surprising_solutions:\n",
    "    pretty_print(sol.simplify_full())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "82dd0f29-6862-4e00-942b-8732bc57bd68",
   "metadata": {},
   "source": [
    "This follows from [Viete's formula](https://en.wikipedia.org/wiki/Cubic_equation?utm_source=chatgpt.com#Trigonometric_solution_for_three_real_roots), which asserts that, if the roots of \n",
    "$$x^3 + px + q = 0$$\n",
    "are all real, then they can be expressed as\n",
    "$$\\frac{2 \\sqrt{-p} \\cos \\left(\\frac{2 \\pi  k}{3}-\\frac{1}{3} \\arccos\\left(\\frac{3 \\sqrt{3} \\sqrt{-\\frac{1}{p}} q}{2p}\\right)\\right)}{\\sqrt{3}}.$$\n",
    "There's a different root for each choice of $k \\, \\text{ mod } 3$ so that there are three roots.\n",
    "\n",
    "This is not an algebraic expression, though, as it explicitly involves the transcendental cosine and arccosine functions. More to the point, this is not the type of formula that $16^{\\text{th}}$ century mathematicians were looking for when they were forced by cubics to accept complex numbers."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "7a7c3ba7-a541-4d3d-a6be-50a9c78c74d5",
   "metadata": {},
   "source": [
    "## The depressed cubic"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "89c9ab7d-12db-4a9c-beae-cd4c62e16b6f",
   "metadata": {},
   "source": [
    "A *depressed cubic* is a monic cubic with no square term. Thus, we can write it in the form\n",
    "$$g(y) = y^3 + p \\, y + q.$$\n",
    "Or, in Sage code:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "id": "46c23450-9046-40eb-8397-0ce50b680245",
   "metadata": {},
   "outputs": [],
   "source": [
    "var('y,p,q')\n",
    "g(y) = y^3 + p*y + q"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6e3bf7a9-b56c-478a-9959-ba53927681b3",
   "metadata": {},
   "source": [
    "We're going to derive a formula for the general solution to this using Sage to do a bit of the algebraic work for us. We begin by plugging in $u+v$ and expanding:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "id": "f7e0c52b-bef0-4148-afc5-66491f975837",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "u^3 + 3*u^2*v + 3*u*v^2 + v^3 + p*u + p*v + q"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "var('u,v')\n",
    "expand(g(u+v))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "93b7a29b-8bb2-4e85-8744-31f0530f3556",
   "metadata": {},
   "source": [
    "Note that we've got \n",
    "$$(u^3 + v^3 + q) + (3\\,u^2\\,v + 3\\,u\\,v^2 + p\\,u + p\\,v).$$\n",
    "Let's factor that more complicated looking second bunch of stuff:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "id": "a7788e5c-db0b-4b94-b98f-9d9ec02b5ed7",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "(3*u*v + p)*(u + v)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "factor(3*u^2*v + 3*u*v^2 + p*u + p*v)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "92a267d1-2db8-4f9d-99ea-ca81144fbca6",
   "metadata": {},
   "source": [
    "This might all seem a bit mysterious. Our objective, though, is to express $g$ in a simpler form in terms of these new variables $u$ and $v$ - which will be related. Note, for example, that if we choose \n",
    "$$v=-\\frac{p}{3u},$$\n",
    "then that second bunch is zero and we should have a simpler expression for our cubic. In fact, let's just plug in\n",
    "$$u + v = u - \\frac{p}{3u}$$\n",
    "straight away:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "id": "45a5a761-eab1-417b-9361-afd2239b9242",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "u^3 + q - 1/27*p^3/u^3"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "g(u - p/(3*u)).collect(u)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "65e10567-ab47-4a69-938e-dcf31baf9822",
   "metadata": {},
   "source": [
    "Note that if we multiply that through by $u^3$, we obtain a quadratic in $u^3$, which we can solve for $u^3$. We can then take the real cube root and plug back in to $u-p/(3u)$ to get the roots of $g$.\n",
    "\n",
    "Since we're using SageMath, we can just solve for $u^3$ directly:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "id": "41c8d684-4617-4f50-83d9-c371ec81ee45",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle -\\frac{1}{2} \\, q - \\frac{1}{18} \\, \\sqrt{12 \\, p^{3} + 81 \\, q^{2}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle -\\frac{1}{2} \\, q - \\frac{1}{18} \\, \\sqrt{12 \\, p^{3} + 81 \\, q^{2}}$"
      ],
      "text/plain": [
       "-1/2*q - 1/18*sqrt(12*p^3 + 81*q^2)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle -\\frac{1}{2} \\, q + \\frac{1}{18} \\, \\sqrt{12 \\, p^{3} + 81 \\, q^{2}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle -\\frac{1}{2} \\, q + \\frac{1}{18} \\, \\sqrt{12 \\, p^{3} + 81 \\, q^{2}}$"
      ],
      "text/plain": [
       "-1/2*q + 1/18*sqrt(12*p^3 + 81*q^2)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# Define a single variable to represent u^3:\n",
    "var('u3')\n",
    "\n",
    "# Solve and pretty print the solutions:\n",
    "h(u3) = u3 + q - 1/27*p^3/u3\n",
    "sols = solve(h(u3), u3, solution_dict=True)\n",
    "v3 = sols[1][u3]\n",
    "u3 = sols[0][u3]\n",
    "pretty_print(u3)\n",
    "pretty_print(v3)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a8a7196b-ab79-4ef4-8cc6-b3aaa39d3a69",
   "metadata": {},
   "source": [
    "The two solutions look exactly as if they might arise from the quadratic formula.\n",
    "\n",
    "Now, remember that `u3` was shorthand for $u^3$ and that `v` was already defined via `v=-p/(3u)`.  Thus, I guess that `v3` should be shorthand for $v^3=-p^3/(27u^3)$.  If we had solved for `v3` in the first place, the roles would've been switched and $u^3=-p^3/(27v^3)$.\n",
    "\n",
    "The point here is that the function $z \\to -p^3/(27z)$ maps from one root of $z + q - 1/27*p^3/z$ to the conjugate root and those two roots are produced by the two pretty print statements above. We can double check that as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "id": "a2adae8e-dfb6-41eb-bde6-c91af5500d12",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "(-p^3/(27*u3) - v3).simplify_full()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "32070cf7-e4f8-490d-899f-b4d4b0aa3749",
   "metadata": {},
   "source": [
    "Putting this all together, here's our formula for a root of the depressed cubic\n",
    "$$y^3 + p \\, y + q = 0:$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "id": "b3cff6e0-c74b-4a71-87e9-328ae8189bb1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle {\\left(-\\frac{1}{2} \\, q + \\frac{1}{18} \\, \\sqrt{12 \\, p^{3} + 81 \\, q^{2}}\\right)}^{\\frac{1}{3}} + {\\left(-\\frac{1}{2} \\, q - \\frac{1}{18} \\, \\sqrt{12 \\, p^{3} + 81 \\, q^{2}}\\right)}^{\\frac{1}{3}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle {\\left(-\\frac{1}{2} \\, q + \\frac{1}{18} \\, \\sqrt{12 \\, p^{3} + 81 \\, q^{2}}\\right)}^{\\frac{1}{3}} + {\\left(-\\frac{1}{2} \\, q - \\frac{1}{18} \\, \\sqrt{12 \\, p^{3} + 81 \\, q^{2}}\\right)}^{\\frac{1}{3}}$"
      ],
      "text/plain": [
       "(-1/2*q + 1/18*sqrt(12*p^3 + 81*q^2))^(1/3) + (-1/2*q - 1/18*sqrt(12*p^3 + 81*q^2))^(1/3)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "pretty_print(u3^(1/3) + v3^(1/3))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1ff80a60-6768-42c0-bdee-3a06cb879bc2",
   "metadata": {},
   "source": [
    "Once you have one root, you can (in principle) use long division to obtain a quadratic whose roots tell you the other two.\n",
    "\n",
    "More better is to use the other cube roots; that's a little down the road, though."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e76ae71f-d61e-4551-839f-e3ccbb55b523",
   "metadata": {},
   "source": [
    "## Examples"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "65a6ff60-1189-4d49-b5ad-a653225c9edb",
   "metadata": {},
   "source": [
    "Let's apply this formula to a couple of examples."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6c76568f-62e9-4d65-b6e2-766e60707a39",
   "metadata": {},
   "source": [
    "### Example 1"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "ac227efa-821f-4620-9e8e-28c05666db1d",
   "metadata": {},
   "source": [
    "Consider $x^3 - 15x - 4$ so that $p=-15$ and $q=-4$. This is exactly the [example shown here in our text](https://complexanalysis.org/web/sec_origin.html#sec_origin-2-10).\n",
    "\n",
    "Asking Sage to do the arithmetic for us, we get:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "id": "51d12519-a998-4271-8046-a11d1e795bef",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "a = 2\n",
      "b = 11*I\n"
     ]
    }
   ],
   "source": [
    "p = -15\n",
    "q = -4\n",
    "a = -q/2\n",
    "b = sqrt(12*p^3 + 81*q^2)/18\n",
    "print(f\"a = {a}\")\n",
    "print(f\"b = {b}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f0f9aacf-997e-4fac-b7e8-f28814f2073e",
   "metadata": {},
   "source": [
    "Thus, one root should have the form\n",
    "$$\\sqrt[3]{2+11i} + \\sqrt[3]{2-11i}.$$\n",
    "If we assume that the root is real, then these must be complex conjugates of one another. That is, there are real numbers $u$ and $v$ with \n",
    "$$(u + v \\, i)^3 = 2+11i \\text{ and } (u - v \\, i)^3 = 2-11i.$$\n",
    "If we expand these equations we obtain four real equations in the two real unknowns $u$ and $v$. As often happens in complex variables, these equations are actually consistent and redundant. That is, we can expand the left hand side of one, like so:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "id": "4c41345b-acfe-4129-9a73-14c81a00672a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "u^3 + 3*I*u^2*v - 3*u*v^2 - I*v^3 == (11*I + 2)"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "var('u,v')\n",
    "expand((u + v*I)^3 == (2+11*I))"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3096679f-c358-4cd2-88f8-604bf8931ae7",
   "metadata": {},
   "source": [
    "We can then set real and imaginary parts equal to one another and solve:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "id": "c65a54df-c0dd-483f-8e7b-c939227fadd2",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[[u == -1/2*sqrt(-4*sqrt(3) + 7), v == -sqrt(3) - 1/2],\n",
       " [u == -1/2*sqrt(4*sqrt(3) + 7), v == sqrt(3) - 1/2],\n",
       " [u == -I*sqrt(3) - 1, v == -1/2*I*sqrt(3) - 1/2],\n",
       " [u == I*sqrt(3) - 1, v == 1/2*I*sqrt(3) - 1/2],\n",
       " [u == 2, v == 1],\n",
       " [u == (0.9330127018922191 - 1.6160254037844384*I),\n",
       "  v == (-0.6160254037844388 + 1.0669872981077808*I)],\n",
       " [u == (0.9330127018922193 + 1.6160254037844386*I),\n",
       "  v == (-0.6160254037844388 - 1.0669872981077808*I)],\n",
       " [u == (0.06698729810778227 + 0.1160254037844424*I),\n",
       "  v == (1.1160254037844388 + 1.9330127018922196*I)],\n",
       " [u == (0.06698729810778227 - 0.11602540378444239*I),\n",
       "  v == (1.1160254037844388 - 1.9330127018922196*I)]]"
      ]
     },
     "execution_count": 13,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "solve([u^3 - 3*u*v^2 == 2, 3*u^2*v - v^3 == 11], [u,v])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cc399f62-54f4-416c-abea-84119e5fcfde",
   "metadata": {},
   "source": [
    "If we're looking for a real root, perhaps, we should select $u=2$ and $v=1$. That is,\n",
    "\n",
    "$$\\sqrt[3]{2+11i} + \\sqrt[3]{2-11i} = (2-i) + (2+i) = 4.$$\n",
    "\n",
    "This exactly the [the solution shown in our text](https://complexanalysis.org/web/sec_origin.html#Depressed_Cubic). It's easy enough to check that $4$ is, indeed a solution:\n",
    "\n",
    "$$4^3 - 15\\times4 - 4 = 64-60-4=0.$$\n",
    "With that root in hand, it's not hard to factor the polynomial to get\n",
    "\n",
    "$$x^3 - 15x - 4 = (x^2 + 4*x + 1)(x-4).$$\n",
    "The remaining quadratic doesn't factor further over the reals but it's easy to apply the quadratic formula."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3e5c849d-dbcf-4d1b-b009-9182ff46db13",
   "metadata": {},
   "source": [
    "## Example 2"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e9eaeef2-4ad8-4704-b29b-3f2536a22ae4",
   "metadata": {},
   "source": [
    "Consider $x^3 - 3x - 1$, which we already know from a graph to have three real roots but for which Sage produced these solutions involving $i$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "id": "ca6dbb72-9afa-408e-9179-423bc0730ef3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = -\\frac{1}{2} \\, {\\left(i \\, \\sqrt{3} + 1\\right)} {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} - \\frac{-i \\, \\sqrt{3} + 1}{2 \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = -\\frac{1}{2} \\, {\\left(i \\, \\sqrt{3} + 1\\right)} {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} - \\frac{-i \\, \\sqrt{3} + 1}{2 \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}$"
      ],
      "text/plain": [
       "x == -1/2*(I*sqrt(3) + 1)*(1/2*I*sqrt(3) + 1/2)^(1/3) - 1/2*(-I*sqrt(3) + 1)/(1/2*I*sqrt(3) + 1/2)^(1/3)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = -\\frac{1}{2} \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} {\\left(-i \\, \\sqrt{3} + 1\\right)} - \\frac{i \\, \\sqrt{3} + 1}{2 \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = -\\frac{1}{2} \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} {\\left(-i \\, \\sqrt{3} + 1\\right)} - \\frac{i \\, \\sqrt{3} + 1}{2 \\, {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}$"
      ],
      "text/plain": [
       "x == -1/2*(1/2*I*sqrt(3) + 1/2)^(1/3)*(-I*sqrt(3) + 1) - 1/2*(I*sqrt(3) + 1)/(1/2*I*sqrt(3) + 1/2)^(1/3)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    },
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} + \\frac{1}{{\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} + \\frac{1}{{\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}$"
      ],
      "text/plain": [
       "x == (1/2*I*sqrt(3) + 1/2)^(1/3) + 1/(1/2*I*sqrt(3) + 1/2)^(1/3)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "for sol in surprising_solutions:\n",
    "    pretty_print(sol)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "89e4a570-8380-40c0-b9d9-dd638c0772bb",
   "metadata": {},
   "source": [
    "We can now apply what we've learned to see how these complex expressions arise:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "id": "91638162-a25c-4ce9-848d-527ef8fc9ec1",
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "a = 1/2\n",
      "b = 1/2*sqrt(-3)\n"
     ]
    }
   ],
   "source": [
    "p = -3\n",
    "q = -1\n",
    "a = -q/2\n",
    "b = sqrt(12*p^3 + 81*q^2)/18\n",
    "print(f\"a = {a}\")\n",
    "print(f\"b = {b}\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "3ace7ea8-0ac1-4639-819b-1562eeb9aeff",
   "metadata": {},
   "source": [
    "Thus, one solution can be written as"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "id": "708d9cc4-58a6-4361-a49f-1071b7e5bc09",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle {\\left(\\frac{1}{2} \\, \\sqrt{-3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} + {\\left(-\\frac{1}{2} \\, \\sqrt{-3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle {\\left(\\frac{1}{2} \\, \\sqrt{-3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} + {\\left(-\\frac{1}{2} \\, \\sqrt{-3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}$"
      ],
      "text/plain": [
       "(1/2*sqrt(-3) + 1/2)^(1/3) + (-1/2*sqrt(-3) + 1/2)^(1/3)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "z0 = (a + b)^(1/3) + (a - b)^(1/3)\n",
    "pretty_print(z0)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a2348fbe-d05c-4ebd-86bf-370d1f98590e",
   "metadata": {},
   "source": [
    "This looks a touch different from the solution provided by Sage:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "id": "474a879f-a930-4639-bc02-5feefa8291d1",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/html": [
       "<html>\\(\\displaystyle x = {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} + \\frac{1}{{\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}\\)</html>"
      ],
      "text/latex": [
       "$\\displaystyle x = {\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}} + \\frac{1}{{\\left(\\frac{1}{2} i \\, \\sqrt{3} + \\frac{1}{2}\\right)}^{\\frac{1}{3}}}$"
      ],
      "text/plain": [
       "x == (1/2*I*sqrt(3) + 1/2)^(1/3) + 1/(1/2*I*sqrt(3) + 1/2)^(1/3)"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "pretty_print(surprising_solutions[2])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d1e1e0ab-d200-40aa-bfa1-e184f26d39f1",
   "metadata": {},
   "source": [
    "We can use numerical evaluation, though, to show that these are real and equal:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "id": "18e02e58-4c1b-40e1-891d-84129c80d34b",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.87938524157182"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "N(z0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "id": "7d4dff78-1d2c-411e-a6fd-f74a2a49f87a",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "1.87938524157182"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "N(surprising_solutions[2].rhs())"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f5d3a423-c1de-41a6-a9d3-e102d830a4f7",
   "metadata": {},
   "source": [
    "While it might be nice to push this farther, it's not actually feasible in this case. In fact, it can be proved that\n",
    "\n",
    "> If $f(x) is an irreducible cubic polynomial with real roots, then any algebraic representation of those roots *must* contain the imaginary unit."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cd834f25-14b9-40d1-ac05-441a71c2eb3e",
   "metadata": {},
   "source": [
    "## The general cubic"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "f46a7504-aa58-4d3e-8e84-ccb0bd8e1993",
   "metadata": {},
   "source": [
    "Finally, it's worth mentioning that our solution for the depressed cubic can be applied to solve the general cubic:\n",
    "$$z^3 + a\\,z^2 + b\\,z + c = 0.$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "id": "0fff8ebe-58ba-457e-b897-a9019daa8e0d",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "a*t^2 + t^3 + (a + 3*t)*z^2 + z^3 + b*t + (2*a*t + 3*t^2 + b)*z + c"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "var('a,b,c, y,z,t')\n",
    "p(z) = z^3 + a*z^2 + b*z + c\n",
    "expand(p(z+t)).collect(z)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "33a76b6a-5bef-4372-9e7c-45649167c128",
   "metadata": {},
   "source": [
    "Note that, if we set $t = -a/3$, then the $z^2$ term is gone and we obtain depressed cubic:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "id": "517c293e-50f2-421a-8e4b-ba6c2c5c26c3",
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "2/27*a^3 + y^3 - 1/3*a*b - 1/3*(a^2 - 3*b)*y + c"
      ]
     },
     "execution_count": 21,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "expand(p(y-a/3)).collect(y)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d9a4301b-5114-450d-b3b4-6cf6d42c1b53",
   "metadata": {},
   "source": [
    "Since $y$ and $z$ are related via $z=y-a/3$, it's easy to obtain the solution to the original cubic."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "SageMath 10.8",
   "language": "sage",
   "name": "sagemath-10.8"
  },
  "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.13.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}
