LCD Quick User Guide.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "id": "3096494a",
   "metadata": {},
   "source": [
    "## *LCD* - Quick User Guide"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1d98de98",
   "metadata": {},
   "source": [
    "Welcome to the *LCD* Quick User Guide! Throughout this *Jupyter Notebook* you'll be able to test the full functionality of *LCD* implemented algorithms in a very simple way. Just follow this file, modify the code cells as you wish, and then select and execute a method to obtain the results. Recall that you can check the source code for any specific details, and/or the article *Algorithms for computing basins of attraction associated with a rational self-map of the Hopf fibration based on Lyapunov exponents* for a more in-depth mathematical undestanding of what the algorithms do, and how. "
   ]
  },
  {
   "cell_type": "markdown",
   "id": "6b1aee12",
   "metadata": {},
   "source": [
    "First of all, we'll have to install the *LCD* module in order to work. To do this, just follow the instructions in the *How to install* section of the *README.txt* file."
   ]
  },
  {
   "cell_type": "markdown",
   "id": "56baefc6",
   "metadata": {},
   "source": [
    "Regarding the dependencies of the *LCD* module, will plot the graphics using *PyPlot*, and the package *Polynomials* will simplify a few lines of code when we have to work with polynomials.\n",
    "\n",
    "Also, one IMPORTANT NOTE: most *LCD* methods benefit from some multi-threading Julia macros, so make sure you are running Julia with 4 threads! If you don't, the code will work aswell, but it will be much slower. You can easily check the number of threads you are working with simply by executing the next cell."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e83fe9c4",
   "metadata": {},
   "outputs": [],
   "source": [
    "Threads.nthreads()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "19829271",
   "metadata": {},
   "source": [
    "To start julia with 4 threads just execute the command line *julia --threads 4*. You can also change the value of the environment variable *JULIA_NUM_THREADS*, or modify the correspondent setting if you are working with an IDE like Visual Studio Code.\n",
    "\n",
    "Finally, we can start using the *LCD* algorithms! First of all, let's make a quick summary of what they do:\n",
    "The module *LCD* has many different algorithms implemented, and all of them share the same aim; to compute the basins of attraction of a complex rational map. If you are unfamiliar with the theory of iteration of rational maps, or with complex dynamics, this is how some cool-looking fractals like Julia sets are obtained! These basins of attraction are also a very rich and interesting mathematical objects to study. In particular, given a rational map, *LCD* computes the basin of attraction of every *n*-cycle that appears when such rational map is iterated; that is, fixed points (*1*-cycles), *2*-cycles, etc... This is a nice way of visualizing the Julia and Fatou sets associated with the rational map, and provides more information about the basins of attraction than other traditional algorithms. Also, unlike more traditional algorithms, *LCD* does not need to compute the list of fixed points of the rational map (or to have any previous knowledge of the fixed points) before computing the basins of attraction. In addition, *LCD* algorithms are implemented in a way that avoids many computational problems that often appear in numerical analysis, such as mathematical indeterminations or overflows.\n",
    "\n",
    "So, let's start!"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "031c5b62",
   "metadata": {},
   "source": [
    "First, we have to choose a rational map, given by 2 arrays of complex numbers, each one of which represent the coefficients of the numerator and denominator respectively. Here's an example; let's represent f(z)=(z^3-2z+2)(z^2+1):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a832dbca",
   "metadata": {},
   "outputs": [],
   "source": [
    "coefficientlistnum=complex([2.0,-2.0,0.0,1.0])\n",
    "coefficientlistden=complex([1.0,0.0,1.0])"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "e4312f56",
   "metadata": {},
   "source": [
    "It is also very common to choose a complex polynomial and study its dynamical behaviour considering the rational map induced by applying some iterative method (such as Newton-Raphson method, Chebyshev's,...) to the chosen polynomial.\n",
    "With *LCD* you can define a polynomial by its coefficients, and apply an iterative method of your choice to obtain the induced rational map. Currently, only Newton, Chebyshev, and Halley methods are supported."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "f7b07b94",
   "metadata": {},
   "outputs": [],
   "source": [
    "poly=complex([-1.0,0.0,0.0,1.0])\n",
    "\n",
    "pair=LCD.applyIterative_method(cubicPol2cycle,\"newton\")\n",
    "coefficientlistnum=pair[1]\n",
    "coefficientlistden=pair[2]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a0e3e49c",
   "metadata": {},
   "source": [
    "You can of course choose the polynomial you prefer, but here are some interesting polynomials for you to play with:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "e5e4fcfb",
   "metadata": {},
   "outputs": [],
   "source": [
    "tradPol3=complex([-1.0,0.0,0.0,1.0]) # z^3-1\n",
    "tradPol5=complex([-1.0,0.0,0.0,0.0,0.0,1.0]) # z^5-1\n",
    "interestingFatou=complex([-1.0,-20.0,200.0,-550.0,600.0]) # 600z^4-550z^3+200z^2-20z-1\n",
    "interestingJulia=complex([-4.0,0.0,-3.0,0.0,1.0]) # z^4-3z^2-4\n",
    "# With 2-cycles (under Newton's method):\n",
    "pol2cycle=complex([2.0,-2.0,0.0,1.0]) # z^3-2z+2\n",
    "anotherPol2cycle=complex([-5.29,0.0,-4.29,0.0,1.0]) # (z-im)(z+im)(z-2.3)(z+2.3)\n",
    "cubicPol2cycle=complex([-0.845-0.81*im,1.0,-0.845-0.81*im,1.0]) # (z-0.845-0.81*im)(z^2+1)\n",
    "# With 3-cycles (under Newton's method):\n",
    "pol3cycle=complex([-1.49175,1.0,0.0,-1.0]) # -z^3+z-1.49175\n",
    "# With 4-cycles (under Newton's method):\n",
    "pol4cycle=complex([-0.832-0.82*im,1.0,-0.832-0.82*im,1.0]) # (z-0.832-0.82*im)(z^2+1)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "1f0dc206",
   "metadata": {},
   "source": [
    "Once you have chosen a rational map, the following method will compute its list of fixed points. Recall that doing this might be useful, but is not necessary to compute the basin of attraction. This method is completely optional.\n",
    "Note that you also have to choose a rectangular region of the complex plane in which *LCD* is going to look for cycles ((-1.5,1.5)x(-1.5,1.5), in this case), a maximum number of iterations (200), the maximum length of the cycles that are going to be detected (1, because in this method we are only interested in finding the fixed points), and a given tolerance (8, which means that we are going to consider 2 numbers to be equal if they differ in less than 1/10^8). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "6a857147",
   "metadata": {},
   "outputs": [],
   "source": [
    "fixedPoints=LCD.positionIteration_superAttractingFixedPoints(coefficientlistnum,coefficientlistden,(-1.5,1.5),(-1.5,1.5),200,1,8)[2]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "608c8068",
   "metadata": {},
   "source": [
    "Thanks to the next method, we'll be able to obtain our first graphic! It just computes and plots the basins of attraction of a given rational map in a chosen rectangle, also given the rest of the parameters mentioned before. Note that you can choose a coloring strategy; you can color the basins with a different color for each one (\"position\"), color the rectangle taking into account the number of iterations it took for each point to converge (\"iteration\"), taking into account both things (\"positionIteration\"), etc...\n",
    "The results of the method are the plot itself, and the number of colors used."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d15a66c0",
   "metadata": {},
   "outputs": [],
   "source": [
    "graphic=LCD.plotBasinsOfAttraction_Lyapunov(coefficientlistnum,coefficientlistden,(-1.5,1.5),(-1.5,1.5),200,1,8;coloringStrat=\"position\")\n",
    "img=graphic[1]\n",
    "numColors=graphic[2]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "32d38f6e",
   "metadata": {},
   "source": [
    "In order to visualize the results, just execute the following cell. You can uncomment the first line to show the colorbar."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "5352813c",
   "metadata": {},
   "outputs": [],
   "source": [
    "#PyPlot.colorbar(img,ticks=0:numColors)\n",
    "\n",
    "gcf()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2d016504",
   "metadata": {},
   "source": [
    "Alternatively, you can also save the generated graphic on your computer instead of showing it, by execute the following cell."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c484e3b9",
   "metadata": {},
   "outputs": [],
   "source": [
    "savefig(fname=\"cool_JuliaSet.png\")\n",
    "println(\"Figure saved succesfully!\")"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "d4fb7403",
   "metadata": {},
   "source": [
    "If you add \"inverted=true\" in the call for \"positionIteration_superAttractingFixedPoints\", the chosen rectangle will be taken as a neighborhood of infinity in the augmented complex plane; the Riemann sphere (one-point compactification of the complex plane)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "a44843d2",
   "metadata": {},
   "outputs": [],
   "source": [
    "graphic=LCD.plotBasinsOfAttraction_Lyapunov(coefficientlistnum,coefficientlistden,(-1.5,1.5),(-1.5,1.5),200,1,8;coloringStrat=\"position\",inverted=true)\n",
    "img=graphic[1]\n",
    "numColors=graphic[2]\n",
    "\n",
    "gcf()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "c36a00b8",
   "metadata": {},
   "source": [
    "You can also plot ONLY the *n*-cycles (for a given *n*) by executing the following cell. Each cycle is represented as a polygonal in the complex plane whose vertices are the elements of the cycle."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "13e6666e",
   "metadata": {},
   "outputs": [],
   "source": [
    "LCD.plot_nCycles_Lyapunov(coefficientlistnum,coefficientlistden,(-1.5,1.5),(-1.5,1.5),200,2,8)\n",
    "\n",
    "gcf()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "132b8f1d",
   "metadata": {},
   "source": [
    "If you want to plot several *n*-cycles with different lenghts, the following cell plots the cycle which have length between 2 and 4."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "08fd787c",
   "metadata": {},
   "outputs": [],
   "source": [
    "LCD.plot_manyCycles_Lyapunov(coefficientlistnum,coefficientlistden,(-1.5,1.5),(-1.5,1.5),200,2,4,8)"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "a580890a",
   "metadata": {},
   "source": [
    "It might be useful to plot the basins of attraction along with the detected *n*-cycles itselves. In order to do this, you can execute:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c1dce2ba",
   "metadata": {},
   "outputs": [],
   "source": [
    "graphic=LCD.plotWithCycles_BasinsOfAttraction_Lyapunov(coefficientlistnum,coefficientlistden,(-0.2,0.7),(0.0,1.2),200,2,8)\n",
    "img=graphic[1]\n",
    "numColors=graphic[2]\n",
    "\n",
    "gcf()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "2defd818",
   "metadata": {},
   "source": [
    "Sometimes it's interesting to plot a neighborhood of the origin next to a neighborhood of infinity; we can just make a subplot."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "d65598f2",
   "metadata": {},
   "outputs": [],
   "source": [
    "PyPlot.subplot(121)\n",
    "graphic=LCD.plotBasinsOfAttraction_Lyapunov(coefficientlistnum,coefficientlistden,(-1.5,1.5),(-1.5,1.5),25,1,8;coloringStrat=\"position\")\n",
    "numColors=graphic[2]\n",
    "ax1=graphic[1]\n",
    "PyPlot.subplot(122)\n",
    "ax2=LCD.plotBasinsOfAttraction_Lyapunov(coefficientlistnum,coefficientlistden,(-1.5,1.5),(-1.5,1.5),25,1,8;coloringStrat=\"position\",inverted=true)[1]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "9f3bd972",
   "metadata": {},
   "source": [
    "These methods make up the basic functionality of *LCD*. However, there are a couple more, a bit more technical, but very interesting.\n",
    "These methods are specifically designed to study which cubic polynomials present attracting *n*-cycles in their dynamical behaviour under the iteration of complex Newton-Raphson's method. If you are interested, we recommend you to check out the article *Algorithms for computing basins of attraction associated with a rational self-map of the Hopf fibration based on Lyapunov exponents*, but roughly speaking, a one-parameter plane is constructed in a way that every point alpha corresponds to the polynomial (z-alpha)(z^2+1). Since we can fix two of the roots of any cubic polynomial to i and -i via a Möbius transformation, studying the dynamical behaviour of this family of polynomials is sort of \"dynamically equivalent\" to studying every cubic polynomial. The following method plots this parameter plane, and colors each point based on if its corresponding polynomial present some kind of *n*-cyclic behaviour. In particular, 4 different colors can be associated with given point alpha; 3 of them in case the arithmetic mean of the 3 roots of the corresponding polynomial converges to one of the roots under iteration of Newton's method, and another color in case this alpha/3 point diverges or converges to an attracting *n*-cycle. This specific criteria for the coloring of the points is a consequence of Fatou's Theorem.\n",
    "This parameter plane gives us very valuable information, since when an attracting *n*-cycle appears, some iterative methods like Newton's might fail to \"approximate\" the roots of the polynomial."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "17886a4f",
   "metadata": {},
   "outputs": [],
   "source": [
    "graphic=LCD.detectCathastrophicBehaviour_distinguishCycles_cubicPolynomials_Newton((-2.0,2.0),(-2.0,2.0),100,12,8)\n",
    "img=graphic[1]\n",
    "numColors=graphic[2]\n",
    "\n",
    "gcf()"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "29f4ecc7",
   "metadata": {},
   "source": [
    "Also, the next method allows us to plot this parameter plane in a way that it provides even more information about the dynamical behaviour of the cubic polynomials it represents; we can choose the maximum length of the attracting *n*-cycles we want to detect, and a different color is assigned to each detected *n*-cycle. This way, when the arithmetic mean of the 3 roots of the polynomial corresponding to the point alpha converges to a *2*-cycle it uses one color, when it converges to a *3*-cycle it uses another different color, and so on.\n",
    "Graphics colored following this strategy allow us to see the distribution of the *n*-cyclic dynamical behaviour (for an specific *n*) along the parameter plane."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "327f900b",
   "metadata": {},
   "outputs": [],
   "source": [
    "graphic=LCD.detectCathastrophicBehaviour_distinguishCycles_cubicPolynomials_Newton((0.81,0.87),(0.79,0.84),100,12,8)\n",
    "img=graphic[1]\n",
    "numColors=graphic[2]"
   ]
  },
  {
   "cell_type": "markdown",
   "id": "cf940187",
   "metadata": {},
   "source": [
    "The last method is also about this parameter plane. \n",
    "Recall that the algorithms we present calculate a finite set of real non-negative constants (which we call *Lyapunov constants*, since they are closely related to Lyapunov exponents), each one of which is associated with a different basin of attraction of the given rational map. This Lyapunov constants provide us extremely useful information, not only about whether the end point of each basin of attraction is super-attracting, attracting, indiferent or repulsive, but also about how much attractive it is. Lyapunov constants between 0.0 and 1.0 correspond to attracting end points, and the lower the constant is, the more attractive the basins of attraction will be.\n",
    "So, with that in mind, this last algoritms plots the said parameter plane, and colors it based on the value of the Lyapunov constant associated with the orbit of the arithmetic mean of the 3 roots of the corresponding polynomial of each point in the plane. If the constant is 0.0 (super-attracting) one color is assigned, if it is between 0.0 and 1.0 another different color is assigned, and so on for 1.0 (indiferent) and constant values greater than 1.0 (repulsive). This way, the parameter plane provides information not only about which polynomials present *n*-cyclic dynamical behaviour, but also about how fast does the arithmetic root of its 3 roots converges, either to a root or to an attracting *n*-cycle, under iteration of Newton's method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "id": "c3fc3887",
   "metadata": {},
   "outputs": [],
   "source": [
    "graphic=LCD.plotLyapunovConstants_cubicPolynomials_Newton((0.81,0.87),(0.79,0.84),500,10,10)\n",
    "img=graphic[1]\n",
    "numColors=graphic[2]"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Julia 1.6.4",
   "language": "julia",
   "name": "julia-1.6"
  },
  "language_info": {
   "file_extension": ".jl",
   "mimetype": "application/julia",
   "name": "julia",
   "version": "1.6.4"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}