diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml
index a63ba8b87..81eaed341 100644
--- a/.github/workflows/ci.yml
+++ b/.github/workflows/ci.yml
@@ -262,7 +262,7 @@ jobs:
         run: |
           pacman -S --needed --noconfirm mingw-w64-x86_64-gmp mingw-w64-x86_64-llvm
           nimble refresh --verbose -y
-          nimble install --verbose -y gmp@#head jsony asynctools yaml@1.1.0
+          nimble install --verbose -y gmp@#head jsony asynctools yaml@1.1.0 cliche
 
           cd constantine
           go mod download -modfile=go_test.mod
@@ -272,7 +272,7 @@ jobs:
         shell: bash
         run: |
           nimble refresh --verbose -y
-          nimble install --verbose -y gmp@#head jsony asynctools yaml@1.1.0
+          nimble install --verbose -y gmp@#head jsony asynctools yaml@1.1.0 cliche
 
           cd constantine
           go mod download -modfile=go_test.mod
@@ -337,6 +337,24 @@ jobs:
           cd constantine
           CTT_ASM=0 cargo test -- --nocapture
 
+      - name: Compile Constantine Zkalc benchmark (no assembly)
+        # Skip 32-bit as that would need clang-multilib or -m32
+        # Skip Windows as clang throws fatal error LNK1107
+        if: matrix.target.BACKEND == 'NO_ASM' && matrix.target.cpu != 'i386' && runner.os != 'Windows'
+        shell: bash
+        run: |
+          cd constantine
+          CTT_ASM=0 nimble make_zkalc
+
+      - name: Compile Constantine Zkalc benchmark (with assembly)
+        # Skip 32-bit as that would need clang-multilib or -m32
+        # Skip Windows as clang throws fatal error LNK1107
+        if: matrix.target.BACKEND == 'ASM' && matrix.target.cpu != 'i386' && runner.os != 'Windows'
+        shell: bash
+        run: |
+          cd constantine
+          nimble make_zkalc
+
       - name: Run Constantine in-depth tests (Unix - with GMP, with Assembly)
         if: runner.os != 'Windows' && matrix.target.BACKEND == 'ASM'
         shell: bash
diff --git a/.gitignore b/.gitignore
index b17171963..e3e179a31 100644
--- a/.gitignore
+++ b/.gitignore
@@ -1,4 +1,4 @@
-# Executables shall be put in an ignored build/ directory
+# Build artifacts shall be put in an ignored build/ directory
 # Ignore dynamic, static libs and libtool archive files
 # -----------------------------------------------------------------------------------------
 build/
@@ -54,4 +54,5 @@ Cargo.lock
 # Perf artifacts
 # -----------------------------------------------------------------------------------------
 perf.data
-perf.data.old
\ No newline at end of file
+perf.data.old
+*.bench.json
diff --git a/benchmarks/bench_elliptic_parallel_template.nim b/benchmarks/bench_elliptic_parallel_template.nim
index d90154061..afc0a649c 100644
--- a/benchmarks/bench_elliptic_parallel_template.nim
+++ b/benchmarks/bench_elliptic_parallel_template.nim
@@ -56,7 +56,7 @@ proc multiAddParallelBench*(EC: typedesc, numInputs: int, iters: int) =
 type BenchMsmContext*[EC] = object
   tp: Threadpool
   numInputs: int
-  coefs: seq[getBigInt(EC.F.Name, kScalarField)]
+  coefs: seq[getBigInt(EC.getName(), kScalarField)]
   points: seq[affine(EC)]
 
 proc createBenchMsmContext*(EC: typedesc, inputSizes: openArray[int]): BenchMsmContext[EC] =
@@ -66,6 +66,7 @@ proc createBenchMsmContext*(EC: typedesc, inputSizes: openArray[int]): BenchMsmC
   const bits = EC.getScalarField().bits()
   type ECaff = affine(EC)
 
+  result.numInputs = maxNumInputs
   result.points = newSeq[ECaff](maxNumInputs)
   result.coefs = newSeq[BigInt[bits]](maxNumInputs)
 
@@ -81,7 +82,7 @@ proc createBenchMsmContext*(EC: typedesc, inputSizes: openArray[int]): BenchMsmC
       var tmp = threadRng.random_unsafe(EC)
       tmp.clearCofactor()
       points[i].affine(tmp)
-      coefs[i] = rng.random_unsafe(BigInt[bits])
+      coefs[i] = threadRng.random_unsafe(BigInt[bits])
 
   let chunks = balancedChunksPrioNumber(0, maxNumInputs, result.tp.numThreads)
 
diff --git a/benchmarks/zkalc.nim b/benchmarks/zkalc.nim
new file mode 100644
index 000000000..fb452d66b
--- /dev/null
+++ b/benchmarks/zkalc.nim
@@ -0,0 +1,536 @@
+# Constantine
+# Copyright (c) 2018-2019    Status Research & Development GmbH
+# Copyright (c) 2020-Present Mamy André-Ratsimbazafy
+# Licensed and distributed under either of
+#   * MIT license (license terms in the root directory or at http://opensource.org/licenses/MIT).
+#   * Apache v2 license (license terms in the root directory or at http://www.apache.org/licenses/LICENSE-2.0).
+# at your option. This file may not be copied, modified, or distributed except according to those terms.
+
+# ############################################################
+#
+#             Benchmark results for zka.lc
+#
+# ############################################################
+
+# https://zka.lc/
+# https://github.com/mmaker/zkalc
+
+import
+  constantine/threadpool,
+  constantine/hashes,
+  constantine/lowlevel_fields,
+  # constantine/lowlevel_extension_fields,
+  constantine/math/extension_fields,
+  constantine/lowlevel_elliptic_curves,
+  constantine/lowlevel_elliptic_curves_parallel,
+  constantine/lowlevel_pairing_curves,
+
+  # Helpers
+  helpers/prng_unsafe,
+  # Standard library
+  std/[stats, monotimes, times, strformat, strutils, os, macros],
+  # Third-party
+  jsony, cliche
+
+# Workarounds generic sandwich bug in 1.6.x
+from constantine/named/algebras import matchingBigInt, matchingOrderBigInt, getLimbs2x
+
+type
+  ZkalcBenchDetails = object
+    `range`: seq[int]
+    results: seq[float64]
+    stddev: seq[float64]
+
+  ZkalcBenchResult = object
+    add_ff, mul_ff, invert: ZkalcBenchDetails
+    ip_ff: ZkalcBenchDetails # ip: inner-product
+    fft: ZkalcBenchDetails
+
+    add_G1, mul_G1, msm_G1: ZkalcBenchDetails
+    is_in_sub_G1: ZkalcBenchDetails
+    hash_G1: ZkalcBenchDetails
+
+    add_G2, mul_G2, msm_G2: ZkalcBenchDetails
+    is_in_sub_G2: ZkalcBenchDetails
+    hash_G2: ZkalcBenchDetails
+
+    mul_Gt: ZkalcBenchDetails
+    multiexp_Gt: ZkalcBenchDetails
+
+    pairing: ZkalcBenchDetails
+    multipairing: ZkalcBenchDetails
+
+type AggStats = tuple[rs: RunningStat, batchSize: int]
+
+# Utilities
+# -------------------------------------------------------------------------------------
+
+template bench(body: untyped): AggStats =
+  const warmupMs = 100
+  const batchMs = 10
+  const benchMs = 5000
+
+  block:
+    var stats: RunningStat
+    stats.clear()
+
+    proc warmup(warmupMs: int): tuple[num_iters: int, elapsedNs: int64] =
+      ## Warmup for the specified time and returns the number of iterations and time used
+      let start = getMonotime().ticks()
+      let stop = start + 1_000_000'i64*int64(warmupMs)
+
+      var num_iters = 0
+
+      while true:
+        body
+
+        let cur = getMonotime().ticks()
+        num_iters += 1
+
+        if cur >= stop:
+          return (num_iters, cur - start)
+
+    # Warmup and measure how many iterations are done during warmup
+    let (candidateIters, elapsedNs) = warmup(warmupMs)
+
+    # Deduce batch size for bench iterations so that each batch is atleast 10ms to amortize clock overhead
+    # See https://gms.tf/on-the-costs-of-syscalls.html on clock and syscall latencies and vDSO.
+    let batchSize = max(1, int(candidateIters.float64 * batchMs.float64 / warmupMs.float64))
+    # Compute the number of iterations for ~5s of benchmarks
+    let iters = int(
+      (candidateIters.float64 / batchSize.float64) *       # Divide the computed number of iterations by the size of the batch
+      max(1, benchMs.float64 / (elapsedNs.float64 * 1e-6)) # Scale by the ratio of bench time / warmup time
+    )
+
+    for _ in 0 ..< iters:
+      let start = getMonotime()
+
+      for _ in 0 ..< batchSize:
+        body
+
+      let stop = getMonotime()
+      let elapsedNs = (stop.ticks() - start.ticks()) div batchSize
+
+      # We can store integers up to 2⁵³ in a float64 without loss of precision (see also ulp)
+      # 1 billion is ~ 2³⁰, so you would need 2²³ seconds = 8388608s = 13 weeks 6 days 2 hours 10 minutes 8 seconds
+      stats.push(elapsedNs.int)
+
+    (stats, batchSize)
+
+proc report(op: string, curve: Algebra, aggStats: AggStats) =
+  let avg = aggStats.rs.mean()
+  let stddev = aggStats.rs.standardDeviationS() # Sample standard deviation (and not population)
+  let coefvar = stddev / avg * 100 # coefficient of variation
+  let throughput = 1e9 / float64(avg)
+  let iters = aggStats.rs.n
+  let batchSize = aggStats.batchSize
+  echo &"{op:<50} {$curve:<10} {throughput:>15.3f} ops/s {avg:>15.1f} ns/op (avg)    ±{coefvar:>4.1f}% (coef var)    {iters:>4} iterations of {batchSize:>6} operations"
+
+proc separator(length: int) =
+  echo "-".repeat(length)
+
+proc separator() = separator(174)
+
+proc toZkalc(stats: AggStats, size = 1): ZkalcBenchDetails =
+  ZkalcBenchDetails(
+    `range`: @[size],
+    results: @[stats.rs.mean()],
+    stddev: @[stats.rs.standardDeviationS()] # Sample standard deviation (and not population)
+  )
+
+proc append(details: var ZkalcBenchDetails, stats: AggStats, size: int) =
+  details.`range`.add size
+  details.results.add stats.rs.mean()
+  details.stddev.add  stats.rs.standardDeviationS() # Sample standard deviation (and not population)
+
+# Prevent compiler optimizing benchmark away
+# -------------------------------------------------------------------------------------
+# This doesn't always work unfortunately ...
+
+proc volatilize(x: ptr byte) {.codegenDecl: "$# $#(char const volatile *x)", inline.} =
+  discard
+
+template preventOptimAway*[T](x: var T) =
+  volatilize(cast[ptr byte](addr x))
+
+template preventOptimAway*[T](x: T) =
+  volatilize(cast[ptr byte](unsafeAddr x))
+
+# Field benches
+# -------------------------------------------------------------------------------------
+
+proc benchFrAdd(rng: var RngState, curve: static Algebra): ZkalcBenchDetails =
+  var x = rng.random_unsafe(Fr[curve])
+  let y = rng.random_unsafe(Fr[curve])
+
+  preventOptimAway(x)
+  preventOptimAway(y)
+
+  let stats = bench():
+    x += y
+
+  report("𝔽r Addition", curve, stats)
+  stats.toZkalc()
+
+proc benchFrMul(rng: var RngState, curve: static Algebra): ZkalcBenchDetails =
+  var x = rng.random_unsafe(Fr[curve])
+  let y = rng.random_unsafe(Fr[curve])
+
+  preventOptimAway(x)
+  preventOptimAway(y)
+
+  let stats = bench():
+    x *= y
+
+  report("𝔽r Multiplication", curve, stats)
+  stats.toZkalc()
+
+proc benchFrInv(rng: var RngState, curve: static Algebra, useVartime: bool): ZkalcBenchDetails =
+  var x = rng.random_unsafe(Fr[curve])
+
+  if useVartime:
+    let stats = bench():
+      x.inv_vartime()
+
+    report("𝔽r Inversion " & align("| vartime", 28), curve, stats)
+    stats.toZkalc()
+  else:
+    let stats = bench():
+      x.inv()
+
+    report("𝔽r Inversion " & align("| constant-time", 28), curve, stats)
+    stats.toZkalc()
+
+proc benchFrIP(rng: var RngState, curve: static Algebra): ZkalcBenchDetails =
+
+  var r: Fr[curve]
+  let a = rng.random_unsafe(Fr[curve])
+  let b = rng.random_unsafe(Fr[curve])
+  let u = rng.random_unsafe(Fr[curve])
+  let v = rng.random_unsafe(Fr[curve])
+
+  preventOptimAway(r)
+  preventOptimAway(a)
+  preventOptimAway(b)
+  preventOptimAway(u)
+  preventOptimAway(v)
+
+  let stats = bench():
+    r.sumprod([a, b], [u, v])
+
+  report("𝔽r Sum of products of size 2", curve, stats)
+  stats.toZkalc(2)
+
+# EC benches
+# -------------------------------------------------------------------------------------
+
+proc benchEcAdd(rng: var RngState, EC: typedesc, useVartime: bool): ZkalcBenchDetails =
+  const G =
+    when EC.G == G1: "𝔾₁"
+    else: "𝔾₂"
+  const curve = EC.getName()
+
+  var r {.noInit.}: EC
+  let P = rng.random_unsafe(EC)
+  let Q = rng.random_unsafe(EC)
+
+  preventOptimAway(r)
+  preventOptimAway(P)
+  preventOptimAway(Q)
+
+  if useVartime:
+    let stats = bench():
+      r.sum_vartime(P, Q)
+
+    report(G & " Addition " & align("| vartime", 29), curve, stats)
+    stats.toZkalc()
+  else:
+    let stats = bench():
+      r.sum(P, Q)
+
+    report(G & " Addition " & align("| constant-time", 29), curve, stats)
+    stats.toZkalc()
+
+proc benchEcMul(rng: var RngState, EC: typedesc, useVartime: bool): ZkalcBenchDetails =
+  const G =
+    when EC.G == G1: "𝔾₁"
+    else: "𝔾₂"
+  const curve = EC.getName()
+
+  var r {.noInit.}: EC
+  var P = rng.random_unsafe(EC)
+  P.clearCofactor()
+  let k = rng.random_unsafe(Fr[curve].getBigInt())
+
+  preventOptimAway(r)
+  preventOptimAway(P)
+
+  if useVartime:
+    let stats = bench():
+      r.scalarMul_vartime(k, P)
+
+    report(G & " Scalar Multiplication " & align("| vartime", 16), curve, stats)
+    stats.toZkalc()
+  else:
+    let stats = bench():
+      r.scalarMul(k, P)
+
+    report(G & " Scalar Multiplication " & align("| constant-time", 16), curve, stats)
+    stats.toZkalc()
+
+# EC Msm benches
+# -------------------------------------------------------------------------------------
+
+type BenchMsmContext[EC] = object
+  numInputs: int
+  coefs: seq[getBigInt(EC.getName(), kScalarField)]
+  points: seq[affine(EC)]
+
+proc createBenchMsmContext*(rng: var RngState, EC: typedesc, maxNumInputs: int): BenchMsmContext[EC] =
+  let tp = Threadpool.new()
+
+  type Big = typeof(result.coefs[0])
+  type ECaff = affine(EC)
+
+  result.numInputs = maxNumInputs
+  result.points = newSeq[ECaff](maxNumInputs)
+  result.coefs = newSeq[Big](maxNumInputs)
+
+  proc genCoefPointPairsChunk[EC, ECaff](rngSeed: uint64, start, len: int, points: ptr ECaff, coefs: ptr Big) {.nimcall.} =
+    let points = cast[ptr UncheckedArray[ECaff]](points)
+    let coefs = cast[ptr UncheckedArray[Big]](coefs)
+
+    # RNGs are not threadsafe, create a threadlocal one seeded from the global RNG
+    var threadRng: RngState
+    threadRng.seed(rngSeed)
+
+    for i in start ..< start + len:
+      var tmp = threadRng.random_unsafe(EC)
+      tmp.clearCofactor()
+      points[i].affine(tmp)
+      coefs[i] = threadRng.random_unsafe(Big)
+
+  let chunks = balancedChunksPrioNumber(0, maxNumInputs, tp.numThreads)
+
+  stdout.write &"Generating {maxNumInputs} (coefs, points) pairs ... "
+  stdout.flushFile()
+
+  let start = getMonotime()
+
+  syncScope:
+    for (id, start, size) in items(chunks):
+      tp.spawn genCoefPointPairsChunk[EC, ECaff](rng.next(), start, size, result.points[0].addr, result.coefs[0].addr)
+
+  # Even if child threads are sleeping, it seems like perf is lower when there are threads around
+  # maybe because the kernel has more overhead or time quantum to keep track off so shut them down.
+  tp.shutdown()
+
+  let stop = getMonotime()
+  stdout.write &"in {float64(inNanoSeconds(stop-start)) / 1e6:6.3f} ms\n"
+
+proc benchEcMsm[EC](ctx: BenchMsmContext[EC]): ZkalcBenchDetails =
+  const G =
+    when EC.G == G1: "𝔾₁"
+    else: "𝔾₂"
+  const curve = EC.getName()
+
+  let tp = Threadpool.new()
+  var size = 2
+  while size <= ctx.numInputs:
+    var r{.noInit.}: EC
+    template coefs: untyped = ctx.coefs.toOpenArray(0, size-1)
+    template points: untyped = ctx.points.toOpenArray(0, size-1)
+
+    let stats = bench():
+      tp.multiScalarMul_vartime_parallel(r, coefs, points)
+
+    report(G & " MSM " & align($size, 9) & ", " & align($tp.numThreads & " threads", 11) & align("| vartime", 12), curve, stats)
+    result.append(stats, size)
+
+    size *= 2
+
+  tp.shutdown()
+
+# EC serialization benches
+# -------------------------------------------------------------------------------------
+
+proc benchEcIsInSubgroup(rng: var RngState, EC: type): ZkalcBenchDetails =
+  const G =
+    when EC.G == G1: "𝔾₁"
+    else: "𝔾₂"
+  const curve = EC.getName()
+
+  var P = rng.random_unsafe(EC)
+  P.clearCofactor()
+  preventOptimAway(P)
+
+  let stats = bench():
+    discard P.isInSubgroup()
+
+  report(G & " Subgroup Check", curve, stats)
+  stats.toZkalc()
+
+proc benchEcHashToCurve(rng: var RngState, EC: type): ZkalcBenchDetails =
+  const G =
+    when EC.G == G1: "𝔾₁"
+    else: "𝔾₂"
+  const curve = EC.getName()
+
+  const dst = "Constantine_Zkalc_Bench_HashToCurve"
+  # Gnark uses a message of size 54, probably to not spill over padding with SHA256
+  let msg = "Privacy is necessary for an open society [...]"
+
+  var P {.noInit.}: EC
+
+  let stats = bench():
+    sha256.hashToCurve(
+      k = 128,
+      output = P,
+      augmentation = "",
+      message = msg,
+      domainSepTag = dst
+    )
+
+  report(G & " Hash-to-Curve", curve, stats)
+  stats.toZkalc()
+
+# Pairing benches
+# -------------------------------------------------------------------------------------
+
+func clearCofactor[F; G: static Subgroup](
+       ec: var EC_ShortW_Aff[F, G]) =
+  # For now we don't have any affine operation defined
+  var t {.noInit.}: EC_ShortW_Prj[F, G]
+  t.fromAffine(ec)
+  t.clearCofactor()
+  ec.affine(t)
+
+func random_point*(rng: var RngState, EC: typedesc): EC {.noInit.} =
+  result = rng.random_unsafe(EC)
+  result.clearCofactor()
+
+proc benchPairing*(rng: var RngState, curve: static Algebra): ZkalcBenchDetails =
+  let
+    P = rng.random_point(EC_ShortW_Aff[Fp[curve], G1])
+    Q = rng.random_point(EC_ShortW_Aff[Fp2[curve], G2])
+
+  var f: Fp12[curve]
+  let stats = bench():
+    f.pairing(P, Q)
+
+  report("Pairing", curve, stats)
+  stats.toZkalc()
+
+proc benchMultiPairing*(rng: var RngState, curve: static Algebra, maxNumInputs: int): ZkalcBenchDetails =
+  var
+    Ps = newSeq[EC_ShortW_Aff[Fp[curve], G1]](maxNumInputs)
+    Qs = newSeq[EC_ShortW_Aff[Fp2[curve], G2]](maxNumInputs)
+
+  stdout.write &"Generating {maxNumInputs} (𝔾₁, 𝔾₂) pairs ... "
+  stdout.flushFile()
+
+  let start = getMonotime()
+
+  for i in 0 ..< maxNumInputs:
+    Ps[i] = rng.random_point(typeof(Ps[0]))
+    Qs[i] = rng.random_point(typeof(Qs[0]))
+
+  let stop = getMonotime()
+  stdout.write &"in {float64(inNanoSeconds(stop-start)) / 1e6:6.3f} ms\n"
+  separator()
+
+  var size = 2
+  while size <= maxNumInputs:
+    var f{.noInit.}: Fp12[curve]
+    let stats = bench():
+      f.pairing(Ps.toOpenArray(0, size-1), Qs.toOpenArray(0, size-1))
+
+    report("Multipairing " & align($size, 5), curve, stats)
+    result.append(stats, size)
+
+    size *= 2
+
+# Run benches
+# -------------------------------------------------------------------------------------
+
+proc runBenches(curve: static Algebra, useVartime: bool): ZkalcBenchResult =
+  var rng: RngState
+  rng.seed(42)
+
+  var zkalc: ZkalcBenchResult
+
+  # Fields
+  # --------------------------------------------------------------------
+  separator()
+  zkalc.add_ff = rng.benchFrAdd(curve)
+  zkalc.mul_ff = rng.benchFrMul(curve)
+  zkalc.invert = rng.benchFrInv(curve, useVartime)
+  zkalc.ip_ff  = rng.benchFrIP(curve)
+  separator()
+
+  # Elliptic curve
+  # --------------------------------------------------------------------
+  type EcG1 = EC_ShortW_Jac[Fp[curve], G1]
+
+  zkalc.add_g1 = rng.benchEcAdd(EcG1, useVartime)
+  zkalc.mul_g1 = rng.benchEcMul(EcG1, useVartime)
+  separator()
+  let ctxG1    = rng.createBenchMsmContext(EcG1, maxNumInputs = 2097152)
+  separator()
+  zkalc.msm_g1 = benchEcMsm(ctxG1)
+  separator()
+  zkalc.is_in_sub_G1 = rng.benchEcIsInSubgroup(EcG1)
+  when curve in {BN254_Snarks, BLS12_381}:
+    zkalc.hash_G1 = rng.benchEcHashToCurve(EcG1)
+  separator()
+
+  # Pairing-friendly curve only
+  # --------------------------------------------------------------------
+
+  when curve.isPairingFriendly():
+
+    # Elliptic curve 𝔾2
+    # --------------------------------------------------------------------
+
+    type EcG2 = EC_ShortW_Jac[Fp2[curve], G2] # For now we only supports G2 on Fp2 (not Fp like BW6 or Fp4 like BLS24)
+
+    zkalc.add_g2 = rng.benchEcAdd(EcG2, useVartime)
+    zkalc.mul_g2 = rng.benchEcMul(EcG2, useVartime)
+    separator()
+    let ctxG2    = rng.createBenchMsmContext(EcG2, maxNumInputs = 2097152)
+    separator()
+    zkalc.msm_g2 = benchEcMsm(ctxG2)
+    separator()
+    zkalc.is_in_sub_G2 = rng.benchEcIsInSubgroup(EcG2)
+    when curve in {BN254_Snarks, BLS12_381}:
+      zkalc.hash_G2 = rng.benchEcHashToCurve(EcG2)
+    separator()
+
+    # Pairings
+    # --------------------------------------------------------------------
+
+    zkalc.pairing = rng.benchPairing(curve)
+    separator()
+    zkalc.multipairing = rng.benchMultiPairing(curve, maxNumInputs = 1024)
+    separator()
+
+  return zkalc
+
+proc main() =
+  let cmd = commandLineParams()
+  cmd.getOpt (curve: BN254_Snarks, vartime: true, o: "constantine-bench-zkalc-" & $curve & "-" & now().format("yyyy-MM-dd--HH-mm-ss") & ".bench.json")
+
+  let results =
+    case curve
+    of BN254_Snarks: BN254_Snarks.runBenches(vartime)
+    of Pallas:       Pallas      .runBenches(vartime)
+    of Vesta:        Vesta       .runBenches(vartime)
+    of BLS12_377:    BLS12_377   .runBenches(vartime)
+    of BLS12_381:    BLS12_381   .runBenches(vartime)
+    else:
+      raise newException(ValueError, "This curve '" & $curve & "' is not configured for benchmarking at the moment.")
+
+  writeFile(o, results.toJSON())
+
+when isMainModule:
+  main()
diff --git a/bin/.gitignore b/bin/.gitignore
new file mode 100644
index 000000000..03fff8904
--- /dev/null
+++ b/bin/.gitignore
@@ -0,0 +1,12 @@
+# Ignore all
+*
+
+# Unignore all with extensions
+!*.*
+
+# Unignore all dirs
+!*/
+
+# Reignore executables
+*.exe
+*.out
diff --git a/constantine.nimble b/constantine.nimble
index dd05e12a2..dbdccb3f8 100644
--- a/constantine.nimble
+++ b/constantine.nimble
@@ -9,6 +9,24 @@ license       = "MIT or Apache License 2.0"
 
 requires "nim >= 1.6.12"
 
+when (NimMajor, NimMinor) >= (2, 0): # Task-level dependencies
+  taskRequires "make_zkalc", "jsony"
+  taskRequires "make_zkalc", "cliche"
+
+  taskRequires "test", "jsony"
+  taskRequires "test", "yaml"
+  taskRequires "test", "gmp#head"
+
+  taskRequires "test_parallel", "jsony"
+  taskRequires "test_parallel", "yaml"
+  taskRequires "test_parallel", "gmp#head"
+
+  taskRequires "test_no_gmp", "jsony"
+  taskRequires "test_no_gmp", "yaml"
+
+  taskRequires "test_parallel_no_gmp", "jsony"
+  taskRequires "test_parallel_no_gmp", "yaml"
+
 # Nimscript imports
 # ----------------------------------------------------------------
 
@@ -277,6 +295,14 @@ task make_lib_rust, "Build Constantine library (use within a Rust build.rs scrip
                  else: "--passC:-fPIC"
   genStaticLib(rustOutDir, rustOutDir/"nimcache", extflags)
 
+task make_zkalc, "Build a benchmark executable for zkalc (with Clang)":
+  exec "nim c --cc:clang " &
+       releaseBuildOptions(bmBinary) &
+       " --threads:on " &
+       " --out:bin/constantine-bench-zkalc " &
+       " --nimcache:nimcache/bench_zkalc " &
+       " benchmarks/zkalc.nim"
+
 proc testLib(path, testName: string, useGMP: bool) =
   let dynlibName = if defined(windows): "constantine.dll"
                    elif defined(macosx) or defined(macos): "libconstantine.dylib"
@@ -640,7 +666,8 @@ const benchDesc = [
   "bench_eth_eip2537_subgroup_checks_impact",
   "bench_verkle_primitives",
   "bench_eth_evm_precompiles",
-  "bench_multilinear_extensions"
+  "bench_multilinear_extensions",
+  # "zkalc", # Already tested through make_zkalc
 ]
 
 # For temporary (hopefully) investigation that can only be reproduced in CI
diff --git a/constantine/lowlevel_bigints.nim b/constantine/lowlevel_bigints.nim
index fd67dfb88..50e58eaaf 100644
--- a/constantine/lowlevel_bigints.nim
+++ b/constantine/lowlevel_bigints.nim
@@ -27,7 +27,10 @@ import
 # Base types
 # ------------------------------------------------------------
 
-export abstractions
+export
+  abstractions.SecretBool,
+  abstractions.SecretWord,
+  abstractions.BigInt
 
 # BigInt
 # ------------------------------------------------------------
diff --git a/constantine/lowlevel_elliptic_curves.nim b/constantine/lowlevel_elliptic_curves.nim
index d1430498f..fb95bf64a 100644
--- a/constantine/lowlevel_elliptic_curves.nim
+++ b/constantine/lowlevel_elliptic_curves.nim
@@ -14,7 +14,6 @@ import
     ./math/elliptic/[
       ec_scalar_mul_vartime,
       ec_multi_scalar_mul],
-    ./math/io/io_ec,
     ./hash_to_curve/hash_to_curve
 
 # ############################################################
@@ -35,8 +34,16 @@ import
 # ------------------------------------------------------------
 
 export
-  abstractions,
-  algebras.Algebra
+  abstractions.SecretBool,
+  abstractions.SecretWord,
+  abstractions.BigInt,
+  algebras.Algebra,
+  algebras.getBigInt,
+  algebras.FieldKind,
+  algebras.isPairingFriendly
+
+# Generic sandwich
+export abstractions
 
 # Elliptic curve
 # ------------------------------------------------------------
@@ -46,7 +53,8 @@ export
   ec_shortweierstrass.EC_ShortW_Aff,
   ec_shortweierstrass.EC_ShortW_Jac,
   ec_shortweierstrass.EC_ShortW_Prj,
-  ec_shortweierstrass.EC_ShortW
+  ec_shortweierstrass.EC_ShortW,
+  ec_shortweierstrass.getName
 
 export ec_shortweierstrass.`==`
 export ec_shortweierstrass.isNeutral
diff --git a/constantine/lowlevel_extension_fields.nim b/constantine/lowlevel_extension_fields.nim
index 23edabe5e..554550851 100644
--- a/constantine/lowlevel_extension_fields.nim
+++ b/constantine/lowlevel_extension_fields.nim
@@ -31,18 +31,25 @@ import
 # ------------------------------------------------------------
 
 export
-  abstractions,
-  algebras.Algebra
+  abstractions.SecretBool,
+  abstractions.SecretWord,
+  algebras.Algebra,
+  algebras.getBigInt
+
+export
+  algebras.Fp,
+  algebras.Fr,
+  algebras.FF
 
 # Extension fields
 # ------------------------------------------------------------
 
 export
-  extension_fields.Fp2
+  extension_fields.Fp2,
   # TODO: deal with Fp2->Fp6 vs Fp3->Fp6 and Fp2->Fp6->Fp12 vs Fp2->Fp4->Fp12
   # extension_fields.Fp4,
   # extension_fields.Fp6,
-  # extension_fields.Fp12
+  extension_fields.Fp12
 
 # Generic sandwich - https://github.com/nim-lang/Nim/issues/11225
 export extension_fields.c0, extension_fields.`c0=`
@@ -86,4 +93,3 @@ export extension_fields.sqrt_if_square
 export extension_fields.sqrt
 
 export frobenius.frobenius_map
-
diff --git a/constantine/lowlevel_fields.nim b/constantine/lowlevel_fields.nim
index c5f87fc85..14dc1c6a0 100644
--- a/constantine/lowlevel_fields.nim
+++ b/constantine/lowlevel_fields.nim
@@ -30,8 +30,11 @@ import
 # ------------------------------------------------------------
 
 export
-  abstractions,
-  algebras.Algebra
+  abstractions.SecretBool,
+  abstractions.SecretWord,
+  abstractions.BigInt,
+  algebras.Algebra,
+  algebras.getBigInt
 
 # Scalar field Fr and Prime Field Fp
 # ------------------------------------------------------------
@@ -88,6 +91,7 @@ export arithmetic.prod
 export arithmetic.`*=`
 export arithmetic.square
 export arithmetic.square_repeated
+export arithmetic.sumprod
 
 export arithmetic.csetZero
 export arithmetic.csetOne
@@ -97,6 +101,7 @@ export arithmetic.csub
 
 export arithmetic.div2
 export arithmetic.inv
+export arithmetic.inv_vartime
 
 export arithmetic.isSquare
 export arithmetic.invsqrt
diff --git a/constantine/lowlevel_pairing_curves.nim b/constantine/lowlevel_pairing_curves.nim
index d4ef869c0..abbcc2aae 100644
--- a/constantine/lowlevel_pairing_curves.nim
+++ b/constantine/lowlevel_pairing_curves.nim
@@ -33,8 +33,10 @@ import
 # ------------------------------------------------------------
 
 export
-  abstractions,
-  algebras.Algebra
+  abstractions.SecretBool,
+  abstractions.SecretWord,
+  algebras.Algebra,
+  algebras.getBigInt
 
 # Pairings
 # ------------------------------------------------------------
diff --git a/constantine/math/elliptic/ec_multi_scalar_mul.nim b/constantine/math/elliptic/ec_multi_scalar_mul.nim
index 4d2e3d8fd..a95416010 100644
--- a/constantine/math/elliptic/ec_multi_scalar_mul.nim
+++ b/constantine/math/elliptic/ec_multi_scalar_mul.nim
@@ -429,13 +429,13 @@ template withEndo[coefsBits: static int, EC, ECaff](
            coefs: ptr UncheckedArray[BigInt[coefsBits]],
            points: ptr UncheckedArray[ECaff],
            N: int, c: static int) =
-  when hasEndomorphismAcceleration(EC.F.Name) and
+  when hasEndomorphismAcceleration(EC.getName()) and
         EndomorphismThreshold <= coefsBits and
         coefsBits <= EC.getScalarField().bits() and
         # computeEndomorphism assumes they can be applied to affine repr
         # but this is not the case for Bandersnatch/wagon
         # instead Twisted Edwards MSM should be overloaded for Projective/ProjectiveExtended
-        EC.F.Name notin {Bandersnatch, Banderwagon}:
+        EC.getName() notin {Bandersnatch, Banderwagon}:
     let (endoCoefs, endoPoints, endoN) = applyEndomorphism(coefs, points, N)
     # Given that bits and N changed, we are able to use a bigger `c`
     # but it has no significant impact on performance
diff --git a/constantine/math/elliptic/ec_multi_scalar_mul_parallel.nim b/constantine/math/elliptic/ec_multi_scalar_mul_parallel.nim
index 3c9824cb5..a132caff5 100644
--- a/constantine/math/elliptic/ec_multi_scalar_mul_parallel.nim
+++ b/constantine/math/elliptic/ec_multi_scalar_mul_parallel.nim
@@ -495,13 +495,13 @@ template withEndo[coefsBits: static int, EC, ECaff](
            coefs: ptr UncheckedArray[BigInt[coefsBits]],
            points: ptr UncheckedArray[ECaff],
            N: int, c: static int) =
-  when hasEndomorphismAcceleration(EC.F.Name) and
+  when hasEndomorphismAcceleration(EC.getName()) and
         EndomorphismThreshold <= coefsBits and
         coefsBits <= EC.getScalarField().bits() and
         # computeEndomorphism assumes they can be applied to affine repr
         # but this is not the case for Bandersnatch/wagon
         # instead Twisted Edwards MSM should be overloaded for Projective/ProjectiveExtended
-        EC.F.Name notin {Bandersnatch, Banderwagon}:
+        EC.getName() notin {Bandersnatch, Banderwagon}:
     let (endoCoefs, endoPoints, endoN) = applyEndomorphism_parallel(tp, coefs, points, N)
     # Given that bits and N changed, we are able to use a bigger `c`
     # but it has no significant impact on performance
@@ -518,7 +518,7 @@ template withEndo[coefsBits: static int, EC, ECaff](
            coefs: ptr UncheckedArray[BigInt[coefsBits]],
            points: ptr UncheckedArray[ECaff],
            N: int, c: static int, useParallelBuckets: static bool) =
-  when coefsBits <= EC.getScalarField().bits() and hasEndomorphismAcceleration(EC.F.Name):
+  when coefsBits <= EC.getScalarField().bits() and hasEndomorphismAcceleration(EC.getName()):
     let (endoCoefs, endoPoints, endoN) = applyEndomorphism_parallel(tp, coefs, points, N)
     # Given that bits and N changed, we are able to use a bigger `c`
     # but it has no significant impact on performance
diff --git a/constantine/math/elliptic/ec_multi_scalar_mul_scheduler.nim b/constantine/math/elliptic/ec_multi_scalar_mul_scheduler.nim
index ea3eb083f..f18af92ac 100644
--- a/constantine/math/elliptic/ec_multi_scalar_mul_scheduler.nim
+++ b/constantine/math/elliptic/ec_multi_scalar_mul_scheduler.nim
@@ -8,7 +8,7 @@
 
 import
   constantine/platforms/abstractions,
-  constantine/math/arithmetic,
+  constantine/math/[arithmetic, extension_fields],
   ./ec_shortweierstrass_affine,
   ./ec_shortweierstrass_jacobian,
   ./ec_shortweierstrass_projective,
@@ -16,7 +16,7 @@ import
   ./ec_twistededwards_projective,
   ./ec_twistededwards_affine
 
-export abstractions, arithmetic,
+export abstractions, arithmetic, extension_fields,
        ec_shortweierstrass_affine,
        ec_shortweierstrass_jacobian,
        ec_shortweierstrass_projective,
diff --git a/constantine/math/elliptic/ec_scalar_mul.nim b/constantine/math/elliptic/ec_scalar_mul.nim
index c837864ca..eaa57c6e7 100644
--- a/constantine/math/elliptic/ec_scalar_mul.nim
+++ b/constantine/math/elliptic/ec_scalar_mul.nim
@@ -242,7 +242,7 @@ func scalarMul*[EC](P: var EC, scalar: BigInt) {.inline, meter.} =
   ## - Cofactor to be cleared
   ## - 0 <= scalar < curve order
   ## Those will be assumed to maintain constant-time property
-  when EC.F.Name.hasEndomorphismAcceleration() and
+  when EC.getName().hasEndomorphismAcceleration() and
        BigInt.bits >= EndomorphismThreshold:
     when EC.F is Fp:
       P.scalarMulGLV_m2w2(scalar)
diff --git a/constantine/math/elliptic/ec_scalar_mul_vartime.nim b/constantine/math/elliptic/ec_scalar_mul_vartime.nim
index 438b403de..825a4394b 100644
--- a/constantine/math/elliptic/ec_scalar_mul_vartime.nim
+++ b/constantine/math/elliptic/ec_scalar_mul_vartime.nim
@@ -21,7 +21,8 @@ import
   constantine/math/io/io_bigints,
   constantine/platforms/abstractions,
   constantine/math_arbitrary_precision/arithmetic/limbs_views,
-  constantine/named/zoo_endomorphisms
+  constantine/named/zoo_endomorphisms,
+  constantine/named/algebras
 
 {.push raises: [].} # No exceptions allowed in core cryptographic operations
 {.push checks: off.} # No defects due to array bound checking or signed integer overflow allowed
@@ -336,11 +337,9 @@ func scalarMul_vartime*[scalBits; EC](P: var EC, scalar: BigInt[scalBits]) {.met
   else:
     {.error: "Unconfigured".}
 
-  const L = scalBits.ceilDiv_vartime(M) + 1
-
   let usedBits = scalar.limbs.getBits_LE_vartime()
 
-  when EC.F.Name.hasEndomorphismAcceleration():
+  when EC.getName().hasEndomorphismAcceleration():
     when scalBits >= EndomorphismThreshold: # Skip static: doAssert when multiplying by intentionally small scalars.
       if usedBits >= EndomorphismThreshold:
         when EC.F is Fp:
@@ -352,8 +351,8 @@ func scalarMul_vartime*[scalBits; EC](P: var EC, scalar: BigInt[scalBits]) {.met
         return
 
   if 64 < usedBits:
-    # With a window of 5, we precompute 2^3 = 8 points
-    P.scalarMul_minHammingWeight_windowed_vartime(scalar, window = 5)
+    # With a window of 4, we precompute 2^4 = 4 points
+    P.scalarMul_minHammingWeight_windowed_vartime(scalar, window = 4)
   elif 16 < usedBits:
     # With a window of 3, we precompute 2^1 = 2 points
     P.scalarMul_minHammingWeight_windowed_vartime(scalar, window = 3)
diff --git a/constantine/math/elliptic/ec_shortweierstrass_affine.nim b/constantine/math/elliptic/ec_shortweierstrass_affine.nim
index 53b9c9d0a..fc4ccde4f 100644
--- a/constantine/math/elliptic/ec_shortweierstrass_affine.nim
+++ b/constantine/math/elliptic/ec_shortweierstrass_affine.nim
@@ -38,6 +38,9 @@ type
 
   SexticNonResidue* = NonResidue
 
+template getName*(EC: type EC_ShortW_Aff): untyped =
+  EC.F.Name
+
 template getScalarField*(EC: type EC_ShortW_Aff): untyped =
   Fr[EC.F.Name]
 
diff --git a/constantine/math/elliptic/ec_shortweierstrass_batch_ops.nim b/constantine/math/elliptic/ec_shortweierstrass_batch_ops.nim
index 8ca885c19..9475e965e 100644
--- a/constantine/math/elliptic/ec_shortweierstrass_batch_ops.nim
+++ b/constantine/math/elliptic/ec_shortweierstrass_batch_ops.nim
@@ -7,6 +7,7 @@
 # at your option. This file may not be copied, modified, or distributed except according to those terms.
 
 import
+  constantine/named/algebras,
   constantine/platforms/abstractions,
   constantine/math/arithmetic,
   constantine/math/extension_fields,
diff --git a/constantine/math/elliptic/ec_shortweierstrass_jacobian.nim b/constantine/math/elliptic/ec_shortweierstrass_jacobian.nim
index 19936e0c2..682b5711a 100644
--- a/constantine/math/elliptic/ec_shortweierstrass_jacobian.nim
+++ b/constantine/math/elliptic/ec_shortweierstrass_jacobian.nim
@@ -37,6 +37,9 @@ type EC_ShortW_Jac*[F; G: static Subgroup] = object
   ## Note that jacobian coordinates are not unique
   x*, y*, z*: F
 
+template getName*(EC: type EC_ShortW_Jac): untyped =
+  EC.F.Name
+
 template getScalarField*(EC: type EC_ShortW_Jac): untyped =
   Fr[EC.F.Name]
 
@@ -213,6 +216,8 @@ template sumImpl[F; G: static Subgroup](
   # | Y₃ = R*(V-X₃)-S₁*HHH          | Y₃ = M*(S-X₃)-YY*YY                      |                   |                |
   # | Z₃ = Z₁*Z₂*H                  | Z₃ = Y₁*Z₁                               |                   |                |
 
+  bind mulCheckSparse
+
   # "when" static evaluation doesn't shortcut booleans :/
   # which causes issues when CoefA isn't an int but Fp or Fp2
   when CoefA is int:
@@ -1028,7 +1033,7 @@ func `~-`*(a: EC_ShortW_Jac, b: EC_ShortW_Aff): EC_ShortW_Jac {.noInit, inline.}
   ## This MUST NOT be used with secret data.
   ##
   ## This is highly VULNERABLE to timing attacks and power analysis attacks.]
-  ## 
+  ##
   ## Out-of-place functions SHOULD NOT be used in performance-critical subroutines as compilers
   ## tend to generate useless memory moves or have difficulties to minimize stack allocation
   ## and our types might be large (Fp12 ...)
diff --git a/constantine/math/elliptic/ec_shortweierstrass_jacobian_extended.nim b/constantine/math/elliptic/ec_shortweierstrass_jacobian_extended.nim
index f73584b13..8c37b4929 100644
--- a/constantine/math/elliptic/ec_shortweierstrass_jacobian_extended.nim
+++ b/constantine/math/elliptic/ec_shortweierstrass_jacobian_extended.nim
@@ -39,6 +39,9 @@ type EC_ShortW_JacExt*[F; G: static Subgroup] = object
   ## Note that extended jacobian coordinates are not unique
   x*, y*, zz*, zzz*: F
 
+template getName*(EC: type EC_ShortW_JacExt): untyped =
+  EC.F.Name
+
 func fromAffine*[F; G](jacext: var EC_ShortW_JacExt[F, G], aff: EC_ShortW_Aff[F, G]) {.inline.}
 
 func isNeutral*(P: EC_ShortW_JacExt): SecretBool {.inline, meter.} =
diff --git a/constantine/math/elliptic/ec_shortweierstrass_projective.nim b/constantine/math/elliptic/ec_shortweierstrass_projective.nim
index d4c68f1f6..912959f6c 100644
--- a/constantine/math/elliptic/ec_shortweierstrass_projective.nim
+++ b/constantine/math/elliptic/ec_shortweierstrass_projective.nim
@@ -37,6 +37,9 @@ type EC_ShortW_Prj*[F; G: static Subgroup] = object
   ## Note that projective coordinates are not unique
   x*, y*, z*: F
 
+template getName*(EC: type EC_ShortW_Prj): untyped =
+  EC.F.Name
+
 template getScalarField*(EC: type EC_ShortW_Prj): untyped =
   Fr[EC.F.Name]
 
diff --git a/constantine/math/elliptic/ec_twistededwards_affine.nim b/constantine/math/elliptic/ec_twistededwards_affine.nim
index 309641602..eac828867 100644
--- a/constantine/math/elliptic/ec_twistededwards_affine.nim
+++ b/constantine/math/elliptic/ec_twistededwards_affine.nim
@@ -28,6 +28,9 @@ type EC_TwEdw_Aff*[F] = object
   ## over a field F
   x*, y*: F
 
+template getName*(EC: type EC_TwEdw_Aff): untyped =
+  EC.F.Name
+
 template getScalarField*(EC: type EC_TwEdw_Aff): untyped =
   Fr[EC.F.Name]
 
diff --git a/constantine/math/elliptic/ec_twistededwards_projective.nim b/constantine/math/elliptic/ec_twistededwards_projective.nim
index 7fedf1d49..d43393859 100644
--- a/constantine/math/elliptic/ec_twistededwards_projective.nim
+++ b/constantine/math/elliptic/ec_twistededwards_projective.nim
@@ -33,6 +33,9 @@ type EC_TwEdw_Prj*[F] = object
   ## hence (aX² + Y²)Z² = Z⁴ + dX²Y²
   x*, y*, z*: F
 
+template getName*(EC: type EC_TwEdw_Prj): untyped =
+  EC.F.Name
+
 template getScalarField*(EC: type EC_TwEdw_Prj): untyped =
   Fr[EC.F.Name]
 
diff --git a/constantine/math/pairings/pairings_bls12.nim b/constantine/math/pairings/pairings_bls12.nim
index f831e7f1a..d10264330 100644
--- a/constantine/math/pairings/pairings_bls12.nim
+++ b/constantine/math/pairings/pairings_bls12.nim
@@ -163,15 +163,28 @@ func pairing_bls12*[Name](
   gt.finalExpEasy()
   gt.finalExpHard_BLS12()
 
-func pairing_bls12*[N: static int, Name](
+func pairing_bls12*[Name: static Algebra](
        gt: var Fp12[Name],
-       Ps: array[N, EC_ShortW_Aff[Fp[Name], G1]],
-       Qs: array[N, EC_ShortW_Aff[Fp2[Name], G2]]) {.meter.} =
+       Ps: ptr UncheckedArray[EC_ShortW_Aff[Fp[Name], G1]],
+       Qs: ptr UncheckedArray[EC_ShortW_Aff[Fp2[Name], G2]],
+       len: int) {.meter.} =
   ## Compute the optimal Ate Pairing for BLS12 curves
   ## Input: an array of Ps ∈ G1 and Qs ∈ G2
   ## Output:
   ##   The product of pairings
   ##   e(P₀, Q₀) * e(P₁, Q₁) * e(P₂, Q₂) * ... * e(Pₙ, Qₙ) ∈ Gt
-  gt.millerLoopAddchain(Qs.asUnchecked(), Ps.asUnchecked(), N)
+  gt.millerLoopAddchain(Qs, Ps, len)
   gt.finalExpEasy()
   gt.finalExpHard_BLS12()
+
+func pairing_bls12*[Name: static Algebra](
+       gt: var Fp12[Name],
+       Ps: openArray[EC_ShortW_Aff[Fp[Name], G1]],
+       Qs: openArray[EC_ShortW_Aff[Fp2[Name], G2]]) {.inline.} =
+  ## Compute the optimal Ate Pairing for BLS12 curves
+  ## Input: an array of Ps ∈ G1 and Qs ∈ G2
+  ## Output:
+  ##   The product of pairings
+  ##   e(P₀, Q₀) * e(P₁, Q₁) * e(P₂, Q₂) * ... * e(Pₙ, Qₙ) ∈ Gt
+  debug: doAssert Ps.len == Qs.len
+  gt.pairing_bls12(Ps.asUnchecked(), Qs.asUnchecked(), Ps.len)
diff --git a/constantine/math/pairings/pairings_bn.nim b/constantine/math/pairings/pairings_bn.nim
index a732b4756..2c7f479eb 100644
--- a/constantine/math/pairings/pairings_bn.nim
+++ b/constantine/math/pairings/pairings_bn.nim
@@ -175,18 +175,31 @@ func pairing_bn*[Name](
   gt.finalExpEasy()
   gt.finalExpHard_BN()
 
-func pairing_bn*[N: static int, Name](
+func pairing_bn*[Name: static Algebra](
        gt: var Fp12[Name],
-       Ps: array[N, EC_ShortW_Aff[Fp[Name], G1]],
-       Qs: array[N, EC_ShortW_Aff[Fp2[Name], G2]]) {.meter.} =
+       Ps: ptr UncheckedArray[EC_ShortW_Aff[Fp[Name], G1]],
+       Qs: ptr UncheckedArray[EC_ShortW_Aff[Fp2[Name], G2]],
+       len: int) {.meter.} =
   ## Compute the optimal Ate Pairing for BLS12 curves
   ## Input: an array of Ps ∈ G1 and Qs ∈ G2
   ## Output:
   ##   The product of pairings
   ##   e(P₀, Q₀) * e(P₁, Q₁) * e(P₂, Q₂) * ... * e(Pₙ, Qₙ) ∈ Gt
   when Name == BN254_Nogami:
-    gt.millerLoopAddChain(Qs.asUnchecked(), Ps.asUnchecked(), N)
+    gt.millerLoopAddChain(Qs, Ps, len)
   else:
-    gt.millerLoopGenericBN(Qs.asUnchecked(), Ps.asUnchecked(), N)
+    gt.millerLoopGenericBN(Qs, Ps, len)
   gt.finalExpEasy()
   gt.finalExpHard_BN()
+
+func pairing_bn*[Name: static Algebra](
+       gt: var Fp12[Name],
+       Ps: openArray[EC_ShortW_Aff[Fp[Name], G1]],
+       Qs: openArray[EC_ShortW_Aff[Fp2[Name], G2]]) {.inline.} =
+  ## Compute the optimal Ate Pairing for BLS12 curves
+  ## Input: an array of Ps ∈ G1 and Qs ∈ G2
+  ## Output:
+  ##   The product of pairings
+  ##   e(P₀, Q₀) * e(P₁, Q₁) * e(P₂, Q₂) * ... * e(Pₙ, Qₙ) ∈ Gt
+  debug: doAssert Ps.len == Qs.len
+  gt.pairing_bn(Ps.asUnchecked(), Qs.asUnchecked(), Ps.len)
diff --git a/constantine/math/pairings/pairings_generic.nim b/constantine/math/pairings/pairings_generic.nim
index 893891c18..a8a79f696 100644
--- a/constantine/math/pairings/pairings_generic.nim
+++ b/constantine/math/pairings/pairings_generic.nim
@@ -13,7 +13,7 @@ import
   constantine/math/extension_fields,
   constantine/named/zoo_pairings
 
-func pairing*[Name](gt: var Fp12[Name], P, Q: auto) {.inline.} =
+func pairing*[Name: static Algebra](gt: var Fp12[Name], P, Q: auto) {.inline.} =
   when family(Name) == BarretoNaehrig:
     pairing_bn(gt, P, Q)
   elif family(Name) == BarretoLynnScott:
@@ -21,13 +21,13 @@ func pairing*[Name](gt: var Fp12[Name], P, Q: auto) {.inline.} =
   else:
     {.error: "Pairing not implemented for " & $Name.}
 
-func millerLoop*[Name](gt: var Fp12[Name], Q, P: auto, n: int) {.inline.} =
+func millerLoop*[Name: static Algebra](gt: var Fp12[Name], Q, P: auto, n: int) {.inline.} =
   when Name == BN254_Snarks:
     gt.millerLoopGenericBN(Q, P, n)
   else:
     gt.millerLoopAddchain(Q, P, n)
 
-func finalExp*[Name](gt: var Fp12[Name]){.inline.} =
+func finalExp*[Name: static Algebra](gt: var Fp12[Name]){.inline.} =
   gt.finalExpEasy()
   when family(Name) == BarretoNaehrig:
     gt.finalExpHard_BN()
diff --git a/constantine/math/polynomials/fft.nim b/constantine/math/polynomials/fft.nim
index 8af0af95a..6efc8e118 100644
--- a/constantine/math/polynomials/fft.nim
+++ b/constantine/math/polynomials/fft.nim
@@ -32,12 +32,12 @@ type
   ECFFT_Descriptor*[EC] = object
     ## Metadata for FFT on Elliptic Curve
     order*: int
-    rootsOfUnity*: ptr UncheckedArray[getBigInt(EC.F.Name, kScalarField)]
+    rootsOfUnity*: ptr UncheckedArray[getBigInt(EC.getName(), kScalarField)]
       ## domain, starting and ending with 1, length is cardinality+1
       ## This allows FFT and inverse FFT to use the same buffer for roots.
 
 func computeRootsOfUnity[EC](ctx: var ECFFT_Descriptor[EC], generatorRootOfUnity: auto) =
-  static: doAssert typeof(generatorRootOfUnity) is Fr[EC.F.Name]
+  static: doAssert typeof(generatorRootOfUnity) is Fr[EC.getName()]
 
   ctx.rootsOfUnity[0].setOne()
 
diff --git a/constantine/named/properties_curves.nim b/constantine/named/properties_curves.nim
index 0ad8621c3..30c303f69 100644
--- a/constantine/named/properties_curves.nim
+++ b/constantine/named/properties_curves.nim
@@ -10,7 +10,8 @@ import
   std/macros,
   # Internal
   ./config_fields_and_curves,
-  ./deriv/parser_curves
+  ./deriv/parser_curves,
+  ./properties_fields
 
 export Algebra, CurveFamily, SexticTwist
 
@@ -46,6 +47,9 @@ template getField*(Name: static Algebra, kind: static FieldKind): untyped =
 template family*(Name: Algebra): CurveFamily =
   CurveFamilies[Name]
 
+template isPairingFriendly*(Name: Algebra): bool =
+  family(Name) in {BarretoNaehrig, BarretoLynnScott, BrezingWeng}
+
 macro getEquationForm*(Name: static Algebra): untyped =
   ## Returns the equation form
   ## (ShortWeierstrass, Montgomery, Twisted Edwards, Weierstrass, ...)
diff --git a/constantine/named/zoo_subgroups.nim b/constantine/named/zoo_subgroups.nim
index 1be5a22d5..2942ade07 100644
--- a/constantine/named/zoo_subgroups.nim
+++ b/constantine/named/zoo_subgroups.nim
@@ -26,12 +26,14 @@ export
   bn254_nogami_subgroups,
   bn254_snarks_subgroups,
   bw6_761_subgroups,
+  pallas_subgroups,
+  vesta_subgroups,
   secp256k1_subgroups
 
-func clearCofactor*[ECP](P: var ECP) {.inline.} =
+func clearCofactor*[EC](P: var EC) {.inline.} =
   ## Clear the cofactor of a point on the curve
   ## From a point on the curve, returns a point on the subgroup of order r
-  when ECP.F.Name in {BN254_Nogami, BN254_Snarks, BLS12_377, BLS12_381}:
+  when EC.F.Name in {BN254_Nogami, BN254_Snarks, BLS12_377, BLS12_381}:
     P.clearCofactorFast()
   else:
     P.clearCofactorReference()
diff --git a/research/kzg/fft_g1.nim b/research/kzg/fft_g1.nim
index 19f20fad0..2b807450b 100644
--- a/research/kzg/fft_g1.nim
+++ b/research/kzg/fft_g1.nim
@@ -176,7 +176,7 @@ func ifft_vartime*[EC](
   var voutput = output.toView()
   fft_internal(voutput, vals.toView(), rootz)
 
-  var invLen {.noInit.}: Fr[EC.F.Name]
+  var invLen {.noInit.}: Fr[EC.getName()]
   invLen.fromUint(vals.len.uint64)
   invLen.inv_vartime()
   let inv = invLen.toBig()
@@ -192,7 +192,7 @@ func ifft_vartime*[EC](
 proc init*(T: type FFTDescriptor, maxScale: uint8): T =
   result.maxWidth = 1 shl maxScale
 
-  let root = scaleToRootOfUnity(T.EC.F.Name)[maxScale]
+  let root = scaleToRootOfUnity(T.EC.getName())[maxScale]
   result.rootOfUnity = root.toBig()
   result.expandedRootsOfUnity = root.expandRootOfUnity()
     # Aren't you tired of reading about unity?
diff --git a/tests/math_elliptic_curves/t_ec_frobenius.nim b/tests/math_elliptic_curves/t_ec_frobenius.nim
index ee0c06495..a17e0bd35 100644
--- a/tests/math_elliptic_curves/t_ec_frobenius.nim
+++ b/tests/math_elliptic_curves/t_ec_frobenius.nim
@@ -278,7 +278,7 @@ suite "ψ²(P) - [t]ψ(P) + [p]P = Inf" & " [" & $WordBitWidth & "-bit words]":
       {.error: "Not implemented".}
 
   proc test(EC: typedesc, randZ: static bool, gen: static RandomGen) =
-    let trace = trace(EC.F.Name)
+    let trace = trace(EC.getName())
 
     for i in 0 ..< Iters:
       let P = rng.random_point(EC, randZ, gen)
diff --git a/tests/math_elliptic_curves/t_ec_sage_template.nim b/tests/math_elliptic_curves/t_ec_sage_template.nim
index 965fdc417..4753738f1 100644
--- a/tests/math_elliptic_curves/t_ec_sage_template.nim
+++ b/tests/math_elliptic_curves/t_ec_sage_template.nim
@@ -146,7 +146,7 @@ proc run_scalar_mul_test_vs_sage*(
   const coord = when EC is EC_ShortW_Prj: " Projective coordinates "
                 elif EC is EC_ShortW_Jac: " Jacobian coordinates "
 
-  const testSuiteDesc = "Scalar Multiplication " & $EC.F.Name & " " & G1_or_G2 & " " & coord & " vs SageMath - " & $bits & "-bit scalar"
+  const testSuiteDesc = "Scalar Multiplication " & $EC.getName() & " " & G1_or_G2 & " " & coord & " vs SageMath - " & $bits & "-bit scalar"
 
   suite testSuiteDesc & " [" & $WordBitWidth & "-bit words]":
     for i in 0 ..< vec.vectors.len:
diff --git a/tests/math_elliptic_curves/t_ec_template.nim b/tests/math_elliptic_curves/t_ec_template.nim
index b16e5d597..11b63c54a 100644
--- a/tests/math_elliptic_curves/t_ec_template.nim
+++ b/tests/math_elliptic_curves/t_ec_template.nim
@@ -105,7 +105,7 @@ proc run_EC_addition_tests*(
   echo "\n------------------------------------------------------\n"
   echo moduleName, " xoshiro512** seed: ", seed
 
-  const testSuiteDesc = "Elliptic curve in " & $ec.F.Name.getEquationForm() & " form"
+  const testSuiteDesc = "Elliptic curve in " & $ec.getName().getEquationForm() & " form"
 
   suite testSuiteDesc & " - " & $ec & " - [" & $WordBitWidth & "-bit mode]":
     test "The infinity point is the neutral element w.r.t. to EC " & " addition":
@@ -292,7 +292,7 @@ proc run_EC_addition_vartime_tests*(
   echo "\n------------------------------------------------------\n"
   echo moduleName, " xoshiro512** seed: ", seed
 
-  const testSuiteDesc = "Elliptic curve in " & $ec.F.Name.getEquationForm() & " form"
+  const testSuiteDesc = "Elliptic curve in " & $ec.getName().getEquationForm() & " form"
 
   suite testSuiteDesc & " - " & $ec & " (vartime) - [" & $WordBitWidth & "-bit mode]":
     test "The infinity point is the neutral element w.r.t. to EC " & $ec.G & " addition (vartime)":
@@ -479,7 +479,7 @@ proc run_EC_mul_sanity_tests*(
   echo "\n------------------------------------------------------\n"
   echo moduleName, " xoshiro512** seed: ", seed
 
-  const testSuiteDesc = "Elliptic curve in " & $ec.F.Name.getEquationForm() & " form"
+  const testSuiteDesc = "Elliptic curve in " & $ec.getName().getEquationForm() & " form"
 
   suite testSuiteDesc & " - " & $ec & " - [" & $WordBitWidth & "-bit mode]":
     test "EC " & " mul [0]P == Inf":
@@ -584,7 +584,7 @@ proc run_EC_mul_distributive_tests*(
   echo "\n------------------------------------------------------\n"
   echo moduleName, " xoshiro512** seed: ", seed
 
-  const testSuiteDesc = "Elliptic curve in " & $ec.F.Name.getEquationForm() & " form"
+  const testSuiteDesc = "Elliptic curve in " & $ec.getName().getEquationForm() & " form"
 
   suite testSuiteDesc & " - " & $ec & " - [" & $WordBitWidth & "-bit mode]":
 
@@ -647,7 +647,7 @@ proc run_EC_mul_vs_ref_impl*(
   echo "\n------------------------------------------------------\n"
   echo moduleName, " xoshiro512** seed: ", seed
 
-  const testSuiteDesc = "Elliptic curve in " & $ec.F.Name.getEquationForm() & " form"
+  const testSuiteDesc = "Elliptic curve in " & $ec.getName().getEquationForm() & " form"
 
   suite testSuiteDesc & " - " & $ec & " - [" & $WordBitWidth & "-bit mode]":
     test "EC " & $ec.G & " mul constant-time is equivalent to a simple double-and-add and recoded algorithms":
@@ -727,7 +727,7 @@ proc run_EC_mul_endomorphism_impl*(
   echo "\n------------------------------------------------------\n"
   echo moduleName, " xoshiro512** seed: ", seed
 
-  const testSuiteDesc = "Elliptic curve in " & $ec.F.Name.getEquationForm() & " form"
+  const testSuiteDesc = "Elliptic curve in " & $ec.getName().getEquationForm() & " form"
 
   suite testSuiteDesc & " - " & $ec & " - [" & $WordBitWidth & "-bit mode]":
     test "EC " & $ec.G & " multiplication with endomorphism":
@@ -740,7 +740,7 @@ proc run_EC_mul_endomorphism_impl*(
           let scalarUnreduced = rng.random_long01Seq(BigInt[bits])
           # Ensure scalar is smaller than curve order
           var scalar {.noInit.}: BigInt[bits]
-          discard scalar.limbs.reduce_vartime(scalarUnreduced.limbs, EC.F.Name.scalarFieldModulus().limbs)
+          discard scalar.limbs.reduce_vartime(scalarUnreduced.limbs, EC.getName().scalarFieldModulus().limbs)
 
           proc diagnostic(expected, computed: EC): string {.used.} =
             return "\n" &
@@ -1126,7 +1126,7 @@ proc run_EC_affine_conversion*(
   echo "\n------------------------------------------------------\n"
   echo moduleName, " xoshiro512** seed: ", seed
 
-  const testSuiteDesc = "Elliptic curve in " & $ec.F.Name.getEquationForm() & " form"
+  const testSuiteDesc = "Elliptic curve in " & $ec.getName().getEquationForm() & " form"
 
   suite testSuiteDesc & " - " & $ec & " - [" & $WordBitWidth & "-bit mode]":
     test "EC " & $ec.G & " batchAffine is consistent with single affine conversion":
diff --git a/tests/math_pairings/t_pairing_bls12_377_multi.nim b/tests/math_pairings/t_pairing_bls12_377_multi.nim
index b7924598d..6c9a08096 100644
--- a/tests/math_pairings/t_pairing_bls12_377_multi.nim
+++ b/tests/math_pairings/t_pairing_bls12_377_multi.nim
@@ -14,7 +14,7 @@ import
   constantine/math/[arithmetic, extension_fields, ec_shortweierstrass],
   constantine/math/io/io_extfields,
   constantine/named/algebras,
-  constantine/math/pairings/pairings_bls12,
+  constantine/math/pairings/pairings_generic,
   # Test utilities
   helpers/prng_unsafe
 
@@ -42,7 +42,7 @@ proc testMultiPairing(rng: var RngState, N: static int) =
   let clockSimpleStart = cpuTime()
   var GTsimple {.noInit.}: Fp12[BLS12_381]
   for i in 0 ..< N:
-    GTs[i].pairing_bls12(Ps[i], Qs[i])
+    GTs[i].pairing(Ps[i], Qs[i])
 
   GTsimple = GTs[0]
   for i in 1 ..< N:
@@ -52,7 +52,7 @@ proc testMultiPairing(rng: var RngState, N: static int) =
   # Multipairing
   let clockMultiStart = cpuTime()
   var GTmulti {.noInit.}: Fp12[BLS12_381]
-  GTmulti.pairing_bls12(Ps, Qs)
+  GTmulti.pairing(Ps, Qs)
   let clockMultiStop = cpuTime()
 
   echo &"N={N}, Simple: {clockSimpleStop - clockSimpleStart:>4.4f}s, Multi: {clockMultiStop - clockMultiStart:>4.4f}s"
diff --git a/tests/math_pairings/t_pairing_bls12_381_multi.nim b/tests/math_pairings/t_pairing_bls12_381_multi.nim
index 5be83c9ce..45c765ca7 100644
--- a/tests/math_pairings/t_pairing_bls12_381_multi.nim
+++ b/tests/math_pairings/t_pairing_bls12_381_multi.nim
@@ -14,7 +14,7 @@ import
   constantine/math/[arithmetic, extension_fields, ec_shortweierstrass],
   constantine/math/io/io_extfields,
   constantine/named/algebras,
-  constantine/math/pairings/pairings_bls12,
+  constantine/math/pairings/pairings_generic,
   # Test utilities
   helpers/prng_unsafe
 
@@ -42,7 +42,7 @@ proc testMultiPairing(rng: var RngState, N: static int) =
   let clockSimpleStart = cpuTime()
   var GTsimple {.noInit.}: Fp12[BLS12_381]
   for i in 0 ..< N:
-    GTs[i].pairing_bls12(Ps[i], Qs[i])
+    GTs[i].pairing(Ps[i], Qs[i])
 
   GTsimple = GTs[0]
   for i in 1 ..< N:
@@ -52,7 +52,7 @@ proc testMultiPairing(rng: var RngState, N: static int) =
   # Multipairing
   let clockMultiStart = cpuTime()
   var GTmulti {.noInit.}: Fp12[BLS12_381]
-  GTmulti.pairing_bls12(Ps, Qs)
+  GTmulti.pairing(Ps, Qs)
   let clockMultiStop = cpuTime()
 
   echo &"N={N}, Simple: {clockSimpleStop - clockSimpleStart:>4.4f}s, Multi: {clockMultiStop - clockMultiStart:>4.4f}s"
diff --git a/tests/math_pairings/t_pairing_bn254_nogami_multi.nim b/tests/math_pairings/t_pairing_bn254_nogami_multi.nim
index dbf3e26ea..10045ef22 100644
--- a/tests/math_pairings/t_pairing_bn254_nogami_multi.nim
+++ b/tests/math_pairings/t_pairing_bn254_nogami_multi.nim
@@ -12,9 +12,8 @@ import
   # Internals
   constantine/platforms/abstractions,
   constantine/math/[arithmetic, extension_fields, ec_shortweierstrass],
-  constantine/math/io/io_extfields,
   constantine/named/algebras,
-  constantine/math/pairings/pairings_bn,
+  constantine/math/pairings/pairings_generic,
   # Test utilities
   helpers/prng_unsafe
 
@@ -42,7 +41,7 @@ proc testMultiPairing(rng: var RngState, N: static int) =
   let clockSimpleStart = cpuTime()
   var GTsimple {.noInit.}: Fp12[BN254_Nogami]
   for i in 0 ..< N:
-    GTs[i].pairing_bn(Ps[i], Qs[i])
+    GTs[i].pairing(Ps[i], Qs[i])
 
   GTsimple = GTs[0]
   for i in 1 ..< N:
@@ -52,7 +51,7 @@ proc testMultiPairing(rng: var RngState, N: static int) =
   # Multipairing
   let clockMultiStart = cpuTime()
   var GTmulti {.noInit.}: Fp12[BN254_Nogami]
-  GTmulti.pairing_bn(Ps, Qs)
+  GTmulti.pairing(Ps, Qs)
   let clockMultiStop = cpuTime()
 
   echo &"N={N}, Simple: {clockSimpleStop - clockSimpleStart:>4.4f}s, Multi: {clockMultiStop - clockMultiStart:>4.4f}s"
diff --git a/tests/math_pairings/t_pairing_bn254_snarks_multi.nim b/tests/math_pairings/t_pairing_bn254_snarks_multi.nim
index 33c65c050..e85323b5a 100644
--- a/tests/math_pairings/t_pairing_bn254_snarks_multi.nim
+++ b/tests/math_pairings/t_pairing_bn254_snarks_multi.nim
@@ -14,7 +14,7 @@ import
   constantine/math/[arithmetic, extension_fields, ec_shortweierstrass],
   constantine/math/io/io_extfields,
   constantine/named/algebras,
-  constantine/math/pairings/pairings_bn,
+  constantine/math/pairings/pairings_generic,
   # Test utilities
   helpers/prng_unsafe
 
@@ -42,7 +42,7 @@ proc testMultiPairing(rng: var RngState, N: static int) =
   let clockSimpleStart = cpuTime()
   var GTsimple {.noInit.}: Fp12[BN254_Snarks]
   for i in 0 ..< N:
-    GTs[i].pairing_bn(Ps[i], Qs[i])
+    GTs[i].pairing(Ps[i], Qs[i])
 
   GTsimple = GTs[0]
   for i in 1 ..< N:
@@ -52,7 +52,7 @@ proc testMultiPairing(rng: var RngState, N: static int) =
   # Multipairing
   let clockMultiStart = cpuTime()
   var GTmulti {.noInit.}: Fp12[BN254_Snarks]
-  GTmulti.pairing_bn(Ps, Qs)
+  GTmulti.pairing(Ps, Qs)
   let clockMultiStop = cpuTime()
 
   echo &"N={N}, Simple: {clockSimpleStop - clockSimpleStart:>4.4f}s, Multi: {clockMultiStop - clockMultiStart:>4.4f}s"
diff --git a/tests/t_hash_to_curve.nim b/tests/t_hash_to_curve.nim
index edf6b917e..94327080a 100644
--- a/tests/t_hash_to_curve.nim
+++ b/tests/t_hash_to_curve.nim
@@ -113,7 +113,7 @@ proc run_hash_to_curve_test(
     const G1_or_G2 = "G2"
   let vec = loadVectors(HashToCurveTest[EC_ShortW_Aff[EC.F, EC.G]], filename)
 
-  let testSuiteDesc = "Hash to Curve " & $EC.F.Name & " " & G1_or_G2 & " - official specs " & spec_version & " test vectors"
+  let testSuiteDesc = "Hash to Curve " & $EC.getName() & " " & G1_or_G2 & " - official specs " & spec_version & " test vectors"
 
   suite testSuiteDesc & " [" & $WordBitWidth & "-bit words]":
 
@@ -148,7 +148,7 @@ proc run_hash_to_curve_svdw_test(
     const G1_or_G2 = "G2"
   let vec = loadVectors(HashToCurveTest[EC_ShortW_Aff[EC.F, EC.G]], filename)
 
-  let testSuiteDesc = "Hash to Curve " & $EC.F.Name & " " & G1_or_G2 & " - official specs " & spec_version & " test vectors"
+  let testSuiteDesc = "Hash to Curve " & $EC.getName() & " " & G1_or_G2 & " - official specs " & spec_version & " test vectors"
 
   suite testSuiteDesc & " [" & $WordBitWidth & "-bit words]":