3M1 Examples Paper 1 Crib.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 3M1 Examples Paper crib (linear algebra)\n",
    "\n",
    "\n",
    "Garth N. Wells (gnw20@cam.ac.uk)\n",
    "\n",
    "(C) 2016-2020"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from scipy import linalg"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 1\n",
    "\n",
    "The usual rules for real matrices extend without modification to complex valued matrices.\n",
    "\n",
    "a. $\\det(\\boldsymbol{A}) = (4 + 4i)(4 + i) - (2-i)(-3 + 2i) = 16 + 13 i$\n",
    "\n",
    "b. $\\det(\\boldsymbol{A}^{H}) = (4 - 4i)(4 - i) - (-3 - 2i)(2 + i) = 16 - 13i$\n",
    "\n",
    "c. Inverse:\n",
    "\n",
    "   $$\n",
    "   \\boldsymbol{A}^{-1} = \\frac{1}{16 + 13 i } \n",
    "   \\begin{bmatrix}\n",
    "      4 + i   & -2 + i  \\\\\n",
    "      3 - 2i  & 4 + 4i  \n",
    "   \\end{bmatrix}\n",
    "   $$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Extra:** [SymPy](https://www.sympy.org/) can be used for symbolic computations:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAJEAAAAzCAYAAACXDOrGAAAACXBIWXMAAA7EAAAOxAGVKw4bAAAG/klEQVR4Ae2di3XUOBiFZ3K2gEAJQwdh6WDoIJytYKEDclLBHraDTQtJBwkVQOgASoB0kL2fkSa2xvZYT3tAOkfxS4/7X139kh+arFer1aniXnh8fHzYO1lP/NYMrNfrXq2ciJV7xR9O/PBbs1WNH2Lg2tEJurlGRITX8jzrVnz383T9Wxl4YkD66OhEVxqdWBE9pax7lQFPBpKLSOPmmeJ/njiyJT+ER9ffLwGvMGzAYeK9treKZ9mImViwMBzkJ7mIhI1x8/lEjCWSHcKDS571JkINtRGGDxou3pn4UsdfFBHTVts5w0F+/kiJDtXGlmdI24jMqwRlHcSjel7E1pMgPzcyf7fLEa4LcfFW5+gEz9rXSu5P4SeZJ5LBuF56dGyv5jaSGBUS4onCMTEz3uaHMLt23+n8qc7jqRYbkolIFv6VwnskZGppeMZMQyzfxN9QB3TFNVZW8WtJhjP1FIaNJU2mR/GYHs9QwVD2VY33b3HmWxWq/jetw/ZuM7HWdeZHxYIvP9GeyLjaBxn6rZiVIxVNxHNphMN843KkuNkuyQ4ExDB2MQMIL36iRSQDuaOIngQnJGoUjxHZJ1MfHmAR4u+xH4FfGbH3XM5zKoSfqOFMFTIkBA1jyks+JpRuaB4P6HrfU/MvInXI9a8m4vmuMm5MpeCfo6e7NneODTd3wtnHgZuW+dJHRZ950xuVPTREevMTLCIZiqs9FZignjxEkMo9V7nc4nvNU6biUbnN5NXUo6pWS/KitiOshvgBcDsYe3iulCSE8BMsIiFGRK/UGLjddmjGcnOeO45SPd0XD3OhG0ta24C59o2wXwjTzgPpHHYhqqDOGmGLFz9fVdFWIAEaHVVW82Y3tCzlxxO9D83v5uvDo3O4/kfFc9Ir8LQ42vaYMoSBzreHQ+e408TjF8NHfYoH+VEapgO3MZ5I+XsDAIhLCX14/gScGubG9HQ60mzBYMCj32m/PcdkfngmnF5DewJDvPhJJiJjfFO5jNjq+Fbb2xkIaDgcwyNMNBaRXs7jibnnRQiIYYue7YahCbCbLtmxLz/JRKSKd+N4MmsiCjqER9dfRxSfNKuwJJsYpwLmw89JqkoTlvOgsog1HAkDyTxRKnvVA3iPVMMRMbBET3RE9FWoMFBFVHUQzUAVUTSFtYAqoqqBaAaqiKIprAVUEVUNRDNQRRRNYS2giqhqIJqB7A8b9X5qI5T2lQgvQznmbfUsDxUNHvt5Cu/6viteCE/vOyrzfq3zeYbSZw+ql7f6zTq07JWNVDDF/qwiEgBEQwNZEdmPrljdydd19gvDETPSXVKdVsC7ryN1jjVfLBJknXmfsMFeFKexmJeyvcI210ttJtmf9HsiNcTuuxdZyVvy3XcpXFNAWOzct9NO3Ve+reLbqenb6ZSPhtn7Nkfnml9Faaedc1944I12uQ7FobzBPE2tU3U03xPlnhPRk3iZSmyCAO727TnPLSIkhgSIXfQiQXlEhjE4mpMnL26zikiC4WPzZ2wtKpF0bvbbH1/Zy7m34Fj6IsFjWnTZtFfWOZGrCAkIT8AcZHRZj5sv1bHEvJsLOWXS+xmGmzmIcOLpcNVFFzeqXoaxOTqXqn0KvvZn9UQWlkDxczMQxCSNhvpsr829BZswMOG2d2xA8lq8R4bYIBxg4CvL0h/k90H3sr+IJzI93PZyhjPuhorfnfWxpXNMtneLBE1jfjJpSy5uxDu3hTwAN+/pEPtHRaQCcesfFdlODWML4+zH8Uwar1U+8yX294Ku4dYZ/tzwnBO6vnts0Eowurixla7ZNXW4iwS9F++55foeCwdDZ9AwloEnb/tHRWQaOPj7XxnYmWu0yGU4QyDE3mcwqrtPJIgHT+a9uFF5OsE0HKLu1GNFbeohT9aP+FUPwxiPHYKGMRe/NTKUpxD7R0VkAYVsZQTei1+mpeEHPU5I2bF5DMGdp9A6R2MiKtuYlzossbiRepe0CBQaCJPtzyYiFK2GYahiuHCHLF43EHa3/j8P8/8VJrzjK2Fy5x94uMbrKA0dgHT/KNIJeE3jpudSdFC5cLDHg+rkASjcDd1RRtc9VICv/dlEZADuES+ANBaNxETSFdeQXUnOq256PRPpQ4sEG5EL35yLG+GIOEfwsj+riNQIV2o4FjK2J4005NB7qtyEISDqH10kKNyzLW40XDWNKJyzLAL1tT+riFAEgLTZc9dcKx2EZfJNgtLOsrhR9XYm+qU5svX52H9iMx3RliGw6DB4RNy0oRbjKbsnaluVYl89ZBFeLYUtOcsoydMxeqKc3NeyAxioIgogrWbpMlBF1OWjHgUwYOdE/HOSzjMJjal18hpA6K+cxdWIbG00s9YOn2FuHON5q72IW00HVz2ckQGJiB8u2zoQ7v4HBHO8iZT6LMMAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}4 + 4 i & 2 - i\\\\-3 + 2 i & 4 + i\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "[4 + 4*I   2 - I]\n",
       "[               ]\n",
       "[-3 + 2*I  4 + I]"
      ]
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "import sympy\n",
    "sympy.interactive.printing.init_printing(use_unicode=False, wrap_line=False)\n",
    "A = sympy.Matrix([[4 + sympy.I*4 , 2- sympy.I], [-3 + sympy.I*2 , 4 + sympy.I]])\n",
    "A"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAEgAAAAQCAYAAAC1MDndAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACwElEQVRYCc2XgVEbMRBF7zIpAFICdACkgpAOgFQAdJAMFWSgA5IKiN0BdEDsDkgHSdyB854i3ciXO1sOGG5n1tLtrlarL+1KrufzeZWorus9+iN4H/ksydstdpct2Q3205bs2T5XxY1+h2DOY0BbtH5fEvNdlC002H9EsIs+jHGAoFzDE1jEtlBWbUau4wf4MOnoO/YhfT9Xy5xFcUe76zwuZGew6zzK5amP3DUKIFYZEAhFbhlAAtie7BbZbe6ntM+4Q/is1L7PDh+9cWe6Bgxkgus6J30+k/w1RkXEsTvC0BR8lw/A0fv8e82+gcqbJFN/FjnMQ8wz1lM0ZzFAeLtwEp0XeR6IEfFaZ7bzcOJmK7KsLKV1APL0THFu+wH+Ce/CoxgE3eET8ZvWXjLnxP0lRYzck2xtck3W1Ct1RQDFwdpLBwz+9LdbVeh+w6fIxkk2xJYY3VjBeQubdt/hnC5cF3baePEEgF7lFgX9PZw0qEf7b7RfcbzpWlIQXr8JcU/hK/gYqxt4QszWVTfZ2/nePqT+R+jxU3SCcJqKWjMwOaD1ZvNoHsB97wpz3Z1p0xsFBJjeKLneBRnskxN+x8xpLR3RWp9+KYsTuZYmQ4oAigMFZ1mBdhc6icm7ABAYd3AHfTjOnYMfKWQOU8vnTPsha4q5ab7pAjgxHs2bLFkXoF4QcNp1upzsxYgFm/aecDdjGyCWbbBm3tTj3G6dGmTl7wJoH7nXf2d6oXsxigsVlIVFx4AsCVKIO4LpabM+CajrrYoBigDcpYHRiTt0Ap/6PVCynvjab4g1mNrG7lWfTlUAjG/rkwfBvxtVjUC0vNYcoJGt6WLO3qNfqA8RIG0ki+xnbNr5HZSrfmKg/12DSuPGzlqTF3wB+OfPKnbhbxM6MyLUoQDQqoVsSv9YgDYVV+63OMXyQU/Y93jLg6U/+bwzMvE+y7YAAAAASUVORK5CYII=\n",
      "text/latex": [
       "$\\displaystyle 16 + 13 i$"
      ],
      "text/plain": [
       "16 + 13*I"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "det = A.det()\n",
    "sympy.simplify(det)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAEgAAAAOCAYAAACM7Fo2AAAACXBIWXMAAA7EAAAOxAGVKw4bAAAChUlEQVRIDc2X4VHcMBBG7QwFQEoIHQCpINABhAoIHSRDBRnSAemAuesAUkECHUAHCdfB5T0heXSKz4gf53hn9iTtruRPn7UrX7tcLpskbdvu0Z+h+9gXyV62xF0Wtmvi7wvbaMOXcON/B5jzCGib1vElmG+jbaUh/jOGXfxhjhMk5Qq9Q2VsG2dTKnYXfkAPk4++cx/SeKyWZ1bhjnFXOS5sn1D3eZzbUx+7e5RAojIiMMrcEEESWD7sBttNvs7Y/SHcma8jA5vkus+7l7BuEVQlHLtjAk3BD/kEHnCUjyfYN/UXUQM8MC/YTxXUaoJY7cKHuHjVyhMJAq91ZieHE1+2JsvKoLyGIE/PPYvbnqK/0V10FkHQnb6A/xCUXjLn4P6eEGM37axN7sma+k1fFUFxsvHKAZO/PHebBt8TeoZtnmxTbMHoi5Wc96hp9wvN5cJ9EWeMF08gqKpIE5yKGmus3m74PKZPxpS+scY8e/ByKXEQbz3tbjH63s6hiNOG2zzNqTpBBKei9sgCpXizeTQP0HXfFRL8A7WtlROeu5FvK9adc1IWAJnRWp/+aIvA3EuXIVUExYmS46LrxLfQKxKMY7/XuWFjTC0zpSTbFDOd/KYL5BDryVK62vRagtaSwKJ9pys87X/9sGFPrCfcWrkTX9QQHG/qeR73Zii68Fn5+wjyZHj996ZXscaow7hRT+/KpiMIS4IScEcyLeTXGhmHv1PVBEUCbtPEuIhv6CN65niiYj3xa78T9mAqid2rXgKVQBhj65MHwb8bTYtBtrzWnGCQrelizv7E/3zdMVAiQcYob9GvxJT5HZyb/qnFTZy15iTDIwH//FklLvxtwmdGhDr0FysfvJ8quWxFAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$\\displaystyle 16 - 13 i$"
      ],
      "text/plain": [
       "16 - 13*I"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "detH = A.H.det()\n",
    "sympy.simplify(detH)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAANsAAAAzCAYAAAAelH2vAAAACXBIWXMAAA7EAAAOxAGVKw4bAAALi0lEQVR4Ae2c2ZEUORCGZ4h52EcCHvcNPODwADzgsADwAAILCPBgWQtY8ACwYAEPwAOW8WD2/9RKjaq6uqS6utWiFKHSrcpDqcxSSTo9OTm5Kr/lLi4uzrcyo4zT09MbqvMjympEU+WNyhUlhDf0vNZHm4rQ/W1R6eOzL9uizalyvsvfaJW81WB51sprJNXhX74dYUPo1PZbqrzRWQEJwQsNHsgzySAwd+VfgYvC4FTvtRI/fca5yt+GQkWsXPlv4vxjjQsfaPFI/plwut3GQ+XPlXddHppcV50X7TolpwX/1xgvj+9LD7PJxRPVYVwEp3puHCh/i88q+6iK90LlTeQTAcJ2T41Ohni1ocOLDv8X/aTKh7xrH3UF73v51/YuxRlEvyxNKPdV/paP31Jc0WF0O6b64CoPHfDf27Arr00zJqv37Xqlpj1eAu+Sh8pz49fySMtv4W7lOaHaP5X/eEWPse6HXnQae3UUa8RU+dj3LtXulTp+F3XObB00tp/Jvghfp+l8eD+qX10UHOWZuQMdWkgiXDHNmL0feO3QqlpWUjCiebrweurLDGA02A3lMfFMclOEDSCCEzCkYxMiVR7alhDxAys2GRlIzGrmmJ3eQ3R5N6DUxpkGVFDePfnn8sz21TvhaSZWMK9ED4vfKZkAgh3TGOEJ/Ivg5fPpS5RuRNV2NJ9HC5sIG2YFAcDA/BgRG/uqt7yBQUEJ4YLAYC5iTrjvMaVhDt40GUx6rXwE0Bxl1D+Kmd2AHht6/iJcJnRMODb7h7yx/S/c7qng3/rW4p3wXB68zCF8WGk2EY/m85n1ODb0A5GPZxuIja5S5Y3KBSSExxvBjNCgxSA+TLnmQcN2d0RXGVr8l0KbHTFDnWmiOjGzfNMqg4fC6oVoYJrANNp/pWIrWFEMH3LgU10mD+q7hSGlmURG83m0ZouAfaR43+BKlUddlRH1woIwob1YEDEtbSFCaDizuMTMR5pZMDY9y0BoISiEMxMNAsckg6D9I48zLbBJTXyKB5husRUxqkcvLFm/ZXxdPoVuG6+n8nkOYYPYYRB2UCFV3tFkv1ki7FV5tJSZQQBgwmSzNTh2mUcOd/pQOYPubasfZdXpDE8Nwg9e8KAViyp942EMMaAtfqqjj5uCm0kUQXLrCj6NBnNOafjMyvR9+XOl+U53vFc4ms9zCBsE5vfBLpcq39Vub/kQVC/D9IkHCYJDvs3WMCaYyjCA+mprZiT1GWi0eSz/OzhoguVijvQTS5QWijfw54V5weeETSH/U51p6YUK6+QVPPZ8xmIx03g0n8/UyVQHEGazd/WVKu9qc4g8NPBLEdd+WN9VOjYh0FhoQBhEnZuUy5tD6B6rnH9SVZiSwoXZnIHGAGP5G7z4brXFBWhx3+PMrxJ+/n5TWLwTzJil8Bz3t9LvvMCxOIb2IgxOZdABN5rPtoOEBQ6boTddrs+VAoVQQIKAicf2QBPyQiDLA8MEew4zMu+Na62VAr85BZJmpKTS7VNakk6asdCwO51gQK1/lifMdQ/Vb9KkKQG/XIR21ZuDPiXQQTBgpmKytp379aJyM+Xicr7DzByM80O8BNwAJilsKUEIGC0YEQwsOsTfR7O9rQT8piIzB31KoINg6BImfpZPMiNLwA0er2bk1JG+tl8pkEmBpGbr6kczzaBjCapvq1qcFMCRjpfNXWYpj9rxy6VzzXQ4FG6DjtiIUSxtSzPnH0tQfezwX7ST533hKEvcTwlxwVY1frk0LokOggUz8nku7Kl6+8ZN7+M3g1M02cKmBggNiAufhrAhROFMnOJoLvLs7BftXDxuV1pcMFaNXy69S6OD4JlN2A6Bm97phC37m01ql5VAdk10/Y/jw7bvx7aKwzGUInfF146fY0DGo2Y6lIBblmYTn5waV4jQaZK81GztuMrZWRBOtyqOxuCcFxqP9sQftNsdMi14qsYvl7Yl0kEwMX44FsPAG+0PhRuwy+eZkapoS6+K9gubytF+CPDVXYTxdVQ8nnBzthU8VeOXS6ua6XBI3PTuPDNSqhdtNPpYgtrafxKi5mynPDPWQV3t+OUSt2Y6lILbWQYzMPvcsQRflzQChKn4r2bNeLe0O5bgy9Fw3D6FYHEQk029tqPD7QhQfrzLnmaHcLXjl0vTmulQDG5Z32xmiohzttIYzESfh02KgJln642ro7Cx1K80y+tF3sIkuKrGz/iYCmumw75x0/ucGTlo17+0E43Yh4b5h0ZzxxKUzz80p/EUBieGuj2PKqfspTxHUziKwUdufDkQWQd1ghEhA6+p+P2pPv4Qfp1bjw6KZObLK+fz3sew0XOQsMW8igSIbAYqLlxmmSrfVC/nKXjRxOBB2DBvJTjO/PU4PVJ554Wlyjfz+qfaFHkcJIVDqjzGkbgcnwqNi2o32Yd79uHgy5j4cVvjdpM9Lx/1TifgOd9s9v52iGkYZm91yCDlwB2HKnGp8k2tcp6m2dBusWMgPRN+mMdWtqXFrYFoUpTGNrgIUzikyn0f8NgdEvX1SRcjbBk4ZI3Lpfg46JtNQLj/HCJwatdIb7n1U0oofBo34QLXjjx+E4R/iKXAPwQOwd+Lw65y5bMo1qCT0mHn0BAYlq7bg8Pex6Vgcd9sUzRbatdIqlwwFOUYSMFphiSdraVUH62H9rurgdR7viq85Pgizhzy2gNL4JNw7dpRVDJmveNyST5eGUsVEbn3MstU+dj3LtVO8IbvNBGcmb9x6WzGe0df3pnR98GriCaYzni79Agha19Ue3A4UwBkjMvF+DhFswW8/EzHAO084JkqDx0VEPGDauels10gqg2z/Ct5tJvdIdlV9ZjzrnngOy+qjSerY0GyPS6X5uNozWYE9QBicoWbqKyMMFUe1y0kzmrj+RBYGGjytMFEYaGoOhcJU7AAPM7gyiQz2mmMzHIJ6xAAusbl0nycJGwGsIDcuswSxFPlQ4izx7p8b4UBlfte4YqJxaCr+ZJW6IIWb7vB9Gp1AO3we3F943JJPo4WNg8ws3jnZZap8r1QddxL7qgZK7RDHYJW+yWtWDD2zWa/Eoo9cd/FwIxxuSgfxy7928lrllKDl5azXwO95VavtFC4QI+tg67KY0ZnsPFfCXyZaMLpYcWZmTk6xFY07jgcfRRkqbbAJd+HQ285cMmBH30QQoOwbW8s3OqD7/1Ay7H9ePh6cdB7escl+MjPykf155b+R+8gUQerWykwCwWkbRC2o72ENUUE4ed+mVxJVVzLVwqsFJiHAmepbiSVznZI1ZtSLvXvNizv6kMwoNo/yw/5iP5tLmndRbc4vxA+YnbyTdR27reCYGQ1t+2O5hLWNuDtdFLYUoLQ7nCJtGA4V7+d//Cmvq8E/KbikNO+BDwFQ5cwsdAyyYwsAbccHqxmZA6V1jorBWagQFKzdb1DM9F6SWvzOBErYMza7lIXhaSLXxKvmY+l4jZo6V+DiCVfae7LpW2l2zvBsc3DznjFsdNtyZX3ccyhuKVxYJKrGj+je0l4CpbZlv5L5KHwc0v/Ct1/paxjEqqL0EAY4XQpLKQpszzFmdnJWy9pjehk9Dl0CK/ki+Gjh2Wu/2xF4Qav5ZywZX+zSS2zEsgRkq4jFZhQ6yWtmz1+RV5CK/44VzMfjwG3LDNSnKr6EtPa8TNtWiKeggltdLSXsBptd4XgJr9e0urVvC09iybuX56yL83kOK5ytDsT1M5tSr7Ozj7i/vYZF1zV4lkyboItz4yUaub7a72kVUTwtGBfYOM4kfIZxLFzu+CVz4xdhKuZj8eC21nGSOBbbb2kdTPpsIrqdr2LwWg4bpZCsEq+hNZYXDMfjwY3TCLUHAAHv8u8UR1baQxmlM9bL2n1pqfowe+DIi+hNb7WzMdD46b3BznyccbDR9v1jwDFjvtFtrbWaDZHKKdeYlrkJa0gXzt+xuCa8SwBN8GA4ml/Qnz6H1R6mj0i8nBKAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}\\frac{77}{425} - \\frac{36 i}{425} & - \\frac{19}{425} + \\frac{42 i}{425}\\\\\\frac{22}{425} - \\frac{71 i}{425} & \\frac{116}{425} + \\frac{12 i}{425}\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "[ 77   36*I     19   42*I]\n",
       "[--- - ----  - --- + ----]\n",
       "[425   425     425   425 ]\n",
       "[                        ]\n",
       "[ 22   71*I   116   12*I ]\n",
       "[--- - ----   --- + ---- ]\n",
       "[425   425    425   425  ]"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Ainv = A.inv()\n",
    "sympy.simplify(Ainv)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAADUAAAAzCAYAAAA3v0BtAAAACXBIWXMAAA7EAAAOxAGVKw4bAAACyklEQVRoBe1ai03kMBBN0BWwgg72OuBTAUsHLB3AlXCiAgQdABUg6AA6QGwHXAeH6GDvPc4TmTC2E++uM0IZydgZx37veZzEO6KuqmqC8sWWy+X7F6cxR13XKnfSfEVZtso1RFXWCzg/tnhTx+MW/tCOIKD2yq//btt/wfcTb7D94C2ibLPvye5H1/uxf3dx7z3KHlao2PMG3EvH8S/qnyiXwP/jfGoVFeUexFuMfEPZR5mqs2zICfwXTH0BEQ+EcHxeUHPbBYVFtx8jgjJH4V6948SlDMTPgDUB9ocg4pIPKl5f8zpkUVGhQYX8c+AsFKxn+GYuakp3VVkWNQNjbvu2yfPMftVMiopFwVOx7bU/NU2KAkMhLFHxSUv0Jr7Tb1sV5XMMtXdCHVZFSTQ03hJFfrdUMynKvbpJWNti4sv7TqnLUM75BKipAieRYr9qJiPlmPJIxlNM2/bgWHjRbPfb/U6B9A3YvuH1fiys3av+BNen4tPq6NmPAzARV4z7WFaNZy9+6Z8BfIV6k8aonANPtuEBrg+BS/ygJUVhAh5XBjFgvwP4d19wy89UXy3N/aOoZimMN8ZIGQ9QQ2+MVLMUxhtjpIwHqKH3LSOVPCZRPs5evROKzbKtoQH8XUzTOZGaFIUJsxKKq2oBLg/RWYnU6PbDxNkJxVVF8TDLwzRK70RqVBSIZScUVxW1yviUqBkm15Ig/ElAY785C4pyezpFWPIFqfuK9gdFgYUQlqj4xCR6fJjNWUxUF7LBhGKXwZu6JyZKoqFhSxSDCUVtUClfUJTLD5CHtsXEF0wolhKg4QRFuZuzE4oaWClfSlR2QrGUAA0nKgpbMDuhqIGV8iXPfiCSlVBchwB8K7MSqUlR7oXRO6G4DlHAzkqkRrffOogNMccoaohVz8EcI5WzakOMkbfftP1TwzsmDcGrE2abMwbJ8e37/RPjP2U+4S7Q5V7uAAAAAElFTkSuQmCC\n",
      "text/latex": [
       "$\\displaystyle \\left[\\begin{matrix}1 & 0\\\\0 & 1\\end{matrix}\\right]$"
      ],
      "text/plain": [
       "[1  0]\n",
       "[    ]\n",
       "[0  1]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "sympy.simplify(A.multiply(A.inv()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 2 \n",
    "\n",
    "The diagonal entries must be real since $A_{ii} = \\bar{A}_{ii}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 3\n",
    "\n",
    "\\begin{align}\n",
    "\\det \\left(\\boldsymbol{Q} - \\lambda \\boldsymbol{I} \\right) &=\n",
    "\\det\\left(\n",
    "\\begin{bmatrix}\n",
    "\\cos \\theta - \\lambda & - \\sin \\theta \\\\ \n",
    "\\sin \\theta           & \\cos \\theta - \\lambda\n",
    "\\end{bmatrix} \\right) \n",
    "\\\\\n",
    "&= (\\cos \\theta - \\lambda)^{2} + \\sin^{2}\\theta \\\\\n",
    "&= \\lambda^{2} - (2 \\cos \\theta) \\lambda  + 1 \\\\\n",
    "&= 0\n",
    "\\end{align}\n",
    "\n",
    "Computing roots,\n",
    "\n",
    "\\begin{align}\n",
    "\\lambda &= \\cos \\theta \\pm \\sqrt{\\cos^{2}\\theta - 1} \\\\\n",
    "        &= \\cos \\theta \\pm i \\sin\\theta\n",
    "\\end{align}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 4\n",
    "\n",
    "a. If $\\boldsymbol{x}$ is an eigenvector and $\\lambda$ is an eigenvalue,\n",
    "\n",
    "   $$\n",
    "   \\boldsymbol{M} \\boldsymbol{x} = \\lambda \\boldsymbol{x}\n",
    "   $$\n",
    "\n",
    "   Premultiplying by $\\boldsymbol{x}^{H}$, \n",
    "\n",
    "   $$\n",
    "   \\boldsymbol{x}^{H} \\boldsymbol{M} \\boldsymbol{x} = \\lambda \\boldsymbol{x}^{H} \\boldsymbol{x}\n",
    "   $$\n",
    "\n",
    "   Taking the complex conjugate of both sides and noting that $\\boldsymbol{x}^{H} \\boldsymbol{x}$ is real,\n",
    "\n",
    "   $$\n",
    "   \\left(\\boldsymbol{x}^{H} \\boldsymbol{M} \\boldsymbol{x}\\right)^{H} = \n",
    "   \\boldsymbol{x}^{H} \\boldsymbol{M}^{H} \\boldsymbol{x} = \n",
    "   \\bar{\\lambda} \\boldsymbol{x}^{H} \\boldsymbol{x}\n",
    "   $$\n",
    "\n",
    "   Since $\\boldsymbol{M}$ is Hermitian ($\\boldsymbol{M} =\\boldsymbol{M}^{H}$), \n",
    "   comparing the above two equations they can hold ony if $\\lambda = \\bar{\\lambda}$, \n",
    "   which is true only if $\\lambda$ is *real*. \n",
    "\n",
    "b. Consider two eigenpairs $(\\lambda_{1}, \\boldsymbol{x}_{1})$ and $( \\lambda_{2}, \\boldsymbol{x}_{2})$, where $\\lambda_{1} \\ne \\lambda_{2}$. We have:\n",
    "\n",
    "   \\begin{align}\n",
    "   \\boldsymbol{M} \\boldsymbol{x}_{1} &= \\lambda_{1} \\boldsymbol{x}_{1} \\\\\n",
    "   \\boldsymbol{M} \\boldsymbol{x}_{2} &= \\lambda_{2} \\boldsymbol{x}_{2}\n",
    "   \\end{align}\n",
    "\n",
    "   Multiplying by $\\boldsymbol{x}_{2}$ and $\\boldsymbol{x}_{1}$, respectively,\n",
    "\n",
    "   \\begin{align}\n",
    "   \\boldsymbol{x}_{2}^{H} \\boldsymbol{M} \\boldsymbol{x}_{1} &= \\lambda_{1} \\boldsymbol{x}_{2}^{H}   \n",
    "   \\boldsymbol{x}_{1} \\\\\n",
    "   \\boldsymbol{x}_{1}^{H} \\boldsymbol{M} \\boldsymbol{x}_{2} &= \\lambda_{2} \\boldsymbol{x}_{1}^{H} \n",
    "   \\boldsymbol{x}_{2}\n",
    "   \\end{align}\n",
    "\n",
    "   Taking the complex conjugate transpose of the first equation and exploiting that $\\boldsymbol{M} = \n",
    "   \\boldsymbol{M}^{H}$,\n",
    "\n",
    "   \\begin{equation}\n",
    "   \\left(\\boldsymbol{x}_{2}^{H} \\boldsymbol{M} \\boldsymbol{x}_{1}\\right)^{H} =  \\boldsymbol{x}_{1}^{H}    \n",
    "   \\boldsymbol{M} \\boldsymbol{x}_{2} = \\lambda_{1} \\boldsymbol{x}_{1}^{H} \\boldsymbol{x}_{2} \n",
    "   \\end{equation}\n",
    "\n",
    "   Comparing this to $\\boldsymbol{x}_{1}^{H}\\boldsymbol{M} \\boldsymbol{x}_{2} = \\lambda_{2} \n",
    "   \\boldsymbol{x}_{1}^{H}\\boldsymbol{x}_{2}$, since $\\lambda_{1} \\ne   \n",
    "   \\lambda_{2}$ both equations can hold only if $\\boldsymbol{x}_{1}^{H} \\boldsymbol{x}_{2} = 0$, i.e. the   eigenvectors are *orthogonal*."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 5\n",
    "\n",
    "To show that $\\boldsymbol{A}^{H} \\boldsymbol{A}$ is positive semi-definite:\n",
    "\\begin{equation}\n",
    "\\boldsymbol{x}^{H} \\boldsymbol{A}^{H} \\boldsymbol{A} \\boldsymbol{x}  =\n",
    "(\\boldsymbol{A}\\boldsymbol{x})^{H} (\\boldsymbol{A} \\boldsymbol{x})\n",
    "\\ge 0\n",
    "\\end{equation}\n",
    "If $(\\lambda, \\boldsymbol{x})$ is an eigenpair of $\\boldsymbol{A}^{H} \\boldsymbol{A}$,\n",
    "\n",
    "\\begin{equation}\n",
    "  \\boldsymbol{A}^{H} \\boldsymbol{A} \\boldsymbol{x} \n",
    "  = \\lambda \\boldsymbol{x}\n",
    "\\end{equation}\n",
    "\n",
    "Premultiplying both sides by $\\boldsymbol{x}^{H}$:\n",
    "\n",
    "\\begin{equation}\n",
    "  \\boldsymbol{x}^{H}  \\boldsymbol{A}^{H} \\boldsymbol{A} \\boldsymbol{x} \n",
    "  = \\lambda \\boldsymbol{x}^{H} \\boldsymbol{x}\n",
    "\\end{equation}\n",
    "\n",
    "Since $\\boldsymbol{x}^{H} \\boldsymbol{A}^{H} \\boldsymbol{A} \\boldsymbol{x} \\ge 0$ and $\\boldsymbol{x}^{H} \\boldsymbol{x} \\ge 0$, all eigenvalues $\\lambda$ must be greater than or equal to zero."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 6\n",
    "\n",
    "a. Eigenvalues of $\\boldsymbol{A}$ satisfy $\\det(\\boldsymbol{A} - \\lambda\\boldsymbol{I}) = 0$. \n",
    "   For the transpose, we have \n",
    "   \\begin{equation}\n",
    "     \\det\\left(\\boldsymbol{A}^{T} - \\lambda\\boldsymbol{I}\\right) = \n",
    "     \\det\\left(\\left(\\boldsymbol{A} - \\bar{\\lambda} \\boldsymbol{I}\\right)^{T}\\right) = \n",
    "     \\det(\\boldsymbol{A} - \\lambda \\boldsymbol{I}),\n",
    "   \\end{equation}\n",
    "   since the $\\det\\boldsymbol{A} = \\det \\boldsymbol{A}^{T}$. Hence eigenvalues of $\\det\\boldsymbol{A}$ and \n",
    "   $\\det \\boldsymbol{A}^{T}$ are the same.\n",
    "  \n",
    "b. Noting that $\\det\\left(\\boldsymbol{A}^{H}\\right) = \\det\\left(\\overline{\\boldsymbol{A}}^{T}\\right) = \\overline{\\det\\left(\\boldsymbol{A}^{T}\\right)} = \\overline{\\det\\left(\\boldsymbol{A}\\right)}$,\n",
    "   \\begin{equation}\n",
    "   \\det\\left(\\boldsymbol{A}^{H} - \\lambda\\boldsymbol{I}\\right) = \n",
    "   \\det\\left(\\left(\\boldsymbol{A} - \\bar{\\lambda} \\boldsymbol{I}\\right)^{H}\\right) = \n",
    "   \\det\\left( \\overline{(\\boldsymbol{A} - \\bar{\\lambda}\\boldsymbol{I})^{T}}\\right) = \n",
    "   \\overline{\\det(\\boldsymbol{A} - \\bar{\\lambda} \\boldsymbol{I}}).\n",
    "   \\end{equation}\n",
    "   \n",
    "c. If $\\lambda$ and $\\boldsymbol{x}$ are an eigenvalue and eigenvector, respectively, of \n",
    "   $\\boldsymbol{A}\\boldsymbol{B}$, then\n",
    "   \\begin{equation}\n",
    "    \\boldsymbol{A}\\boldsymbol{B} \\boldsymbol{x} = \\lambda \\boldsymbol{x}\n",
    "   \\end{equation}\n",
    "   Multiplying $\\boldsymbol{x}$ by $\\boldsymbol{B}$ on both sides,\n",
    "   \\begin{equation}\n",
    "     (\\boldsymbol{B}\\boldsymbol{A}) (\\boldsymbol{B} \\boldsymbol{x}) \n",
    "     = \\lambda (\\boldsymbol{B}\\boldsymbol{x})\n",
    "   \\end{equation}\n",
    "   We see that $\\lambda$ is an eigenvalue of $\\boldsymbol{B}\\boldsymbol{A}$, and the eigenvector is now    \n",
    "   $\\boldsymbol{B} \\boldsymbol{x}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 7\n",
    "\n",
    "Follows directly from the definition of the norms: \n",
    "\n",
    "\\begin{equation}\n",
    "  \\max_{i=1}^{n} |x_{i}| \\le \\sum_{i=1}^{n} |x_{i}| \\le n \\max_{i=1}^{n} |x_{i}|.\n",
    "\\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 8\n",
    "\n",
    "Recall that from the defintion of a matrix operator norm it follows that $\\| \\boldsymbol{A}\\boldsymbol{x} \\| \\le \\| \\boldsymbol{A} \\| \\| \\boldsymbol{x} \\|$ \n",
    "\n",
    "a. From the definition of the matrix norm:\n",
    "\n",
    "   \\begin{equation}\n",
    "      \\| \\boldsymbol{A} \\boldsymbol{B} \\boldsymbol{x} \\| \n",
    "      \\le \\| \\boldsymbol{A} \\| \\| \\boldsymbol{B} \\boldsymbol{x}\\| \n",
    "      \\le \\| \\boldsymbol{A} \\| \\| \\boldsymbol{B}\\| \\| \\boldsymbol{x}\\|\n",
    "   \\end{equation}\n",
    "   \n",
    "   for all $\\boldsymbol{x}$. Re-arranging,\n",
    "\n",
    "   \\begin{equation}\n",
    "      \\frac{\\| \\boldsymbol{A} \\boldsymbol{B} \\boldsymbol{x} \\|}{\\|\\boldsymbol{x}\\|} \n",
    "      \\le \\| \\boldsymbol{A} \\| \\| \\boldsymbol{B}\\|\n",
    "   \\end{equation}\n",
    "\n",
    "   for all $\\boldsymbol{x} \\ne \\boldsymbol{0}$. From the definition of the norm, $\\|    \n",
    "   \\boldsymbol{A} \\boldsymbol{B} \\|$ is the largest possible value of $\\| \\boldsymbol{A} \n",
    "   \\boldsymbol{B} \\boldsymbol{x} \\| / \\|\\boldsymbol{x}\\|$ (over all\n",
    "   $\\boldsymbol{x}$, exluding the zero vector), hence\n",
    "\n",
    "   \\begin{equation}\n",
    "      \\| \\boldsymbol{A} \\boldsymbol{B} \\|  \\le \\| \\boldsymbol{A} \\| \\| \\boldsymbol{B}\\|\n",
    "   \\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "b. Consider\n",
    "   \\begin{equation}\n",
    "      \\| (\\boldsymbol{A} + \\boldsymbol{B}) \\boldsymbol{x} \\| \n",
    "      = \\| \\boldsymbol{A}\\boldsymbol{x} + \\boldsymbol{B} \\boldsymbol{x} \\| \n",
    "   \\end{equation}\n",
    "\n",
    "   From the triagle inequality for vectors,\n",
    "\n",
    "   \\begin{align}\n",
    "      \\| \\boldsymbol{A}\\boldsymbol{x} + \\boldsymbol{B} \\boldsymbol{x} \\| \n",
    "      &\\le\n",
    "      \\| \\boldsymbol{A}\\boldsymbol{x} \\| + \\| \\boldsymbol{B} \\boldsymbol{x} \\|  \\\\\n",
    "      &\\le\n",
    "      \\| \\boldsymbol{A} \\| \\| \\boldsymbol{x} \\| + \\| \\boldsymbol{B} \\| |\\boldsymbol{x} \\| \n",
    "   \\end{align}\n",
    "\n",
    "   Re-arranging\n",
    "\n",
    "   \\begin{equation}\n",
    "      \\frac{\\| (\\boldsymbol{A} + \\boldsymbol{B}) \\boldsymbol{x} \\|}{\\| \\boldsymbol{x} \\|}\n",
    "      \\le\n",
    "      \\| \\boldsymbol{A} \\| + \\| \\boldsymbol{B} \\|  \n",
    "   \\end{equation}\n",
    "\n",
    "   Using same argument as from part a), we get $\\| \\boldsymbol{A} + \\boldsymbol{B} \\| \\le  \n",
    "   \\| \\boldsymbol{A} \\| + \\| \\boldsymbol{B} \\|$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 9\n",
    "\n",
    "Recall that the smallest possible value of $C$ for which $\\|{\\boldsymbol{A}\\boldsymbol{x}}\\| \\le C \\| \\boldsymbol{x}\\|$ holds is the operator norm $\\|\\boldsymbol{A}\\|$. The task is to find $C$ for the different norms. It simplifies the proofs if we consider all vectors for which $\\|\\boldsymbol{x}\\| = 1$ (this is possible since $\\|\\alpha \\boldsymbol{x}\\| = |\\alpha| \\|\\boldsymbol{x}\\|$), in which case the task is to find the smallest possible $C$ such that $\\|\\boldsymbol{A}\\boldsymbol{x}\\| \\le C$.\n",
    "\n",
    "What we need to do is (i) prove the inequality, and then (ii) find an equality to show that is is a weak inequality. \n",
    "\n",
    "a. For the 1-norm, denoting the $j$th column of $\\boldsymbol{A}$ by\n",
    "   $\\boldsymbol{a}_{j}$, and for a vector $\\|\\boldsymbol{x}\\|_{1} = 1$:\n",
    "\n",
    "   \\begin{align*}\n",
    "        \\| \\boldsymbol{A} \\boldsymbol{x}\\|_{1}\n",
    "        = \\|\\sum_{j=1}^{n} \\boldsymbol{a}_{j} x_{j}\\|_{1}\n",
    "        \\le\n",
    "        \\sum_{j=1}^{n} \\|\\boldsymbol{a}_{j}\\|_{1} |x_{j}|\n",
    "        \\le \\max_{j=1}^{n} \\|\\boldsymbol{a}_{j}\\|_{1}\n",
    "   \\end{align*}\n",
    "\n",
    "   This is the maximum column sum.\n",
    "\n",
    "   The term on the right could possibly be larger than $\\| \\boldsymbol{A} \\|_{1}$, whereas the norm is the smallest \n",
    "   possible value for the RHS that still satisfies the inequalities.  If we can show a case for which equality is\n",
    "   reached, $\\|\\boldsymbol{A} \\boldsymbol{x}\\|_{1} = \\max_{j=1}^{n} \\|\\boldsymbol{a}_{j}\\|_{1}$, we have the norm.\n",
    "   For the vector the $\\boldsymbol{e}_{j}$ with $e_{j} = 1$, where $j$ is the column with the greatest 1-norm,\n",
    "   and $e_{i \\ne j} = 0$, we have\n",
    "   \\begin{equation}\n",
    "        \\|\\boldsymbol{A} \\boldsymbol{e}_{j}\\|_{1} = \\max_{j=1}^{n} \\|\\boldsymbol{a}_{j}\\|_{1}\n",
    "   \\end{equation}\n",
    "   Therefore,\n",
    "   \\begin{equation}\n",
    "      \\|\\boldsymbol{A}\\|_{1} = \\max_{j=1}^{n} \\|\\boldsymbol{a}_{j}\\|_{1}\n",
    "   \\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "b. For a vector $\\|\\boldsymbol{x}\\|_{\\infty} = 1$:\n",
    "\n",
    "   \\begin{align*}\n",
    "        \\|\\boldsymbol{A} \\boldsymbol{x}\\|_{\\infty}\n",
    "          = \\max_{i=1}^{m} |\\sum_{j=1}^{n} a_{ij} x_{j}|\n",
    "          &\\le \\max_{i=1}^{m} \\sum_{j=1}^{n} |a_{ij}| |x_{j}|\n",
    "          \\\\\n",
    "          &\\le \\max_{i=1}^{m} \\sum_{j=1}^{n} |a_{ij}|\n",
    "   \\end{align*}\n",
    "\n",
    "   The RHS is the maximum row sum.\n",
    "\n",
    "   As before, we need to find a case with equality.  If the row\n",
    "   with the maximum sum is row $k$, we choose a vector $\\boldsymbol{x}$\n",
    "   where $x_{i} = \\pm 1$ such that the sign of $x_{i}$ is the\n",
    "   same as the sign of the entry $a_{ki}$. We then have\n",
    "\n",
    "   \\begin{equation}\n",
    "      \\|\\boldsymbol{A}\\|_{\\infty}\n",
    "      = \\max_{i=1}^{m} \\sum_{j=1}^{n} |a_{ij}|\n",
    "   \\end{equation}"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 10\n",
    "\n",
    "Recall that the $l_{2}$ norm is the square root of the maximum eigenvalue of $\\boldsymbol{A}^{T}\\boldsymbol{A}$\n",
    "(matrix is real in this question). We have\n",
    "$$\n",
    "    \\boldsymbol{A}^{T}\\boldsymbol{A}\n",
    "    =\n",
    "     \\begin{bmatrix}\n",
    "       26 & 5 \\\\ 5 & 1\n",
    "     \\end{bmatrix}\n",
    "$$\n",
    "For this matrix, $\\lambda = (27 \\pm \\sqrt{27^{2} -4})/2$. \n",
    "Hence $\\|\\boldsymbol{A}\\|_{2} \\approx 5.1926$.\n",
    "\n",
    "Checking for a vector of length $1$:\n",
    "$$\n",
    "    \\boldsymbol{A}\n",
    "    \\begin{bmatrix} \\cos \\theta \\\\ \\sin \\theta   \\end{bmatrix}\n",
    "    =    \n",
    "    \\begin{bmatrix} \\cos \\theta \\\\ 5\\cos \\theta + \\sin \\theta  \\end{bmatrix}\n",
    "$$\n",
    "Comptuing the norm,\n",
    "$$\n",
    "  F(\\theta) = \\|\\boldsymbol{A}\\boldsymbol{x}\\|_{2}^{2}\n",
    "    = 25 \\cos^{2} \\theta + 10 \\cos\\theta\\sin\\theta  + 1\n",
    "$$\n",
    "To find the extreme points of the function, we differentiate $F$ with respect to $\\theta$ \n",
    "and set the derivative equal to zero:\n",
    "\\begin{align}\n",
    "\\frac{d F}{d \\theta} &= -50 \\cos\\theta \\sin\\theta + 10 (\\cos^{2} \\theta - \\sin^{2} \\theta) \\\\\n",
    "&= -25 \\sin 2\\theta + 10 \\cos 2\\theta \\\\\n",
    "&= 0\n",
    "\\end{align}\n",
    "The norm is maximised when $\\theta = \\arctan(2/5) /2$. Evaluating:\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "θ = 0.19025318855618245\n",
      "F = 5.192582403567252\n"
     ]
    }
   ],
   "source": [
    "θ = np.arctan(2/5)/2\n",
    "print(\"θ =\", θ)\n",
    "\n",
    "F = 25*np.cos(θ)*np.cos(θ) + 10*np.cos(θ)*np.sin(θ) + 1\n",
    "print(\"F =\", np.sqrt(F))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 11 (condition number)\n",
    "\n",
    "Recall that $\\kappa(\\boldsymbol{A}) = \\| \\boldsymbol{A} \\| \\| \\boldsymbol{A}^{-1} \\|$.\n",
    "\n",
    "We have $\\| \\boldsymbol{A}\\|_{1} = 3$ and $\\| \\boldsymbol{A}\\|_{\\infty} \\approx 2$.\n",
    "Note that\n",
    "$$\n",
    "\\boldsymbol{A}^{-1} \\approx -\\frac{1}{2}\n",
    "\\begin{bmatrix}\n",
    "1 & -2 \\\\\n",
    "-1 & 10^{-4} \n",
    "\\end{bmatrix}\n",
    "$$ \n",
    "so $\\| \\boldsymbol{A}^{-1} \\|_{1} \\approx 1$ and $\\| \\boldsymbol{A}^{-1} \\|_{\\infty} = 3/2$.\n",
    "Therefore $\\kappa_{1}(\\boldsymbol{A}) \\approx 3$ and $\\kappa_{\\infty}(\\boldsymbol{A}) \\approx 3$.\n",
    "\n",
    "Recall that $\\kappa_{2}(\\boldsymbol{A}) = \\sqrt{\\lambda_{\\max}(\\boldsymbol{A}^{T}\\boldsymbol{A}})/\\sqrt{\\lambda_{\\min}(\\boldsymbol{A}^{T}\\boldsymbol{A}})$. As an approximation we ignore the $(1, 1)$ entry:\n",
    "$$\n",
    "  \\boldsymbol{A}^{T}\\boldsymbol{A} \\approx\n",
    "  \\begin{bmatrix}\n",
    "     1 & 1 \\\\ 1 & 5\n",
    "  \\end{bmatrix}\n",
    "$$\n",
    "Computing the eigenvalues, $\\lambda \\approx 3 \\pm \\sqrt{5}$, therefore $\\kappa_{2} \\approx \\sqrt{3 + \\sqrt{5}}/\\sqrt{3 - \\sqrt{5}} \\approx 2.618$. \n",
    "\n",
    "The matrix is very well conditioned, however LU factorisation may require pivoting. This is an issue with LU factorisation rather than a pathological problem with the matrix."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 12 (least squares)\n",
    "\n",
    "A projection matrix $\\boldsymbol{P}$  has the property $\\boldsymbol{P} = \\boldsymbol{P}\\boldsymbol{P}$, and $\\boldsymbol{P}^{H} = \\boldsymbol{P}$.\n",
    " \n",
    "a. The solution to the least squares problem is $\\hat{\\boldsymbol{x}} = \\boldsymbol{A}\n",
    "   (\\boldsymbol{A}^{H} \\boldsymbol{A})^{-1} \\boldsymbol{A}^{H} \\boldsymbol{b}$.\n",
    "   Therefore $\\boldsymbol{r} = \\boldsymbol{A} \\hat{\\boldsymbol{x}} - \\boldsymbol{b}\n",
    "   = \\boldsymbol{A}(\\boldsymbol{A}^{H} \\boldsymbol{A})^{-1} \\boldsymbol{A}^{H}\n",
    "   \\boldsymbol{b} - \\boldsymbol{b}$. Insert this expression for $\\boldsymbol{r}$ into the expression in the\n",
    "   question and re-arrange. to show the result.\n",
    "\n",
    "   Vectors $\\boldsymbol{A}\\boldsymbol{z}$ lie in the *column space* of $\\boldsymbol{A}$,\n",
    "   hence the expression says that the least-squares residual is *orthogonal* to the column\n",
    "   space of $\\boldsymbol{A}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "b.\n",
    "   $$\n",
    "     \\boldsymbol{P}\\boldsymbol{P}\n",
    "        = \\boldsymbol{A} (\\boldsymbol{A}^{H}\\boldsymbol{A})^{-1}\\boldsymbol{A}^{H}\n",
    "    \\boldsymbol{A}(\\boldsymbol{A}^{H}\\boldsymbol{A})^{-1}\\boldsymbol{A}^{H} \n",
    "     = \\boldsymbol{A} (\\boldsymbol{A}^{H}\\boldsymbol{A})^{-1}\\boldsymbol{A}^{H}\n",
    "    = \\boldsymbol{P}\n",
    "   $$\n",
    "   and\n",
    "   $$\n",
    "        \\boldsymbol{P}^{H}\n",
    "        = \\boldsymbol{A}(\\boldsymbol{A}^{H}\\boldsymbol{A})^{-H}\\boldsymbol{A}^{H}\n",
    "        = \\boldsymbol{A}(\\boldsymbol{A}^{H}\\boldsymbol{A})^{-1}\\boldsymbol{A}^{H}\n",
    "   $$\n",
    "   by $\\boldsymbol{A}^{H}\\boldsymbol{A}$ being Hermitian\n",
    "\n",
    "c. We can phrase a least squares problem as\n",
    "   $$\n",
    "        \\boldsymbol{A} \\hat{\\boldsymbol{x}}\n",
    "        = \\boldsymbol{A}(\\boldsymbol{A}^{H}\\boldsymbol{A})^{-1}\\boldsymbol{A}^{H} \n",
    "        \\boldsymbol{b}\n",
    "        = \\boldsymbol{P} \\boldsymbol{b}\n",
    "   $$\n",
    "   which says that $\\boldsymbol{P}$ projects $\\boldsymbol{b}$ into the column space of $\\boldsymbol{A}$. If\n",
    "   $\\boldsymbol{b}^{\\prime} = \\boldsymbol{P}\\boldsymbol{b}$ is in the column space of $\\boldsymbol{A}$, then\n",
    "   $\\boldsymbol{A} \\hat{\\boldsymbol{x}} = \\boldsymbol{b}^{\\prime}$ has a solution."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 13 (pseudo inverse)\n",
    "\n",
    "a.\n",
    "\n",
    "   - Firstly, if the $m \\times n$ matrix $\\boldsymbol{A}$ has linearly independent\n",
    "     rows, then the rank of $\\boldsymbol{A}$ is $m$ and the column space of\n",
    "     $\\boldsymbol{A}$ spans $\\mathbb{C}^{m}$, and the nullspace space\n",
    "     of $\\boldsymbol{A}^{H}$ contains the zero vector only.\n",
    "\n",
    "   - Now, consider the nullspace of $\\boldsymbol{A}\\boldsymbol{A}^{H}$:\n",
    "\n",
    "     \\begin{equation}\n",
    "      \\boldsymbol{A}\\boldsymbol{A}^{H} \\boldsymbol{x} = \\boldsymbol{0} \\ \\rightarrow\n",
    "      \\ \\boldsymbol{x}^{H}\\boldsymbol{A} \\boldsymbol{A}^{H} \\boldsymbol{x} = 0 \\ \\rightarrow\n",
    "      \\ (\\boldsymbol{A}^{H} \\boldsymbol{x})^{H}\\boldsymbol{A}^{H} \\boldsymbol{x} = 0\n",
    "     \\end{equation}\n",
    "\n",
    "     The above holds only if $\\boldsymbol{A}^{H} \\boldsymbol{x} = \\boldsymbol{0}$, which\n",
    "     says that $\\boldsymbol{x}$ must come from the nullspace of\n",
    "     $\\boldsymbol{A}^{H}$.  We have already determined that the nullspace of\n",
    "     $\\boldsymbol{A}^{H}$ contains only the zero vector, therefore\n",
    "     $\\boldsymbol{A}\\boldsymbol{A}^{H}$ is full rank (the nullspace of\n",
    "     $\\boldsymbol{A}\\boldsymbol{A}^{H}$ contains the zero vector only) and can be inverted.\n",
    "\n",
    "b. Since $\\boldsymbol{A}\\boldsymbol{A}^{H}$ is square and full rank, it can be inverted,\n",
    "   $$\n",
    "     \\boldsymbol{A} \\boldsymbol{A}^{+} = \\boldsymbol{A} \\boldsymbol{A}^{H}\n",
    "     (\\boldsymbol{A}\\boldsymbol{A}^{H})^{-1} = \\boldsymbol{I}.\n",
    "   $$\n",
    "   hence $\\boldsymbol{A}^{+}$ is a right-inverse of $\\boldsymbol{A}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 14 (stationary iterative methods)\n",
    "\n",
    "Define matrix $\\boldsymbol{A}$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 2 -1]\n",
      " [-1  2]]\n"
     ]
    }
   ],
   "source": [
    "A = np.array([[2, -1], [-1, 2]])\n",
    "print(A)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We split $\\boldsymbol{A}$ such that $\\boldsymbol{A} = \\boldsymbol{N} - \\boldsymbol{P}$. A method will converge if the largest absolute eigenvalue of $\\boldsymbol{N}^{-1}\\boldsymbol{P}$ is less the one.\n",
    "\n",
    "For the Richardson method $\\boldsymbol{N} = \\boldsymbol{I}$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 0. -2.]\n"
     ]
    }
   ],
   "source": [
    "# Richardson\n",
    "N = np.identity(2)\n",
    "P = N - A\n",
    "M = np.linalg.inv(N).dot(P) \n",
    "print(np.linalg.eigvals(M))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Largest eigenvalue (absolute value) is greater than 1, therefore method will not converge.\n",
    "\n",
    "For the Jacobi method, $\\boldsymbol{N} = \\text{diag}(\\boldsymbol{A})$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[ 0.5 -0.5]\n"
     ]
    }
   ],
   "source": [
    "# Jacobi\n",
    "N = np.diag(np.diag(A))\n",
    "P = N - A\n",
    "M = np.linalg.inv(N).dot(P) \n",
    "print(np.linalg.eigvals(M))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Largest eigenvalue (absolute value) is less than 1, therefore method will converge.\n",
    "\n",
    "For Gauss-Seidel, $\\boldsymbol{N}$ is the lower triangular part of $\\boldsymbol{A}$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "scrolled": true
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0.   0.25]\n"
     ]
    }
   ],
   "source": [
    "# Gauss-Seidel\n",
    "N = np.tril(A)\n",
    "P = N - A\n",
    "M = np.linalg.inv(N).dot(P) \n",
    "print(np.linalg.eigvals(M))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Gauss-Seidel will converge because largest eigenvalue is less than one."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 15 (SVD)\n",
    "\n",
    "\\begin{align}\n",
    "  \\boldsymbol{A}^{-1} &= \\left(\\boldsymbol{U} \\boldsymbol{\\Sigma} \\boldsymbol{V}^{H}\\right)^{-1} \\\\\n",
    "  &=\\boldsymbol{V}^{-H} \\boldsymbol{\\Sigma}^{-1} \\boldsymbol{U}^{-1} \\\\\n",
    "  &=\\boldsymbol{V} \\boldsymbol{\\Sigma}^{-1} \\boldsymbol{U}^{H}\n",
    "\\end{align}\n",
    "Non-singular matrix cannot have any zero singular values. In fact, smallest singular values is a measure of the 'distance' to a singular matrix."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 16 (SVD)\n",
    "\n",
    "Define matrix $\\boldsymbol{A}$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[1 3]\n",
      " [2 2]\n",
      " [3 1]]\n"
     ]
    }
   ],
   "source": [
    "A = np.array([[1, 3], [2, 2], [3, 1]])\n",
    "print(A)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute the reduced SVD (recall that NumPy uses $\\boldsymbol{A} = \\boldsymbol{U} \\boldsymbol{\\Sigma} \\boldsymbol{V}$ rather than  $\\boldsymbol{A} = \\boldsymbol{U} \\boldsymbol{\\Sigma} \\boldsymbol{V}^{T}$):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[4.89897949 2.        ]\n",
      "[[-5.77350269e-01  7.07106781e-01]\n",
      " [-5.77350269e-01  8.98662938e-17]\n",
      " [-5.77350269e-01 -7.07106781e-01]]\n",
      "[[-0.70710678 -0.70710678]\n",
      " [-0.70710678  0.70710678]]\n"
     ]
    }
   ],
   "source": [
    "U, s, V = np.linalg.svd(A, full_matrices=False)\n",
    "print(s)\n",
    "print(U)\n",
    "print(V)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute pseudoinverse by creating $\\boldsymbol{\\Sigma}^{+} =  \\boldsymbol{\\Sigma}_{1}^{-1}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-0.16666667  0.08333333  0.33333333]\n",
      " [ 0.33333333  0.08333333 -0.16666667]]\n"
     ]
    }
   ],
   "source": [
    "# Pseudoinverse\n",
    "Sigma_p = np.diag(1.0/s)\n",
    "Ap = (V.T).dot(Sigma_p.dot(U.T))\n",
    "print(Ap)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute $\\boldsymbol{A}^{+}\\boldsymbol{A}$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 1.00000000e+00 -3.88578059e-16]\n",
      " [ 2.49800181e-16  1.00000000e+00]]\n"
     ]
    }
   ],
   "source": [
    "# Check that A^{+}A = I\n",
    "print(Ap.dot(A))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "which is the identity. Compute now $\\boldsymbol{A}\\boldsymbol{A}^{+}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 0.83333333  0.33333333 -0.16666667]\n",
      " [ 0.33333333  0.33333333  0.33333333]\n",
      " [-0.16666667  0.33333333  0.83333333]]\n"
     ]
    }
   ],
   "source": [
    "# Check that AA^{+} \\ne I\n",
    "print(A.dot(Ap))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Which is clearly not the identity."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "(e) Recall that from $\\boldsymbol{A}^{T} \\boldsymbol{A} \\hat{\\boldsymbol{x}} = \\boldsymbol{A}^{T} \\boldsymbol{b}$ we have \n",
    "\\begin{align}\n",
    "\\hat{\\boldsymbol{x}} &= (\\boldsymbol{A}^{T} \\boldsymbol{A})^{-1}\\boldsymbol{A}^{T} \\boldsymbol{b} \\\\\n",
    "&= \\boldsymbol{A}^{+} \\boldsymbol{b}\n",
    "\\end{align}\n",
    "Multiplying both sides by $\\boldsymbol{A}$, \n",
    "$$\n",
    "\\boldsymbol{A}\\hat{\\boldsymbol{x}} = \\underbrace{\\boldsymbol{A} \\boldsymbol{A}^{+}}_{\\boldsymbol{P}} \\boldsymbol{b}\n",
    "$$\n",
    "where $\\boldsymbol{P}$ is the projection matrix from an earlier questions. Recall that $\\boldsymbol{P}$ projects a vector into the column space of $\\boldsymbol{A}$.\n",
    "Since $\\boldsymbol{P}$ is a projector, it does nothing if $\\boldsymbol{b}$ is already in column space. Therefore any vector $\\boldsymbol{b}$ in column space of $\\boldsymbol{A}$ is a solution."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 17 (pseudo inverse)\n",
    "\n",
    "Define matrix $\\boldsymbol{A}$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0 3 0]\n",
      " [2 0 0]]\n"
     ]
    }
   ],
   "source": [
    "A = np.array([[0, 3, 0], [2, 0, 0]])\n",
    "print(A)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute SVD:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[1. 0.]\n",
      " [0. 1.]]\n",
      "[3. 2.]\n",
      "[[0. 1. 0.]\n",
      " [1. 0. 0.]]\n"
     ]
    }
   ],
   "source": [
    "U, s, V = np.linalg.svd(A, full_matrices=False)\n",
    "print(U)\n",
    "print(s)\n",
    "print(V)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute $\\boldsymbol{\\Sigma}^{+} =  \\boldsymbol{\\Sigma}_{1}^{-1}$, and then $\\boldsymbol{A}^{+}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0.         0.5       ]\n",
      " [0.33333333 0.        ]\n",
      " [0.         0.        ]]\n"
     ]
    }
   ],
   "source": [
    "# \\Sigma^{+}\n",
    "Sigma_p = np.diag(1.0/s)\n",
    "\n",
    "# A^{+}\n",
    "Ap = (V.T).dot(Sigma_p).dot(U)\n",
    "print(Ap)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute $\\boldsymbol{A}^{+}\\boldsymbol{A}$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[1. 0. 0.]\n",
      " [0. 1. 0.]\n",
      " [0. 0. 0.]]\n"
     ]
    }
   ],
   "source": [
    "ApA = Ap.dot(A)\n",
    "print(ApA)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For this matrix, any $\\boldsymbol{x} = [x_{1} \\ \\ x_{2} \\ \\ 0]$ (which is from the row space of $\\boldsymbol{A}$) satisfies $\\boldsymbol{A}^{+} \\boldsymbol{A} \\boldsymbol{x} = \\boldsymbol{x}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Question 18 (rank deficient least squares)\n",
    "\n",
    "Define matrix and RHS vector"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[1 0 0]\n",
      " [1 0 0]\n",
      " [1 1 1]]\n",
      "[0 2 2]\n"
     ]
    }
   ],
   "source": [
    "A = np.array([[1, 0, 0], [1, 0, 0], [1, 1, 1]])\n",
    "print(A)\n",
    "\n",
    "b = np.array([0, 2, 2])\n",
    "print(b)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compute the SVD and print singular values:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[2. 1. 0.]\n"
     ]
    }
   ],
   "source": [
    "# Compute SVD\n",
    "U, s, V = np.linalg.svd(A)\n",
    "print(s)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There is one zero singular value, so we need to 'trim' the last column from $\\boldsymbol{U}$ and the last row from $\\boldsymbol{U}$, and compute $\\boldsymbol{\\sigma}^{+}$ "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[-0.40824829 -0.57735027]\n",
      " [-0.40824829 -0.57735027]\n",
      " [-0.81649658  0.57735027]]\n",
      "[[-0.81649658 -0.40824829 -0.40824829]\n",
      " [-0.57735027  0.57735027  0.57735027]]\n",
      "[[0.5 0. ]\n",
      " [0.  1. ]]\n"
     ]
    }
   ],
   "source": [
    "# Create view of U with last columns \n",
    "U1 = U[:, :2]\n",
    "print(U1)\n",
    "V1 = V[:2,:]\n",
    "print(V1)\n",
    "\n",
    "# Create Sigma^{+}\n",
    "S1 = np.diag(1.0/s[:-1])\n",
    "print(S1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Solve problem $\\boldsymbol{x} = \\boldsymbol{V}_{1} \\boldsymbol{\\Sigma}_{1}^{-1} \\boldsymbol{U}_{1}^{T}$ "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1.  0.5 0.5]\n"
     ]
    }
   ],
   "source": [
    "x = np.transpose(V1).dot(S1.dot(U1.T).dot(b))\n",
    "print(x)"
   ]
  }
 ],
 "metadata": {
  "celltoolbar": "Raw Cell Format",
  "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.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}