sieve_shortened.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Basic working code"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import math\n",
    "import numpy as np\n",
    "import sympy\n",
    "import time "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calc_time(start, end):\n",
    "    seconds = end - start\n",
    "    mins = seconds / 60\n",
    "    print('Time:', str(seconds), '    ', str(mins))\n",
    "    return seconds, mins"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 21,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "3276.67097683093\n"
     ]
    }
   ],
   "source": [
    "## I have kept the generator cell, but I generate more numbers in the next one. The prime number is set in this cell.\n",
    "## bound is commented out, as it is already set, and the array of primes is already created\n",
    "\n",
    "# n = 16921456439215439701\n",
    "# n = 12110279\n",
    "n = 62615533\n",
    "the_number = n\n",
    "s_const = int(math.sqrt(n)) + 1\n",
    "A_const = (math.sqrt(2)-1)*math.sqrt(n)-1      \n",
    "\n",
    "print(A_const)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "def f(x):\n",
    "    return (x+s_const)**2 - the_number"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Final optimized version \n",
    "def get_primes_less_thanB(bound):                      # eratosthenes way\n",
    "    check_until = int(bound**0.5)\n",
    "    unmarked = np.array([x for x in range(2, bound)])\n",
    "    q = 2\n",
    "    i = 0\n",
    "    while True:\n",
    "        u = unmarked[i]\n",
    "        if u> check_until:               # only check numbers up to sqrt(bound)\n",
    "            break\n",
    "        u_idx = np.where(unmarked == u)[0][0]\n",
    "        max_idx = bound//u\n",
    "        listof_multiples = [u*i for i in range(q, max_idx+1)]\n",
    "        u_multiples = np.array(listof_multiples)\n",
    "        unmarked = np.setdiff1d(unmarked, u_multiples)   # numpy set diff without changing order of elements in list\n",
    "        q = u\n",
    "        i +=1\n",
    "    # print(\"final unmarked elements (for eye balling!): \", unmarked)\n",
    "    print(\"Number of elements less than bound = \", bound, \" is: \", len(unmarked))\n",
    "    return (unmarked)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of elements less than bound =  10000  is:  1229\n",
      "Time: 0.02193737030029297      0.00036562283833821614\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(0.02193737030029297, 0.00036562283833821614)"
      ]
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "start = time.time()\n",
    "bound = 10000\n",
    "primes_less_thanB = get_primes_less_thanB(bound)\n",
    "end = time.time()\n",
    "\n",
    "calc_time(start, end)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Correct optimized version\n",
    "def check_Bsmooth(primes, num):\n",
    "    powers_of_primes= []\n",
    "    num_remaining= num                          # keep dividing by primes in B smooth primes, until == 1\n",
    "    for prime in primes: \n",
    "        #print(\"current prime: \", prime)\n",
    "        i = 0                                   # i = power of prime\n",
    "        if (num_remaining%prime == 0):\n",
    "            while (num_remaining%prime ==0):   #keep dividing for each power of the prime\n",
    "                quotient = num_remaining/prime\n",
    "                num_remaining = quotient\n",
    "                i +=1\n",
    "        powers_of_primes.append(i)\n",
    "        #print(\"num remaining: \", num_remaining)\n",
    "        if num_remaining == 1:     # ---------------(**) INDENT PROPERLY!\n",
    "            diff = len(primes)-len(powers_of_primes)\n",
    "            if diff>0:\n",
    "                zeros = [0 for k in range(diff)]\n",
    "                return powers_of_primes + zeros\n",
    "            return powers_of_primes\n",
    "    #print(\"Number is not B smooth\")\n",
    "    return []"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Still more checks, on B-smooth, on Sieve; "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of elements less than bound =  100  is:  25\n",
      "Time: 0.0      0.0\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(0.0, 0.0)"
      ]
     },
     "execution_count": 23,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "# Was checking check_Bsmooth function\n",
    "num = 100\n",
    "bound = 100\n",
    "primes = np.array(get_primes_less_thanB(bound))\n",
    "\n",
    "start = time.time()\n",
    "result = check_Bsmooth(primes, num)\n",
    "end = time.time()\n",
    "\n",
    "calc_time(start, end)\n",
    "#print(results)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {},
   "outputs": [],
   "source": [
    "def sieve(A_const, bounded_primes):\n",
    "    B_smooth_squares = []\n",
    "    powers_of_primes = []\n",
    "    for i in range(0, int(A_const)):\n",
    "        k = f(i)\n",
    "#         if i%100 ==0:\n",
    "#             print(\"in this iteration i=\", i, \" candidate k is: \", k)\n",
    "        exponents = check_Bsmooth(bounded_primes, k)\n",
    "        if len(exponents)>0:                # same as if exponents != []\n",
    "            B_smooth_squares.append(k)\n",
    "            powers_of_primes.append(exponents)\n",
    "            #print(\"found a B-smooth number(!) = \", k)\n",
    "            #print(\"with prime powers: \", exponents)\n",
    "        else:\n",
    "            continue\n",
    "    return (B_smooth_squares, powers_of_primes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 25,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Time: 0.06623697280883789      0.0011039495468139649\n",
      "[121, 14045, 153725, 174746, 237917, 613210, 763465, 799370, 1197625, 1743005, 1772797, 2535650, 2658370, 3279650, 3342482, 3690346, 3985865, 4066205, 7674425, 8428745, 8884445, 9280346, 9298850, 10234250, 11519042, 11782265]\n"
     ]
    }
   ],
   "source": [
    "# primes defined above\n",
    "start = time.time()\n",
    "results = sieve(A_const, primes)\n",
    "end = time.time()\n",
    "\n",
    "calc_time(start, end)\n",
    "print(results[0])\n",
    "#print(results[1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of B smooth candidates with bound:  100  are  26\n"
     ]
    }
   ],
   "source": [
    "print(\"Number of B smooth candidates with bound: \", bound, \" are \", len(results[0]))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### More tests, ignore to run main prog; start after Markdown cell again"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of elements less than bound =  120  is:  30\n"
     ]
    }
   ],
   "source": [
    "bound = 120\n",
    "primes = np.array(get_primes_less_thanB(bound))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [],
   "source": [
    "start = time.time()\n",
    "result = sieve(A_const, primes)\n",
    "end = time.time()\n",
    "\n",
    "calc_time(start, end)\n",
    "print(\"Number of B smooth candidates with bound: \", bound, \" are \", len(result[0]))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "B_smooth_squares = results[0]\n",
    "# print(B_smooth_squares)\n",
    "powers_of_primes = results[1]\n",
    "# print(powers_of_primes)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "num rows:  53 num columns  46\n",
      "num rows:  46 num columns  53\n"
     ]
    }
   ],
   "source": [
    "exp_array = (np.array(powers_of_primes))%2\n",
    "rownum = len(exp_array)\n",
    "colnum = len(exp_array[0])\n",
    "print(\"num rows: \", rownum, \"num columns \", colnum)\n",
    "\n",
    "exp_array = np.transpose(exp_array)\n",
    "rownum = len(exp_array)\n",
    "colnum = len(exp_array[0])\n",
    "print(\"num rows: \", rownum, \"num columns \", colnum)\n",
    "# print(exp_array)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 37,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Just checking how long this alternative you had takes; seems comparable in terms of time\n",
    "# But this method actually doesn't find ANY candidates k, because the check at (**) means you check \n",
    "# whether EVERY one of the primes in'B smooth primes' divide the candidate k! But that's too strict a requirement!\n",
    "# It also does not check if one prime repeatedly divides k - we need to do that; e.g. 11 goes into 2658370 3 times!\n",
    "\n",
    "def sieve2(A_const, primes):\n",
    "    sieve_array = []\n",
    "\n",
    "    for i in range(0, int(A_const)):\n",
    "        k = f(i)\n",
    "        #print(k)\n",
    "        trigger = 0\n",
    "        flag_cont = 0\n",
    "\n",
    "        for j in range(0, len(primes)):\n",
    "            if(k%primes[j]!=0):   # if any prime does not divide k ----------------(**)\n",
    "                flag_cont = 1     # then set flag to 1\n",
    "\n",
    "        if flag_cont == 1:\n",
    "            flag_cont = 0\n",
    "            continue\n",
    "        else:\n",
    "            sieve_array.append(k)\n",
    "\n",
    "    print(len(sieve_array))\n",
    "    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 40,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n",
      "Time: 0.05785202980041504      0.000964200496673584\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "(0.05785202980041504, 0.000964200496673584)"
      ]
     },
     "execution_count": 40,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "start = time.time()\n",
    "result = sieve2(A_const, primes)\n",
    "end = time.time()\n",
    "calc_time(start, end)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Start here (but load sieve, checkBsmooth, primes less than B smooth)\n",
    "### Checking whether we still get x is  3480 , y is  11.0 GIVEN the new mod2 matrix is DIFFERENT!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Number of elements less than bound =  40  is:  12\n"
     ]
    }
   ],
   "source": [
    "# A first simulation\n",
    "bound = 40\n",
    "primes = np.array(get_primes_less_thanB(bound))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "iterations,  11368\n",
      "B smooth numbers found:  [36, 937251, 22521996, 34781628, 119350, 455455, 859375, 1691690, 41204735, 144, 3749004, 35360, 16344835, 27287936, 47580416, 58679296] total:  16\n",
      "Time: 0.17159318923950195      0.0028598864873250326\n",
      "[[2, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 5, 0, 1, 0, 0, 0, 1, 0, 1, 0, 0], [2, 3, 0, 1, 0, 0, 0, 0, 0, 0, 3, 0], [2, 1, 0, 1, 0, 0, 0, 2, 0, 0, 1, 1], [1, 0, 2, 1, 1, 0, 0, 0, 0, 0, 1, 0], [0, 0, 1, 2, 1, 2, 0, 0, 0, 0, 0, 0], [0, 0, 7, 0, 1, 0, 0, 0, 0, 0, 0, 0], [1, 0, 1, 1, 1, 3, 0, 0, 0, 0, 0, 0], [0, 0, 1, 0, 2, 3, 0, 0, 0, 0, 1, 0], [4, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [2, 5, 0, 1, 0, 0, 0, 1, 0, 1, 0, 0], [5, 0, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0], [0, 0, 1, 0, 0, 2, 0, 0, 1, 2, 0, 0], [7, 0, 0, 0, 0, 1, 0, 0, 2, 0, 1, 0], [8, 0, 0, 0, 0, 1, 1, 0, 0, 2, 0, 0], [13, 0, 0, 0, 0, 1, 0, 1, 0, 1, 0, 0]]\n"
     ]
    }
   ],
   "source": [
    "## this is the generator cell, it works for the small tests I am doing, but we will have to update it later\n",
    "nums = []\n",
    "results=[]\n",
    "size_of_bound = len(primes)\n",
    "\n",
    "numbers = []\n",
    "count_found = 0\n",
    "begin = int(math.sqrt(n)) + 1\n",
    "i = 0\n",
    "start = time.time()\n",
    "while(True):\n",
    "    temp = check_Bsmooth(primes, (begin**2)%n)\n",
    "    if temp!=[]:\n",
    "        #print(\"exponents \", temp)\n",
    "        exponents = np.array(temp)\n",
    "        result = np.power(primes, exponents)\n",
    "        the_num = np.prod(result)\n",
    "        #print(\"the number itself: \", the_num)\n",
    "        numbers.append(the_num)\n",
    "        \n",
    "        nums.append(begin)\n",
    "        results.append(temp)\n",
    "        count_found = count_found +1\n",
    "            \n",
    "    if count_found> size_of_bound +3:\n",
    "        break\n",
    "    begin +=1\n",
    "    i +=1\n",
    "end = time.time()\n",
    "print(\"iterations, \", i)\n",
    "print(\"B smooth numbers found: \", numbers, \"total: \", len(numbers))\n",
    "#print(\"Should equal nums: \", len(nums))\n",
    "calc_time(start, end)\n",
    "print (results)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0 0 0 0 1 0 0 1 0 0 0 1 0 1 0 1]\n",
      " [0 1 1 1 0 0 0 0 0 0 1 0 0 0 0 0]\n",
      " [0 0 0 0 0 1 1 1 1 0 0 1 1 0 0 0]\n",
      " [0 1 1 1 1 0 0 1 0 0 1 0 0 0 0 0]\n",
      " [0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 0]\n",
      " [0 0 0 0 0 0 0 1 1 0 0 1 0 1 1 1]\n",
      " [0 0 0 0 0 0 0 0 0 0 0 1 0 0 1 0]\n",
      " [0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1]\n",
      " [0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0]\n",
      " [0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1]\n",
      " [0 0 1 1 1 0 0 0 1 0 0 0 0 1 0 0]\n",
      " [0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0]]\n"
     ]
    }
   ],
   "source": [
    "## store_pows is the original matrix of powers, which we use later in the computation of y\n",
    "# Got rid of the 'matrix' calculation from 'results', since results already has no empty list\n",
    "# produces same array; checked using 'np.array_equal(matrix, results)'\n",
    "\n",
    "np_array = np.array(results)\n",
    "store_pows = np_array.copy()\n",
    "\n",
    "np_array = np.transpose(np_array)\n",
    "print(np_array %2)\n",
    "# print(store_pows)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[0 0 0 0 1 0 0 1 0 0 0 1 0 1 0 1 0]\n",
      " [0 1 1 1 0 0 0 0 0 0 1 0 0 0 0 0 0]\n",
      " [0 0 0 0 0 1 1 1 1 0 0 1 1 0 0 0 0]\n",
      " [0 1 1 1 1 0 0 1 0 0 1 0 0 0 0 0 0]\n",
      " [0 0 0 0 1 1 1 1 0 0 0 0 0 0 0 0 0]\n",
      " [0 0 0 0 0 0 0 1 1 0 0 1 0 1 1 1 0]\n",
      " [0 0 0 0 0 0 0 0 0 0 0 1 0 0 1 0 0]\n",
      " [0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0]\n",
      " [0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0]\n",
      " [0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0]\n",
      " [0 0 1 1 1 0 0 0 1 0 0 0 0 1 0 0 0]\n",
      " [0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0]]\n",
      "12\n"
     ]
    }
   ],
   "source": [
    "a = (np_array%2)\n",
    "n = len(a)            # row number + number of equations we have\n",
    "m = len(a[0])        # column number + number of primes (in factor base)\n",
    "\n",
    "b =  np.array([[0 for i in range(n)]])     # all zeros\n",
    "a = np.concatenate((a, b.T), axis=1)\n",
    "print(a)\n",
    "print(len(a))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "# changed it to remove all zero rows from A\n",
    "def drop_zero_rows(A): \n",
    "    return A[~np.all(A == 0, axis=1)]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# I got rid of the upp_triangular function since it wasn't necessary to use\n",
    "# aug_matrix = upp_triangular(row_num, col_num,  aug_matrix)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Reduced row echelon form is: \n",
      "[[0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 1 0]\n",
      " [0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 -1 0]\n",
      " [0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0]\n",
      " [0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0]\n",
      " [0 0 0 0 0 1 1 0 0 0 0 0 0 0 0 0 0]\n",
      " [0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0]\n",
      " [0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0]\n",
      " [0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0]\n",
      " [0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0]\n",
      " [0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 1 0]\n",
      " [0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0]]\n",
      "12\n"
     ]
    }
   ],
   "source": [
    "from sympy import *\n",
    "\n",
    "# aug_matrix = a\n",
    "# row_num = len(a)\n",
    "col_num = len(a[0])\n",
    "\n",
    "L = Matrix(a)\n",
    "# print(L.rref())\n",
    "aug_matrix = np.array(L.rref()[0])\n",
    "aug_matrix  = drop_zero_rows(aug_matrix)\n",
    "\n",
    "print(\"Reduced row echelon form is: \")\n",
    "print(aug_matrix)\n",
    "print(len(a))\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1 0 0 0 0 0 1 0 0 1 1 0 0 0 0 1 1]\n"
     ]
    }
   ],
   "source": [
    "# v: 0s are non-free variables, 1s are fixed variables which will be 1 (free or non-free),\n",
    "# 2s are fixed variables which will be 0 (free or non-free)\n",
    "# - Just changed the inner loop to start from 'i', not 0\n",
    "\n",
    "v = np.ones(col_num, dtype = 'int')\n",
    "\n",
    "row_num = len(aug_matrix)   # the row_num has changed after dropping the 0 rows\n",
    "\n",
    "# find the leading entries, aka the non-free variables, set them to 0 in the v vector\n",
    "for i in range(0, row_num):\n",
    "    check_main = 1\n",
    "    for j in range(i, col_num):\n",
    "        if(aug_matrix[i][j]==1 and check_main ==1):\n",
    "            v[j]=0\n",
    "            break\n",
    "\n",
    "print(v)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1 0 0 0 0 0 2 0 0 2 2 0 0 0 0 2 2]\n"
     ]
    }
   ],
   "source": [
    "# set one of the free variable to 1, all the rest of the free variables to 2.\n",
    "## THIS IS A HUGE PLACE FOR OPTIMIZATION, IF WE PICK THE RIGHT FREE VARIABLE LIFE WILL BE EASIER LATER\n",
    "\n",
    "delete_ones = False\n",
    "for i in range(0, col_num):\n",
    "    if(delete_ones==True and v[i]==1):\n",
    "        v[i]=2\n",
    "    if(v[i]==1):\n",
    "        delete_ones = True\n",
    "\n",
    "print(v)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### I actually don't understand what's happening in the cell below\n",
    "#### Do you think you can: 1. Give the action in this cell a name (for easy reference)\n",
    "#### 2. Explain what it does in a section of our paper?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 99,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "final v vector is  [1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]\n"
     ]
    }
   ],
   "source": [
    "# Optimized with numpy\n",
    "\n",
    "for i in reversed(range(0, row_num)):\n",
    "    sum = np.multiply(v, aug_matrix[i]).sum() # equal to sum = sum + sum + v[j]*aug_matrix[i][j] for all j\n",
    "    sum = sum%2;\n",
    "    zer =(np.argwhere(v == 0).flatten())\n",
    "    if(len(zer)!=0):\n",
    "        index = zer[len(zer)-1]   # last entry\n",
    "        if(sum==0):\n",
    "            v[index] = 2\n",
    "        else:\n",
    "            v[index] = 1\n",
    "        \n",
    "print(\"final v vector is \", v%2)\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      " final v vector is  [1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]\n"
     ]
    }
   ],
   "source": [
    "# calculate the other variables according to the values of the free variables we have:\n",
    "\n",
    "for i in reversed(range(0, row_num)):\n",
    "    sum =0\n",
    "    for j in range(0, col_num):\n",
    "        sum = sum + v[j]*aug_matrix[i][j];\n",
    "    zer =(np.argwhere(v == 0).flatten())     # an\n",
    "    sum = sum%2;\n",
    "    if(len(zer)!=0):\n",
    "        index = zer[len(zer)-1]   # why are we using the length of the array?\n",
    "    \n",
    "        if(sum==0):\n",
    "            v[index] = 2\n",
    "        else:\n",
    "            v[index] = 1\n",
    "        \n",
    "print(\" final v vector is \", v%2)\n",
    "        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 102,
   "metadata": {},
   "outputs": [],
   "source": [
    "# QUESTION IS: Do we need the augmented matrix for the process using which you are doing this?\n",
    "\n",
    "v = v[:-1] "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 112,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "16\n",
      "[3480    0    0    0    0    0    0    0    0    0    0    0    0    0\n",
      "    0    0]\n",
      "3480\n"
     ]
    }
   ],
   "source": [
    "# optimized with numpy\n",
    "v = v%2\n",
    "print(len(v))\n",
    "x = np.multiply(nums, v%2)     # element-wise multiplication of two arrays\n",
    "print(x)\n",
    "x[x==0] = 1\n",
    "x = np.prod(x)      # prod is essentially multiplying elements of an array together\n",
    "print(x)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 113,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[0]\n",
      "[[0 0 0 0 2 0 0 0 0 0 0 0]]\n"
     ]
    }
   ],
   "source": [
    "y = 1\n",
    "y_v  = np.zeros(size_of_bound)\n",
    "indices = np.argwhere(v==1)[0]\n",
    "print(indices)\n",
    "\n",
    "print(store_pows)\n",
    "store_pows = store_pows[indices]\n",
    "#print(store_pows)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 114,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "y_v:  [0. 0. 0. 0. 1. 0. 0. 0. 0. 0. 0. 0.]\n",
      "[ 1.  1.  1.  1. 11.  1.  1.  1.  1.  1.  1.  1.]\n",
      "x is  3480 , y is  11.0\n"
     ]
    }
   ],
   "source": [
    "for row in store_pows:\n",
    "    y_v = np.add(y_v, row)\n",
    "\n",
    "y_v = y_v/2\n",
    "print(\"y_v: \", y_v)\n",
    "\n",
    "z = np.power(primes, y_v)\n",
    "print(z)\n",
    "y = y*z\n",
    "\n",
    "\n",
    "y = np.prod(y)   # no need to check y_v == 0, since anything raised to 0 is 1, and on multiplication, 1 has no effect\n",
    "    \n",
    "print(\"x is \", x, \", y is \", y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "x:  7913\n",
      "y_v:  [2. 2. 0. 0. 0. 0. 0. 0. 0. 0. 0. 0.]\n",
      "x is  7913 , y is  6.0\n"
     ]
    }
   ],
   "source": [
    "x = 1;\n",
    "y = 1;\n",
    "# just for recall: nums is a list of 'begin' values, where, begin = int(math.sqrt(n))+1+i, (i= 0, 1, 2,..)\n",
    "\n",
    "\n",
    "\n",
    "for i in range(0, col_num):\n",
    "    if v[i] ==1:\n",
    "        x = x * nums[i]\n",
    "        #for j in range(0, len(primes)):\n",
    "        #    if(store_pows[i][j]!=0):\n",
    "        #        #print(int(y*(primes[j]**store_pows[i][j])))\n",
    "        #        y = (y* ((primes[j]**store_pows[i][j])))\n",
    "\n",
    "print(\"x: \", x)\n",
    "\n",
    "y_v  = np.zeros(size_of_bound)\n",
    "\n",
    "for i in range(0, col_num):\n",
    "    if v[i]==1:\n",
    "        for j in range(0, (size_of_bound)):\n",
    "            y_v[j] = y_v[j] + store_pows[i][j]\n",
    "\n",
    "print(\"y_v: \", y_v)\n",
    "y_v = y_v/2\n",
    "\n",
    "for i in range(0, size_of_bound):\n",
    "    if(y_v is not 0):\n",
    "        y = y*(primes[i]**y_v[i])\n",
    "        \n",
    "\n",
    "print(\"x is \", x, \", y is \", y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[[ 2  2  0  0  0  0  0  0  0  0  0  0]\n",
      " [ 0  5  0  1  0  0  0  1  0  1  0  0]\n",
      " [ 2  3  0  1  0  0  0  0  0  0  3  0]\n",
      " [ 2  1  0  1  0  0  0  2  0  0  1  1]\n",
      " [ 1  0  2  1  1  0  0  0  0  0  1  0]\n",
      " [ 0  0  1  2  1  2  0  0  0  0  0  0]\n",
      " [ 0  0  7  0  1  0  0  0  0  0  0  0]\n",
      " [ 1  0  1  1  1  3  0  0  0  0  0  0]\n",
      " [ 0  0  1  0  2  3  0  0  0  0  1  0]\n",
      " [ 4  2  0  0  0  0  0  0  0  0  0  0]\n",
      " [ 2  5  0  1  0  0  0  1  0  1  0  0]\n",
      " [ 5  0  1  0  0  1  1  0  0  0  0  0]\n",
      " [ 0  0  1  0  0  2  0  0  1  2  0  0]\n",
      " [ 7  0  0  0  0  1  0  0  2  0  1  0]\n",
      " [ 8  0  0  0  0  1  1  0  0  2  0  0]\n",
      " [13  0  0  0  0  1  0  1  0  1  0  0]]\n"
     ]
    }
   ],
   "source": [
    "print(store_pows)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "7907\n",
      "7919\n"
     ]
    }
   ],
   "source": [
    "k = math.gcd(x -int(y), the_number)\n",
    "print(k)\n",
    "print (int(the_number/k)) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 41,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "### works for 5429 = 61*89\n",
    "### works for 8633 = 89 * 97 \n",
    "### works for 8051 = 83 *97 \n",
    "### works for 2257 = 61*37\n",
    "### 1 817 429 = 1579*1151 - works\n",
    "### 12 110 279 = 3469*3461 - works\n",
    "### 62 615 533 = 7907*7919 - works "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "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.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}