Basic Polynomial Algebra Subprograms (BPAS)  v. 1.791
SmallPrimeField.hpp
1 
2 #ifndef _SMALL_PRIME_FIELD_H
3 #define _SMALL_PRIME_FIELD_H
4 
5 #include "BPASFiniteField.hpp"
6 #include <iostream>
7 #include <sstream>
8 #include <math.h>
9 
10 using std::cout;
11 using std::endl;
12 
13 //forward declarations
14 class Integer;
15 class RationalNumber;
17 class BigPrimeField;
21 template <class Ring>
23 
24 //NOTE the below code is not properly formatted with the corresponding
25 //implementation file. Only the montegomery version is setup properly.
26 #if ( 0 )
27 // simple small prime field using long long int to represent the number
28 class SmallPrimeField : public BPASField{
29 
30 private:
31 
32  static long int prime;
33  long long int a;
34 
35 public:
36 
37  static mpz_class characteristic;
38  static bool isPrimeField;
39  static bool isSmallPrimeField;
40  static bool isComplexField;
41 
42 
43  SmallPrimeField () : a(0) {}
44  SmallPrimeField (long long int _a) {
45  while (_a < 0){
46  _a = _a + prime;
47  }
48  a = _a%prime;
49  }
50  SmallPrimeField (const SmallPrimeField& c) {
51  a = c.a;
52  }
53  //SmallPrimeField (Integer c);
57  //SmallPrimeField (DenseUnivariateIntegerPolynomial c);
58  //SmallPrimeField (DenseUnivariateRationalPolynomial c);
59  //SmallPrimeField (SparseUnivariatePolynomial<Integer> c);
60  //SmallPrimeField (SparseUnivariatePolynomial<RationalNumber> c);
61  //SmallPrimeField (SparseUnivariatePolynomial<ComplexRationalNumber> c);
62  //template <class Ring>
63  //SmallPrimeField (SparseUnivariatePolynomial<Ring> c);
64 
65  SmallPrimeField* SPFpointer(SmallPrimeField* b) {
66  return b;
67  }
68  SmallPrimeField* SPFpointer(RationalNumber* a) {
69  std::cout << "BPAS error, try to cast pointer to Rational Number to pointer to SmallPrimeField" << std::endl;
70  exit(1);
71  }
72  SmallPrimeField* SPFpointer(BigPrimeField* a) {
73  std::cout << "BPAS error, try to cast pointer to BigPrimeField to pointer to SmallPrimeField" << std::endl;
74  exit(1);
75  }
77  std::cout << "BPAS error, try to cast pointer to GeneralizedFermatPrimeField to pointer to SmallPrimeField" << std::endl;
78  exit(1);
79  }
80 
81  static void setPrime(long int p){
82 
83  prime = p;
84  //using string for mpz_class assignment
85  std::ostringstream ss;
86  ss << p;
87  std::string sp = ss.str();
88 
89  //mpz_t z;
90  //mpz_init(z
91 
92  mpz_class w(sp);
93  // mpz_class w(p);
94  characteristic = w;
95  }
96 
97  long int Prime(){
98  return prime;
99  }
100 
101 
102  void whichprimefield(){
103  cout << "SmallPrimeField" << endl;
104  }
105 
106  long long int number(){
107  return a;
108  }
109 
110 //fing the nth primitive root of unity
111  static SmallPrimeField findPrimitiveRootofUnity(long int n){
112  if ( ((prime - 1) % n != 0)){
113  cout << "ERROR: n does not divide prime - 1." << endl;
114  return -1;
115  }
116  bool flag = false;
117  long int q = (prime - 1) / n;
118  SmallPrimeField p1 (prime - 1);
119  SmallPrimeField c;
120  int i = 0;
121  long int test = q * n / 2;
122  srand (time(NULL));
123  while(i < 20){
124  c = rand();
125  if ((c^test) == p1) {
126  flag = true;
127  return (c^q);
128  }
129  i ++;
130  }
131  if (!flag ){
132  cout << "No primitive root found!"<< endl;
133  return 0;
134  }
135  // If no primitive root found
136  }
137 
139  if (this != &c) {
140  a = c.a;
141  }
142  return *this;
143  }
144  SmallPrimeField& operator= (long long int k) {
145  while (k < 0){
146  k = k + prime;
147  }
148  a = k%prime;
149  return *this;
150  }
151 
152  inline bool isZero() {
153  return (a == 0);
154  }
155  inline void zero() {
156  a = 0;
157  }
158  inline bool isOne() {
159  return (a == 1);
160  }
161  inline void one() {
162  a = 1;
163  }
164 
165  inline bool isNegativeOne() {
166  return (a == (prime - 1));
167  }
168  inline void negativeOne() {
169  a = prime - 1;
170  }
171 
172  inline int isConstant() {
173  if (a >= 0)
174  return 1;
175  else { return -1; }
176  }
177 
178 
179  inline SmallPrimeField gcd (SmallPrimeField c) {
180  SmallPrimeField q (0);
181  SmallPrimeField r (0);
182  SmallPrimeField b (*this);
183  if(b.a < c.a){
184  return c.gcd(b);
185  }
186  while (c.a > 0) {
187  q.a = b.a / c.a;
188  r.a = b.a % c.a;
189  b = c;
190  c = r;
191  }
192  return b;
193  }
194  inline SmallPrimeField gcd (long long int c){
195  SmallPrimeField b (c);
196  SmallPrimeField r (*this);
197  return r.gcd(b);
198  }
199 
200  inline bool operator== (SmallPrimeField& c) {
201  if (a == c.a)
202  return 1;
203  else { return 0; }
204  }
205 
206  inline bool operator== (long long int k){
207  SmallPrimeField r (*this);
208  SmallPrimeField b (k);
209  if (b == r){
210  return 1;
211  }
212  else {
213  return 0;
214  }
215  }
216  inline bool operator!= (SmallPrimeField& c) {
217  if (a == c.a)
218  return 0;
219  else { return 1; }
220  }
221 
222  inline bool operator!= (long long int k){
223  SmallPrimeField r (*this);
224  SmallPrimeField b (k);
225  if (b == r){
226  return 0;
227  }
228  else {
229  return 1;
230  }
231  }
232 
233 
235  SmallPrimeField r (*this);
236  return (r += c);
237  }
238  inline SmallPrimeField operator+ (long long int& c) {
239  SmallPrimeField r (*this);
240  SmallPrimeField b (c);
241  return (r += b);
242  }
243  inline SmallPrimeField operator+ (int& c) {
244  SmallPrimeField r (*this);
245  SmallPrimeField b (c);
246  return (r += b);
247  }
248  inline SmallPrimeField operator+= (long long int& c) {
249  SmallPrimeField r (*this);
250  SmallPrimeField b (c);
251  return (r += b);
252  }
253  inline SmallPrimeField& operator+= (SmallPrimeField c) {
254  a = (a + c.a)%prime;
255  return *this;
256  }
257 
259  SmallPrimeField r (*this);
260  return (r -= c);
261  }
262  inline SmallPrimeField operator- (long long int& c) {
263  SmallPrimeField r (*this);
264  SmallPrimeField b (c);
265  return (r -= b);
266  }
267  inline SmallPrimeField operator-= (long long int& c) {
268  SmallPrimeField r (*this);
269  SmallPrimeField b (c);
270  return (r -= b);
271  }
272  inline SmallPrimeField& operator-= (SmallPrimeField c) {
273  if ((a - c.a)<0){
274  a = prime+(a - c.a);
275  }
276  else{
277  a = a - c.a;
278  }
279  return *this;
280  }
281  inline SmallPrimeField operator- () {
282  a = prime - a;
283  return *this;
284  }
286  SmallPrimeField r (*this);
287  return (r *= c);
288  }
289  inline SmallPrimeField operator* (long long int c) {
290  SmallPrimeField r (*this);
291  return (r *= c);
292  }
293 
295  a = (a * c.a)%prime;
296  return *this;
297  }
298  inline SmallPrimeField operator*= (long long int& c) {
299  a = a * c;
300  SmallPrimeField r (a);
301  return r;
302  }
303  void egcd (long long int x, long long int y, long long int *ao, long long int *bo, long long int *vo, long long int P){
304  long long int t, A, B, C, D, u, v, q;
305 
306  u = y;
307  v = x;
308  A = 1;
309  B = 0;
310  C = 0;
311  D = 1;
312 
313  do
314  {
315  q = u / v;
316  t = u;
317  u = v;
318  v = t - q * v;
319  t = A;
320  A = B;
321  B = t - q * B;
322  t = C;
323  C = D;
324  D = t - q * D;
325  }
326  while (v != 0);
327 
328  *ao = A;
329  *bo = C;
330  *vo = u;
331  }
332 
333  inline SmallPrimeField inverse2(){
334  long long int n, b, v;
335  egcd (a, prime, &n, &b, &v, prime);
336  if (b < 0)
337  b += prime;
338  SmallPrimeField r (b);
339  return r;
340  }
341 
342 
343 //Binary inverse
344 
346  long long int u = a;
347  long long int v = prime;
348  long long int x1 = 1;
349  long long int x2 = 0;
350  SmallPrimeField r(0);
351  while(u!=1 && v!=1){
352  while (u%2 == 0){
353  u = u >> 1;
354  if(x1%2==0)
355  x1 = x1 >> 1;
356  else
357  x1 = (x1 + prime) >> 1;
358  }
359  while (v%2 ==0){
360  v = v >> 1;
361  if(x2%2 == 0)
362  x2 = x2 >> 1;
363  else
364  x2 = (x2 + prime) >> 1;
365  }
366  if(u>=v){
367  u = u - v;
368  x1 = x1 - x2;
369  }
370  else{
371  v = v - u;
372  x2 = x2 - x1;
373  }
374  }
375  if(u == 1)
376  r = x1;
377  else
378  r = x2;
379  return r;
380  }
381 
383  SmallPrimeField r (*this);
384  return (r /= c);
385  }
386  inline SmallPrimeField operator/ (long long int c) {
387  SmallPrimeField r (*this);
388  SmallPrimeField b (c);
389  return (r /= b);
390  }
391  inline SmallPrimeField operator/= (long long int c) {
392  SmallPrimeField r (*this);
393  SmallPrimeField b (c);
394  return (r /= b);
395  }
397  if (c.isZero()) {
398  std::cout << "BPAS: error, dividend is zero from SmallPrimeField."<< std::endl;
399  exit(1);
400  }
401  c=c.inverse();
402  SmallPrimeField r (*this);
403  r*=c;
404  a = r.a;
405  return *this;
406  }
407 
408  inline SmallPrimeField operator^ (long long int e){
409  SmallPrimeField r;
410  if (isZero() || isOne() || e == 1)
411  r = *this;
412  else if (e == 2) {
413  r = *this * *this;
414  }
415  else if (e > 2) {
416  SmallPrimeField x (*this);
417  r.one();
418 
419  while (e != 0) {
420  if (e % 2)
421  r *= x;
422  x = x * x;
423  e >>= 1;
424  }
425  }
426  else if (e == 0) {
427  r.one();
428  }
429  else {
430  r = *this ^ (-e);
431  r=r.inverse();
432  }
433  return r;
434  }
435 
436 
437  inline friend std::ostream& operator<< (std::ostream &out, SmallPrimeField c){
438  out << c.a;
439  return out;
440  }
441 };
442 
443 
444 #else
445 
446 /**
447  * A prime field whose prime is 32 bits or less. Elements of this field
448  * are encoded using montgomery trick.
449  */
450 class SmallPrimeField : public BPASFiniteField<SmallPrimeField> {
451 
452 private:
453 
454  long long int a;
455  static long long int prime;
456  //static long long int R;
457  static unsigned long long int Pp;
458 
459 //R = 2^Wt;
460 // static long long int Rp;
461 
462 //Pp = -p^-1 mod R;
463 // static long long int N;
464 
465 public:
466 
467  static mpz_class characteristic;
468 // static bool isPrimeField;
469 // static bool isSmallPrimeField;
470 // static bool isComplexField;
471 
472  SmallPrimeField ();
473 
474 // using aR mod p to represent the number
475  SmallPrimeField (long long int _a);
476 
477  SmallPrimeField (const SmallPrimeField& c);
478 
479  explicit SmallPrimeField (const Integer& c);
480 
481  explicit SmallPrimeField (const RationalNumber& c);
482 
483  explicit SmallPrimeField (const ComplexRationalNumber& c);
484 
485  explicit SmallPrimeField (const BigPrimeField& c);
486 
487  explicit SmallPrimeField (const GeneralizedFermatPrimeField& c);
488 
490 
492 
494 
496 
498 
499  template <class Ring>
501 
502  SmallPrimeField* SPFpointer(SmallPrimeField* b);
503 
504  SmallPrimeField* SPFpointer(RationalNumber* a);
505 
506  SmallPrimeField* SPFpointer(BigPrimeField* a);
507 
509 
510 // return the actual number
511  long long int number() const;
512 
513  void whichprimefield();
514 
515  static void setPrime(long long int p){
516 
517  prime = p;
518 //using string for mpz_class assignment
519  std::ostringstream ss;
520  ss << p;
521  std::string sp = ss.str();
522 
523 //mpz_t z;
524 //mpz_init(z
525 
526  mpz_class w(sp);
527 // mpz_class w(p);
528  SmallPrimeField::characteristic = w;
529 
530  // long long int p1 = R - prime;
531 
532  // for(int i = 0; i < R; i ++){
533  // if ((i * p1)%R == 1){
534  // Pp = i;
535  // break;
536  // }
537  // }
538 //cout << Pp << endl;
539  long long int t, A, B, C, D, u, v,q;
540  __asm__( //compute the first step of egcd where R = q*prime + v
541  // "movq %%rax,%%rsi\n\t"
542  "xor %%rax,%%rax\n\t"
543  "movq $1,%%rdx\n\t"
544  "divq %2\n\t"
545  // "movq %%rsi,%%rax\n\t"
546  // "movq %%rdx,%0\n\t"
547  // "movq %%rax,%1\n\t"
548  : "=&d" (v),"=&a" (q)
549  : "b"(prime)
550  :"rsi","rdi");
551  A = 1;
552  B = 0;
553  C = 0;
554  D = 1;
555  // q = R/prime;
556  u = prime;
557  //v = R - R mod prime;
558  t = A;
559  A = B;
560  B = t - q * B;
561  t = C;
562  C = D;
563  D = t - q * D;
564 
565  while (v != 0){
566  q = u / v;
567  t = u;
568  u = v;
569  v = t - q * v;
570  t = A;
571  A = B;
572  B = t - q * B;
573  t = C;
574  C = D;
575  D = t - q * D;
576  // printf("%llu,%llu\n",A,C);
577  }
578  //unsigned long long int res;
579  if(C < 0){
580  C = 0-C;
581  Pp = (unsigned long long int)C;
582  }
583  else{
584  __asm__( //compute -C mod R;
585  // "movq %%rax,%%rsi\n\t"
586  "xor %%rax,%%rax\n\t"
587  "movq $1,%%rdx\n\t"
588  "sub %1,%%rax\n\t"
589  "sbb $0,%%rdx\n\t"
590  // "movq %%rsi,%%rax\n\t"
591  // "movq %%rdx,%0\n\t"
592  "movq %%rax,%%rdx\n\t"
593  : "=&d" (Pp)
594  : "b"((unsigned long long int)C)
595  :"rsi","rdi");
596  }
597 
598  }
599 
600  mpz_class getCharacteristic() const override {
601  return characteristic;
602  }
603 
604  long long int Prime();
605 
607 
608  SmallPrimeField& operator= (long long int k);
609 
610  SmallPrimeField findPrimitiveRootOfUnity(long int n) const {
611  return SmallPrimeField::findPrimitiveRootofUnity(n);
612  }
613 
614 
615  static SmallPrimeField findPrimitiveRootofUnity(long long int n){
616  if ( ((prime - 1) % n != 0)){
617  cout << "ERROR: n does not divide prime - 1." << endl;
618  return -1;
619  }
620  bool flag = false;
621  long long int q = (prime - 1) / n;
622  SmallPrimeField p1 (prime - 1);
623  SmallPrimeField c;
624  int i = 0;
625  long long int test = q * n / 2;
626  srand (time(NULL));
627  while(i < 20){
628  c = rand();
629  if ((c^test) == p1) {
630  flag = true;
631  return (c^q);
632  }
633  i ++;
634  }
635  if (!flag ){
636  cout << "No primitive root found!"<< endl;
637  return 0;
638  }
639 // If no primitive root found
640  return 0;
641  }
642 
643  inline bool isZero() const {
644  return (a == 0);
645  }
646 
647  inline void zero() {
648  a = 0;
649  }
650 
651  inline bool isOne() const {
652  SmallPrimeField b(1);
653  return (a == b.a);
654  }
655 
656  inline void one() {
657  SmallPrimeField b(1);
658  a = b.a;// R % prime;
659  }
660 
661  inline bool isNegativeOne() {
662  SmallPrimeField b(-1);
663  return (a == b.a);
664  }
665  inline void negativeOne() {
666  SmallPrimeField b(-1);
667  a = b.a;
668  }
669 
670  inline int isConstant() {
671  if (a >= 0)
672  return 1;
673  else { return -1; }
674  }
675 
676  SmallPrimeField unitCanonical(SmallPrimeField* u = NULL, SmallPrimeField* v = NULL) const ;
677 
678  inline SmallPrimeField operator+ (const SmallPrimeField& c) const {
679  SmallPrimeField r (*this);
680  return (r += c);
681  }
682 
683  inline SmallPrimeField operator+ (const long long int& c) const {
684  SmallPrimeField r (*this);
685  SmallPrimeField b (c);
686  return (r += b);
687  }
688  inline SmallPrimeField operator+ (const int& c) const {
689  SmallPrimeField r (*this);
690  SmallPrimeField b (c);
691  return (r += b);
692  }
693  inline SmallPrimeField operator+= (const long long int& c) {
694  SmallPrimeField r (*this);
695  SmallPrimeField b (c);
696  return (r += b);
697  }
698 
699  inline SmallPrimeField& operator+= (const SmallPrimeField& c) {
700  //a = (a + c.a)%prime;
701  //return *this;
702  a = a + c.a;
703  a -= prime;
704  a += (a >> 63) & prime;
705  return *this;
706  }
707 
708  inline SmallPrimeField operator- (const SmallPrimeField& c) const {
709  SmallPrimeField r (*this);
710  return (r -= c);
711  }
712 
713  inline SmallPrimeField operator- (const long long int& c) const {
714  SmallPrimeField r (*this);
715  SmallPrimeField b (c);
716  return (r -= b);
717  }
718  inline SmallPrimeField operator-= (const long long int& c) {
719  SmallPrimeField r (*this);
720  SmallPrimeField b (c);
721  return (r -= b);
722  }
723 
724  inline SmallPrimeField& operator-= (const SmallPrimeField& c) {
725  // if ((a - c.a)<0){
726  // a = prime+(a - c.a);
727  // }
728  // else{
729  // a = a - c.a;
730  // }
731  // return *this;
732  a = a - c.a;
733  a += (a >> 63) & prime;
734  return *this;
735  }
736  inline SmallPrimeField operator- () const {
737  SmallPrimeField ret (*this);
738  ret.a = prime - a;
739  return ret;
740  }
741 
742  inline SmallPrimeField operator* (const SmallPrimeField& c) const {
743  SmallPrimeField r (*this);
744 
745  r*=c;
746 //cout << r << endl;
747  return r;
748  }
749 
750  inline SmallPrimeField operator* (long long int c) const {
751  SmallPrimeField r (*this);
752  SmallPrimeField b(c);
753  return (r *= b);
754  }
755 // Mont(a,b) = ab/R mod p
756 // Mont(aR,bR) = abR mod p
758  // long long int x = (c.a * a);
759  // long long int w = x*Pp;
760  // long long int y = x + prime*(w&(R-1));
761  // long long int z = y>>32 ;
762  // if(z >= prime){
763  // z = z - prime;
764  // }
765  // a = z;
766  // return *this;
767  long long int _a = a;
768  __asm__(
769  "mulq %2\n\t"
770  "movq %%rax,%%rsi\n\t"
771  "movq %%rdx,%%rdi\n\t"
772  "imulq %3,%%rax\n\t"
773  "mulq %4\n\t"
774  "add %%rsi,%%rax\n\t"
775  "adc %%rdi,%%rdx\n\t"
776  "subq %4,%%rdx\n\t"
777  "mov %%rdx,%%rax\n\t"
778  "sar $63,%%rax\n\t"
779  "andq %4,%%rax\n\t"
780  "addq %%rax,%%rdx\n\t"
781  : "=&d" (_a)
782  : "a"(_a),"rm"(c.a),"b"((unsigned long long int)Pp),"c"(prime)
783  :"rsi","rdi");
784  a = _a;
785  return *this;
786  }
787 
788 //partial Montgomery inversion
789  long long int* pinverse();
790 // long long int MontMul(long long int a, long long int b);
791 
792 //return bc/R mod p;
793  static long long int Mont(long long int b, long long int c);
794  static long long int getRsquare();
795 //Montgomery inversion
796  inline SmallPrimeField inverse() const {
797  SmallPrimeField r(*this);
798  // SmallPrimeField rsquare(R);
799  long long int rsquare = getRsquare();
800  // long long int* result = r.pinverse();
801  // if(result[0] < 32){
802  // result[0] = Mont(result[0], rsquare.a);
803  // result[1] += 32;
804  // }
805  // result[0] = Mont(result[0], rsquare.a);
806  // result[0] = Mont(result[0], pow(2,64-result[1]));
807  // r.a = result[0];
808  // return r;
809 
810  if (r.a == 0)
811  {
812  printf("Not invertible!\n");
813  return 0;
814  }
815  // long long int rsquare = getRsquare(prime);
816  long long int* result = r.pinverse();
817  if(result[1] < 64){
818  result[0] = Mont(result[0], rsquare);
819  result[1] += 64;
820  }
821  result[0] = Mont(result[0], rsquare);
822  long long int tmp;
823  if(result[1] != 64){
824  r.a = 1L << (128 - result[1]);
825  result[0] = Mont(result[0], r.a);
826  }
827  r.a = result[0];
828  free(result);
829  return r;
830  }
831 
832 //Binary inverse
833  SmallPrimeField inverse2();
834 
835  inline SmallPrimeField operator^ (long long int e) const {
836  SmallPrimeField r;
837  if (isZero() || isOne() || e == 1)
838  r = *this;
839  else if (e == 2) {
840  r = *this * *this;
841  }
842  else if (e > 2) {
843  SmallPrimeField x (*this);
844  r.one();
845 
846  while (e != 0) {
847  if (e % 2)
848  r *= x;
849  x = x * x;
850  e >>= 1;
851  }
852  }
853  else if (e == 0) {
854  r.one();
855  }
856  else {
857  r = *this ^ (-e);
858  r=r.inverse();
859  }
860  return r;
861 
862  }
863 
864  inline SmallPrimeField& operator^= (long long int e) {
865  *this = *this ^ e;
866  return *this;
867  }
868 
869  inline bool operator== (const SmallPrimeField& c) const {
870  if ((*this).number() == c.number())
871  return 1;
872  else { return 0; }
873  }
874 
875  inline bool operator== (long long int k) const {
876  SmallPrimeField r (*this);
877  SmallPrimeField b (k);
878  if (b.number() == r.number()){
879  return 1;
880  }
881  else {
882  return 0;
883  }
884  }
885 
886  inline bool operator!= (const SmallPrimeField& c) const {
887  if ((*this).number() == c.number())
888  return 0;
889  else { return 1; }
890  }
891 
892 
893  inline bool operator!= (long long int k) const {
894  SmallPrimeField r (*this);
895  SmallPrimeField b (k);
896  if (b.number() == r.number()){
897  return 0;
898  }
899  else {
900  return 1;
901  }
902  }
903 
905  return ExpressionTree(new ExprTreeNode(this->number()));
906  }
907 
908  inline SmallPrimeField operator/ (const SmallPrimeField& c) const {
909  SmallPrimeField r (*this);
910  r /= c;
911  return r;
912  }
913 
914  inline SmallPrimeField operator/ (long long int c) const {
915  SmallPrimeField r (*this);
916  SmallPrimeField b(c);
917  return (r /= b);
918  }
919 
921  if (c.isZero()) {
922  std::cout << "BPAS: error, dividend is zero from SmallPrimeField."<< std::endl;
923  exit(1);
924  }
925  //cout << "division" << endl;
926  SmallPrimeField inv = c.inverse();
927  SmallPrimeField r (*this);
928  r *= inv;
929  (*this).a = r.a;
930  return *this;
931  }
932 
933  inline SmallPrimeField operator% (const SmallPrimeField& c) const {
934  return *this % c;
935  }
936 
938  (*this).a = (*this).a % c.a;
939  return *this;
940  }
941 
942  inline SmallPrimeField gcd (const SmallPrimeField& other) const {
943  SmallPrimeField q (0);
944  SmallPrimeField r (0);
945  SmallPrimeField b (*this);
946  SmallPrimeField c (other);
947  if(b.a < c.a){
948  return c.gcd(b);
949  }
950  while (c.a > 0) {
951  q.a = b.a / c.a;
952  r.a = b.a % c.a;
953  b = c;
954  c = r;
955  }
956  return b;
957  }
958 
959  /**
960  * Compute squarefree factorization of *this
961  */
963  std::vector<SmallPrimeField> ret;
964  ret.push_back(*this);
965  return ret;
966  }
967 
968  /**
969  * Get the euclidean size of *this.
970  */
971  inline Integer euclideanSize() const {
972  return Integer(1);
973  }
974 
975  /**
976  * Perform the eucldiean division of *this and b. Returns the
977  * remainder. If q is not NULL, then returns the quotient in q.
978  */
980 
981  /**
982  * Perform the extended euclidean division on *this and b.
983  * Returns the GCD. If s and t are not NULL, returns the bezout coefficients in them.
984  */
986 
987  /**
988  * Get the quotient of *this and b.
989  */
990  SmallPrimeField quotient(const SmallPrimeField& b) const;
991 
992  /**
993  * Get the remainder of *this and b.
994  */
995  SmallPrimeField remainder(const SmallPrimeField& b) const;
996 
997 };
998 
999 #endif //montgomery switch
1000 
1001 #endif //include guard
Factors< SmallPrimeField > squareFree() const
Compute squarefree factorization of *this.
Definition: SmallPrimeField.hpp:962
SmallPrimeField & operator=(const SmallPrimeField &c)
Copy assignment.
A sparsely represented univariate polynomial over an arbitrary ring.
Definition: BigPrimeField.hpp:21
SmallPrimeField operator*(const SmallPrimeField &c) const
Multiplication.
Definition: SmallPrimeField.hpp:742
bool isOne() const
Determine if *this ring element is one, that is the multiplication identity.
Definition: SmallPrimeField.hpp:651
SmallPrimeField operator/(const SmallPrimeField &c) const
Exact division.
Definition: SmallPrimeField.hpp:908
SmallPrimeField inverse() const
Get the inverse of *this.
Definition: SmallPrimeField.hpp:796
SmallPrimeField & operator*=(const SmallPrimeField &c)
Multiplication assignment.
Definition: SmallPrimeField.hpp:757
void zero()
Make *this ring element zero.
Definition: SmallPrimeField.hpp:647
An arbitrary-precision complex rational number.
Definition: ComplexRationalNumber.hpp:23
SmallPrimeField operator^(long long int e) const
Exponentiation.
Definition: SmallPrimeField.hpp:835
SmallPrimeField euclideanDivision(const SmallPrimeField &b, SmallPrimeField *q=NULL) const
Perform the eucldiean division of *this and b.
SmallPrimeField operator-() const
Negation.
Definition: SmallPrimeField.hpp:736
An ExpressionTree encompasses various forms of data that can be expressed generically as a binary tre...
Definition: ExpressionTree.hpp:17
A finite field whose prime should be a generalized fermat number.
Definition: GeneralizedFermatPrimeField.hpp:36
friend std::ostream & operator<<(std::ostream &ostream, const SmallPrimeField &d)
Output operator.
Definition: BPASRing.hpp:155
ExpressionTree convertToExpressionTree() const
Convert this to an expression tree.
Definition: SmallPrimeField.hpp:904
void one()
Make *this ring element one.
Definition: SmallPrimeField.hpp:656
A prime field whose prime is 32 bits or less.
Definition: SmallPrimeField.hpp:450
A univariate polynomial with Integer coefficients using a dense representation.
Definition: uzpolynomial.h:14
SmallPrimeField gcd(const SmallPrimeField &other) const
Get GCD of *this and other.
Definition: SmallPrimeField.hpp:942
A univariate polynomial with RationalNumber coefficients represented densely.
Definition: urpolynomial.h:16
A prime field whose prime can be arbitrarily large.
Definition: BigPrimeField.hpp:27
bool isZero() const
Determine if *this ring element is zero, that is the additive identity.
Definition: SmallPrimeField.hpp:643
A simple data structure for encapsulating a collection of Factor elements.
Definition: Factors.hpp:95
bool operator!=(const SmallPrimeField &c) const
Inequality test,.
Definition: SmallPrimeField.hpp:886
SmallPrimeField operator%(const SmallPrimeField &c) const
Get the remainder of *this and b;.
Definition: SmallPrimeField.hpp:933
An arbitrary-precision Integer.
Definition: Integer.hpp:22
SmallPrimeField & operator^=(long long int e)
Exponentiation assignment.
Definition: SmallPrimeField.hpp:864
mpz_class getCharacteristic() const override
The characteristic of this ring class.
Definition: SmallPrimeField.hpp:600
An arbitrary-precision rational number.
Definition: RationalNumber.hpp:24
SmallPrimeField unitCanonical(SmallPrimeField *u=NULL, SmallPrimeField *v=NULL) const
Obtain the unit normal (a.k.a canonical associate) of an element.
SmallPrimeField remainder(const SmallPrimeField &b) const
Get the remainder of *this and b.
An abstract class defining the interface of a field.
Definition: BPASField.hpp:11
SmallPrimeField & operator/=(const SmallPrimeField &c)
Exact division assignment.
Definition: SmallPrimeField.hpp:920
SmallPrimeField quotient(const SmallPrimeField &b) const
Get the quotient of *this and b.
ExprTreeNode is a single node in the bianry tree of an ExpressionTree.
Definition: ExprTreeNode.hpp:76
An abstract class defining the interface of a prime field.
Definition: BPASFiniteField.hpp:12
Integer euclideanSize() const
Get the euclidean size of *this.
Definition: SmallPrimeField.hpp:971
SmallPrimeField operator+(const SmallPrimeField &c) const
Addition.
Definition: SmallPrimeField.hpp:678
bool operator==(const SmallPrimeField &c) const
Equality test,.
Definition: SmallPrimeField.hpp:869
SmallPrimeField extendedEuclidean(const SmallPrimeField &b, SmallPrimeField *s=NULL, SmallPrimeField *t=NULL) const
Perform the extended euclidean division on *this and b.
SmallPrimeField & operator%=(const SmallPrimeField &c)
Assign *this to be the remainder of *this and b.
Definition: SmallPrimeField.hpp:937