Square root of 2x2 matrix.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Square root of 2x2 matrix using basic algebra\n",
    "=======================\n",
    "\n",
    "We want to find some matrix $B$ where $ B = \\sqrt{A} $ or, equivalently, $ A = B^2 $. We take\n",
    "\n",
    "$$ A = \\begin{pmatrix}\n",
    "  a & b \\\\\n",
    "  c & d\n",
    " \\end{pmatrix} $$\n",
    "\n",
    "$$ B = \\begin{pmatrix}\n",
    "  a' & b' \\\\\n",
    "  c' & d'\n",
    " \\end{pmatrix} $$\n",
    "\n",
    "Expand $B$ and equate terms."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 71,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/latex": [
       "$$\\left[\\begin{matrix}a & b\\\\c & d\\end{matrix}\\right] = \\left[\\begin{matrix}a'^{2} + b' c' & a' b' + b' d'\\\\a' c' + c' d' & b' c' + d'^{2}\\end{matrix}\\right]$$"
      ],
      "text/plain": [
       "         ⎡   2                        ⎤\n",
       "⎡a  b⎤ = ⎢ a'  + b'⋅c'   a'⋅b' + b'⋅d'⎥\n",
       "⎢    ⎥   ⎢                            ⎥\n",
       "⎣c  d⎦   ⎢                          2 ⎥\n",
       "         ⎣a'⋅c' + c'⋅d'   b'⋅c' + d'  ⎦"
      ]
     },
     "execution_count": 71,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from sympy import *\n",
    "init_printing()\n",
    "A,a,b,c,d,e,f,g,h = symbols('A a b c d a\\' b\\' c\\' d\\'')\n",
    "\n",
    "A = Matrix([[a,b], [c,d]])\n",
    "B = Matrix([[e,f], [g,h]])\n",
    "\n",
    "e1 = Eq(A,B**2)\n",
    "e1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 72,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkcAAAAZCAYAAAA7WIpZAAAABHNCSVQICAgIfAhkiAAACcZJREFU\neJztnW2sHUUZx3/F2pfQ2qKpiE3TSszFCkotKpJo09QXqoSoMZqIShaDMQaRD34oaEivGl68voAa\nDdZAFnxBtLxZMAFDJNEgeCkFUVGIpBhrhRaVGnIVaK8fntmcc8+Z2Z3ZnZ2zJ/f5JScn2Z2Z53nm\n/Gdmd3Z2DiiKoiiKoijzlouAaeAQcADYBZw0Uo8URVEURWmCju0NuQM4B6m01wE3A/8AXjpKpxRF\nURRFqU3ysf0s4F7gGeAwsKEtQy1xD/CFkvPLkLjOTOPOEJ8CZoGPB+Spiqnt/DGpE39qe676uhZ4\nCjg6gl9d+k0K1iB1ddPA8UFfTzHpzk3k17ijeppLv6/jrCVXfHUZxdgQq4zY+PZF/ZSN7ZPA88jF\n0y3AhK2Ao0ocei3wA+BE4EbgUlPYuDABnAb8rSTNcqQO/pXEo2FOMd+7PdP7xNRm/tiExp/anqu+\n3gR8DLgceLahT137TQqKunqg75jN191IB/MlpENS3KieyvU0zlqyxRejvFRjQ6wy2sC3L+qnbGy/\nG/ga8CDwXuQ6J4jzkKu1z4Rm7AjnA08DS0vS/ATYA7woiUfD7AFmgIWe6X1iajN/bELjT23PVV93\nAv+2HK9D136Tgi8i7f/dfcdcvr7ZpP1cGtfGFtVTtZ7GVUu2+JqQemyIVUYbhPRFBb5j+27gCIEX\n4xcbh94ZkqlD3Ibcibn4OvB34Pia5WdI/WyumX8x8BxwX0Ceqpjazh+TOvGntmerrwmkMe2I5Ffd\n3ySjmf6quN2Uf2zfsTJfHwGeoHw2er6jevLT0zhqyRZfXUYxNjQpI6M72oGwsf2Hpuw1IQ5N4hfw\nWcbAo8B/kGmsaWRxVGx8bS1C7sRcAV8B7Ade08CXjGaCeKPJ/x1kLdeNyJXws8id5OBK+6qYNgM3\nINOM/0PiuwN4n2f+1ITGn9qeq74uN+W8vcRWLJ2WkdG8Q1oIXAD8DvgvMiBtAxYg+tnXl7bK1+3G\nn9Mb+NMGmylvF6nsqJ789dRVLYXE14TYYwO0Oz5kdEc7oWN7jvi+LsTZSaoDXg68APwGuBq4zHwf\nNHm3eTqWefgTYmsL8FNHOd9A1k6t97BZRkYzQXzS5L8dEf1OYAr4lTn+FLCyL70rpgXAlSbPAaRO\nLwN+hNTNVEX+UREav40cf/3Equ/7ER26Fs7G0mkVGc30twjpaGeR6eevANcgU/nfM8d/1pe+ytd3\nmDxfrelPbHzbhY0cf101bX+qp2G6piUIj89Gjp+uYvVVkGZ8yOiGduqM7Tk1Lo6KCj21JM0y4BWW\n48chdzh/8nQs8/AnxNaFwFstab+N7IOwxZRVfOos/stoJojvmvwHgTcMnLuO4efurpguNWl3MtzB\nLgPWVuQfFaHx28jx10+M+j4aGageLrETQ6c+ZDTTX9HpXIx0oAWbzPFZ5AapoMrXFSbPb2v6Exvf\ndmEjx19XTdqf6slO17QE4fHZyPHTVayxAdKMDxmj107dsf0qU/7Jvs4uoLdQydYwffgzcoVbxnHI\n9NeKmjZCbEGvogc/kzVsZjQTxP24X9U81Zy7vqKMjcjritPI1fc4ESP+EP3EsDdh0t3pYc+Gr059\nyKivv2LR6y7H+T+a86FbXMzQjTdam7YLX101taN6ctMVLUG8+Hx1FaOvgnTjQ8botVN3bL/QpLug\nytGtyCtu0ybDVRXpj0Gu9u5FnoMfHnAs1muNqW0Nshd35ds+eUV5L0ae++7F/jbCWlPOzRXlfN+k\nK1uv0EVixZ/a3mkm3Q0ladrQ6V7i6q/QzVsc54vp+9WBfu5DZkJGTap20dSO6slNV7QE7bUXGzH7\nxjbawV66r50QViMxPY+84XYJ5hHbYOVvpXcF9XvgWyWFvh654zkWmf78MbJo7AXgVcDZwEMRnE9t\ny8aVDK9/2YDskXAtUrn9PFhR3knIlfwu7B1AMdX5REU5pyMd5i8r0nWNWPGntjdjvpc4zrel09j6\ne5fxy/U2zPHIjEToAtOl9OpolKRqF03tqJ7cdEVL0F57sRGzb2yjHYyDdkLYh6zBugL4oDn2C4bj\nAOAlSKAzyMp21z4BDyGryjdbzhX7EpxX0+FR2vIlo/5U4rkm72cd579szr+npIwlJs2eGvZHTYz4\nR2HvlSbdrx3nU+o0o57+Ct24Zhw2mvM/Dyz3KOQx/F8C88UmVbuIYUf1ZKcrWoL22ouLWH1VyvEh\nozvaCeUMem3wRGTmDrDvI3EIuBV5fXA19pXfa5C7mrvNp5+V9F4xjbHzcUpbqSh2/LT978sq4BPA\nY5SvRSgWrr08ol+piBH/KOztR974OMFyblx0eth8XLr5vPkOfVxzAqLJqjvFtknVLmLYUT3Z6YqW\noL324iJWXzUO40PqurWxxXxPAn9AHq9VcglyRbXJcm6VOfcofVdawMuQvRNmjRHXdHGBzwK1WLZi\nk1F/5ug+k/dx5r5BsBy4CxHMVo9yHjblfMhybgL/nb9zU07mmb4pseL3XeAYyx7IWx+zwKsHjqfW\naUZ9/T1i8g4uctxGb63A+wPLPMfk+/TA8Zy02oLm7cJXVzHan+ppGJeWYDR6ihWfj65i9lWxxocq\nMrqjnVCuwd7+SpmkPOC76E2JTSELqw4ijf0wfs/Cc/yEHsNWbDLqCWIh8sjyAeQ13MeQfR2+iTzG\nPIL/X7acgcR/BNkTYwrZaXcauSv1pXg99KMBeeoSM/6cav3EtAfwYdyPM1LqNKN+h3S2yfsc8ttP\nIWtaZoC/mnNrnbntXI+skRjcmC2ltgqatoscv34pRvtTPQ3j0hKMRk+x4ssp11XsvirW+FBFRne0\nE0pubKwLybSd8oBXIcEcQPbbuAcJtHhOeHWAY1lFuhi2YpNRTxAnm3w7kMeWNyGPMg8hd4Sh5b0N\nEf4/kTvJ/ciU60cCythj7B8TaLsOMePPqdZP7PpeBDyJfQFhSp1m1O+QQF68eBzRzJPIY/SNyELP\ng4FlrUA6s1ss51Jqq58m7SLHf3aiaftTPc2lTEswOj3FiC+nXFex+yqIMz5UkdEN7dQhp8YFWPH+\nf5tTWsroWYncXZTtGKzM5SKkbQxu0DZfOR+pj8GN2VRbfqieeri0BKonJT63InoL+k+8D5hMtyHT\nm+P0J4CKP2cib8TU3exzPrIEuatxbVw2n1iK/MnjTss51ZYfqiehTEugelLisRjZgPIZ8wm6vllK\nb0FX8dkQ2UFFGVc2IY+eXf+JNV9Yj6xPXDdaN8Ye1ZNqSUnDJHOva7bbEi2wHexjMfKq23rkP0p2\n0J0t3RVFURRFUULYjDy2fRr5Q+cubBmhKIqiKIqiKIqiKIqiKIqiKIqiKOPI/wHIskFSgWVz5QAA\nAABJRU5ErkJggg==\n",
      "text/latex": [
       "$$\\left \\{ a : a'^{2} + b' c', \\quad b : b' \\left(a' + d'\\right), \\quad c : c' \\left(a' + d'\\right), \\quad d : b' c' + d'^{2}\\right \\}$$"
      ],
      "text/plain": [
       "⎧     2                                                         2⎫\n",
       "⎨a: a'  + b'⋅c', b: b'⋅(a' + d'), c: c'⋅(a' + d'), d: b'⋅c' + d' ⎬\n",
       "⎩                                                                ⎭"
      ]
     },
     "execution_count": 72,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "solve(e1, (a,b,c,d))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can calculate the determinant $A$ and $B^2$ (confirming that $det(B^2) = det(B)^2$)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 73,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAANsAAAAYCAYAAABz/4haAAAABHNCSVQICAgIfAhkiAAABcJJREFU\neJztmnuIFVUcxz9r67qLK1qy2QORpJYEe629hIrNIgsR6p+gIsEIIsr0v80gMnvJLXpBRYZhWZm1\nq9EmpSAFgWW7miVkGdW1UsvdLAzZcl/98TvDzs6eM3fOzLlzr3E+cBk4Z36PM9+Z87zg8Xg8Hs8J\nxnKgCzgK9ACdwOyKZuTx/E/ZAixGPrDzgE3Ab8ApeQSfDgwDG/MIFuJuFfeOnOO65DXgMDDRga/t\nwMOWNibt0vhyQRpNs+YatZ+jcrgzoX0jMAgsBBiXIZEkzFHXXWWOY4q7M+e4rrgEuB1YBRzL6KsZ\nmAv8ammn0y6tLxfYapo1V539TuA94BHkQyrFJOQb+zNlDlasRHqCG/IIFuJLoA+ozTmuK7YCfwEN\nDnwtAf5I4UunXVpfLrDVNGuuJvtLkefyQAIf7yB5n5QyBys2I4lNyyOYYgJwHNiRY0yXNANDwGpH\n/j5ARkhbdNql9ZWVNJpmzTXOfi+wn/iZ4dPAQWBmXJBbgTeBfcDfyBDYhSz8dNQCS4GvgX9UEm1A\nDXAIOBAXrAxcjLwkLwIXAh1ID3UMGTHidodagQ3I1OFfJP8twI3lS3cMq5D8r4m5J6lGdcgIOd3g\nx0a7Ur7Kia2mpXJtJV7nUvYPqXzmG+qfUT7PNTdJ5pgDwGfAGuAJde1Vztsi99chjR1GhssngVeR\n4f4VVf5+XMAycJeKuxkRox0oAJ+q8sPAlIhNDfCsqu8B1iJtfwtpeyGHvAO6EQ1MGyM2Gs0D3jX4\nsdUuzle5sdXUlGtSnUu19Vrl4ylN3XPIDuSsUo1qBE7TlJ+O9KDfRsoDUR5UDQm4SpUPAytKBXXM\nyypuL3BRpO519PPtx1V5O2Nf8kZgRky8ZUgbk/7iRsmJyIe0J+YeG43uB64w+LHVLs5XubHV1JRr\nUp1LtXWy8vNFpPwF5IxtHqJR8EuymTKK75AeJCBYKHYa7v9G1S+0DZSRbsxbxJepuvWhshZke7YL\n6e1tKTLycib5rY3x1azu2ZoiDxirkYlq1c6EraY6suocpQ8ZwcKYNF+hc3Ay0tN9jqwDBiNG4W3g\ndarsckMywRB/ZkzCxZgEdb83YnwBjEfm4EX0u1YzlJ9NmnbErZHyYi6Sy4aYe2w0MuFCOxNFKq+p\nDtc6H0BmIYkJJ38+0qNOQ4bHt5FF6ABwFrAI+Cp0/3Wq3rRDNBPpZeM2SH5AFuZJOViifjbSa3Wi\nfxDBNGF/qGw+8tJ+bJFHuehT13pDva1GJlxoZ6IaNNXhWucGRvRKRPhjW4csMq8GPonct1Jdu9W1\nHjgVWVgPa/y2AGcAH5aI73o0CQ4+i4b6BeoaTNPqgSZgN7LdnoZljN1wiWM3cjCqI5gCTjXU22hk\nwpV2JiqtqQ4XOocZh+jwUxrj4K85H2nqpgC/MHraMR7pZUyn8x3q/kfTJJOBl1TcxzR1TcARZLs8\n6GQa1P1ZjieKuFuz1SAfXI+mzlYjE9WqnQlbTXW40DnMLOWvI41xkzLeh4gRMBU5fxgG+hk9vdmL\nfhHdxsiLdVOaZDKwQ8X9kdG7TZOAbcj65vqIzR5lc7PGXzM5nf6HaEfyOTtSnkYjE9WonYk0mupw\nqfNi5eteC5tRbGNkgV1Apiy9iPiDjF0LLFL3H0e2XwvIOqIP+FnVxW2Zu6ZWxd6FbH9/j5wdPY/0\n4kPAfRq7BUj7hpBznALy740u5GAyb25Bnt09mjpbjUxUm3Ym0mqqw6XO65HZQeoD/ibkwfcg5zXb\nEVFakIe/RmOzFOlx+oHfkWG1BVms9qZNJCUXIHmuRnbRNiJnHkeRnr81xvZKRIAjSFsOIWuA28qX\nrpE65FnqNi/SaGSimrQzkUVTHS50nox0AKZ1t+cEYznykkUPcD2VZwmiTaUO+D2OqUdGGNOhs6cy\nNCBHFe1pjPNe/HuSMYCsUyYgW/n9lU3HozgHOWAvIH9U9ng8Ho/H4/F4ys1/BBZRXZF7uVsAAAAA\nSUVORK5CYII=\n",
      "text/latex": [
       "$$a d - b c = \\left(a' d' - b' c'\\right)^{2}$$"
      ],
      "text/plain": [
       "                           2\n",
       "a⋅d - b⋅c = (a'⋅d' - b'⋅c') "
      ]
     },
     "execution_count": 73,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Eq(A.det(), factor((B**2).det()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also calculate the trace of $A$ and $B^2$ $$ trace(A) = trace(B^2)$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 74,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAOMAAAAWCAYAAAAy0SvmAAAABHNCSVQICAgIfAhkiAAABZVJREFU\neJztmntollUcxz/V1FkOrDXNaowkpD+0bIVWFKihFRFdhKKi9lj9UVG6/6zIWKJZEGmEBYGxlOhC\nU2lKabeFJdlW77rNGhGv3WbN1Kzoomv98T0Pe3l97rd30+cDL8/Lc875nd/zO79zzu9cICcnJ+co\n5n6gEzgA9APtwNSKapSTNnmbD1O2AAtQY0wDNgC7gZMqqVROqhzxbV4PDALrK61IALYDD7ukjQMG\ngKuyUyc2dyHb3xaijJcNsiifBl4+GKrNj01ctWw5zzw/qagW/kwBLgR+cEmvQW2xLzONRC1wBxql\nvwH+An4D3gdux9s/bNt/HLAuPxukXT4t3HxwuLZ5aixFo9IVlVbEh3uBX4GxLumvAAXguMw0Enci\n+/0EvACsAJ4D9pv3rwLHuJQtoM5bFbAuPxukXT4t3HxwWLa5hZSdlYLszUb2xBRkJ8km4FGXtCdQ\nZ5gcQ75FNBvPQWFS+Qx4CvCdkTnfodwY4F9gR4i6vGyQdnmL7H0wkTa/CY2SvcDvaBrtRIvPKFjE\nM0QVsAj4DPgb2AUsRiN2H/BjRLlxCWqn0WimqXeQsRJ9w1kxdbFI3tkeMDKfckg736Q9DUwH2tAs\n8CewlcN3Cb1sYDMLeBmFdf8gu2wBrglY3guLbH0wUpuXhxg1wFrkVNuAX4AJwNUofJkAPBblayIy\nGo0wc4Fu5Bi1QAtwJhrB2zPUxyaMnS4G3gS+L5PxJHADMBv4Kn2VQ3PQPA85pNnrpAbgA+B1YA1a\nI80F3kFrpv0mn5sNQA69Ejn7HjTT9BnZ84CL0HGAW/m0ieKDkdq8vDMOAqejLddSHkQzwAKy7Yyr\nkREeApYZ/QBagffMf7/Nm2ZgfIg6u4GNPnnC2GkGaoRSVgO3oFF/H2pQgD/Mr9JUAbea/284pDea\n50zkeIWStLXo2+4GHjHvnGxgsxx1xDagCc2uNuOQ49/oUT5tovhg6m3+NZoBwmIRLUSYYcq5zXw9\nJt3vOKBo8gX9tYbUs5wgdnKruyVinRbJhqmPG3mbXdK7cD/WmGnSXgxQTyPa3u9EM1BaWFTWByFC\nm58ILAE+RL13oKyg3yxU9Kg0rOOvM3kucEnfZtJP89EpDeLaKQ5FkrOxEwtNuZ04H0iPQmu6Is47\nqQ2m/IYAddltfGlIHb0oMkJ9sNSYZ6PF90TgI+AltCg/BJyBwpZPfeSt4vCQcDpaSz2PDFVKt4es\neaZ+tx27yWgGynoDJwk7xSFJG5dzDwqvelAH2euQZyqaxdpxXk82mOeuAPVdhgazd0Po6MeI9cHS\nzrgOfcRsoKMs31Lz7PKRt8rhnYUM0eog141qtAlSYChGL6UROBVtHPiR9JoxCTvFISkbl9OMNlK+\nQB3RLdS2N2+KLulXmudWn/qqgTpk7/8Ca+nPcPTBUNhXepwW6+PRrpDXdO2FRfh4fRQadd1uL7QZ\nmcsCyCqSXNiSpp3iYBFvzbjYlC8AJ/vkfcbkXe6QVodm0178LwOMNXKyiGwsKuuDoagzgnuNEja1\n6KxnEG11V0eQbRHNUXbivDi2HWcQuDaCPnFI005xsIjeGZeYsl0Eu7S8w+T/Fjih5H0N8DZaP18e\nsO7PjazrHdKmkNztFIsR4IP26NWPzobmIGO/BUxCV3w6UBjRgw48s2IFivHb0LpsNzLmNDQD1ZP9\nndThaKc4NKHQegBtRix0yFNkKFqoQmvmAnA8QyH9GOA6FLY14xw5OHEf8Bpq3ybgSxRhnIuOjiaF\n+5zEqZgP1qEzon50o2Q72oxoRCPAmohyLaKP2ovQCHwQ+BkZpRFtDuyJqE9c0rJTHCyi2bgF/7C9\noyT/Oebds2gHcT06kD+AIoOw9QNcgo5Q9qJ27kPrzZsjyHLD4sjywZycnJycnJycnJycnJyjjP8B\nXeIpfJ4LBpYAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$$a + d = a'^{2} + 2 b' c' + d'^{2}$$"
      ],
      "text/plain": [
       "          2               2\n",
       "a + d = a'  + 2⋅b'⋅c' + d' "
      ]
     },
     "execution_count": 74,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Eq(A.trace(), (B**2).trace())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "But from the equation for the determinant we know that $$ ad-bc=(a'd'-b'c')^2 $$ or $$b'c' = a'd'-\\sqrt{ad-bc}$$ so $$ a+d = a'^2 + 2a'd' \\pm 2\\sqrt{ad-bc} +d'^2$$ $$a+d = (a'+d')^2 \\pm 2\\sqrt{ad-bc}$$ so $$ trace(B) = a'+d' = \\pm \\sqrt{a+d \\pm 2\\sqrt{ad-bc}} = \\pm \\sqrt{trace(A) \\pm 2\\sqrt{det(A)}}$$\n",
    "\n",
    "Giving 4 solutions. Using these two relations we can calculate the expressions for $a'$,$b'$, $c'$ and $d'$ in terms of $a$, $b$, $c$ and $d$. For example, $$b=b'(a'+d')$$ $$b=b'(\\pm \\sqrt{a+d \\pm 2\\sqrt{ad-bc}})$$ $$b'=\\pm \\frac{b}{\\sqrt{a+d \\pm 2\\sqrt{ad-bc}}}$$ etc..."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}