in src/shamir.c [25:96]
static void SSS_lagrange_coefficients(int k, const octet* X, BIG_256_56* lc, const BIG_256_56 q)
{
int i;
BIG_256_56 x2[k];
DBIG_256_56 w;
for(i = 0; i < k; i++)
{
BIG_256_56_fromBytesLen(x2[i], X[i].val, X[i].len);
}
// Compute numerators in place using partial products
// to achieve it in O(n)
// c_i = x_0 * ... * x_(i-1) * x_(i+1) * ... * x_(k-1)
// Compute partial left products
// leave c_0 alone since it only has a right partial product
BIG_256_56_copy(lc[1], x2[0]);
for(i = 2; i < k; i++)
{
// lp_i = x_0 * ... * x_(i-1) = lp_(i-1) * x_(i-1)
BIG_256_56_mul(w, lc[i-1], x2[i-1]);
BIG_256_56_dmod(lc[i], w, q);
}
// Compute partial right products and combine
// Store partial right products in c_0 so at the end
// of the procedure c_0 = x_1 * ... x_(k-1)
BIG_256_56_copy(lc[0], x2[k-1]);
for(i = k-2; i > 0; i--)
{
// c_i = lp_i * rp_i
BIG_256_56_mul(w, lc[i], lc[0]);
BIG_256_56_dmod(lc[i], w, q);
// rp_(i-1) = x_i * ... * x_k = x_i * rp_i
BIG_256_56_mul(w, lc[0], x2[i]);
BIG_256_56_dmod(lc[0], w, q);
}
BIG_256_56 cneg;
BIG_256_56 denominator;
BIG_256_56 s;
for(i = 0; i < k; i++)
{
BIG_256_56_one(denominator);
// cneg = -x_i mod r
BIG_256_56_sub(cneg, q, x2[i]);
for(int j = 0; j < k; j++)
{
if (i == j) continue;
// denominator = denominator * (x_j - x_i)
BIG_256_56_add(s, x2[j], cneg);
BIG_256_56_norm(s);
BIG_256_56_mul(w, denominator, s);
BIG_256_56_dmod(denominator, w, q);
}
BIG_256_56_invmodp(denominator, denominator, q);
BIG_256_56_mul(w, lc[i], denominator);
BIG_256_56_dmod(lc[i], w, q);
}
}