Loading...
Searching...
No Matches
fields.cpp
1#include <symengine/fields.h>
2#include <symengine/add.h>
4#include <symengine/mul.h>
5#include <symengine/pow.h>
6#include <symengine/symengine_exception.h>
7
8namespace SymEngine
9{
10GaloisField::GaloisField(const RCP<const Basic> &var, GaloisFieldDict &&dict)
11 : UIntPolyBase(var, std::move(dict))
12{
13 SYMENGINE_ASSIGN_TYPEID()
14 SYMENGINE_ASSERT(is_canonical(get_poly()))
15}
16
18{
19 // Check if dictionary contains terms with coeffienct 0
20 if (dict.modulo_ <= integer_class(0))
21 return false;
22 if (not dict.empty())
23 if (dict.dict_[dict.dict_.size() - 1] == integer_class(0))
24 return false;
25 return true;
26}
27
29{
30 hash_t seed = SYMENGINE_GALOISFIELD;
31
32 seed += get_var()->hash();
33 for (const auto &it : get_poly().dict_) {
34 hash_t temp = SYMENGINE_GALOISFIELD;
35 hash_combine<hash_t>(temp, mp_get_si(it));
36 seed += temp;
37 }
38 return seed;
39}
40
41int GaloisField::compare(const Basic &o) const
42{
43 const GaloisField &s = down_cast<const GaloisField &>(o);
44
45 if (get_poly().size() != s.get_poly().size())
46 return (get_poly().size() < s.get_poly().size()) ? -1 : 1;
47
48 int cmp = unified_compare(get_var(), s.get_var());
49 if (cmp != 0)
50 return cmp;
51
52 cmp = unified_compare(get_poly().modulo_, s.get_poly().modulo_);
53 if (cmp != 0)
54 return cmp;
55
56 return unified_compare(get_poly().dict_, s.get_poly().dict_);
57}
58
59RCP<const GaloisField> GaloisField::from_dict(const RCP<const Basic> &var,
61{
62 return make_rcp<const GaloisField>(var, std::move(d));
63}
64
65RCP<const GaloisField>
66GaloisField::from_vec(const RCP<const Basic> &var,
68 const integer_class &modulo)
69{
70 return make_rcp<const GaloisField>(var,
71 GaloisFieldDict::from_vec(v, modulo));
72}
73
74RCP<const GaloisField> GaloisField::from_uintpoly(const UIntPoly &a,
75 const integer_class &modulo)
76{
77 GaloisFieldDict wrapper(a.get_poly().get_dict(), modulo);
78 return GaloisField::from_dict(a.get_var(), std::move(wrapper));
79}
80
82{
83 vec_basic args;
84 if (get_poly().dict_.empty())
85 args.push_back(zero);
86 else {
87 for (unsigned i = 0; i < get_poly().dict_.size(); i++) {
88 if (get_poly().dict_[i] == integer_class(0))
89 continue;
90 if (i == 0) {
91 args.push_back(integer(get_poly().dict_[i]));
92 } else if (i == 1) {
93 if (get_poly().dict_[i] == 1) {
94 args.push_back(get_var());
95 } else {
96 args.push_back(Mul::from_dict(integer(get_poly().dict_[i]),
97 {{get_var(), one}}));
98 }
99 } else {
100 if (get_poly().dict_[i] == 1) {
101 args.push_back(pow(get_var(), integer(i)));
102 } else {
103 args.push_back(Mul::from_dict(integer(get_poly().dict_[i]),
104 {{get_var(), integer(i)}}));
105 }
106 }
107 }
108 }
109 return args;
110}
111
112GaloisFieldDict::GaloisFieldDict(const int &i, const integer_class &mod)
113 : modulo_(mod)
114{
115 integer_class temp;
116 mp_fdiv_r(temp, integer_class(i), modulo_);
117 if (temp != integer_class(0))
118 dict_.insert(dict_.begin(), temp);
119}
120
121GaloisFieldDict::GaloisFieldDict(const map_uint_mpz &p,
122 const integer_class &mod)
123 : modulo_(mod)
124{
125 if (p.size() != 0) {
126 dict_.resize(p.rbegin()->first + 1, integer_class(0));
127 for (auto &iter : p) {
128 integer_class temp;
129 mp_fdiv_r(temp, iter.second, modulo_);
130 dict_[iter.first] = temp;
131 }
132 gf_istrip();
133 }
134}
135
136GaloisFieldDict::GaloisFieldDict(const integer_class &i,
137 const integer_class &mod)
138 : modulo_(mod)
139{
140 integer_class temp;
141 mp_fdiv_r(temp, i, modulo_);
142 if (temp != integer_class(0))
143 dict_.insert(dict_.begin(), temp);
144}
145
146GaloisFieldDict GaloisFieldDict::from_vec(const std::vector<integer_class> &v,
147 const integer_class &modulo)
148{
149 GaloisFieldDict x;
150 x.modulo_ = modulo;
151 x.dict_.resize(v.size());
152 for (unsigned int i = 0; i < v.size(); ++i) {
153 integer_class a;
154 mp_fdiv_r(a, v[i], modulo);
155 x.dict_[i] = a;
156 }
157 x.gf_istrip();
158 return x;
159}
160
161GaloisFieldDict &GaloisFieldDict::negate()
162{
163 for (auto &a : dict_) {
164 a *= -1;
165 if (a != 0_z)
166 a += modulo_;
167 }
168 return down_cast<GaloisFieldDict &>(*this);
169}
170
171void GaloisFieldDict::gf_istrip()
172{
173 for (auto i = dict_.size(); i-- != 0;) {
174 if (dict_[i] == integer_class(0))
175 dict_.pop_back();
176 else
177 break;
178 }
179}
180
181GaloisFieldDict GaloisFieldDict::mul(const GaloisFieldDict &a,
182 const GaloisFieldDict &b)
183{
184 if (a.modulo_ != b.modulo_)
185 throw std::runtime_error("Error: field must be same.");
186 if (a.get_dict().empty())
187 return a;
188 if (b.get_dict().empty())
189 return b;
190
191 GaloisFieldDict p;
192 p.dict_.resize(a.degree() + b.degree() + 1, integer_class(0));
193 p.modulo_ = a.modulo_;
194 for (unsigned int i = 0; i <= a.degree(); i++)
195 for (unsigned int j = 0; j <= b.degree(); j++) {
196 auto temp = a.dict_[i];
197 temp *= b.dict_[j];
198 if (temp != integer_class(0)) {
199 auto t = p.dict_[i + j];
200 t += temp;
201 mp_fdiv_r(t, t, a.modulo_);
202 p.dict_[i + j] = t;
203 }
204 }
205 p.gf_istrip();
206 return p;
207}
208
209void GaloisFieldDict::gf_div(const GaloisFieldDict &o,
210 const Ptr<GaloisFieldDict> &quo,
211 const Ptr<GaloisFieldDict> &rem) const
212{
213 if (modulo_ != o.modulo_)
214 throw SymEngineException("Error: field must be same.");
215 if (o.dict_.empty())
216 throw DivisionByZeroError("ZeroDivisionError");
218 if (dict_.empty()) {
219 *quo = GaloisFieldDict::from_vec(dict_out, modulo_);
220 *rem = GaloisFieldDict::from_vec(dict_, modulo_);
221 return;
222 }
223 auto dict_divisor = o.dict_;
224 auto deg_dividend = this->degree();
225 auto deg_divisor = o.degree();
226 if (deg_dividend < deg_divisor) {
227 *quo = GaloisFieldDict::from_vec(dict_out, modulo_);
228 *rem = GaloisFieldDict::from_vec(dict_, modulo_);
229 } else {
230 dict_out = dict_;
231 integer_class inv;
232 mp_invert(inv, *(dict_divisor.rbegin()), modulo_);
233 integer_class coeff;
234 for (auto it = deg_dividend + 1; it-- != 0;) {
235 coeff = dict_out[it];
236 auto lb = deg_divisor + it > deg_dividend
237 ? deg_divisor + it - deg_dividend
238 : 0;
239 auto ub = std::min(it + 1, deg_divisor);
240 for (size_t j = lb; j < ub; ++j) {
241 mp_addmul(coeff, dict_out[it - j + deg_divisor],
242 -dict_divisor[j]);
243 }
244 if (it >= deg_divisor)
245 coeff *= inv;
246 mp_fdiv_r(coeff, coeff, modulo_);
247 dict_out[it] = coeff;
248 }
249 std::vector<integer_class> dict_rem, dict_quo;
250 dict_rem.resize(deg_divisor);
251 dict_quo.resize(deg_dividend - deg_divisor + 1);
252 for (unsigned it = 0; it < dict_out.size(); it++) {
253 if (it < deg_divisor)
254 dict_rem[it] = dict_out[it];
255 else
256 dict_quo[it - deg_divisor] = dict_out[it];
257 }
258 *quo = GaloisFieldDict::from_vec(dict_quo, modulo_);
259 *rem = GaloisFieldDict::from_vec(dict_rem, modulo_);
260 }
261}
262
263GaloisFieldDict GaloisFieldDict::gf_lshift(const integer_class n) const
264{
266 auto to_ret = GaloisFieldDict::from_vec(dict_out, modulo_);
267 if (!dict_.empty()) {
268 auto n_val = mp_get_ui(n);
269 to_ret.dict_.resize(n_val, integer_class(0));
270 to_ret.dict_.insert(to_ret.dict_.end(), dict_.begin(), dict_.end());
271 }
272 return to_ret;
273}
274
275void GaloisFieldDict::gf_rshift(const integer_class n,
276 const Ptr<GaloisFieldDict> &quo,
277 const Ptr<GaloisFieldDict> &rem) const
278{
280 *quo = GaloisFieldDict::from_vec(dict_quo, modulo_);
281 auto n_val = mp_get_ui(n);
282 if (n_val < dict_.size()) {
283 quo->dict_.insert(quo->dict_.end(), dict_.begin() + n_val, dict_.end());
284 std::vector<integer_class> dict_rem(dict_.begin(),
285 dict_.begin() + n_val);
286 *rem = GaloisFieldDict::from_vec(dict_rem, modulo_);
287 } else {
288 *rem = down_cast<const GaloisFieldDict &>(*this);
289 }
290}
291
292GaloisFieldDict GaloisFieldDict::gf_sqr() const
293{
294 return (*this * *this);
295}
296
297GaloisFieldDict GaloisFieldDict::gf_pow(const unsigned long n) const
298{
299 if (n == 0) {
300 return GaloisFieldDict({integer_class(1)}, modulo_);
301 }
302 if (n == 1)
303 return down_cast<const GaloisFieldDict &>(*this);
304 if (n == 2)
305 return gf_sqr();
306 auto num = n;
307 GaloisFieldDict to_sq = down_cast<const GaloisFieldDict &>(*this);
308 GaloisFieldDict to_ret = GaloisFieldDict({integer_class(1)}, modulo_);
309 while (1) {
310 if (num & 1) {
311 to_ret *= to_sq;
312 }
313 num >>= 1;
314 if (num == 0)
315 return to_ret;
316 to_sq = to_sq.gf_sqr();
317 }
318}
319
320void GaloisFieldDict::gf_monic(integer_class &res,
321 const Ptr<GaloisFieldDict> &monic) const
322{
323 *monic = down_cast<const GaloisFieldDict &>(*this);
324 if (dict_.empty()) {
325 res = integer_class(0);
326 } else {
327 res = *dict_.rbegin();
328 if (res != integer_class(1)) {
329 integer_class inv, temp;
330 mp_invert(inv, res, modulo_);
331 for (auto &iter : monic->dict_) {
332 temp = inv;
333 temp *= iter;
334 mp_fdiv_r(iter, temp, modulo_);
335 }
336 }
337 }
338}
339
340GaloisFieldDict GaloisFieldDict::gf_gcd(const GaloisFieldDict &o) const
341{
342 if (modulo_ != o.modulo_)
343 throw SymEngineException("Error: field must be same.");
344 GaloisFieldDict f = down_cast<const GaloisFieldDict &>(*this);
345 GaloisFieldDict g = o;
346 GaloisFieldDict temp_out;
347 while (not g.dict_.empty()) {
348 f %= g; // f, g = f % g, g
349 f.dict_.swap(g.dict_);
350 }
351 integer_class temp_LC;
352 f.gf_monic(temp_LC, outArg(f));
353 return f;
354}
355
356GaloisFieldDict GaloisFieldDict::gf_lcm(const GaloisFieldDict &o) const
357{
358 if (modulo_ != o.modulo_)
359 throw SymEngineException("Error: field must be same.");
360 if (dict_.empty())
361 return down_cast<const GaloisFieldDict &>(*this);
362 if (o.dict_.empty())
363 return o;
364 GaloisFieldDict out, temp_out;
365 out = o * (*this);
366 out /= gf_gcd(o);
367 integer_class temp_LC;
368 out.gf_monic(temp_LC, outArg(out));
369 return out;
370}
371
372GaloisFieldDict GaloisFieldDict::gf_diff() const
373{
374 auto df = degree();
375 GaloisFieldDict out = GaloisFieldDict({}, modulo_);
376 out.dict_.resize(df, integer_class(0));
377 for (unsigned i = 1; i <= df; i++) {
378 if (dict_[i] != integer_class(0)) {
379 out.dict_[i - 1] = i * dict_[i];
380 mp_fdiv_r(out.dict_[i - 1], out.dict_[i - 1], modulo_);
381 }
382 }
383 out.gf_istrip();
384 return out;
385}
386
387integer_class GaloisFieldDict::gf_eval(const integer_class &a) const
388{
389 integer_class res = 0_z;
390 for (auto rit = dict_.rbegin(); rit != dict_.rend(); ++rit) {
391 res *= a;
392 res += (*rit);
393 res %= modulo_;
394 }
395 return res;
396}
397
398vec_integer_class
399GaloisFieldDict::gf_multi_eval(const vec_integer_class &v) const
400{
401 vec_integer_class res(v.size());
402 for (unsigned int i = 0; i < v.size(); ++i)
403 res[i] = gf_eval(v[i]);
404 return res;
405}
406
407bool GaloisFieldDict::gf_is_sqf() const
408{
409 if (dict_.empty())
410 return true;
411 integer_class LC;
412 GaloisFieldDict monic;
413 gf_monic(LC, outArg(monic));
414 monic = monic.gf_gcd(monic.gf_diff());
415 return monic.is_one();
416}
417
419GaloisFieldDict::gf_sqf_list() const
420{
422 if (degree() < 1)
423 return vec_out;
424 unsigned n = 1;
425 // This cast is okay, because the multiplicities are unsigned
426 unsigned r = numeric_cast<unsigned>(mp_get_ui(modulo_));
427 bool sqf = false;
428 integer_class LC;
429 GaloisFieldDict f;
430 gf_monic(LC, outArg(f));
431 while (true) {
432 GaloisFieldDict F = f.gf_diff();
433 if (not F.dict_.empty()) {
434 GaloisFieldDict g = f.gf_gcd(F);
435 GaloisFieldDict h = f / g;
436
437 unsigned i = 1;
438
439 while (not h.is_one()) {
440 GaloisFieldDict G = h.gf_gcd(g);
441 GaloisFieldDict H = h / G;
442
443 if (H.degree() > 0)
444 vec_out.push_back({H, i * n});
445
446 ++i;
447 g /= G;
448 h = G;
449 }
450 if (g.is_one())
451 sqf = true;
452 else
453 f = g;
454 }
455 if (not sqf) {
456 auto deg = f.degree();
457 auto d = deg / r;
458 GaloisFieldDict temp = f;
459 for (unsigned int i = 0; i <= d; ++i) {
460 f.dict_[d - i] = temp.dict_[deg - i * r];
461 }
462 n *= r;
463 f.dict_.resize(d + 1);
464 f.gf_istrip();
465 } else
466 break;
467 }
468 return vec_out;
469}
470
471GaloisFieldDict GaloisFieldDict::gf_sqf_part() const
472{
473 auto sqf = gf_sqf_list();
474 GaloisFieldDict g = GaloisFieldDict::from_vec({1_z}, modulo_);
475
476 for (auto &f : sqf)
477 g *= f.first;
478
479 return g;
480}
481
482GaloisFieldDict GaloisFieldDict::gf_compose_mod(const GaloisFieldDict &g,
483 const GaloisFieldDict &h) const
484{
485 if (g.modulo_ != h.modulo_)
486 throw SymEngineException("Error: field must be same.");
487 if (g.modulo_ != modulo_)
488 throw SymEngineException("Error: field must be same.");
489 if (g.dict_.size() == 0)
490 return g;
491 GaloisFieldDict out
492 = GaloisFieldDict::from_vec({*(g.dict_.rbegin())}, modulo_);
493 if (g.dict_.size() >= 2) {
494 for (auto i = g.dict_.size() - 2;; --i) {
495 out *= h;
496 out += g.dict_[i];
497 out %= (*this);
498 if (i == 0)
499 break;
500 }
501 }
502 return out;
503}
504
505GaloisFieldDict GaloisFieldDict::gf_pow_mod(const GaloisFieldDict &f,
506 const unsigned long &n) const
507{
508 if (modulo_ != f.modulo_)
509 throw SymEngineException("Error: field must be same.");
510 if (n == 0)
511 return GaloisFieldDict::from_vec({1_z}, modulo_);
512 GaloisFieldDict in = f;
513 if (n == 1) {
514 return f % (*this);
515 }
516 if (n == 2) {
517 return f.gf_sqr() % (*this);
518 }
519 GaloisFieldDict h = GaloisFieldDict::from_vec({1_z}, modulo_);
520 auto mod = n;
521 while (true) {
522 if (mod & 1) {
523 h *= in;
524 h %= *this;
525 }
526 mod >>= 1;
527
528 if (mod == 0)
529 break;
530
531 in = in.gf_sqr() % *this;
532 }
533 return h;
534}
535
536std::vector<GaloisFieldDict> GaloisFieldDict::gf_frobenius_monomial_base() const
537{
538 auto n = degree();
540 if (n == 0)
541 return b;
542 b.resize(n);
543 b[0] = GaloisFieldDict::from_vec({1_z}, modulo_);
544 GaloisFieldDict temp_out;
545 if (mp_get_ui(modulo_) < n) {
546 for (unsigned i = 1; i < n; ++i) {
547 b[i] = b[i - 1].gf_lshift(modulo_);
548 b[i] %= (*this);
549 }
550 } else if (n > 1) {
551 b[1] = gf_pow_mod(GaloisFieldDict::from_vec({0_z, 1_z}, modulo_),
552 mp_get_ui(modulo_));
553 for (unsigned i = 2; i < n; ++i) {
554 b[i] = b[i - 1] * b[1];
555 b[i] %= (*this);
556 }
557 }
558 return b;
559}
560
561GaloisFieldDict
562GaloisFieldDict::gf_frobenius_map(const GaloisFieldDict &g,
563 const std::vector<GaloisFieldDict> &b) const
564{
565 if (modulo_ != g.modulo_)
566 throw SymEngineException("Error: field must be same.");
567 auto m = g.degree();
568 GaloisFieldDict temp_out(*this), out;
569 if (this->degree() >= m) {
570 temp_out %= g;
571 }
572 if (temp_out.empty()) {
573 return temp_out;
574 }
575 m = temp_out.degree();
576 out = GaloisFieldDict::from_vec({temp_out.dict_[0]}, modulo_);
577 for (unsigned i = 1; i <= m; ++i) {
578 auto v = b[i];
579 v *= temp_out.dict_[i];
580 out += v;
581 }
582 out.gf_istrip();
583 return out;
584}
585
586std::pair<GaloisFieldDict, GaloisFieldDict> GaloisFieldDict::gf_trace_map(
587 const GaloisFieldDict &a, const GaloisFieldDict &b,
588 const GaloisFieldDict &c, const unsigned long &n) const
589{
590 unsigned long n_val(n);
591 auto u = this->gf_compose_mod(a, b);
592 GaloisFieldDict v(b), U, V;
593 if (n_val & 1) {
594 U = a + u;
595 V = b;
596 } else {
597 U = a;
598 V = c;
599 }
600 n_val >>= 1;
601 while (n_val) {
602 u += this->gf_compose_mod(u, v);
603 v = gf_compose_mod(v, v);
604
605 if (n_val & 1) {
606 auto temp = gf_compose_mod(u, V);
607 U += temp;
608 V = gf_compose_mod(v, V);
609 }
610 n_val >>= 1;
611 }
612 return std::make_pair(gf_compose_mod(a, V), U);
613}
614
615GaloisFieldDict
616GaloisFieldDict::_gf_trace_map(const GaloisFieldDict &f, const unsigned long &n,
617 const std::vector<GaloisFieldDict> &b) const
618{
619 GaloisFieldDict x = f % (*this);
620 auto h = f;
621 auto r = f;
622 for (unsigned i = 1; i < n; ++i) {
623 h = gf_frobenius_map(h, b);
624 r += h;
625 r %= (*this);
626 }
627 return r;
628}
629
631GaloisFieldDict::gf_ddf_zassenhaus() const
632{
633 unsigned i = 1;
634 GaloisFieldDict f(*this);
635 GaloisFieldDict g = GaloisFieldDict::from_vec({0_z, 1_z}, modulo_);
636 GaloisFieldDict to_sub(g);
638
639 auto b = f.gf_frobenius_monomial_base();
640 while (2 * i <= f.degree()) {
641 g = g.gf_frobenius_map(f, b);
642
643 GaloisFieldDict h = f.gf_gcd(g - to_sub);
644
645 if (not h.is_one()) {
646 factors.push_back({h, i});
647 f /= h;
648 g %= f;
649 b = f.gf_frobenius_monomial_base();
650 }
651 ++i;
652 }
653 if (not(f.is_one() || f.empty())) {
654 factors.push_back({f, f.degree()});
655 }
656 return factors;
657}
658
659GaloisFieldDict
660GaloisFieldDict::_gf_pow_pnm1d2(const GaloisFieldDict &f, const unsigned &n,
661 const std::vector<GaloisFieldDict> &b) const
662{
663 GaloisFieldDict f_in(f);
664 f_in %= *this;
665 GaloisFieldDict h, r;
666 h = r = f_in;
667 for (unsigned i = 1; i < n; ++i) {
668 h = h.gf_frobenius_map(*this, b);
669 r *= h;
670 r %= *this;
671 }
672 auto res = gf_pow_mod(r, (mp_get_ui(modulo_) - 1) / 2);
673 return res;
674}
675
676GaloisFieldDict GaloisFieldDict::gf_random(const unsigned int &n_val,
677 mp_randstate &state) const
678{
679 std::vector<integer_class> v(n_val + 1);
680 for (unsigned i = 0; i < n_val; ++i) {
681 state.urandomint(v[i], modulo_);
682 }
683 v[n_val] = 1_z;
684 return GaloisFieldDict::from_vec(v, modulo_);
685}
686
688GaloisFieldDict::gf_edf_zassenhaus(const unsigned &n) const
689{
691 factors.insert(*this);
692 if (this->degree() <= n)
693 return factors;
694
695 auto N = this->degree() / n;
696
698 if (modulo_ != 2_z)
699 b = this->gf_frobenius_monomial_base();
700 mp_randstate state;
701 while (factors.size() < N) {
702 auto r = gf_random(2 * n - 1, state);
703 GaloisFieldDict g;
704 if (modulo_ == 2_z) {
705 GaloisFieldDict h = r;
706 unsigned ub = 1 << (n * N - 1);
707 for (unsigned i = 0; i < ub; ++i) {
708 r = gf_pow_mod(r, 2);
709 h += r;
710 }
711 g = this->gf_gcd(h);
712 } else {
713 GaloisFieldDict h = _gf_pow_pnm1d2(r, n, b);
714 h -= 1_z;
715 g = this->gf_gcd(h);
716 }
717
718 if (!g.is_one() and g != (*this)) {
719 factors = g.gf_edf_zassenhaus(n);
720 auto to_add = ((*this) / g).gf_edf_zassenhaus(n);
721 if (not to_add.empty())
722 factors.insert(to_add.begin(), to_add.end());
723 }
724 }
725 return factors;
726}
727
729GaloisFieldDict::gf_ddf_shoup() const
730{
732 if (dict_.empty())
733 return factors;
734 GaloisFieldDict f(*this);
735 auto n = this->degree();
736 auto k = static_cast<unsigned>(std::ceil(std::sqrt(n / 2)));
737 auto b = gf_frobenius_monomial_base();
738 auto x = GaloisFieldDict::from_vec({0_z, 1_z}, modulo_);
739 auto h = x.gf_frobenius_map(f, b);
740
742 U.push_back(x);
743 U.push_back(h);
744 U.resize(k + 1);
745 for (unsigned i = 2; i <= k; ++i)
746 U[i] = U[i - 1].gf_frobenius_map(*this, b);
747 h = U[k];
748 U.resize(k);
750 V.push_back(h);
751 V.resize(k);
752 for (unsigned i = 1; i + 1 <= k; ++i)
753 V[i] = this->gf_compose_mod(V[i - 1], h);
754 for (unsigned i = 0; i < V.size(); i++) {
755 h = GaloisFieldDict::from_vec({1_z}, modulo_);
756 auto j = k - 1;
757 GaloisFieldDict g;
758 for (auto &u : U) {
759 g = V[i] - u;
760 h *= g;
761 h %= f;
762 }
763 g = f.gf_gcd(h);
764 f /= g;
765 for (auto rit = U.rbegin(); rit != U.rend(); ++rit) {
766 h = V[i] - (*rit);
767 auto F = g.gf_gcd(h);
768 if (not F.is_one()) {
769 unsigned temp = k * (i + 1) - j;
770 factors.push_back({F, temp});
771 }
772 g /= F;
773 --j;
774 }
775 }
776 if (not f.is_one())
777 factors.push_back({f, f.degree()});
778 return factors;
779}
780
782GaloisFieldDict::gf_edf_shoup(const unsigned &n) const
783{
784 auto N = this->degree();
786 if (N <= n) {
787 if (N != 0)
788 factors.insert(*this);
789 return factors;
790 }
791 auto x = GaloisFieldDict::from_vec({0_z, 1_z}, modulo_);
792 mp_randstate state;
793 auto r = gf_random(N - 1, state);
794 if (modulo_ == 2_z) {
795 auto h = gf_pow_mod(x, mp_get_ui(modulo_));
796 auto H = gf_trace_map(r, h, x, n - 1).second;
797 auto h1 = gf_gcd(H);
798 auto h2 = (*this) / h1;
799 factors = h1.gf_edf_shoup(n);
800 auto temp = h2.gf_edf_shoup(n);
801 factors.insert(temp.begin(), temp.end());
802 } else {
803 auto b = gf_frobenius_monomial_base();
804 auto H = _gf_trace_map(r, n, b);
805 auto h = gf_pow_mod(H, (mp_get_ui(modulo_) - 1) / 2);
806 auto h1 = gf_gcd(h);
807 auto h2 = gf_gcd(h - 1_z);
808 auto h3 = (*this) / (h1 * h2);
809 factors = h1.gf_edf_shoup(n);
810 auto temp = h2.gf_edf_shoup(n);
811 factors.insert(temp.begin(), temp.end());
812 temp = h3.gf_edf_shoup(n);
813 factors.insert(temp.begin(), temp.end());
814 }
815 return factors;
816}
817
819GaloisFieldDict::gf_zassenhaus() const
820{
822 auto temp1 = gf_ddf_zassenhaus();
823 for (auto &f : temp1) {
824 auto temp2 = f.first.gf_edf_zassenhaus(f.second);
825 factors.insert(temp2.begin(), temp2.end());
826 }
827 return factors;
828}
829
831GaloisFieldDict::gf_shoup() const
832{
834 auto temp1 = gf_ddf_shoup();
835 for (auto &f : temp1) {
836 auto temp2 = f.first.gf_edf_shoup(f.second);
837 factors.insert(temp2.begin(), temp2.end());
838 }
839 return factors;
840}
841
843 GaloisFieldDict::DictLess>>
844GaloisFieldDict::gf_factor() const
845{
846 integer_class lc;
848 GaloisFieldDict monic;
849 gf_monic(lc, outArg(monic));
850 if (monic.degree() < 1)
851 return std::make_pair(lc, factors);
853 = monic.gf_sqf_list();
854 for (auto a : sqf_list) {
855 auto temp = (a.first).gf_zassenhaus();
856 for (auto f : temp)
857 factors.insert({f, a.second});
858 }
859 return std::make_pair(lc, factors);
860}
861} // namespace SymEngine
Classes and functions relating to the binary operation of addition.
T ceil(T... args)
The lowest unit of symbolic representation.
Definition: basic.h:97
int compare(const Basic &o) const
Definition: fields.cpp:41
virtual vec_basic get_args() const
Returns the list of arguments.
Definition: fields.cpp:81
GaloisField(const RCP< const Basic > &var, GaloisFieldDict &&dict)
Constructor of GaloisField class.
Definition: fields.cpp:10
bool is_canonical(const GaloisFieldDict &dict) const
Definition: fields.cpp:17
hash_t __hash__() const
Definition: fields.cpp:28
static RCP< const Basic > from_dict(const RCP< const Number > &coef, map_basic_basic &&d)
Create a Mul from a dict.
Definition: mul.cpp:115
T insert(T... args)
T make_pair(T... args)
T min(T... args)
T move(T... args)
Main namespace for SymEngine package.
Definition: add.cpp:19
RCP< const Integer > mod(const Integer &n, const Integer &d)
modulo round toward zero
Definition: ntheory.cpp:66
int unified_compare(const T &a, const T &b)
Definition: dict.h:205
std::enable_if< std::is_integral< T >::value, RCP< constInteger > >::type integer(T i)
Definition: integer.h:200
STL namespace.
T push_back(T... args)
T resize(T... args)
T size(T... args)
T sqrt(T... args)