Skip to content

Commit

Permalink
feat: add sage implementation for ph23
Browse files Browse the repository at this point in the history
  • Loading branch information
@surfer05
surfer05 committed Dec 24, 2024
1 parent 29e6f71 commit 56e1065
Showing 1 changed file with 325 additions and 0 deletions.
325 changes: 325 additions & 0 deletions ph23/ph23.ipynb
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,325 @@
{
"cells": [
{
"cell_type": "code",
"execution_count": 50,
"id": "b0286ed1",
"metadata": {},
"outputs": [
{
"data": {
"text/html": [
"<html>\\(\\displaystyle -a_{0} {\\left(u_{0} - 1\\right)} {\\left(u_{1} - 1\\right)} {\\left(u_{2} - 1\\right)} + a_{4} u_{0} {\\left(u_{1} - 1\\right)} {\\left(u_{2} - 1\\right)} + a_{2} {\\left(u_{0} - 1\\right)} u_{1} {\\left(u_{2} - 1\\right)} - a_{6} u_{0} u_{1} {\\left(u_{2} - 1\\right)} + a_{1} {\\left(u_{0} - 1\\right)} {\\left(u_{1} - 1\\right)} u_{2} - a_{5} u_{0} {\\left(u_{1} - 1\\right)} u_{2} - a_{3} {\\left(u_{0} - 1\\right)} u_{1} u_{2} + a_{7} u_{0} u_{1} u_{2}\\)</html>"
],
"text/latex": [
"$\\displaystyle -a_{0} {\\left(u_{0} - 1\\right)} {\\left(u_{1} - 1\\right)} {\\left(u_{2} - 1\\right)} + a_{4} u_{0} {\\left(u_{1} - 1\\right)} {\\left(u_{2} - 1\\right)} + a_{2} {\\left(u_{0} - 1\\right)} u_{1} {\\left(u_{2} - 1\\right)} - a_{6} u_{0} u_{1} {\\left(u_{2} - 1\\right)} + a_{1} {\\left(u_{0} - 1\\right)} {\\left(u_{1} - 1\\right)} u_{2} - a_{5} u_{0} {\\left(u_{1} - 1\\right)} u_{2} - a_{3} {\\left(u_{0} - 1\\right)} u_{1} u_{2} + a_{7} u_{0} u_{1} u_{2}$"
],
"text/plain": [
"-a0*(u0 - 1)*(u1 - 1)*(u2 - 1) + a4*u0*(u1 - 1)*(u2 - 1) + a2*(u0 - 1)*u1*(u2 - 1) - a6*u0*u1*(u2 - 1) + a1*(u0 - 1)*(u1 - 1)*u2 - a5*u0*(u1 - 1)*u2 - a3*(u0 - 1)*u1*u2 + a7*u0*u1*u2"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<html>\\(\\displaystyle \\left[-{\\left(u_{0} - 1\\right)} {\\left(u_{1} - 1\\right)} {\\left(u_{2} - 1\\right)}, {\\left(u_{0} - 1\\right)} {\\left(u_{1} - 1\\right)} u_{2}, {\\left(u_{0} - 1\\right)} u_{1} {\\left(u_{2} - 1\\right)}, -{\\left(u_{0} - 1\\right)} u_{1} u_{2}, u_{0} {\\left(u_{1} - 1\\right)} {\\left(u_{2} - 1\\right)}, -u_{0} {\\left(u_{1} - 1\\right)} u_{2}, -u_{0} u_{1} {\\left(u_{2} - 1\\right)}, u_{0} u_{1} u_{2}\\right]\\)</html>"
],
"text/latex": [
"$\\displaystyle \\left[-{\\left(u_{0} - 1\\right)} {\\left(u_{1} - 1\\right)} {\\left(u_{2} - 1\\right)}, {\\left(u_{0} - 1\\right)} {\\left(u_{1} - 1\\right)} u_{2}, {\\left(u_{0} - 1\\right)} u_{1} {\\left(u_{2} - 1\\right)}, -{\\left(u_{0} - 1\\right)} u_{1} u_{2}, u_{0} {\\left(u_{1} - 1\\right)} {\\left(u_{2} - 1\\right)}, -u_{0} {\\left(u_{1} - 1\\right)} u_{2}, -u_{0} u_{1} {\\left(u_{2} - 1\\right)}, u_{0} u_{1} u_{2}\\right]$"
],
"text/plain": [
"[-(u0 - 1)*(u1 - 1)*(u2 - 1),\n",
" (u0 - 1)*(u1 - 1)*u2,\n",
" (u0 - 1)*u1*(u2 - 1),\n",
" -(u0 - 1)*u1*u2,\n",
" u0*(u1 - 1)*(u2 - 1),\n",
" -u0*(u1 - 1)*u2,\n",
" -u0*u1*(u2 - 1),\n",
" u0*u1*u2]"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<html>\\(\\displaystyle \\left[14 X^{7} + 11 X^{6} + X^{5} + 5 X^{4} + 3 X^{3} + 13 X^{2} + 10 X + 12, X^{7} + 4 X^{6} + 4 X^{5} + 5 X^{4} + 10 X^{3} + 4 X^{2} + 13 X + 10, 9 X^{7} + 5 X^{6} + 7 X^{5} + 6 X^{3} + 7 X^{2} + 4 X + 13, X^{7} + 11 X^{6} + 10 X^{4} + 10 X^{3} + 6 X^{2} + 10 X + 3, 5 X^{7} + 2 X^{6} + 3 X^{5} + 4 X^{4} + 10 X^{3} + 5 X + 5, 3 X^{7} + X^{6} + 15 X^{5} + 3 X^{4} + 7 X^{2} + 4 X + 1, 14 X^{7} + 3 X^{6} + X^{5} + 2 X^{4} + 11 X^{3} + 5 X^{2} + 4 X + 11, 4 X^{7} + 14 X^{6} + 3 X^{5} + 5 X^{4} + X^{3} + 9 X^{2} + X + 14\\right]\\)</html>"
],
"text/latex": [
"$\\displaystyle \\left[14 X^{7} + 11 X^{6} + X^{5} + 5 X^{4} + 3 X^{3} + 13 X^{2} + 10 X + 12, X^{7} + 4 X^{6} + 4 X^{5} + 5 X^{4} + 10 X^{3} + 4 X^{2} + 13 X + 10, 9 X^{7} + 5 X^{6} + 7 X^{5} + 6 X^{3} + 7 X^{2} + 4 X + 13, X^{7} + 11 X^{6} + 10 X^{4} + 10 X^{3} + 6 X^{2} + 10 X + 3, 5 X^{7} + 2 X^{6} + 3 X^{5} + 4 X^{4} + 10 X^{3} + 5 X + 5, 3 X^{7} + X^{6} + 15 X^{5} + 3 X^{4} + 7 X^{2} + 4 X + 1, 14 X^{7} + 3 X^{6} + X^{5} + 2 X^{4} + 11 X^{3} + 5 X^{2} + 4 X + 11, 4 X^{7} + 14 X^{6} + 3 X^{5} + 5 X^{4} + X^{3} + 9 X^{2} + X + 14\\right]$"
],
"text/plain": [
"[14*X^7 + 11*X^6 + X^5 + 5*X^4 + 3*X^3 + 13*X^2 + 10*X + 12,\n",
" X^7 + 4*X^6 + 4*X^5 + 5*X^4 + 10*X^3 + 4*X^2 + 13*X + 10,\n",
" 9*X^7 + 5*X^6 + 7*X^5 + 6*X^3 + 7*X^2 + 4*X + 13,\n",
" X^7 + 11*X^6 + 10*X^4 + 10*X^3 + 6*X^2 + 10*X + 3,\n",
" 5*X^7 + 2*X^6 + 3*X^5 + 4*X^4 + 10*X^3 + 5*X + 5,\n",
" 3*X^7 + X^6 + 15*X^5 + 3*X^4 + 7*X^2 + 4*X + 1,\n",
" 14*X^7 + 3*X^6 + X^5 + 2*X^4 + 11*X^3 + 5*X^2 + 4*X + 11,\n",
" 4*X^7 + 14*X^6 + 3*X^5 + 5*X^4 + X^3 + 9*X^2 + X + 14]"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<html>\\(\\displaystyle \\left[1, 2, 3, 4, 5, 6, 7, 8\\right]\\)</html>"
],
"text/latex": [
"$\\displaystyle \\left[1, 2, 3, 4, 5, 6, 7, 8\\right]$"
],
"text/plain": [
"[1, 2, 3, 4, 5, 6, 7, 8]"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<html>\\(\\displaystyle 16 X^{7} + 6 X^{6} + 13 X^{5} + 11 X^{4} + 12 X^{3} + 11 X^{2} + 3 X + 14\\)</html>"
],
"text/latex": [
"$\\displaystyle 16 X^{7} + 6 X^{6} + 13 X^{5} + 11 X^{4} + 12 X^{3} + 11 X^{2} + 3 X + 14$"
],
"text/plain": [
"16*X^7 + 6*X^6 + 13*X^5 + 11*X^4 + 12*X^3 + 11*X^2 + 3*X + 14"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<html>\\(\\displaystyle \\left[\\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}\\right]\\)</html>"
],
"text/latex": [
"$\\displaystyle \\left[\\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}, \\mathrm{True}\\right]$"
],
"text/plain": [
"[True, True, True, True, True, True, True, True]"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<html>\\(\\displaystyle X^{7} + 2 X^{6} + 11 X^{5} + 4 X^{4} + 16 X^{3} + 7 X^{2} + 8 X + 13\\)</html>"
],
"text/latex": [
"$\\displaystyle X^{7} + 2 X^{6} + 11 X^{5} + 4 X^{4} + 16 X^{3} + 7 X^{2} + 8 X + 13$"
],
"text/plain": [
"X^7 + 2*X^6 + 11*X^5 + 4*X^4 + 16*X^3 + 7*X^2 + 8*X + 13"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<html>\\(\\displaystyle X^{6} + 15 X^{5} + 2 X^{4} + 13 X^{3} + 15 X^{2} + 15 X + 16\\)</html>"
],
"text/latex": [
"$\\displaystyle X^{6} + 15 X^{5} + 2 X^{4} + 13 X^{3} + 15 X^{2} + 15 X + 16$"
],
"text/plain": [
"X^6 + 15*X^5 + 2*X^4 + 13*X^3 + 15*X^2 + 15*X + 16"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"data": {
"text/html": [
"<html>\\(\\displaystyle X^{4} + 5 X^{3} + 4 X^{2} + 12 X + 1\\)</html>"
],
"text/latex": [
"$\\displaystyle X^{4} + 5 X^{3} + 4 X^{2} + 12 X + 1$"
],
"text/plain": [
"X^4 + 5*X^3 + 4*X^2 + 12*X + 1"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"# Define variables\n",
"n = 3\n",
"N = 2^n\n",
"X = var(['X{}'.format(i) for i in range(n)]) # Define the variables X_0, X-1,.....\n",
"u = var(['u{}'.format(i) for i in range(n)]) # Define the variables u_0, u-1,.....\n",
"\n",
"# Function to get binary representation as a list of bits\n",
"def bits(i,n):\n",
" return list(map(int, format(i,'0{}b'.format(n))))\n",
"\n",
"# Define the eq_tilde function\n",
"def eq_tilde(bits_i, u_vector):\n",
" result=1\n",
" for bit,u in zip(bits_i,u_vector):\n",
" result *= (1-bit)*(1-u) + bit*u\n",
" return result\n",
"\n",
"# Coefficients of the polynomial\n",
"a = [var('a{}'.format(i)) for i in range(N)] # Coefficients a_0, a_1, ...., a_(N-1)\n",
"\n",
"# MLE polynomial\n",
"f_tilde = sum(a[i]*eq_tilde(bits(i,n), u) for i in range(N))\n",
"show(f_tilde)\n",
"\n",
"# Generate all combinations of (1-u[i]) and u[i] based on binary representation\n",
"def generate_c_vector(n, u):\n",
" c_vector = []\n",
" for i in range(2^n): # Loop over all binary numbers from 0 to 2^n - 1\n",
" binary = list(map(int, format(i, f'0{n}b'))) # Binary representation of i\n",
" product = 1\n",
" for j, bit in enumerate(binary):\n",
" if bit == 0:\n",
" product *= (1 - u[j]) # Use (1 - u[j]) for 0\n",
" else:\n",
" product *= u[j] # Use u[j] for 1\n",
" c_vector.append(product)\n",
" return c_vector\n",
"\n",
"# Compute the c vector\n",
"c = generate_c_vector(n, u)\n",
"\n",
"# Display the c vector\n",
"show(c)\n",
"\n",
"# SageMath Implementation for Polynomial Encoding using Lagrange Basis\n",
"\n",
"# Define the finite field and subgroup H\n",
"p = 17 # Example prime (p must be chosen appropriately)\n",
"F = GF(p) # Finite field F_p\n",
"omega = F(3) # Example primitive 8th root of unity in F_p\n",
"H = [omega^i for i in range(8)] # Multiplicative subgroup H of size 8\n",
"\n",
"# Define Lagrange Basis Polynomials\n",
"def lagrange_basis(H, X):\n",
" basis = []\n",
" for i in range(len(H)):\n",
" Li = 1\n",
" for j in range(len(H)):\n",
" if i != j:\n",
" Li *= (X - H[j]) / (H[i] - H[j]) # Lagrange basis polynomial\n",
" basis.append(Li)\n",
" return basis\n",
"\n",
"# Symbolic variable for the polynomial\n",
"X = polygen(F, 'X')\n",
"\n",
"# Compute the Lagrange basis polynomials\n",
"L = lagrange_basis(H, X)\n",
"show(L)\n",
"\n",
"# Vector c as computed earlier\n",
"c = [F(c_val) for c_val in [1, 2, 3, 4, 5, 6, 7, 8]] # Example values for c vector\n",
"show(c)\n",
"# Compute the polynomial encoding c(X)\n",
"c_X = sum(c[i] * L[i] for i in range(len(c)))\n",
"\n",
"# Compute the polynomial encoding c(X) step by step\n",
"# c_X = 0 # Start with an empty polynomial\n",
"# print(\"Building c(X) step by step:\")\n",
"# for i in range(len(c)):\n",
"# show(c[i])\n",
"# show(L[i])\n",
"# term = c[i] * L[i] # Compute the current term\n",
"# show(term)\n",
"# c_X += term # Add the current term to the polynomial\n",
"# show(term)\n",
"# '\\n'\n",
"# # print(f\"Step {i+1}: Adding term {term}\")\n",
"# # print(f\"Partial sum: {c_X}\\n\")\n",
"\n",
"# # Final polynomial\n",
"# print(\"Final polynomial c(X):\")\n",
"# show(c_X)\n",
"\n",
"# Display the resulting polynomial c(X)\n",
"show(c_X)\n",
"\n",
"# Verification: Check that c(omega^i) = c_i\n",
"verification = [c_X(H[i]) == c[i] for i in range(len(H))]\n",
"show(verification)\n",
"\n",
"\n",
"# Define subsets H0, H1, H2 based on the Group Tower relationship\n",
"H0 = [H[0]] # {1}\n",
"H1 = [H[0], H[4]] # {1, ω^4}\n",
"H2 = [H[0], H[2], H[4], H[6]] # {1, ω^2, ω^4, ω^6}\n",
"H3 = H # Full set {1, ω, ω^2, ..., ω^7}\n",
"\n",
"# Define the vanishing polynomials v_H(X) and v_Hi(X)\n",
"X = polygen(F, 'X') # Define X as a polynomial variable\n",
"\n",
"def vanishing_polynomial(domain):\n",
" \"\"\"Compute the vanishing polynomial for a given domain.\"\"\"\n",
" v = 1\n",
" for alpha in domain:\n",
" v *= (X - alpha)\n",
" return v\n",
"\n",
"# Compute vanishing polynomials\n",
"v_H = vanishing_polynomial(H3) # Full vanishing polynomial for H\n",
"v_H0 = vanishing_polynomial(H0)\n",
"v_H1 = vanishing_polynomial(H1)\n",
"v_H2 = vanishing_polynomial(H2)\n",
"\n",
"# Compute the selector polynomials s_0(X), s_1(X), s_2(X)\n",
"s0 = v_H/v_H0\n",
"s1 = v_H/v_H1\n",
"s2 = v_H/v_H2\n",
"\n",
"# Display the selector polynomials\n",
"show(s0)\n",
"show(s1)\n",
"show(s2)\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "65879a63",
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "SageMath 10.0",
"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.9.19"
}
},
"nbformat": 4,
"nbformat_minor": 5
}

0 comments on commit 56e1065

Please sign in to comment.