generateMds<F extends PastaFieldElement<F>> static method

Implementation

static MdsGenerateResult<F> generateMds<F extends PastaFieldElement<F>>(
  Grain<F> grain,
  int t,
  int select,
  F one,
  F zero,
) {
  List<F> xs = [];
  List<F> ys = [];
  List<List<F>> mds = [];
  int counter = select;
  while (true) {
    // Generate two lists of unique field elements
    while (true) {
      List<F> vals = List.generate(
        2 * t,
        (_) => grain.nextFieldElementWithoutRejection(),
      );

      // Check uniqueness
      List<F> unique = List.from(vals);
      unique.sort((a, b) => a.compareTo(b));
      unique = unique.toSet().toList();

      if (vals.length == unique.length) {
        xs = vals.sublist(0, t);
        ys = vals.sublist(t);
        break;
      }
    }

    // Skip if select counter not yet zero
    if (counter != 0) {
      counter -= 1;
      continue;
    }

    // Generate MDS matrix: a_ij = 1 / (x_i + y_j)
    mds = List.generate(
      t,
      (i) => List.generate(t, (j) {
        final sum = xs[i] + ys[j];
        final s = sum.invert();
        if (s == null) {
          throw PoseidonException("Division by zero.");
        }
        return s;
      }),
    );

    break;
  }

  // Compute inverse MDS
  List<List<F>> mdsInv = List.generate(
    t,
    (_) => List.generate(t, (_) => zero),
  );

  F lagrange(List<F> arr, int j, F x) {
    final xj = arr[j];
    F acc = one;
    for (int m = 0; m < arr.length; m++) {
      if (m == j) continue;
      final denom = xj - arr[m];
      final denomInv = denom.invert();
      if (denomInv == null) {
        throw PoseidonException("Division by zero.");
      }
      acc = acc * (x - arr[m]) * denomInv;
    }
    return acc;
  }

  // Negate ys for positive Cauchy formulation
  final negYs = ys.map((y) => -y).toList();

  for (int i = 0; i < t; i++) {
    for (int j = 0; j < t; j++) {
      mdsInv[i][j] =
          (xs[j] - negYs[i]) *
          lagrange(xs, j, negYs[i]) *
          lagrange(negYs, i, xs[j]);
    }
  }
  return MdsGenerateResult(mds: mds, mdsInv: mdsInv);
}