in FourQ_ARM_side_channel/eccp2.c [875:1006]
bool ecc_mul_double(digit_t* k, point_t Q, digit_t* l, point_t R)
{ // Double scalar multiplication R = k*G + l*Q, where the G is the generator. Uses DOUBLE_SCALAR_TABLE, which contains multiples of G, Phi(G), Psi(G) and Phi(Psi(G)).
// Inputs: point Q in affine coordinates,
// scalars "k" and "l" in [0, 2^256-1].
// Output: R = k*G + l*Q in affine coordinates (x,y).
// The function uses wNAF with interleaving.
// SECURITY NOTE: this function is intended for a non-constant-time operation such as signature verification.
unsigned int position;
int i, digits_k1[65] = {0}, digits_k2[65] = {0}, digits_k3[65] = {0}, digits_k4[65] = {0};
int digits_l1[65] = {0}, digits_l2[65] = {0}, digits_l3[65] = {0}, digits_l4[65] = {0};
point_precomp_t V;
point_extedwards_t Q1, Q2, Q3, Q4, T;
point_extedwards_t U, Q_table1[NPOINTS_DOUBLEMUL_WQ], Q_table2[NPOINTS_DOUBLEMUL_WQ], Q_table3[NPOINTS_DOUBLEMUL_WQ], Q_table4[NPOINTS_DOUBLEMUL_WQ];
f2elm_t t0, t1;
uint64_t k_scalars[4], l_scalars[4];
point_setup(Q, Q1); // Convert to representation (X,Y,1,Ta,Tb)
if (ecc_point_validate(Q1) == false) { // Check if point lies on the curve
return false;
}
// Computing endomorphisms over point Q
ecccopy(Q1, Q2);
ecc_phi(Q2);
ecccopy(Q1, Q3);
ecc_psi(Q3);
ecccopy(Q2, Q4);
ecc_psi(Q4);
decompose((uint64_t*)k, k_scalars); // Scalar decomposition
decompose((uint64_t*)l, l_scalars);
wNAF_recode(k_scalars[0], WP_DOUBLEBASE, digits_k1); // Scalar recoding
wNAF_recode(k_scalars[1], WP_DOUBLEBASE, digits_k2);
wNAF_recode(k_scalars[2], WP_DOUBLEBASE, digits_k3);
wNAF_recode(k_scalars[3], WP_DOUBLEBASE, digits_k4);
wNAF_recode(l_scalars[0], WQ_DOUBLEBASE, digits_l1);
wNAF_recode(l_scalars[1], WQ_DOUBLEBASE, digits_l2);
wNAF_recode(l_scalars[2], WQ_DOUBLEBASE, digits_l3);
wNAF_recode(l_scalars[3], WQ_DOUBLEBASE, digits_l4);
ecc_precomp_double(Q1, Q_table1, NPOINTS_DOUBLEMUL_WQ); // Precomputation
ecc_precomp_double(Q2, Q_table2, NPOINTS_DOUBLEMUL_WQ);
ecc_precomp_double(Q3, Q_table3, NPOINTS_DOUBLEMUL_WQ);
ecc_precomp_double(Q4, Q_table4, NPOINTS_DOUBLEMUL_WQ);
fp2zero1271(T->x); // Initialize T as the neutral point (0:1:1)
fp2zero1271(T->y); T->y[0][0] = 1;
fp2zero1271(T->z); T->z[0][0] = 1;
for (i = 64; i >= 0; i--)
{
eccdouble(T); // Double (X_T,Y_T,Z_T,T_T) = 2(X_T,Y_T,Z_T,T_T)
if (digits_l1[i] < 0) {
position = (-digits_l1[i])/2;
eccneg_extproj_precomp(Q_table1[position], U); // Load and negate U = (X_U,Y_U,Z_U,Td_U) <- -(X+Y,X-Y,2Z,2dT) from a point in the precomputed table
eccadd_core(T, U, T, t0, t1); // T = T+U = (X_T,Y_T,Z_T,T_T) = (X_T,Y_T,Z_T,T_T) + (X_U,Y_U,Z_U,Td_U)
fp2mul1271(t0, t1, T->t);
} else if (digits_l1[i] > 0) {
position = (digits_l1[i])/2; // Take U = (X_U,Y_U,Z_U,Td_U) <- (X+Y,X-Y,2Z,2dT) from a point in the precomputed table
eccadd_core(T, Q_table1[position], T, t0, t1); // T = T+U = (X_T,Y_T,Z_T,T_T) = (X_T,Y_T,Z_T,T_T) + (X_U,Y_U,Z_U,Td_U)
fp2mul1271(t0, t1, T->t);
}
if (digits_l2[i] < 0) {
position = (-digits_l2[i])/2;
eccneg_extproj_precomp(Q_table2[position], U);
eccadd_core(T, U, T, t0, t1);
fp2mul1271(t0, t1, T->t);
} else if (digits_l2[i] > 0) {
position = (digits_l2[i])/2;
eccadd_core(T, Q_table2[position], T, t0, t1);
fp2mul1271(t0, t1, T->t);
}
if (digits_l3[i] < 0) {
position = (-digits_l3[i])/2;
eccneg_extproj_precomp(Q_table3[position], U);
eccadd_core(T, U, T, t0, t1);
fp2mul1271(t0, t1, T->t);
} else if (digits_l3[i] > 0) {
position = (digits_l3[i])/2;
eccadd_core(T, Q_table3[position], T, t0, t1);
fp2mul1271(t0, t1, T->t);
}
if (digits_l4[i] < 0) {
position = (-digits_l4[i])/2;
eccneg_extproj_precomp(Q_table4[position], U);
eccadd_core(T, U, T, t0, t1);
fp2mul1271(t0, t1, T->t);
} else if (digits_l4[i] > 0) {
position = (digits_l4[i])/2;
eccadd_core(T, Q_table4[position], T, t0, t1);
fp2mul1271(t0, t1, T->t);
}
if (digits_k1[i] < 0) {
position = (-digits_k1[i])/2;
eccneg_precomp(((point_precomp_t*)&DOUBLE_SCALAR_TABLE)[position], V); // Load and negate V = (X_V,Y_V,Z_V,Td_V) <- -(x+y,y-x,2dt) from a point in the precomputed table
eccmadd(V, T); // T = T+V = (X_T,Y_T,Z_T,Ta_T,Tb_T) = (X_T,Y_T,Z_T,Ta_T,Tb_T) + (X_V,Y_V,Z_V,Td_V)
} else if (digits_k1[i] > 0) {
position = (digits_k1[i])/2; // Take V = (X_V,Y_V,Z_V,Td_V) <- (x+y,y-x,2dt) from a point in the precomputed table
eccmadd(((point_precomp_t*)&DOUBLE_SCALAR_TABLE)[position], T); // T = T+V = (X_T,Y_T,Z_T,Ta_T,Tb_T) = (X_T,Y_T,Z_T,Ta_T,Tb_T) + (X_V,Y_V,Z_V,Td_V)
}
if (digits_k2[i] < 0) {
position = (-digits_k2[i])/2;
eccneg_precomp(((point_precomp_t*)&DOUBLE_SCALAR_TABLE)[NPOINTS_DOUBLEMUL_WP+position], V);
eccmadd(V, T);
} else if (digits_k2[i] > 0) {
position = (digits_k2[i])/2;
eccmadd(((point_precomp_t*)&DOUBLE_SCALAR_TABLE)[NPOINTS_DOUBLEMUL_WP+position], T);
}
if (digits_k3[i] < 0) {
position = (-digits_k3[i])/2;
eccneg_precomp(((point_precomp_t*)&DOUBLE_SCALAR_TABLE)[2*NPOINTS_DOUBLEMUL_WP+position], V);
eccmadd(V, T);
} else if (digits_k3[i] > 0) {
position = (digits_k3[i])/2;
eccmadd(((point_precomp_t*)&DOUBLE_SCALAR_TABLE)[2*NPOINTS_DOUBLEMUL_WP+position], T);
}
if (digits_k4[i] < 0) {
position = (-digits_k4[i])/2;
eccneg_precomp(((point_precomp_t*)&DOUBLE_SCALAR_TABLE)[3*NPOINTS_DOUBLEMUL_WP+position], V);
eccmadd(V, T);
} else if (digits_k4[i] > 0) {
position = (digits_k4[i])/2;
eccmadd(((point_precomp_t*)&DOUBLE_SCALAR_TABLE)[3*NPOINTS_DOUBLEMUL_WP+position], T);
}
}
eccnorm(T, R); // Output R = (x,y)
return true;
}