21 constexpr limb_t MAX_PRIME_DIFF = 1103717;
24 inline void extract3(limb_t& c0, limb_t& c1, limb_t& c2, limb_t& n)
33 inline void mul(limb_t& c0, limb_t& c1,
const limb_t& a,
const limb_t& b)
35 double_limb_t t = (double_limb_t)a * b;
41 inline void mulnadd3(limb_t& c0, limb_t& c1, limb_t& c2, limb_t& d0, limb_t& d1, limb_t& d2,
const limb_t& n)
43 double_limb_t t = (double_limb_t)d0 * n + c0;
46 t += (double_limb_t)d1 * n + c1;
53 inline void muln2(limb_t& c0, limb_t& c1,
const limb_t& n)
55 double_limb_t t = (double_limb_t)c0 * n;
58 t += (double_limb_t)c1 * n;
63 inline void muladd3(limb_t& c0, limb_t& c1, limb_t& c2,
const limb_t& a,
const limb_t& b)
65 double_limb_t t = (double_limb_t)a * b;
66 limb_t th = t >> LIMB_SIZE;
70 th += (c0 < tl) ? 1 : 0;
72 c2 += (c1 < th) ? 1 : 0;
76 inline void muldbladd3(limb_t& c0, limb_t& c1, limb_t& c2,
const limb_t& a,
const limb_t& b)
78 double_limb_t t = (double_limb_t)a * b;
79 limb_t th = t >> LIMB_SIZE;
83 limb_t tt = th + ((c0 < tl) ? 1 : 0);
85 c2 += (c1 < tt) ? 1 : 0;
87 th += (c0 < tl) ? 1 : 0;
89 c2 += (c1 < th) ? 1 : 0;
96 inline void addnextract2(limb_t& c0, limb_t& c1,
const limb_t& a, limb_t& n)
117 inline void square_n_mul(
Num3072& in_out,
const int sq,
const Num3072& mul)
119 for (
int j = 0; j < sq; ++j) in_out.
Square();
128 if (this->
limbs[0] <= std::numeric_limits<limb_t>::max() - MAX_PRIME_DIFF)
return false;
129 for (
int i = 1; i <
LIMBS; ++i) {
130 if (this->
limbs[i] != std::numeric_limits<limb_t>::max())
return false;
137 limb_t c0 = MAX_PRIME_DIFF;
139 for (
int i = 0; i <
LIMBS; ++i) {
140 addnextract2(c0, c1, this->
limbs[i], this->
limbs[i]);
155 for (
int i = 0; i < 11; ++i) {
157 for (
int j = 0; j < (1 << i); ++j) p[i + 1].
Square();
163 square_n_mul(out, 512, p[9]);
164 square_n_mul(out, 256, p[8]);
165 square_n_mul(out, 128, p[7]);
166 square_n_mul(out, 64, p[6]);
167 square_n_mul(out, 32, p[5]);
168 square_n_mul(out, 8, p[3]);
169 square_n_mul(out, 2, p[1]);
170 square_n_mul(out, 1, p[0]);
171 square_n_mul(out, 5, p[2]);
172 square_n_mul(out, 3, p[0]);
173 square_n_mul(out, 2, p[0]);
174 square_n_mul(out, 4, p[0]);
175 square_n_mul(out, 4, p[1]);
176 square_n_mul(out, 3, p[0]);
183 limb_t c0 = 0, c1 = 0, c2 = 0;
187 for (
int j = 0; j <
LIMBS - 1; ++j) {
188 limb_t d0 = 0, d1 = 0, d2 = 0;
190 for (
int i = 2 + j; i <
LIMBS; ++i) muladd3(d0, d1, d2, this->
limbs[i], a.
limbs[
LIMBS + j - i]);
191 mulnadd3(c0, c1, c2, d0, d1, d2, MAX_PRIME_DIFF);
192 for (
int i = 0; i < j + 1; ++i) muladd3(c0, c1, c2, this->
limbs[i], a.
limbs[j - i]);
193 extract3(c0, c1, c2, tmp.
limbs[j]);
202 muln2(c0, c1, MAX_PRIME_DIFF);
203 for (
int j = 0; j <
LIMBS; ++j) {
204 addnextract2(c0, c1, tmp.
limbs[j], this->limbs[j]);
208 assert(c0 == 0 || c0 == 1);
220 limb_t c0 = 0, c1 = 0, c2 = 0;
224 for (
int j = 0; j <
LIMBS - 1; ++j) {
225 limb_t d0 = 0, d1 = 0, d2 = 0;
226 for (
int i = 0; i < (
LIMBS - 1 - j) / 2; ++i) muldbladd3(d0, d1, d2, this->
limbs[i + j + 1], this->
limbs[LIMBS - 1 - i]);
227 if ((j + 1) & 1) muladd3(d0, d1, d2, this->
limbs[(LIMBS - 1 - j) / 2 + j + 1], this->
limbs[LIMBS - 1 - (
LIMBS - 1 - j) / 2]);
228 mulnadd3(c0, c1, c2, d0, d1, d2, MAX_PRIME_DIFF);
229 for (
int i = 0; i < (j + 1) / 2; ++i) muldbladd3(c0, c1, c2, this->
limbs[i], this->
limbs[j - i]);
230 if ((j + 1) & 1) muladd3(c0, c1, c2, this->
limbs[(j + 1) / 2], this->
limbs[j - (j + 1) / 2]);
231 extract3(c0, c1, c2, tmp.
limbs[j]);
235 for (
int i = 0; i <
LIMBS / 2; ++i) muldbladd3(c0, c1, c2, this->
limbs[i], this->
limbs[LIMBS - 1 - i]);
239 muln2(c0, c1, MAX_PRIME_DIFF);
240 for (
int j = 0; j <
LIMBS; ++j) {
241 addnextract2(c0, c1, tmp.
limbs[j], this->limbs[j]);
245 assert(c0 == 0 || c0 == 1);
258 for (
int i = 1; i <
LIMBS; ++i) this->
limbs[i] = 0;
279 for (
int i = 0; i <
LIMBS; ++i) {
280 if (
sizeof(
limb_t) == 4) {
282 }
else if (
sizeof(
limb_t) == 8) {
289 for (
int i = 0; i <
LIMBS; ++i) {
290 if (
sizeof(
limb_t) == 4) {
292 }
else if (
sizeof(
limb_t) == 8) {
310 m_numerator = ToNum3072(in);
315 m_numerator.Divide(m_denominator);
316 m_denominator.SetToOne();
319 m_numerator.ToBytes(data);
326 m_numerator.Multiply(mul.m_numerator);
327 m_denominator.Multiply(mul.m_denominator);
333 m_numerator.Multiply(div.m_denominator);
334 m_denominator.Multiply(div.m_numerator);
339 m_numerator.Multiply(ToNum3072(in));
344 m_denominator.Multiply(ToNum3072(in));