Lab-11 gaussian,jacobi,sediel.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# **Lab 11**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 59,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Libraries\n",
    "from numpy import array, zeros, diag, diagflat, dot, tril , allclose\n",
    "from numpy.linalg import norm , inv\n",
    "from math import sqrt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 60,
   "metadata": {},
   "outputs": [],
   "source": [
    "#general jacobi method\n",
    "#Ax=b\n",
    "#this requires A and b as seperate inputs\n",
    "def Jacobi(A,b,x,tol):\n",
    "    dummy_x=x\n",
    "    #1st run of Jacobi\n",
    "    D = diag(A)\n",
    "    R = A - diagflat(D)\n",
    "    x= (b - dot(R,x))/D\n",
    "    temp_array=x-dummy_x\n",
    "    #loop to continue\n",
    "    while (norm(temp_array)/norm(x)) > tol:\n",
    "        x= (b - dot(R,x))/D\n",
    "        temp_array=x-dummy_x\n",
    "        dummy_x=x\n",
    "    return x"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 61,
   "metadata": {},
   "outputs": [],
   "source": [
    "#general jacobi without libraries\n",
    "#solving for a system of linear equations\n",
    "def norm(x):\n",
    "    return max(x)\n",
    "def subtract_list_elem(a,b):\n",
    "    ans=[]\n",
    "    for i in range(len(a)):\n",
    "        ans.append(a[i]-b[i])\n",
    "    return ans\n",
    "def summation(A,i,x):\n",
    "    final_sum=0\n",
    "    for j in range(len(A[i])):\n",
    "        if j!=i:\n",
    "            final_sum+=x[j]*A[i][j]\n",
    "    return final_sum\n",
    "\n",
    "def Jacobi_alt(A, b, x, tol):\n",
    "    prev_x = x\n",
    "    curr_x = x\n",
    "    while True:\n",
    "        new_x = []\n",
    "        for i in range(len(A)):\n",
    "            x_i = (b[i] - summation(A, i, curr_x)) / A[i][i]\n",
    "            new_x.append(x_i)\n",
    "        if norm(subtract_list_elem(new_x, curr_x)) / norm(new_x) < tol:\n",
    "            return new_x\n",
    "        curr_x = new_x\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 62,
   "metadata": {},
   "outputs": [],
   "source": [
    "#general gauss-seidel method\n",
    "def gauss_seidel(A, b, x0, TOL, max_iterations=1000000):\n",
    "    n = len(A)\n",
    "    x = x0.copy()\n",
    "    for i in range(max_iterations):\n",
    "        x_new = zeros(n)\n",
    "        for j in range(n):\n",
    "            s1 = dot(A[j, :j], x_new[:j])\n",
    "            s2 = dot(A[j, j + 1:], x[j + 1:])\n",
    "            x_new[j] = (b[j] - s1 - s2) / A[j, j]\n",
    "        if allclose(x, x_new, rtol=TOL):\n",
    "            return x_new\n",
    "        x = x_new\n",
    "    return x\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#general gauss seidel method without numpy\n",
    "#solving for a system of linear equations\n",
    "def norm(x):\n",
    "    return max(x)\n",
    "def subtract_list_elem(a,b):\n",
    "    ans=[]\n",
    "    for i in range(len(a)):\n",
    "        ans.append(a[i]-b[i])\n",
    "    return ans\n",
    "def summation_seidel(A, i, curr_x, new_x):\n",
    "    final_sum = 0\n",
    "    for j in range(0, i):\n",
    "        final_sum += new_x[j] * A[i][j]\n",
    "    for j in range(i+1, len(A[i])):\n",
    "        final_sum += curr_x[j] * A[i][j]\n",
    "    return final_sum\n",
    "\n",
    "def Seidel_alt(A, b, x, tol):\n",
    "    prev_x = x\n",
    "    curr_x = x\n",
    "    while True:\n",
    "        new_x = []\n",
    "        for i in range(len(A)):\n",
    "            x_i = (b[i] - summation_seidel(A, i, curr_x, new_x)) / A[i][i]\n",
    "            new_x.append(x_i)\n",
    "        print(new_x)\n",
    "        if norm(subtract_list_elem(new_x, curr_x)) / norm(new_x) < tol:\n",
    "            return new_x\n",
    "        curr_x = new_x\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 77,
   "metadata": {},
   "outputs": [],
   "source": [
    "#general gaussian elimination method\n",
    "#note that it is a matriz that includes b matrix\n",
    "def gaussian_elimination(a):\n",
    "    m_21=a[1][0]/a[0][0]\n",
    "    m_31=a[2][0]/a[0][0]\n",
    "    m_32=a[2][1]/a[1][1]\n",
    "    #eliminating a_21\n",
    "    temp_arr=[]\n",
    "    for i in a[0]:\n",
    "        temp_arr.append(m_21*i)\n",
    "    for i in range(len(temp_arr)):\n",
    "        a[1][i]=a[1][i]-temp_arr[i]\n",
    "    #eliminating a_31  \n",
    "    temp_arr=[]\n",
    "    for i in a[0]:\n",
    "        temp_arr.append(m_31*i)\n",
    "    for i in range(len(temp_arr)):\n",
    "        a[2][i]=a[2][i]-temp_arr[i]\n",
    "                    \n",
    "    #eliminating a_32\n",
    "    temp_arr=[]\n",
    "    for i in a[1]:\n",
    "        temp_arr.append(m_32*i)\n",
    "    for i in range(len(temp_arr)):\n",
    "        a[2][i]=a[2][i]-temp_arr[i]\n",
    "\n",
    "    #finding solution\n",
    "    x_3=a[2][3]/a[2][2]\n",
    "    x_2=(a[1][3]-a[1][2]*x_3)/a[1][1]\n",
    "    x_1=(a[0][3]-a[0][2]*x_3-a[0][1]*x_2)/a[0][0]\n",
    "\n",
    "    return (round(x_1,2),round(x_2,2),round(x_3,2))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q1.**Use Gaussian elimination with backward substitution with tolerance $10^{-2}$\n",
    "to solve the following linear system of equations.\n",
    "\n",
    "$4x_1 − x_2 + x_3 = 8$,\n",
    "\n",
    "$2x_1 + 5x_2 + 2x_3 = 3$,\n",
    "\n",
    "$x_1 + 2x_2 + 4x_3 = 11$.\n",
    "\n",
    "The exact solution of the system is $x_1 = 1, x_2 = −1, x_3 = 3$\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "matrix=[[4,-1,1,8],[2,5,2,3],[1,2,4,11]]\n",
    "print(gaussian_elimination(matrix))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q2.**\n",
    "The following linear system\n",
    "\n",
    "$10x_1 - x_2 +2x_3 = 6$\n",
    "\n",
    "$-x_1 + 11x_2 -x_3 +3x_4=25$\n",
    "\n",
    "$2x_1-x_2 +10x_3-x_4= -11$\n",
    "\n",
    "$3x_2 -x_3 +8x_4=15$\n",
    "\n",
    "has the unique solution  $x = (1, 2, −1, 1)^T$. Use Gauss Jacobi’s iterative technique to find the approximations  $x(k)$ to x with $x_0 = (0, 0, 0, 0)^T$ until TOL = $10^{-3}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A=array([[10.,-1.,2.,0.],[-1.,11.,-1.,3.],[2.,-1.,10.,-1.],[0.,3.,-1.,8.]])\n",
    "B=array([6.,25.,-11.,15.])\n",
    "x=array([0.,0.,0.,0.])\n",
    "tol=0.001\n",
    "sol=Jacobi_alt(A,B,x,tol)\n",
    "print(\"The solution of the system is: \",sol)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q3.**\n",
    "Solve Problem 2 by Gauss Seidel iterative Technique.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A=array([[10.,-1.,2.,0.],[-1.,11.,-1.,3.],[2.,-1.,10.,-1.],[0.,3.,-1.,8.]])\n",
    "B=array([6.,25.,-11.,15.])\n",
    "x=array([0.,0.,0.,0.])\n",
    "tol=0.001\n",
    "sol=Seidel_alt(A,B,x,tol)\n",
    "print(\"The solution of the system is: \",sol)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Q4.** Use Gauss-Jacobi Iterations to attempt solving the linear system\n",
    "\n",
    "$x_1 + 2x_2 + 3x_3 = 5$\n",
    "\n",
    "$2x_1 − x_2 + 2x_3 = 1$\n",
    "\n",
    "$3x_1 + x_2 − 2x_3 = −1$\n",
    "\n",
    "(Whether the method converges?)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A=([[1.,2.,3.],[2.,-1.,2.],[3.,1.,-2.]])\n",
    "B=array([5.,1.,-1.])\n",
    "x=array([0.,0.,0.])\n",
    "tol_1=0.1\n",
    "tol_2=0.0001\n",
    "tol_3=0.00001\n",
    "sol_1=Jacobi_alt(A,B,x,tol_1)\n",
    "sol_2=Jacobi_alt(A,B,x,tol_2)\n",
    "sol_3=Jacobi_alt(A,B,x,tol_3)\n",
    "print(\"The solution of the system is: \",sol_1,\" with tolerance: \",tol_1)\n",
    "print(\"The solution of the system is: \",sol_2,\" with tolerance: \",tol_2)\n",
    "print(\"The solution of the system is: \",sol_3,\" with tolerance: \",tol_3)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The method does not converge.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "** Q5.** Use Gauss-Seidel Iterations to attempt solving the linear system\n",
    "\n",
    "$2x_1 + 8x_2 + 3x_3 + x_4 = −2$\n",
    "\n",
    "$2x_2 − x_3 + 4x_4 = 4$\n",
    "\n",
    "$7x_1 − 2x_2 + x_3 + 2x_4 = 3$\n",
    "\n",
    "$−x_1 + 5x_3 + 2x_2 = 5$.\n",
    "\n",
    "(Whether the method converges?)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A=[[2,8,3,1],[0,2,-1,4],[7,-2,1,2],[-1,0,5,2]]\n",
    "B=[-2,4,3,5]\n",
    "x=[0,0,0,0]\n",
    "tol_1=0.1\n",
    "tol_2=0.0001\n",
    "tol_3=0.00001\n",
    "sol_1=Seidel_alt(A,B,x,tol_1)\n",
    "sol_2=Seidel_alt(A,B,x,tol_2)\n",
    "sol_3=Seidel_alt(A,B,x,tol_3)\n",
    "print(\"The solution of the system is: \",sol_1,\" with tolerance: \",tol_1)\n",
    "print(\"The solution of the system is: \",sol_2,\" with tolerance: \",tol_2)\n",
    "print(\"The solution of the system is: \",sol_3,\" with tolerance: \",tol_3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As we can see from the matrix, the method does not converge."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.10.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}