Maths for transforming time series data of walking into a linear least squares form: Ax=b that can solve for the gains.

walking_least_squares

{
 "metadata": {
  "name": ""
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "This notebook describes the details of forming the coefficient matrix and ordinate vector from walking data for a linear least squares solution to the optimal gains from a simple full state feedback controller employed by the human while walking."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "!pip freeze"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Cython==0.19.1\r\n",
        "-e git+git@github.com:moorepants/DynamicistToolKit.git@83822a2634756ad8a0a1701a0cb3cca3c07877d4#egg=DynamicistToolKit-dev\r\n",
        "Jinja2==2.7.1\r\n",
        "MarkupSafe==0.18\r\n",
        "argparse==1.2.1\r\n",
        "coverage==3.6\r\n",
        "ipython==1.0.0\r\n",
        "matplotlib==1.3.0\r\n",
        "nose==1.3.0\r\n",
        "numpy==1.7.1\r\n",
        "pandas==0.12.0\r\n",
        "pyparsing==2.0.1\r\n",
        "python-dateutil==2.1\r\n",
        "pytz==2013b\r\n",
        "pyzmq==13.1.0\r\n",
        "scipy==0.12.0\r\n",
        "six==1.3.0\r\n",
        "sympy==0.7.3\r\n",
        "tornado==3.1\r\n",
        "wsgiref==0.1.2\r\n"
       ]
      }
     ],
     "prompt_number": 1
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from sympy import init_printing\n",
      "init_printing()"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 2
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "from sympy import symbols, Matrix, Symbol, zeros, eye"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 3
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "n = 3 # number of time samples in each cycle (step)\n",
      "m = 2 # number of cycles (steps)\n",
      "p = 2 # number of sensors\n",
      "q = 2 # number of controls"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 4
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "First form the equation that generates the controls given the error in the sensors at a single time step during a single foot step.\n",
      "\n",
      "$$m_m(t) = m(t)_{lc} + K(t) s_e(t) = m(t)_{lc} + K(t) (s_{lc}(t) - s_m(t))$$\n",
      "\n",
      "Now rearrange the equations such that we have one linear in both the gains, $K(t)$, and in $m^*(t)$:\n",
      "\n",
      "$$m_m(t) = m(t)_{lc} + K(t) s_{lc}(t) - K(t) s_m(t) = m^*(t) - K(t) s_m(t)$$"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The \"measured\" controls, $m_m(t)$:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "mm = Matrix(q, 1, lambda i, j: Symbol('m_m_q={}'.format(i)))\n",
      "mm"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}m_{m q=0}\\\\m_{m q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAEYAAAAZBAMAAAB+2EWVAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAu90iZpkydkTvqxCJ\nVM1ghewcAAABHElEQVQoFWMQUjJgYAlJL2UvTGDABjiUdBgUgBIs/Q6cq9kXoCvhPnIaJGQFVmNd\nxMDfwHoAXQ3jho0INQzLGTIZmC6gq+Fn4NwAFIOYwzCNoZDBno3B9Ypr7zEGBvcag1OrVq26YMnA\nlABXw/wBaJIjD8P1B5wX7jDwBnBCTERRw7GA4TGD1WkGroRNDGcZ+BISIWpQ7YKIMTAxXGVQuMC0\n4TgD2C5UN0PV7GKoYFhuwHFEDcLnDjkMYkDdDFUDoybAGGAaqxrWv4TVoKjAZReqIohd+OIUZg6O\nOGXgAkUF/jjlOGmAUIMjThmYkNVgj1MUNTjiFEUNjjhFUQMLD7Q4xaoGLU7Z6mOR/AUzCEqTE6dK\nwPyFD3AoaQEAe8NjjECjxxQAAAAASUVORK5CYII=\n",
       "prompt_number": 5,
       "text": [
        "\u23a1m_m_q=0\u23a4\n",
        "\u23a2       \u23a5\n",
        "\u23a3m_m_q=1\u23a6"
       ]
      }
     ],
     "prompt_number": 5
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The unknown gains, $K$:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "K = Matrix(q, p, lambda i, j: Symbol('k_q={},p={}'.format(i, j)))\n",
      "K"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}k_{q=0,p=0} & k_{q=0,p=1}\\\\k_{q=1,p=0} & k_{q=1,p=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAJoAAAAmBAMAAADO2Z8qAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEIondMqvN\nZlSyKd/8AAACQUlEQVRIDe2UP2gTURzHvxfT5PKapKE66NRLRXQshS6CbYR2cZAiKRSlWJdat+Ag\ndDK2SAdJaaF/wOkodhAcdHAKiJsgCAUR4iIHorOhu3jv/X53757EvGAjOPSGy+/P933yfu8lXwyP\nnkd/nvHRMob7g1IUpg1csjNv1K2aaG+7ViVyi1YN05y2VYmRklXDtMymVYl7dgnT3LGDlzbxmrhl\nkzAt/+WOntVpPDKWnW2pGT98T8zKtVj3SkVMG7oc14F06UwiA646FZkf+fLND9eidPmiipi2Xw2i\nBjCElHHgP1AOm+LCdS0B1XRhT4VMqwrPx6f52pNmsxncR84HChOtwA3TN1kPT0PtQGVfoDC3nqx9\nXNiAWgKDNiM268WxFH0X0YrzDl0009Kzz118hac0XPtWP0dLTNo7vPbz/jNq8aTXxFvKaap8bfAm\n3otKspb2X1Bq0mQtV3pM25a3IADPPVRTRbcgNVvurK75SKHRaVKpzDTk8YSPc2UdD0vZnQVKcbpV\n4wjLn30Kw1q2jZW5Q0oPjgIZ8C1QzaMP+S6iUNeZjvhoVSFASzdUlKSd+qmbRaxUdKajqUDHQXZb\nJypK0n5r/UX6T2gnbpm8if/CLTFoGBI6uSXIHuO9/9ktM6vxP0mpO7kl22NE6+KWyJm0Tm7JphHR\nTA8x3JJoXd0StLwXtyRaV7dkWi9uyZN2dUvaWy9uKWk2t5S03txS3J60umVojyduGf8GjxuUMTI9\nc1wIr69Ot5FZutsn2oOlxV/mf9/t7EoKewAAAABJRU5ErkJggg==\n",
       "prompt_number": 6,
       "text": [
        "\u23a1k_q=0,p=0  k_q=0,p=1\u23a4\n",
        "\u23a2                    \u23a5\n",
        "\u23a3k_q=1,p=0  k_q=1,p=1\u23a6"
       ]
      }
     ],
     "prompt_number": 6
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The unknown \"limit cycle\" controls, $m^*(t)$:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "ms = Matrix(q, 1, lambda i, j: Symbol('m^*_q={}'.format(i)))\n",
      "ms"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}m^{*}_{q=0}\\\\m^{*}_{q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAADwAAAAnBAMAAAChwzYZAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEImYyq92J\nVM0hoUD5AAABRElEQVQ4EWNgQAYcAsg8EFtISRkuxHs0C84GMoyUFBmEkAVuL0fmMTCgSrOl1qCa\njirNwIaqGaa7uTPUINQTTQ7Ihepu+8GgJRCFS5pddAPDNwYbXNIMXAGcCgx7EdKMJw6COVDDmSaw\nOrB/mACXZxEQQZauB+rncDgAl+ZnYAL7EKo7i4GvgH2NAENrTMGpVatWXShi4AYbhepvXgMmiAHY\npfkmTIFIoxkOtZRb4CgD2HA0p0Gl2U7shLAYLY6DGah2MzAoQBVCKTRp5r94pVEl4VGCLgzjoxkO\nE4bRaNLoSRFVGi0pYthNXlKEpTUcSZG6aY2BC5IdoB5DT2tshwrAIQOVRk9rDNwo0uhpDU0aGsbw\ntIZdGp7WsEvD0xp2aXha44y3RXY5LH4x05q8iytckgE1rYW4fGRgS0tHSKOyytISAPZgUwRC8z7J\nAAAAAElFTkSuQmCC\n",
       "prompt_number": 7,
       "text": [
        "\u23a1m^*_q=0\u23a4\n",
        "\u23a2       \u23a5\n",
        "\u23a3m^*_q=1\u23a6"
       ]
      }
     ],
     "prompt_number": 7
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The measured sensors, $s_m(t)$:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "sm = Matrix(p, 1, lambda i, j: Symbol('s_m_p={}'.format(i)))\n",
      "sm"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}s_{m p=0}\\\\s_{m p=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAEAAAAAZBAMAAABzxjXSAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAu90iiXZUmUQyEKtm\nze+zsTonAAABBUlEQVQoFWMQUjJgYGAJTWDABjiUVBgUQBIVnALo8lw+XiAha4gCwawF6AoYFyxE\nUjD7/wR0BfwMfCBNEBPYDLgDeMrOVKRMYDge48bgvXv37gQ7BiaQJogCpgnsE9hiGZ6xOTCkOiyC\nmIWigCcmgoFhJ5cCSwID44S5EAUoVoCFFNgKmNgY+Bh8IFagOhKogmsDawN3MINR6AGICVw1ziAG\n1JsQMRB5BMEEs9AVMCsSUIAmjcUKdBVQK3BFFtwELJHFwL0AZBrUBCyRxeEHTAZwBVgii4EJSQG2\nyEJRgC2yUBRgiywUBSDXMKBFFoYC9MjijC8DaUOPC3IiSwmUL3ABDiVtAAJzSZ85FRodAAAAAElF\nTkSuQmCC\n",
       "prompt_number": 8,
       "text": [
        "\u23a1s_m_p=0\u23a4\n",
        "\u23a2       \u23a5\n",
        "\u23a3s_m_p=1\u23a6"
       ]
      }
     ],
     "prompt_number": 8
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "And the equation:\n",
      "\n",
      "$$0 = m(t) - m^*(t) + K(t) s(t)$$"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "zero = mm - ms + K * sm\n",
      "zero"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}k_{q=0,p=0} s_{m p=0} + k_{q=0,p=1} s_{m p=1} - m^{*}_{q=0} + m_{m q=0}\\\\k_{q=1,p=0} s_{m p=0} + k_{q=1,p=1} s_{m p=1} - m^{*}_{q=1} + m_{m q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAAnBAMAAAAGFxWOAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEIondMqvN\nZlSyKd/8AAAEwklEQVRYCe1YS4gcVRQ9/e+unq6UUcHgIu2IKCjSG4MEjBWYbLKIg6RjUEfbhWMk\nII0rs0qrDAFlQgtqYBZShGQhuoiLQXFEgxtBOjAwjIyCUiLRuBASAi4EP+//blVXZ7q6R6cV36Lq\nvnvPvffUe6+aOg1MzKh5qansnL4TKNyTOo8mPNqhsxFtd+F4mszc9LSHnSLj7TR5FFvkk0qLuka1\nf/4wkrlpzawin7kWyUsxyXPs7vQb3t+i+OxLkTJ+PyTq0eSL3ah/+Jkg//zw+BsgxSbauG/NZEuT\nLzfOX0hGKG+5GSbHBflXnSeSo9q78JS2Eu5r3zQbzYN9AZ95BoQkVpOf+u4Ze3Ayiy9HKt224WGP\nw2tFhoQJ8pd+JAeHw+moesi1zlGPtD/Qrq9+x71e/9P7LJ4cUgQ1+R17dSV2z3u3kBmwP+Nj5vso\nIwMT5K8HJIHDySi+0kb+0E3EI81jdytX6dYufsWDfQAfGBBSBDX5s4dDm7wDzE3GVdyB9/4MiEeY\nHPZJ7/Neb9W56xAJcjgdlTacz36jHmmfgTvDxn5UZ5063rEAvrDVXu++3qW+kAQpgpr8YaceYP1o\n+/WVlZXwRVQCoPbARlhm009LdSy57epsbd+3e64EWJs7DQLjK1/wzzqoNU8ZOPMZGBj5n7CMtV/W\nXnjNVgXOSCZANij4pauspRpqYX02VSHVmrUgnfkSi9/5A0634zayMluSd49mumLOyVeCXOA+ifvd\nDi53dlEYJ5+ffbeMH1A3cGYYGCe/hL24fDEbXoGpSsifZAtc8zsim1/UwvrMVCHVWreQBA35L7Ac\nTAXneKpJPuhcFFOwc1CbOwJ8nOmWQ+SD96Vb9uDkp9rVx/GlfqXlsTEwTn59vYVqsAunAVPVrvxx\nVqD0kRfbefissg7J1rpF7NhwOhVvQZ4Hvm0OUC+vinMgXlgO6LqNiossFi2MrTqP8PFGedbCAwvj\n5PmosL3phraqJS/CQGznBXkVUq1ZC3FsYi8sBxUXlyQ289ApnPRKb82p1Js3RPvMhUKr+hhONFct\nDMgpEI59HUiTwUvXLMx5ep/w344jWG7bquevhxKvr/Gdb+gAoFrrFpwgG/rMS1xd3vjVRa1jZ9Ta\noBNqqzdGuEIMgA2qyrIiO08LW5u2iJHP/WFhLk74dkas0ptkEjEfDu00HAQbUJVnRnbelqIWbREj\nT2HbYtfTdY0em3S5W46mOz9M8YkiPwxhivmfPF2Nf9LWKz8JGhZpBbgmjwnQsHEB3tpsFzX5SdCw\niAlwf1jyxe5myEFx8W2zJRp2ZAE+ERo2QYCn17BgqpOOlBoWRpzKIluvYRVBfeYjGparTjpSalgj\nTmWNv0HDKoKafETD6i9w/QApNWz8Q519z6fTsHxhb6hhFUFNPqJhZWxUDav03RgaVi2szxYvWcPG\nyEc0rIwZtZlSwyryY2hYxY2TT9awMfIRDatiRm2m1LDy2IyhYQn5ZA0bI88e0ioZHhtDw3LywTga\nlpDnvPiQ8llr2ATyRslw1TmGhmXidDwNq7gN1rBKFusXVj5fXd74dTs1rCJvyRjrX6Bh9f8NhrMx\n/lsadvfMAfNk22uk07CFmUc8FOef217OI3Yvzc/jL9lTJETMjCn1AAAAAElFTkSuQmCC\n",
       "prompt_number": 9,
       "text": [
        "\u23a1k_q=0,p=0\u22c5s_m_p=0 + k_q=0,p=1\u22c5s_m_p=1 - m^*_q=0 + m_m_q=0\u23a4\n",
        "\u23a2                                                         \u23a5\n",
        "\u23a3k_q=1,p=0\u22c5s_m_p=0 + k_q=1,p=1\u22c5s_m_p=1 - m^*_q=1 + m_m_q=1\u23a6"
       ]
      }
     ],
     "prompt_number": 9
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now for a simple linear least squares the equations can be massage into this form:\n",
      "\n",
      "$$Ax=b$$\n",
      "\n",
      "which gives and over determined system of linear equations.\n",
      "\n",
      "We define $x$ as:\n"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x = K.reshape(q * p, 1).col_join(ms)\n",
      "x"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}k_{q=0,p=0}\\\\k_{q=0,p=1}\\\\k_{q=1,p=0}\\\\k_{q=1,p=1}\\\\m^{*}_{q=0}\\\\m^{*}_{q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAFgAAABsCAMAAADZoIW9AAAANlBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABHL6OuAAAAEXRSTlMAMquZdlQQ\nQN0iRO+JzWa7fNDgyowAAAOnSURBVGgF7ZrpeqwgDIajonVhsN7/zZ6AiiQsCsPTsw0/WtbXgPjF\n4AC4SSzSLVbMj3NFmIMatsEpVczOYz1Ys5nUGqLsoBF7tvQK484DaLZOYOoN6KX6Obwcg1QysUxX\ns9K0dtPg5jJu65zCVY25dYBhJTWkwJoVAzeLWHfTyShdwJ6wxE1mzRwsVmgW0EbjhOQEst1TD9OC\nleMxHdUJJdGASLMezixeBTQj4DPSd2gfdrCJgHshYSD7hzTjIA5+TQCtQLteDTR0RclcV2x+2Yvq\nDGn2wbbvOIAU7lztzTOrgZORKth8ELjFFjxLdS7oUTd0ChcdJv1oTqMQyHXT0XxWRcHYgSzxOQBX\nH7PqVqkS4G4L7zusnUa6/NdVbS4Btn2KMj8HriX0vsWVhN4Dh5Utf5k9cC2h98AVhN5MzwMXCz1q\njKupHrhY6Ge54ANvEweXCz2KdQpcLPRoahJcLPR34HKhv7HYLn620D8FZwt9I7bO0Vm+K6zFmPkI\nvbsa5fnUGpdTcaQH/ngQu57/kgdhLsLO8chcoYKp0Lti6B2fEvcgzEVwMAsVJtQNNwZJeBAmXxx8\n834c9yAIulxEzVCBgKuGCq7FUDNUOMHVQ4XDRXxCBb5t78qeHt8NeNrugf8Kof/avq4J1goVvlE7\nXBGC/0jomR/wdkWx0DM/wMHlQn9q1LEPOLhY6JF3+QEscHA8VCAjfaEnzQFwIlQgJvlCT5p9i/Fa\ne/JDhXNkROjP5mM8X4qTC36ocCP0PxkqDMJ5YbEmw3tnQhOerTVkn1zkt3I9niMSEXqL5gyO3jyn\nT1E2C6yC9yJ83QxwL8X8/Iw9Awz92LonEmFLz9oM8DyLKXVIfyL3/ymw6lQjRAdCoDLolPPpJQme\nNzzKXfAMqeD7RQI8DUIflm9oZ3d7qmtm5P5JgPWZOj49+nBhdA4inNEsnnFaErKpey24tOhS8JDb\nfGcgA3WBxTOkPWUx8nar1aq/M/iJxTOkQwpsHgftUqbWbLPklxFCxUIKzPqmv4ywzjngmy8jjJxh\nse8FS28es8H3guzTB+mfYTGOixyYEeJRyALHvODb4BAgVpdlcQwSqv8zwB+fd9ybj8/7HT4v9GCc\ndckH5M7nsZjxZJr/KfCtz2Mx42Mw93lk4F6Iv7KnLOYg7kyxvQrYc6YPwM5PTbiZV9l3phGL7U9N\nBv3TECHu3iV9ZxoBm5+aCHw9e5h8Z8pixoecULcMZxoaHq3LcaYI+QUUbjPuTbLH3gAAAABJRU5E\nrkJggg==\n",
       "prompt_number": 10,
       "text": [
        "\u23a1k_q=0,p=0\u23a4\n",
        "\u23a2         \u23a5\n",
        "\u23a2k_q=0,p=1\u23a5\n",
        "\u23a2         \u23a5\n",
        "\u23a2k_q=1,p=0\u23a5\n",
        "\u23a2         \u23a5\n",
        "\u23a2k_q=1,p=1\u23a5\n",
        "\u23a2         \u23a5\n",
        "\u23a2 m^*_q=0 \u23a5\n",
        "\u23a2         \u23a5\n",
        "\u23a3 m^*_q=1 \u23a6"
       ]
      }
     ],
     "prompt_number": 10
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "To form the A matrix for a single time step in a single foot step the equation can be rearranged:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A = zeros(q, q * p)\n",
      "for row in range(q):\n",
      "    A[row, row * p:(row + 1) * p] =  sm.T\n",
      "A = A.row_join(-eye(q))\n",
      "A"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}s_{m p=0} & s_{m p=1} & 0 & 0 & -1 & 0\\\\0 & 0 & s_{m p=0} & s_{m p=1} & 0 & -1\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAARYAAAAmBAMAAADzbqllAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEIolUMqvd\nZs0eYxJ8AAADf0lEQVRYCc2WzWvUQBjGn2w/dpPN0rV6ECz0Q0R6q4LF+gEtrAc9FbS1Hqp7qnpb\nFDz1sFqKVCnEj2JBhYB6Emw99KAV3INVKIiFUtiL0D/BRVCoB80kM8lkPlhsU9pcZt73nfeZ3yaZ\n7IPWroPYDdfRrk607gYQnyEZlpER3Q9qmiwrSy/DbNSbCEuubM6H2vFJFU/jiSB68YFlud5EWGwX\nNaYtjMtoyQspP3zGklxvIiztJXQz7fho/EbKiaeCKGThehnL5CVVB8llhpRi/PK5PL7ycTS3arBX\nojCahSxcL2VpKEYvE2vYX/Xvbq/VzzLhmPUrYejpLYVBbEJYBmMZGvAsrJeyNJ7dI3UMGP0kV3gQ\n39hLNd8pkUp4zZVwKAxiE8t7RnVYol7KYn3ciGmQ4Ac6yfD6r0uG2GXHWdpd/IrVw8DYQEsxjLyJ\nWfCuASC8L1wvZXmCBaw+X71+D7ljVSezuLj4Lt2BZdJcyg7m+r73zrpYHZ3Cfa/kQGCxi9pz9Lnu\nOYp6KcsyTmBmPeXMwhwxKv4PoSy22+Cal9FtljFTbvMrEFmaerTfl/Pq7wseUiVwvZRlba2IrNuG\nKeCctR6sC55RbnQYeG9UMg4a3TdqFkxMSO8UXZk5foTOYsPNbwssjnopCynYeISKg47Miv+MQN9d\nUqqYPbaJFKbVz4gsSeDiWA5gGAul9NIold1bZW+oMd9UzF7E+NBKULKu9CWwtSTBsQS1XFlawxJV\nNtmmUWIZ79ftlP6kqySUl1gS0t2MjMQS2QlRTmdFxHV1Y0mI/suILJydEDU1VkRcVj8WhayTFb9J\nZOHshCiqsyLiurqxJGRW/B6RhbMTgqbWigjr6oaykIaFsxOCqMKK6E2PvgLIQnoWZicULIOxnML0\nGNO3yRJFBXhLm2VPo2OJ7ARtZYNsRRSmpzG/j6xXVHDtMFWShTQsnJ1gEHSUrAgUpqcFKfIvqagA\nX2Shx8TN5KFhsSM7IbBAsiIK03MD3kEEFJWIRRbSsHB2QmSRrIjC9FAWRYVjkYTMDn8v8Uz/jxVR\nmB76jBQVjkX0NOlTf/xDIbGId6NOLJge+u6SJqHCsegkt8oimB7j9F22k1DBq58OK2nGrbIEsnrT\no6/IQMmw6E2PvrJdLLLuZjLJ3JfN7Cz37C6W9sIZGXEHMhcKNTSPXd2BneUtb40V/wH69wfWJ8Iq\nCQAAAABJRU5ErkJggg==\n",
       "prompt_number": 11,
       "text": [
        "\u23a1s_m_p=0  s_m_p=1     0        0     -1  0 \u23a4\n",
        "\u23a2                                          \u23a5\n",
        "\u23a3   0        0     s_m_p=0  s_m_p=1  0   -1\u23a6"
       ]
      }
     ],
     "prompt_number": 11
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "b = -mm\n",
      "b"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}- m_{m q=0}\\\\- m_{m q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAFUAAAAmBAMAAABHbcoGAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEIlRmMqvd\nic1EIzl8AAABPklEQVQ4EWMQUlJmIAIwKykJMAgRoRCshGlIq+W0XBrMGrIAh2dR/caZP4FpI+sG\ndLWMNy+ChFDVFsYw8CcwX0BXyyIggqy2xQUIBBi2Mixi4GlAV8sPMhLdXAYthhCGei6GyZ2TM24x\nMEyNLbi9e/fuhiIGngUYatk/AE2ewsvQ+oCpoZeBw4AJYgNWtWwbGJ4zVF9k4F4gw3CDgW/BEoha\n7G6AyDHwMDQzHGjgEbjKAHYDmt+giqCUNEMYw9YCtpunIHxGy+sgBmqYQZXCKQU4C8TAq5b5B/Fq\nUVQSMHfQq6VC+mXgBiUztPDFkX7ZLhUgq8Wbfhl4UNRCgxJ7+sWqFkf6xaoWR/rFqhbqEvT0i1ct\nWvrlirdADgeYiai0Agp3CKZfeRdXFC/g4LC6+AswsKWl45BGEWZPS2MAANK2Y7+yx1MrAAAAAElF\nTkSuQmCC\n",
       "prompt_number": 12,
       "text": [
        "\u23a1-m_m_q=0\u23a4\n",
        "\u23a2        \u23a5\n",
        "\u23a3-m_m_q=1\u23a6"
       ]
      }
     ],
     "prompt_number": 12
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now we can check to make sure that our new form gives the same answer as the previous form:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "zero2 = A * x - b \n",
      "zero2"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}k_{q=0,p=0} s_{m p=0} + k_{q=0,p=1} s_{m p=1} - m^{*}_{q=0} + m_{m q=0}\\\\k_{q=1,p=0} s_{m p=0} + k_{q=1,p=1} s_{m p=1} - m^{*}_{q=1} + m_{m q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAAnBAMAAAAGFxWOAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEIondMqvN\nZlSyKd/8AAAEwklEQVRYCe1YS4gcVRQ9/e+unq6UUcHgIu2IKCjSG4MEjBWYbLKIg6RjUEfbhWMk\nII0rs0qrDAFlQgtqYBZShGQhuoiLQXFEgxtBOjAwjIyCUiLRuBASAi4EP+//blVXZ7q6R6cV36Lq\nvnvPvffUe6+aOg1MzKh5qansnL4TKNyTOo8mPNqhsxFtd+F4mszc9LSHnSLj7TR5FFvkk0qLuka1\nf/4wkrlpzawin7kWyUsxyXPs7vQb3t+i+OxLkTJ+PyTq0eSL3ah/+Jkg//zw+BsgxSbauG/NZEuT\nLzfOX0hGKG+5GSbHBflXnSeSo9q78JS2Eu5r3zQbzYN9AZ95BoQkVpOf+u4Ze3Ayiy9HKt224WGP\nw2tFhoQJ8pd+JAeHw+moesi1zlGPtD/Qrq9+x71e/9P7LJ4cUgQ1+R17dSV2z3u3kBmwP+Nj5vso\nIwMT5K8HJIHDySi+0kb+0E3EI81jdytX6dYufsWDfQAfGBBSBDX5s4dDm7wDzE3GVdyB9/4MiEeY\nHPZJ7/Neb9W56xAJcjgdlTacz36jHmmfgTvDxn5UZ5063rEAvrDVXu++3qW+kAQpgpr8YaceYP1o\n+/WVlZXwRVQCoPbARlhm009LdSy57epsbd+3e64EWJs7DQLjK1/wzzqoNU8ZOPMZGBj5n7CMtV/W\nXnjNVgXOSCZANij4pauspRpqYX02VSHVmrUgnfkSi9/5A0634zayMluSd49mumLOyVeCXOA+ifvd\nDi53dlEYJ5+ffbeMH1A3cGYYGCe/hL24fDEbXoGpSsifZAtc8zsim1/UwvrMVCHVWreQBA35L7Ac\nTAXneKpJPuhcFFOwc1CbOwJ8nOmWQ+SD96Vb9uDkp9rVx/GlfqXlsTEwTn59vYVqsAunAVPVrvxx\nVqD0kRfbefissg7J1rpF7NhwOhVvQZ4Hvm0OUC+vinMgXlgO6LqNiossFi2MrTqP8PFGedbCAwvj\n5PmosL3phraqJS/CQGznBXkVUq1ZC3FsYi8sBxUXlyQ289ApnPRKb82p1Js3RPvMhUKr+hhONFct\nDMgpEI59HUiTwUvXLMx5ep/w344jWG7bquevhxKvr/Gdb+gAoFrrFpwgG/rMS1xd3vjVRa1jZ9Ta\noBNqqzdGuEIMgA2qyrIiO08LW5u2iJHP/WFhLk74dkas0ptkEjEfDu00HAQbUJVnRnbelqIWbREj\nT2HbYtfTdY0em3S5W46mOz9M8YkiPwxhivmfPF2Nf9LWKz8JGhZpBbgmjwnQsHEB3tpsFzX5SdCw\niAlwf1jyxe5myEFx8W2zJRp2ZAE+ERo2QYCn17BgqpOOlBoWRpzKIluvYRVBfeYjGparTjpSalgj\nTmWNv0HDKoKafETD6i9w/QApNWz8Q519z6fTsHxhb6hhFUFNPqJhZWxUDav03RgaVi2szxYvWcPG\nyEc0rIwZtZlSwyryY2hYxY2TT9awMfIRDatiRm2m1LDy2IyhYQn5ZA0bI88e0ioZHhtDw3LywTga\nlpDnvPiQ8llr2ATyRslw1TmGhmXidDwNq7gN1rBKFusXVj5fXd74dTs1rCJvyRjrX6Bh9f8NhrMx\n/lsadvfMAfNk22uk07CFmUc8FOef217OI3Yvzc/jL9lTJETMjCn1AAAAAElFTkSuQmCC\n",
       "prompt_number": 13,
       "text": [
        "\u23a1k_q=0,p=0\u22c5s_m_p=0 + k_q=0,p=1\u22c5s_m_p=1 - m^*_q=0 + m_m_q=0\u23a4\n",
        "\u23a2                                                         \u23a5\n",
        "\u23a3k_q=1,p=0\u22c5s_m_p=0 + k_q=1,p=1\u22c5s_m_p=1 - m^*_q=1 + m_m_q=1\u23a6"
       ]
      }
     ],
     "prompt_number": 13
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "assert zero == zero2"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 14
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now for n timesteps we want to form the A matrix such that there are different gains for each time step."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "gain_matrices = []\n",
      "control_vectors = []\n",
      "control_star_vectors = []\n",
      "sensor_vectors = []\n",
      "for time_step in range(n):\n",
      "    gain_matrices.append(Matrix(q, p, lambda i, j: Symbol('k^{}_q={},p={}'.format(time_step, i, j))))\n",
      "    control_vectors.append(Matrix(q, 1, lambda i, j: Symbol('m^{}_m_q={}'.format(time_step, i))))\n",
      "    control_star_vectors.append(Matrix(q, 1, lambda i, j: Symbol('m^*{}_q={}'.format(time_step, i))))\n",
      "    sensor_vectors.append(Matrix(p, 1, lambda i, j: Symbol('s^{}_m_p={}'.format(time_step, i))))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 15
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "gain_matrices"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\begin{bmatrix}\\left[\\begin{smallmatrix}{}k^{0}_{q=0,p=0} & k^{0}_{q=0,p=1}\\\\k^{0}_{q=1,p=0} & k^{0}_{q=1,p=1}\\end{smallmatrix}\\right], & \\left[\\begin{smallmatrix}{}k^{1}_{q=0,p=0} & k^{1}_{q=0,p=1}\\\\k^{1}_{q=1,p=0} & k^{1}_{q=1,p=1}\\end{smallmatrix}\\right], & \\left[\\begin{smallmatrix}{}k^{2}_{q=0,p=0} & k^{2}_{q=0,p=1}\\\\k^{2}_{q=1,p=0} & k^{2}_{q=1,p=1}\\end{smallmatrix}\\right]\\end{bmatrix}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAi4AAAAqBAMAAABmas3TAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEIondMqvN\nZlSyKd/8AAAFjklEQVRoBe1aTYgcVRCuXje7M5PdyZJ4iChsuyJ6XAQRAtERNggeZJENBCUYL2u8\njR6EnJxEJYGwsoI/4MUhGFDwoAdPC7K3gCAERFgvMiAK3gzx4sXYr96b7tfVr6re9vQsiziH2X6v\nq75X8/Xr6ur6FqD0SbYul8ahwcm9pdB03NzRCN9v4qBCVp1rz4WmS3NfvNsvjWMGs0v3q2bPJj3V\nhjOYu6LHdPExzlufn4XXVaPe3LZqQw2OwYx6Qf+Eh6lb/Lit8wKfxMNRywVYVuPfTe5Qt8D4+Moj\nkKysrAIceRzgTWgPA0bjKWMyn8Kn43H17xMrIdLMIogPCi8vDjJIiRcNH/6oxuTNIH5325shh4hv\nvo5nJ5I0+wL4WOXFmMi8QHAzmUUQX+OlfSGzk3jR8JNbZiX2g/gPiXsWr6vPi9leyn2EO1C5j4L7\nxQSKzsp+wbugJi+I/6B8Gxn8RM7rFV5MOlLyLmYsJe+yvKCzwgtmzZq8GPzkqZPsXjEnDP6JNwbm\nkPtUeGmt3vw6OX2VszfzxgRO7PUlG5YX49x55WnJF97pvAw3744kGxG/fe9vyRfxP7uXEch/Krws\n/PKqlqgjTPj8EuP8w+/q84TlpUl8P78cO8WT6M5EmPC8xDjfHaoxsLw0ie/zcmNjVARFC19b5Iom\nYPMZG7fobPE7j75QhAC0sJ42vq3FK/fRRicd5lHRBGyTrWji6lSWF9HZ4h/p3ejkIQBJ8NPGd7V4\nhZczne1BHhR9YNuHs2ji6g6WF9HZ4s+uf9nKQwBaENjn1PTw7bOywsst+HYIP53rX9/Z2RnZwnfx\nyb1RKxt+54o5yWRcj7FxS84Of6F/9CVYPHvVX/PH8+8DhqTxPjE+w4u5UN3VGXu9LC/dc4l9pnlF\nLmei8hKJD79CijG4NX8bPGBD0niZGF/gZWH4uQ3C3UfPd3btuNjTrElM3KxzgQ/fd3r+mrPDr+xw\n6vgCL+2l9+ymNXk3S4Fp6zbuaS8HsiYxcbPOXo79oLVerDmEGdiKu48MexPhC7zMbbl3ZVP4vr00\n/9F5d62KIpc1cXUqm18MEutc4MPFn4d20Wxu/g5cOnvbDqeN72rxSt61q6f2j/nuwuKgGBVHaXEY\nMBF5yXag6IwnXYrD4xHsFQ54dCD4fl2Hq973TxFGFy71ilF+pJjIcSvOuMgzo3wtGM1/WAzw6EDw\nK7yQIGoM5bhrABKXA8H/nxfCuut7UV7oa1HFK5ugry3ERr6ek+sBIn4DegDiU17oaxH5zTj0Hqmh\n02LcDegBIn4DekCQF/paFPrhXgkWOi3GrfV3EbBmv874NqAHVHhpRA9w9+fyLv7C8Rf2vQ2+xgv2\n6yN4YfEn1QNC+aVBPaBd7vw1rgew+IddD3htvFXwL/LSqB7A4R92PcDV7o4d5KVRPYDBP+x6QBKI\nu0k9gMM/7HrAYuk2QlETmuzXTxvfr18i+ukRJv9xPQBIYRrQA6jJvvQA6hzSAxxgvvX2owfUwY/R\nA2hhaovcUkufmuyrX0+dLX5ZDyD/CDNtfBdSpa4r9etpAWaLXNFE69eJzha/rAc4wHy/7EcPCMcv\n4zP9ulI/3eKKegAxUfveEfglPcAB1tMDSHAhvaGCz/Birkve7Le4oh5ATFReIvFzPcAB1tMDSHCe\nnsHjC7zk/XqLC5IeQE20+8jwEoFf6AEWsJ4eQIPzXnZzvYHiC7zk/XSDq+gB1CSGlwh8yPUAA1hX\nD6DBec0RFr/Ey/LaGUj+Wls3V9Pr15vmuKIHUBPs12+sld8YLaxZBD+5HkCdQ3pABljRA6aIj3oA\n4puvuc3sFWxzs2/jjurXp2PbgGTw1uaF4mx+hIvYUZpPsnqDpAccCH5gkYh+fYRJ8dvpUYyzrAdQ\nxPK4Kfx/ARSk9V4KMxp4AAAAAElFTkSuQmCC\n",
       "prompt_number": 16,
       "text": [
        "\u23a1\u23a1k^0_q=0,p=0  k^0_q=0,p=1\u23a4, \u23a1k^1_q=0,p=0  k^1_q=0,p=1\u23a4, \u23a1k^2_q=0,p=0  k^2_q=0\n",
        "\u23a2\u23a2                        \u23a5  \u23a2                        \u23a5  \u23a2                    \n",
        "\u23a3\u23a3k^0_q=1,p=0  k^0_q=1,p=1\u23a6  \u23a3k^1_q=1,p=0  k^1_q=1,p=1\u23a6  \u23a3k^2_q=1,p=0  k^2_q=1\n",
        "\n",
        ",p=1\u23a4\u23a4\n",
        "    \u23a5\u23a5\n",
        ",p=1\u23a6\u23a6"
       ]
      }
     ],
     "prompt_number": 16
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "control_vectors"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\begin{bmatrix}\\left[\\begin{smallmatrix}{}m^{0}_{m q=0}\\\\m^{0}_{m q=1}\\end{smallmatrix}\\right], & \\left[\\begin{smallmatrix}{}m^{1}_{m q=0}\\\\m^{1}_{m q=1}\\end{smallmatrix}\\right], & \\left[\\begin{smallmatrix}{}m^{2}_{m q=0}\\\\m^{2}_{m q=1}\\end{smallmatrix}\\right]\\end{bmatrix}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAATgAAAAqBAMAAADVOFp9AAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEImYyq92J\nVM0hoUD5AAAC6ElEQVRYCdWYP2gUQRTG3x5J7s/exSNqI4G7xEIQhGBhLAyeEBG7FB6ChaaKoCkO\nC1sFQRRTXJNeSw+FsxILIdrEwiJgqYWd2EUEESz0dnZG3s6+NzMs+wSvmNmZed/Ob9/M7N13APYn\nevfG7gpqN9pBYVZQ/f0Fq8fZnGofco4zgzNvB8yIs3sKbjrHrcEDUCmUg7gQ3Cx0HLPNLR4FUIWG\nvAXx0OIFOLm4kOvLyoCCC5DBx/yNlSwp5pIxVeggEg6AhMP3peACZNEXfBNzrWYj4JhllYI7Qq4q\nB8ccCCG4aHneZAvXObh7D/tL/YsQnd7FYeaahdOy+rUVE4prnyz+/QuHm+sc3P2fcLx9xQzbNTuL\nnAztuerhMfyAMzaUaXNwgjIEB421ehdeG5i/9fyW2rIcnKAMw1WG073qfu4ddy7qJaAsnJwMw92Z\nJKHW20lI8Gc/5WLh5GQY7gbMDqrP2rWVB6c+DaHW322ORqPn1S68SlhZODkZhtPpal2FE60d+ABd\n1eGFk5MRcPA0Gjf34Gu9l07rW1YNJyCj4MatpbgFL5pralnBdyAMXPkyAi56PL3euAybj/S5Pbg1\nSOZn95yGE5ARcCYTFXOR1j44E12izAF3ds/Mp+pQuBJlDrgMmn9ZrXDTDH0mE69rJfvf4Bj35U0B\n7b58MsZ9cZkr+GOTcV8+OMZ9cXD/9mc64744OJfB6XzObFzsi0j3pc+RS0a5Ly0jDoQLLv5WCM4l\nK9F9XS8EBw5Zie7rbjE4Xlai+4qIWQLcl0NWovuqZRKX/lEQ4L4KytCBELRR2UdSzxQyG4LjbBQ8\nUfdm36aMaStBhuEYG7V5zA0nJ8NwjI2Cl244ORmGo20U+ODkZBhOb1vLfXnh5GQEnG2jfJnTcAIy\nCs6yUaFw5csWoLN6HtJCp8C2UdvfJ2bi0mr2Kz+NVVopWTLlzMbkW1kVZn10nbFRtzfWreGkKSqj\np9QUlo0i2MiuMmV/ACadTQRg0XHNAAAAAElFTkSuQmCC\n",
       "prompt_number": 17,
       "text": [
        "\u23a1\u23a1m^0_m_q=0\u23a4, \u23a1m^1_m_q=0\u23a4, \u23a1m^2_m_q=0\u23a4\u23a4\n",
        "\u23a2\u23a2         \u23a5  \u23a2         \u23a5  \u23a2         \u23a5\u23a5\n",
        "\u23a3\u23a3m^0_m_q=1\u23a6  \u23a3m^1_m_q=1\u23a6  \u23a3m^2_m_q=1\u23a6\u23a6"
       ]
      }
     ],
     "prompt_number": 17
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "control_star_vectors"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\begin{bmatrix}\\left[\\begin{smallmatrix}{}m^{*0}_{q=0}\\\\m^{*0}_{q=1}\\end{smallmatrix}\\right], & \\left[\\begin{smallmatrix}{}m^{*1}_{q=0}\\\\m^{*1}_{q=1}\\end{smallmatrix}\\right], & \\left[\\begin{smallmatrix}{}m^{*2}_{q=0}\\\\m^{*2}_{q=1}\\end{smallmatrix}\\right]\\end{bmatrix}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAARIAAAAqBAMAAACNR8lkAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEImYyq92J\nVM0hoUD5AAAEdUlEQVRYCc1YO2hUQRS9k5jdjZvoww9YiI5rIVoFKwvFFWJhZcBEQQXTuIIIrhYW\nIhjxA6LFItpZJIUIQSFWoiAmoihWAUst7MQuooiNn/ndyX0zd95bY+MU87v33DmZOzPvbAC4It7M\ncdMFc1MFNs5092Wbm47mlmVrornCiVNbCs2xsVmZiSeZmZXQkzHTBVOPC2yc6aP4Ek6vamwGMBVa\nahmchXrHDnsbDYaTaDSG0B3gkup6Jjwgt4ABDM64ANsbmwB0tUpPmMpZBl+dJEyA3R0hnbNuPhyc\nJEx4AF3AAja0fQjFBIBhAu/vA80OlychNdiVwacqqN+TLphogCBHPMWkcvxCRk9sKZPWobm/Y6IB\nq8/M4l+itkMVuieXr48Nje0DqACIHa+9X5JJbcedA32jHYBxgDtf5wsBJjsU8Pz3jAdETK79gG3Z\nYW/HTprJ6Ymeqb5J9PMtB7BMEoCQSXXtDHyHnT4gdrjAQirrmSOwcrx3VvXyhQMYJilAyASWj/RL\neLYY1T20XGAhtds9uAX1xazA+puZnuYA9u4EAHzJIyY9nb5mdaGjg5niji0XWEjtsRVG4WK/7tmy\nRzR1hwNYJgEA70XE5KLalVpzVgczxV1lLrCQyqO6oHblyoDxNdWCvQUcwDAJAfhWRExOwop29UEG\nV4+2305PT8+7h5YLLKRaujIJn+D8HNTGXg8o/4dVCU80Iw5gmBAAXSC+xTqKKoNDPbZTxsR6AbwD\nabqlTEIAflOiPXGOKzq3ba8sOxj4c3/TdsuyEwKS2XGO9ewVmM0rPbEO8GhgxGQHSk/sIoAukMxO\n5Y3Jtn9oubQL6YLq5tQNd99W32zrMQewd0dbVUEAvuSp7IBLu0XxgYVEq27ducKpUiYhIMmk9yfG\nNC0XWEjqsnuejsr3BAJAkkkubDd70gUgl53AP3lOQr/yPQkQHGDJTPCLoNfgAgupLbRMLQ44QMSE\navui7OAXQYfnAgupLaRQbc8BIiZU2xcxwTdHr8UFFlJbaPkb9ahxVNsnmeS0PTLZ+JGuKyQdXVKD\nkIn4RT1ykt3+GPDaXgtHVah6dFBe29dzP0+EdM664bX984y45JkYwNK1/QkaWEgy4rV9xby36JY7\nJ/+o7ScwqG6FVBUWXtvX0GzaHJN/0/ZigkYWUo2oVI+1/Trqb7NDAUvX9vk/0TJJSHV72SYYJglA\neGL/Y20Py83xT78ngVQHmDL7wAHsOQkBdoH4FofavvLCHH8usMlOpO3dQ8sBLJNA27sFYiahtod6\nMZNQquPzlmQSAewCMZNQ25cxIVLdqkf30CaZEIBRj26BmIk77F7blzHBy4HavnRPIkBqT5yj1/bd\nMvHavmxPkIkHlDDx2r5bJl7bd8tEAbjsbBzeC2AqR9lr+/5juwD6hvebu4x/j23Ft+ERP4NS3f4T\nhQfQBby2NwvA6LD6uOqq0lLfNlP50NL31I/fVouMfLfVMvfKjnPangfkF8gBzrXGAUzlo2Mn0PY4\nnW4DqZ52RAsL+AN8jEJslCydcAAAAABJRU5ErkJggg==\n",
       "prompt_number": 18,
       "text": [
        "\u23a1\u23a1m^*0_q=0\u23a4, \u23a1m^*1_q=0\u23a4, \u23a1m^*2_q=0\u23a4\u23a4\n",
        "\u23a2\u23a2        \u23a5  \u23a2        \u23a5  \u23a2        \u23a5\u23a5\n",
        "\u23a3\u23a3m^*0_q=1\u23a6  \u23a3m^*1_q=1\u23a6  \u23a3m^*2_q=1\u23a6\u23a6"
       ]
      }
     ],
     "prompt_number": 18
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "sensor_vectors"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\begin{bmatrix}\\left[\\begin{smallmatrix}{}s^{0}_{m p=0}\\\\s^{0}_{m p=1}\\end{smallmatrix}\\right], & \\left[\\begin{smallmatrix}{}s^{1}_{m p=0}\\\\s^{1}_{m p=1}\\end{smallmatrix}\\right], & \\left[\\begin{smallmatrix}{}s^{2}_{m p=0}\\\\s^{2}_{m p=1}\\end{smallmatrix}\\right]\\end{bmatrix}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAASQAAAAqBAMAAADoqjsbAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEIolUMqvd\nZs0eYxJ8AAAEAklEQVRYCd2YPWgUQRiGv8nv3eVOY7QIJJLxLMQuCEYiCBFORBtTmKhFNNWpYHFY\nWKWIgmiUwIkEUyhcoZUgWqTQFKZQhDQGQiCNkM5W0UrEOD87m53Zd2c3cgg6xfzcfO83z84sk31D\nZJXe9W5rnGXQtX1J4e7JLJlVzHE2kjk2COy4VduuhNroWmbNF9qXOdYEFrePtIMGPFvbU95PrFwe\nVAt0clo2K4Vta7kM9KGGEBLWyLWC8th0ttpDZbEdsuoRvzEuKlkgElELQuJKISqEhDVyLV3YZ9OL\ntuqEHCTCB9d8pD7wkKRfGhcJv95NR2JHeqO7Y/rxXbp9gXav18x8pE1Gyo3VC5eGI6FhF2n0wQlN\ncfNnGBjpxJBaJ59GpqNdlJ5xGTFUGJENKEijkZI1MaS2U7tAavkTSs+4nKk8gC+FmEEajZSsiSEV\n3v6Qq1hF3+YoPeMy8PlmQzbRwmZvyiHSaKRkTQzpES1EU6u+ft1ResZFQL7WNariIlVb9x45QhqF\n5NHEkJbpaCSz7upLAaVnXEQUG60N0Vhlp6ZBGoXk0cSQ1tYmafXJ6tV7VDq8Xs8tLi6+Dq5OlJ5x\nAVKaGKfS8Keh+QatTszQfaGpXxegYgppFJJHE0MSeWhuo6U+T/nzbEmOzG2O0jMuA0TJX6SD+Wma\nm+7X4zQkjwYidTX6aYbodGFD5085OBX0hi3l6tTWeKElaQfn0UCkIj2kpTrx3Io6OEp5vVX6pfxg\nMU8tNKsPLu319mggUh+N00Kt8/2EfuLgNvceHHvZPtl1jqbGVrSGHbsjO0ijLwExmaCBSDpraVq3\npkbpGTezul23h36kINbVeJCmRuz8GZA639mSLEgxjQfJyQ7TM+5G2WP0GOHB2aHh6F9Bwg4FPTHj\n5uGwQ0GacJcSHArcpT/4hEtwKF6kBIcCkf7Sh26CQ0FIPjvAfpmTUi3jZuixA44mPDgi4FDwt7cP\niT5aX2uZkBzNFlKTHEpHzeyLbBkXlSqeXSJbs4XUJIdSMgyqZdwMfUi2JkRqlkPZaxhUy7hs0hyK\nrdFITXQowZ9VSSIK46JKdSi2RiP9fw6FXsm9QNce43IGuA3StznS6F1K1sTuJeRQrhyQC6P0jIsJ\n5DaC2xxpFJJHE0NCDoU+iIWTkZDbCP6Tkojk0cSQgEOhFCTkNtKQPJoYktwPx6GkIUmJ61DSkDwa\niOQ6lJRdkunJcShZkBI0EMl1KFmQHIeSCQlrAqSBygli3yuj6oHJdSjPvtWpvXKmW89G61Djug31\n/yaskWupAjVnK1+JZNVRvUxUrdaiy9kOpbNajU6avqOx3QbWqLWMXrS25kZ1kkhVkRDTdR2K+d3T\nxtyGJ9ZMJWh+A5QnM2g1BTWuAAAAAElFTkSuQmCC\n",
       "prompt_number": 19,
       "text": [
        "\u23a1\u23a1s^0_m_p=0\u23a4, \u23a1s^1_m_p=0\u23a4, \u23a1s^2_m_p=0\u23a4\u23a4\n",
        "\u23a2\u23a2         \u23a5  \u23a2         \u23a5  \u23a2         \u23a5\u23a5\n",
        "\u23a3\u23a3s^0_m_p=1\u23a6  \u23a3s^1_m_p=1\u23a6  \u23a3s^2_m_p=1\u23a6\u23a6"
       ]
      }
     ],
     "prompt_number": 19
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Picking off a trio we should get our familiar equation back:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "zero = control_vectors[0] - control_star_vectors[0] + gain_matrices[0] * sensor_vectors[0]\n",
      "zero"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}k^{0}_{q=0,p=0} s^{0}_{m p=0} + k^{0}_{q=0,p=1} s^{0}_{m p=1} - m^{*0}_{q=0} + m^{0}_{m q=0}\\\\k^{0}_{q=1,p=0} s^{0}_{m p=0} + k^{0}_{q=1,p=1} s^{0}_{m p=1} - m^{*0}_{q=1} + m^{0}_{m q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXcAAAAqBAMAAAC6iQZQAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAEO+Zu3ZEIondMqvN\nZlSyKd/8AAAFc0lEQVRYCe1YTYgcRRR+8/+302mjgsFDxhVRUGQuBhGMHdhccoiLZGIwro4H10hA\nBk/mlFFZAsqGEfwBD9KE5CB6SA6L4ooGL4JMYGFZWQVDS4jGg5AQ8CBoUn+v6nVN1bbjKLiydZh+\nP99775ua6qqpB6DHHeuhlh1CbvFVh1WbMqIhyy8SeWp4zLo2E/bkIqracjG8zTZRPSM6K7tM5anh\nMdPyV+EuqtryNshv9MtkREOWX5Tz1PCYAbZP3w1Qug+g0oL3bcJUfxlqMdWJ/EQ/MzorO8vWDMFT\nw2UuTE+HsF2QePfvki/z8Fr3HyAfLBxNke8KYuLDRR74QhDkc9cYKOOHdf92RZ59J1tPGdGZfoBf\nPgFaI+KZ5aBmtBny5QGzZbxy7rdGkH8xOzorO0D5+VdCWiMyNKnZWHHmq+0zZ+HW9Z7xWFK1k+Qe\nPWEZuSrIv14/vGE0LDyzgX/1+067s4+xB6A1IpZduaiZF5UDyU9dfI4vHP/YVY/cTkH+wk984fhH\noXva74Rv/4D7w8MjgIhZPC6JRfLbHhmJTRtmfvSwE+Svx2m0rRX332KbjF65fQC/wcPGoKSI7SEe\nl0Qg+VMHEhXCHvZ5xs/Gj27EBiAlDvt8+NVwuFK/Zz9x2kdpI4T6l78TgBLPQTDDxh5ozNZb8IEB\n8MSN4fCB4YURlwQpgkj+QL0V62j79WAvctBrzGq/EiSMz3wpOlU3Xuu9L7/Wg59hyfiVdOReNOXj\nUlS5GqMKqn7EDLZLYhQAye+tD/rSwT7tjYltgrW4EGu/EiSMky/Oflg1XnvPrPXY4edYlu9hyHE2\nwc2ojyrWj5jBdkmMIojkv4alGNYO9d5cXl5O5JHQfGg9qTL1C342NucOQnP3D7uuxLA6dxIIjJOf\n6jWegmbnhIYzm4YBI7+21oXVX1dfegN0VgBN/ihLUPk0tOpDxLIolyrNSpDK+pBiOAjaef4AdcwF\nh3IDoXLyXAiehgeDPlzu7xBmBePkxbgELfFUcA3j5Pm4fD6fXAGdlZAX7pH6gjy6ZGksoQ5cnHkO\nmopPS6z6VfbVz0sd18FnuUE1gWL8MYUVpALwDW6mEq5hSL4R74CTADqrmXmVwa7fxszsKUtjCWvZ\ncFwtXJDrgb8P7AVsVVfEOtBn4yBo1wLIw6KBkfxvVWcNPDYwJF+DSzBITNYR8qn6JDEXZWlWQiwb\n64XlgPKiWB3ymDseVt6Z41Y21MmbO1vqNp6EY50VabaOvSPfxdLO4JVrBlZ/drew3wkHYalnsp65\nnkg8fqbqo1E+VWksoSrTZcPm2oQE0OwbjUrrVKGyemOEKQEPzJdVRLVoOpdMS5g/ZgJZ+NMEBHAs\nMhqRKm8TJSU+lhg18cE8WUUkrW9SUYmWsMhT2GaQ08tmMzAmHFPk7T8lBMdF+y+P5c6I/he6Byny\n1p8Sixz+5bDNqGdE6w0X8c6n2gRtn9ucIo+HkR2qdHU0eLyZ17yM7DKtp4bbjOQ3YfeA7DZb3QPf\nihL27GWz1T1wTeBW98A1K8zmaRFEfpdMhLvN/6d7AOy+T8eY3QM4R4OBZxuve8AjNuweKII486nu\nAb/v0zFm98C0BUSS8bsHqn7Ewt3dAwVA8qnuAd598AuwTW6c7oF9RRq7e6Dqc/Lu7oECIPlU90D6\n9D1/zO6BullP0D0g5N3dA4s8n2TdPZA+fc8fs3ugyE/QPSDkOS82VOMCuwcO8vr2rnz6no9n41/s\nHshlM0H3QNX3dw8c5PXtnfsm6B5w8vEk3QPFTU66/Ex3Dwz5nTN7JUDf3vl9f4LuAWsLTNY9wH6D\nYW91DwSgNPN4COX5FxDVQoGtsv9896AyP2/oAtDb+2boHgDcBM7Haf5wIcMkAAAAAElFTkSuQmCC\n",
       "prompt_number": 20,
       "text": [
        "\u23a1k^0_q=0,p=0\u22c5s^0_m_p=0 + k^0_q=0,p=1\u22c5s^0_m_p=1 - m^*0_q=0 + m^0_m_q=0\u23a4\n",
        "\u23a2                                                                    \u23a5\n",
        "\u23a3k^0_q=1,p=0\u22c5s^0_m_p=0 + k^0_q=1,p=1\u22c5s^0_m_p=1 - m^*0_q=1 + m^0_m_q=1\u23a6"
       ]
      }
     ],
     "prompt_number": 20
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The larger A matrix for all of the time steps can be formed by stacking the A matrices from a single time step diagonally:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "Am = zeros(n * q, n * (q + q * p))\n",
      "for time_step in range(n):\n",
      "    An = zeros(q, q * p)\n",
      "    for row in range(q):\n",
      "        An[row, row * p:(row + 1) * p] =  sensor_vectors[time_step].T\n",
      "    An = An.row_join(-eye(q))\n",
      "    num_rows, num_cols = An.shape\n",
      "    Am[time_step * num_rows:(time_step + 1) * num_rows, time_step * num_cols:(time_step + 1) * num_cols] = An"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 21
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "Am"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}s^{0}_{m p=0} & s^{0}_{m p=1} & 0 & 0 & -1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & s^{0}_{m p=0} & s^{0}_{m p=1} & 0 & -1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & s^{1}_{m p=0} & s^{1}_{m p=1} & 0 & 0 & -1 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & s^{1}_{m p=0} & s^{1}_{m p=1} & 0 & -1 & 0 & 0 & 0 & 0 & 0 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & s^{2}_{m p=0} & s^{2}_{m p=1} & 0 & 0 & -1 & 0\\\\0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & s^{2}_{m p=0} & s^{2}_{m p=1} & 0 & -1\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAywAAACECAMAAABF91+xAAAAP1BMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADFBd4eAAAAFHRS\nTlMAMquZdlQQQO0wRIm7It1m7818bInf38oAABKLSURBVHgB7V1rY+s6CnSax7nbtM25u/n/v3Ul\nu3bBwjAiVaK05EulwMCAzcnLnjMMwuN0vpxPwvPqUy6QGjGM0YH+O/B2Gk5vtTRdoNok4R8d6KwD\n18Tnb+1LiwvUWeFBJzpQ2YH3vwnwsatDuUB1KcI7OvDYDuyu4+OFsHCd9y4QSToMh/PhzJ5ANpWg\nnfDuEg2B+lHaHkzbRiiMFBOtycUPjX1zIpidwmht+phmZBh21/0hPY6UpesdlQtEsl7SpBzpzBLb\n5rIOtNu/fBSh0BCoH03gwQxNQUpwxURrcvFDY9+cCGanMCpMlzwhL+kE313L91uuz+ouEGnO+L7v\nWvlJqRZ0LIcFDYH6kZKm97KtaxoTouwUP8VEa3IVhca+ORHMTmEkmy4bw3LaX87vjDewcYG+4p6u\nOeNf9gr3Zd1YVYPKYUFDoH6UqQczNAUpwRUTrcnFD419cyKYncJow1QOy3n1hoyx39q4QOtgu/Ef\n4I/L+nl1Xw0qhwUNgfpRvh5MeqnPr66NGqEEV0y0Jhc/NPbNiWB2CqMNUzEs+9dB+gjMaig2LlAR\n5ZPioTBoT1SDNofFzFudKvH2YNqCFEaKiR0C1I+CPJjeGlEMy9vb8VT5qWEYXCDaynG9m96GmSct\nA1aDpGEZ3/6ZeatTJaIeTFuQwkgx3dbzH9KIYlh2L1fh8z7rVd7wK1tA0DC8aHN4GhNfX4tk2hM2\n6Pg2P8bk5bDYISYCqB+l68Gkd935G5dvb4RZBZoX9ftpjVgPyzH9I/ti/hubXkvo5TAo6PX8V/3S\n4CPPSfU3R5WgclgGNC/qR08SDwYmRBPBIIWRYnKloiA0NsXANblACiPZtB6W/J3ZvvwmmZHJG/aT\nCgpKH1rVYTnsHb+z1IKO+eIE/kBDoH40ugczNAUpwRUTrcnFD419cyKYncJINq2H5XK4HPI/1Zf9\n8XLZvx7PeXAu+8PlnL7RPb9Mj+PAf6wHQSmSPizDeX/ea2/UWCPnTRXoff9xLV850RCo30wt//Vg\n2oIURoqJ1uTih8a+ORHMTmEkmtbDMjM9Hv4df/I45stPjukilBP7LY8PCwhKbsawzIHib3Sgww5s\nDcvwshvnI7+gpE8ou2GXhoc82Nuw5XkDFMOydCoWT9iBzWFJLyTH9BHiY8jTkt7mn9NvhV9vw/gH\n/KVuAxTDsnQqFk/Yga1hOaXrGfOHl5dD+kj+/nE4rH5XF69ssUC7w3VfdzHLE3Y0KP/YDmwNCys4\nX4NZ/XCBqrMEIDpwtw4gw/L+IdwCYjF0gaygYY8OPLADyLA8kF6kjg7004EYln6ORTDpvAOPHxZ+\nlRnWLg8Gixxe0YHNDjx+WNhVZps8ucGD4RFiFx2o7sDjh0X+eVMvxIPRI4Y1OmB2ABmWtdKFGTQ7\noCD5whk9hQeDMDLvekOLovQ9GLh7NBEMUhgpJlcqCkJjUwxckwukMBJNwLAUSheM2MYGBnlOfA/G\n1vwQhV9YeXBRBOXB2FxJgmWJZlL8FNOSJi9QPwryYFyJYJDCSDYBwyIrXdA+CGsc5HlL5cEAjIRb\nXVhpQAjmnzceTFuQwkgxscJQPwryYDprhD0sG0oXtA/lugLk+bDuwCCMjGFBQqxb4cHA+iQsGZpJ\n8VNMrlQUhMammN4aYQ/Lp9LA6towVlO5gUCTJIx4lVkZ8fMZD2aEIoyMYUFCrIl7MLNOQ4OWJ3oK\nI8XE6kL9KMiD0bjS2HyNZlL8NkzwsAC3GhPOn9lUkEcSxoOZaCGMsGFRiyI9GJdI2jVmPkeqEsEg\nhZFiYhxRPwryYOCaaCIYpDDaMCHDkm8F/lt75GyQRxLGg5k6iSiXmMNiF8UOW9ogadeYtiCFkWJi\nHFE/CvJgemuEPSweKQ9InwSThOE/1mOYdJAKGRmxDFP4hR5uqCgGSBsx7dppvW8KUoIrJkYR9aMg\nD6Zt9xRGGyZ7WFopbICSMOyzPIgZJBkZWbCDHs/BeGVxdQJIyziMm6YgJbhiYiRRPwryYFwdh0EK\nI9kEDIusdEH7IKxtUP4qEdCRYd8Sg5jEp7zZ32Y0CMIvrDQgBPPPGw+mLUhhpJhYYagfBXkwnTUC\nGBZYLIO2xgZhkjD890cMk3mUw2IykoVf6opi3tNG1AkR/NhTTUFKcMV0Mz809s2JzCM9Z1AYiSZk\nWObYbf7qOjJ8WGYGOiZ7CcMyg+NvdMDXgccPi6Ejw96GLTWGjMzSiljcrQMdDIsuCcM+4C9t0THJ\nLV5Zll7F4rs68PhhMSRhxB/4DcwQMjLfdX5EHNKBPCynY7VkKonwrct8tWftw4OpzRH+0YEkCZZ0\nvKT/U/IxvfFIwngwj6kusj55Bx7/NuzJGxj0f08HYlh+z7GOSm/sQAyL2cDiMjMTkRxcICRw+Dyu\nAzEsRu+ly8wMyCBem2aCwqH3DsSwmEfI9ZONC2RSCYeHdgAZFlHpwmLdH0hhpJicv2/eOCwqoa3W\nV4IkLRs0BOpHqXowoe5CO8jWsqYGcyk3KEjxU0wpn+u8d4GW4nRCixtf1IFELRs0BOpHCXowsFAL\nTQSDFEayCXhl6U+Wg7UGpaf4KaaUynXeu0BLXTqhxY0vakHC3TtoCNSPEvRgQt2FdpCtmwqAKMEV\nU+bnOu9doLkbBqHZjf+tBpXDgoZA/ShBDybdKTneyF33H2ChIMVvw2S/snzevN9GaoS2M11JkC+7\n+WiTSQm+YZqEZCovM3OBWBfSZoPQ2o3vq0HlsKAhUD9K0INp2wiF0YYJHpZawYrxvO8I9Fm/xEg2\nuYRkXCB6To1rmVDhxp+oBm0Oi9Qjlqs6VUJ7MG1BCqMNUx6WP9c/rBd8058sB+OH0lP8ZJNLSMYF\nYvXkjUyocONPVIOkYcHka6pT3a2m3BKUneK3YfpvGhb9QsoNpQt+pNa7/kAKI9kECsl8l/oM7aBM\niHoIaxtkatnYIaa8qB9l6cGEugvt4Gota2qsnNZbFKT4SSZUSIbdnIaCrMsCJELruot9Lah8ZfkO\nXZSC1vxELb0R1xSkBJdN9meWtgobcyvHv00FQJTgkil/0wmIzwzstmcUZH3FJhFinZI2tSBBywYN\ngfpRmh5M23NPYSSbgGGBxTJoazoEiYIdE2fBBArJcEENEJSSGl8tC4RYd8VNFUjWskFDoH6UpwfT\n9jRSGIkmZFhoxb9ubQjJ8GGZu2OAkpsxLHOg+NtTB2JYrKNhCMmwt2FLLAMUw7J06pkWMSzW0TKE\nZNgH/CWWAYphWTr1TIsYFuNoWUIyoT5jNPAHmWNYag5mvhi1+uECVWcJQPsOxLBU9NglJOMCVZAK\n17t1IIblbq2ORM/egTwsp0M3InvP3s7g/5M78J4uMdWvDfvJ1XdWm0cSxoPprOynoRNvw7o5VNYF\nYxJRD0aKE88hHYhhQbp0Jx/Pz/oezJ3K+XFpkGHpT5aDHQaUnuKnmFypKAiNnTCeE9+BsRlJwi+0\nJpfoip2WpZg2TUFKcNEEDIusdCEURp/qD6QwUky0Jlg2hILQ2BnjOPEdGJORKPxCa2reiCWZyXXx\nJAsUpPjJJmBY+pPlII1JZ1i6ln4Y795nTxcbxU8xsTCoHwXVYO4zLAAj4VYXWhPccwoC0lL3ad0U\npASXTfawbChdlIXRZ/oDKYwUE63JpTWCxh4T3WVYEEbGsCAhWOfSxoNpC1IYbZjsYfm8eb+N5gpr\nadNMSnDFdDM/KLZHEsaDGatBGBnDgoRgnUsbD6YtSGG0YYKHxZT8YN35zNYRSGGkmG4uContkYTx\nYKZiEEbYsFQdXSQt63beNAUpwTdMyLBgkh+s0g15DOZTbJqClOCKiXFE/SgIwXgkYTyYiRfCyByW\n+lMCSUsbN66bgpTgGyZ7WPqT5WA9RekpforJlYqCkNgeHRkQI/xHMSIjU/iF1uQSXRHTsqjCpilI\nCb5hsocFlvxgxcryGMyl3DQFKcEVEyOJ+lGQjUElYehtZihG+oXfZjQYryyuUwJIS/s2rZuClOCy\nKQ+LLrLXVmGD9UfW1GAu5QYFKX6KieVD/SjIxuSvKWt1ZFCM9NuNzWgQhF9oTa5TAkjLcoybpiAl\nuGyyRfaGtgobrEGipgbzEDYoSPFTTCwh6kdBJgaUhGHSGCAm8RC+jrYYycIvtCbXKWGlZQnmTVOQ\nElw0AW/DZt7xt20HDEkYNiwzEwOT3IRhmcHxt7IDMSyVDWvobkjCiDoyBiaG5TuPVwzLd3bztliG\nJAz9gL8kMjAxLEunvmERw/INTfyeEB4dGQtT+b/LfE8hPzZKDEufhzZf9lr78GBqc/xq/xiWLg+/\nRxLGg+my+G5JxbB0e2iCWG8diGHp7YgEn247EMPS7aHpjNjufH7L109WPVygqgz3dI5huWe3nznX\n2zC8pi+q6x4uUF2KO3rnYfnzz3/umDFSPWcH/s13O9aqMbpA3fbnf//YInui0oVVUX8ghZFiYnWi\nfhTkwbjUU2CQwkgxpZqsi5Fp2cvaBVrQOqHFjS8qQZKWjRgCeBsmK11wesWuP5DCSDGxwlA/CvJg\nXOopMEhhpJhSTad/qz+zOEFLA3VCixtf1IFELRs5BDAsstIF51fs+gMpjBQTKwz1oyAPxqWeAoMU\nRoopnfYepVgX6Kt/KqEvN76qBQmvfXIIe1g2lC44v/WuP5DCSDGxulA/CvJgehM1meo57d9P1ZcV\nuEBf/btP98ph2chrD8vnzfuh7jLLJ1R1or/uKYw2TJOQzMf1eq34NswF+hqTabVBaO3G99Wgclg2\nQsDDUiXlMZ9XHYE+65cYKSZ2GFA/CvJg2nZPYSSbXEIyLhBt3LiWCRVu/Ilq0OawrE8WZFjqpTzS\n4e4NpDBSTOwwoH4U5MG07Z7CSDa5hGRcINq4cS0TKtz4E9UgaVjE09celg2lC05wvesPpDBSTKwu\n1I+CPBiXegoMUhjJJlhIhtbtAtEA41omVLjxJ2yQqWWzEcIeFpeUR4cgWbBj7LNiYscB9aMgD6Zt\n9xRGkgkVkjmdL+flN0sUJEg10ea1bcScqXxl2cgLDIusdDFn2vjbH0hhpJhYeagfBXkwLvUUGKQw\nkkz5W1REfIbdyImCJKmmu3dP0LKRGjEMwLC4pDw6BImCHdORUUz00LmKQmPfnAhmpzASTLOQDOMn\nbJhEAAoyb3sWCAm5V09VgWQtGzEEMiwrKrH9lR247I+Xy/71eM7/xcdlf7icj0kl62V6HAdRfCb5\nqaAU6ZnUZ2JYfuWZ7yjaUF2Sh8UAxbA4DkRAnqADhuoSexu2lGOA4pVl6VQsflIHDNUl9gF/qdsA\nxbAsnYrFD+qApbp02l/OxVXJFui5pJriM8sPOp/vV0q+hL364QJVZ2kIiGFp2NwfG9qluuQCddXC\nGJauDkeQ6bkDMSw9H53fyc0lCeMCVfY3hqWyYeHevAMuSRgXqLKUGJbKhoV78w64JGFcoMpSYlgq\nGxbu9+iAcCGwndYFssN+eSDDIsrCfIWQV/2BFEaKiVWH+lGQBwOrGtFEMEhhpJhcqSgIjZ0wd9KR\nsRlxlSRgWGRZGNoGYd0fSGGkmFhpqB8FeTCwqhFNBIMURorJlYqC0NgJ45KEqQeZjNYqScCwyLIw\ntA/Cuj+QwkgxsdJQPwryYGBVI5oIBimMFJMrFQWhsdOs3ElHBmDE39nZw7IhC0P7UK77AymMFBOr\nDPWjIA+mTykkUpWnqApMpY7MRKwehDCqHZZPsYwqAaBZn6QjkFKGYiJnSJa5zbfNflQV5cG4EsEg\nhZFiuk8jJv0klsveuEApLFKtc1jWsjB6EZ9EOgIpjBQTqxL1oyAPZj6MbbqnMFJMtCYXPyS2Sz/J\nBcrlIIzqh0WUhWHdKzbVcjQT+4aZFEaKiRWG+lGQB5MO4+9shEs/yQVCT7faYdmQhaFnRLnuD6Qw\nUkysMtSPgjwYWNWIJoJBCiPF5EpFQUhsUD+J6cgMIKjUkREZqSpJ9gf8DVkY2ghhLYnqCG78qaYg\nJbhiYgRRPwryYNq2XGGkmGhNLn52bFQ/id1mhoIEHRmb0eo/2QCGRZaFYc0rN/2BFEaKiVWG+lGQ\nBwOrGtFEMEhhpJhcqSjIjp2/ykVEl9gNzChIuCfTZjRwlSRgWGCJHdqaDkGius3EWTHdXBQa++ZE\ncMsVRorpZn5m7Fk/yZCE4dIYICixL3VkLEZrlSRkWFibYhMdaNwBQxKGD8vMxQAlt3JYZiz6dxqW\n9L8JXK8vKCb8ogNNO2BIwrC3YQsRA3TjsOTfPNMjXV5wGB+vS95YRAce2QFDEoZ9wF94GqAbh+Uy\nzciSLRbRgS46YEnC/AIdmS6OQ5B4tg7kS4OrHy5QdZYARAe66oBLEsYFwsr+P9uzDyDJ0+lIAAAA\nAElFTkSuQmCC\n",
       "prompt_number": 22,
       "text": [
        "\u23a1s^0_m_p=0  s^0_m_p=1      0          0      -1  0       0          0         \n",
        "\u23a2                                                                             \n",
        "\u23a2    0          0      s^0_m_p=0  s^0_m_p=1  0   -1      0          0         \n",
        "\u23a2                                                                             \n",
        "\u23a2    0          0          0          0      0   0   s^1_m_p=0  s^1_m_p=1     \n",
        "\u23a2                                                                             \n",
        "\u23a2    0          0          0          0      0   0       0          0      s^1\n",
        "\u23a2                                                                             \n",
        "\u23a2    0          0          0          0      0   0       0          0         \n",
        "\u23a2                                                                             \n",
        "\u23a3    0          0          0          0      0   0       0          0         \n",
        "\n",
        " 0          0      0   0       0          0          0          0      0   0 \u23a4\n",
        "                                                                             \u23a5\n",
        " 0          0      0   0       0          0          0          0      0   0 \u23a5\n",
        "                                                                             \u23a5\n",
        " 0          0      -1  0       0          0          0          0      0   0 \u23a5\n",
        "                                                                             \u23a5\n",
        "_m_p=0  s^1_m_p=1  0   -1      0          0          0          0      0   0 \u23a5\n",
        "                                                                             \u23a5\n",
        " 0          0      0   0   s^2_m_p=0  s^2_m_p=1      0          0      -1  0 \u23a5\n",
        "                                                                             \u23a5\n",
        " 0          0      0   0       0          0      s^2_m_p=0  s^2_m_p=1  0   -1\u23a6"
       ]
      }
     ],
     "prompt_number": 22
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The control vector and the gain vector are simple column stacks:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "bm = -control_vectors[0]\n",
      "for control in control_vectors[1:]:\n",
      "    bm = bm.col_join(-control)\n",
      "bm"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}- m^{0}_{m q=0}\\\\- m^{0}_{m q=1}\\\\- m^{1}_{m q=0}\\\\- m^{1}_{m q=1}\\\\- m^{2}_{m q=0}\\\\- m^{2}_{m q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAFkAAACDCAMAAADRcEfUAAAANlBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABHL6OuAAAAEXRSTlMAMquZdlQQ\nQN0iRLtm74nNfDyb/jsAAAOQSURBVGgF7Vpht6QgCKUy31SOs/3/P7vYO1NChTDuntPbzS+jDt0I\nwQsZgLL5GKJXytrEJg9+sl2ilJ5R7vk3lB6eiPxolGosYs28tLZwjR758Q0I0Mydw9YXkEFtjZDg\nWhRvZt0j2lYwSMihC41zHTg3pXv7LsSh9GDr/zLyOAd0h9aDKy3CCrh2JOTBuwcKziNAF9crtB2O\n3E/vlhy3dbgKi7MVlxeNRUOUI1ONnmheh2E3zFBeZLbAIjICfusdJiibgzmliDymdXuhOwxt2opC\n14fQjX1c9EenRU+J7XfrgQeSiExNA717oaf00Cf9+w77uUANMrQNePSWmJbz1UBDNz6LNXKVUh9R\nk6oPQOiHh4ies1kDLCvIkD2a3aFztw5NP8bA9j0WohY7sxtRM/M/oQK5m8XwqUDeaUknroHMNgaq\n4m5k0Zm51Q6LTliQWShQoN3IgMzDd4dFJy6BrOfuRXmDznxjoA+/G4nIN3e/qRvDwsjdKJ4nqqKd\nwcTdY3zmGZSIbONu5Bk1som70Xh6ZOajMnfXIBe4u0LnAndXIBe4u3FzlzGj6BvMznxIUiT+58Ld\nX/PXbl4xUeDuX0gU2jpFcbdMpMYaGcxB14ZMNoYDtHzKgsw2hhzmoG9BZg57gJZP3cgXsAbjbrYx\n5Aoe9MUVDLV1t3c5M263r6u7h1RZZ+837ro7t8ZmZeyZuJun16JvkNtg/XrX3ZlFCtx9192ZrXjX\n4nX8WnmckPWZjIW7LZnMzd10lUi5QP/ajWy+8QOQL8zdZ5xSx918BSu4m3kHR6Z/m+rug4r+dN+w\n1d0s+MV9g3E3fZ6DEXF32Rr0avltPMp+jCxnBTXId0VP15COyhX9WXRTnN2oUNFbvG6HLU5YkJsY\np/ztlghsOpvAI58xHRjqWtL5rALiCHhq5efjYomL4nipgEhUHghtU71e5ySq9w2fjvSUzbKC4C3Z\nlwXZd4MvH22+H8mCnA719XYWkVlW8NZG9ysj11b0Z75xmayAHZiJ1rC9jWcHZiKyLStgB2YiMssK\nZO6+K/r7JH3dNX4Ed6/aqjpipDCEK2UFZ7soUxlsWYHFGngnQ1ZgQ/6vswLG3X8yK/jnufug7j6P\nFFNFz+pu0Z9t3G2pjk3cjRFKMnxRZ7ZxXKvuVn09i08rfgW3WWP9etanj11d+r5NboW38duB2fL1\nrMPTpU+bXHd/iorXFbgbJX4D5s49FG12xmIAAAAASUVORK5CYII=\n",
       "prompt_number": 23,
       "text": [
        "\u23a1-m^0_m_q=0\u23a4\n",
        "\u23a2          \u23a5\n",
        "\u23a2-m^0_m_q=1\u23a5\n",
        "\u23a2          \u23a5\n",
        "\u23a2-m^1_m_q=0\u23a5\n",
        "\u23a2          \u23a5\n",
        "\u23a2-m^1_m_q=1\u23a5\n",
        "\u23a2          \u23a5\n",
        "\u23a2-m^2_m_q=0\u23a5\n",
        "\u23a2          \u23a5\n",
        "\u23a3-m^2_m_q=1\u23a6"
       ]
      }
     ],
     "prompt_number": 23
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "xm = gain_matrices[0].reshape(q * p, 1).col_join(control_star_vectors[0])\n",
      "for gain, control_star in zip(gain_matrices[1:], control_star_vectors[1:]):\n",
      "    xm = xm.col_join(gain.reshape(q * p, 1)).col_join(control_star)\n",
      "xm"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}k^{0}_{q=0,p=0}\\\\k^{0}_{q=0,p=1}\\\\k^{0}_{q=1,p=0}\\\\k^{0}_{q=1,p=1}\\\\m^{*0}_{q=0}\\\\m^{*0}_{q=1}\\\\k^{1}_{q=0,p=0}\\\\k^{1}_{q=0,p=1}\\\\k^{1}_{q=1,p=0}\\\\k^{1}_{q=1,p=1}\\\\m^{*1}_{q=0}\\\\m^{*1}_{q=1}\\\\k^{2}_{q=0,p=0}\\\\k^{2}_{q=0,p=1}\\\\k^{2}_{q=1,p=0}\\\\k^{2}_{q=1,p=1}\\\\m^{*2}_{q=0}\\\\m^{*2}_{q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAFgAAAGQCAMAAAAdh1drAAAAP1BMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADFBd4eAAAAFHRS\nTlMAMquZdlQQQO0wRO+J3c1muyJ8bBZ/NCcAAAqQSURBVHgB7V3tYqM4DHQamt1LSEn28v7PeiJg\nG42Q/IH32qbwY2uwPQhZ8kjGYZ1bO0790J/WKgquHR7P4411uZ3c6cauZJ9cJzznDo9jR8c76/mg\ns0udyMOI9kYAh8chYnaXfjy5X+if6+J6bJFXGhDYXc9jz/bAp8f0/BtU8XwkIfH5Oj3qhsFbB+6P\n7tCRgZyOQ3+f7lH1r5D4Y3g/e3VUIc6dBPDjuMEUFpIg8OHS3bhBLxqXFBG4u7nDxa0KnXByqEbg\nW+cOV/f0EZQvYSdQjcAfZAlv3arECcuGagSOYpLHk731b9Pxzn1xOHZDT2OhVBOKCvx+pEko3gac\n/L3r3Wl2pakRzgEq8MfBHfi0yZ71RtUfy/s6Vm1IfD25vls+q/Oj8xwAeph+WK2e76ZKfO4HmDVn\nJ7+Ps9T92nWEuzxgDlCBqQ9TccQY/WdYNcjYxlCFc8fHugvS1fuVq3+JOJctiVea519SgcFD8xFT\nEnsbKEZMAYNZFuNrqkBH2g78t+kfPXS7xJ7vmg+ep3/w0O0Se/ovRoIOwipeif6TTv62jHpRFQb9\nJ+zk3F+WIRkCG/SfdPKrBRzpv5ClyShM4GA0pSy9Bvzr8SvghUIxS4PEf0hzLAfxyMUsDcA4eB7X\nFbP0oXscFySpAtMddpYOat5UUHWcnBkSt1WBEzNDAtYIsZIzQwJaSPy3WbpJXHHqFjO0Z+mtqrhT\n1M7mCs/SzQfPs3Rz+n8llqaUdTHIwozBVdGODZYGGkZkGG0ENlgaKAKBwT4RWGdpAoo03DKXZsBN\nc+mlxK5lLu2Bm+fSMw0X5NKrAQupdBNLqwHL98ul0f5Lz9FBQn9w/XA9t6ACg+vn4oV2KjC4fuiQ\nWxDA34KlGZm2Ymmhim/B0ixgacXSImD5lKV0k6WBxHHw6lkaSByB61naE8zsmghczdKEF0mcThA4\nTAUil2Y9JUuzagtY5tJMJMnSrNqQWObSvqfC0r56fmZVFTKXTrD0nkvPKhV/Rh1rAYtoXHLBCFg0\nmDz+Vq1Cg3XhRYveYqypAM7j7yLgYUxuMpPhAuD3vjvT69/2wO79+jYuMjZXxfnc3cedEXnBl6WK\n4Tgcuu7ouo7msvF47grIfBVuAp8f9GrtQkw3vlovPAzg+6kbX14+SM5j8i2buK0BPL7upcWMMRO5\nLtYWBcT6hRFYnSsupFqibXd/lL8Ct+YKwpukHm7lr8AtVYzu4Ebavr+tb8CxLNoCRuWZi+uiMd2V\nxcfYIJ6LgMB07gKJZUDQSBUyILCcu0BiGRBYS2gFwKRtJW2P4xBLRcDaK/AIF0tFwLFbulQBnM/S\nLLlJy2LZQuwtkptYpZUs6419ilSxs/SsuHYsHUfCQz//Wo4cu1iDJ+g/dkuXTOBXpH89bbdU4VL0\nD3k5U7wFnKR/yBkFsBawIP0jSxMQy8sZsBWwsIbOCZY2gS1VALBkaUviAmDJ0o2AJUtDXs6esEBi\n6rezNFNe5omhY30eyMFWga15YBOwZaO5wCxgmZfSqa8+D+QAy4Bl2pW+GVjo2C+lb5VYAPul9ObA\n/oW3NQ/k6FhI3GopXQB/wlJ6IuCGapTYWEpPBNxQjcDGUnoi4IZqBNaX0lmC2/KFNwOWS+msmuwR\nJQ4mKjmZPatcSmfVE/BqwCI52Y+OspTuq2fR1IBFcvIccBe88NbWKxROHheeftTm8WAelQXV3OAt\nWjG8CvylWZqZWyuWlqr4DizNApZWLC0Clp2lvWcZvx0DGvY94l8WUKO51bM0uCoC17M0BHsIXM3S\npBEWUCNw0JjInBkNS5bOBi5m6VyJi1k6F7iYpSGgVnVMOtvM0qsBy9ZtaWrAEsyjsmCpohJy6mYA\nM9cvvokKDK7fDhjMcjvwz2ZpFrC0YmlhFd+CpVnA0iqXxoDFYOkyU0YdGyydCMUhOkBgg6UTTg55\nHgJHll55ckbDWL+SmTI7Du3NhfOWSTrdMUos6Z9FB9QWVREElvS/BH7VXenrg9eA/hVg7fX2F/iF\ntxKwbGNpNWBJOHAwS62g2jGzWq23cf3/A97pPwyDp38I0EN9bkEMnqf/XACtnQB+paX0BEtD9aiK\nZcBi0H/CyaEaAxaD/pNOHimRRhQHL9J/IUsTlgkcrEck6aynZGlWTScocQAuZulciWWS7nsqS+m+\nehZNlVgm6bOTa0vpMAeowHTjzUn692NpReJgIFUFS8cqYE5gUAEMk4Jy+xFYCViUHnSZ+a7STA1Y\n1to/97tlAheoYt6V3h7Y70rPU0WBxH6/23ZgsS1t2pUOk8LaYNA1S+Lh9balQfrMlGKpIrktDdLn\nbODktjRIn7OBU9vSMMsVwJlzheBtE7hgrpC8XasK9mjkb+r306Dh89S0Ct5B8ra157QAmG6j8DYX\nYDorAtay683AawDatSKJNZC161XAn87SxQFLDksXqaKUpZfJzdo4zNfKWBqTGwN4Z+kX++2Y7oOm\nHX/ZzeOaS6dYmmxed25TFdxbBEu3ApYs3UhiydIp4MyARbK0FSkXBCyk0p2ludnlnal2fOj72/K7\ng3lwsZUKTD8v9eunsXVBSQWmL/rGzdMFgL6pAI4r3sRxvlXFXwE8f9KcoE7j4l71IYDD85/0STzn\nbgLYj9jpeD+V/5g53nIEZgGLX/Eev1a/RcciYHmlFW8rcybdQjUOnrHibWXOBAzVCGyseFvpIgFD\nNQLrK94swW35XpoByxVvVk0PgBLTpemQnMyeNWvz+GrAIjnZj46y4u2rZ8FUiSUn/4TN45kBix/W\nvL9lAUse5rOVOng7Sy+12JalI/03ZulA/61Z2tN/C5Zmk5Cn/60sLezY0/9yFGvKArgl/bOA5RP2\nuwG/S/2wYA9VUU//8G4LgevpH7JJBK6mf9IMy9gROGhOpOSM3yX9ZwMX03+uxMX0nwtcTP+Qsas6\nJp0pKfmP++1YC5ZeDVi25tJqwNI8lw6uRx9T2JLnCXPbWTqo9mezNAtYWrE0mpu5K53RcBiVuQAk\njsAGSwMNI/BKnrcMWAyWBopAYJa40hfW4Av6OksTUKThlkk6A5YszUic2qKOwwNKll5KvG8efyrq\nx7F0MI/KgmpuLeifzW5ewK30r0rcgv5ZwBJZehv9y4Al/MK7cZLuc2nXOkn3LN08Sfe5dPMkvWUu\nzez467H0yubxpcT1LA0kjp5Xz9KMuSTnVbM0zTGRxOkEJfZzUOLbaJKl14CXAUsALmZpkBgDlgAs\nc2nfU1lK99UzgqoKmUvPyfK+eTwoP7Og6jizv9qsAjiPvyuA8/i7AjiPv4uA/ba0HP4uAI7b0nL4\nuwCY5Hx+0jyPv0fgX7//Ua1mUeE3j+fx97+/9U9BK5vH8/jbUsW+eTyOV2pbGqTPsSOVLB1/3c3j\ny4CFPQ4/sdfAeVtbFdBW8jasJLD2lo5ZQyLhffN41Mi+eTzq4lnaWZrU8PosnTm7lbG0NQntLB0d\nbet/PKIOXoqlIS+PEo02T/SiASdZGvLybGDcPM46TicsfWb1lsSsIZ2Yv/kWjQ1VQFtB/1TfRGJJ\n/42AJf03AlbTdlDZdFoyeOor8DXkCZj2wT4e6f+IKJP+x1yCjnHT7fOYfp66dvfCa2SW4+H+AzSQ\n3Ah21oRUAAAAAElFTkSuQmCC\n",
       "prompt_number": 24,
       "text": [
        "\u23a1k^0_q=0,p=0\u23a4\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^0_q=0,p=1\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^0_q=1,p=0\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^0_q=1,p=1\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2 m^*0_q=0  \u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2 m^*0_q=1  \u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^1_q=0,p=0\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^1_q=0,p=1\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^1_q=1,p=0\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^1_q=1,p=1\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2 m^*1_q=0  \u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2 m^*1_q=1  \u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^2_q=0,p=0\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^2_q=0,p=1\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^2_q=1,p=0\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2k^2_q=1,p=1\u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a2 m^*2_q=0  \u23a5\n",
        "\u23a2           \u23a5\n",
        "\u23a3 m^*2_q=1  \u23a6"
       ]
      }
     ],
     "prompt_number": 24
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Now the original equations can be regenerated with:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "Am * xm - bm"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}k^{0}_{q=0,p=0} s^{0}_{m p=0} + k^{0}_{q=0,p=1} s^{0}_{m p=1} - m^{*0}_{q=0} + m^{0}_{m q=0}\\\\k^{0}_{q=1,p=0} s^{0}_{m p=0} + k^{0}_{q=1,p=1} s^{0}_{m p=1} - m^{*0}_{q=1} + m^{0}_{m q=1}\\\\k^{1}_{q=0,p=0} s^{1}_{m p=0} + k^{1}_{q=0,p=1} s^{1}_{m p=1} - m^{*1}_{q=0} + m^{1}_{m q=0}\\\\k^{1}_{q=1,p=0} s^{1}_{m p=0} + k^{1}_{q=1,p=1} s^{1}_{m p=1} - m^{*1}_{q=1} + m^{1}_{m q=1}\\\\k^{2}_{q=0,p=0} s^{2}_{m p=0} + k^{2}_{q=0,p=1} s^{2}_{m p=1} - m^{*2}_{q=0} + m^{2}_{m q=0}\\\\k^{2}_{q=1,p=0} s^{2}_{m p=0} + k^{2}_{q=1,p=1} s^{2}_{m p=1} - m^{*2}_{q=1} + m^{2}_{m q=1}\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAXsAAACDCAMAAABbVkv5AAAANlBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAABHL6OuAAAAEXRSTlMAMquZdlQQ\nQN0iRO+JzWa7fNDgyowAAA/USURBVHgB7V3tlqo4EERF15nx6vr+L7sdkk76g07SRMc958KPgUBX\nV6cEBFKD07Q6XW7322V1S9fKIfgQGMobxZMe+lL5ogkNX3xcpsuDr/K0huBDYKhyFE866kvliyY0\nfPEJzZ/tO/4QfAgMdY/iiRK+VL5oQrMszj+3Zf71A7PrQW7ubW+DJ/Jt4FLaKD5muocdz5fKFx1o\nDs9lOkbK6/cy96eJ8PTXBY+MARnJXWBGGxuj+JDlfJu/QRErVSmZ8lvRNGZZvkbJg/anGabzsvby\nTOeZscPHddjnrxUkH+R2kStd0orz9fgFi0YpuWQON6J5ELTuQfAjRB+e5NTyfU2Bg18bHnjuCJJ7\nwKpfsGIUDym+v+evcJ1npMolc3ojmgdh6y61v52mwxxOP5fT/RY+ePd0Ww4hDzx3BMk9YF6fnxzx\n99P9MM+naZ4fcVeMpxWjFCxZoIxo5OBzpf2f+/kbD30e2tk6fU8HrKwTksM/QY413r+fd7iyO16m\nZdfD1cYcS/aheDKl/fNETkA8tq/1eJwvzktT7Mj0CfLUq6/LHM62T9jbT/FKr9rdVLITxVNK7Q8/\n8yN+6fI4aDXu2XDz4fik3x8qTVqRwo9/YPoJf2Cn204+TbC/wtRJnsLPD5wC+DiHCmDhui4AdpCV\n3ESFqnDCDKkttZ8fgX911298j6TNZ/iOCL1oTTQb7vfbyb9vP+G7qZc8hfMaf6DXUMH09ezoP5Y8\ntVCUg/YZ1kvtH/N0uE6rB13j+iltDrdjPactmg07sp0cdtagfS95Cqe6BMnjXnN/dPQfS26iKAft\nM6yX2v9Zdtuw38MlKFzn3I5xOvO7jPtpvt/g0NSb7/N9DlcI99P5fj99n29LLjN8KQ07spkcsiza\nd5Jj+EKOf8Kd1BQq+DouzxCr/Z+wZIFivY4SaomQUWqP66fzCb70cwsW2D3beb5NF7wTWKLYZlhz\nnv9AgjPcHy6Najh2ZEkU/njJAbJojwka5DIcYXReLyFrTyGwzImFhFIiU/s/B3mpyI4YuAY+gLhk\nYpth/fGwfDrh6IA7lGo41MgnN7nQvkUuwzn70mqUoEpOKWivVQohkan99TLd4CuzHDH5Di+cRMIh\ncYMrk5XNqYgJjoqw71ynIH47HGHL3E0uxWyQy3BGHht2CSvBZRUljilsiUztv2938RQz3bN9hQc+\nX9d5BunpJG7pLnD6DCf+4xxOoc1wmglu6L3kh/l5IleGLXIRzsljyyxhLTivY8QqhZDI1B7SsdN9\nTr/sx/fVK6ESI5ac4QH9UfJYvlGC6Fu1WU1R0f70JHsSYYC1X1f19UgC1KIzPOA/Sh47YJWgumev\nqKeoaG+n3Le8RIFd+5fIuCmJrb14+ODNPgQfAkOlo3jSWV8qV7StvXj4QMrpWhyCD4GhvFE86aEv\nlSva1l7cCJByuhaH4ENgKG8UT3roS+WKVtrvPgUU/u0+BaV9sgrwxzdYTv9cPruoIuPoXAj5W3wK\nixxK+zxe6Dp8tLQeeL5ZQHIPWFO/5pzzVp/CuvZoFRj9vvJ862TtkdwDXtN+FA85d5/CmrDVdbtP\nIe/IVZ3yxhy++xSyJpsWdp9Cr2zyu7ZiFYDnwTXvB97SdVoFUjgb9N9OjrV1kqdwn08BOVjJLp8C\nZkgfjtS+YhVYHdwvn3H6fuu1CtCvQzznbCf/DZ8C6z+W7PIpsAx6rLxiFWgM9qTrwl6rAL2MxI5s\nJ0+19ZKvdaXDcUDGhLHkDlTZPQWt3O9tqwCkKNy/7VNokOPmd/oUkAPmMKH2Hp8CwIqC0JDah7zL\npMbpYW1B/rZPoUHONkOD2wWUp0KGQ1tN9f5n7QWOEwufAsQWBaFhaq8G2TmybjwIFdER+3GfAitb\nk7PNbXIZHuoVU6P/m3wKgtbUXo/TI/IjPoUGOW5GAaldYMUkIcMRRuZ2/0mQXqTEyqcgaE3t1SA7\n2PIXL8BnfAoNcmE8YHaBFZOECNcSrlklOkATI1YSigym9lCPMcgeRtCdxgNneNDio+ShALOEuLHv\nr9GLCA7a//P8Zy2TNci++xTW1FpdZ0kYg/8F7dn/W63m2Fe+Q4HaOecdfHvOokBN++rjm5LCWhqC\nD4GholE86ZQvlSPa1l48fCDFdC0OwYfAUN4onvTQl8oVbWsvLkZJOZ2L7B6uE5PDhsCQZRSfC/Gm\nchAr7dGnMF6/owrsaSHfAMYky3wUD0kWnwLMfakc0Up79Cm4SZcu0z/9VUifwji3UzBadlpO71Ow\nSiklC2x/r/XzHLQKWKSCqtLsrwIfCk6FvB+8XsEoHp7HxfcpGPt9LlnSO4jVfo9WAcjpSCMrCO1+\neO5IIe8Hr1F7yNfx2adgpMolS7ijcKU9vtIAH9/I3O12tAqIZxdVWO4IknvAPLOfHPHizQhwsbRs\nMUrBkgXKiEYOPlfaf8IqgB2ZPkGOevjejIAl+1DIFedB+8tMxsA/8UoD7Mhf9j6Fr/ASAfL+g4pV\nAI0I/MPLLdzcaRVw+hQwe6aTC/F+spM83ff6fApYwnafAmZItctzTsUqQJ0FsufQfq9PoUGe7id7\nTRKrt5+tNyOwEvKh2kJRpVgGPWZYsQpQZwHNmJbT5l6rAM2GHdlOnq5GesnXLl6ajgNacR6vbaKo\nUiyD1t72KTBX92/7FBrk0MHl2u6dPgVWQtbe41PgGbT2+WNS4/QM+ds+hQY5ao/Vc7vAa3wKrISs\nPTKmOScWPgWeoaK9HqdnR4y2CrDNUMprfQosuyaXJ5EGuQwXCoZmo/+bfAqsExXt9Tg9flN8xKfQ\nIJdiUrvAq3wKWMLKB1VWUWLlUxAZ5HVOzqIG2fHFgZ/xKTRe5iDuJ5ld4FU+BfE2hKwUXWDESkKR\nwdQeMhqD7LtPgapdXzYkjKCK9tYg++5TqOtNtloSNrUnOfbFdyhQ2e/fQbfnJAoE7Q1v1OhYv2PE\nnhSUFofAkGMUTyrypXJE296o1WcepKTG4hB8CAyFjeJJ33ypXNG1c45jCIYUmxeH4ENgKGEUn3vh\nTeUgVtoXq8Bo/Y4qsKeFfAMYkyzzUTwk2X0KTNH+xrD2u0+hX2wROaz97lMQivY3h7XH9ykY53sc\nclAVOYjV+R6tArtPYfcpqB2rumLDj05gPp/jAPd7Hwq54lzt97tP4f/4uw/8QzNb3VaBJQMb9K+Y\nJEw6vsFH7vMpIBMr2fk+BcyxzOV+X/EpNG7UkwGi1ypA/RJ4AG8nR9PH8vZ+1kGjsXLr33YcEBCW\n7HqfglAwaE+f51SsAo375TQo02sVoGM42JHt5Km2XvK1rjQdBwyEJTdR9LNnGaZJPs8pPgUKSsvV\n66c0GIlWgRU4W0XHLrEj28nTlWAv+dqFo3AcsGJ1/7HkDhTNxBSU55wSKH9zArYUpO0RwQTs1yca\n4dgRxML9vPjBjwY52wyNBrkMh7aeqiVYPgVGLH9zAziKgtAwtVceEYZ8rUdEDfp7yVlt0OBWjdd4\nRDiHKhk2w8SJhUcEtvdprz0SDKltGvQkEspo2DRkeIDkyU3OaoM0DXIZnonLQqOEEsiXKPFKij7t\ntUcEC/6IR6RBjptRCGrVeJVHRHIgF5tTYuURERnMc44yOOBDhs94RJIJxCL/DY+I4GCSY6PuEREZ\nTO0hm2Fw2D0iqHR7bkgYgRXtLYPD7hFpa54iLAmb2ncz7IHbFLD3+8Pt9oBb9K3TEHwIDBWP4kmn\nfalc0bb2cMNT/uWSFNO5OAQfAkOBo3jSR18qV3TQnv2vW6aFX8wr/2qc13YvDMGHwFDiKJ700pfK\nFb38rxu74ie8Z7haHZiG4ENgKHoUX3wKzlT9xCGS/p/hVGwa0yU829o+bYAX8g1gVukoHn5/dF6e\nk8HR75LBEa3P9/FnL6AjF/K4mvWrr+GAlxdDILkDvFrNKB6S4vsU1lOVkjn/ejSPwZbSPp/jL6ev\ni+9XCzHnMvfA4RsqTkjuASOWzkfxkAt9CkaqXDKlhf3VI5rSPl/bXJ/P58D53gPPHUFyD5h3PrZG\n8UuWuGsbqXLJnN6I5kHYUtpnjwhGuOfxdRIeWO7IJ8ixUPFWClxtzLFkH4onU9p/4nUS2JG/8F0W\n7Dpn94j8Hz0i1FnAD56lhZs7bRopnBkuKh4RzL5CHFfFi7JOcrQ1PHAKL1I5htd6hAeP1/CkVk9Y\nAiu5iaJ5MENaJ885FZsGdRbQjGk5be71iNBseM7ZTp4cAL3kwjAQe9DyiNCKy3htC0WVYhn0eG3F\nplEd5cu/pdZr06DZUPvt5GlEqJdcDCAt+jTdHrTirH0TRbVnGbT2xaYhx+nZ2wBs4wHaNNiIvR2+\nlIbabyaHLMtzkU5yDKe6TNLtUe1/1l6gWK+lT4EpCNzynJPLUVYBhnytTwG130wOQPZMitsFXuNT\nYP3P2ueS4wInFj4FnqGivR5kZ0fMS30KynDhJhfav8WnwPqvSk4fQtWnwDJUtNc+Bfym+IhPoUEu\ntad2gVf5FLAEsavzJiVWPgWRwTznaJ9C43US4kUNbMR+5XUSIpz3wE0uHAAtchHOyWPLLGEtOK9j\nxCqF6LOpPaQzBtl3n0KWurlgSBhxFe2tQfbdp9CUHAMsCZvaY4J9/iYFKvv9mxj3tKiArb3L7oDp\nynwIPgSGGkbxpRvOVC5iW3u44cGhDFJK9+IQfAgMJY7iSS99qVzRtvYuuwMpNi0OwYfAUMAonnTH\nl8oVrbQvVoFhn0W/WwL7Ssj9YEwS56N4uO3EN3T7UvVHK+3Jb2447A6837HVDy+D/uhTcBozNH0/\nucYuaxoekVKywDuIlfZoFfhNjwicJeOUyT1WCwTT+SgectU9IrlkygrLHmKlff5+ddkdRAWhCIdb\nIncEyT1gTe0jX8PDuk94RLJVwGV3oB2IPgUPPGuP5B4wpZ4mPznileOgyyMiUK7C1X6/+xSm+Yif\nhz3H3WV/l0XW6PE4X/DyJK/sWvi6zMEJ9oS9/dRhx0vaO1G8ErnfV6wCaXCf40sLB+E7rQIpnA36\nbyfH2jrJN/kUkIOV7PIpYIYkmtS+YhVYHdwv2qeBgV6rAB1HwAN4O/lv+BRY/7Fk17ssWAY9blWx\nCsixoSL7spQGxHqtAnT8DDuynTzV1ku+1pUOxwEZE8aSO1BUJ5JBa29bBSBFQdrGg06rABs3xo5s\nJsfaOskxnKrS9ClwEJbs8SnwDFr7XI7yKTDkb/sUGuRsMzS4XeA1PgXOgdrDWjpxYuFTgMCy90JD\nnu9zIm0VYMjf9ik0yNlm6AO1C0y6Vhmee10WGv3f5FMQtFF7sNo/n/yqVvsUEPkRn0KDHDejdtQu\n8CqfguRALjanxMqnUDKEezCY4A58Xib+eEgNsu/vU5ifp3WXbFG/7lMo9ggwvYWpAOWSMci++xSk\nUHbbkNAGpC3WIPvuU2hKhwGWhLh9+g+7kPST5OF8bQAAAABJRU5ErkJggg==\n",
       "prompt_number": 25,
       "text": [
        "\u23a1k^0_q=0,p=0\u22c5s^0_m_p=0 + k^0_q=0,p=1\u22c5s^0_m_p=1 - m^*0_q=0 + m^0_m_q=0\u23a4\n",
        "\u23a2                                                                    \u23a5\n",
        "\u23a2k^0_q=1,p=0\u22c5s^0_m_p=0 + k^0_q=1,p=1\u22c5s^0_m_p=1 - m^*0_q=1 + m^0_m_q=1\u23a5\n",
        "\u23a2                                                                    \u23a5\n",
        "\u23a2k^1_q=0,p=0\u22c5s^1_m_p=0 + k^1_q=0,p=1\u22c5s^1_m_p=1 - m^*1_q=0 + m^1_m_q=0\u23a5\n",
        "\u23a2                                                                    \u23a5\n",
        "\u23a2k^1_q=1,p=0\u22c5s^1_m_p=0 + k^1_q=1,p=1\u22c5s^1_m_p=1 - m^*1_q=1 + m^1_m_q=1\u23a5\n",
        "\u23a2                                                                    \u23a5\n",
        "\u23a2k^2_q=0,p=0\u22c5s^2_m_p=0 + k^2_q=0,p=1\u22c5s^2_m_p=1 - m^*2_q=0 + m^2_m_q=0\u23a5\n",
        "\u23a2                                                                    \u23a5\n",
        "\u23a3k^2_q=1,p=0\u22c5s^2_m_p=0 + k^2_q=1,p=1\u22c5s^2_m_p=1 - m^*2_q=1 + m^2_m_q=1\u23a6"
       ]
      }
     ],
     "prompt_number": 25
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "We have $n(p q + q)$ unknowns and $nq$ equations so we need at least $1 + p$ cycles to solve for the unknowns."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "n * (p*q + q)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$18$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAABMAAAAPBAMAAAD0aukfAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAzRAiu5mrdu/dZolE\nVDLjuNgcAAAAeklEQVQIHWNgYBBiYGAq33qAgYHB5BMDAysDQy4Dg0oYkHmLgcEQKMoOZG5hYNCE\nMs87MeyCMtn+r14AZTLEf1SAMrni5HugzGUMXD+AwiATmhgY3l8AM5k+MDCwLYCINjMwsBhAmAcV\nGJYDrc/4mcXAlwp2DtAQCAAASJ4i3e599ZEAAAAASUVORK5CYII=\n",
       "prompt_number": 26,
       "text": [
        "18"
       ]
      }
     ],
     "prompt_number": 26
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "1 + p"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$3$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAAAoAAAAOBAMAAADkjZCYAAAALVBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADAOrOgAAAADnRSTlMAIom7VJlmdt1E7xDNMpCR\nWcAAAABSSURBVAgdY2AQUjJhYGAMYPBLYGB/wsDXwMC5kmHfAQYgAIoAwVEg5tUIADEZtC6ASK7V\nDIwCDMxAxa8ZmJ8xcBkwcDxkYEtg8CtgYJgaaskAAFKHDvy4QzOnAAAAAElFTkSuQmCC\n",
       "prompt_number": 27,
       "text": [
        "3"
       ]
      }
     ],
     "prompt_number": 27
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "Am.shape, bm.shape"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\begin{pmatrix}\\begin{pmatrix}6, & 18\\end{pmatrix}, & \\begin{pmatrix}6, & 1\\end{pmatrix}\\end{pmatrix}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAALoAAAAaBAMAAADyPFssAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMkS7zRCZdiKJ71Rm\nq90icBAQAAACxElEQVRIDZ2VPWhTURTH/y9NXl6blxh0c2kGV2lAHNziIFJQDNhJxMRvHMROik6B\nulsQKgWHIAiCy5tcdHibVChEBAd10EEnwVpaKFbU/znvI++93NvBA7nn3HPO79zvFwD7ATgtNnZZ\nNIeENEpVvAo5bVoHjEmp83BqZQ0ls47ULo1oKlTu0jqRBozGdGhyK2kKAI/pVug6DVdWsXq3WUzl\noJU7Ky3gazEifSHNEBoyZULeNvUDdhZabodmVpYYc4GPwO2sO7aFtECo9plEaFr0Cn/zmGpTZeTs\nTfKvgPscPuNOTCEtELwNZhGaCaifA/5OQo11ldUfAheARjj2JhZJGwR8ZxahWW6KwypTw4Qaa6k+\nd1SOyC0sS3JI2iDgMhMIvaHyh8C+448u0cyJVK//XQuB6jAX0A5JGwTcYwqhK1RynL1T0SYpGDdS\nHbe2uLqSWAUhaYUwF0FHyMhx9rbhPi3wUt2/eOYJ9+5nIcQuSRvEbWsqdJp55REX+Rn+ZqGEVF+H\n/6sLrxhiJkkbxPvS4mXfxDfmNfq8Ox34v2lnRapzPbMBvN2sP7JJ2iBWD1h9V2agc+fxm+Ze2eDB\nhra5WyCWjOeeVK9x34fRpNJW5v6MxzKAY3gNJG0Q1zvgvu/gHWvV22yO4VxQ2UoriyHVD3Xx2nxn\nhLRAfCZdheKLyU/CtbfADzpTcT/8+QTnvXzF5EWgMLaQFgjJM4oHYaZIaaBqooneam5snZ4mmqD4\nrc60mMEbHYmfGAVdXhRHvoyS4jZBJ+kn5C5T35AkkfORmmh7umP5MkpaoApPDIScIfVBSRIJtZ1s\n1tWVH1tJ8YcazDX+iF2B+N+AepCLTXa+qCvMB4S0SKPJgEAL/FU6bPYQf9kUFNIia/QrVOvTemHJ\nit1LxrCSxojXpjuCXtKq7b01V401IKRR3C7dFsgI/JfzH5YmuXfLGAVgAAAAAElFTkSuQmCC\n",
       "prompt_number": 28,
       "text": [
        "((6, 18), (6, 1))"
       ]
      }
     ],
     "prompt_number": 28
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "b = bm\n",
      "A = Am\n",
      "for step in range(p):\n",
      "    b = b.col_join(bm)\n",
      "    A = A.col_join(Am)"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 29
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "A.shape, b.shape"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\begin{pmatrix}\\begin{pmatrix}18, & 18\\end{pmatrix}, & \\begin{pmatrix}18, & 1\\end{pmatrix}\\end{pmatrix}$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAANAAAAAaBAMAAAA0+c+VAAAAMFBMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAv3aB7AAAAD3RSTlMAMkS7zRCZdiKJ71Rm\nq90icBAQAAADC0lEQVRIDbWVP2gUQRTG316ye5fL5QhqJchdYWEjORDFzmuClbigYCPm/IekCF6l\naGMg9icIqSyCIAgWprLRYislQjAgWJhAtNBKMWoEScT4fW//zM7cXjof5Jtvfm/mzcztzkZE9ghC\nhWZwlItTXrOYJ7TLVsVrwakkqUFNaaUws7eQZvAwncpwCKeSJQeYB4V8spBmcCSCVblGRuFi525d\nYteKjNe5JTeCLggW82/ON90cuXwihVR+oaVUTqBOR8YWmciF4eVODqe2h1lzmB6IrKYsbcnlBnuQ\nkQ5aiH9sPZTRSCpt9HOR45WNHE/tvMiZGRR8IXI3ZUmrXHrsQao8gcqdUEZa4n9gIh+Gf8njxD9B\nW8ZC90UuuGlyqUexNHB0UUHB0s9QnxuQCcMvG5g6j7VYcOKo9L0sulDQwgjIa05RQUFZ/9M3Wgy/\nzbF21BbQZ8GxnaXITsVcyjpiQa4wq8KCwQ5fFDsMn+Dp7Qja6OvOr2/2ZZWXeGbIEU5UYcGX+7aa\nBPkwvDGe5+qHWmhYsHbx9EM3qwt534Ehp5hVQcHairx55A43vNd0czK8AsSCy1Lbwk6t0IUqP8Ag\nn5lSQcFqKN62NRYdw3uLbk7qHSAWxP4abjpeiAUr2/aJpsDW8GcFFkr4LifyN/A+RNa85NmlJ7J+\nuuqsyHFntJ4o5g00TmQ/3WORITetJ/J+YwrkLWeqYOelSal1/U0iEynHXXEfAk7RwkAWPBjKK5lq\nm2kJ5wunb53eDcrM06+RHJjGR/WbVc/w+LrZxXhLgvd/18R7h49qnUWzUC560yC6S3urpdlsrG30\ny2AX0/2aUYeMTV36Zag2QVTSjNQy55iT2reL2ZcncmagO9yNJbiHViUbcz5ztvHjHyay6HS+59s5\nTU3xQUC8BbQqiilR5myD64xwiu3Pjyn6LZY5gLJKo0KzW9THmXWK9f2jdCt8JKCcpVGh2S2WipJ+\nu4gaVuOjURntwKmYbKGrtArxs0KawTk6FXlOq0IzOIKwMDe6WIhTeJVGJSX/uf0HSmHWidRWRpoA\nAAAASUVORK5CYII=\n",
       "prompt_number": 30,
       "text": [
        "((18, 18), (18, 1))"
       ]
      }
     ],
     "prompt_number": 30
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "The following should compute the least squares solution but since I didn't actually create new equations for the steps, then you get NaN as a solution."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "x = A.LDLsolve(b)\n",
      "x"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "latex": [
        "$$\\left[\\begin{smallmatrix}{}\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\\\\\bot\\end{smallmatrix}\\right]$$"
       ],
       "metadata": {},
       "output_type": "pyout",
       "png": "iVBORw0KGgoAAAANSUhEUgAAACQAAADrCAMAAAA2YMTLAAAAM1BMVEX///8AAAAAAAAAAAAAAAAA\nAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAAADxgEwMAAAAEHRSTlMAMquZdlQQQN0i\nRLuJ72Z8fOynKQAAAZtJREFUaAXtmuF2gjAMhYtWxOmU939aIT0XanoV1sFWMf1B0MS25ssNHMW5\nrLGb86nz/wdVrYxxu9GemuBzrmr3vhuHYbtR0LX37No+qBr8chIFyetrTtDl3I3v/nAZZ6czje5w\ntnLQmLEXyy2SApl/szMZ4FA/BjjkITqurOBoJTmly2XL/NgeHxfQgG857VCm1DPRjRcW9FYyr339\nmt3Jr3W9M8CPmccrWuNwwtIgpmAT589v10pRcCJODdjE6agQpq4IpQBOxKkBf9aVU5rbFLtZQX8A\nOBEnWjMsFSecsLR8dRXQoEXyJNvY7Ey/q4JEnGAGS8UJJyxlZ4DdIkX3FPCK4gRYWAOcdw2W/E1V\nwVPAK4oTYGEN8MKAY3YUcH5jBTNYyg5OWBpkLZq26KSxagXrux4KmGYcOGBLCPq8KojFKSQ0YCpO\nHaTZFVsFmwI8T5zNCQoLVrM7NBm/+BlgnkwtBCqpdwgywAY4ZCA50vLVUflBRbToGY8O1P3jAd5/\nDd88aqzy6IDv/kpJRhSU+IY3FgmiCu6WuAMhFTqkOvNfLgAAAABJRU5ErkJggg==\n",
       "prompt_number": 31,
       "text": [
        "\u23a1nan\u23a4\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a2nan\u23a5\n",
        "\u23a2   \u23a5\n",
        "\u23a3nan\u23a6"
       ]
      }
     ],
     "prompt_number": 31
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "type(x[0])"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "metadata": {},
       "output_type": "pyout",
       "prompt_number": 32,
       "text": [
        "sympy.core.numbers.NaN"
       ]
      }
     ],
     "prompt_number": 32
    }
   ],
   "metadata": {}
  }
 ]
}