System of Linear Equations Symbolic Solver.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "def solve_linear_system(F, A, B, n):\n",
    "    solution = []\n",
    "    P = PolynomialRing(F, n, names='x', order='lex')\n",
    "    x = P.gens()\n",
    "    X = matrix([[x[i]] for i in range(n)]); X\n",
    "    equations = A*X-B\n",
    "    I = Ideal([equations[i][0] for i in range(n)])\n",
    "    print(\"The Ideal Corresponding to the equations is I = \", I)\n",
    "    gb = I.groebner_basis()\n",
    "    print(\"Groebner Basis of I is \", gb)\n",
    "\n",
    "    rref = [list(gb[i].dict().items()) for i in range(len(gb))];\n",
    "    \n",
    "    if len(rref) == 1 and sum(rref[0][0][0]) == 0 and rref[0][0][1] == 1:\n",
    "        print(\"No Solution!, Exiting\")\n",
    "        return solution\n",
    "    \n",
    "    is_pivot_variable =  [0 for i in range(n)] # get alll pivot variables\n",
    "    for li in rref:\n",
    "        is_pivot_variable[[i for i, x in enumerate(li[0][0]) if x][0]] = 1\n",
    "        \n",
    "    answer = []  # computes B part of rref(A:B)\n",
    "    for i in range(n):\n",
    "        if is_pivot_variable[i] == 0:\n",
    "            answer.append(var('c'+str(i), domain = CC))\n",
    "        else:\n",
    "            answer.append(rref[i][-1][-1])\n",
    "    \n",
    "    for i in range(len(rref)): # back-substitution\n",
    "        for j in rref[i]:\n",
    "            indices = [k for k, x in enumerate(j[0]) if x]\n",
    "            if len(indices) != 0:\n",
    "                index = indices[0]\n",
    "                if index != i and index < n:\n",
    "                    answer[i] = answer[i] + answer[index]*j[1]\n",
    "                    \n",
    "    for i in range(n):\n",
    "        solution.append(answer[i])\n",
    "    return solution"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def print_solution(solution):\n",
    "    if len(solution) == 0:\n",
    "        print(\"No Solution!\")\n",
    "    for i, sol in enumerate(solution):\n",
    "        print(\"x[\"+str(i)+\"]\", \" = \", sol)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 4"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "F = QQ"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[  1   1   1   1]\n",
       "[  1   2   3   4]\n",
       "[  1   4   9  16]\n",
       "[ -3  -7 -13 -21]"
      ]
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = matrix([[1,1,1,1],[1,2,3,4],[1,4,9,16],[-3,-7,-13,-21]]); A"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[ 1]\n",
       "[ 2]\n",
       "[ 3]\n",
       "[-6]"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "B = matrix([[1],[2],[3],[-6]]); B"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The Ideal Corresponding to the equations is I =  Ideal (x0 + x1 + x2 + x3 - 1, x0 + 2*x1 + 3*x2 + 4*x3 - 2, x0 + 4*x1 + 9*x2 + 16*x3 - 3, -3*x0 - 7*x1 - 13*x2 - 21*x3 + 6) of Multivariate Polynomial Ring in x0, x1, x2, x3 over Rational Field\n",
      "Groebner Basis of I is  [x0 + x3 + 1/2, x1 - 3*x3 - 2, x2 + 3*x3 + 1/2]\n",
      "x[0]  =  c3 + 1/2\n",
      "x[1]  =  -3*c3 - 2\n",
      "x[2]  =  3*c3 + 1/2\n",
      "x[3]  =  c3\n"
     ]
    }
   ],
   "source": [
    "s = solve_linear_system(F, A, B, n)\n",
    "print_solution(s)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [],
   "source": [
    "set_random_seed(1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "n = 5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "F = QQ # Field "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[  0   7   1   0  -1]\n",
       "[  1 -12   0   4   3]\n",
       "[ -1  -2   1  -2  -1]\n",
       "[  0   8   3  -2  -1]\n",
       "[  1   0   2   0  -1]"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "A = random_matrix(ZZ,n,n); A"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "[ 0]\n",
       "[-3]\n",
       "[ 1]\n",
       "[ 1]\n",
       "[ 1]"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "B = random_matrix(ZZ,n,1); B"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "The Ideal Corresponding to the equations is I =  Ideal (7*x1 + x2 - x4, x0 - 12*x1 + 4*x3 + 3*x4 + 3, -x0 - 2*x1 + x2 - 2*x3 - x4 - 1, 8*x1 + 3*x2 - 2*x3 - x4 - 1, x0 + 2*x2 - x4 - 1) of Multivariate Polynomial Ring in x0, x1, x2, x3, x4 over Rational Field\n",
      "Groebner Basis of I is  [x0 - 6/7, x1 + 1/21, x2 + 4/21, x3 + 5/7, x4 + 11/21]\n",
      "x[0]  =  -6/7\n",
      "x[1]  =  1/21\n",
      "x[2]  =  4/21\n",
      "x[3]  =  5/7\n",
      "x[4]  =  11/21\n"
     ]
    }
   ],
   "source": [
    "s = solve_linear_system(F, A, B, n)\n",
    "print_solution(s)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "SageMath 9.3",
   "language": "sage",
   "name": "sagemath"
  },
  "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.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}