diff --git a/sage/.gitignore b/sage/.gitignore new file mode 100644 index 000000000..a2c29b1d0 --- /dev/null +++ b/sage/.gitignore @@ -0,0 +1,313 @@ +# Byte-compiled / optimized / DLL files +__pycache__/ +*.py[cod] +*$py.class + +# C extensions +*.so + +# Distribution / packaging +.Python +build/ +develop-eggs/ +dist/ +downloads/ +eggs/ +.eggs/ +lib/ +lib64/ +parts/ +sdist/ +var/ +wheels/ +share/python-wheels/ +*.egg-info/ +.installed.cfg +*.egg +MANIFEST + +# PyInstaller +# Usually these files are written by a python script from a template +# before PyInstaller builds the exe, so as to inject date/other infos into it. +*.manifest +*.spec + +# Installer logs +pip-log.txt +pip-delete-this-directory.txt + +# Unit test / coverage reports +htmlcov/ +.tox/ +.nox/ +.coverage +.coverage.* +.cache +nosetests.xml +coverage.xml +*.cover +*.py,cover +.hypothesis/ +.pytest_cache/ +cover/ + +# Translations +*.mo +*.pot + +# Django stuff: +*.log +local_settings.py +db.sqlite3 +db.sqlite3-journal + +# Flask stuff: +instance/ +.webassets-cache + +# Scrapy stuff: +.scrapy + +# Sphinx documentation +docs/_build/ + +# PyBuilder +.pybuilder/ +target/ + +# Jupyter Notebook +.ipynb_checkpoints + +# IPython +profile_default/ +ipython_config.py + +# pyenv +# For a library or package, you might want to ignore these files since the code is +# intended to run in multiple environments; otherwise, check them in: +# .python-version + +# pipenv +# According to pypa/pipenv#598, it is recommended to include Pipfile.lock in version control. +# However, in case of collaboration, if having platform-specific dependencies or dependencies +# having no cross-platform support, pipenv may install dependencies that don't work, or not +# install all needed dependencies. +#Pipfile.lock + +# UV +# Similar to Pipfile.lock, it is generally recommended to include uv.lock in version control. +# This is especially recommended for binary packages to ensure reproducibility, and is more +# commonly ignored for libraries. +#uv.lock + +# poetry +# Similar to Pipfile.lock, it is generally recommended to include poetry.lock in version control. +# This is especially recommended for binary packages to ensure reproducibility, and is more +# commonly ignored for libraries. +# https://python-poetry.org/docs/basic-usage/#commit-your-poetrylock-file-to-version-control +#poetry.lock + +# pdm +# Similar to Pipfile.lock, it is generally recommended to include pdm.lock in version control. +#pdm.lock +# pdm stores project-wide configurations in .pdm.toml, but it is recommended to not include it +# in version control. +# https://pdm.fming.dev/latest/usage/project/#working-with-version-control +.pdm.toml +.pdm-python +.pdm-build/ + +# PEP 582; used by e.g. github.com/David-OConnor/pyflow and github.com/pdm-project/pdm +__pypackages__/ + +# Celery stuff +celerybeat-schedule +celerybeat.pid + +# SageMath parsed files +*.sage.py + +# Environments +.env +.venv +env/ +venv/ +ENV/ +env.bak/ +venv.bak/ + +# Spyder project settings +.spyderproject +.spyproject + +# Rope project settings +.ropeproject + +# mkdocs documentation +/site + +# mypy +.mypy_cache/ +.dmypy.json +dmypy.json + +# Pyre type checker +.pyre/ + +# pytype static type analyzer +.pytype/ + +# Cython debug symbols +cython_debug/ + +# PyCharm +# JetBrains specific template is maintained in a separate JetBrains.gitignore that can +# be found at https://github.com/github/gitignore/blob/main/Global/JetBrains.gitignore +# and can be added to the global gitignore or merged into this file. For a more nuclear +# option (not recommended) you can uncomment the following to ignore the entire idea folder. +#.idea/ + +# PyPI configuration file +.pypirc + +*.bundle.* +lib/ +node_modules/ +*.egg-info/ +.ipynb_checkpoints +*.tsbuildinfo + +# Created by https://www.gitignore.io/api/python +# Edit at https://www.gitignore.io/?templates=python + +### Python ### +# Byte-compiled / optimized / DLL files +__pycache__/ +*.py[cod] +*$py.class + +# C extensions +*.so + +# Distribution / packaging +.Python +build/ +develop-eggs/ +dist/ +downloads/ +eggs/ +.eggs/ +lib/ +lib64/ +parts/ +sdist/ +var/ +wheels/ +pip-wheel-metadata/ +share/python-wheels/ +.installed.cfg +*.egg +MANIFEST + +# PyInstaller +# Usually these files are written by a python script from a template +# before PyInstaller builds the exe, so as to inject date/other infos into it. +*.manifest +*.spec + +# Installer logs +pip-log.txt +pip-delete-this-directory.txt + +# Unit test / coverage reports +htmlcov/ +.tox/ +.nox/ +.coverage +.coverage.* +.cache +nosetests.xml +coverage.xml +*.cover +.hypothesis/ +.pytest_cache/ + +# Translations +*.mo +*.pot + +# Scrapy stuff: +.scrapy + +# Sphinx documentation +docs/_build/ + +# PyBuilder +target/ + +# pyenv +.python-version + +# celery beat schedule file +celerybeat-schedule + +# SageMath parsed files +*.sage.py + +# Spyder project settings +.spyderproject +.spyproject + +# Rope project settings +.ropeproject + +# Mr Developer +.mr.developer.cfg +.project +.pydevproject + +# mkdocs documentation +/site + +# mypy +.mypy_cache/ +.dmypy.json +dmypy.json + +# Pyre type checker +.pyre/ + +# OS X stuff +*.DS_Store + +# End of https://www.gitignore.io/api/python + +_temp_extension +junit.xml +[uU]ntitled* +notebook/static/* +!notebook/static/favicons +notebook/labextension +notebook/schemas +docs/source/changelog.md +docs/source/contributing.md + +# playwright +ui-tests/test-results +ui-tests/playwright-report + +# VSCode +.vscode + +# RTC +.jupyter_ystore.db + +# yarn >=2.x local files +.yarn/* +.pnp.* +ui-tests/.yarn/* +ui-tests/.pnp.* + +# generated html +notebook/templates/*.html + diff --git a/sage/fri-and-friends/README.md b/sage/fri-and-friends/README.md new file mode 100644 index 000000000..0f7ce85f6 --- /dev/null +++ b/sage/fri-and-friends/README.md @@ -0,0 +1,9 @@ +# Working with FRI, DEEP-FRI, and WHIR + +The basic IOPP for FRI, DEEP-FRI and WHIR are in the [src](./src/) folder. The Jupter notebooks are meant for using these IOPPs/PCS for experimentation. + +## Sage Resources + +* [installation guide](https://doc.sagemath.org/html/en/installation/index.html). +* Running Sage as [Jupyter Notebook](https://doc.sagemath.org/html/en/installation/launching.html) +* VSCode [Jupter Extension](https://marketplace.visualstudio.com/items?itemName=ms-toolsai.jupyter). (Supports inline python debugging!!!) diff --git a/sage/fri-and-friends/Zero Knowledge for WHIR.md b/sage/fri-and-friends/Zero Knowledge for WHIR.md new file mode 100644 index 000000000..46ffab1d4 --- /dev/null +++ b/sage/fri-and-friends/Zero Knowledge for WHIR.md @@ -0,0 +1,141 @@ +# Zero Knowledge for WHIR + +### Zero-Knowledge WHIR as a PCS + +The following scheme describes a Zero-Knowledge variant of WHIR as a PCS. The original ideas came from other people. The scheme uses a single run of non-ZK WHIR as a subprotocol. The scheme is zero knowledge in the bounded query model. + +Given a witness polynomial $f(x_1,\cdots, x_\mu)$, the protocol works as follows: + +1. Prover ($\mathsf{P}$) masks $f$ by adding a new variable $y$ to create a new randomized *multilinear* polynomial + + $$ + \hat{f}(x_1,\cdots,x_{\mu},y) := f(x_1,\cdots,x_{\mu}) + y\cdot \mathsf{msk}(x_1, \cdots, x_\mu) + $$ + + where $\mathsf{msk}(\vec{x})$ should be random with at least as many non-zero coefficients as the *total* query bound. Furthermore, it must also not evaluate to zero on the NTT evaluation domain. (Looking ahead, the oracle corresponding to $\hat{f}$ will only be queried on NTT evaluation domain and never as part of queries to the Sumcheck protocol.) + +2. $\mathsf{P}$ additionally samples a random $\mu + 1$-variate *multilinear polynomial* $g(x_1,\cdots,x_\mu, y)$. +3. **Commitment**: $\mathsf{P}$ commits to $f$ by sending two FRI oracles $([[\hat{f}]], [[g]])$, to the verifier $\mathsf{V}$. +4. **Opening**: Verifier sends the evaluation point $\vec{a} := (a_1, \cdots, a_\mu)$ to the prover who opens $f$ by computing the following: + - Honest Prover computes $F = f(a_1, \cdots, a_\mu) = \hat{f}(a_1, \cdots, a_\mu, 0)$. + - Honest Prover computes $G = g(a_1, \cdots, a_\mu, 0)$. (Notice $y$ is set to zero.) + - $\mathsf{P} \longrightarrow \mathsf{V}$: Prover sends the pair $(F, G)$ to the verifier. +5. **Proof Verification**: After receiving $(F, G)$, the protocol starts verification as follows: + - $\mathsf{V} \longrightarrow \mathsf{P}$: Verifier sends a random challenge $\rho$ to the prover. Verifier expects $\mathsf{P}$ to prove the claim + + $$ + + \rho\cdot F + G = \rho\cdot f(\vec{a}) + g(\vec{a}, 0) + + $$ + + - Prover and Verifier proceed with non-ZK WHIR opening proof of the previous claim as follows: + - Verifier: + - Creates the **virtual oracle** $L := \rho \cdot [[\hat{f}]] + [[g]]$ and sets it as the first WHIR FRI oracle. + - Verifier creates the WHIR weight polynomial $w(z,x_1,\cdots,x_\mu, y) := z\cdot\mathsf{eq}(\vec{x},y;\; \vec{a},0)$. + - Verifier runs the standard WHIR Sumcheck rounds receiving new folded FRI oracles as usual. If the prover is honest, the folded FRI oracles will consist of the + folding of the polynomial + + $$ + + h(\vec{x}, y) := \rho\cdot \hat{f}(x_1,\cdots,x_\mu, y) + g(x_1,\cdots,x_\mu, y) + + $$ + + And an honest prover's Sumcheck rounds will prove the claim: + + $$ + + \rho\cdot F + G = \rho\cdot \hat{f}(\vec{a},0) + g(\vec{a},0) = \sum_{\vec{x} \in \{0,1\}^{\mu+1}} \left [\rho\cdot\hat{f}(\vec{x}) + g(\vec{x}) \right ]\mathsf{eq}(\vec{x};\; \vec{a}, 0) + + $$ + + Notice that since $f(\vec{a}) = \hat{f}(\vec{a}, 0)$, the Sumcheck opening proof is indeed a proof of evaluation of $f$ and is therefore binding, in spite of the mask! + + - **OOD Queries**: The verifier should make its out-of-domain (OOD) queries as usual. + - **Shift Queries**: At the end of the first Sumcheck round and on receiving the first *folded* oracle corresponding to $h(\vec{x},y)$, the Verifier should compute it's shift-query response using the virtual oracle $L$. For subsequent rounds, the protocol should follow normal WHIR operations. + - **Termination Check**: The verifier should do termination check as normal WHIR terminal check, i.e., locally compute the Sumcheck over final $w(z,\vec{x},y)$, with $z$ replaced by the last round of FRI oracle $h(\vec{r},y)$ now sent in clear. Verifier should also check that the final FRI oracle is not a constant, i.e., it must have a non-zero $y$-coefficient. (Degree of $y$ should also be one.) + + Note that at the $\mu$-th Sumcheck round, an honest prover would have sent + + $$ + + h(\vec{r}, y)\cdot \mathsf{eq}(\vec{r}, y; \vec{a},0) = (\rho\cdot[f(\vec{r}) + y\cdot \mathsf{msk}(\vec{r})] + g(\vec{r},y) )\cdot \mathsf{eq}(\vec{r}, y; \vec{a},0) + + $$ + + and the corresponding FRI oracle, sent in clear to the Verifier will be + + $$ + + \rho\cdot[f(\vec{r}) + y\cdot \mathsf{msk}(\vec{r})] + g(\vec{r},y) + + $$ + + where $\vec{r}$ is the randomness sent during Sumcheck. Since the verifier had set $\vec{w}(z,\vec{x},y) = z\cdot\mathsf{eq}(\vec{x},y;\; \vec{a},0)$ at the start of the protocol, the result of locally computing + + $$ + + \sum_{y\in\{0,1\}} h(\vec{r},y)\cdot\mathsf{eq}(\vec{r}, y; \vec{a},0)) = h(\vec{r}, 0) = \rho\cdot f(\vec{r}) + g(\vec{r},0) + + $$ + + - Honest Prover + - On receiving verifier's challenge $\rho$, the prover starts WHIR PCS evaluation proof for $(\vec{a}, 0)$ on the polynomial $h(\vec{x}, y) = \rho\cdot \hat{f}(x_1,\cdots,x_\mu, y) + g(x_1,\cdots,x_\mu, y)$. However, unlike a standard proof, the Prover should not compute the first proof oracle. + - **OOD Queries**: For an out-of-domain (OOD) query $q_o$, the prover should respond with the evaluation of $h\left(\textsf{pow2}(q_o), q_o^{2^\mu}\right )$ (i.e., without setting y to zero) to the verifier. Notice that out-of-domain queries are meant for $h(\vec{x})$, which by construction is blinded by $g$, and therefore doesn't leak information about $f$. + +### Zero-Knowledge WHIR for $\Sigma$-IOP + +The following Zero-Knowledge variant of WHIR is applicable to $\Sigma$-IOP. + +Given: + +- A private *multilinear* witness polynomial $f(x_1,\cdots,x_\mu)$ +- A public WHIR weight polynomial $w$ of the form $w(z,x_1,\cdots,x_\mu) = z\cdot W(x_1,\cdots,x_\mu)$, and +- A $\Sigma$-IOP claim $F = \sum_{\vec{b}\in \{0,1\}^\mu} f(\vec{b})\cdot W(\vec{b})$ + +The following protocol proves the above claim in zero-knowledge using non-ZK WHIR as a subprotocol. The protocol proceeds as follows: + +1. As in the case of PCS, the Prover ($\mathsf{P}$) masks $f$ by adding a new variable $y$ to create a new randomized multilinear polynomial with same requirements as in the case of PCS: + +$$ +\hat{f}(x_1,\cdots,x_{\mu},y) := f(x_1,\cdots,x_{\mu}) + y\cdot \mathsf{msk}(x_1, \cdots, x_\mu) +$$ + +1. $\mathsf{P}$ additionally samples a random $\mu + 1$-variate multilinear polynomial $g(x_1,\cdots,x_\mu, y)$. +2. **Commitment**: $\mathsf{P}$ commits to the claim $F = \sum_{\vec{b}\in \{0,1\}^\mu} f(\vec{b})\cdot W(\vec{b})$ as follows: + - Prover computes a new weight polynomial + + $$ + \hat{w}(z,\vec{x}, y) := z\cdot W(\vec{x})\cdot \textsf{eq}(0,y) = z\cdot (1-y) \cdot W(\vec{x}) + $$ + + - Prover computes + + $$ + G = \sum_{\vec{b}\in \{0,1\}^{\mu}} g(\vec{b}, 0)\cdot W(\vec{b}) + $$ + + - Prover computes the FRI oracles $[[\hat{f}]]$ and $[[g]]$ + - $\mathsf{P} \longrightarrow \mathsf{V}$: Prover sends the quadruple $(F, G, [[\hat{f}]], [[g]])$ to the verifier. +3. **Proof Verification**: After receiving $(F,G, [[\hat{f}]], [[g]])$, the protocol starts verification as follows: + - $\mathsf{V} \longrightarrow \mathsf{P}$: Verifier sends a random challenge $\rho$ to the prover. Verifier expects $\mathsf{P}$ to prove the claim + + $$ + \rho\cdot F + G = \sum_{\vec{b}\in \{0,1\}^\mu} [\rho\cdot f(\vec{b}) + g(\vec{b},0)]\cdot W(\vec{b}) + $$ + + - Prover and Verifier proceed with non-ZK WHIR opening proof of the previous claim as follows: + - Verifier: + - Creates the **virtual oracle** $L := \rho \cdot [[\hat{f}]] + [[g]]$ and sets it as the first WHIR FRI oracle. + - Verifier creates a new WHIR weight polynomial + + $$ + \hat{w}(z,x_1,\cdots,x_\mu, y) := z\cdot (1-y) \cdot W(x_1,\cdots, x_\mu) + $$ + + + and runs non-ZK WHIR as usual, expecting $\mathsf{P}$ to prove the $\Sigma$-IOP claim on $\rho\cdot F + G$. + + - Rest of the protocol should follow the same algorithm for OOD-queries, Shift-queries, and Termination checks as in the case of WHIR PCS. + - Honest Prover should use the weight polynomial computed during the commitment phase, and follow rest of the protocol as in case of WHIR PCS. \ No newline at end of file diff --git a/sage/fri-and-friends/coding_theory.ipynb b/sage/fri-and-friends/coding_theory.ipynb new file mode 100644 index 000000000..d9ab85213 --- /dev/null +++ b/sage/fri-and-friends/coding_theory.ipynb @@ -0,0 +1,196 @@ +{ + "cells": [ + { + "cell_type": "code", + "execution_count": 3, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [], + "source": [ + "from sage.coding.all import *\n", + "from sage.coding.grs_code import ReedSolomonCode\n", + "from sage.coding.guruswami_sudan.gs_decoder import GRSGuruswamiSudanDecoder\n", + "from sage.rings.polynomial.all import *\n", + "from sage.rings.finite_rings.all import *\n", + "from sage.modules.free_module import VectorSpace\n", + "from sage.matrix.all import *\n", + "from sage.misc.functional import sqrt\n", + "from sage.functions.other import floor\n", + "from copy import deepcopy\n", + "from sage.misc.prandom import random" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Minimum Distance: 769\n", + "Dimension: 256\n", + "Rate: 1/4\n" + ] + } + ], + "source": [ + "P = 15*(2**27) + 1\n", + "Fq = GF(P);\n", + "R = PolynomialRing(Fq, \"X\");\n", + "X = R.gen()\n", + "\n", + "π = Fq.multiplicative_generator() # Multiplicative Generator 31\n", + "\n", + "omega = π**(15*2**(27 - 10));\n", + "omega_order = omega.multiplicative_order()\n", + "\n", + "k = omega_order // 4\n", + "\n", + "RS = ReedSolomonCode(Fq, omega_order, k, primitive_root=omega)\n", + "evaluation_points = RS.evaluation_points(); # print(f\"Evaluation points: {evaluation_points}\")\n", + "minimum_distance = RS.minimum_distance(); print(f\"Minimum Distance: {minimum_distance}\")\n", + "dimension = RS.dimension(); print(f\"Dimension: {dimension}\")\n", + "rate = RS.rate(); print(f\"Rate: {rate}\")\n", + "E = RS.encoder(\"EvaluationPolynomial\", polynomial_ring=R)" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [], + "source": [ + "def enumerate_points_in_radius(base_point, hamming_distance):\n", + " from itertools import combinations\n", + " base_len = len(base_point)\n", + " result = list()\n", + " for r in range(hamming_distance):\n", + " for entries in combinations(range(base_len), r+1):\n", + " value = list(base_point)\n", + " for j in entries:\n", + " value[j] = 0\n", + " result.append(value)\n", + " return result" + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [], + "source": [ + "def berlekamp_welch(code : ReedSolomonCode, received_vector : VectorSpace):\n", + " n = code.length()\n", + " dim = code.dimension()\n", + " d = code.minimum_distance()\n", + " field = code.base_ring()\n", + " eval_points = code.evaluation_points()\n", + " correctable = d // 2\n", + "\n", + " if len(received_vector) != n:\n", + " raise ValueError(f\"Invalid received vector\")\n", + " \n", + " row_entry = lambda eval_point, y_val : [ eval_point**i for i in range(dim+correctable)] + [field(-1)*y_val*(eval_point**i) for i in range(correctable)]\n", + "\n", + " m = Matrix(field, [ row_entry(x,y) for (x,y) in zip(eval_points, received_vector)])\n", + " # print(m)\n", + " echelon = m.echelon_form()\n", + " # print(echelon)\n", + " if echelon.rank() == n:\n", + " raise ValueError(\"Failed to decode: Possibly too many errors\")\n", + " \n", + " result_vector = field(-1)*echelon.column(echelon.rank())\n", + " a_poly = R(list(result_vector[0:dim+correctable]))\n", + " b_poly = result_vector[dim+correctable:]\n", + " b_poly[echelon.rank() - dim - correctable] = field(1)\n", + " b_poly = R(list(b_poly))\n", + " # print(a_poly)\n", + " # print(b_poly)\n", + " if a_poly % b_poly == 0:\n", + " decoded = a_poly // b_poly\n", + " # print(factor(b_poly))\n", + " return decoded\n", + " else:\n", + " raise ValueError(\"failed to decode!\")\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [], + "source": [ + "fx = R.random_element(k-1); print(fx)\n", + "code = E.encode(fx)\n", + "error_rate = 1 - sqrt(rate)\n", + "error_count = floor(error_rate*omega_order) - 1\n", + "D = GRSGuruswamiSudanDecoder(RS, tau=error_count)\n", + "\n", + "corrupted_code = deepcopy(code)\n", + "\n", + "for i in range(omega_order):\n", + " if error_count > 0:\n", + " if random() < error_rate:\n", + " corrupted_code[i] = Fq.random_element()\n", + " error_count = error_count - 1\n", + " else:\n", + " break\n", + "\n", + "print(D.decode_to_message(corrupted_code))\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "SageMath 10.5", + "language": "sage", + "name": "SageMath-10.5" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "sage", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.12.5" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/sage/fri-and-friends/fri.ipynb b/sage/fri-and-friends/fri.ipynb new file mode 100644 index 000000000..00dc64430 --- /dev/null +++ b/sage/fri-and-friends/fri.ipynb @@ -0,0 +1,208 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# FRI: Fast Reed Solomon IOP of Proximity" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "---\n", + "## FRI Low Degree Tests\n", + "---" + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [], + "source": [ + "%load_ext autoreload\n", + "%autoreload 2\n", + "\n", + "from src.proth_primes import proth_in_range\n", + "from src.fri import FRIOPP, FRIPCS_Prover, FRIPCS_Verifier\n", + "from sage.rings.finite_rings.all import GF\n", + "from sage.rings.polynomial.all import PolynomialRing\n", + "from sage.misc.functional import log\n", + "from sage.misc.prandom import random" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Using prime 2013265921 = 15⋅2^27 + 1\n", + "Maximal 2-sylow multiplicative generator: 440564289\n", + "Maximal 2-sylow multiplicative power: 27\n" + ] + } + ], + "source": [ + "(p, k, n) = proth_in_range(15, 15, 27,27)\n", + "print(f\"Using prime {p} = {k}⋅2^{n} + 1\")\n", + "Fq = GF(p)\n", + "omega = Fq.multiplicative_generator() \n", + "omega = omega**k\n", + "omega_2sylow_expo = log(omega.multiplicative_order(), 2)\n", + "print(f\"Maximal 2-sylow multiplicative generator: {omega}\")\n", + "print(f\"Maximal 2-sylow multiplicative power: {omega_2sylow_expo}\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Setup" + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Using polynomial: 189446473*X^511 + 1887525604*X^510 + 1995045324*X^509 + 92468193*X^508 + 37250970*X^507 + 1782909639*X^506 + 1097782517*X^505 + 1792228040*X^504 + 1508883015*X^503 + 1604780547*X^502 + 831912909*X^501 + 170663884*X^500 + 208740638*X^499 + 1096769607*X^498 + 1320297757*X^497 + 1210580977*X^496 + 686744765*X^495 + 1912317870*X^494 + 583015882*X^493 + 609794307*X^492 + 1143539288*X^491 + 238884464*X^490 + 414376753*X^489 + 1577859556*X^488 + 1478831601*X^487 + 320044573*X^486 + 1068683403*X^485 + 1643769744*X^484 + 968826857*X^483 + 713524189*X^482 + 1559677840*X^481 + 1808414572*X^480 + 616753196*X^479 + 369442567*X^478 + 1281367581*X^477 + 1674756031*X^476 + 1174814323*X^475 + 1200755310*X^474 + 1788093022*X^473 + 1396156594*X^472 + 952892639*X^471 + 235193249*X^470 + 1027910567*X^469 + 1525610627*X^468 + 353148555*X^467 + 94058086*X^466 + 1030223912*X^465 + 1156090947*X^464 + 191786119*X^463 + 628689551*X^462 + 1183025583*X^461 + 527879846*X^460 + 1544754304*X^459 + 1328048902*X^458 + 622098138*X^457 + 423247087*X^456 + 799719438*X^455 + 1219633323*X^454 + 551666405*X^453 + 1815338475*X^452 + 1772679509*X^451 + 21535809*X^450 + 824109832*X^449 + 1502506084*X^448 + 170048639*X^447 + 507606550*X^446 + 1756552490*X^445 + 710132313*X^444 + 289442767*X^443 + 1760375410*X^442 + 852195305*X^441 + 393498842*X^440 + 878539995*X^439 + 451093279*X^438 + 508204195*X^437 + 1851434606*X^436 + 1195804336*X^435 + 872368488*X^434 + 1590197919*X^433 + 756435632*X^432 + 1041236287*X^431 + 1763409252*X^430 + 1467691185*X^429 + 1053850359*X^428 + 789124605*X^427 + 1479151161*X^426 + 1824474173*X^425 + 1247212630*X^424 + 647808870*X^423 + 1467985303*X^422 + 1664515386*X^421 + 429526509*X^420 + 1543889385*X^419 + 1369260402*X^418 + 1202550921*X^417 + 493273609*X^416 + 662114152*X^415 + 1070635791*X^414 + 905264134*X^413 + 202264978*X^412 + 702885828*X^411 + 57677289*X^410 + 1763856263*X^409 + 1310974177*X^408 + 1371097547*X^407 + 1074768231*X^406 + 982645670*X^405 + 493726941*X^404 + 1963735751*X^403 + 1524255046*X^402 + 838324038*X^401 + 82190277*X^400 + 812292128*X^399 + 263913418*X^398 + 1554577155*X^397 + 1934693504*X^396 + 38261105*X^395 + 57343381*X^394 + 1665131977*X^393 + 326767498*X^392 + 1451465397*X^391 + 1552539030*X^390 + 513972749*X^389 + 544183355*X^388 + 1282049580*X^387 + 562276406*X^386 + 645608829*X^385 + 906953496*X^384 + 293187992*X^383 + 554661935*X^382 + 481201063*X^381 + 327790756*X^380 + 419774431*X^379 + 738411995*X^378 + 1919647879*X^377 + 869163814*X^376 + 1903071183*X^375 + 1530682679*X^374 + 1291735824*X^373 + 971000210*X^372 + 1483526357*X^371 + 394462341*X^370 + 632976391*X^369 + 466753642*X^368 + 1902704037*X^367 + 1526492507*X^366 + 512238581*X^365 + 601476050*X^364 + 1096562802*X^363 + 570826774*X^362 + 1544095127*X^361 + 321330125*X^360 + 1578627030*X^359 + 1410338350*X^358 + 618133825*X^357 + 1867538556*X^356 + 649913093*X^355 + 1930232583*X^354 + 1008438008*X^353 + 4324096*X^352 + 1870531586*X^351 + 1101150682*X^350 + 1630963240*X^349 + 469808401*X^348 + 542572408*X^347 + 1642983791*X^346 + 1612824465*X^345 + 313994970*X^344 + 375394635*X^343 + 1840622847*X^342 + 1977684249*X^341 + 1596179631*X^340 + 1256180032*X^339 + 1248429538*X^338 + 652719307*X^337 + 819936998*X^336 + 1905654833*X^335 + 1677948375*X^334 + 675342560*X^333 + 152300371*X^332 + 223334844*X^331 + 1002875053*X^330 + 1360033022*X^329 + 1415470556*X^328 + 1117575678*X^327 + 1424371347*X^326 + 403803173*X^325 + 572966429*X^324 + 946348645*X^323 + 642473629*X^322 + 1949442948*X^321 + 1052848292*X^320 + 1195652273*X^319 + 588431547*X^318 + 1819304866*X^317 + 1207067677*X^316 + 1551931028*X^315 + 953974828*X^314 + 1271327096*X^313 + 127840297*X^312 + 892368977*X^311 + 1752167578*X^310 + 1774065778*X^309 + 674485001*X^308 + 736968797*X^307 + 1807416573*X^306 + 1226651450*X^305 + 402138558*X^304 + 1459925223*X^303 + 1258253478*X^302 + 572449722*X^301 + 1923680340*X^300 + 502677815*X^299 + 253938084*X^298 + 1641881859*X^297 + 1685371939*X^296 + 581154436*X^295 + 715950725*X^294 + 1842149265*X^293 + 924838203*X^292 + 337194851*X^291 + 525818198*X^290 + 1929951982*X^289 + 1208341879*X^288 + 445307976*X^287 + 511032931*X^286 + 1285665474*X^285 + 241078385*X^284 + 1580790359*X^283 + 1985159054*X^282 + 684018134*X^281 + 75237368*X^280 + 1313635714*X^279 + 636743997*X^278 + 1908220390*X^277 + 244374412*X^276 + 1765374566*X^275 + 403233835*X^274 + 165465309*X^273 + 1050353545*X^272 + 1717828450*X^271 + 1634172243*X^270 + 135120435*X^269 + 838574449*X^268 + 1448020483*X^267 + 690186045*X^266 + 372481552*X^265 + 104897040*X^264 + 1787924855*X^263 + 325279998*X^262 + 1846571308*X^261 + 885934331*X^260 + 1995868438*X^259 + 283202203*X^258 + 1431960667*X^257 + 1757449215*X^256 + 1298943657*X^255 + 766449678*X^254 + 1534638189*X^253 + 1977696337*X^252 + 453299350*X^251 + 1788204933*X^250 + 52448742*X^249 + 52979105*X^248 + 1155301902*X^247 + 445982544*X^246 + 745087493*X^245 + 381118281*X^244 + 1622951530*X^243 + 211986480*X^242 + 646157812*X^241 + 1401285583*X^240 + 23966606*X^239 + 1795608601*X^238 + 1482505840*X^237 + 854480431*X^236 + 235954805*X^235 + 1640966116*X^234 + 1363425327*X^233 + 1762214456*X^232 + 841376356*X^231 + 890244364*X^230 + 231673656*X^229 + 1810304453*X^228 + 1656066355*X^227 + 1628109488*X^226 + 1122240141*X^225 + 62379506*X^224 + 198937824*X^223 + 1959519879*X^222 + 617754439*X^221 + 182156090*X^220 + 2004328789*X^219 + 1525003804*X^218 + 1565876698*X^217 + 391069527*X^216 + 84161492*X^215 + 1046041108*X^214 + 208383541*X^213 + 196590431*X^212 + 1070509051*X^211 + 235203173*X^210 + 1571433373*X^209 + 768443247*X^208 + 1439492642*X^207 + 615325143*X^206 + 698809175*X^205 + 80624599*X^204 + 1478841383*X^203 + 1180788455*X^202 + 450178630*X^201 + 746547550*X^200 + 792104307*X^199 + 620806935*X^198 + 1749865336*X^197 + 1951507053*X^196 + 152821486*X^195 + 1509849176*X^194 + 1508853729*X^193 + 1938881556*X^192 + 855588988*X^191 + 1717212881*X^190 + 768345866*X^189 + 793877478*X^188 + 551494720*X^187 + 1525849705*X^186 + 354821665*X^185 + 1405747621*X^184 + 463843666*X^183 + 1917755864*X^182 + 20286255*X^181 + 90479982*X^180 + 1890256469*X^179 + 504264870*X^178 + 283448314*X^177 + 923458140*X^176 + 149377195*X^175 + 512320803*X^174 + 1498716774*X^173 + 439089927*X^172 + 438675833*X^171 + 452684349*X^170 + 1719495284*X^169 + 112432536*X^168 + 1184290568*X^167 + 446685624*X^166 + 1365608428*X^165 + 677606367*X^164 + 1793429113*X^163 + 912585583*X^162 + 636506700*X^161 + 1945131138*X^160 + 1639587267*X^159 + 530923664*X^158 + 1343559637*X^157 + 1522491447*X^156 + 164010869*X^155 + 315784022*X^154 + 1414724736*X^153 + 146903209*X^152 + 1526118659*X^151 + 1377154394*X^150 + 605268500*X^149 + 1662067295*X^148 + 1362588607*X^147 + 1970988203*X^146 + 993447883*X^145 + 156395187*X^144 + 226194282*X^143 + 277613679*X^142 + 1349050218*X^141 + 467585630*X^140 + 886010502*X^139 + 1553702058*X^138 + 1214727876*X^137 + 1827526145*X^136 + 1506276331*X^135 + 120098725*X^134 + 1527555382*X^133 + 884260912*X^132 + 824900166*X^131 + 1711477792*X^130 + 232183281*X^129 + 1655212868*X^128 + 1641598259*X^127 + 1844163621*X^126 + 1251917886*X^125 + 865160444*X^124 + 917798698*X^123 + 436070582*X^122 + 1921152674*X^121 + 1221538369*X^120 + 1507778342*X^119 + 1438089210*X^118 + 33768497*X^117 + 770763602*X^116 + 1337544937*X^115 + 283112724*X^114 + 58132628*X^113 + 826798370*X^112 + 144298771*X^111 + 101394975*X^110 + 1527762484*X^109 + 1179999892*X^108 + 1291335500*X^107 + 662759961*X^106 + 945531676*X^105 + 1919979274*X^104 + 1978016917*X^103 + 95151601*X^102 + 266551751*X^101 + 798150064*X^100 + 628450501*X^99 + 447255949*X^98 + 1975933810*X^97 + 1595492254*X^96 + 1043164156*X^95 + 1182599848*X^94 + 427307373*X^93 + 693976618*X^92 + 959124752*X^91 + 915881406*X^90 + 1710036081*X^89 + 188171293*X^88 + 1647597052*X^87 + 1180565402*X^86 + 1221780901*X^85 + 806374051*X^84 + 310729438*X^83 + 146509083*X^82 + 1125825197*X^81 + 1927254837*X^80 + 56031384*X^79 + 2002449607*X^78 + 1726383262*X^77 + 1081138937*X^76 + 1317990132*X^75 + 77903253*X^74 + 266813979*X^73 + 1471170030*X^72 + 1444124641*X^71 + 1256862904*X^70 + 1825435042*X^69 + 1145109615*X^68 + 1230283053*X^67 + 663005676*X^66 + 35264336*X^65 + 73617278*X^64 + 1284452553*X^63 + 424283812*X^62 + 1716707861*X^61 + 1459849820*X^60 + 1774484883*X^59 + 1797862846*X^58 + 730700496*X^57 + 418592510*X^56 + 1416533445*X^55 + 681761189*X^54 + 597417997*X^53 + 1735180475*X^52 + 434870218*X^51 + 85145610*X^50 + 394520519*X^49 + 1005497894*X^48 + 1160586236*X^47 + 1841999005*X^46 + 928268084*X^45 + 16116123*X^44 + 1627147831*X^43 + 409843777*X^42 + 807804328*X^41 + 1368392557*X^40 + 1279821125*X^39 + 1009982203*X^38 + 1920671559*X^37 + 214081064*X^36 + 1797991119*X^35 + 922213725*X^34 + 760514223*X^33 + 1441198790*X^32 + 1334138967*X^31 + 1772984401*X^30 + 1236831835*X^29 + 1358771284*X^28 + 500321057*X^27 + 280875445*X^26 + 290795411*X^25 + 1762448696*X^24 + 1806818813*X^23 + 761393155*X^22 + 535205564*X^21 + 732571398*X^20 + 65252649*X^19 + 1547717018*X^18 + 1796563894*X^17 + 989207484*X^16 + 326667473*X^15 + 316644654*X^14 + 96551328*X^13 + 950155262*X^12 + 1445348680*X^11 + 936403683*X^10 + 1247751751*X^9 + 459976706*X^8 + 1580525340*X^7 + 1726129736*X^6 + 268590302*X^5 + 178556091*X^4 + 229961231*X^3 + 33549262*X^2 + 379345611*X + 1272304792\n" + ] + } + ], + "source": [ + "polynomial_degree = 512\n", + "folding_factor = 4\n", + "rate_factor = 4\n", + "rs_code_length = rate_factor*polynomial_degree\n", + "rs_code_log_len = int(log(rs_code_length, 2))\n", + "PolyRing = PolynomialRing(Fq, \"X\")\n", + "ntt_omega = omega**(2**(omega_2sylow_expo - rs_code_log_len))\n", + "assert ntt_omega.multiplicative_order() == rs_code_length\n", + "\n", + "input_polynomial = PolyRing.random_element(polynomial_degree-1)\n", + "\n", + "print(f\"Using polynomial: {input_polynomial}\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Low degree test" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "✓ Positive case of low degree test passed\n", + "✓ Negative case case of low degree test passed\n" + ] + } + ], + "source": [ + "iopp = FRIOPP(input_polynomial=input_polynomial, \n", + " ntt_omega=ntt_omega, \n", + " rate_factor=rate_factor, \n", + " folding_factor=folding_factor)\n", + "\n", + "max_degree = polynomial_degree\n", + "\n", + "if iopp.verify_proximity(max_degree=max_degree) is True:\n", + " print(\"✓ Positive case of low degree test passed\")\n", + "else:\n", + " print(\"❌ERROR: Positive case of low dergee test failed!\")\n", + "\n", + "if iopp.verify_proximity(max_degree // folding_factor ) is False:\n", + " print(\"✓ Negative case case of low degree test passed\")\n", + "else:\n", + " print(\"❌ ERROR: Negative case of low dergee test failed!\")\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "---\n", + "## FRI Polynomial Commitment Scheme\n", + "---" + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "✓ A valid proof opening succeeded\n", + "✓ An invalid proof opening failed as expected\n", + "✓ Another invalid proof opening failed as expected\n" + ] + } + ], + "source": [ + "evaluation_point = PolyRing.base_ring().random_element()\n", + "\n", + "prover = FRIPCS_Prover(ntt_omega, rate_factor, folding_factor)\n", + "commitment = prover.commit(input_poly=input_polynomial)\n", + "(expected_value, proof) = prover.open(evaluation_point)\n", + "\n", + "verifier = FRIPCS_Verifier(polynomial_degree - 1)\n", + "\n", + "if verifier.verify_proof(commitment, proof, evaluation_point, expected_value) == True:\n", + " print(f\"✓ A valid proof opening succeeded\")\n", + "else:\n", + " print(f\"❌A valid opening should succeed, but it didn't\")\n", + "\n", + "if verifier.verify_proof(commitment, proof, evaluation_point + 1, expected_value) == False:\n", + " print(\"✓ An invalid proof opening failed as expected\")\n", + "else:\n", + " print(f\"❌Invalid opening of an evaluation point should not succeed, but it did!\")\n", + "\n", + "if verifier.verify_proof(commitment, proof, evaluation_point, expected_value - 1) == False:\n", + " print(\"✓ Another invalid proof opening failed as expected\")\n", + "else:\n", + " print(f\"❌Invalid opening of an evaluation point should not succeed, but it did!\")" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "SageMath 10.5", + "language": "sage", + "name": "sagemath-10.5" + }, + "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.12.5" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +} diff --git a/sage/fri-and-friends/src/chat_gpt_generated.py b/sage/fri-and-friends/src/chat_gpt_generated.py new file mode 100644 index 000000000..755d43659 --- /dev/null +++ b/sage/fri-and-friends/src/chat_gpt_generated.py @@ -0,0 +1,34 @@ +### +### WARNING: All code in this file were generated by ChatGPT!!! +### + +def to_subscript(number): + """ + Convert a decimal number to its subscript representation. + + Parameters: + number (int, float, or str): The number to convert. + + Returns: + str: A string with each digit (and '-' sign) converted to a subscript. + """ + # Mapping for digits and the minus sign to their subscript equivalents. + subscripts = { + '0': '₀', + '1': '₁', + '2': '₂', + '3': '₃', + '4': '₄', + '5': '₅', + '6': '₆', + '7': '₇', + '8': '₈', + '9': '₉', + '-': '₋' + } + + # Convert the input to a string to iterate over each character. + num_str = str(number) + + # Replace each character with its subscript equivalent if available. + return ''.join(subscripts.get(char, char) for char in num_str) diff --git a/sage/fri-and-friends/src/fri.py b/sage/fri-and-friends/src/fri.py new file mode 100644 index 000000000..f3d78a190 --- /dev/null +++ b/sage/fri-and-friends/src/fri.py @@ -0,0 +1,400 @@ +from sage.arith.misc import is_power_of_two +from sage.misc.functional import log, sqrt +from sage.misc.misc_c import prod +from sage.misc.prandom import random +from sage.rings.finite_rings.all import * +from sage.rings.polynomial.all import * +from sage.coding.grs_code import ReedSolomonCode +from sage.matrix.all import * +from sage.modules.all import vector; +from typing import * + +def is_majority(lst): + votes = sum([1 if v else 0 for v in lst]) + return votes > len(lst)//2 + +def is_overwhelming_majority(lst): + votes = sum([1 if v else 0 for v in lst]) + return votes > (2*len(lst))//3 + +## Split a polynomial into powers of `split_degree` +def poly_ksplit(poly, split_degree): + chunks = poly.degree()//split_degree + splits = [ poly.parent()(0) for _ in range(split_degree) ] + indeterminate = poly.parent().gen() + + for (i,coef) in enumerate(poly.list()): + ndx = i % split_degree + deg = i // split_degree + splits[ndx] = splits[ndx] + coef*(indeterminate**deg) + + pp = poly.parent()(0) + + for (j,mini_poly) in enumerate(splits): + pp = pp + mini_poly.subs({ indeterminate : indeterminate** split_degree})*(indeterminate**j) + assert pp == poly, "Reconstruction works!!!" + + return splits + + +class FRIRound: + def __init__(self, + input_polynomial : Polynomial, + ntt_omega : int, # must have order 2^k + rate_factor : int, # Rate of Reed Solomon Code + folding_factor: int # folding factor of FRI + ): + + self.poly = input_polynomial + self.PolyRing = self.poly.parent() + self.Fq = self.PolyRing.base_ring() + self.omega = self.Fq(ntt_omega) + self.rate_factor = rate_factor + + # Assuming small field size, typically less than 64-bits + self.mutl_order = self.omega.multiplicative_order() + + if not is_power_of_two(self.mutl_order): + raise ValueError(f"The multiplicative order of RS generator {self.omega} must be a power of two!") + + if not is_power_of_two(rate_factor) or rate_factor > self.mutl_order: + raise ValueError(f"Invalid rate factor: {rate_factor}. Must be a \ + power of 2 and smaller than the multiplicative \ + order of the RS generator") + + if not is_power_of_two(folding_factor): + raise ValueError(f"Folding factor must be a power of 2") + + self.dimension = self.mutl_order // rate_factor + self.folding_factor = folding_factor + self.folding_root = self.omega**(self.mutl_order // self.folding_factor) + + assert self.folding_root**folding_factor == 1 + + self.RS_CODE = ReedSolomonCode(self.Fq, self.mutl_order, self.dimension, self.omega) + self.encoder = self.RS_CODE.encoder("EvaluationPolynomial", polynomial_ring=input_polynomial.parent()) + self._evaluation_points = None + + def evaluation_points(self): + if not self._evaluation_points: + self._evaluation_points = self.RS_CODE.evaluation_points() + return self._evaluation_points + + def ndx_with_folding_roots(self, ndx): + order = self.mutl_order + step = self.mutl_order // self.folding_factor + indexes = [ (ndx + j*step) % order for j in range(self.folding_factor) ] + if __debug__: + reference = self.evaluation_points()[ndx]**self.folding_factor + entries = [self.evaluation_points()[j]**self.folding_factor for j in indexes] + assert all([reference == e for e in entries]) + + return indexes + + def ud_radius(self): + return (1 - float(1/self.rate_factor))/2.0 + + def jb_radius(self): + return 1 - sqrt(float(1/self.rate_factor)) + + def gen_proof(self, corrupt_fn = None) -> list[FiniteField]: + encoded = self.encoder.encode(self.poly) + + if corrupt_fn: + for i in range(len(encoded)): + if corrupt_fn(i): + encoded[i] = self.Fq.random_element() + return encoded + + def fold(self, fold_challenge, ood_challenges = None): + splits = poly_ksplit(self.poly, self.folding_factor) + poly_tbf = self.PolyRing(0) + + for i in range(self.folding_factor): + poly_tbf += (fold_challenge**i)*splits[i] + + ood_response = None + + if ood_challenges: + if not isinstance(ood_challenges, (list,)): + ood_challenges = [ood_challenges] + + num,den = FRIRound.deep_rational_poly(poly_tbf, ood_challenges=ood_challenges) + + new_num = poly_tbf - num + assert new_num % den == 0 + poly_tbf = new_num // den + ood_response = num / den + + new_omega = self.omega ** self.folding_factor + return (FRIRound(poly_tbf, new_omega, self.rate_factor, self.folding_factor), ood_response) + + def decode_upto_ud(self, proof_code): + decoder = self.RS_CODE.decoder("KeyEquationSyndrome") + decoded = self.PolyRing(list(decoder.decode_to_message(proof_code))) + return decoded + + def decode_upto_jb(self, proof_code): + pass + + @staticmethod + def deep_rational_poly(eval_poly, ood_challenges): + PolyRing = eval_poly.parent() + poly_evals = [(v, eval_poly(v)) for v in ood_challenges] + numerator = PolyRing.lagrange_polynomial(poly_evals) + x = PolyRing.gen() + denomonator = prod([(x - v) for v in ood_challenges]) + return (numerator, denomonator) + + @staticmethod + def consistancy_cross_check(this_round: Self, prev_round : Self, folding_challenge, queries = 1, ood_response = None): + previous_round_proof = prev_round.gen_proof() + this_round_proof = this_round.gen_proof() + alphas = [prev_round.folding_root**i for i in range(prev_round.folding_factor)] + success = list() + + for _ in range(queries): + ndx = int(random()*len(previous_round_proof)) + check_indices = prev_round.ndx_with_folding_roots(ndx) + new_ndx = ndx % this_round.mutl_order + if __debug__: + this_elem = this_round.evaluation_points()[new_ndx] + prev_elem = prev_round.evaluation_points()[check_indices[0]] + expected = prev_elem**prev_round.folding_factor + assert this_elem == expected + + values = [previous_round_proof[i] for i in check_indices] + eval_point = prev_round.evaluation_points()[ndx] + + vander = [[(eval_point*alpha)**k for k in range(this_round.folding_factor)] for alpha in alphas ] + M = Matrix(this_round.Fq, vander) + val = vector(this_round.Fq, values) + coeff = M.inverse()*val + expected = sum([c*(folding_challenge**i) for (i,c) in enumerate(coeff)]) + + value_this_round = this_round_proof[new_ndx] + + check_passed = value_this_round == expected + + if ood_response: + this_round_omega = this_round.evaluation_points()[new_ndx] + num = ood_response.numerator()(this_round_omega) + den = ood_response.denominator()(this_round_omega) + value_this_round = den*value_this_round + num + check_passed = value_this_round == expected + + success.append(check_passed) + + return is_overwhelming_majority(success) + +class FRIOPP: + + class RoundData: + def __init__(self, fri_round, fold_random, deep_poly = None): + self.fri_round = fri_round + self.fold_randomness = fold_random + self.deep_poly = deep_poly + + def __init__(self, + input_polynomial : Polynomial, + ntt_omega : int, # must have order 2^k + rate_factor : int, # Rate of Reed Solomon Code + folding_factor: int, # folding factor of FRI + ood_count: int = 0, # Number of OOD Random samples + ): + self.Rx = input_polynomial.parent() + self.Fq = input_polynomial.base_ring() + self.rounds = list() + first_round = FRIRound(input_polynomial, ntt_omega, rate_factor, folding_factor) + + self.oracle = first_round + + self.rounds.append(FRIOPP.RoundData(first_round, None, None)) + + this_round = first_round + + while this_round.poly.degree() > folding_factor: + fold_challenge = this_round.Fq.random_element() + ood_challenge = [self.Fq.random_element() for _ in range(ood_count) ] if ood_count > 0 else None + (this_round, ood_frac) = this_round.fold(fold_challenge, ood_challenge) + self.rounds.append(FRIOPP.RoundData(this_round, fold_challenge, ood_frac)) + + def verify_proximity(self, max_degree, query_count = 10): + prev_round = self.rounds[0].fri_round + + round_max_degree = max_degree + + if prev_round.dimension < max_degree: + # Can't test degree if max_degree is higher than code dimension + return False + + for this_round in self.rounds[1:]: + this_fri = this_round.fri_round + fold_random = this_round.fold_randomness + deep_frac = this_round.deep_poly + is_good = FRIRound.consistancy_cross_check(this_fri, + prev_round, + fold_random, + query_count, + deep_frac) + if not is_good: + return False + + prev_round = this_fri + + round_max_degree = round_max_degree // this_fri.folding_factor + 1 + + if deep_frac: + round_max_degree -= deep_frac.denominator().degree() + + if round_max_degree <= this_fri.folding_factor : + break + + last_round_poly = self.Rx.lagrange_polynomial(zip(prev_round.evaluation_points(), + prev_round.gen_proof())) + + if last_round_poly.degree() > prev_round.folding_factor: + return False + else: + return True + + def proof_oracle(self): + return self.oracle + +class FRIPCS_Prover: + def __init__(self, + ntt_omega : int, # must have order 2^k + rate_factor : int, # Rate of Reed Solomon Code + folding_factor: int, # folding factor of FRI + ood_count: int = 0, # Number of OOD Random samples + ): + self.ntt_omega = ntt_omega + self.rate_factor = rate_factor + self.folding_factor = folding_factor + self.ood_count = ood_count + + def commit(self, input_poly : Polynomial): + self.input_poly = input_poly + self.input_oracle = FRIOPP(self.input_poly, + self.ntt_omega, + self.rate_factor, + self.folding_factor, + self.ood_count) + return self.input_oracle + + def open(self, evaluation_point): + evaluated_value = self.input_poly(evaluation_point) + quot = (self.input_poly - evaluated_value) // (self.input_poly.variables()[0] - evaluation_point) + quot_oracle = FRIOPP(quot, + self.ntt_omega, + self.rate_factor, + self.folding_factor, + self.ood_count) + return (evaluated_value, quot_oracle) + + +class FRIPCS_Verifier: + def __init__(self, poly_degree, query_count = 10): + self.poly_degree = poly_degree + self.query_count = query_count + pass + + def verify_proof(self, + original_commit : FRIOPP, + quot_commit : FRIOPP, + evaluation_point : FiniteField, + expected_value : FiniteField): + + if original_commit.verify_proximity(self.poly_degree, self.query_count) == False: + return False + + if quot_commit.verify_proximity(self.poly_degree - 1, self.query_count) == False: + return False + + # TODO: Check that RS data rates are compatible + + original = original_commit.proof_oracle() + quot = quot_commit.proof_oracle() + eval_points = original.evaluation_points() + + if eval_points != quot.evaluation_points(): + # Reed Solomon evaluation domains are incompatible! + return False + + rs_orignal = original.gen_proof() + rs_quot = quot.gen_proof() + votes = 0 + + for _ in range(self.query_count): + ndx = int(random()*len(rs_orignal)) + reconstructed = expected_value + (rs_quot[ndx] * (eval_points[ndx] - evaluation_point)) + if reconstructed == rs_orignal[ndx]: + votes += 1 + + # Overwhelming majority might be too weak! + if votes > (2*self.query_count) // 3: + return True + else: + return False + + +if __name__ == '__main__': + from proth_primes import proth_in_range + + (p, cofactor, twodicity) = proth_in_range(15, 15, 27,27) + Fq = GF(p) + + PolyRing = PolynomialRing(Fq, "X") + omega = Fq.multiplicative_generator() + + def run_test_friopp(code_dimension=2048, rate_factor=4, folding_factor=4, ood_count = 5): + ntt_omega = omega**(cofactor*2**(twodicity - log(code_dimension, 2))) + poly_degree = code_dimension//rate_factor - int(random()*(code_dimension//10)) + input = PolyRing.random_element(poly_degree - 1) + + iopp = FRIOPP(input_polynomial=input, + ntt_omega=ntt_omega, + rate_factor=rate_factor, + folding_factor=folding_factor, + ood_count=ood_count) + + assert iopp.verify_proximity(poly_degree), "FRI IOPP varification should work" + assert not iopp.verify_proximity(poly_degree // folding_factor - 1), "FRI IOPP varification should fail if the max degree is 1/folding factor away" + + + def run_test_fri_round(code_dimension=2048, rate_factor=4, folding_factor=4, ood_count = 5): + ntt_omega = omega**(cofactor*2**(twodicity - log(code_dimension, 2))) + poly_degree = code_dimension//rate_factor + input = PolyRing.random_element(poly_degree-1) + + first_oracle = FRIRound(input, ntt_omega, rate_factor, 2) + folding_randomness = Fq.random_element() + ood_randomness = [Fq.random_element() for _ in range(ood_count)] + (second_oracle, ood_response) = first_oracle.fold(folding_randomness, ood_randomness) + assert FRIRound.consistancy_cross_check(second_oracle, first_oracle, folding_challenge=folding_randomness, queries=5, ood_response=ood_response) + + def run_test_pcs(code_dimension=2048, rate_factor=4, folding_factor=4, ood_count = 5): + ntt_omega = omega**(cofactor*2**(twodicity - log(code_dimension, 2))) + poly_degree = code_dimension//rate_factor - int(random()*(code_dimension//10)) + input_poly = PolyRing.random_element(poly_degree - 1) + evaluation_point = PolyRing.base_ring().random_element() + + prover = FRIPCS_Prover(ntt_omega, rate_factor, folding_factor, ood_count) + commitment = prover.commit(input_poly=input_poly) + (expected_value, proof) = prover.open(evaluation_point) + + verifier = FRIPCS_Verifier(poly_degree) + + assert verifier.verify_proof(commitment, proof, evaluation_point, expected_value), \ + "A valid opening should succeed, but it didn't" + + assert not verifier.verify_proof(commitment, proof, evaluation_point + 1, expected_value), \ + "An invalid opening of evaluation point should not succeed, but it did" + + assert not verifier.verify_proof(commitment, proof, evaluation_point, expected_value + 1), \ + "An invalid opening of evaluation value shoud not succeed, but it did" + + run_test_fri_round() + run_test_friopp() + run_test_pcs() + \ No newline at end of file diff --git a/sage/fri-and-friends/src/proth_primes.py b/sage/fri-and-friends/src/proth_primes.py new file mode 100644 index 000000000..58b73d83a --- /dev/null +++ b/sage/fri-and-friends/src/proth_primes.py @@ -0,0 +1,110 @@ +# +# Taken from http://www.prothsearch.com/riesel1.html +# +# This is a list of PROTH primes, which are of the form (k.2^n + 1) +# This dictionary has keys that k and values that are n +# +# + +from sage.arith.misc import is_prime, is_power_of_two, factor + +PROTH_PRIMES_300 = { +1 : [2, 4, 8, 16], + +3 : [1, 2, 5, 6, 8, 12, 18, 30, 36, 41, 66, 189, 201, 209, 276, 353, 408, 438, 534, 2208, 2816, 3168, 3189, 3912, 20909, 34350, 42294, 42665, 44685, 48150, 54792, 55182, 59973, 80190, 157169, 213321, 303093, 362765, 382449, 709968, 801978, 916773, 1832496, 2145353, 2291610, 2478785, 5082306, 7033641, 10829346, 16408818 ], + +5 : [1, 3, 7, 13, 15, 25, 39, 55, 75, 85, 127, 1947, 3313, 4687, 5947, 13165, 23473, 26607, 125413, 209787, 240937, 819739, 1282755, 1320487, 1777515 ], + +7 : [2, 4, 6, 14, 20, 26, 50, 52, 92, 120, 174, 180, 190, 290, 320, 390, 432, 616, 830, 1804, 2256, 6614, 13496, 15494, 16696, 22386, 54486, 88066, 95330, 207084, 283034, 561816, 804534, 811230, 1491852, 2139912, 2167800, 2915954, 3015762, 3511774, 5775996, 18233956, 20267500 ], + +9 : [1, 2, 3, 6, 7, 11, 14, 17, 33, 42, 43, 63, 65, 67, 81, 134, 162, 206, 211, 366, 663, 782, 1305, 1411, 1494, 2297, 2826, 3230, 3354, 3417, 3690, 4842, 5802, 6937, 7967, 9431, 13903, 22603, 24422, 39186, 43963, 47003, 49902, 67943, 114854, 127003, 145247, 147073, 149143, 304607, 384990, 412034, 435743, 461081, 834810, 1051026, 1807574, 2543551, 3497442, 5642513, 9778263, 11158963, 11366286, 11500843, 12406887, 13334487], + +11 : [1, 3, 5, 7, 19, 21, 43, 81, 125, 127, 209, 211, 3225, 4543, 10179, 15329, 18759, 28277, 93279, 105741, 268009, 412447, 525589, 644677, 886071, 960901, 1343347, 2230369, 2476839, 2691961, 2897409, 3771821, 8103463, 9381365, 10797109, 10803449, 15502315], + +13 : [2, 8, 10, 20, 28, 82, 188, 308, 316, 1000, 28280, 38008, 43856, 88018, 109258, 114296, 521306, 562456, 684560, 1038896, 1499876, 1861732, 4998362, 5523860, 6481780, 8989858, 15294536, 16828072], + +15 : [1, 2, 4, 9, 10, 12, 27, 37, 38, 44, 48, 78, 112, 168, 229, 297, 339, 517, 522, 654, 900, 1518, 2808, 2875, 3128, 3888, 4410, 6804, 7050, 7392, 19219, 21445, 21550, 24105, 24995, 34224, 34260, 43388, 48444, 61758, 184290, 294894, 300488, 403929, 483098, 635989, 701902, 734400, 1229600, 1244377, 1276177, 1368428, 1418605, 1597510, 1667744, 2393365, 2785940, 2988834, 3162659, 4246384, 4800315, 7300254, 7619838, 9312889, 9830108], + +17 : [3, 15, 27, 51, 147, 243, 267, 347, 471, 747, 2163, 3087, 5355, 6539, 7311, 99231, 824451, 1388355, 1990299, 8636199], + +19 : [6, 10, 46, 366, 1246, 2038, 4386, 4438, 6838, 7498, 7998, 9450, 11890, 17034, 23290, 73338, 149146, 323638, 853546, 2206266, 6833086], + +21: [1, 4, 5, 7, 9, 12, 16, 17, 41, 124, 128, 129, 187, 209, 276, 313, 397, 899, 1532, 1613, 1969, 2245, 2733, 4585, 4644, 6712, 6981, 13344, 17524, 27124, 29769, 47337, 55828, 91008, 94801, 164901, 179457, 191337, 200568, 220992, 686632, 856865, 1022168, 1240067, 1421741, 1830919, 3065701, 6048861, 9000000], + +23: [1, 9, 13, 29, 41, 49, 69, 73, 341, 381, 389, 649, 1961, 3929, 977541, 1399841, 1448461, 2425641, 4300741, 9000000], + +25: [2, 4, 6, 10, 20, 22, 52, 64, 78, 184, 232, 268, 340, 448, 554, 664, 740, 748, 1280, 1328, 1640, 3314, 3904, 3938, 5152, 9522, 57488, 66872, 148060, 154254, 216092, 302194, 367294, 627710, 966414, 1211488, 1258562, 1298186, 2141884, 2583690, 2927222, 3733144, 4481024, 8456828, 8788628, 9000000, 13719266], + +27: [2, 4, 7, 16, 19, 20, 22, 26, 40, 44, 46, 47, 50, 56, 59, 64, 92, 175, 215, 275, 407, 455, 1076, 1090, 3080, 3322, 6419, 7639, 19360, 30500, 38770, 44164, 50696, 114760, 117784, 153847, 160874, 178739, 379880, 405514, 672007, 974752, 1164664, 1476347, 2218064, 5213635, 7046834, 7963247, 8437000, 12184319], + +29: [1, 3, 5, 11, 27, 43, 57, 75, 77, 93, 103, 143, 185, 231, 245, 391, 1053, 1175, 2027, 3627, 4727, 5443, 7927, 8533, 9067, 14185, 25723, 88117, 96947, 144937, 219317, 252227, 261031, 339123, 343601, 627167, 742191, 800191, 1152765, 1416873, 1574753, 3964697, 4532463, 7374577, 7899985, 9000000], + +31: [8, 60, 68, 140, 216, 416, 1808, 1944, 9096, 76596, 93168, 4673544, 5560820, 8348000, 9000000], + +33: [1, 6, 13, 18, 21, 22, 25, 28, 66, 93, 118, 289, 412, 453, 525, 726, 828, 1420, 1630, 3076, 3118, 4452, 4941, 5236, 6346, 9133, 13401, 14214, 18766, 37249, 42685, 47805, 61372, 74178, 600270, 922782, 1130884, 2345001, 3176269, 3570132, 3649810, 9000000], + +35: [1, 3, 7, 9, 13, 15, 31, 45, 47, 49, 55, 147, 245, 327, 355, 663, 1423, 1443, 2493, 3627, 4737, 6541, 9667, 22645, 27321, 61963, 68281, 144397, 151741, 187859, 191509, 273251, 494607, 696935, 768063, 831411, 1040005, 1357881, 1765449, 2099769, 2241049, 3570777, 3587843, 9000000], + +37: [2, 4, 8, 10, 12, 16, 22, 26, 68, 82, 84, 106, 110, 166, 236, 254, 286, 290, 712, 1240, 1706, 1804, 1904, 2240, 2632, 3104, 5870, 8750, 11514, 11972, 14434, 16322, 18044, 19050, 23014, 32424, 59354, 62350, 120196, 157128, 167326, 211056, 218550, 230716, 233968, 239902, 256954, 331144, 548012, 551830, 939364, 2128328, 4046360, 9000000, 10599476, 11855148, 14166940, 15474010] + +} + +def proth_prime(cofactor, twoadicity): + try: + cofactor_list = PROTH_PRIMES_300[cofactor] + for v in cofactor_list: + if v == twoadicity: + return cofactor*(2**twoadicity) + 1 + raise IndexError + except IndexError as ie: + raise ValueError(f"proth prime with multiplicative cofactor {cofactor} and to-adic valuation of multiplicative order {twoadicity} not found") + +def proth_in_range(cofactor_min, cofactor_max, twodicity_min, twodicity_max) -> tuple[int, int, int]: + for key,values in PROTH_PRIMES_300.items(): + if key < cofactor_min: + continue + if key > cofactor_max: + continue + + for v in values: + if v < twodicity_min: + continue + elif v > twodicity_max: + break + else: + return (key*(2**v) + 1, key, v) + + raise ValueError(f"No proth prime found with multiplicative co-factor in the range ({cofactor_min}, {cofactor_max}) and multiplicative two-adic valuation in the range ({twodicity_min}, {twodicity_max}) found") + + +if __name__=='__main__': + + def test_range_primes(): + k_min = 22 + k_max = 26 + two_min = 100 + two_max = 500 + + (p, cofactor, twodic_val) = proth_in_range(k_min, k_max, two_min, two_max) + assert p == cofactor*(2**twodic_val) + 1 + assert two_min <= twodic_val and twodic_val <= two_max + assert k_min <= cofactor and cofactor <= k_max + + + def test_proth_primes(): + for key,values in PROTH_PRIMES_300.items(): + for v in values: + if v < 256: + should_find = proth_prime(key, v) + assert is_prime(should_find) + assert is_power_of_two((should_find - 1)//key) + else: + continue + + try: + should_not_find = proth_prime(15, 11) + except ValueError as e: + pass + + test_range_primes() + test_proth_primes() \ No newline at end of file diff --git a/sage/fri-and-friends/src/whir.py b/sage/fri-and-friends/src/whir.py new file mode 100644 index 000000000..19e20273e --- /dev/null +++ b/sage/fri-and-friends/src/whir.py @@ -0,0 +1,673 @@ +import os +import sys +src_dir = os.path.dirname(os.path.abspath(__file__)) +sys.path.insert(0, src_dir) + +from sage.rings.polynomial.all import PolynomialRing, Polynomial +from sage.arith.misc import is_power_of_two +from sage.rings.finite_rings.all import GF, FiniteField +from sage.misc.functional import log +from sage.misc.misc_c import prod +from sage.rings.integer import Integer +from chat_gpt_generated import to_subscript +from sage.coding.grs_code import ReedSolomonCode +from sage.matrix.all import * +from sage.modules.all import vector; +from sage.misc.prandom import random +from typing import * + + +def compute_hypercube_sum(sumcheck_poly, skip_variables = [], hypercube = [0,1]): + sc_val = sumcheck_poly + variables = sumcheck_poly.variables() + skip_variables = skip_variables or [] + + for xvar in variables: + if xvar in skip_variables: + continue + + this_round = 0 + + for hval in hypercube: + this_round += sc_val.subs({xvar : hval}) + + sc_val = this_round + + return sc_val + + +def eq_list(poly_vars, eq_values): + if len(poly_vars) != len(eq_values): + raise ValueError(f"The number of variables in the base ring does not match eq_values") + + one_term = lambda x_var, r_val : (x_var*r_val + (1-x_var)*(1-r_val)) + + result = poly_vars[0].parent()(1) + + for (x_var, r_val) in zip(poly_vars, eq_values): + result *= one_term(x_var, r_val) + + return result + +def eq_polynomial(poly_ring, eq_values): + Fq = poly_ring.base_ring() + poly_vars = poly_ring.gens() + return eq_list(poly_vars, eq_values) + +def monomial_basis(poly_ring): + poly_vars = list(poly_ring.gens()) + basis = list() + for i in range(2**len(poly_vars)): + b = Integer(i).bits() + this_round = prod([x**b for (x,b) in zip(poly_vars, b)]) + basis.append(this_round) + return basis + +def random_multilinear_poly(poly_ring): + basis = monomial_basis(poly_ring) + result = poly_ring(0) + + for b in basis: + c = poly_ring.base_ring().random_element() + result = result + c*b + + return result + +def eq_weight_polynomial(WeightRing, evaluation_point): + base_poly = WeightRing.base_ring() + ep = eq_polynomial(base_poly, evaluation_point) + z = WeightRing.gen() + return z* WeightRing(ep) + +def gen_mvariate_ring(base_ring, variables_count : int, var_name = "X"): + variable_names = tuple([f"{var_name}{to_subscript(i)}" for i in range(variables_count)]) + return PolynomialRing(base_ring=base_ring, names = variable_names ) + +def tensor_map(sumcheck_challenges): + """ + Given the list of Sumcheck challenges, say + [a_0, a_1, •••, a_{ℓ-1}], it computes the 2^ℓ tensor products + of (1, a_0)⊗(1, a_1)⊗ ••• ⊗(1, a_{ℓ-1}) and stores + them in a dictionary of size 2^ℓ. + + Since there's a bijection between univariate polynomials and multilinear + polynomial, given a polynomial written as a bilinear form + + f(x) = <(f_0, f_1, •••, f_{2^{ℓ}-1}), (1, x, x^2, •••, x^{2^{ℓ}-1})> + + one can write + + fold(f, [a_0, a_1, •••, a_ℓ] ) = + <(f_0, f_1, •••, f_{2^{ℓ}-1}), (1, a_0)⊗(1, a_0)⊗ ••• ⊗(1, a_{ℓ-1})> + + This function therefore computes (1, a_0)⊗(1, a_0)⊗ ••• ⊗(1, a_{ℓ-1}) where a_i's + are the sumcheck challenges. + """ + degree_map = dict() + + for i in range(2**len(sumcheck_challenges)): + bits = Integer(i).bits() + value = 1 + + for (j,v) in enumerate(bits): + if v == 1: + value *= sumcheck_challenges[j] + + degree_map[i] = value + + return degree_map + + +def fold_virtually(oracle_ndxes : list[int], # list of indices in proof oracle + proof_oracle : vector, # Reed-Solomon code word + fold_challenges : list[FiniteField], # fold_factor = 2^|fold_challenge| + ntt_omega : FiniteField, # Generator of Reed-Solomon code word + ): + """ + This function computes the folded value of `ndx`-th element from + proof oracle. There might be more efficient ways to achieve this + that does not involve Vandermonde Matrix, but this is the + cleanest algorithm I could think of. Here's how the algorithm + works. + + 1. First compute the fold_factor ℓ := 2^{len(fold_challenges)} + 2. Compute ℓ-th roots of unity `fold_omega` from `ntt_omega`. + 2. Given `ndx`, compute other indicies that are ℓ-th root + conjugates corresponding to `ndx`. + 3. Given a polynomial f(x), it can be written as + f(x) = f_0(x^ℓ) + x•f_1(x^ℓ) + ••• + x^{ℓ}•f_{ℓ-1}(x^ℓ) + Therefore for all i, (ntt_omega•fold_omega^i) ^ ℓ = ntt_omega^ℓ. + This gives us a way to compute f_j(x^ℓ) for all j. + + 4. Once the values of f_j(ntt_omega^{4•ndx}) are computed for all j, + uses the bijection between the multilinear polynomial and univariate + polynomial to compute the fold. See function `tensor_map` for details. + + """ + + # Compute the Sumcheck tensor and cache it + schk_tensor = tensor_map(fold_challenges) + + folding_factor = 2**len(fold_challenges) + order = len(proof_oracle) + + assert order > folding_factor and order % folding_factor == 0 + + fold_omega = ntt_omega**(order // folding_factor) + + assert fold_omega**folding_factor == 1 + + Fq = proof_oracle[0].base_ring() + step = order // folding_factor + + result = list() + + for ndx in oracle_ndxes: + # Compute the query indices. + query_ndx = [ (ndx + j*step) % order for j in range(folding_factor) ] + query_resp = vector([proof_oracle[i] for i in query_ndx]) + q_omega = ntt_omega**ndx + + # + # Compute the Vandermonde matrix for the fold_factor root + # conjugates of q_omega + # + vander = [[(q_omega*(fold_omega**j))**k for k in range(folding_factor)] \ + for j in range(folding_factor)] + + M = Matrix(Fq, vander) + coefficients = M.inverse()*query_resp + + fold_value = Fq(0) + for (i, coeff) in enumerate(coefficients): + fold_value = fold_value + schk_tensor[i]*coeff + + result.append((q_omega**folding_factor, fold_value)) + + return result + +class WhirRound: + + class FoldData: + def __init__(self, + last_whir_round: Self, + whir_round : Self, + ood_data: dict[FiniteField, FiniteField] , + shift_qeries : list[FiniteField], + ): + self.last_round = last_whir_round + self.this_round = whir_round + self.ood_data = ood_data + self.shift_queries = shift_qeries + + + def __init__(self, + input_polynomial : Polynomial, # Witness multilinear polynomial + weight_polynomial : Polynomial, # Compatible with input Ring + ntt_omega : FiniteField, # must have order 2^k + claimed_sum : FiniteField = None, # Claimed value of sum + ntt_order_omega : int = None, # Multiplicative order of generator + hypercube = [0,1] # Hypercube to sum over + ): + + assert weight_polynomial.parent().base_ring() == input_polynomial.parent(), \ + "Weight polynomial should be a polynomial ring over input polynomial" + + self.poly = input_polynomial + self.weight_poly = weight_polynomial + + self.Fq = input_polynomial.base_ring() + self.PolyRing = input_polynomial.parent() + + self.polyvars = input_polynomial.variables() + self.weight_var = weight_polynomial.variables()[0] + + self.omega = self.Fq(ntt_omega) + self.multi_order = ntt_order_omega or ntt_omega.multiplicative_order() + self.dimension = 2**len(self.polyvars) + + if not is_power_of_two(self.multi_order): + raise ValueError(f"The order of multiplicative generator must be a power of two") + + if self.dimension > self.multi_order // 2: + raise ValueError(f"The multiplicative order of RS is too small compared to the input variables count") + + self.RS_CODE = ReedSolomonCode(self.Fq, self.multi_order, self.dimension, self.omega) + self.RS_RING = PolynomialRing(self.Fq, "h") + + self.encoder = self.RS_CODE.encoder("EvaluationPolynomial", polynomial_ring=self.RS_RING) + self._evaluation_points = None + + self.hypercube = hypercube + self.sumcheck_challenges = list() + self.sumcheck_poly = self.weight_poly.subs({ self.weight_var : self.poly }) + self.claimed_sum = claimed_sum or self.do_hypercube_sum(None) + + if claimed_sum and __debug__: + check_match = self.do_hypercube_sum(None) + assert check_match == claimed_sum + + temp_var = self.RS_RING.gen() + univar = WhirRound.to_univariate(self.poly, temp_var) + self._code_value = self.encoder.encode(univar) + + + def weight_coefficient(self): + return self.weight_poly.coefficient(1) + + + def evaluation_points(self): + if not self._evaluation_points: + self._evaluation_points = self.RS_CODE.evaluation_points() + + return self._evaluation_points + + def fri_oracle(self): + assert self._code_value is not None + return self._code_value + + def set_fri_oracle(self, new_oracle): + self.set_fri_oracle = new_oracle + + def folded_poly(self): + temp_var = self.RS_RING.gen() + return WhirRound.to_univariate(self.poly, temp_var) + + def ntt_value(self, index): + if index > self.RS_CODE.length(): + raise ValueError("Invalid index for RS Code") + + + @staticmethod + def eval_as_univarite(multi_variate : Polynomial, value : FiniteField) -> FiniteField: + from copy import deepcopy + + result = deepcopy(multi_variate) + variables = result.variables() + + for (i,v) in enumerate(variables): + subst_value = value**(2**i) + result = result.subs( { v : subst_value }) + + return result + + @staticmethod + def to_univariate(mvar_poly : Polynomial, univariate_gen : Polynomial) -> FiniteField: + from copy import deepcopy + result = deepcopy(mvar_poly).change_ring(univariate_gen.parent()) + variables = result.variables() + for (ndx,var) in enumerate(variables): + result = result.subs( {var : univariate_gen**(2**ndx) }) + return univariate_gen.parent()(result) + + + def do_hypercube_sum(self, skip_variables = []): + return compute_hypercube_sum(self.sumcheck_poly, skip_variables, self.hypercube) + + def sumcheck_prover_round(self, sumcheck_challenge=None): + if sumcheck_challenge is None and len(self.sumcheck_challenges) > 0: + raise ValueError(f"Only first round of sumcheck is allowed to have None as challenge") + + if len(self.sumcheck_challenges) == len(self.polyvars): + return self.sumcheck_poly + + + if sumcheck_challenge is not None: + this_var = self.polyvars[len(self.sumcheck_challenges)] + self.sumcheck_challenges.append(sumcheck_challenge) + + # Update the Sumcheck polynomial with new randomness + self.sumcheck_poly = self.sumcheck_poly.subs({this_var : sumcheck_challenge }) + + # Update the polynomial in tandem with Sumcheck polynomial + self.poly = self.poly.subs({this_var : sumcheck_challenge }) + + # Update the weight polynomial so that at the end we only need to add + # swift combination rounds + weight_update = self.weight_coefficient().subs( {this_var : sumcheck_challenge }) + self.weight_poly = self.weight_var*weight_update + + skip_var = self.polyvars[len(self.sumcheck_challenges)] + + # Get the Sumcheck round polynomial by skipping the first Variable from updated Sumcheck polynomial + result = self.do_hypercube_sum([skip_var]) + interminate = result.variables()[0] + self.claimed_sum = result.subs({ interminate : 0}) + result.subs({ interminate : 1}) + + return result + + def update_weight(self, + gamma: FiniteField, + ood_sample: FiniteField, + shift_query: list[FiniteField], + trace : bool = True + ): + + folded_poly = self.folded_poly() + ood_reply = folded_poly(ood_sample) + + weight_eq_vars = self.poly.variables() + + frobenious_char2 = lambda inp : [ inp**(2**i) for i,_ in enumerate(weight_eq_vars) ] + weight_updt = lambda point : eq_list(weight_eq_vars, frobenious_char2(point) ) + + weight = weight_updt(ood_sample) + + if trace: + print(f"\n P> Weight update for OOD challenge {ood_sample} => {ood_reply} : {self.weight_var}*({weight})") + + update_sum = gamma*weight_updt(ood_sample) + + for (i,query) in enumerate(shift_query): + weight = weight_updt(query) + if trace: + print(f" P> Weight update for shift query {query} : {self.weight_var}*({weight})") + update_sum += weight*gamma**(i+2) + + new_weight = self.weight_poly + self.weight_var * self.weight_var.parent()(update_sum) + + if trace: + print(f" P> Updated new weight: {new_weight}") + + h_terms = [ood_reply] + [folded_poly(q) for q in shift_query] + h_update = sum([ sresp * gamma**(i+1) for (i,sresp) in enumerate(h_terms)]) + + if trace: + print(f" P> Sumcheck claim combination randomness additional sum : {h_update}") + + return (ood_reply, new_weight, h_update) + + def fold_round(self, + gamma : FiniteField, + ood_sample : FiniteField, # OOD randomness + shift_queries : list[FiniteField], # Verifier's shift queries + trace = True + ) -> FoldData : + + generator = self.omega + (ood_resp, new_weight, h_update) = self.update_weight(gamma, ood_sample, shift_queries, trace=trace) + claimed_sum = h_update + self.claimed_sum + + if trace: + print(f" P> Updated Sumcheck claim: {claimed_sum}") + + next_whir_round = WhirRound(self.poly, + new_weight, + ntt_omega=generator**2, + claimed_sum=claimed_sum) + + last_whir_round = self + + return WhirRound.FoldData(last_whir_round, + next_whir_round, + { ood_sample : ood_resp}, + shift_queries) + +class WhirVerifier: + + def __init__(self, + weight_polynomial : Polynomial, + ntt_omega : FiniteField, # must have order 2^k + claimed_sum : FiniteField, # Claimed value of sum + fold_factor : int = 2, # Determines Sumcheck rounds + shift_query_count : int = 5, # Number of shift queries to make + ntt_omega_order : int = None, # Multiplicative order of generator + hypercube = [0,1] # Hypercube to sum over + ): + self.weight_polynomial = weight_polynomial + self.verifier_omega = ntt_omega + self.claimed_sum = claimed_sum + self.hypercube = hypercube + self.schk_rounds = int(log(fold_factor, 2)) + self.shift_query_count = shift_query_count + + self.Fq = weight_polynomial.base_ring().base_ring() + self.weight_var = weight_polynomial.variables()[0] + + # + # Warning: If the order is not given and the multiplicative group is + # not smooth for what ever reason (not our use case), Sage/Pari + # falls back to factoring and can take a long time for large primes. + # + self.mult_order = ntt_omega_order or ntt_omega.multiplicative_order() + assert is_power_of_two(self.mult_order) and \ + self.verifier_omega**self.mult_order == 1 and \ + self.verifier_omega**(self.mult_order // 2) != 1 + + self.schk_challenges = [] + + def weight_coefficient(self): + return self.weight_polynomial.coefficient(1) + + def polyvars(self): + variables = self.weight_coefficient().variables() + return variables + + def do_sumcheck(self, + prover : WhirRound, + print_round_poly : bool = True): + prev_poly = None + prev_var = None + claimed_sum = self.claimed_sum + + if print_round_poly: + print(f"Current Sumcheck claim: {self.claimed_sum}") + + polyvars = self.polyvars() + + for (i,v) in zip(range(self.schk_rounds+1), polyvars): + if i > 0: + challenge = self.Fq.random_element() + claimed_sum = prev_poly.subs({prev_var : challenge}) + self.schk_challenges.append(challenge) + weight_update = self.weight_coefficient().subs({prev_var : challenge }) + self.weight_polynomial = self.weight_var*weight_update + else: + challenge = None + + h_poly = prover.sumcheck_prover_round(challenge) + + if print_round_poly: + print(f"Round-{i} poly: {h_poly}, challenge: {challenge}") + + hyper_sum = sum( [ h_poly.subs( {v : h} ) for h in self.hypercube] ) + if hyper_sum != claimed_sum: + return False + elif print_round_poly: + print(f"Hypercube sum {hyper_sum} matches claimed sum {claimed_sum}") + + prev_poly = h_poly + prev_var = v + + self.claimed_sum = claimed_sum + + if print_round_poly: + print(f"Updated Sumcheck claim: {self.claimed_sum}") + + def virtual_shift_response(self, + last_round : WhirRound, + shift_queries : list[FiniteField], + shift_queries_ndxes : list[int]): + + assert len(shift_queries) == len(shift_queries_ndxes) + schk_challenges = self.schk_challenges[-self.schk_rounds:] + proof_oracle = last_round.fri_oracle() + ntt_omega = last_round.omega + + result = fold_virtually(oracle_ndxes=shift_queries_ndxes , + proof_oracle=proof_oracle, + fold_challenges=schk_challenges, + ntt_omega=ntt_omega) + + if __debug__: + folded_poly = last_round.folded_poly() + virtual_responses = [(query, folded_poly(query)) for query in shift_queries] + + for (a,b) in zip(virtual_responses, result): + assert a == b + + return result + + def update_weight(self, + gamma: FiniteField, + ood_query: dict[FiniteField, FiniteField], + shift_query: list[(FiniteField, FiniteField)], + trace : bool = True + ): + ood_challenge = list(ood_query.keys())[0] + ood_reponse = ood_query[ood_challenge] + + weight_eq_vars = self.polyvars() + + frobenious_char2 = lambda inp : [ inp**(2**i) for i,_ in enumerate(weight_eq_vars) ] + weight_updt = lambda point : eq_list(weight_eq_vars, frobenious_char2(point) ) + + weight = weight_updt(ood_challenge) + + if trace: + print(f"\n V> Weight update for OOD challenge {ood_challenge} => {ood_reponse} : {self.weight_var}*({weight})") + + update_sum = gamma*weight + + for (i,(query,resp)) in enumerate(shift_query): + weight = weight_updt(query) + if trace: + print(f" V> Weight update for shift query {query} => {resp} : {self.weight_var}*({weight})") + update_sum += weight*gamma**(i+2) + + self.weight_polynomial += self.weight_var * self.weight_var.parent()(update_sum) + + if trace: + print(f" V> Updated new weight: {self.weight_polynomial}") + + h_terms = [ood_reponse] + [r for (_,r) in shift_query] + h_update = sum([ sresp * gamma**(i+1) for (i,sresp) in enumerate(h_terms)]) + + if trace: + print(f" V> Sumcheck claim combination randomness additional sum : {h_update}") + + self.claimed_sum += h_update + + if trace: + print(f" V> Updated Sumcheck claim: {self.claimed_sum}") + + def fold_check(self, prover: WhirRound, shift_query_count : int, trace=True): + gamma = self.Fq.random_element() # This should be sent after OOD response + ood_challenge = self.Fq.random_element() # TODO: Check it's not in evaluation domain + + fold_factor = 2**self.schk_rounds + assert prover.multi_order % fold_factor == 0 + + shift_omega = prover.omega**fold_factor + shift_query_ndxes = [int(prover.multi_order*random()) for _ in range(shift_query_count)] + + shift_queries = [shift_omega**i for i in shift_query_ndxes] + + if trace: + print(f"Current Sumcheck claim: {self.claimed_sum}") + print(f"OOD challenge: {ood_challenge}") + print(f"Shift queries: {shift_queries}") + print(f"Combination randomness 𝛾: {gamma}") + + prover_fold = prover.fold_round(gamma, ood_challenge, shift_queries, trace=trace) + + last_round = prover_fold.last_round + new_round = prover_fold.this_round + ood_data = prover_fold.ood_data + + shift_responses = self.virtual_shift_response(last_round, + shift_queries, + shift_query_ndxes) + + if trace: + print(f"\nProver OOD response: {next(iter(ood_data.values()))}") + print(f"Verifier computed shift response: {shift_responses}") + + self.update_weight(gamma, ood_query=ood_data, shift_query=shift_responses, trace=trace) + return new_round + + def validate_pcs_claim(self, first_round : WhirRound, trace = True): + folding_factor = 2**self.schk_rounds + current_round = first_round + + while current_round.RS_CODE.dimension() > folding_factor: + if trace: + print("\n----- New Sumcheck Round -----\n") + + if self.do_sumcheck(current_round, trace) is False: + print("ERROR: Sumcheck claim failed") + return False + + if trace: + print("\n----- Next OOD + STIR Round -----\n") + current_round = self.fold_check(current_round, self.shift_query_count, trace=trace) + + if trace: + print("\n----- Final Round -----\n") + + last_round_poly = current_round.poly + + if trace: + print(f"Final Sumcheck claim: {self.claimed_sum}") + print(f"Final Folded Poly: {last_round_poly}") + + sumcheck_poly = self.weight_polynomial.subs({self.weight_var : last_round_poly }) + + expected_sum = sumcheck_poly + + for v in sumcheck_poly.variables(): + this_round = self.Fq(0) + for h in self.hypercube: + this_round += expected_sum.subs( { v : h}) + expected_sum = this_round + + if trace: + print(f"Final evaluated sum: {expected_sum}. Expected claimed sum: {self.claimed_sum}") + + if expected_sum != self.claimed_sum: + print(f"PCS claim invalid") + return False + else: + print(f"It's a match! PCS claim is valid.") + return True + + +if __name__ == '__main__': + from proth_primes import proth_in_range + + def test_whir(variables_count = 10, rate_factor=8, fold_factor = 4, shift_query_count = 10, trace=True): + (p, cofactor, twodicity) = proth_in_range(15, 15, 27,27) + Fq = GF(p) + + + code_length = rate_factor*(2**variables_count) + + omega = Fq.multiplicative_generator() + omega = omega**(cofactor*2**(twodicity - log(code_length, 2))) + mult_order = omega.multiplicative_order() + mult_order_log = int(log(mult_order, 2)) + + MultilinearRing = gen_mvariate_ring(base_ring=Fq, variables_count=variables_count) + WeightRing = PolynomialRing(MultilinearRing, "Z") + input_poly = random_multilinear_poly(MultilinearRing) + + evaluation_point = [Fq.random_element() for _ in range(variables_count)] + weight_poly = eq_weight_polynomial(WeightRing, evaluation_point) + pcs_claim = input_poly(evaluation_point) + + print(f"TEST SETUP:") + print(f" Finite Field: {Fq} ({cofactor}⋅2^{twodicity} + 1)") + print(f" RS Generator: {omega}, multiplicative order: {mult_order} (2^{mult_order_log})") + print(f" Polynomial Variable: {input_poly.variables()}") + print(f" PCS Evaluation point: {evaluation_point}") + print(f" PCS Evaluation claim: {pcs_claim}") + print(f" Input polynomial : {input_poly}") + print(f" Weight Polynomial: {weight_poly}\n") + + wr = WhirRound(input_poly, weight_poly, omega, pcs_claim) + wv = WhirVerifier(weight_poly, omega, pcs_claim, fold_factor, shift_query_count) + + assert wv.validate_pcs_claim(wr, trace) + + + test_whir(variables_count=10) diff --git a/sage/fri-and-friends/src/zkwhir.py b/sage/fri-and-friends/src/zkwhir.py new file mode 100644 index 000000000..26940704c --- /dev/null +++ b/sage/fri-and-friends/src/zkwhir.py @@ -0,0 +1,235 @@ +import os +import sys +src_dir = os.path.dirname(os.path.abspath(__file__)) +sys.path.insert(0, src_dir) + +from sage.functions.other import ceil +from sage.misc.functional import log +from sage.rings.integer import Integer +from sage.coding.grs_code import ReedSolomonCode +from whir import * + + +def bit_decom(i, min_len): + b = Integer(i).bits() + pad_len = min_len - len(b) + if pad_len > 0: + b.extend([0 for _ in range(pad_len)]) + + return b + +def uni_to_multi(univar, MVarPoly): + poly_vars = MVarPoly.gens() + results_dict = {} + + for (i,b) in enumerate(univar.list()): + decompose = bit_decom(i, len(poly_vars)) + results_dict[tuple(decompose)] = b + + mpoly = MVarPoly(results_dict) + return mpoly + +def upgrade_ring(input_poly, OutputRing): + new_dict = {} + + for (k,v) in input_poly.dict().items(): + vv = list(k) + [0] + new_dict[tuple(vv)] = v + + return OutputRing(new_dict) + +# +# WARNING: This way of masking is probably not uniformly random. +# Finding a _multilinear_ polynomial mask such that its evaluation +# outside the boolean hypercube is uniformly random and zero within +# the boolean hypercube is WIP!!! +# +def mask_polynomial(input_poly, query_threshold, folding_factor): + mask_degree = query_threshold*folding_factor + input_degree = 2**len(input_poly.variables()) + nzc = 2**ceil(log(input_degree + mask_degree, 2)) - input_degree - 1 + + MaskPolyRing = PolynomialRing(input_poly.base_ring(), "T") + + # create a random univariate polynomial of degree d-1 + random_poly = MaskPolyRing.random_element(nzc) + + multivars = input_poly.parent().gens() + + FinalPolyRing = gen_mvariate_ring(input_poly.base_ring(), len(multivars) + 1) + last_var = FinalPolyRing.gens()[-1] + random_multivar = uni_to_multi(random_poly, FinalPolyRing) + + input_poly = upgrade_ring(input_poly, FinalPolyRing) + + output_poly = input_poly + random_multivar*last_var # Multiplied to the last variable + + if __debug__: + BR = output_poly.base_ring() + # print(f"Input polynomial: {input_poly}\n") + # print(f"Masked polynomial: {output_poly}\n") + # print(f"Delta: {output_poly - input_poly}\n") + + random_eq = eq_polynomial(FinalPolyRing, [BR.random_element() for _ in range(len(FinalPolyRing.gens()) -1)] + [0]) + old_sum = compute_hypercube_sum(input_poly * random_eq) + new_sum = compute_hypercube_sum(output_poly * random_eq) + assert old_sum == new_sum + assert input_poly == output_poly.subs( { last_var : 0 }) + + return output_poly + +class CommitmentData: + def __init__(self, + committed_poly : ReedSolomonCode, + scheck_mask : ReedSolomonCode, + PolyRing : PolynomialRing): + + self.PolyRing = PolyRing + self.func = committed_poly + self.scheck_mask = scheck_mask # The mask during the sumcheck protocol + + self.pcs_evaluation = None # Value of evaluation + self.scheck_evaluation = None # Value of evaluation of schek randomness + self.whir_first_round = None # First FRI WHIR oracle + + +class ZKWhirPCS: + + def __init__(self, + input_polynomial, + ntt_omega, + query_threshold = 3, + folding_factor = 4, + ntt_order_omega = None): + + self.input_poly = input_polynomial + self.Fq = input_polynomial.base_ring() + # Mask to hide information about the input polynomial. + self.masked_poly = mask_polynomial(self.input_poly, + query_threshold, + folding_factor) + print(f" Masked input: {self.masked_poly}\n") + self.MaskingRing = self.masked_poly.parent() + self.ProverWeightRing = PolynomialRing(self.masked_poly.parent(), "Z") + self.RS_RING = PolynomialRing(self.Fq, "h") + self.folding_factor = folding_factor + + # Mask used during the sumcheck round + self.scheck_mask = random_multilinear_poly(self.MaskingRing) + + print(f" Sumcheck Mask Polynomial: {self.masked_poly}") + + self.omega = self.Fq(ntt_omega) + self.multi_order = ntt_order_omega or ntt_omega.multiplicative_order() + self.dimension = 2**len(self.MaskingRing.gens()) + + + def commit(self) -> CommitmentData: + rs_code = ReedSolomonCode(self.Fq, self.multi_order, self.dimension, self.omega) + encoder = rs_code.encoder("EvaluationPolynomial", polynomial_ring=self.RS_RING) + fn_masked = WhirRound.to_univariate(self.masked_poly, self.RS_RING.gen()) + fn_scheck = WhirRound.to_univariate(self.scheck_mask, self.RS_RING.gen()) + + masked_oracle = encoder.encode(fn_masked) + scheck_oracle = encoder.encode(fn_scheck) + + return CommitmentData(masked_oracle, scheck_oracle, self.input_poly.parent()) + + def open(self, + evaluation_point : list[FiniteField], + verifier_challenge : FiniteField, + commitment : CommitmentData): + + commitment.pcs_evaluation = self.input_poly(evaluation_point) + + # + # Set the last variable in the the weight polynomial to zero. This undoes + # the effect of masking. Note that the challenger does not get to select + # this last co-ordinate. + # + + extended_eval_point = evaluation_point + [0] + + assert self.masked_poly(extended_eval_point) == commitment.pcs_evaluation + commitment.scheck_evaluation = self.scheck_mask(extended_eval_point) + + # + # Sumcheck weight polynomial has last co-ordinate set to zero. + # + weight_poly = eq_weight_polynomial(self.ProverWeightRing, extended_eval_point) + + # + # First round of oracle happens over the virtual oracle + virtual_oracle = verifier_challenge*commitment.func + commitment.scheck_mask + + # + # Virtual sumcheck polynomial is + # verifier_challenge * masked_polynomial + sumcheck_polynomial + # + + sumcheck_poly = verifier_challenge*self.masked_poly + self.scheck_mask + + claimed_sum = verifier_challenge*commitment.pcs_evaluation + \ + commitment.scheck_evaluation + + # print(f"P: Prover's claimed sum: {claimed_sum}") + + assert sumcheck_poly(extended_eval_point) == claimed_sum + + first_round = WhirRound(sumcheck_poly, weight_poly, self.omega, claimed_sum) + + assert first_round.fri_oracle() == virtual_oracle, "The virtual oracle doesn't match with expected codeword of sumcheck polynomial" + + # first_round.set_fri_oracle(virtual_oracle) + commitment.whir_first_round = first_round + + return commitment + + def verify(self, evaluation_point, verifier_challenge, commitment_data): + claimed_value = commitment_data.pcs_evaluation + random_sum = commitment_data.scheck_evaluation + whir_prover = commitment_data.whir_first_round + + combined_sum = verifier_challenge*claimed_value + random_sum + + # The verifier weight polynomial is just as before + weight_poly = eq_weight_polynomial(self.ProverWeightRing, evaluation_point + [0]) + + verifier = WhirVerifier(weight_poly, + self.omega, + combined_sum, + self.folding_factor) + + return verifier.validate_pcs_claim(whir_prover) + + +if __name__ == '__main__': + from proth_primes import proth_in_range + + def test_zk_whir(variables_count = 10, rate_factor=8, fold_factor = 4, shift_query_count = 10, trace=True): + (p, cofactor, twodicity) = proth_in_range(15, 15, 27,27) + Fq = GF(p) + + code_length = rate_factor*(2**(variables_count+1)) # Accommodate for masking + + omega = Fq.multiplicative_generator() + omega = omega**(cofactor*2**(twodicity - log(code_length, 2))) + + MultilinearRing = gen_mvariate_ring(base_ring=Fq, variables_count=variables_count) + + input_poly = random_multilinear_poly(MultilinearRing) + evaluation_point = [Fq.random_element() for _ in range(variables_count)] + verifier_challenge = Fq.random_element() + + print(f"Test setup:") + print(f" Input polynomial: {input_poly}") + print(f" Evaluation point: {evaluation_point}") + print(f" Verifier's Sigma ZK Challenge: {verifier_challenge}") + + zkWHIR = ZKWhirPCS(input_poly, ntt_omega=omega) + commitment = zkWHIR.commit() + commitment = zkWHIR.open(evaluation_point, verifier_challenge, commitment) + zkWHIR.verify(evaluation_point, verifier_challenge, commitment) + + test_zk_whir(variables_count=10, fold_factor = 4) + diff --git a/sage/fri-and-friends/whir.ipynb b/sage/fri-and-friends/whir.ipynb new file mode 100644 index 000000000..29dda7260 --- /dev/null +++ b/sage/fri-and-friends/whir.ipynb @@ -0,0 +1,427 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "# WHIR: Reed–Solomon Proximity Testing with Super-Fast Verification" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "[WHIR](https://eprint.iacr.org/2024/1586) is a multivariate Polynomial Commitment Scheme (PCS) with three main ideas combined into one:\n", + "\n", + "1. Using Sumcheck as a _proof_ of multilinear polynomial commitment. \n", + " Given, a multilinear polynomial $f(x_1, \\cdots, x_\\mu)$, it's \n", + " evaluation at a point $\\vec{a}$ can be written as \n", + " $f(\\vec{a}) = \\sum_{\\vec{b}\\in \\{0,1\\}^\\mu} f(\\vec{b})\\cdot \\mathsf{eq}(\\vec{b},\\vec{a})$.\n", + " The Sumcheck protocol _reduces_ a claim of this form to a \n", + " claim about evaluating $f(\\vec{x})$ at a random point \n", + " $\\vec{r}$ of verifier's choosing. But this causes a problem: \n", + " If the verifier could not evaluate $f(\\vec{a})$ in the first place, \n", + " why would it be able to evaluate $f(\\vec{r})$? \n", + "\n", + " The first main idea in WHIR (which first appeared in \n", + " [BaseFold](https://eprint.iacr.org/2023/1705.pdf)) \n", + " is that instead of reducing the sumcheck claim to a _random evaluation_, \n", + " one can **reduce the sumcheck claim to a Proof of Proximity** and a final evaluation. \n", + " To achieve this, WHIR intermixes sumcheck rounds with FRI rounds, and relies on \n", + " consistency checks across FRI proof oracles to ensure soundness \n", + " across rounds. By the last round, the FRI proof oracle is small enough that\n", + " the sumcheck value can be directly computed and checked against expected value. \n", + " This is the gist of BaseFold/WHIR reduction.\n", + "\n", + "3. The second idea in WHIR is to use lower rate codes to minimize FRI query \n", + " complexity (this idea first appeared in [STIR](https://eprint.iacr.org/2024/390.pdf) paper).\n", + " Unlike in STIR, however, WHIR elegantly avoids division operations during inter-round \n", + " consistency checks by simply updating the sumcheck claim on prover and verifier side. \n", + "\n", + "4. The third idea is to use out of domain (OOD) queries to eliminate \"pretenders.\" This allows \n", + " WHIR to achieve better soundness even in in list-decoding region. This idea first \n", + " appeared in [DEEP-FRI](https://www.math.toronto.edu/swastik/deep-fri.pdf) paper." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The following code describes a run of the protocol. All the SageMath WHIR code is located in [src/whir.py](./src/whir.py).\n", + "\n", + "File [src/proth_primes.py](./src/proth_primes.py) has a long list of [proth prime](https://en.wikipedia.org/wiki/Proth_prime), which are prime numbers of the form $k\\cdot2^n + 1$ for different values of $k$ and $n$. These primes are convenient for NTT operations." + ] + }, + { + "cell_type": "code", + "execution_count": 13, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "The autoreload extension is already loaded. To reload it, use:\n", + " %reload_ext autoreload\n" + ] + } + ], + "source": [ + "%load_ext autoreload\n", + "%autoreload 2\n", + "\n", + "from src.proth_primes import proth_in_range\n", + "from src.whir import *\n", + "from sage.rings.finite_rings.all import GF\n", + "from sage.rings.polynomial.all import PolynomialRing\n", + "from sage.misc.functional import log" + ] + }, + { + "cell_type": "code", + "execution_count": 14, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Using prime 2013265921 = 15⋅2^27 + 1\n", + "Maximal 2-sylow multiplicative generator: 440564289\n", + "Maximal 2-sylow multiplicative power: 27\n" + ] + } + ], + "source": [ + "(p, k, n) = proth_in_range(15, 15, 27,27)\n", + "print(f\"Using prime {p} = {k}⋅2^{n} + 1\")\n", + "Fq = GF(p)\n", + "omega = Fq.multiplicative_generator() \n", + "omega = omega**k\n", + "omega_2sylow_expo = log(omega.multiplicative_order(), 2)\n", + "print(f\"Maximal 2-sylow multiplicative generator: {omega}\")\n", + "print(f\"Maximal 2-sylow multiplicative power: {omega_2sylow_expo}\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "---\n", + "### Setup\n", + "\n", + "Creating WHIR sumcheck with three variables and using a random multilinear polynomial." + ] + }, + { + "cell_type": "code", + "execution_count": 15, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Using polynomial: X₀*X₁*X₂ - 2*X₀*X₁ + 3*X₀*X₂ - 4*X₁*X₂ + 5*X₀ - 6*X₂ + 7\n", + "Evaluation point: [5, 6, 7]\n", + "Evaluated value: 77\n" + ] + } + ], + "source": [ + "variables_count = 3\n", + "code_rate_factor = 8\n", + "sumcheck_rounds = 2\n", + "shift_queries = 5\n", + "\n", + "PolyRing = gen_mvariate_ring(Fq, variables_count=variables_count)\n", + "WeightRing = PolynomialRing(PolyRing, \"Z\") \n", + "\n", + "X0, X1, X2 = PolyRing.gens()\n", + "Z = WeightRing.gen()\n", + "\n", + "input_poly = X0*X1*X2 - 2*X0*X1 + 3*X0*X2 - 4*X1*X2 + 5*X0 - 6*X2 + 7 # random_multilinear_poly(PolyRing);\n", + "evaluation_point = [5, 6, 7] # [Fq.random_element() for _ in range(variables_count)]\n", + "\n", + "pcs_evaluation_claim = input_poly(evaluation_point)\n", + "\n", + "print(f\"Using polynomial: {input_poly}\")\n", + "print(f\"Evaluation point: {evaluation_point}\")\n", + "print(f\"Evaluated value: {input_poly(evaluation_point)}\")" + ] + }, + { + "cell_type": "code", + "execution_count": 16, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Using NTT group generator: 1721589904, order: 64\n" + ] + } + ], + "source": [ + "ntt_order = code_rate_factor*(2**variables_count)\n", + "omega_order = omega.multiplicative_order()\n", + "assert omega_order % ntt_order == 0\n", + "coset_exponent = omega_order // ntt_order\n", + "ntt_omega = omega**coset_exponent\n", + "print(f\"Using NTT group generator: {ntt_omega}, order: {ntt_omega.multiplicative_order()}\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Create the weight polynomial for the given evaluation point. This evaluation point will typically come from the verifier, and is known to the verifier. " + ] + }, + { + "cell_type": "code", + "execution_count": 17, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Prover's Weight polynomial: (1287*X₀*X₁*X₂ - 594*X₀*X₁ - 585*X₀*X₂ - 572*X₁*X₂ + 270*X₀ + 264*X₁ + 260*X₂ - 120)*Z\n" + ] + } + ], + "source": [ + "weight_poly = eq_weight_polynomial(WeightRing, evaluation_point)\n", + "print(f\"Prover's Weight polynomial: {weight_poly}\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Verifier (positive case)" + ] + }, + { + "cell_type": "code", + "execution_count": 18, + "metadata": { + "vscode": { + "languageId": "python" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "Sumcheck polynomial: 1287*X₀^2*X₁^2*X₂^2 - 3168*X₀^2*X₁^2*X₂ + 3276*X₀^2*X₁*X₂^2 - 5720*X₀*X₁^2*X₂^2 + 1188*X₀^2*X₁^2 + 6093*X₀^2*X₁*X₂ + 3784*X₀*X₁^2*X₂ - 1755*X₀^2*X₂^2 - 6838*X₀*X₁*X₂^2 + 2288*X₁^2*X₂^2 - 3510*X₀^2*X₁ - 528*X₀*X₁^2 - 2115*X₀^2*X₂ + 8785*X₀*X₁*X₂ - 1056*X₁^2*X₂ + 4290*X₀*X₂^2 + 2392*X₁*X₂^2 + 1350*X₀^2 - 2598*X₀*X₁ - 4775*X₀*X₂ - 5108*X₁*X₂ - 1560*X₂^2 + 1290*X₀ + 1848*X₁ + 2540*X₂ - 840\n" + ] + } + ], + "source": [ + "prover = WhirRound(input_poly, weight_poly, ntt_omega, pcs_evaluation_claim)\n", + "print(f\"Sumcheck polynomial: {prover.sumcheck_poly}\")" + ] + }, + { + "cell_type": "code", + "execution_count": 19, + "metadata": { + "editable": true, + "scrolled": true, + "slideshow": { + "slide_type": "" + }, + "tags": [], + "vscode": { + "languageId": "python" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "\n", + "----- New Sumcheck Round -----\n", + "\n", + "Current Sumcheck claim: 77\n", + "Round-0 poly: 504*X₀^2 - 2051*X₀ + 812, challenge: None\n", + "Hypercube sum 77 matches claimed sum 77\n", + "Round-1 poly: 710417221*X₁^2 + 933868125*X₁ - 122024773, challenge: 913189950\n", + "Hypercube sum -613030121 matches claimed sum -613030121\n", + "Updated Sumcheck claim: -613030121\n", + "\n", + "----- Next OOD + STIR Round -----\n", + "\n", + "Current Sumcheck claim: -613030121\n", + "OOD challenge: 1961264585\n", + "Shift queries: [1835875777, 1253260071]\n", + "Combination randomness 𝛾: 1576224895\n", + "\n", + " P> Weight update for OOD challenge 1961264585 => 424934549 : Z*(-70261355*X₁*X₂ - 16870659*X₁ + 726465213*X₂ - 337231938)\n", + " P> Weight update for shift query 1835875777 : Z*(633836656*X₁*X₂ + 512324488*X₁ + 76728129*X₂ + 50331008)\n", + " P> Weight update for shift query 1253260071 : Z*(191116099*X₁*X₂ - 855563900*X₁ + 100838210*X₂ - 677049140)\n", + " P> Updated new weight: (417065636*X₁*X₂ - 463474140*X₁ + 745632449*X₂ - 226734683)*Z\n", + " P> Sumcheck claim combination randomness additional sum : 1291366824\n", + " P> Updated Sumcheck claim: 678336703\n", + "\n", + "Prover OOD response: 424934549\n", + "Verifier computed shift response: [(1835875777, 2002511067), (1253260071, 1846892102)]\n", + "\n", + " V> Weight update for OOD challenge 1961264585 => 424934549 : Z*(-70261355*X₁*X₂ - 16870659*X₁ + 726465213*X₂ - 337231938)\n", + " V> Weight update for shift query 1835875777 => 2002511067 : Z*(633836656*X₁*X₂ + 512324488*X₁ + 76728129*X₂ + 50331008)\n", + " V> Weight update for shift query 1253260071 => 1846892102 : Z*(191116099*X₁*X₂ - 855563900*X₁ + 100838210*X₂ - 677049140)\n", + " V> Updated new weight: (417065636*X₁*X₂ - 463474140*X₁ + 745632449*X₂ - 226734683)*Z\n", + " V> Sumcheck claim combination randomness additional sum : 1291366824\n", + " V> Updated Sumcheck claim: 678336703\n", + "\n", + "----- New Sumcheck Round -----\n", + "\n", + "Current Sumcheck claim: 678336703\n", + "Round-0 poly: 1001931861*X₁^2 - 937798515*X₁ - 699531282, challenge: None\n", + "Hypercube sum 678336703 matches claimed sum 678336703\n", + "Round-1 poly: -675808816*X₂^2 - 827779045*X₂ - 367190350, challenge: 312420118\n", + "Hypercube sum -224702640 matches claimed sum -224702640\n", + "Updated Sumcheck claim: -224702640\n", + "\n", + "----- Next OOD + STIR Round -----\n", + "\n", + "Current Sumcheck claim: -224702640\n", + "OOD challenge: 323854840\n", + "Shift queries: [211723194, 78945800]\n", + "Combination randomness 𝛾: 1947755793\n", + "\n", + " P> Weight update for OOD challenge 323854840 => 824441233 : Z*(647709679*X₂ - 323854839)\n", + " P> Weight update for shift query 211723194 : Z*(423446387*X₂ - 211723193)\n", + " P> Weight update for shift query 78945800 : Z*(157891599*X₂ - 78945799)\n", + " P> Updated new weight: (-460563677*X₂ - 405589246)*Z\n", + " P> Sumcheck claim combination randomness additional sum : 1140985209\n", + " P> Updated Sumcheck claim: 916282569\n", + "\n", + "Prover OOD response: 824441233\n", + "Verifier computed shift response: [(211723194, 1085228624), (78945800, 763160990)]\n", + "\n", + " V> Weight update for OOD challenge 323854840 => 824441233 : Z*(647709679*X₂ - 323854839)\n", + " V> Weight update for shift query 211723194 => 1085228624 : Z*(423446387*X₂ - 211723193)\n", + " V> Weight update for shift query 78945800 => 763160990 : Z*(157891599*X₂ - 78945799)\n", + " V> Updated new weight: (-460563677*X₂ - 405589246)*Z\n", + " V> Sumcheck claim combination randomness additional sum : 1140985209\n", + " V> Updated Sumcheck claim: 916282569\n", + "\n", + "----- Final Round -----\n", + "\n", + "Final Sumcheck claim: 916282569\n", + "Final Folded Poly: 353556209*X₂ + 798818320\n", + "Final evaluated sum: 916282569. Expected claimed sum: 916282569\n", + "It's a match! PCS claim is valid.\n", + "PCS verification succeeded\n" + ] + } + ], + "source": [ + "verifier = WhirVerifier(weight_poly, ntt_omega, pcs_evaluation_claim)\n", + "succes = verifier.validate_pcs_claim(prover)\n", + "if succes == True:\n", + " print(f\"PCS verification succeeded\")\n", + "else:\n", + " print(f\"PCS verification failed\")" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Verifier (negative case)" + ] + }, + { + "cell_type": "code", + "execution_count": 20, + "metadata": { + "collapsed": true, + "jupyter": { + "outputs_hidden": true + }, + "vscode": { + "languageId": "python" + } + }, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "\n", + "----- New Sumcheck Round -----\n", + "\n", + "Current Sumcheck claim: 78\n", + "Round-0 poly: 504*X₀^2 - 2051*X₀ + 812, challenge: None\n", + "ERROR: Sumcheck claim failed\n", + "PCS verification failed\n" + ] + } + ], + "source": [ + "prover = WhirRound(input_poly, weight_poly, ntt_omega, pcs_evaluation_claim)\n", + "\n", + "verifier = WhirVerifier(weight_poly, ntt_omega, pcs_evaluation_claim + 1)\n", + "succes = verifier.validate_pcs_claim(prover)\n", + "if succes == True:\n", + " print(f\"PCS verification succeeded\")\n", + "else:\n", + " print(f\"PCS verification failed\")" + ] + } + ], + "metadata": { + "kernelspec": { + "display_name": "SageMath 10.5", + "language": "sage", + "name": "sagemath-10.5" + }, + "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.12.5" + } + }, + "nbformat": 4, + "nbformat_minor": 4 +}