Skip to content

Instantly share code, notes, and snippets.

@davehorner

davehorner/simple.py

Forked from cswiercz/simple.py
Created November 9, 2021 00:31
Show Gist options
  • Save davehorner/adc5fdc28e7b2e1b37bb88db60a63f9d to your computer and use it in GitHub Desktop.
Save davehorner/adc5fdc28e7b2e1b37bb88db60a63f9d to your computer and use it in GitHub Desktop.
Download ZIP

Revisions

  1. @cswiercz cswiercz revised this gist Dec 4, 2014. 1 changed file with 10 additions and 3 deletions.
    13 changes: 10 additions & 3 deletions simple.py
    View file
    Original file line number Diff line number Diff line change
    @@ -3,7 +3,7 @@

    from mpl_toolkits.mplot3d import Axes3D
    import matplotlib
    from matplotlib import cm
    from matplotlib import cm, colors
    from matplotlib import pyplot as plt

    branching_number = 2
    @@ -21,12 +21,19 @@
    # compute w^2 = z. THE KEY IDEA is to pass the exponentiation by 1/2 into exp().
    w = numpy.sqrt(R)*numpy.exp(1j*Theta/2)

    # color by argument
    arguments = numpy.angle(w)
    norm = colors.Normalize(arguments.min(), arguments.max())
    color = cm.jet(norm(arguments))

    fig = plt.figure(figsize=(16,8))

    # plot the real part
    ax_real = fig.add_subplot(1,2,1,projection='3d')
    ax_real.plot_surface(z.real, z.imag, w.real, rstride=1, cstride=1, cmap=cm.jet, alpha=0.5)
    ax_real.plot_surface(z.real, z.imag, w.real,
    rstride=1, cstride=1, alpha=0.5, facecolors=color)

    # plot the imaginary part
    ax_imag = fig.add_subplot(1,2,2,projection='3d')
    ax_imag.plot_surface(z.real, z.imag, w.imag, rstride=1, cstride=1, cmap=cm.jet, alpha=0.5)
    ax_imag.plot_surface(z.real, z.imag, w.imag,
    rstride=1, cstride=1, alpha=0.5, facecolors=color)
  2. @cswiercz cswiercz revised this gist Dec 4, 2014. 1 changed file with 32 additions and 0 deletions.
    32 changes: 32 additions & 0 deletions simple.py
    View file
    Original file line number Diff line number Diff line change
    @@ -0,0 +1,32 @@
    import numpy
    import sympy

    from mpl_toolkits.mplot3d import Axes3D
    import matplotlib
    from matplotlib import cm
    from matplotlib import pyplot as plt

    branching_number = 2

    Nr = 16
    Ntheta = 32

    # compute the theta,R domain
    theta = numpy.linspace(0,2*numpy.pi*branching_number, Ntheta)
    r = numpy.linspace(0,1,Nr)
    Theta, R = numpy.meshgrid(theta,r)

    z = R*numpy.exp(1j*Theta)

    # compute w^2 = z. THE KEY IDEA is to pass the exponentiation by 1/2 into exp().
    w = numpy.sqrt(R)*numpy.exp(1j*Theta/2)

    fig = plt.figure(figsize=(16,8))

    # plot the real part
    ax_real = fig.add_subplot(1,2,1,projection='3d')
    ax_real.plot_surface(z.real, z.imag, w.real, rstride=1, cstride=1, cmap=cm.jet, alpha=0.5)

    # plot the imaginary part
    ax_imag = fig.add_subplot(1,2,2,projection='3d')
    ax_imag.plot_surface(z.real, z.imag, w.imag, rstride=1, cstride=1, cmap=cm.jet, alpha=0.5)
  3. @cswiercz cswiercz created this gist Dec 4, 2014.
    182 changes: 182 additions & 0 deletions simple_notebook.ipynb
    View file
    Original file line number Diff line number Diff line change
    @@ -0,0 +1,182 @@
    {
    "metadata": {
    "name": "",
    "signature": "sha256:e64d7262a527a7b1d4f96247f34f5a60ddf1879cb8a56485959d028e3168dab3"
    },
    "nbformat": 3,
    "nbformat_minor": 0,
    "worksheets": [
    {
    "cells": [
    {
    "cell_type": "code",
    "collapsed": false,
    "input": [
    "%matplotlib inline"
    ],
    "language": "python",
    "metadata": {},
    "outputs": []
    },
    {
    "cell_type": "code",
    "collapsed": false,
    "input": [
    "import numpy\n",
    "import sympy\n",
    "\n",
    "from mpl_toolkits.mplot3d import Axes3D\n",
    "import matplotlib\n",
    "from matplotlib import cm\n",
    "from matplotlib import pyplot as plt"
    ],
    "language": "python",
    "metadata": {},
    "outputs": []
    },
    {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
    "A Naive Approach\n",
    "================\n",
    "\n",
    "...consisting of creating a complex grid and simply taking the square root of the grid."
    ]
    },
    {
    "cell_type": "code",
    "collapsed": false,
    "input": [
    "a,b = -1,1\n",
    "n = 16\n",
    "x,y = numpy.mgrid[a:b:(1j*n),a:b:(1j*n)]\n",
    "z = x + 1j*y\n",
    "w = z**2"
    ],
    "language": "python",
    "metadata": {},
    "outputs": []
    },
    {
    "cell_type": "code",
    "collapsed": false,
    "input": [
    "fig = plt.figure(figsize=(16,8))\n",
    "\n",
    "# plot the real part\n",
    "ax_real = fig.add_subplot(1,2,1,projection='3d')\n",
    "ax_real.plot_surface(z.real, z.imag, w.real,\n",
    " rstride=1, cstride=1, cmap=cm.jet)\n",
    "\n",
    "# plot the imaginary part\n",
    "ax_imag = fig.add_subplot(1,2,2,projection='3d')\n",
    "ax_imag.plot_surface(z.real, z.imag, w.imag,\n",
    " rstride=1, cstride=1, cmap=cm.jet)"
    ],
    "language": "python",
    "metadata": {},
    "outputs": []
    },
    {
    "cell_type": "code",
    "collapsed": false,
    "input": [
    "a,b = -1,1\n",
    "n = 16\n",
    "x,y = numpy.mgrid[a:b:(1j*n),a:b:(1j*n)]\n",
    "z = x + 1j*y\n",
    "w = numpy.sqrt(z)"
    ],
    "language": "python",
    "metadata": {},
    "outputs": []
    },
    {
    "cell_type": "code",
    "collapsed": false,
    "input": [
    "fig = plt.figure(figsize=(16,8))\n",
    "\n",
    "# plot the real part\n",
    "ax_real = fig.add_subplot(1,2,1,projection='3d')\n",
    "ax_real.plot_surface(z.real, z.imag, w.real,\n",
    " rstride=1, cstride=1, cmap=cm.jet)\n",
    "\n",
    "# plot the imaginary part\n",
    "ax_imag = fig.add_subplot(1,2,2,projection='3d')\n",
    "ax_imag.plot_surface(z.real, z.imag, w.imag,\n",
    " rstride=1, cstride=1, cmap=cm.jet)"
    ],
    "language": "python",
    "metadata": {},
    "outputs": []
    },
    {
    "cell_type": "markdown",
    "metadata": {},
    "source": [
    "The Issue\n",
    "=========\n",
    "\n",
    "In many (all) complex square-root functions a branch cut is chosen (usually the negative real axis) and only one branch is computed in order to make square root single-valued.\n",
    "\n",
    "From complex analysis we learn that we can use polar coordinates to easily verify that $w = \\sqrt{z}$, which I will read as $w^2 = z$ since the latter makes it clear that I won't be chosing a branch of $w$ / making a branch cut in the $z$ plane. Because the branch point $z=0$ is 2-ramified / has branching number two I will need to make two rotations to capture all of the behavior near the branch point."
    ]
    },
    {
    "cell_type": "code",
    "collapsed": false,
    "input": [
    "branching_number = 2\n",
    "\n",
    "Nr = 16\n",
    "Ntheta = 32\n",
    "\n",
    "# compute the theta,R domain\n",
    "theta = numpy.linspace(0,2*numpy.pi*branching_number, Ntheta)\n",
    "r = numpy.linspace(0,1,Nr)\n",
    "Theta, R = numpy.meshgrid(theta,r)\n",
    "\n",
    "z = R*numpy.exp(1j*Theta)\n",
    "\n",
    "# compute w^2 = z. THE KEY IDEA is to pass the exponentiation by 1/2 into exp().\n",
    "w = numpy.sqrt(R)*numpy.exp(1j*Theta/2)"
    ],
    "language": "python",
    "metadata": {},
    "outputs": []
    },
    {
    "cell_type": "code",
    "collapsed": false,
    "input": [
    "fig = plt.figure(figsize=(16,8))\n",
    "\n",
    "# plot the real part\n",
    "ax_real = fig.add_subplot(1,2,1,projection='3d')\n",
    "ax_real.plot_surface(z.real, z.imag, w.real,\n",
    " rstride=1, cstride=1, cmap=cm.jet, alpha=0.5)\n",
    "\n",
    "# plot the imaginary part\n",
    "ax_imag = fig.add_subplot(1,2,2,projection='3d')\n",
    "ax_imag.plot_surface(z.real, z.imag, w.imag,\n",
    " rstride=1, cstride=1, cmap=cm.jet, alpha=0.5)"
    ],
    "language": "python",
    "metadata": {},
    "outputs": []
    },
    {
    "cell_type": "code",
    "collapsed": false,
    "input": [],
    "language": "python",
    "metadata": {},
    "outputs": []
    }
    ],
    "metadata": {}
    }
    ]
    }