extra_2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Tutoraggio 2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "addpath(\"./functions\");\n",
    "pkg load symbolic"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Esercizio 1\n",
    "\n",
    "Si consideri la matrice:\n",
    "\n",
    "$$ A =\n",
    " \\begin{pmatrix}\n",
    "  4 & 0 & 1 & 0 \\\\\n",
    "  0 & \\alpha - 1 & 0 & 0 \\\\\n",
    "  1 & 0 & 1 & 0 \\\\\n",
    "  0 & 0 & 0 & 2\n",
    " \\end{pmatrix}\n",
    "$$\n",
    "\n",
    "in cui $\\alpha$ è parametro reale."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Punto a\n",
    "\n",
    "Scrivere una funzione che, preso in input $\\alpha$, calcoli il numero di condizionamento in norma infinito di A, utilizzando per il calcolo di $A^{−1}$ la fattorizzazione LU di A. \n",
    "\n",
    "Testare la function con $\\alpha = 1.2 : 0.05 : 3.2$ e produrre un grafico di $K_{inf}(A)$ in funzione di $\\alpha$, confrontando i valori ottenuti con l’output del comando built-in `cond` per il calcolo del numero di condizionamento in norma infinito."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Soluzione\n",
    "\n",
    "Per il calcolo dell'inversa della matrice $A$ di ordine $n$ si risolvono $n$ sistemi lineari sfruttando la fattorizzazione LU, senza strategie di pivoting.\n",
    "\n",
    "Per il calcolo dell'indice di condizionamento si impiega:\n",
    "\n",
    "$$ K(A) := \\|A\\| \\|A^{-1}\\| $$\n",
    "\n",
    "In base a quanto visto nello studio del condizionamento dei sistemi lineari."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "function [K, K_built] = cond_inf(alpha)\n",
    "\n",
    "A = [4 0 1 0; 0 alpha-1 0 0; 1 0 1 0; 0 0 0 2];\n",
    "[L, U, err] = gauss_simple(A);\n",
    "\n",
    "% Per la risoluzione si sfrutta l'equazione matriciale A * A_inv = I,\n",
    "% in particolare la forma A * A_inv(:, j) = I(:,j)\n",
    "\n",
    "A_inv = zeros(4,4);\n",
    "for j = 1:4\n",
    "    % aggiornamento del termine noto per la j-esima colonna\n",
    "    b = zeros(4,1);\n",
    "    b(j) = 1;\n",
    "\n",
    "    % in base al valore del termine noto si ottengono i valori\n",
    "    % della j-esima colonna di A_inv.\n",
    "    [x, err] = lusolve(L, U, [], [], b);\n",
    "    A_inv(:,j) = x;\n",
    "end\n",
    "\n",
    "% calcolo indice di condizionamento\n",
    "K = norm(A, inf) * norm(A_inv, inf);\n",
    "K_built = cond(A, inf);\n",
    "\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAJMmlDQ1BkZWZhdWx0X3JnYi5pY2MAAEiJlZVnUJNZF8fv8zzphUASQodQQ5EqJYCUEFoo0quoQOidUEVsiLgCK4qINEWQRQEXXJUia0UUC4uCAhZ0gywCyrpxFVFBWXDfGZ33HT+8/5l7z2/+c+bec8/5cAEgiINlwct7YlK6wNvJjhkYFMwE3yiMn5bC8fR0A9/VuxEArcR7ut/P+a4IEZFp/OW4uLxy+SmCdACg7GXWzEpPWeGjy0wPj//CZ1dYsFzgMt9Y4eh/eexLzr8s+pLj681dfhUKABwp+hsO/4b/c++KVDiC9NioyGymT3JUelaYIJKZttIJHpfL9BQkR8UmRH5T8P+V/B2lR2anr0RucsomQWx0TDrzfw41MjA0BF9n8cbrS48hRv9/z2dFX73kegDYcwAg+7564ZUAdO4CQPrRV09tua+UfAA67vAzBJn/eqiVDQ0IgALoQAYoAlWgCXSBETADlsAWOAAX4AF8QRDYAPggBiQCAcgCuWAHKABFYB84CKpALWgATaAVnAad4Dy4Aq6D2+AuGAaPgRBMgpdABN6BBQiCsBAZokEykBKkDulARhAbsoYcIDfIGwqCQqFoKAnKgHKhnVARVApVQXVQE/QLdA66At2EBqGH0Dg0A/0NfYQRmATTYQVYA9aH2TAHdoV94fVwNJwK58D58F64Aq6HT8Id8BX4NjwMC+GX8BwCECLCQJQRXYSNcBEPJBiJQgTIVqQQKUfqkVakG+lD7iFCZBb5gMKgaCgmShdliXJG+aH4qFTUVlQxqgp1AtWB6kXdQ42jRKjPaDJaHq2DtkDz0IHoaHQWugBdjm5Et6OvoYfRk+h3GAyGgWFhzDDOmCBMHGYzphhzGNOGuYwZxExg5rBYrAxWB2uF9cCGYdOxBdhK7EnsJewQdhL7HkfEKeGMcI64YFwSLg9XjmvGXcQN4aZwC3hxvDreAu+Bj8BvwpfgG/Dd+Dv4SfwCQYLAIlgRfAlxhB2ECkIr4RphjPCGSCSqEM2JXsRY4nZiBfEU8QZxnPiBRCVpk7ikEFIGaS/pOOky6SHpDZlM1iDbkoPJ6eS95CbyVfJT8nsxmpieGE8sQmybWLVYh9iQ2CsKnqJO4VA2UHIo5ZQzlDuUWXG8uIY4VzxMfKt4tfg58VHxOQmahKGEh0SiRLFEs8RNiWkqlqpBdaBGUPOpx6hXqRM0hKZK49L4tJ20Bto12iQdQ2fRefQ4ehH9Z/oAXSRJlTSW9JfMlqyWvCApZCAMDQaPkcAoYZxmjDA+SilIcaQipfZItUoNSc1Ly0nbSkdKF0q3SQ9Lf5RhyjjIxMvsl+mUeSKLktWW9ZLNkj0ie012Vo4uZynHlyuUOy33SB6W15b3lt8sf0y+X35OQVHBSSFFoVLhqsKsIkPRVjFOsUzxouKMEk3JWilWqUzpktILpiSTw0xgVjB7mSJleWVn5QzlOuUB5QUVloqfSp5Km8oTVYIqWzVKtUy1R1WkpqTmrpar1qL2SB2vzlaPUT+k3qc+r8HSCNDYrdGpMc2SZvFYOawW1pgmWdNGM1WzXvO+FkaLrRWvdVjrrjasbaIdo12tfUcH1jHVidU5rDO4Cr3KfFXSqvpVo7okXY5upm6L7rgeQ89NL0+vU++Vvpp+sP5+/T79zwYmBgkGDQaPDamGLoZ5ht2GfxtpG/GNqo3uryavdly9bXXX6tfGOsaRxkeMH5jQTNxNdpv0mHwyNTMVmLaazpipmYWa1ZiNsulsT3Yx+4Y52tzOfJv5efMPFqYW6RanLf6y1LWMt2y2nF7DWhO5pmHNhJWKVZhVnZXQmmkdan3UWmijbBNmU2/zzFbVNsK20XaKo8WJ45zkvLIzsBPYtdvNcy24W7iX7RF7J/tC+wEHqoOfQ5XDU0cVx2jHFkeRk4nTZqfLzmhnV+f9zqM8BR6f18QTuZi5bHHpdSW5+rhWuT5z03YTuHW7w+4u7gfcx9aqr01a2+kBPHgeBzyeeLI8Uz1/9cJ4eXpVez33NvTO9e7zofls9Gn2eedr51vi+9hP0y/Dr8ef4h/i3+Q/H2AfUBogDNQP3BJ4O0g2KDaoKxgb7B/cGDy3zmHdwXWTISYhBSEj61nrs9ff3CC7IWHDhY2UjWEbz4SiQwNCm0MXwzzC6sPmwnnhNeEiPpd/iP8ywjaiLGIm0iqyNHIqyiqqNGo62ir6QPRMjE1MecxsLDe2KvZ1nHNcbdx8vEf88filhICEtkRcYmjiuSRqUnxSb7JicnbyYIpOSkGKMNUi9WCqSOAqaEyD0tandaXTlz/F/gzNjF0Z45nWmdWZ77P8s85kS2QnZfdv0t60Z9NUjmPOT5tRm/mbe3KVc3fkjm/hbKnbCm0N39qzTXVb/rbJ7U7bT+wg7Ijf8VueQV5p3tudATu78xXyt+dP7HLa1VIgViAoGN1tubv2B9QPsT8M7Fm9p3LP58KIwltFBkXlRYvF/OJbPxr+WPHj0t6ovQMlpiVH9mH2Je0b2W+z/0SpRGlO6cQB9wMdZcyywrK3BzcevFluXF57iHAo45Cwwq2iq1Ktcl/lYlVM1XC1XXVbjXzNnpr5wxGHh47YHmmtVagtqv14NPbogzqnuo56jfryY5hjmceeN/g39P3E/qmpUbaxqPHT8aTjwhPeJ3qbzJqamuWbS1rgloyWmZMhJ+/+bP9zV6tua10bo63oFDiVcerFL6G/jJx2Pd1zhn2m9az62Zp2WnthB9SxqUPUGdMp7ArqGjzncq6n27K7/Ve9X4+fVz5ffUHyQslFwsX8i0uXci7NXU65PHsl+spEz8aex1cDr97v9eoduOZ67cZ1x+tX+zh9l25Y3Th/0+LmuVvsW523TW939Jv0t/9m8lv7gOlAxx2zO113ze92D64ZvDhkM3Tlnv296/d5928Prx0eHPEbeTAaMip8EPFg+mHCw9ePMh8tPN4+hh4rfCL+pPyp/NP637V+bxOaCi+M24/3P/N59niCP/Hyj7Q/Fifzn5Ofl08pTTVNG02fn3Gcufti3YvJlykvF2YL/pT4s+aV5quzf9n+1S8KFE2+Frxe+rv4jcyb42+N3/bMec49fZf4bmG+8L3M+xMf2B/6PgZ8nFrIWsQuVnzS+tT92fXz2FLi0tI/QiyQvpTNDAsAAAAJcEhZcwAACxMAAAsTAQCanBgAAAAfdEVYdFNvZnR3YXJlAEdQTCBHaG9zdHNjcmlwdCA5LjUzLjNvnKwnAAAalUlEQVR4nO3dMWgc2eE/8Df//4+UknWQFIn3Ch05FSs4gs+F5SbcWaAtQoiI7SMEYkGkK5XCRk6hEGJCzpGLuF2nMCSFJScmIYUMFiHNrYtYhINVoQRv4TUp8gPr1mWa/RVjy/JKu5a8s7tPO58Pwszum515926tr9+bN2+SZrMZAGDQ/t+gKwAAIQgkACIhkACIgkACIAoZBFKlUmk0Grsvq9Xq3pcAcBhdBVKj0Zidnf38889nZ2fX19dDCJcuXVpdXZ2bm6tUKhnVEIBc+J9uPry2tlYqlebn5z/55JPl5eUQQqFQuHbtWr1eX15enpqayqiSAAy/JJP7kNbX1yuVysjIyLvvvnvx4sUQwsTExPb2dvdHBiAnuuohpVZWVu7fv7+wsPDkyZNCoZC+eebMme6PnJWf/uFb6cavvv+PwdYEgHYymNRw5cqVe/fu3bhxI4RQr9fTNx8+fHjgzhMTE8ke3Z8dgOHQVSAtLy+ncxlGR0dDCJOTk0+ePAkh1Ov1YrHY7lPNPbo5++GlHSPdI4CYdXUNqVqtLi4uzszMbG1tnT17dn5+fnZ29syZMw8fPpyfny+VSvs/4toSAAfqdlJDo9HY2toqFAq7V48qlcrely0EEhAhVxC6kdVwVzaz7A5PIAERSpJ+/zIcGhk2naWDAIiCQAIgCgIJgCgIJACiIJAAYler1TY2Ng75slwu97Vy2RFIAEeXJK9+eq9Wqz148OCQL9fW1vpQpV7IYC07gHxJkrB3onPLyzfZ2dlZWloaGxsLIVy9enVpaSmEMDY2dv369c3Nzc3NzcePH9dqtbSjs7CwMD4+vrOzk+4/3PSQALrTbB6pn1Qulz/88MPr16+HED777LP33nsvzZ5yubyzs/Po0aPr169PT0+Xy+WlpaXp6enr16/nIY2CQAI4wN4RuZZxuXS7XWnnz4YQQnjw4MHCwkIIIe0SnT9/PoQwPT39+PHjEEKaPePj48+ePavVaume09PTPfyPjUYuAulnH/wm/Rl0RYBjotls/dlb1LLD7ptv/GwIIYTx8fHNzc0QQrlc/sY3vpFu7+zs7K/F2NhYrVYLIaR/trO72xv3jJxrSAB99emnny4tLZ06dWpzc7NcLp8/f/7Ro0d3797d2NhoiZOrV68uLCyke546dWpv0d27d9MkO3XqVDqmt7CwMDY2trGxcXxn2eViLbu0b/SLL37S5/MCx8WRF2TbOxB39N+iOzs7m5ub586dS19ubGzsbu/fs1artaTRgWq1Wq1Wa3ec3slwLbtcBBJAZxZXfWsWVwVg2AgkAKIgkACIgkACIAoCCYAoCCQAoiCQAIiCQAIgCgIJgCgIJACiIJAAiEIuVvt+9ufxdOOd7x7jhdkBhpseEkDsarXaxsbG7sv9D6rYKy09cJ/OHxw4gQQQu1qt9uDBg92XDx486JAr4+PjY2Nj6T4bGxtLS0stRb2taxdyMWRnpA6Ix87OztLSUhoMV69eTQNjbGwsfaL55ubm48ePa7Va+py9hYWF8fHxnZ2dliBZW1tbW1vb2dkpl8t37949derUqVOnyuVyy9P81tbWNjY2Lly4kL6zsbGRbuw9SzwRpYcE0CpJWn+yKg0hlMvlDz/88Pr16yGEzz777L333kuzp1wu7+zsPHr06Pr169PT0+VyeWlpaXp6On0gbMtBxsbG0uN89tlnjx8/Tp+Anm7svgwhXLhw4fz587v5lBa1nCWbJsuCQAJo1Wy2/mRVGkJ48ODBwsJCCCHtEp0/fz6EMD09/fjx4xBCmj3j4+PPnj2r1WrpntPT0y0HSd9ZWlpKH2TeWTpwt/cq1N6zHKZB+kMgAfTV+Ph4miLlcvkb3/hGur3bp9lrbGwsvVa0/4pR+s7m5ub4+Pgbz3ju3Lnr16/3/+nmR5WLa0gA8fj000+XlpbSiz3lcvn8+fOPHj26e/fu/ilwV69eXVhYaLkslEpH6jY2NtKBvqWlpXPnzm1ubrb0pcbGxjY3N2u12mFya+D6/Rj5iYmJ7e3tfp4R4I2SpK+/DHd2djY3N3e7LBsbG+26Lzs7O7VarSWNdj916tSpdPCtVqvt7OwcuFvai+rdzIUMm04gAfQ7kIZJhk3nGhIAURBIAERBIAEQBYEEQBRM+wYIIYRk/5oK9FcuAumf57+ebrx/99+DrQkQJ1PsYmDIDoAoCCQAopCLITsjdQDx00MCIAoCCYAoCCQAoiCQAIiCQAIgCgIJgCgIJACiIJAAiIJAAiAKAgmAKAgkAKKQi7Xskst/TTeaNz4abE0AaEcPCYAoCCQAopCLITsjdQDx00MCIAoCCYAoCCQAoiCQAIiCQAIgCgIJgCgIJACiIJAAiEIGgVStVuv1evfHASDPulqpodFozM3NFYvFer1eLBavXLly+vTpYrEYQkhfZlTJriXJi41mc6D1AKCtrgJpbW3tzJkzafB8/PHHn3zySbFYvH37djZVy0qSvMqhvdsAxKSrQJqZmUk3Go1GCKFer4+Oji4vL4+MjCwsLIyOjmZQwS4lSWg2f//VD9JXP2w2ZRJAnLq6hlQoFAqFQqVSmZubW1hYaDQaJ0+eLJVKJ06cWFxcbPepZI9uzg7AMOl2te+VlZWnT5/evHmzUCiEEEqlUghhampqfX293UeaOigA7NNVIK2urqZplL68detWsVicmprKomJZ+uH/fvFioM54HUCsugqkdML3pUuX0peXL19eXFycmZnZ2tq6ePFiBrXrXnrRKIQQXD0CiFqS+QBapVJJry0dWDoxMbG9vZ3tGQ9DGAFELvsnxkY4ZBde9pRkEkC0crR0kDQCiFmOAgmAmAkkAKIgkACIQo4CyboQADHLUSABEDOBBEAUchRIrxZtACA+2d8YG6Hv/O3Zy813BlkPANrLUQ8JgJjlK5As1gAQrVwM2f3l20bqAGKXrx4SANESSABEIXeBZOY3QJxyF0gAxEkgARAFgQRAFAQSAFEQSABEIXeBZLEGgDjlLpAAiJNAAiAKuQikn/7hWz/9w7cGXQsAOslFIKV2M8liDQARylEgARCzpNnfaWcTExPb29v9PON+SWKuHUB09JAAiIJAAiAKeQwk43UAEcpjIAEQIYEEQBQEEgBREEgARCGngWSxBoDY5DSQAIiNQAIgCjkNpGbTqB1AXPISSD/74DeDrgIAneQikNI0kkkAMctFIB3IAkIAUclFIP3ii5/s/glAnPL4PCQAIpSLHhIA8ct1IJn5DRCPXAcSAPEQSABEQSABEAWBBEAUch1I7o0FiEeuAwmAeAgkAKIgkACIgkACIAr/M+gK9MOzP4+nG+98t9ZSlCSmNgBEQQ8JgCjkKJD2d48AiEcuhuxEEUD8ctRDAiBmeQ8kMxoAIpH3QAIgEgIJgCgIJACiIJAAiIJACkky6BoAIJAAiEQGgVStVuv1+t6XjUaj+8MCkCtdrdTQaDTm5uaKxWK9Xi8Wi1euXLl06VKhUNja2rp8+fLU1FRWteyhJGmGkCTNZkhCcF8SwMB0FUhra2tnzpy5cuVKCOHjjz+enJwsFArXrl2r1+vLy8vHIJB21/pOXkaR1b8BBqSrQJqZmUk30jG6arU6OTkZQigUCg8fPuy+cln55/mvpxvv3/33gTvIIICB6+oaUqFQKBQKlUplbm5uYWEhfSctOnPmTLtPJXt0c/YsxVMTgLzqdrXvlZWVp0+f3rx5s1AorKys7M5u6NBDakbYH4mwSgA501Ugra6upmmUvpycnKxWqyGEdI5DBrXLSLuRutBsvnbRyAUkgMHpKpDSCd+XLl1KX96+ffvWrVsrKysPHz6cn5/PoHZ9kGbS7jYAA5JkPoBWqVTSa0sHlk5MTGxvb2d7xkzoHQEMVvZPjD0Gs70BiI+lgwCIgkACIAoC6YW9kxsA6D+B9IpJDQADJJAAiIJAAiAKAgmAKAgkAKKQ/Y2xEUou/zXdaN746A17Wq8BYED0kACIgkACIAq5GLJ740gdAAOnh/Qa6zUADIpAamVSA8BACCQAoiCQAIiCQAIgCgfMsltdXa1Wq/V6PYQwOjo6OTk5MzPT7pHkQ8ntsQD991ogra+v379//+TJk6VSqVgsjo6OViqVer1+48aNkydPLiwsjI6ODqqiAAy3pLmnL1CpVKamplr2aDQao6OjaYep+37SxMTE9vZ2lwfpNT0kgP57LZBarK+vVyqVtbW1DCNEIAFwoAOuIaU5dP/+/ZGRkU8++WRjY6P/1QIgb17rIa2srNy/f390dLRUKo2MjDx58uTKlSvZnu9Y9JCCThJA373WQ/rtb387MzMzPz8/OTlZqVSePHkyqGpl6/df/SDd+OH/fnHIj0gjgD57LZA2Njbu3LmzuLg4OjpaLBZHRkYGVS0A8ubgSQ3VanV9fX1tba1QKJRKpfn5+azON5Ahu7foIQHQZ51m2YWXExyuXbuW1fmOyzUkAPrsDYGUuWMUSOY1APTTa2vZLS8v37p1K70Hdq/V1dXFxcX97wNAVlp7SOnqQdVqdXdRhufPn5dKpQsXLmSybpAeEgAHajtkV61Wnz9/nq5ol+H5BBIAB2r7+Il05rdhOgD644Clg1L3799/8uRJpVKZmZnZfTPzhRtipnsE0E9tA2l+fr5Sqfz3v/89ceJEPysEQD51emLs1NTUV77ylS/36Fu1AMibtj2kVL1ev337dl9qAkCuvSGQSqXS4uLizMxMOtdu/+P7ACATbwikYrHY/VNijzWTvwH6o20gDdOE7+/87Vm68ZdvvzPYmgDQzhsCqV6v730qkiE7AHqkbSDtZk+j0ch2sQYA2K91cdWW4mq1Ojs728f69MRfvv1O+vOWn0+SkCSZ1giAVq33Ie3NpJWVlbm5ucuXL/e3SjFJkmYzhGYzNJsyCaCnXguk9EF8y8vL9Xp9dnb26dOnGxsbpVJpQHUbtJYJdjIJoJdae0hpJp07d+7ixYs3b9509SiEIIcA+uCApYOuXbt24cKFarXa/9pEyo1IAL332vOQJiYmDtwpwycYHaPnIYWwb9TOXbIAPfPatO/jFBX98fK6URKazSCNAHroDUsH8SKEEgN3AL3V6fETANA3AulQTPkG6DWBBEAUBNJhuYQE0FO5mNTw0z98K9341ff/MdiaANCOHhIAURBIAEThtZUa+uCYrdTwOgs1APSOHtIRmPwN0DsCCYAoCCQAoiCQjsaoHUCPCCQAoiCQjsxEO4BeEEgAREEgARCFDAKp0WjU6/XujwNAnmWwuOra2tqXX3555cqVEMLp06eLxWIIoVgspu8MJUs2AGSu20CanZ3d2tr68Y9/HEKo1+vFYvH27dsZ1Ctu6eRvmQSQoW4D6d69eysrK+l2vV4fHR1dXl4eGRlZWFgYHR3tunrZ+NkHv0k3fvHFTwZbEwDayXJSQ6PROHnyZKlUOnHixOLiYrvdkj0yPDsAx1qWD+grlUqlUimEMDU1tb6+3m63Pq8vDsCxkGUg3bp1q1gsTk1NZXjMTPRipE6qAmQry0A6c+bM4uLizMzM1tbWxYsXMzwyAEMv+wf0VSqVQqFQKBQOLD3WD+gDoHey7CGlIhyyAyB+lg7qinmCAFkRSABEQSB1LUl0lAC6J5C6kCTNZkhC88VSQmIJoAsC6W21LGbnviSA7ggkAKIgkLqlawSQCYH0ttLrRgBkJC+B9OzP49kfdHcuQ/qjrwTQhVwEUppGvcqk3R8AupCLQOoDo3cAXcp+LbsIvfPd2qCrAMAb6CFlwxQHgC4JJACiIJAyo5ME0A2BlCVT7QDemkACIAoCCYAoCCQAoiCQAIiCQMqeuXYAb0EgARAFgZQ9NyQBvAWBBEAUcrG46j/Pfz3deP/uvwdbEwDa0UPqCUs2ABxVLgIp7RjpHgHELBeBFKQRQPTyEkgARE4g9VCShJAk5oADHIZA6pkkaTZDEpov7ksSSwAdCaTeSJLXZtqZdQfwJgIJgCgIpN6yjBDAIQmk3tgTRIbrAA4jF0sHDUZL50guAXQkkHpJCAEcmiE7AKIgkPrE1AaAznIxZJdc/mu60bzx0aDqkF5RMoYH0I4eEgBRyFEgDbB79KIC7kkCaC8XQ3YDjyIA3ihHPaQYuIYE0I5AAiAKAgmAKAgkAKIgkAbAXDuA/QTSAJj/DbCfQAIgCrm4DylCL1YSCsmr1wD5JpAGJElCaL7KIevcAbknkAYhSUJT/gC8xjWkgUr2DNmZ5wDkWy56SL//6gfpxg//94vB1qSVITuAl/SQBkF/CGCfXPSQYrQnk5LQdEUJIBeBFN1IXeplCDWN2AEYsgMgEgIpCi4qAQikWBiyA3JOIAEQBYEEQBQEEgBREEjRMbsByCeBFB0z7oB8yiCQGo1GvV7ffVmtVhuNRveHBSBXMgiktbW1O3fupNuXLl1aXV2dm5urVCrdHzm3XnSSkuTFD0AOdBtIs7OzN27cSLfX19cLhcK1a9du3rxZLpe7rluOpSHUbL74kUlADnS7lt29e/dWVlbS7Wq1Ojk5GUIoFAoPHz7stmrZ+c7fnqUbf/n2O4OtyeG5TxbIm4wXVy0UCunGmTNn2u2T7Pn3vlWuD7B3pdV0O+0kaStgqGUcSLuzGzr0kITQG+yNn/RPaQTkQJaBNDk5Wa1WQwj1er1YLGZ45C4do5E6gNzKMpBKpdKtW7dWVlYePnw4Pz+f4ZFz5/WJDPpHQB4kmQ+gVSqVQqGwezGpxcTExPb2drZnzAODdsDQy/6JsVNTU5kfE4ChZ+mg48HNSMDQE0jHhiE7YLgJJACiIJAAiIJAAiAK2c+yow+SJDTDy0kOLi4BQ0EgHUNJEkKzdb07gGPOkN1xkyShZRa4KeHAUMhFD+mnf/hWuvGr7/9jsDXJRrMZghXAgWGjh3Tc7OkPvUgjuQQMBYF0PO2O0RmsA4ZF9ourdmZx1czsPuYcYCjk4hrScBJFwHAxZDckDN0Bx51AGhLmfgPHnUAaHjIJONZcQxoqr25OSjouLHSsSw/coZvSrCoWZ+mBO2jMdqUH7hBtcw0dgTRsmiEJyZ6vb8tdSi1z8/bfw7T3nZ6WhtBasUOWZv4yHIfmeuvG7PPLoLn61VzDyJDdMNr7xd0/kNehtOVL39PSvY5U2qUjNUjn0n4211s3Zk9priPJsLmGlB7S8EqSJKRf6D0rgw9h6au/ty//wnYufbXD/n+tDn/py19zb1N6wJHD6208dKVJEprNtyw94MihJQ8PWZqDrtELAml4Nff8Jn7V8R+60vBqUb+9vybalrbusHtspS//Gd58/Xb5fR9X+qog2be2QJLs/UdUa8N3U5oDhuyGXecRiSOVdjOk0LvSo8qwQTqXZjjo1KKfjam54myuIaWHNHRa/gq1fIm7KW3ZoXNptqc+jqXhiM3VuTSS/6ieNlfoOEqluXIgF4H0sw9+k2784oufDLYmfdL5u9tNaU8PPnylAzz18JUO8NRxlg4jQ3YAREEgARAFj58AIAp6SABEQSABEAWBBEAUBBIAURBIAERBIAEQBYEEQBQEEgBREEgAREEgARAFgQRAFHLx+Ilnfx5PN975bm2wNQGgHT0kAKIgkACIQi6G7IzUAcRPDwmAKAgkAKIgkACIgkACIAoCCYAoCCQAoiCQAIiCQAIgCgIJgCgIJACiIJAAiEIuAumf57/+z/NfH3QtAOgkF4GUkkkAMctRIAEQs1w8fuL9u/8edBUAeAM9JACiIJAAiIJAAiAKAgmAKAgkAKIgkACIgkACIAoCCYAoCKSoJUky6CocJ5rr8LTVkWiuI3nr5sp4pYbTp08Xi8UQQrFYvHLlSrYHB2CIZRlI9Xq9WCzevn07w2MCkBNZDtnV6/XR0dHl5eWVlZVGo5HhkbuXXP7roKsAQCdJs9nM6ljr6+vVavXs2bNbW1uff/75gV2liYmJrE4HQLS2t7eP+pEsA2mv2dnZe/fuHXC+1y929ejsrSd92T1q3vioD6fLUJL06n/QUNJch6etjkRzHclbN1eWQ3a3bt2qVCqd92m+LsOzdzrpjY/CMUyj0K/AHhqa6/C01ZForiN56+bKMvar1eri4uLMzMzW1lapVLp48WJWRwZg6GXfD61UKoVCoVAoZHtYAIabgVEAotDXlRoajUa9Xu/nGY8dTfTWqtWqputME721SqUS260sscmkif7/z3/+8ywqcyi/+93vHj16dPbs2b6d8dg5sIlOnz5dqVT+9Kc//etf/9J6+zUajR/84Af/+c9/NFE7HZrIt6uztOmazeYvf/nLr33ta9/85jcHXaPodGiiI3+7mv3yve997/333//1r3/dtzMeOwc20ZMnT370ox8NqEbHQ7lc3m20jz76aLCViVO7JvLteqNyuVwul5vaqr12TfQWLZbxWnYd3Lt3b2VlpW+nO44ObKLd9S9GRkYWFhZGR0cHUreYzczMpBsGVdpp10S+XW80Pz+fblSrVXO1DtSuid7i22W179g1Go2TJ0+WSqUTJ04sLi4OujoxSmd1ViqVubm5hYWFQVcnRu2ayLfrkFZWVm7cuDE5OTnoisRrfxO9xberr7Ps0n/+WwW8g85N1G79C1ZWVp4+fXr58mX/hm3njU3k29VZo9E4d+7c3//+90FXJF4dmuiQ3y49pNgdZv2LnFtdXX369OnNmzelUTvtmsi3642Wl5fX19dDCMYz22nXRG/x7erfNSSOKh1g+eMf/7i7/oXFLw6Uzma+dOlS+tIDUPbb30S+XYd08eLFxcXFarW6tbVlQPhA+5vorb9dbow9Hqx/Qe/4dnXWaDS2trY0UQcdmuhI3y6BBEAUXEMCIAoCCYAoCCQAoiCQAIiCQAIgCgIJgCgIJACiIJAAiIJAgl6pVqsrKyurq6shhPX19Wq1OugaQdSs1AA9ka4dd+/evefPn9+5c+f+/fv37t2zQCd0YHFV6Ik7d+4UCoXR0dHR0dHnz59PTk5KI+jMkB30SrFY3N2empoaYE3gWBBI0BPvvvvu1tZWCKHRaKRPhXENCTpzDQl65dKlS2fPnv3yyy/ffffdSqUyMjJy7dq1QVcK4iWQoIfq9Xr6JJhGo+EaEnQmkACIgmtIAETh/wCrN/IbJG7vWwAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "alpha = 1.2 : 0.05 : 3.2;\n",
    "len_alpha = length(alpha);\n",
    "K = zeros(len_alpha, 1);\n",
    "K_built = K;\n",
    "\n",
    "for i = 1:len_alpha\n",
    "    [x, y] = cond_inf(alpha(i));\n",
    "    K(i) = x;\n",
    "    K_built(i) = y;\n",
    "end\n",
    "\n",
    "plot( alpha, K, \"ro--\", ...\n",
    "      alpha, K_built, \"b--\", ...\n",
    "      1.2, 5:30);\n",
    "legend(\"cond LU\", \"cond built-in\")\n",
    "ylabel('K_{inf}(A)')\n",
    "xlabel('\\alpha');\n",
    "box off"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Per valori di $\\alpha$ vicini a $1.2$ l'indice di condizionamento cresce notevolmente, ciò accade perché per $\\alpha = 1$ la matrice $A$ è singolare. Ciò si verifica calcolando il determinante di $A$ in funzione di $\\alpha$ e trovando gli zeri del polinomio ottenuto."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Symbolic pkg v2.9.0: Python communication link active, SymPy v1.5.\n",
      "Determinante di A: \n",
      "ans = (sym) 6⋅α - 6\n",
      "Valori che annullano il determinante: \n",
      "ans = (sym) 1\n"
     ]
    }
   ],
   "source": [
    "syms alpha\n",
    "A = [[4 0 1 0]; 0 alpha-1 0 0; [1 0 1 0]; [0 0 0 2]];\n",
    "disp('Determinante di A: ')\n",
    "det(A)\n",
    "disp('Valori che annullano il determinante: ')\n",
    "solve(det(A))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Punto b\n",
    "\n",
    "Quando $\\alpha = 2$, determinare l’autovalore massimo e minimo di $A$, calcolando l’opportuno zero del polinomio caratteristico mediante il metodo di Newton."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Soluzione\n",
    "\n",
    "Il polinomio caratteristico è pari a: $p_{char} = \\det(A - xI)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "p_char = (sym)\n",
      "\n",
      "   4      3       2           \n",
      "  x  - 8⋅x  + 20⋅x  - 19⋅x + 6\n",
      "\n",
      "ans = (sym 4×1 matrix)\n",
      "\n",
      "  ⎡   1   ⎤\n",
      "  ⎢       ⎥\n",
      "  ⎢   2   ⎥\n",
      "  ⎢       ⎥\n",
      "  ⎢5   √13⎥\n",
      "  ⎢─ - ───⎥\n",
      "  ⎢2    2 ⎥\n",
      "  ⎢       ⎥\n",
      "  ⎢√13   5⎥\n",
      "  ⎢─── + ─⎥\n",
      "  ⎣ 2    2⎦\n",
      "\n",
      "x3 =  0.69722\n",
      "x4 =  4.3028\n"
     ]
    }
   ],
   "source": [
    "A = double(subs(A, alpha, 2));\n",
    "syms x\n",
    "p_char = det(A - x*eye(4))\n",
    "solve(det(A - x*eye(4)))\n",
    "x3 = (-sqrt(13) + 5) / 2\n",
    "x4 = (sqrt(13) + 5) / 2"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Dato che $\\frac{5 - \\sqrt(13)}{2} \\approx 0.7$ e $\\frac{5 + \\sqrt{13}}{2} \\approx 4.3$ si ha che l'autovalore massimo corrisponde alla soluzione più grande tra le quattro fornite, pari a $\\frac{5 + \\sqrt(13)}{2}$, mentre la soluzione più piccola è $\\frac{5 - \\sqrt{13}}{2}$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Visto che si vogliono calcolare gli zeri corrispondenti all'autovalore massimo e minimo impiegando il metodo di Newton, è necessario calcolare la derivata del polinomio caratteristico."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "der_p_char = (sym)\n",
      "\n",
      "     3       2            \n",
      "  4⋅x  - 24⋅x  + 40⋅x - 19\n",
      "\n"
     ]
    }
   ],
   "source": [
    "der_p_char = diff(p_char)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Bisogna considerare che l'algoritmo per il calcolo degli zeri di funzione che sfrutta il metodo di Newton funziona solo se l'approssimazione iniziale $x_0$ è scelta abbastanza vicino allo zero da trovare. In quel caso, per funzioni non convesse, la convergenza del metodo di Newton è garantita.\n",
    "\n",
    "Inoltre bisogna assicurarsi di scegliere il punto di innesco in modo da non convergere verso uno zero diverso da quello corrispondente all'autovalore massimo o minimo."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAJMmlDQ1BkZWZhdWx0X3JnYi5pY2MAAEiJlZVnUJNZF8fv8zzphUASQodQQ5EqJYCUEFoo0quoQOidUEVsiLgCK4qINEWQRQEXXJUia0UUC4uCAhZ0gywCyrpxFVFBWXDfGZ33HT+8/5l7z2/+c+bec8/5cAEgiINlwct7YlK6wNvJjhkYFMwE3yiMn5bC8fR0A9/VuxEArcR7ut/P+a4IEZFp/OW4uLxy+SmCdACg7GXWzEpPWeGjy0wPj//CZ1dYsFzgMt9Y4eh/eexLzr8s+pLj681dfhUKABwp+hsO/4b/c++KVDiC9NioyGymT3JUelaYIJKZttIJHpfL9BQkR8UmRH5T8P+V/B2lR2anr0RucsomQWx0TDrzfw41MjA0BF9n8cbrS48hRv9/z2dFX73kegDYcwAg+7564ZUAdO4CQPrRV09tua+UfAA67vAzBJn/eqiVDQ0IgALoQAYoAlWgCXSBETADlsAWOAAX4AF8QRDYAPggBiQCAcgCuWAHKABFYB84CKpALWgATaAVnAad4Dy4Aq6D2+AuGAaPgRBMgpdABN6BBQiCsBAZokEykBKkDulARhAbsoYcIDfIGwqCQqFoKAnKgHKhnVARVApVQXVQE/QLdA66At2EBqGH0Dg0A/0NfYQRmATTYQVYA9aH2TAHdoV94fVwNJwK58D58F64Aq6HT8Id8BX4NjwMC+GX8BwCECLCQJQRXYSNcBEPJBiJQgTIVqQQKUfqkVakG+lD7iFCZBb5gMKgaCgmShdliXJG+aH4qFTUVlQxqgp1AtWB6kXdQ42jRKjPaDJaHq2DtkDz0IHoaHQWugBdjm5Et6OvoYfRk+h3GAyGgWFhzDDOmCBMHGYzphhzGNOGuYwZxExg5rBYrAxWB2uF9cCGYdOxBdhK7EnsJewQdhL7HkfEKeGMcI64YFwSLg9XjmvGXcQN4aZwC3hxvDreAu+Bj8BvwpfgG/Dd+Dv4SfwCQYLAIlgRfAlxhB2ECkIr4RphjPCGSCSqEM2JXsRY4nZiBfEU8QZxnPiBRCVpk7ikEFIGaS/pOOky6SHpDZlM1iDbkoPJ6eS95CbyVfJT8nsxmpieGE8sQmybWLVYh9iQ2CsKnqJO4VA2UHIo5ZQzlDuUWXG8uIY4VzxMfKt4tfg58VHxOQmahKGEh0SiRLFEs8RNiWkqlqpBdaBGUPOpx6hXqRM0hKZK49L4tJ20Bto12iQdQ2fRefQ4ehH9Z/oAXSRJlTSW9JfMlqyWvCApZCAMDQaPkcAoYZxmjDA+SilIcaQipfZItUoNSc1Ly0nbSkdKF0q3SQ9Lf5RhyjjIxMvsl+mUeSKLktWW9ZLNkj0ie012Vo4uZynHlyuUOy33SB6W15b3lt8sf0y+X35OQVHBSSFFoVLhqsKsIkPRVjFOsUzxouKMEk3JWilWqUzpktILpiSTw0xgVjB7mSJleWVn5QzlOuUB5QUVloqfSp5Km8oTVYIqWzVKtUy1R1WkpqTmrpar1qL2SB2vzlaPUT+k3qc+r8HSCNDYrdGpMc2SZvFYOawW1pgmWdNGM1WzXvO+FkaLrRWvdVjrrjasbaIdo12tfUcH1jHVidU5rDO4Cr3KfFXSqvpVo7okXY5upm6L7rgeQ89NL0+vU++Vvpp+sP5+/T79zwYmBgkGDQaPDamGLoZ5ht2GfxtpG/GNqo3uryavdly9bXXX6tfGOsaRxkeMH5jQTNxNdpv0mHwyNTMVmLaazpipmYWa1ZiNsulsT3Yx+4Y52tzOfJv5efMPFqYW6RanLf6y1LWMt2y2nF7DWhO5pmHNhJWKVZhVnZXQmmkdan3UWmijbBNmU2/zzFbVNsK20XaKo8WJ45zkvLIzsBPYtdvNcy24W7iX7RF7J/tC+wEHqoOfQ5XDU0cVx2jHFkeRk4nTZqfLzmhnV+f9zqM8BR6f18QTuZi5bHHpdSW5+rhWuT5z03YTuHW7w+4u7gfcx9aqr01a2+kBPHgeBzyeeLI8Uz1/9cJ4eXpVez33NvTO9e7zofls9Gn2eedr51vi+9hP0y/Dr8ef4h/i3+Q/H2AfUBogDNQP3BJ4O0g2KDaoKxgb7B/cGDy3zmHdwXWTISYhBSEj61nrs9ff3CC7IWHDhY2UjWEbz4SiQwNCm0MXwzzC6sPmwnnhNeEiPpd/iP8ywjaiLGIm0iqyNHIqyiqqNGo62ir6QPRMjE1MecxsLDe2KvZ1nHNcbdx8vEf88filhICEtkRcYmjiuSRqUnxSb7JicnbyYIpOSkGKMNUi9WCqSOAqaEyD0tandaXTlz/F/gzNjF0Z45nWmdWZ77P8s85kS2QnZfdv0t60Z9NUjmPOT5tRm/mbe3KVc3fkjm/hbKnbCm0N39qzTXVb/rbJ7U7bT+wg7Ijf8VueQV5p3tudATu78xXyt+dP7HLa1VIgViAoGN1tubv2B9QPsT8M7Fm9p3LP58KIwltFBkXlRYvF/OJbPxr+WPHj0t6ovQMlpiVH9mH2Je0b2W+z/0SpRGlO6cQB9wMdZcyywrK3BzcevFluXF57iHAo45Cwwq2iq1Ktcl/lYlVM1XC1XXVbjXzNnpr5wxGHh47YHmmtVagtqv14NPbogzqnuo56jfryY5hjmceeN/g39P3E/qmpUbaxqPHT8aTjwhPeJ3qbzJqamuWbS1rgloyWmZMhJ+/+bP9zV6tua10bo63oFDiVcerFL6G/jJx2Pd1zhn2m9az62Zp2WnthB9SxqUPUGdMp7ArqGjzncq6n27K7/Ve9X4+fVz5ffUHyQslFwsX8i0uXci7NXU65PHsl+spEz8aex1cDr97v9eoduOZ67cZ1x+tX+zh9l25Y3Th/0+LmuVvsW523TW939Jv0t/9m8lv7gOlAxx2zO113ze92D64ZvDhkM3Tlnv296/d5928Prx0eHPEbeTAaMip8EPFg+mHCw9ePMh8tPN4+hh4rfCL+pPyp/NP637V+bxOaCi+M24/3P/N59niCP/Hyj7Q/Fifzn5Ofl08pTTVNG02fn3Gcufti3YvJlykvF2YL/pT4s+aV5quzf9n+1S8KFE2+Frxe+rv4jcyb42+N3/bMec49fZf4bmG+8L3M+xMf2B/6PgZ8nFrIWsQuVnzS+tT92fXz2FLi0tI/QiyQvpTNDAsAAAAJcEhZcwAACxMAAAsTAQCanBgAAAAfdEVYdFNvZnR3YXJlAEdQTCBHaG9zdHNjcmlwdCA5LjUzLjNvnKwnAAAbjUlEQVR4nO3dz2sb+d3A8Y+eZ89OnUPhaTMueGkElW5KD3JuGxssejOP7aU8BxnWDhRKCnXqXFQKpmAj7yGXQuRLSntYyyX05oDUva10qAU9jA7OYsHjWbrQw6rTf0DP4dvMo5VkWTOakb7fmffrECxZPyaSpbc/X42lVK/XEwAA5u0/5r0BAACIECQAgCYIEgBACwQJAKCFiYLkuq7jOP3H2Lbtuu5tBwEA8GuiIFWr1c8++8w7WCwWz87OdnZ2Go3G8EEAAAL44M5TbGxstNvtTz75RB28uLiwLOvw8NBxnFKp5Lpu/8GVlZWINxgAEE93T0hv3rzxaiQitm1ns1kRsSyr2WwOHIxuQwEA8RZkpwbLstQX+Xx++CAAAAEECZK3g4MaiQYODkun06k+gbYTABBzvoOUzWZvbm5ExHGcTCYzcPC2c/X6TLO5AIC4ununhgGFQuH09LRcLjebzd3d3YGDUWwiACAJUsFGlkajYVmW9+rRwMEB6XT66uoq+DYCMAor80kT1tJXwCD5QpCAREmlZvHEAk2EeHfz1kEAAC0QJACAFggSAEALBAkAoAWCBCC26vV6p9NR/05zCeFu1ZTXO69NmgGCBCC2arVap9NZXl5eXFwMdgnTnHeker1+cHCg/vV1vd5ZQt8kffj+w1gA0EelUrm+vu52uyJyfHwsIupZe3FxUR0UkXq9nsvlRKTVal1fX3c6nUqlMnDKVqtVrVZFpNvtek/3x8fH6rzLy8vDF9vtdg8ODtSJj4+Pu93u3t7e8vKyuvxOp6Mu8/z8vP/4arVar9e//PLLv/3tb1tbW7lc7uDgwNt+71zqWyLy6tUrEXn06NHl5WW9Xt/a2mq1Wt4mqWt/8eLF8OaZiAkJgMGur69FpFKpfPjhh0dHR+oL1Rv1rzpNt9vtdruXl5fHx8dra2uVSmXglF4Sut3uo0ePjo+P1bfUeUdebKVSUacUEbUwuLW1dXx8vLy83Gq1ut1urVar1+sDx29tbW1ubv7sZz/b3NzM5XKVSuX+/fve5XvnUteryqSiq86Yy+W8TfKu/ejoaHjzTMSEBCByqf3Pw7qo3slHA8esra2JyMHBwerqqrx/Rl5bW6vVagOnVPPE8vJyrVZrtVr9p1xeXlan8RbE/vnPf3pnrNVqwxeryiHvJ7PFxcWjo6NqtdrpdNQmqctcXV0dOL7f5eVlq9Wq1Wrdbldtv7clIrK1tXVwcPDq1autra2BM/Zf++rq6pj/tUEIEoDIDVckROoV/larpVqivlATzxi5XG7CU952YjXxqCknl8tVq1U1sgy8OHR0dDTyeEUtsq2urnY6neFdFVqt1vn5uYisrq7+4he/6P9W/7V///vfn/z/ojOCBMBsr169ur6+rtfrlUplcXFxc3Pz8vLy/Py8Xq+rF2BGevr0af8px++3NnBi78iDgwPVqr29vXq9Xq1W1cQj7+c2Ebl//37/8epFoJWVlVar1el01CWvrq6qiW1gMxYXF9XrT4uLiz/60Y9+//vfeyfov/ZKpTK8eUbqRe/hw4czuBYAmpjNE4vyq1/9qlarXV5efvPNN96RtVptwrNPfsqRJ/7mm2/6j7y8vLy+vlZf9J9s4Hi1tf3bPGYzrq+vvUsb+G8OXLuv/0uIQry7eXNVACGb5ZurHhwcrK2tqVdfMBch3t0ECUDIeLfvROHdvgEAcUOQAABaIEgAAC0QJACAFggSAEALBAkAoAWCBADQAkECAGiBIAEAtECQAABaIEgAAC0QJABxoD6z9bbv3vkBE9ABQQIQB51OZ8yHpXqfAztGvV4f+Rl6mBk+oA+AwbrdrvoIu263q5JzcHCgPjj1+Pi40+m0Wq1qtbq1taU+1PXFixfeh7QuLy+r83Y6nUqlUq1W6/X6+vr67373O+/IOzOGEPHxEwBCNvx5BO82vxfWhT88/3v/wb29vUePHu3t7anh5sMPP+x2uwcHB8fHxyKSy+VevXp1fn6uPjap1Wqpz2BdXV2t1+vq08E3NzfVd0WkVqs9evSo/0g+aelOIX78BBMSgMgNVCREao4RkbW1NfXRsa1Wq1ardbtd1ZLl5WXvxCozuVwul8uJyOrq6tHRUbVa7XQ63ieOjzwSs8FrSAAMtri4qPZWUP+q5Tg1/QznRC3EHR0dPX36VESOjo4ePXp0fn7ePwaNPBKzwYQEwGAvXrzY29vL5XKtViuXyz19+nRzc3N1dbXValUqleE967a2tmq1mhqb7t+/X61W1VClvqUupP9IsjRLvIYEIGQz/gjzbrfb6XTUKpxSr9cnDIl6VWl5eVmlqNVqqSmq/8jINjwmQry7CRKAkM04SJivEO9uXkMCAGiBIAEAtECQAABaIEgAAC2w2zeA8KVSqXlvAsxDkACEjF3sEAxLdgAALRAkAIAWCBIAQAsECQCgBYIEANACQQIAaIEgAQC0QJAAAFogSAAALRAkAIAWCBIAQAsECQCgBYIEANACQQIAaIEgAQC0EDBItm27rnvbQQAA/PIdJNd1NzY2zs7O1L8iUiwWz87OdnZ2Go1GBFsIAEgE358Y+/bt23w+//z5c9d1nz17trCwYFnW4eGh4zilUmllZSWKrQQAmOLd5vcenv89wBl9BymTyZydnTUajXa7nclkbNvOZrMiYllWs9kMsAUAAEiAJTvLshYWFi4uLi4uLpaWltQx6lv5fD7krQMAJIbvCalSqRQKhe3tbRF58uTJ+vq64zjqW2MmpFQq5X3d6/X8bycAwACB1+tk+t2+s9nszc2NiDiOk8lkbjtZr8+U1wgAiCXfE9LHH39cLBZvbm7a7fb6+nqhUDg9PS2Xy81mc3d3N4pNBAAkQSrYyNJoNCzL8l49Gjg4IJ1OX11dBd9GAIAJplmvkwATkjKwezd7ewMApsRbBwEAtECQAAAhmHK9TggSAEATBAkAoAWCBACY1vTrdUKQAACaIEgAAC0QJADAVEJZrxOCBADQBEECAGiBIAEAggtrvU4IEgBAEwQJAKAFggQACCjE9TohSAAATRAkAIAWCBIAIIhw1+uEIAEANEGQAABaIEgAAN9CX68TggQA0ARBAgBogSABAPyJYr1OCBIAQBMECQCgBYIEAPAhovU6IUgAAE0QJACAFggSAGBS0a3XCUECAGiCIAEAtECQAAATiXS9TggSAEATBAkAoAWCBAC4W9TrdUKQAACaIEgAAC0QJADAHWawXicECQCgCYIEANACQQIAjDOb9TohSAAATRAkAIAWCBIA4FYzW68TggQA0ARBAgBogSABAEab5XqdECQAgCYIEgBghBmPR0KQAACaIEgAAC0QJADAoNmv1wlBAgBogiABALRAkAAA3zKX9ToJHCTHcRzH8Q7atu26bkibBABIog8CnKdUKomI4ziFQmF7e7tYLFqW1W639/f3V1ZWwt5CAEAi+A5So9EQkcPDQ9d1f/3rXy8sLFiWdXh46DhOqVQiSABgtHmt10mAJbsvvvhiaWnp7Ozs7du3L1++tG07m82KiGVZzWYzgi0EACRCkNeQPvvsMxG5ubkpFosiYlmWOj6fz4e3YQCAWZvjeCTBXkNaX1/f3t4WkY2NDRHx9m4YMyGlUinv616vF+BKAQDx5ntCWlpa8r52XTebzd7c3IiI4ziZTOa2c/X6BNtQAEC8pQIUYmNjI5/Pt9tttZedOthsNnd3dwuFwvDp0+n01dVVGFsLAIjKfNfrJFiQRKTRaFiW5b16NHBwAEECAP3NPUhBXkMSkYHdu9nbGwAwJd46CAAw//FICBIAQBMECQCSTofxSAgSAEATBAkAoAWCBACJpsl6nRAkAIAmCBIAJJc+45EQJACAJggSAEALBAkAEkqr9TohSAAATRAkAIAWCBIAJJFu63VCkAAAmiBIAJA4Go5HQpAAAJogSAAALRAkAEgWPdfrhCABADRBkAAgQbQdj4QgAQA0QZAAAFogSACQFDqv1wlBAgBogiABQCJoPh4JQQIAaIIgAQC0QJAAIP70X68TggQA0ARBAoCYM2I8EoIEANAEQQKAODNlPBKCBADQBEECAGiBIAFAbBm0XicECQCgCYIEAPFk1ngkBAkAoAmCBAAxZNx4JAQJAKAJggQA0AJBAoC4MXG9TggSAEATBAkAYsXQ8UgIEgBAEwQJAOLD3PFICBIAQBMECQCgBYIEADFh9HqdECQAgCYIEgDEgenjkRAkAIAmCBIAGC8G45EQJACAJoIHyXEc13XV17Zte18DAGYpHuORBA6S67obGxv37t0TkWKxeHZ2trOz02g0Qt02AECCfBDsbCcnJwsLC67rNhoNy7IODw8dxymVSisrK+FuHwBgjNiMRxJsQjo9PV1aWrIs6969e7ZtZ7NZEbEsq9lshr15AICk8B0k27Zt297d3fWOsSxLfZHP50PbLgBAwvhesjs9PX3w4EG5XFZrdKlUynEc9a0xE1IqlfK+7vV6ATYUADAgTut1EiBIu7u7//rXv0Sk2WwWCoWvv/660+mIiOM4mUzmtnMRIQDAeL6DpF4xEpGFhQW1C8PGxka5XG42m/3reACASMVsPBKRVCizi9rXznsxaUA6nb66upr+WgAAnvgFKeBu3wPY2xsAZil+NRLeOggAoAmCBACGieV4JAQJAKAJggQAJonreCQECQCgCYIEAMaI8XgkBAkAoAmCBABmiPd4JAQJABCu1P7nwc5IkADAALEfj4QgAQBClNr/vHfyUbDzEiQA0F0SxiMhSAAATRAkANCaQePRNOt1QpAAAJogSACgr+SMR0KQAACaIEgAoCmDxqNQECQAwLSmX68TggQAekraeCQECQAwpVDGIyFIAKChBI5HQpAAAJogSACgF7PGo7DW64QgAQA0QZAAQCOJHY+EIAEANEGQAEAXSR6PhCABADRBkABAC2aNR1EgSAAwf8bVKPT1OiFIAABNECQAmDPGI4UgAQC0QJAAYJ4YjzwECQCgBYIEAHNj3HgUKYIEAJhUdOt1QpAAYF4YjwYQJACYAxNrFOl4JAQJAKAJggQAs8Z4NBJBAgBogSABwEwxHt2GIAEAtECQAGB2GI/GIEgAMCMm1miWCBIA4FYzG4+EIAHAbDAe3YkgAQC0QJAAIHKGjkezXK8TggQAUTO0RrNHkAAAI8x4PBKCBACRYjyaHEECgKiYW6PZj0cSOEi2bTuO03/Qdd2QNgkAkEQf+D2D67o7OzuZTMZxnEwm8/z582KxaFlWu93e399fWVmJYisBwDiMR375DlK1Ws3n88+fPxeRJ0+eZLNZy7IODw8dxymVSgQJAMTkGs2R7yCtr6+rL9QanW3b2WxWRCzLajab4W4cAGDG5jUeSYDXkCzLsiyr0Wjs7Ozs7e2pY9S38vl8yFsHAAYydzyaY40kwIQkIuVy+auvvnr58qVlWeVy2du7YcyElEqlvK97vV6AKwUAxJvvIJ2dnakaqYPZbNa2bRFR+zjcdi4iBCAhGI8C8x0ktcN3sVhUB1+/fn16eloul5vN5u7ubshbBwBGMbdGOkiFMrs0Gg312tLI76bT6aurq+mvBQA0Z26Q5j4eSbDXkIaxtzcAmFsjTfDWQQAQAqNrpMN4JAQJAKAJggQA02I8CgVBAoCpGF0jrRAkAAjO9BrpMx4JQQIAaIIgAUBAjEfhIkgAEAQ1Ch1BAgBogSABgG+MR1EgSADgj+k10hZBAgAfYlAjPccjIUgAkCja1kgIEgBMLgbjkc4IEgBMJAY10nk8EoIEAJOIQY30R5AAIBE0H4+EIAHAnWIwHulfIyFIADBeDGpkCoIEALeKR42MGI+EIAHAbajRjBEkAIAWCBIAjMB4NHsECQAGxaNGxiFIAPAtsamRWeORECQA6EeN5oggAUDcmFgjIUgA4InNeGQoggQAIjGqkaHjkRAkABBqpAeCBCDpYlMj0xEkAIkWpxoZPR4JQQKQZNRIKwQJQELFqUbxQJAAJFHMahSD8UgIEoAEokZ6IkgAkoUaaYsgAYCp4lQjIUgAEiVm41HMECQASRGzGsVsPBKCBCAhqJH+ZhSk1P7nqf3PZ3NdADCAGhnhg9lcjbrt4nojAtDWu83viQg1MsKMgqT0Tj5Sc1Jcb00AWonZYBR7Mw2S9I1KQpYARCmWNYrxeCSzD5JClgBEihqZaD5BUsgSgChQI0PNM0gKWQIQImpkrvkHSenPklAmAIFQI6PpEiTFu9EZmAD4RY1Mp1eQPAxMAHyhRjGgaZCUgYFJKBOAIfH701claTUSzYPkoUwARorlYCSJrJGYEiQPZQLgoUYxY1iQPMNlEuIEJAk1ip9wgmTbtmVZ9+7dC+XSfOm/52IfJ7VWPo1YPoCRNHF90UiSXSMJJUjFYtGyrHa7vb+/v7KyMv0FBnZbnCT0PqVSIiK9XpiXKSLql74/fS0i7/77v4a/O/0j8M6kGfkgT6WiuC+gp7gNRn1PJgmvkUwfpIuLC8uyDg8PHccplUrzDVK/gft1+NOYAt7x/Sma+nlwIA8P//T1w/cX/jCaJ9k7H8kji6Xv49+7OyL7FQFaiWGN3j+ZpH75l4TXSKYPkm3b2WxWRCzLajabYWxSJIbv6ds+MPDunwnvWU89D/p5Ehws0MBDq//S/F/4SL4/F/EHfxw+7mqoUulRJ4vCuLtj4OYSRqU4i+EyXd+Pa+qXf+l9+kROkv7TG8KSnWVZ6ot8Pj/9pc3Mbc90Y57Be58+Sf3yL9L/LhITZKM/QuMeTv0/neoqAuRkeJvD+Z1rcLOHExXRM8W4u4OdLRMjboORDD7eeycfyUk4v4MaLdWb7v9fLpeXlpa2t7dFJJ1OX11dDZ8mnU6/e/fOOzjlNc7Tt+drddy/KzXk6n//R30x4TDhXU7v0ycm/r5/x/AXhaG7Y+R9QauMFsPBSHn/05va/9zQh3wUpp2Qstmsbdsi4jhOJpO57WQGR6ifNw95/59Uqv//NnIYmvR/ftLrqQs3c2wfeMqYaZ+8u+DTUNdmMW8xHIw8vd7/v250Qo3+bdoJSUQ2Njby+Xyz2dzd3S0UCsMnuG1yMtKonRq8J99pHzyh7jGhlUnXLf2aYqeG0PZzQQRiOxi99+9leWajbwshSCLSaDQsy/JeTBoQqyAp75/+InnYxH2HsfDjFNKDeeQsRaVmL86DkYgM7N4d98e7L+EEabz4BSm0kSjx5vDKk0/R/kEbvi32g5Hwx0ZjESR/kvCAmaOoVvbCw0JfRBLyyKJG4xGkiTASzZ7+cVIYoaYX+zU64UNHJ0OQ7pCQX9w0Z0qchD75lJDHF4PRhAjSrRLyUDGOQXES+nS7hDy+GIx8IUgjJOShEgPGLaXG/g3pJ5GcxxeDkV8E6VuS81CJGbPGJiWBw1NyHl8MRsEQpH9LzkMl9kyMk8R9eErU44vBKDCCJJKMnXySiTjNXdJSJObfZXOU9CAl6tGScDGIk5jzZGfcy3tTIkWhSG6QSFGS6f8OEbfRf3hK4COLNbqwJDRIrNGhXwyGp7k/ISZtJFIYjMKVuCAl8Nc3+EKcfElmh4QURSNZQWIwgi+GruzN4GWnxHZISFGUkhIkBiNMLwZ9mvJplA4JKYpSIoLEYIQomNinAMOToWuYIWIkmpmYB4nBCDMz0Ccx4QfvtuGJCAkj0TzEOUgMRpgvUxI1vJ3pH/wxyc/CjETz8sG8NyAq1AhzN/wTOPzUP+bEURi5AcNX3TP2D3ID02of+sSK4YTEMh1MNKZVIQr8uIjrR+WyLqeVuAWJwQiYDUMTlbTJzyyxWrKjRsDMDD+VDydq5MlmjAIZJCZBYpkOmLuRz/UjKzXm9NPQs4iYXByCxGAEaGtMD8a0KvTrghGMDxI1AgxFPzDgP+a9AVOhRgAQGwYHiRoBQJyYGiRqBAAxY95rSOxQBwCxZFiQGIwAIK5MWrKjRgAQY8YEiRoBQLyZESRqBACxZ0CQqBEAJIHuQaJGAJAQWgeJGgFAcugbJGoEAImiaZCoEQAkjY5BokYAkEDaBYkaAUAy6RUkagQAiaVRkKgRACSZLkGiRgCQcLoECQCQcFoEifEIADD/IFEjAIDMPUjUCACgzDNI1AgA4JlbkKgRAKDffIJEjQAAA+YQJGoEABg2/73sAACQ2QeJ8QgAMNJMg0SNAAC3mV2QqBEAYIwZBYkaAQDGCxgk27Ydx+k/6LpuSJuUdKlUat6bYBhuMb+4xfziFvMl8M31gd8zuK67s7OTyWQcx8lkMs+fPy8Wi5Zltdvt/f39lZWVkediPAIAjOc7SNVqNZ/PP3/+XESePHmSzWYtyzo8PHQcp1Qq3RYkAADG8x2k9fV19YVao7NtO5vNiohlWc1mM9yNAwAkR6rX6wU4W6PRODk52d7evrm5efz4sRqMisXi69evh0+cTqen3EoAgFmurq78nmXSCanRaHzxxRdLS0vb29vlcvmrr756+fKlZVnlctnbu+G2Cendu3f9B4MlMDlSqYC/JSQWt5hf3GJ+cYv5EvjmmjRIKysragw6OztTNVLHZ7NZ27ZFRO3jMPK83JG+cHP5xS3mF7eYX9xivgS+uXy/hqR2+C4Wi+rg69evT09Py+Vys9nc3d0NthEAAIQzhzYaDcuyLMua/qIAAMnEwigAQAuRv3UQb+Lgl+u6/e+CgTsNvG8I7tRoNHhU+uI4DrfYDPznb37zm+guvVgsOo5zenrKgt7k/vCHP1xeXj5+/HjeG2IA13V/+tOf/uMf//jzn//85ZdfcqPdSd1ivV7vt7/97Xe/+90f/vCH894iA7iu+5Of/OTnP//5vDfEDD/+8Y8bjUawh6TvnRomd3FxwZs4+LWxsdFutz/55JN5b4gZBt43RH2BMarVaqFQ2N3d/fjjj0ulUqFQmPcWGeDk5GRhYcF13Xv37s17W3Sndrce+Qepk4gwSLyJQwBv3rwpl8vz3gpjDLxvCO7k7Qpr2zaLFpM4PT1dWlpyHIcaTULdUKVSaWFhYW9vz++NFu1rSN5PfD6fj/SKkExqKbjRaOzs7Ozt7c17c4xRLpdPTk7U74sYw7Zt27b5g5bJua774MGDQqHwne9859mzZ37PHuFeduVyWb2zg4ik0+kAbyORTGpCYvVpQup9Q/b39/l93xfXdVdXV//617/Oe0O09uzZswcPHojI27dvV1ZW9vb2+DGb3MbGxps3b3ydJcIJKZvN3tzcyNg3cQCm4b1vCE8TEyqVShcXFyLCAtQkdnd3Hz9+/Pjx43v37hUKhYWFhXlvke5OT08bjUbgs0f4GlKhUOBNHBCp4fcNmefWmGB7e/vZs2e2bbfbbRY57+Stai4sLLBb1iTy+fyzZ8/W19fb7bZaHvMl8j+M5U0cAK24rttut3lUIjqBn/Z5pwYAgBYif6cGAAAmQZAAAFogSAAALRAkAIAWCBIAQAv/BzmsBHk70/nDAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "p_char = @(x) x.^4 - 8 * x.^3 + 20 * x.^2 - 19 * x + 6;\n",
    "der_p_char = @(x) 4 * x.^3 - 24 * x.^2 + 40 * x - 19;\n",
    "\n",
    "figure\n",
    "t = linspace(0, 5);\n",
    "plot( t, p_char(t), ...\n",
    "      t, der_p_char(t), ...\n",
    "      1, p_char(1), \"ro\", ... % zeri di p_char\n",
    "      2, p_char(2), \"ro\", ...\n",
    "      x3, p_char(x3), \"ro\", ...\n",
    "      x4, p_char(x4), \"ro\");\n",
    "legend(\"polinomio caratteristico\", \"derivata\");\n",
    "box off"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Considerando la convergenza locale e la distribuzione degli zeri del polinomio caratteristico, si scelgono come punti di innesco per trovare gli zeri di funzione $x_3$ e $x_4$ rispettivamente $(0, p_{char}(0))$ e $(5, p_{char}(5))$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x3_newton =    6.972243622680053e-01\n",
      "x3 =    6.972243622680054e-01\n",
      "x4_newton =  4.302775637731996\n",
      "x4 =  4.302775637731995\n"
     ]
    }
   ],
   "source": [
    "format long\n",
    "[x3_newton, appr_x3, nit_x3] = newton(p_char, der_p_char, ...\n",
    "    0, 10^-12, 10^-12, 100);\n",
    "[x4_newton, appr_x4, nit_x4] = newton(p_char, der_p_char, ...\n",
    "    5, 10^-12, 10^-12, 100);\n",
    "x3_newton\n",
    "x3\n",
    "x4_newton\n",
    "x4"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "L'approssimazione è particolarmente buona e la velocità di convergenza è relativamente elevata."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Numero di iterazioni richieste per approssimare x3: 8\n",
      "Numero di iterazioni richieste per approssimare x4: 6\n"
     ]
    }
   ],
   "source": [
    "fprintf(\"Numero di iterazioni richieste per approssimare x3: %.i\\n\", ...\n",
    "    nit_x3);\n",
    "fprintf(\"Numero di iterazioni richieste per approssimare x4: %.i\\n\", ...\n",
    "    nit_x4);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Esercizio 2\n",
    "\n",
    "Si consideri la matrice $n \\times n$: $A = B^T B$\n",
    "\n",
    "dove `B = (1/n) * hilb(n) + diag(1:n)`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Punto a\n",
    "\n",
    "a) si calcola la fattorizzazione $LL^T$ di $A$ per valori di `n = 10 : 1 : 50` mediante l’algoritmo di Cholesky built-in;\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Punto b\n",
    "\n",
    "Per ogni `n = 10 : 1 : 50` si determina la soluzione del sistema lineare $Ax = b$ con termine noto `b = A ∗ (1 : n)'` sfruttando la fattorizzazione calcolata al punto a);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Punto c\n",
    "\n",
    "Per ogni `n = 10 : 1 : 50` si confronta la soluzione ottenuta al punto b) con quella restituita dal metodo di eliminazione gaussiana senza pivoting, mostrando in un grafico gli errori relativi generati dai due metodi al variare di $n$. Commentare i risultati ottenuti nel grafico."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Soluzione\n",
    "\n",
    "I tre punti, per come sono formulati, possono essere risolti simultaneamente."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = (10 : 1 : 50);\n",
    "N = length(n);\n",
    "err_rel_chol = zeros(N, 1);\n",
    "err_rel_gauss = zeros(N, 1);\n",
    "\n",
    "for k = 1:N\n",
    "    B = (1 / n(k)) * hilb(n(k)) + diag(1:n(k));\n",
    "    A = B' * B;\n",
    "    \n",
    "    x_esatta = (1:n(k))';\n",
    "    b = A * x_esatta;\n",
    "    \n",
    "    [R, p] = chol(A, 'lower');\n",
    "    %x_chol = lusolve(R, R', [], [], b);\n",
    "    y = lsolve(R, b);\n",
    "    x_chol = usolve(R', y);\n",
    "    \n",
    "    [L, U, err] = gauss_simple(A);\n",
    "    x_gauss = lusolve(L, U, [], [], b);\n",
    "    \n",
    "    err_rel_chol(k) = norm(x_esatta - x_chol) / ...\n",
    "        norm(x_esatta);\n",
    "    err_rel_gauss(k) = norm(x_esatta - x_gauss) / ...\n",
    "        norm(x_esatta);\n",
    "end"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Le matrici di Hilbert sono il classico esempio di matrici mal condizionate. Tuttavia gli si sta sommando una matrice diagonale, ben condizionata, con termini tutti interi e positivi. Ciò permette di ottenere matrici che non sono particolarmente mal condizionate."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAjAAAAGkCAIAAACgjIjwAAAJMmlDQ1BkZWZhdWx0X3JnYi5pY2MAAEiJlZVnUJNZF8fv8zzphUASQodQQ5EqJYCUEFoo0quoQOidUEVsiLgCK4qINEWQRQEXXJUia0UUC4uCAhZ0gywCyrpxFVFBWXDfGZ33HT+8/5l7z2/+c+bec8/5cAEgiINlwct7YlK6wNvJjhkYFMwE3yiMn5bC8fR0A9/VuxEArcR7ut/P+a4IEZFp/OW4uLxy+SmCdACg7GXWzEpPWeGjy0wPj//CZ1dYsFzgMt9Y4eh/eexLzr8s+pLj681dfhUKABwp+hsO/4b/c++KVDiC9NioyGymT3JUelaYIJKZttIJHpfL9BQkR8UmRH5T8P+V/B2lR2anr0RucsomQWx0TDrzfw41MjA0BF9n8cbrS48hRv9/z2dFX73kegDYcwAg+7564ZUAdO4CQPrRV09tua+UfAA67vAzBJn/eqiVDQ0IgALoQAYoAlWgCXSBETADlsAWOAAX4AF8QRDYAPggBiQCAcgCuWAHKABFYB84CKpALWgATaAVnAad4Dy4Aq6D2+AuGAaPgRBMgpdABN6BBQiCsBAZokEykBKkDulARhAbsoYcIDfIGwqCQqFoKAnKgHKhnVARVApVQXVQE/QLdA66At2EBqGH0Dg0A/0NfYQRmATTYQVYA9aH2TAHdoV94fVwNJwK58D58F64Aq6HT8Id8BX4NjwMC+GX8BwCECLCQJQRXYSNcBEPJBiJQgTIVqQQKUfqkVakG+lD7iFCZBb5gMKgaCgmShdliXJG+aH4qFTUVlQxqgp1AtWB6kXdQ42jRKjPaDJaHq2DtkDz0IHoaHQWugBdjm5Et6OvoYfRk+h3GAyGgWFhzDDOmCBMHGYzphhzGNOGuYwZxExg5rBYrAxWB2uF9cCGYdOxBdhK7EnsJewQdhL7HkfEKeGMcI64YFwSLg9XjmvGXcQN4aZwC3hxvDreAu+Bj8BvwpfgG/Dd+Dv4SfwCQYLAIlgRfAlxhB2ECkIr4RphjPCGSCSqEM2JXsRY4nZiBfEU8QZxnPiBRCVpk7ikEFIGaS/pOOky6SHpDZlM1iDbkoPJ6eS95CbyVfJT8nsxmpieGE8sQmybWLVYh9iQ2CsKnqJO4VA2UHIo5ZQzlDuUWXG8uIY4VzxMfKt4tfg58VHxOQmahKGEh0SiRLFEs8RNiWkqlqpBdaBGUPOpx6hXqRM0hKZK49L4tJ20Bto12iQdQ2fRefQ4ehH9Z/oAXSRJlTSW9JfMlqyWvCApZCAMDQaPkcAoYZxmjDA+SilIcaQipfZItUoNSc1Ly0nbSkdKF0q3SQ9Lf5RhyjjIxMvsl+mUeSKLktWW9ZLNkj0ie012Vo4uZynHlyuUOy33SB6W15b3lt8sf0y+X35OQVHBSSFFoVLhqsKsIkPRVjFOsUzxouKMEk3JWilWqUzpktILpiSTw0xgVjB7mSJleWVn5QzlOuUB5QUVloqfSp5Km8oTVYIqWzVKtUy1R1WkpqTmrpar1qL2SB2vzlaPUT+k3qc+r8HSCNDYrdGpMc2SZvFYOawW1pgmWdNGM1WzXvO+FkaLrRWvdVjrrjasbaIdo12tfUcH1jHVidU5rDO4Cr3KfFXSqvpVo7okXY5upm6L7rgeQ89NL0+vU++Vvpp+sP5+/T79zwYmBgkGDQaPDamGLoZ5ht2GfxtpG/GNqo3uryavdly9bXXX6tfGOsaRxkeMH5jQTNxNdpv0mHwyNTMVmLaazpipmYWa1ZiNsulsT3Yx+4Y52tzOfJv5efMPFqYW6RanLf6y1LWMt2y2nF7DWhO5pmHNhJWKVZhVnZXQmmkdan3UWmijbBNmU2/zzFbVNsK20XaKo8WJ45zkvLIzsBPYtdvNcy24W7iX7RF7J/tC+wEHqoOfQ5XDU0cVx2jHFkeRk4nTZqfLzmhnV+f9zqM8BR6f18QTuZi5bHHpdSW5+rhWuT5z03YTuHW7w+4u7gfcx9aqr01a2+kBPHgeBzyeeLI8Uz1/9cJ4eXpVez33NvTO9e7zofls9Gn2eedr51vi+9hP0y/Dr8ef4h/i3+Q/H2AfUBogDNQP3BJ4O0g2KDaoKxgb7B/cGDy3zmHdwXWTISYhBSEj61nrs9ff3CC7IWHDhY2UjWEbz4SiQwNCm0MXwzzC6sPmwnnhNeEiPpd/iP8ywjaiLGIm0iqyNHIqyiqqNGo62ir6QPRMjE1MecxsLDe2KvZ1nHNcbdx8vEf88filhICEtkRcYmjiuSRqUnxSb7JicnbyYIpOSkGKMNUi9WCqSOAqaEyD0tandaXTlz/F/gzNjF0Z45nWmdWZ77P8s85kS2QnZfdv0t60Z9NUjmPOT5tRm/mbe3KVc3fkjm/hbKnbCm0N39qzTXVb/rbJ7U7bT+wg7Ijf8VueQV5p3tudATu78xXyt+dP7HLa1VIgViAoGN1tubv2B9QPsT8M7Fm9p3LP58KIwltFBkXlRYvF/OJbPxr+WPHj0t6ovQMlpiVH9mH2Je0b2W+z/0SpRGlO6cQB9wMdZcyywrK3BzcevFluXF57iHAo45Cwwq2iq1Ktcl/lYlVM1XC1XXVbjXzNnpr5wxGHh47YHmmtVagtqv14NPbogzqnuo56jfryY5hjmceeN/g39P3E/qmpUbaxqPHT8aTjwhPeJ3qbzJqamuWbS1rgloyWmZMhJ+/+bP9zV6tua10bo63oFDiVcerFL6G/jJx2Pd1zhn2m9az62Zp2WnthB9SxqUPUGdMp7ArqGjzncq6n27K7/Ve9X4+fVz5ffUHyQslFwsX8i0uXci7NXU65PHsl+spEz8aex1cDr97v9eoduOZ67cZ1x+tX+zh9l25Y3Th/0+LmuVvsW523TW939Jv0t/9m8lv7gOlAxx2zO113ze92D64ZvDhkM3Tlnv296/d5928Prx0eHPEbeTAaMip8EPFg+mHCw9ePMh8tPN4+hh4rfCL+pPyp/NP637V+bxOaCi+M24/3P/N59niCP/Hyj7Q/Fifzn5Ofl08pTTVNG02fn3Gcufti3YvJlykvF2YL/pT4s+aV5quzf9n+1S8KFE2+Frxe+rv4jcyb42+N3/bMec49fZf4bmG+8L3M+xMf2B/6PgZ8nFrIWsQuVnzS+tT92fXz2FLi0tI/QiyQvpTNDAsAAAAJcEhZcwAACxMAAAsTAQCanBgAAAAfdEVYdFNvZnR3YXJlAEdQTCBHaG9zdHNjcmlwdCA5LjUzLjNvnKwnAAAgAElEQVR4nO3dPWwj17338aPnBkhSPHrZIsE1JF1AQVbFbqd1sdpuIwEiUoW4Eo0gQGTAoqtAxbUgF6ZhQAhghWpULtUIyC1WXEBIFS0gXSONuEWWHVloYwswuXCANNK4ShGAT3Gu5xnPkKN5PfOf4fcDwyC1Q/LwzPD85pw5MzMxGAwUAABZ+z9ZFwAAAKUIJACAEAQSAEAEAgkAIEKKgWRZVr/ft5/2+33nUwAAnH4QfNF33333wYMHSqkHDx7s7OzcuXyz2by9vdVL1mo1pVS/3y+VSpVKJWppAQCFFTSQ+v3+gwcPjo+PnX+0LKvb7U5OTj58+NC1fLlc7na7H3zwgVKq1Woppfb29izL+vTTTwkkAIDXv3322WdBlut2u19++WWr1Wq32w8fPvzRj35kWVa5XP7hD3/4l7/8pdPpPHnyxLl8pVL55z//qZR68uTJ8+fPf/rTn3755ZdffvnlRx99lMbXAADkXdBjSJZlzc7Olkql6enp7e1tpVSz2axWqzs7O4eHhy9fvvR/+fPnz5VSvV5vc3MzXoEBAMUUdMiuVCqVSiWl1PLy8tnZmVLq9vb27OxMP56bm1NKtVqty8vL+fl576Dc2tqa/mO5XB76/ouLi2/evLGfcv0IABg3QQPp6OjowYMHy8vL9l+mp6crlYqOmaOjI6XU8vKycwHb/Px8r9fTjy3LGvURhBAAjLOggfT48ePt7e21tbVut6tDaGNjY3t7u9frdbtd3XkapVKplMvler3e7Xar1WoCpQYAFM5EqH5Jq9Wam5vTA3Sj/hL8tU6Li4tXV1fBSwIAPiYmJrIuwhhJanwrXCClh0ACkKCJCSmNW+ElWNVcOggAIAKBBAAQgUACAIhAIAEARCCQACBLFxcX19fXQ//p+vr64uIiwgtDLSMHgQRgrE189MWoB2acn5/7BNL5+XmEF4ZaRo4Qt58AgOIZHDyd+OiLwcFT/dT5OLjd3d2bmxul1P7+/vX1dbvdbjab+/v7+oHu5bTbbfupc3nXW93c3Ozu7s7MzCilVldXr6+vd3d3r6+vG42G/iCl1MzMjPOFrk9/9uyZUurRo0f2hQj0e/74xz/+7LPP9GvX19cXFhbCfs200UMCMO50JqmoadRoNO7du9doNH72s581Go2bm5vz8/OLiwv7gV7Mfupa3vtujx490nnz97//XefH6upqo9FwvsR+oevdms3mxsZGo9H46quv7M9dX1//8MMPHzx48OLFC6VUs9kUmEaKHhKAcRB8FC7Ikq7Qev36dbvdPj8/v7m5WVlZUUrZzb2r3ddPvcs72Rm2v79/cXGhu0oLCwvn5+ftdlvn0Orqqj2U53q3jY2N3d3dZ8+ebWxs6AWePXt2fX29sLCwsLBQrVZnZma8HyoEgQSg+O7s99g5FKGHpDsxKysr19fXQQ7YuJZ3HSVaWFhot9tLS0uNRuNf//qX85+Wlpba7fbCwoIeoBv6bu12W3eDVlZWlpaWlFIff/zxzc3N559/rntdeiwx7Hc0g0ACMO70SJ3z/6Fe/uGHH66vr6+srOgezJ2Z5Fre+6+7u7s6ez7++OOvv/7a9cLXr1+/ePHi4uJCHytyvdvNzU21Wl1YWJiZmbH7Z+vr68+ePbu4uFhdXRU7Xqe4lh2AQjJ/LbuLi4tQQ2E+y9/c3LTb7VH/OvSFzj9eX1/f3Nzo7pGLzr9k77qQYFUTSAAKiIureukpDy9evNDHpZJCIAGAHwLJGK72DQAoGgIJACACgQQAEIFAAgCIQCABAEQgkAAAIhBIAAARCCQAgAgEEgBABAIJACACgQQAEIHbTwCAerP+jn5w/8U3hj/64uJC3z3P8OcKRCABGHdv1t+xc0gnk8lYOj8/X11dJZAUgQRgzDnTSCl1/8U3dm8puN3dXX0X1/39fX3bVn1jVv1A35K83W7bT53Lu97q5uZmd3d3Zmbm3r17Kysr+r7jCwsL19fXjUbjxYsXS0tL+n6y+o5H+jZ9jx49Wlpash8ne8cjYwgkAMXnnzHef/Vf3tV/ajQa9+7d29/f39/f1zlxfn5+cXFxcXGhH+jFbm5u9FPX8q43bzQaOlHW19eXlpYuLi42NjbW19d3d3fb7fZXX32l+1L6wfn5+cbGxsrKyu7u7ldffWU/Dlgt0hBIAIrPZwjO1UMa+hd/r1+/brfb5+fnNzc3+rat9vibayBOP/Uu73R+fv7ixQt74ZWVlc8//7zZbF5fX6+urroW3tjY2N3dffbs2cbGxsLCgv04eOFFIZAAjDU9RhfnoNHMzMz+/v7Kysr19fX19XXY5c/Pz53/qkfn7BuQf/75548ePdrf3x/a72m32zq9VlZWyuWy/XhpaSmPB6UIJADjznXcKGw4ffjhh+vr6ysrK+12u9Fo3JlJruW9/7q7u7u0tNRut1dXV+/du9dsNnWnSim1urq6u7urX7u6ujozM6OPMM3MzNy/f99+nMc0UtzCHEAhmb+F+cXFhXf8LcLyFxcXMzMzSqlms7mxsaGTSWdMu91eWlq6vr6+ubmxu1DOp65/MiPBqiaQABSQ+UBKyvX1tZ4jt7S05J2DJxCBBAB+8htIuZNgVXPpIACACAQSAEAEAgkAIALTvgEU08TERNZFQDgEEoACYkZDHjFkBwAQgUACAIhAIAEARCCQAAAiEEgAABEIJACACAQSAEAEAgkAIAKBBAAQgUACAIhAIAEARCCQAAAiEEgAABEIJACACAQSAEAEAgkAIAKBBAAQgUACAIhAIAEARCCQAAAiEEgAABEIJACACAQSAEAEAgkAIAKBBAAQgUACAIhAIAEARCCQAAAihA6kfr9vWVaQJS3L6vf7zhc6nwIA4PSDUEtbllUul//6178GWbjZbN7e3u7s7CilarWaUqrf75dKpUqlEqGgAIBiC9dDOjg4mJyctHtIlmW1Wq1Op+NdslwuHxwc6MetVksptbe3d3h4qB8DAOASood0dHQ0Pz/f7/enpqbUd72ltbW1t2/fzs7O6p6Q7fT0tF6v68eXl5fz8/MnJydKqcPDw+QKDwAojqA9pE6n0+l0tra27L80m81qtbqzs3N4ePjy5Uv/lz9//lwp1ev1Njc3oxYVAFBkQXtIR0dHs7Oz9Xq93+/XarVqtXp7e3t2dnZ2dqaUmpubU0q1Wi3dGfIeJVpbW9N/LJfLoz5iYmLCfjwYDMJ+EwBArgUNpK2trW+//VYp9erVq1KpNDk5OT09XalUdMwcHR0ppZaXl5eXl72vnZ+f7/V6+rHPDD1CCADGWdBAevjwoX4wOTmpU2djY2N7e7vX63W73VKp5PPaSqVSLpfr9Xq3261WqzFLDAAopImY/ZJWqzU3N6eH7OIsubi4eHV1FackAIBcixtISSGQAGDMcekgAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIP8i6AADwPW/W39EP7r/4JtuSwDACCYAgb9bfsXNIJxOxND4YsgMghTONFFE0fggkAOI4R+3sxyg8AgmAXKTRWCGQAEjh7A8xtWEMTQwGg6zLoJRSi4uLV1dXWZcCQPacvSLSaKzQQwIgix1CpNG4IZAAACIQSAAAEQgkALK4zkbC+CCQAAAiEEgAABEIJAAScY2GMUQgAQBE4GrfAHA3Ttc1gEACgDu4Jv4xDzAlDNkBEERgWy+wSEVFIAFAIEyySBuBBAAQgUACAD+um2IwgpceAgmAUHJORXKWhDRKD4EEAIEQRWkjkAAAIhBIAKSQf3hGziiiWPowW7Ra4sRYAEAynLsUOpNC7WHQQwIAJC9CZ5dAAoA7yB9LlMCuJeeMxFBjdwQSALkEHrMxU6Q4R2LkCFt+jiEBgCwxj8RkRUe1q6gcQwKAvHK16bmIIpvdfYx2SQsCCYAIHKdxiXwkJnP3X3yj/wv7QgIJAISyextZF8QQAgkA/Aw9LpJeSAx98zHpOxJIACCL6/riOUqjmKUlkABAHN2s6yMxDNkBgAhj1SLb8tUxSgqBBAChmYzJvERy/BDlxFgAdzBwbzqxHQKxBSskekgA/OgWedwOZkgzJpVPIAEYyXvVgHFoFhFBIl1JAgnA3cghr/TieWjjPg57AwQSAORD4TOJQAIAiEAgARhJ75LHuX5zgsUw/7kCp9jJ7CQlVVEEEgA/Bm6FILDdDyiNeMhvbcRHIAG429g2kQKllIIS7lHLibEAMNZcfbKwXbQEu3T0kAAgZ8zf/8IMAglAIDIPp/szMAyVVbUk9bl2/ybz8TrFkB0Af/k9xu4seYRvkckXjzBcZj+OVlo72KIlXLK1FK6H1Gq1LMsKuLBlWf1+337a7/edTwEgPYW/6JH9BfVlBuO8j3LcfinDyf0qeCBZllUuly8vL8vl8tnZWZCXNJvN58+f68e1Wq3RaNRqtZOTk4glBVBEAZu/yIlSsBwaJU7/xs6hzOsq6JBds9kslUpbW1vvvfderVYrlUpKKcuyut3u5OTkw4cPXcuXy+Vut/vBBx8opVqtllJqb2/PsqxPP/20Uqkk+hUAIEu6Qc/XKcPOjpH6/nBf8DdM/FsHDaStrS39oNPpzM3Nqe/6TGtra2/fvp2dnd3Z2XEuf3p6Wq/X9ePLy8v5+XndNzo8PEys7ADMyrDlDctZ1LyUORrncbKh39R7Oyv/CslwLYeb1FCv11++fFmtVpVSzWazWq3q7s4vfvELVyC5PH/+vFqt9nq9zc3N4+PjGAUGJHLtb0ICZ+8hw3NrUvpQV2w409deQHlmdjj/KfibGxMukHZ2dqrV6srKSqVSub29PTs708eTdJ+p1WrpzpB3UG5tbU3/sVwuj3rziYkJ+/FgMAhVMCBDrragSLFUgL6Fq5kuEmfiegffRn1lySs0aCDVarXl5eVSqTQ1NaX/Mj09XalUdMwcHR0ppZaXl5eXl72vnZ+f7/V6+rHPJD1CCHnkbLKLN5ULAWXVpQgy8uYMreDb553fKI3vGzSQKpXK9vZ2p9Ppdrt6yG5jY2N7e7vX63W7XT3Hwee15XK5Xq/brwWKwfubzNFRlmgK/wWLJ/LKMr+ugwbSw4cPT09Pu93ue++9pwfopqamjo+PW62W/RcX51Gl09NTnyWBvCveYJ0ZqbZ33s6r5BWURvHkf2uXEMeQpqamvCNyQ8fohgq+JJAX3h98vn7/+ZK75lUC1zBd2NobVecprQguHQTE4rrySpGay4J9neCifXGxG0DMIpn8XlxcFUgG0xmEY8qJfAQSkAB9/ZUxafJy8TVldlZyyrXG06tbhuwAIGPy49PMXgg9JCCWwl9VukjSXjvjsPZTPdGYQAIwhIF9dvndgrB0My3hstmJc10YMKXoJZCA6IY2qeOwmyxcnKiL/FrXDYoKuQ2kvQNBIAEILasGN/7nFjUqjEk1kwgkIHm0euPGe1Ht4NtALoYuzWzSzLIDIspFOwIz7LNH7Ya7kJtHnIs+BEEgAakQe95+EPktuQpQ+LRXTVH7xwY2CQIJAP5XnKCKedU4KI4hAdEEabmKuqesxfx2ue6EjWJfsKPYqz49BBKAsSMnMAoZzJERSECK5DR8Y0JO+86qj4BAAkKT0+qNoRw19DkqqhAEEoDvGZO4JS0EIpCAcMK21wVu+Ar21dK7ibjJT8w1AglAQbCvkHcEEpA6Gj6XseoZsPaDI5CAEGJeRrqQNybAncikgLhSA2CO844ykd8k1csBjFXfReX8Ck/FQw8JOZPHfkaCd5W1r+Ap5HIAmZRh6IcKzxVvmYUXOBMEEnJD51CGbXHMFiSNMkvIpGIgHiQgkJAnSfUzsqLLHLPYeewjCmffMCLVT8njFmsYgYT8yfWvWnfv1He54l3AJ2zsO78JGa8rAGfHyECVstb8MakBCCTBGxMoR1fPzpg7H6thcyKSHWgyM2yVxqck9Z5m5jiQSaMQSMgT182hczToP6qoQyNHfX9wz/VaaXvZwRtxu9jSVp+xKnXNkOTAlQuBhNxwjXSZ/CWn3XDYbVPABsv5x7w0anY5vQ8kMFYS51omk1wIJOSPsebMGQ9mGg5X0qT9cXAxtpZZuUMxqQHwY89BSLsRGfr+QRrH8WzdnN86/uE9+z+T/SQzH5Qv9JCA4cyPLHmH7EK9UH4bJ+3QkWayPK41lYu1ZhI9JOSJ4Um65jm/neGPTnYaoc+S8T8u15x9svGsAR/0kICRMtmjj/ZxwjtJzu6m/ovYohowzt/dHz0k5JiZMxn1A8nNvSa2y+jq19qH5RJ/c+QdPSRgCPbok0JgIDh6SIBbqnv0qcqwkzT0o0kjhEIgwYRErgc6tHVLtgkuwKHmaBUi+XI+/h+R6vvDMIbskDqx81xdcxZElS13nBfMTnuGhetTWHGFQQ8J6XI1FmmMKUXuE6iCzsF1nexp4BPt2jN25lbaGxUyQSDBEFFNxpvvbvSnvmvaCpNGTveN3KXC2wMmIRANgYR88N/pjtwC2n2jIjWgBnqlQ5kP9ULuRowzAgnpcl1zTNTImMwr2SQleA7FXCmuy1fHf8OwRG1UiINJDTBBYP9D7FSLBAW5Aqk9GSFHNRD5on8QjkAqGoE/VDOzoSK/v8CwjMlZFc4ui7cfoxJaO+YTQsi2jWQRSIXialNiBkBRL1LgugpDwb6dzZsQzh6h8+/245iZFLGgwHcIpOJItv/hGtFSSbQ4cToxQV4V9v0L3Ib6fzU7j4vXO0SuEUhFk0i3JqtpWgbk62BJ2qgKiFKcWXYTH31h/9/5oEh8TnV0HppO8ESQpN5ESMMnpyQCUTnIXHF6SIODp840Ghw8zbY8iRs6K2zoYeoEPzHBd9MynNpAg+vCXDVIU5xAUt9lUuHTyP6L8jQiCZ4R4s28XCONhqJOIEpBAsk1QGc/LWQyKd9+wKiICt70jOp+ZSVslngrhzQCcqEggaSDxxlLxYsiFbVhDd5VckVXrntIDEYBuVOoSQ12CDmPJ0Vj8krJQbiG+yMk053THJzXG01V2nP2XOeEAsiLgvSQvMeN4vSQJF9UJs4pQc6xLGMXH8vwsmbFu2oqUGwFCSRn/OjuUeRAGnoKjpBMSurUVOdb3ZlwomoglDyWGRhnxRmyU0nP9ha1Z514JGR7SU0DF6ImjYDcKVQgISD7zNnitdrxD7YByEoxAynOpIYEjzp4LxsRrVQptapFbanvF/Gu5MA4KMgxpGTZmRTzoqKuo1kyz9gN0mqnfVHUxJFDQB4Vs4eUiESmDtt9tchplHib7uxAZNiHCNITpYuDcVbIC3L6K2wgxT8VKT59HaNRTzNkX4CV5h4QaNSlZ4b+sUhXlC5sIEkwOHjqmo8etpNEFwEYQ0EuFe3a55Z5RCCs4gSSzPVhl0pCj00azlotDNdOOpt6fPbl0HyatfhHBKQp8qSGmGfIxqc/3bnFhCqPqO5RhHkNosqPVDm37cI0jlkJMl436lV5r/kiB1LmvNfWy/vmYhJ5li9s5ElxXir6zr3Ygt3ZoDhDdglyXp0zqzElmmPkhZ6t4xy1kzN/J6cCXiradVaJocKlqeCBFOHITSHnrojFYaQC8M7WiTB/BzbnUL82qjLt49M+y+RLuEDqdDr9fj/gwpZlORfu9/vBX5sh73qNs6ajjafL7B4RHhjFNVhXsH04w/umzr5R8JcUo86DBpJlWeVy+eTkpFar1ev1IC9pNpvPnz/Xj2u1WqPRqNVqJycnEUvqK/HjqIlPFirMFhOQzExFSlxtaDH21m3eH+9Y/ZZNChpIzWbz8ePHe3t7x8fHL1++1H+0LKvVanU6He/y5XL54OBAP261Wkqpvb29w8ND/diksEngXdjwxjduTfnQjte4VULeJXt7TIGCnBiUoAgfUYxqDxpIa2tr7733nlLKsiz9F91nury8PDo68vaZTk9PP/jgA/348vJyfn7+5OTk5cuXh4eHCZU8LVdf/8bVFBZsdw9IWzEaRxcmtRsQNJDm5ubm5uZardb7779frVaVUs1ms1qt7uzsHB4e2n2mUfTYXa/X29zcjFdgP0nddzyp31KEbZeeAXKn8G300GmE2RZpqALsB4Q4D6ler799+/bw8HBubk4pdXt7e3Z2dnZ2ppTSf2m1WrozVKlUXK9dW1vTfyyXy6Pef2Jiwn48GAzCfAulPPeuVo7rokY4Q9ZwDz0vkwWKehclJC7z09ITFOrEIMQRtId0cnLiTCOl1PT0dKVSOT4+Pj4+fvLkiVJqeXl5Z2fHm0bz8/P2Y3vEz2vgEO5LeMbZEmkx9TZ3/8U3V1//Jv67qdH7L7qJz/zMpwSFDS3uqpdfY9I6uyZhp9oRiVOlee8kBQ0kPeF78ztKqY2NjbOzs3q9vrm5OTk56fPaSqXy6tUrvaQe7kuWcwXY7Vqclt3k4dkMG99CnmjFaWRyRPv5uCa4Slh9Q69AIaFgBTSI5/LystfrxV/y/v37kcug/ut/rv7z3/Xjq//8d/3YfuBc7M63cr3E9eZhSxXkj86S+xQjQXYZ9IMg1eIUpGwRyp/g13d+r7DfDqEEqd4Iq0D91/+IWn2jimGX08zHmXyHrMS9UsPy8rI9iJfUkhEM7Q+F7XkI2eVJu88U8/BYSoOKCb5t8S6BPG7yMp9NX5Ai2XZD/rdOVXEuHeRs0RJs09PeOMwfNPIOiQhJ4kSIvSli8QRsOkM12QIvi3fn1xR42EZgkQLKOJASP//Z577j/ispwx0TO5PM3FPce89AUXtkMWsg/k0RkSFpl8ULFbrFvimUmaOzGQdSgrOr41yiO40rDwW5yaPNLrmBCQ7OI7RpbFJxEiWRzmLxLoEsUKifTPAddt2s5/GiD64fVLQmRfJ4nWtd+Bc18lrLfsgu5m1UElyFo1rSxf/477RH1czMtbN34pw1JvAHH6e2XVdKzlGLliMptZuujTPzM34inLyoxNwoL40tP8jR2ZiDXlneoG9U0dNbl0M38SBjxG/WUyqROT6dtrBvlcbpsfoN7TOaI7+/6yqfknc5x4p/ungb8dytNfOtmWEBv2DMQa8se0j26LCzEYm2/iK3X0VqsEYN8rqGQWxpzBGKxrv64sz1KNI6TVW0Yx6Rq3fUxqbfMNQQtwERukdJtWZJSbb2hh6dtQ+eOf/TCzgbn+CfkvEtzL3jxYmswoC72OmNPGSyIboqMKdTn/N7jYbcVbVzgzFQeO+Os5DRLa9EjgCFfRP52489Hu7ccoaWOfKazX5Sg8/TgIJ3j+yqjDBRJKl9De9uS7KT1J17u+b30aJ9F+dgXZz38TKwi534TFGTwh7BjdNoeg/7B//cXHANF0v4aglu/65V7/PO9pqN8NHZT2qwGfglu8Y3lWOjKcD107y94yBVGm27SeP0qUROIzO/m5ngTFHDJgKc8TNqHDgaZ/c94N6S+VG7RFaikMHwpLiGslwPfJYMWwmCAimCiY++8N6+SBvVYub6wKm/oYO8GZYnFD3l3cyZWImLOVM0K0PHJ1yx5LMPF9ZEHs5ZzmqXItXDB9HGhFxvEnDhmINeGR9DMmxoByJIld1/8c3V+jtKSW8ohw7y3inBo3cxxYyiUXM30vt2uZ5b5RxaGXoAUnNmUpzv5RqyC9XGmTnElexHSPhZufooEcoT7StE/taCekgGOrkx+xBBihd2/SXYIbCbFf20YIMGMg3dhHKRRtqoowIDxyU8EuzWOI+EC9k4Xd8uw5KkxDVMGuq15gNVUCBFMGq8zoeEn0RKHx1wkDc9eRxti887vC5zMMol4FGBZMeBnUf7Q71Per9W1zun0Um6c7FUe/CR9ycy6d7lO5D8DT2MNPQnMZ4tqVOEH7y0ewn6/H5S3flwjcwMvn9+RkofGlOotkbCPlyqcnoIMIihOxA+KzHz6Tk5DqQ36+8s/sd/m/xE+T/IzMesx5n3cK69wSQ7V82k4JN9cyrIVMM4Mq8x5/6EcmyWPjMdMmxGZAVS2JUnrfHN8ABSfJn/cuIQm8TeTFJZlzbapKmsxoGdJUm1m2s/kLkhRTNqf8LVj3fWbbYbp6xAAlKSXnN25w/YNTsgR2lUbBOeaw2nd3DXf5QsvTVy5/6Etx+f7eaR10AKOF6X+HGOXHcj7hT22wW5KWJ+R6sS5B0SyaQS8p5GobZPn+3N7hZ4+0PF/oEPNXTjzKoS8hpI5iWbbWls9zKbGwO/cJlf3CnZuWrRyK+lZHk3PLup9Z+sn0YtSc65gaTLwooLpCBrLtlDL0m9W9gffIQ56zklZIQ6250A55Lm26bCpFGoleja8FSmTe3QkgtZL65pohlmZ/aBpK8WoxK6W6iXtNnJwqU3jK7S+fkJ+Un7MzxXbZxHR12cR+/kbydZ8Y5bZlWSjANJ904i3Hq8YDva8Yltl70TTDOsvazWnWswJO2V5fqaYreNaAKuRO+xkMx/tjKbDqfMt5Pse0he/qtN1FTpQgr1s7nziqjOQZI7x6yE/1xdJDf0zqlTYguZKleXVAlobV3GdtX4yDKQkr1PaLaibVv5an8jc35Nu2m489S8gO/MT9rL2zMYky3NZmxKd1hCiiFWloHkih/7YNJQ+l9j3v9t6NtGeMP4W1Wy/bwM22XntO+hlTxq+D6rU/MktAhpl0HCpL5UBTyzR86hERtnQfgTNGTnPJJ09fVvXFllH2rKsAuV3w5cWAFbTGes+gzZjXpgP42WRpxh6sOuTwkBnB5Xs27v3CgBl5YYyjnJQonfiszL+H5Ijib+/9qtm/7/xEdfXH3X+ttR9P8zyexV7ODlXF/2Dci9Pb8gP7lR43UF/q2amdrgelAYQ7vUeWncnSXPRYFNyr6H5JxlN+qfvMNBolZkIXfVA54Qph/EHH4cOvHBX5xvXexOg8p6gzTAe4JRXr4yh/f8ZR9IQehkitPqibqGkN2NENsyBhnp9q6OOMW81TAAAApiSURBVAfGXKfmKd+JD5mPv+el+SuwnB4kK/zhvZjkBtKoxjrBu04wg3yUgCey2HO+1YjKDN5wDz3+7DPxIdSby5TS7kjeqyWIXPcznEe5clRsM+QGks1u9fyn4Y2n9Fof5zlDPh/hOryX1Id6/+48scZ/4YAfVMi2YBzSSHnOK8pRP8O1goq6HUaW8aSGgNyNXaRVOPSQexz2xmTPnMnLr+JOQ2cZjPp2o6o07QrJqsIT/NwcHYoXxTXAm6M6vHPS6ZjLQQ9JFPtYlLfJDrhh5WKcUOZIt7dUGRZGoBy1yzF5B3jpZxSD6EAaup0J+dUZO5EzKzFHulOqk6TG3+0jUpnPj0hKITfCO9HPKBjRgZSqmFPCvD2kaC1anIY17SsaOP8vob1OcPw98vyIxOtcSN0CEuTjGFJSkjqM5Dzgr9g780jv1KhRD6K9W75OYfFRgK8AqHHuIcVntwLB93BzcQBpqILtyIuaN5zI8GOC5QGyIj2QxLaD3nk+GRbAjIBfMxftY7RZG7n4akB+SQ8ksbyXFci0OIaI3T+IIO/nJ7rOygIKIGeBFH8X9c4bypmX0wZxqLz0IaSdnxi8AK6RRjU2O0MYBzkIpGQbC+f9e6K9g4SbUGTY7mfedidi6PwI/+8lJGtd5ZRQJCApOQikBLnu32M4V0R1y+IoRiZ5Zfu9gh+iU45+Ur4u4wb4G69A0jLMpMKQfM5yHK6zZQXyzsXIe50DtnEMpLwrQLsv2eD7d7DV5NS56/KJWRcHSFI+AinxsRRpQ2d5HARzXTzJe3uIXHN9O9cDAx/qZXfdsj3fAEhPDgLJe8m4yD9C1xidyYM6MT/LfOMYULQL8ORC5MsLReOzZu289w7QkUkokhxcOijZdsGZSUL6Sa5vNPQLOq+xL+d6+64jLhKKlCzntpf2txs6Tujz0VxXFMUzMRgMsi6DUkotLi5eXV35LCCt1QvV49ER6LO8PRTjnzSiKmHUjrmEsiVi6BdM+9s5twQDHweIkoNAyqRduFPAQLKjKHgmOZ/6y7wStGI3ncZ2Agof8MCdchBImqjOgQoWSPYy3ge2UM2QtEpQjhsLKUmlSor97UwOkxa1MoE75WBSg8r/Zcd8OA9T299xVBpJq4TMrzCbqkwuL+QasgPGSg4CSdplxxIXpFmXWQneGV9ZlSQNQy8vlOonFjvggTvlZshOmuBDds4lhw7Z2XPZnQeQCta4A8CdCKSIgk9q8EkjAIAtB+chyRT8bujSTnsCAJkIpBTRJQKA4HIwqSGnSCMACIVAAgCIQCClgu4RAIRFICWPNAKACAgkAIAIBFLC6B4BQDQEUnSu2/0p0ggAYiCQAAAiEEiJoXsEAHEQSMkgjQAgJgIpLtdhJABANFzLLiLnJVPJJACIL1wPybKsfr8fbeF+vx/8tbnAGB0AJChcD6nZbN7e3u7s7IRduFarKaX6/X6pVKpUKhEKKop9xMh5EwoOIwFAHCF6SOVy+eDgwPkXy7JarVan0/FfuNVqKaX29vYODw/14wLgLkcAkKwQPaTT09N6vW4/tSyrXC6vra29fft2dnbW1W1yLnx5eTk/P39ycqKUOjw8TKLYGeOUWABIXPRZds1ms1qt7uzsHB4evnz50n/h58+fK6V6vd7m5mbkTxRFZ5L9H2kEADFFn2V3e3t7dnZ2dnamlJqbm1NKtVot3RnyHiVaW1vTfyyXy6PecGJiwn48GAwiF8wYu59EGgFAfNEDaXp6ulKp6Jg5OjpSSi0vLy8vL3uXnJ+f7/V6+rFlWaPeMBch5EIUAUBSogfSxsbG9vZ2r9frdrulUslnyUqlUi6X6/V6t9utVquRPxEAUGATMfslrVZrbm5OD9nFWXJxcfHq6ipOSQAAuRY3kJJCIAHAmONadgAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAgAoEEABCBQAIAiEAgAQBEIJAAACIQSAAAEQgkAIAIBBIAQAQCCQAggrlA6nQ6lmUZ+zjcaWJiIusijB3q3Dzq3LzIdW4okDY3N09OTt5///1Wq2XmEwEA+fIDA59xdnY2Nze3t7fX7/drtdry8rKBDwUA5IuJHlKn03n48KFSam5u7tWrVwY+EQCQOyZ6SEqpubk5/eDx48ejlllcXDRTGGj379+nzg2jzs2jzs3TdX51dRX2hYYCqd/v6wejekhv3rxxPh0MBqmXaexNTExQz4ZR5+ZR5+ZFrnMTQ3YPHz7s9XpKqX6//+DBg6HLDL7PQKlAPZtHnZtHnZsXuc4N7TuUy+XHjx+/evVqa2urVCoZ+EQAQL6Y68y2Wq25uTn7YBIAAE6MrgIARMjy0kGWZdmTHTSu5mBAp9NxVjt1bkCr1XJVMtVuQL/fd1YydW5e2Dr/t88++yy1wtzhj3/84+vXr588eaKfbm5u9vv9o6MjRvZSYlnWr3/963/84x9/+tOf/va3vz158oQ6T5uu88Fg8Pvf//4nP/nJz3/+c8WmboRlWb/85S9/97vf6afUedrefffdVqtlty0qWp0PMvKrX/3q/v37f/jDH/TTP//5z5988slgMOj1er/97W+zKlWxNRoNu8KfPn1KnRvQaDQajcbAUclUuxmffPLJ06dPb29vB9R5+rwVG63ODZ2H5HV6elqv1+2nXM3BgLW1Nf1Ad6KpcwO2trb0g06no3cSqXYDjo6O5ufn+/3+1NSUos7Tp6u6VqtNTk5Wq9WpqalodS7o9hNBruaAOHTHudVqvf/++9VqVVHnptTr9YODA/37VFR7yjqdTqfTsXcFNOo8VZZlzc7Olkql6enp7e1t/ccIdS4okO68mgPiq9frJycnh4eHlUpFUeem7OzsnJ6eHhwc6KdUe6qOjo5mZ2fr9bq+mrOubeo8VaVSaWdnZ3l5eWtr69tvv9V/jFDnUgIpyNUcENPJycnbt28PDw/1ngt1bkCtVjs7O1NK6bEjRbWnb2tr68mTJ0+ePJmamiqVSpOTk9R52o6Ojly3FopW55kdQ3IplUpHR0f1el1fzSHr4hSTnvC9ubmpnx4fH1PnaatUKtvb251Op9vt6mFSNvW02UOjk5OT+mY31HnaHj9+vL29vba21u129ehLtDqXdWIsV3MwjzpPm2VZ3W7XVclUu3nUedq8NRy2zmUFEgBgbEk5hgQAGHMEEgBABAIJACACgQQAEIFAAgCI8P8Agvz1H8GU+WgAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "figure\n",
    "plot( n, err_rel_chol, \"-*\", ...\n",
    "      n, err_rel_gauss, \"o-\");\n",
    "legend(\"err rel cholesky\", \"err rel gauss\");\n",
    "box off"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Sia l'algoritmo di Gauss senza strategie di pivoting che l'algoritmo di Cholesky risultano stabili. Ciò significa che le matrici fornite in input sono ben formate. L'errore relativo è molto piccolo $\\left(10^{-16}\\right)$ e vicino per ordine di grandezza alla precisione di macchina."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Octave",
   "language": "octave",
   "name": "octave"
  },
  "language_info": {
   "file_extension": ".m",
   "help_links": [
    {
     "text": "GNU Octave",
     "url": "https://www.gnu.org/software/octave/support.html"
    },
    {
     "text": "Octave Kernel",
     "url": "https://github.com/Calysto/octave_kernel"
    },
    {
     "text": "MetaKernel Magics",
     "url": "https://metakernel.readthedocs.io/en/latest/source/README.html"
    }
   ],
   "mimetype": "text/x-octave",
   "name": "octave",
   "version": "5.2.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}