Mat/Wrzodak_Koszarek_Zadania.ipynb

{
 "cells": [
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Zadanie 4.6"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3> Zadanie 6 </h3>\n",
    "\n",
    "Rozwiąż układ równań $Ax=b$ metodą przybliżoną, gdzie\n",
    "\n",
    "$$A=\\left(\\begin{array}{rrr}\n",
    "2 & 4 & 6 \\\\\n",
    "8 & 10 & 12 \\\\\n",
    "14 & 16 & 18 \\\\\n",
    "20 & 22 & 24 \\\\\n",
    "26 & 28 & 31\n",
    "\\end{array}\\right)$$\n",
    "\n",
    "oraz $b=(-1,0,1,0,1)$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 126,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "A =\n",
      "[[ 2  4  6]\n",
      " [ 8 10 12]\n",
      " [14 16 18]\n",
      " [20 22 24]\n",
      " [26 28 31]] \n",
      "\n",
      "b =\n",
      "[-1  0  1  0  1] \n",
      "\n",
      "A^T * A =\n",
      "Matrix([[1340, 1480, 1646], [1480, 1640, 1828], [1646, 1828, 2041]]) \n",
      "\n",
      "Macierz A^T*A jest kwadratowa, więc rozwiązanie istnieje\n",
      "\n",
      "(A^T * A)^(-1)*A.transpose() =\n",
      "[[0.050000 -0.233333 -0.516667 -0.800000 1.000000]\n",
      " [-0.600000 0.216667 1.033333 1.850000 -2.000000]\n",
      " [0.500000 0.000000 -0.500000 -1.000000 1.000000]] \n",
      "\n",
      "Rozwiązanie:\n",
      "\n",
      "x = (A^T * A)^(-1) * A^T * b =\n",
      "[0.433333 -0.366667 -0.000000] \n",
      "\n",
      "\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from sympy import symbols, Matrix, vector\n",
    "from numpy.linalg import eig\n",
    "\n",
    "m = np.array([[2, 4, 6],\n",
    "                [8, 10, 12],\n",
    "                [14, 16, 18],\n",
    "                [20, 22, 24],\n",
    "                [26, 28, 31]])\n",
    "\n",
    "A=Matrix(m)\n",
    "\n",
    "np.set_printoptions(formatter={'float_kind':'{:f}'.format})\n",
    "\n",
    "print('A =')\n",
    "print(m, '\\n')\n",
    "\n",
    "print('b =')\n",
    "b=np.array([-1,0,1,0,1])\n",
    "print(b, '\\n')\n",
    "\n",
    "print('A^T * A =')\n",
    "print(A.transpose()*A, '\\n')\n",
    "print('Macierz A^T*A jest kwadratowa, więc rozwiązanie istnieje\\n')\n",
    "\n",
    "#print('(A^T * A)^(-1) =')\n",
    "#print(np.linalg.inv(m.transpose() @ m), '\\n')\n",
    "\n",
    "print('(A^T * A)^(-1)*A.transpose() =')\n",
    "print(np.linalg.inv(m.transpose() @ m) @ m.transpose(), '\\n')\n",
    "\n",
    "print('Rozwiązanie:\\n')\n",
    "#u=(A.transpose()*A)^(-1)*A.transpose()*b\n",
    "x = np.linalg.inv(m.transpose() @ m) @ m.transpose() @ b\n",
    "print('x = (A^T * A)^(-1) * A^T * b =')\n",
    "print(x, '\\n\\n')\n",
    "\n",
    "# SPRAWDZENIE\n",
    "x, residuals, _, _ = np.linalg.lstsq(m, b, rcond=None)\n",
    "#print(\"Przybliżone rozwiązanie:\")\n",
    "#print(x)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Zadanie 4.7"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3> Zadanie 7</h3>\n",
    "\n",
    "Przybliż funkcją $f(t)=a+be^{t}$ zbiór punktów $(1,1)$, $(2,3)$, $(4,5)$ metodą z zadania 6.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 127,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Przybliżone parametry:\n",
      "a = 1.641485981693081\n",
      "b = 0.06298603382678646\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from scipy.optimize import curve_fit\n",
    "\n",
    "# Definicja funkcji, której chcemy dokonać przybliżenia\n",
    "def func(t, a, b):\n",
    "    return a + b * np.exp(t)\n",
    "\n",
    "# Dane punktów\n",
    "t_data = np.array([1, 2, 4])\n",
    "y_data = np.array([1, 3, 5])\n",
    "\n",
    "# Przybliżanie funkcji do danych punktów\n",
    "params, _ = curve_fit(func, t_data, y_data)\n",
    "\n",
    "# Rozwiązania parametrów a i b\n",
    "a = params[0]\n",
    "b = params[1]\n",
    "\n",
    "# Wyświetlenie wyników\n",
    "print(\"Przybliżone parametry:\")\n",
    "print(\"a =\", a)\n",
    "print(\"b =\", b)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 128,
   "metadata": {},
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "'module' object is not callable",
     "output_type": "error",
     "traceback": [
      "\u001b[1;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[1;31mTypeError\u001b[0m                                 Traceback (most recent call last)",
      "Cell \u001b[1;32mIn[128], line 9\u001b[0m\n\u001b[0;32m      5\u001b[0m m1\u001b[39m=\u001b[39mnp\u001b[39m.\u001b[39marray([[\u001b[39m1\u001b[39m,np\u001b[39m.\u001b[39mexp(\u001b[39m1.0\u001b[39m)],[\u001b[39m1\u001b[39m,np\u001b[39m.\u001b[39mexp(\u001b[39m2.0\u001b[39m)],[\u001b[39m1\u001b[39m,np\u001b[39m.\u001b[39mexp(\u001b[39m4.0\u001b[39m)]])\n\u001b[0;32m      7\u001b[0m \u001b[39m#a,b,t=var('a,b,t')\u001b[39;00m\n\u001b[1;32m----> 9\u001b[0m m1\u001b[39m*\u001b[39mvector([a,b])\u001b[39m-\u001b[39mvector([\u001b[39m1\u001b[39m,\u001b[39m3\u001b[39m,\u001b[39m5\u001b[39m])\n\u001b[0;32m     11\u001b[0m M1\u001b[39m=\u001b[39mm1\u001b[39m.\u001b[39mtranspose()\u001b[39m*\u001b[39mm1\n\u001b[0;32m     12\u001b[0m M1\u001b[39m.\u001b[39mdet()\n",
      "\u001b[1;31mTypeError\u001b[0m: 'module' object is not callable"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "zb1=np.array([[1,1],[2,3],[4,5]])\n",
    "\n",
    "m1=np.array([[1,np.exp(1.0)],[1,np.exp(2.0)],[1,np.exp(4.0)]])\n",
    "\n",
    "#a,b,t=var('a,b,t')\n",
    "\n",
    "m1*vector([a,b])-vector([1,3,5])\n",
    "\n",
    "M1=m1.transpose()*m1\n",
    "M1.det()\n",
    "\n",
    "M1^(-1)*m1.transpose()*vector([1,3,5])\n",
    "\n",
    "##########################################################\n",
    "\n",
    "def func(t, a, b):\n",
    "    return a + b * np.exp(t)\n",
    "\n",
    "x = np.linalg.inv(m.transpose() @ m) @ m.transpose() @ b\n",
    "\n",
    "zbior=[(1,1),(2,3),(4,5)]\n",
    "print('zbior punktów = ', zbior)\n",
    "m=matrix(3,2,[1,exp(1.0),1,exp(2.0),1,exp(4.0)])\n",
    "\n",
    "a,b,t=var('a,b,t')\n",
    "\n",
    "m*vector([a,b])-vector([1,3,5])\n",
    "\n",
    "print('\\n (m^T * m)^-1 * m^T * vector =')\n",
    "z = (m.transpose()*m)^(-1)*m.transpose()*vector([1,3,5])\n",
    "print(z)\n",
    "\n",
    "plot(z[0] +z[1]*exp(t),(t,0,4))+sum([point(x) for x in zbior])"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Zadanie 4.9"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<h3> Zadanie 9 </h3>\n",
    "\n",
    "Znajdź bazę ortonormalnych wektorów własnych dla macierzy\n",
    "\n",
    "$$\\left(\\begin{array}{rrr}\n",
    "1 & 1 & 0 \\\\\n",
    "1 & 2 & 2 \\\\\n",
    "0 & 2 & 3\n",
    "\\end{array}\\right)$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Macierz A=\n",
      "[[1 1 0]\n",
      " [1 2 2]\n",
      " [0 2 3]]\n",
      "Baza ortonormalnych wektorów własnych:\n",
      "[[-0.593233 -0.786436 0.172027]\n",
      " [0.679313 -0.374362 0.631179]\n",
      " [-0.431981 0.491296 0.756320]]\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "from sympy import symbols, Matrix\n",
    "from numpy.linalg import eig\n",
    "\n",
    "#Definicja macierzy A\n",
    "A = np.array([[1,1,0],[1,2,2],[0,2,3]])\n",
    "#A = Matrix([[1,1,0],[1,2,2],[0,2,3]])\n",
    "print(\"Macierz A=\")\n",
    "print(A)\n",
    "\n",
    "def qr_eigval(matrix, epsilon=1e-10, max_iterations=1000):\n",
    "\n",
    "    eigenv = np.diag(matrix)\n",
    "\n",
    "    for _ in range(max_iterations):\n",
    "\n",
    "        q, r = np.linalg.qr(matrix)\n",
    "        matrix = np.dot(r, q)\n",
    "        new_eigenv = np.diag(matrix)\n",
    "\n",
    "        if np.allclose(eigenv, new_eigenv, atol=epsilon):\n",
    "            break\n",
    "        eigenv = new_eigenv\n",
    "\n",
    "    return eigenv\n",
    "\n",
    "#Obliczenie wektorów własnych i wartości własnych\n",
    "#eigenvalues = qr_eigval(matrix)\n",
    "eigenvalues, eigenvectors = np.linalg.eig(A)\n",
    "\n",
    "#Dekompozycja\n",
    "Q, R = np.linalg.qr(eigenvectors)\n",
    "\n",
    "#Normalizacja wektora własnego\n",
    "orthonormal_basis = Q\n",
    "\n",
    "print(\"Baza ortonormalnych wektorów własnych:\")\n",
    "print(orthonormal_basis)"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Zadanie 3.9"
   ]
  },
  {
   "attachments": {},
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Zadanie 9\n",
    "\n",
    "Oblicz metodą ''power iteration'' wartości własne macierzy\n",
    "\n",
    "$$\\left(\\begin{array}{rrr}\n",
    "1 & 2 & 3 \\\\\n",
    "4 & 5 & 6 \\\\\n",
    "7 & 8 & 9\n",
    "\\end{array}\\right)$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Dominująca wartość własna (power iteration):\n",
      "16.116843969807043\n",
      "Dominujący wektor własny:\n",
      "[[0.107761 0.132408 0.157054]\n",
      " [0.244037 0.299852 0.355666]\n",
      " [0.380312 0.467295 0.554278]]\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "A = np.array([[1,2,3],\n",
    "              [4,5,6],\n",
    "              [7,8,9]])\n",
    "\n",
    "S = np.array([[1,1,1]])\n",
    "\n",
    "def power_iteration(m, n, s):\n",
    "    #punkt startowy\n",
    "    st = s\n",
    "    eigv = 0\n",
    "    for i in range(n):\n",
    "        st = np.dot(A, st)\n",
    "        eigv = np.max(st)\n",
    "        st = st/eigv\n",
    "    return eigv\n",
    "\n",
    "def power_iteration_vec(m, n, s):\n",
    "    tolerance = 1e-6\n",
    "    #punkt startowy\n",
    "    x = s\n",
    "    for i in range(n):      \n",
    "        y = np.dot(m, x)\n",
    "        #Normalizacja wektora\n",
    "        x_new = y / np.linalg.norm(y)\n",
    "        #Sprawdzenie warunku zbieżności\n",
    "        if np.linalg.norm(x - x_new) < tolerance:\n",
    "            break\n",
    "        #Aktualizacja\n",
    "        x = x_new\n",
    "    return x\n",
    "\n",
    "y1 = power_iteration(A, 1000, A)\n",
    "print(\"Dominująca wartość własna (power iteration):\")\n",
    "print(y1)\n",
    "\n",
    "y2 = power_iteration_vec(A, 1000, A)\n",
    "print(\"Dominujący wektor własny:\")\n",
    "print(y2)\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Dominująca wartość własna:\n",
      "16.116840810027913\n",
      "Dominujący wektor własny:\n",
      "[0.231970 0.525322 0.818674]\n"
     ]
    }
   ],
   "source": [
    "import numpy as np\n",
    "\n",
    "# Definicja macierzy A\n",
    "A = np.array([[1,2,3],\n",
    "              [4,5,6],\n",
    "              [7,8,9]])\n",
    "              \n",
    "# Definicja wektora startowego\n",
    "x = np.array([1, 1, 1])\n",
    "\n",
    "# Parametry iteracji\n",
    "max_iterations = 1000\n",
    "tolerance = 1e-6\n",
    "\n",
    "# Iteracyjne obliczanie dominującej wartości własnej i wektora własnego\n",
    "for _ in range(max_iterations):\n",
    "    # Mnożenie macierzy przez wektor\n",
    "    y = np.dot(A, x)\n",
    "\n",
    "    # Normalizacja wektora\n",
    "    x_new = y / np.linalg.norm(y)\n",
    "\n",
    "    # Sprawdzenie warunku zbieżności\n",
    "    if np.linalg.norm(x - x_new) < tolerance:\n",
    "        break\n",
    "\n",
    "    # Aktualizacja wektora\n",
    "    x = x_new\n",
    "\n",
    "# Obliczenie dominującej wartości własnej\n",
    "eigenvalue = np.dot(np.dot(A, x), x) / np.dot(x, x)\n",
    "\n",
    "# Wyświetlenie wyników\n",
    "print(\"Dominująca wartość własna:\")\n",
    "print(eigenvalue)\n",
    "print(\"Dominujący wektor własny:\")\n",
    "print(x)\n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "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.10.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}