Add more experiments into EPA notebook.

1 files changed, 190 insertions, 33 deletions

diff --git a/re/epa.ipynb b/re/epa.ipynb
index afafa62..f36d829 100644
--- a/re/epa.ipynb
+++ b/re/epa.ipynb
@@ -17,18 +17,24 @@
    "source": [
     "import io\n",
     "import tabulate\n",
+    "import secrets\n",
+    "from tqdm.notebook import tqdm, trange\n",
     "from itertools import product\n",
     "from IPython.display import HTML, display\n",
+    "from sympy.ntheory import factorint\n",
+    "from sympy.ntheory.modular import crt\n",
     "\n",
     "from pyecsca.ec.model import ShortWeierstrassModel\n",
     "from pyecsca.ec.coordinates import AffineCoordinateModel\n",
     "from pyecsca.ec.curve import EllipticCurve\n",
     "from pyecsca.ec.params import DomainParameters, load_params_ectester\n",
     "from pyecsca.ec.mod import Mod, miller_rabin, gcd\n",
-    "from pyecsca.ec.point import Point\n",
+    "from pyecsca.ec.point import Point, InfinityPoint\n",
     "from pyecsca.ec.error import NonInvertibleError\n",
     "from pyecsca.ec.mult import LTRMultiplier, AccumulationOrder\n",
+    "from pyecsca.ec.formula.fake import FakeAdditionFormula, FakeDoublingFormula, FakePoint\n",
     "from pyecsca.sca.re.rpa import MultipleContext\n",
+    "from pyecsca.sca.re.zvp import unroll_formula_expr\n",
     "from pyecsca.ec.context import local\n",
     "from pyecsca.ec.error import UnsatisfiedAssumptionError"
    ]
@@ -48,8 +54,12 @@
     "    \"cca6f6718a06cad7094962b2a35f067d,67aa9464eb493fbb7b509d29381b9a9d,cafc69aa517b654a6a608644996cc8d1,4c092beb06cc00751eec39675f680cb8,82800378a47dd6f26ff6a50f69e4c4e6,18a22d20b6de3ff6bdc49329c21163f77,0001\",\n",
     "    # 0.8508806646440022 \n",
     "    \"b3755d654bad73114e4191e9f5f36af9,9fe4f88cfbacba71f4b767ace8580c74,4610526fdcfbd69aed453ac2ee6efeef,542d8e0bbafe40dae36f25cbc350add6,68a65f5a5dc304bfd0d8fe963c250206,118a34a1ea295e78b3a3c960b6f680ee1,0001\",\n",
-    "    # 0.9845701775215489    \n",
+    "    \n",
+    "    # 0.9845701775215489 (has a = 0 for a subcurve)\n",
     "    \"de1406450d5d7e91d81907956019c0c1,5fbe46b9f1086011e18f5d823c6110ce,a859c36ceeadb39c7a978f7b1b0563ee,1cba89c3f099c29401ecf3fe1806e822,345d7282a0114070be91f95fe3db1faa,0fcd24d24e57a40547814b6766b9ea735,0001\",\n",
+    "    # 0.980582605794486  (does not have a = 0 for any subcurve)\n",
+    "    \"cab298b495875d4ab2c8ee3eb03016a7,a7c4f56f286d9eae44424c85c8b2fcb9,5e8c439d939273fdcb5503acbda7d3f8,816c9f865c831223067a88046bf00d75,972ce29ed18d5d73f15cef31187659be,0b0e97ff8c3e72e7ae75eb3f5e759fe03,0001\",\n",
+    "    \n",
     "    # 0.9547100843537808\n",
     "    \"f1a8a441b6d0e9600e33ccf16f9b8291,b3f55185bd6a63528e3d560c6a7b729a,c2fee2d65350e870eda0ac5e2b96b810,29b3e793822fad03a3c2ebca3cf62c12,b937d5389b6c5d0212d0f53e26843092,1153442389f9e1da8dd130bc93c6ef42b,0001\",\n",
     "    # 0.7214369438844093\n",
@@ -71,8 +81,9 @@
     "which = \"projective\"\n",
     "coords = model.coordinates[which]\n",
     "\n",
-    "params = load_params_ectester(io.BytesIO(curves[3].encode()), which)\n",
+    "params = load_params_ectester(io.BytesIO(curves[4].encode()), which)\n",
     "curve = params.curve\n",
+    "p = params.curve.prime\n",
     "g = params.generator\n",
     "n = params.order"
    ]
@@ -86,7 +97,116 @@
    "source": [
     "adds = list(filter(lambda formula: formula.name.startswith(\"add\"), coords.formulas.values()))\n",
     "dbls = list(filter(lambda formula: formula.name.startswith(\"dbl\"), coords.formulas.values()))\n",
-    "formula_pairs = list(product(adds, dbls))"
+    "formula_pairs = list(product(adds, dbls))\n",
+    "\n",
+    "fake_add = FakeAdditionFormula(params.curve.coordinate_model)\n",
+    "fake_dbl = FakeDoublingFormula(params.curve.coordinate_model)\n",
+    "fake_mult = LTRMultiplier(fadd, fdbl, None, False, AccumulationOrder.PeqPR, True, True)\n",
+    "fake_mult.init(params, FakePoint(params.curve.coordinate_model))"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "21ad6946-bb69-40b8-9407-d6c78537db47",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def random_scalar(n):\n",
+    "    \"\"\"Generate a random scalar mod n.\"\"\"\n",
+    "    return secrets.randbelow(n)\n",
+    "\n",
+    "def random_scalar_trivial(n):\n",
+    "    \"\"\"Generate a random scalar with trivial gcd mod n.\"\"\"\n",
+    "    scalar = secrets.randbelow(n)\n",
+    "    while gcd(scalar, n) != 1:\n",
+    "        scalar = secrets.randbelow(n)\n",
+    "    return scalar\n",
+    "\n",
+    "def random_scalar_fully_trivial(n, mult):\n",
+    "    \"\"\"Generate a random scalar with trivial gcd mod n, and also ensure that the given mult computes only multiples with trivial gcd mod n.\"\"\"\n",
+    "    scalar = random_scalar_trivial(n)\n",
+    "    while True:\n",
+    "        with local(MultipleContext()) as ctx:\n",
+    "            mult.multiply(scalar)\n",
+    "        if all(map(lambda x: gcd(x, n) == 1, ctx.points.values())):\n",
+    "            return scalar\n",
+    "        scalar = random_scalar_trivial(n)\n",
+    "\n",
+    "def fixed_point(params):\n",
+    "    return params.generator\n",
+    "\n",
+    "def random_point(splitted, top, randomized=False):\n",
+    "    results = {}\n",
+    "    for factor, params in splitted.items():\n",
+    "        results[factor] = params.curve.affine_random()\n",
+    "    factors = list(results.keys())\n",
+    "    xs = list(map(lambda factor: int(results[factor].x), factors))\n",
+    "    ys = list(map(lambda factor: int(results[factor].y), factors))\n",
+    "    res_x = Mod(int(crt(factors, xs)[0]), top.curve.prime)\n",
+    "    res_y = Mod(int(crt(factors, ys)[0]), top.curve.prime)\n",
+    "    res = Point(affine, x=res_x, y=res_y)\n",
+    "    return res.to_model(top.curve.coordinate_model, top.curve, randomized=randomized)"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "b6e888d3-38a0-4028-add5-a7b428b8b6cd",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "def project_down(point, subcurve):\n",
+    "    return Point(subcurve.coordinate_model, **{name: Mod(int(value), subcurve.prime) for name, value in point.coords.items()})\n",
+    "\n",
+    "def lift_up(point, topcurve):\n",
+    "    return Point(topcurve.coordinate_model, **{name: Mod(int(value), topcurve.prime) for name, value in point.coords.items()})\n",
+    "\n",
+    "def split_params(params):\n",
+    "    factors = factorint(params.curve.prime)\n",
+    "    if set(factors.values()) != {1}:\n",
+    "        raise ValueError(\"Not squarefree\")\n",
+    "    results = {}\n",
+    "    # Construct the curves\n",
+    "    for factor in sorted(factors.keys()):\n",
+    "        p_i = factor\n",
+    "        parameters_i = {name: Mod(int(value), p_i) for name, value in params.curve.parameters.items()}\n",
+    "        curve_i = EllipticCurve(params.curve.model, params.curve.coordinate_model, p_i, params.curve.neutral, parameters_i)\n",
+    "        generator_i = project_down(params.generator, curve_i)\n",
+    "        params_i = DomainParameters(curve_i, generator_i, 0, 1)\n",
+    "        results[factor] = params_i\n",
+    "    # Now map the orders to the curves\n",
+    "    orders = list(factorint(params.order).keys())\n",
+    "    orders.sort()\n",
+    "    for factor_i, params_i in results.items():\n",
+    "        for order in orders:\n",
+    "            try:\n",
+    "                params_i.curve.affine_multiply(params_i.generator.to_affine(), order)\n",
+    "            except NonInvertibleError:\n",
+    "                params_i.order = order\n",
+    "                orders.remove(order)\n",
+    "                break\n",
+    "    return results\n",
+    "\n",
+    "def split_scalarmult(splitted, top, point, scalar):\n",
+    "    results = {}\n",
+    "    for factor, params in splitted.items():\n",
+    "        order = params.order\n",
+    "        projected = project_down(point, params.curve)\n",
+    "        partial_scalar = scalar % order\n",
+    "        if partial_scalar == 0:\n",
+    "            result = InfinityPoint(params.curve.coordinate_model)\n",
+    "        else:\n",
+    "            result = params.curve.affine_multiply(projected.to_affine(), partial_scalar)\n",
+    "        results[factor] = result\n",
+    "    if any(map(lambda point: isinstance(point, InfinityPoint), results.values())):\n",
+    "        return InfinityPoint(top.curve.coordinate_model)\n",
+    "    factors = list(results.keys())\n",
+    "    xs = list(map(lambda factor: int(results[factor].x), factors))\n",
+    "    ys = list(map(lambda factor: int(results[factor].y), factors))\n",
+    "    res_x = Mod(int(crt(factors, xs)[0]), top.curve.prime)\n",
+    "    res_y = Mod(int(crt(factors, ys)[0]), top.curve.prime)\n",
+    "    return Point(affine, x=res_x, y=res_y)"
    ]
   },
   {
@@ -96,40 +216,60 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "scalars = [int(Mod.random(n)) for _ in range(200)]\n",
+    "split = split_params(params)\n",
+    "scalars = [random_scalar_trivial(n) for _ in trange(50, desc=\"Generate scalars\")]\n",
+    "points = [random_point(split, params, randomized=True) for _ in trange(50, desc=\"Generate points\")]\n",
+    "results = []\n",
+    "chains = []\n",
+    "\n",
+    "\n",
+    "gcds = []\n",
+    "fgcds = []\n",
+    "for scalar, point in tqdm(zip(scalars, points), desc=\"Precomp\", total=len(scalars)):\n",
+    "    try:\n",
+    "        result = split_scalarmult(split, params, point, scalar)\n",
+    "    except NonInvertibleError:\n",
+    "        result = None\n",
+    "    results.append(result)\n",
+    "    with local(MultipleContext()) as ctx:\n",
+    "        fake_mult.multiply(scalar)\n",
+    "    chains.append(list(ctx.points.values()))\n",
+    "    scalar_trivial_gcd = gcd(scalar, n) == 1\n",
+    "    all_subscalars_trivial_gcd = all(map(lambda x: gcd(x, n) == 1, ctx.points.values()))\n",
+    "    gcds.append(scalar_trivial_gcd)\n",
+    "    fgcds.append(all_subscalars_trivial_gcd)\n",
+    "\n",
     "\n",
-    "table = [[\"Pair\", \"scalars with trivial gcd\", \"scalars with all multiples with trivial gcds\", \"scalars with invertible final zs\", \"scalars with all multiples's zs invertible\"]]\n",
+    "table = [[\"Pair\", \"scalars with trivial gcd\", \"scalars with all multiples with trivial gcds\", \"scalars with invertible final zs\", \"scalars with all multiples's zs invertible\", \"scalars with correct result\"]]\n",
     "pair_table = [[None for _ in dbls] for _ in adds]\n",
-    "for pair in formula_pairs:\n",
+    "for pair in tqdm(formula_pairs):\n",
     "    mult = LTRMultiplier(*pair, None, False, AccumulationOrder.PeqPR, True, True)\n",
-    "    mult.init(params, g)\n",
-    "    inv = 0\n",
-    "    gcds = 0\n",
-    "    fgcds = 0\n",
-    "    zs = 0\n",
-    "    for scalar in scalars:\n",
-    "        if gcd(scalar, n) == 1:\n",
-    "            gcds += 1\n",
+    "    inv = []\n",
+    "    correct = []\n",
+    "    zs = []\n",
+    "    for scalar, point, result in tqdm(zip(scalars, points, results), leave=None, total=len(scalars)):\n",
+    "        mult.init(params, point)\n",
     "        with local(MultipleContext()) as ctx:\n",
     "            res = mult.multiply(scalar)\n",
-    "        if all(map(lambda x: gcd(x, n) == 1, ctx.points.values())):\n",
-    "            fgcds += 1\n",
-    "        if all(map(lambda x: gcd(int(x.Z), n) == 1, ctx.points.keys())):\n",
-    "            zs += 1\n",
-    "    \n",
+    "        \n",
+    "        all_submultiples_invertible_z = all(map(lambda x: gcd(int(x.Z), p) == 1, ctx.points.keys()))\n",
+    "        result_invertible_z = False\n",
+    "        result_correct = False\n",
     "        try:\n",
-    "            res.Z.inverse()\n",
-    "            inv += 1\n",
+    "            res_aff = res.to_affine()\n",
+    "            result_invertible_z = True\n",
+    "            if res_aff == result:\n",
+    "                result_correct = True\n",
     "        except NonInvertibleError as e:\n",
     "            pass\n",
-    "    #print(\"---\")\n",
-    "    #print(pair[0].name, pair[1].name)\n",
-    "    #print(gcds, \"\\tscalars with trivial gcd\")\n",
-    "    #print(fgcds, \"\\tscalars with all multiples with trivial gcds\")\n",
-    "    #print(inv, \"\\tscalars with invertible final zs\")\n",
-    "    #print(zs, \"\\tscalars with all multiples's zs invertible\")\n",
-    "    pair_table[adds.index(pair[0])][dbls.index(pair[1])] = inv\n",
-    "    table.append([f\"{pair[0].name}, {pair[1].name}\", gcds, fgcds, inv, zs])\n",
+    "        zs.append(all_submultiples_invertible_z)\n",
+    "        inv.append(result_invertible_z)\n",
+    "        correct.append(result_correct)\n",
+    "    pair_table[adds.index(pair[0])][dbls.index(pair[1])] = sum(inv)\n",
+    "    for i in inv:\n",
+    "        print(\"x\" if i else \".\", end=\"\")\n",
+    "    print()\n",
+    "    table.append([f\"{pair[0].name}, {pair[1].name}\", sum(gcds), sum(fgcds), sum(inv), sum(zs), sum(correct)])\n",
     "for pl, add in zip(pair_table, adds):\n",
     "    pl.insert(0, add.name)\n",
     "pair_table.insert(0, [None] + [dbl.name for dbl in dbls])\n",
@@ -220,17 +360,34 @@
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "dd3e076c-cd52-47e1-beb0-2cb59812dacc",
+   "id": "0c8a43dc-1a2f-4135-b83c-6ac84929789f",
+   "metadata": {},
+   "outputs": [],
+   "source": [
+    "split = split_params(params)\n",
+    "split"
+   ]
+  },
+  {
+   "cell_type": "code",
+   "execution_count": null,
+   "id": "d6d8bb2a-2e5d-475e-b775-5aff54eeaf47",
    "metadata": {},
    "outputs": [],
    "source": [
-    "c"
+    "for add in adds:\n",
+    "    r = None\n",
+    "    for iv in unroll_formula_expr(add):\n",
+    "        if iv[0] == \"Z3\":\n",
+    "            r = iv[1]\n",
+    "    print(add, r)\n",
+    "    print(\"---\")"
    ]
   },
   {
    "cell_type": "code",
    "execution_count": null,
-   "id": "aa548870-4236-42d8-a4a6-b6b8082b87bb",
+   "id": "64202cc7-a8f8-4d88-a0b5-5cf4d483a547",
    "metadata": {},
    "outputs": [],
    "source": []