My Project  debian-1:4.1.2-p1+ds-2
AEp.cc
Go to the documentation of this file.
1 #include "misc/auxiliary.h"
2 
3 #ifdef SINGULAR_4_2
4 #include "AEp.h"
5 
6 #include <stdio.h>
7 #include <math.h>
8 
9 
10 using namespace std;
11 
12 //Konstruktoren
13 
14 p_poly::p_poly()
15 {
16  //coef = reinterpret_cast<mpz_t*> (omAlloc(100*sizeof(mpz_t)));
17  deg=-1;
18  mod=2;
19  mpz_init_set_ui(coef[0],0);
20 }
21 
22 
23 
24 p_poly::p_poly(int n,int p, mpz_t *a)
25 {
26 
27  deg=n;
28  mod=p;
29 
30  for ( int i=0;i<=n;i++)
31  {
32  mpz_mod_ui(a[i],a[i],mod);
33  mpz_init_set(coef[i], a[i]);
34  }
35 
36 }
37 
38 /*
39 //Destruktor
40 
41 p_poly::~p_poly()
42 {
43  delete[] coef;
44 }
45 
46 */
47 
48 //Reduktion modulo p
49 
50 void p_poly::p_poly_reduce(p_poly f,int p)
51 {
52  if (f.is_zero()==0)
53  {
54 
55  for (int i=f.deg;i>=0;i--)
56  {
57  mpz_mod_ui(coef[i],f.coef[i],p);
58  }
59  }
60  //Hier nötige Grad Korrektur
61  int k=deg;
62  while(mpz_sgn(coef[k])==0 && k>=0)
63  {deg--;k--;}
64 }
65 
66 
67 // Arithmetik
68 
69 
70 //Additionen
71 
72 //Standard - Addition
73 void p_poly::p_poly_add(const p_poly a, const p_poly b)
74 {
75  if (a.deg >=b.deg)
76  {
77 
78  deg=a.deg;
79  mod=a.mod;
80 
81  for ( int i=0;i<=b.deg;i++)
82  {
83  mpz_add(coef[i],a.coef[i],b.coef[i]);
84  }
85 
86  for ( int i=b.deg+1;i<=a.deg;i++)
87  {
88  mpz_init_set(coef[i],a.coef[i]);
89  }
90 
91  //Hier nötige Grad Korrektur
92  int k=deg;
93  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
94  {deg--;k--;}
95 
96 
97  }
98 
99  else {p_poly_add(b,a); }
100 
101 }
102 
103 //Überschreibende Addition
104 
105 void p_poly::p_poly_add_to(const p_poly g)
106 {
107  this->p_poly_add(*this,g);
108 }
109 
110 //Addition einer Konstanten
111 void p_poly::p_poly_add_const(p_poly f,const mpz_t a)
112 {
113  if (f.is_zero()==1 && mpz_divisible_ui_p(a,f.mod)==0)
114  p_poly_set(a,f.mod);
115 
116  else if (mpz_divisible_ui_p(a,f.mod)==0)
117  {
118  p_poly_set(f);
119  mpz_add(coef[0],coef[0],a);
120  /*/Hier nötige Grad Korrektur
121  int k=deg;
122  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
123  {deg--;k--;}
124  */
125  }
126 
127 }
128 
129 
130 //To Variante Addition einer Konstanten
131 
132 void p_poly::p_poly_add_const_to(const mpz_t a)
133 {
134  this->p_poly_add_const(*this,a);
135 }
136 
137 //Monom Addition
138 void p_poly::p_poly_add_mon(const p_poly f, mpz_t a, int i)
139 {
140  p_poly_set(f);
141  if (i<=deg && is_zero()==0)
142  {
143  mpz_add(coef[i],coef[i],a);
144  }
145  else if (is_zero()==1 && mpz_divisible_ui_p(a,f.mod)==0)
146  {
147  deg=i;
148  for(int j=0;j<=i;j++)
149  {
150  mpz_init_set_ui(coef[j],0);
151  }
152  mpz_t temp;
153  mpz_init_set_ui(temp,0);
154  mpz_mod_ui(temp,a,mod);
155  mpz_add(coef[i],coef[i],temp);
156 
157  }
158  else if(i>deg && mpz_divisible_ui_p(a,f.mod)==0)
159  {
160  deg=i;
161  for(int j=i;j>deg;j--)
162  {
163  mpz_init_set_ui(coef[j],0);
164  }
165  mpz_t temp;
166  mpz_init_set_ui(temp,0);
167  mpz_mod_ui(temp,a,mod);
168  mpz_add(coef[i],coef[i],temp);
169 
170  }
171  /*/Hier nötige Grad Korrektur
172  int k=deg;
173  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
174  {deg--;k--;}
175  */
176 
177 }
178 
179 //To Variante Monomaddition
180 void p_poly::p_poly_add_mon_to(mpz_t a, int i)
181 {
182 
183  if (i<=deg && is_zero()==0)
184  {
185  mpz_add(coef[i],coef[i],a);
186  }
187  else if (is_zero()==1 && mpz_divisible_ui_p(a,mod)==0)
188  {
189  deg=i;
190  for(int j=0;j<=i;j++)
191  {
192  mpz_init_set_ui(coef[j],0);
193  }
194  mpz_t temp;
195  mpz_init_set_ui(temp,0);
196  mpz_mod_ui(temp,a,mod);
197  mpz_add(coef[i],coef[i],temp);
198 
199  }
200  else if(i>deg && mpz_divisible_ui_p(a,mod)==0)
201  {
202  deg=i;
203  for(int j=i;j>deg;j--)
204  {
205  mpz_init_set_ui(coef[j],0);
206  }
207  mpz_t temp;
208  mpz_init_set_ui(temp,0);
209  mpz_mod_ui(temp,a,mod);
210  mpz_add(coef[i],coef[i],temp);
211 
212  }
213  /*Hier nötige Grad Korrektur
214  int k=deg;
215  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
216  {deg--;k--;} */
217 }
218 
219 //Subtraktionen
220 
221 void p_poly::p_poly_sub(const p_poly a, const p_poly b)
222 {
223  if (a.deg >=b.deg)
224  {
225 
226  deg=a.deg;
227  mod=a.mod;
228 
229  for ( int i=0;i<=b.deg;i++)
230  {
231  mpz_sub(coef[i],a.coef[i],b.coef[i]);
232  }
233 
234  for ( int i=b.deg+1;i<=a.deg;i++)
235  {
236  mpz_init_set(coef[i],a.coef[i]);
237  }
238 
239  //Hier nötige Grad Korrektur
240  int k=deg;
241  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
242  {deg--;k--;}
243 
244  }
245 
246  else {p_poly_sub(b,a);p_poly_neg(); }
247 
248 }
249 
250 
251 //Überschreibende Subtraktion
252 
253 void p_poly::p_poly_sub_to(const p_poly b)
254 {
255  this->p_poly_sub(*this,b);
256 }
257 
258 //Subtraktion einer Konstanten
259 void p_poly::p_poly_sub_const(p_poly f,const mpz_t a)
260 {
261  if (f.is_zero()==1)
262  {
263  p_poly_set(a,f.mod);
264  p_poly_neg();
265  }
266  else
267  {
268  p_poly_set(f);
269  mpz_sub(coef[0],coef[0],a);
270  }
271  /*/*Hier nötige Grad Korrektur
272  int k=deg;
273  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
274  {deg--;k--;} */
275 
276 }
277 
278 
279 //To Variante Subtraktion einer Konstanten
280 
281 void p_poly::p_poly_sub_const_to(const mpz_t a)
282 {
283  this->p_poly_sub_const(*this,a);
284 }
285 
286 
287 //Monom Subtraktion
288 void p_poly::p_poly_sub_mon(const p_poly f , mpz_t a, int i)
289 {
290  mpz_t temp;
291  mpz_neg(temp,a);
292  p_poly_add_mon(f,temp,i);
293 }
294 
295 //To Variante Monomaddition
296 void p_poly::p_poly_sub_mon_to(mpz_t a, int i)
297 {
298  mpz_t temp;
299  mpz_neg(temp,a);
300  p_poly_add_mon_to(temp,i);
301 }
302 
303 
304 //Multiplikationen
305 
306 //Multiplikation mit Monom
307 void p_poly::p_poly_mon_mult( p_poly f, int n)
308 {
309  if (f.is_zero()==1)
310  {p_poly_set_zero();}
311  else
312  {
313  deg=f.deg+n;
314  mod=f.mod;
315  for (int i=deg;i>=n;i--)
316  {
317  mpz_init_set(coef[i],f.coef[i-n]);
318  }
319  for (int i=n-1;i>=0;i--)
320  {
321  mpz_init_set_ui(coef[i],0);
322  }
323 
324  /*/Hier nötige Grad Korrektur
325  int k=deg;
326  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
327  {deg--;k--;} */
328  }
329 }
330 
331 void p_poly::p_poly_mon_mult_to(const int n)
332 {
333  this->p_poly_mon_mult(*this,n);
334 }
335 
336 
337 //Multiplikation mit Skalar
338 
339 void p_poly::p_poly_scalar_mult(const p_poly g, const mpz_t n)
340 {
341  if (mpz_divisible_ui_p(n,g.mod)!=0)
342  p_poly_set_zero();
343  else
344  {
345  deg=g.deg;
346  mod=g.mod;
347 
348  mpz_t temp;
349  mpz_init_set_ui(temp,0);
350  for(int i=0;i<=deg;i++)
351  {
352  mpz_mul(temp,n,g.coef[i]);
353  mpz_init_set(coef[i],temp);
354  }
355  }
356 }
357 
358 
359 
360 void p_poly::p_poly_scalar_mult(const mpz_t n, const p_poly g)
361 {
362  if (mpz_divisible_ui_p(n,g.mod)!=0)
363  p_poly_set_zero();
364  else
365  {
366  deg=g.deg;
367  mod=g.mod;
368 
369 
370  mpz_t temp;
371  mpz_init_set_ui(temp,0);
372  for(int i=0;i<=deg;i++)
373  {
374  mpz_mul(temp,n,g.coef[i]);
375  mpz_init_set(coef[i],temp);
376  }
377  /*/Hier nötige Grad Korrektur
378  int k=deg;
379  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
380  {deg--;k--;} */
381  }
382 }
383 
384 
385 void p_poly::p_poly_scalar_mult_to(const mpz_t n)
386 {
387  this->p_poly_scalar_mult(*this,n);
388 }
389 
390 
391 
392 // Negation
393 
394 void p_poly::p_poly_neg()
395 {
396  for (int i=0;i<=deg;i++)
397  {
398  mpz_neg(coef[i],coef[i]);
399  }
400 }
401 
402 // Naive Multiplikation
403 void p_poly::p_poly_mult_n(p_poly a,p_poly b)
404 {
405  //Reduktion mod p
406  a.p_poly_reduce(a,a.mod);
407  b.p_poly_reduce(b,b.mod);
408 
409  if (a.is_zero()==1 || b.is_zero()==1)
410  p_poly_set_zero();
411  else
412  {
413  mpz_t temp;
414  mpz_init_set_ui(temp,0);
415  deg = a.deg + b.deg;
416 
417  // Kopien atemp und btemp von a bzw. b, mit Nullen aufgefüllt
418  p_poly atemp, btemp;
419  atemp.p_poly_set(a);
420  btemp.p_poly_set(b);
421  for(int i=a.deg+1;i<=deg;i++)
422  {
423  mpz_init_set_ui(atemp.coef[i],0);
424  }
425  for(int i=b.deg+1;i<=deg;i++)
426  {
427  mpz_init_set_ui(btemp.coef[i],0);
428  }
429  atemp.deg = deg;
430  btemp.deg = deg;
431 
432  // Multiplikationsalgorithmus
433  for (int k=0; k<=deg; k++)
434  {
435  mpz_init_set_ui(coef[k],0); // k-ter Koeffizient zunächst 0
436  for (int i=0; i<=k; i++) // dann schrittweise Summe der a[i]*b[k-i]/
437  {
438  mpz_mul(temp,atemp.coef[i],btemp.coef[k-i]);
439  mpz_add(coef[k],coef[k],temp);
440  }
441  }
442 
443  /*/Hier nötige Grad Korrektur
444  int k=deg;
445  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
446  {deg--;k--;} */
447 
448 
449  }
450 }
451 
452 
453 //Überschreibende Multiplikation
454 
455 void p_poly::p_poly_mult_n_to(const p_poly g)
456 {
457  this->p_poly_mult_n(*this,g);
458 
459 }
460 
461 // Karatsuba-Multiplikation (Weimerskirch/Paar Alg. 1), ACHTUNG VORLÄUFIGE VERSION, macht noch Fehler beim Grad und ist unelegant !!!
462 void p_poly::p_poly_mult_ka( p_poly A, p_poly B)
463 {
464 
465  if (A.is_zero()==1 || B.is_zero()==1)
466  {
467  p_poly_set_zero();
468  }
469  else
470  {
471  //Reduktion mod p
472  A.p_poly_reduce(A,A.mod);
473  B.p_poly_reduce(B,B.mod);
474 
475  // Größeren Grad feststellen
476  int n;
477  if(A.deg>=B.deg){n=A.deg+1;}
478  else{n=B.deg+1;}
479  // n auf nächste 2er-Potenz setzen (VORLÄUFIG!)
480  n = static_cast<int>(ceil(log(n)/log(2)));
481  n = static_cast<int>(pow(2,n));
482 
483  if (n==1)
484  {
485  mpz_t AB;
486  mpz_mul(AB,A.coef[0],B.coef[0]);
487  p_poly_set(AB,A.mod);
488  this->p_poly_reduce(*this,A.mod);
489  }
490  else
491  {
492  // p_polynome A und B aufspalten
493  p_poly Au, Al, Bu, Bl;
494  Au.p_poly_mon_div(A,n/2);
495  Al.p_poly_mon_div_rem(A,n/2);
496  Bu.p_poly_mon_div(B,n/2);
497  Bl.p_poly_mon_div_rem(B,n/2);
498  p_poly Alu,Blu;
499  Alu.p_poly_add(Al,Au);
500  Blu.p_poly_add(Bl,Bu);
501 
502  // Teile rekursiv multiplizieren
503  p_poly D0, D1, D01;
504  D0.p_poly_mult_ka(Al,Bl);
505  D1.p_poly_mult_ka(Au,Bu);
506  D01.p_poly_mult_ka(Alu,Blu);
507 
508  // Ergebnis zusammensetzen
509  p_poly temp;
510  D01.p_poly_sub_to(D0);
511  D01.p_poly_sub_to(D1);
512  D01.p_poly_mon_mult_to(n/2);
513  D1.p_poly_mon_mult_to(n);
514  D1.p_poly_add_to(D01);
515  D1.p_poly_add_to(D0);
516  p_poly_set(D1);
517 
518  }
519 
520  //Hier nötige Grad Korrektur
521  int k=deg;
522  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
523  {deg--;k--;}
524  }
525 }
526 
527 
528 
529 //Skalare Divisionen
530 
531 void p_poly::p_poly_scalar_div( const p_poly g, const mpz_t n)
532 {
533 
534  if (mpz_divisible_ui_p(n,g.mod)==0) // überprüft invertierbarkeit
535  {
536  deg=g.deg;
537  mod=g.mod;
538 
539 
540  mpz_t temp;
541  mpz_t p;
542  mpz_init_set_ui(temp,0);
543  mpz_init_set_ui(p,mod);
544  for(int i=0;i<=deg;i++)
545  {
546  mpz_invert(temp,n,p);
547  mpz_mul(temp,g.coef[i],temp);
548  mpz_init_set(coef[i],temp);
549  }
550 
551  /*/Hier nötige Grad Korrektur
552  int k=deg;
553  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
554  {deg--;k--;} */
555  }
556 }
557 
558 
559 void p_poly::p_poly_scalar_div_to(const mpz_t n)
560 {
561  this->p_poly_scalar_div(*this,n);
562 }
563 
564 // Division durch Monom - Quotient
565 void p_poly::p_poly_mon_div(const p_poly f, const int n)
566 {
567  if (f.deg<n)
568  {
569  p_poly_set_zero();
570  }
571  else
572  {
573  deg=f.deg-n;
574  mod=f.mod;
575 
576  for (int i=0;i<=f.deg-n;i++)
577  {
578  mpz_init_set(coef[i],f.coef[n+i]);
579  }
580  }
581 }
582 
583 // Division durch Monom - Rest
584 void p_poly::p_poly_mon_div_rem(const p_poly f, const int n)
585 {
586  if (f.deg<n)
587  {
588  p_poly_set(f);
589  }
590  else
591  {
592  deg=n-1;
593  mod=f.mod;
594 
595 
596  for (int i=0;i<=n-1;i++)
597  {
598  mpz_init_set(coef[i],f.coef[i]);
599  }
600  }
601 }
602 
603 
604 
605 
606 //Euklidische Division nach Cohen Algo 3.1.1 (degA muss größer gleich deg B sein)!!
607 
608 void p_poly::p_poly_div_rem( p_poly A, p_poly B)
609 {
610 
611  if (B.is_zero()==0)
612  {
613  //Reduktion mod p
614  A.p_poly_reduce(A,A.mod);
615  B.p_poly_reduce(B,B.mod);
616 
617  p_poly R;
618  p_poly temp;
619  R.p_poly_set(A);
620  mpz_t a;
621  mpz_t u;
622  mpz_t p;
623  mpz_init_set_ui(p,A.mod);
624  mpz_init_set_ui(a,0);
625  mpz_init_set_ui(u,0);
626  int i;
627 
628  mpz_invert(u,B.coef[B.deg],p); // Inverse von lc(B)
629 
630 
631  while (R.deg>=B.deg)
632  {
633 
634  mpz_mul(a,R.coef[R.deg],u);
635  i=R.deg-B.deg;
636 
637  temp.p_poly_mon_mult(B,i);
638  temp.p_poly_scalar_mult_to(a);
639 
640  R.p_poly_sub_to(temp);
641 
642  }
643  p_poly_set(R);
644 
645  /*/Hier nötige Grad Korrektur
646  int k=deg;
647  while(mpz_divisible_ui_p(coef[k],mod)!=0 && k>=0)
648  {deg--;k--;}*/
649  }
650 }
651 //To Variante von Algo 3.1.1 im Cohen
652 
653 void p_poly::p_poly_div_rem_to(const p_poly B)
654 {
655  this->p_poly_div_rem(*this,B);
656 
657 
658 }
659 
660 
661 
662 //Exakte Division nach Cohen 3.1.1
663 void p_poly::p_poly_div(p_poly &Q, p_poly &R, p_poly A, p_poly B)
664 {
665  if (B.is_zero()==0)
666  {
667  //Reduktion mod p
668  A.p_poly_reduce(A,A.mod);
669  B.p_poly_reduce(B,B.mod);
670 
671  //Initialisierungen
672  p_poly temp;
673  R.p_poly_set(A);
674  Q.p_poly_set_zero();
675  Q.mod=A.mod;
676 
677  mpz_t a;
678  mpz_t b;
679  mpz_t p;
680  mpz_init_set_ui(p,A.mod);
681  mpz_init_set_ui(a,0);
682  mpz_init_set_ui(b,0);
683  int i;
684  mpz_invert(b,B.coef[B.deg],p);
685 
686  //Algorithmus TO DO: evtl hier mit auch den Rest ausgeben
687  while (R.deg>=B.deg)
688  {
689 
690  mpz_mul(a,R.coef[R.deg],b);
691  i=R.deg-B.deg;
692 
693  temp.p_poly_mon_mult(B,i);
694  temp.p_poly_scalar_mult_to(a);
695 
696  R.p_poly_sub_to(temp);
697 
698  Q.p_poly_add_mon_to(a,i);
699 
700  R.p_poly_reduce(R,R.mod);
701  Q.p_poly_reduce(Q,Q.mod);
702  }
703 
704  /*/Hier nötige Grad Korrektur Q
705  int k=Q.deg;
706  while(mpz_divisible_ui_p(Q.coef[k],Q.mod)!=0 && k>=0)
707  {Q.deg--;k--;}*/
708 
709  /*/Hier nötige Grad Korrektur R
710  k=R.deg;
711  while(mpz_divisible_ui_p(R.coef[k],R.mod)!=0 && k>=0)
712  {R.deg--;k--;} */
713  }
714 }
715 
716 
717 //To Varainte der exakten Division
718 
719 void p_poly::p_poly_div_to(p_poly &Q,p_poly &R, p_poly B)
720 {
721  this->p_poly_div(Q ,R,*this,B);
722 }
723 
724 
725 // Kombinationen
726 
727 // a := a*b + c
728 void p_poly::p_poly_multadd_to(const p_poly b, const p_poly c)
729 {
730  p_poly_mult_n_to(b);
731  p_poly_add_to(c);
732 }
733 
734 //a=a*b-c
735 void p_poly::p_poly_multsub_to(const p_poly b, const p_poly c)
736 {
737  p_poly_mult_n_to(b);
738  p_poly_sub_to(c);
739 }
740 
741 
742 
743 /*
744 // a := (a+b)* c
745 void p_poly::poly_addmult_to(const p_poly b, const p_poly c)
746 {
747  p_poly a(deg,coef);
748  a.poly_add_to(b);
749  a.poly_mult_n_to(c);
750  poly_set(a);
751 }
752 */
753 
754 
755 
756 //Sonstiges
757 void p_poly::p_poly_horner(mpz_t erg, const mpz_t u)
758 {
759  if (is_zero()==0)
760  {
761  mpz_init_set(erg,coef[deg]);
762  for (int i=deg;i>=1;i--)
763  {
764  mpz_mul(erg,erg,u);
765  mpz_add(erg,erg,coef[i-1]);
766  }
767 
768  //Reduktion
769  mpz_mod_ui(erg,erg,mod);
770  }
771  else
772  {
773  mpz_set_ui(erg,0);
774  }
775 }
776 
777 // p_polynom in p_polynom einsetzen(Horner-Schema) KRITISCHE EINGABE x^2, x^2 ....
778 
779 void p_poly::p_poly_horner_p_poly( p_poly A, p_poly B)
780 {
781  //Reduktion mod p
782  A.p_poly_reduce(A,A.mod);
783  B.p_poly_reduce(B,B.mod);
784 
785  p_poly_set(A.coef[A.deg],A.mod);
786  for (int i=A.deg;i>=1;i--)
787  {
788  p_poly_mult_n_to(B);
789  p_poly_add_const_to(A.coef[i-1]);
790  }
791  /*/Hier nötige Grad Korrektur
792  int i=deg;
793  while(mpz_divisible_ui_p(coef[i],mod)!=0 && i>=0)
794  {deg--;i--;} */
795 }
796 
797 
798 
799 //Hilfsfunktionen
800 
801 
802 //setze p_polynom auf p_polynom b
803 void p_poly::p_poly_set(const p_poly b)
804 {
805  deg=b.deg;
806  mod=b.mod;
807 
808 
809  for(int i=0;i<=deg;i++)
810  {
811  mpz_init_set(coef[i],b.coef[i]); // Hier wird init set dringendst benötigt
812  }
813 
814 }
815 
816 // setze p_polynom auf konstantes p_polynom b
817 void p_poly::p_poly_set(const mpz_t b,int p)
818 {
819  deg=0;
820  mod=p;
821 
822  if (mpz_divisible_ui_p(b,mod)!=0)
823  p_poly_set_zero();
824  else
825  {
826  mpz_t temp;
827  mpz_init_set(temp,b);
828  mpz_mod_ui(temp,temp,p);
829  mpz_init_set(coef[0],b);
830  }
831 }
832 
833 
834 //setze p_polynom auf Nullpolynom
835 void p_poly::p_poly_set_zero()
836 {
837  deg = -1;
838 }
839 
840 
841 //Vergleiche ob 2 p_polynome gleich return 1 falls ja sont 0
842 
843 int p_poly::is_equal(const p_poly g) const
844 {
845  if (deg!=g.deg)
846  return 0;
847  else
848 
849  for (int i=deg;i>=0; i--)
850  {
851  if (mpz_cmp(coef[i],g.coef[i])!=0)
852  {return 0;}
853  }
854  return 1;
855 }
856 
857 //Überprüft ob das p_polynom 0 ist
858 
859 int p_poly::is_zero() const
860 {
861  if (deg<0)
862  return 1;
863  else
864  return 0;
865 
866 }
867 
868 int p_poly::is_one() const
869 {
870  if (deg==0)
871  {
872  if (mpz_cmp_ui(coef[0],1)==0) { return 1; }
873  else { return 0; }
874  }
875  else { return 0; }
876 }
877 
878 int p_poly::is_monic() const
879 {
880  if (mpz_cmpabs_ui(coef[deg],1)==0)
881  return 1;
882  else
883  return 0;
884 }
885 
886 // klassischer GGT nach Cohen 3.2.1
887 
888 void p_poly::p_poly_gcd( p_poly A, p_poly B)
889 {
890 
891  //Reduktion mod p
892  A.p_poly_reduce(A,A.mod);
893  B.p_poly_reduce(B,B.mod);
894 
895  if (A.deg<B.deg)
896  p_poly_gcd(B,A);
897  else if (B.is_zero()==1)
898  p_poly_set(A);
899  else
900  {
901  p_poly App;
902  p_poly Bpp;
903  p_poly R;
904  App.p_poly_set(A);
905  Bpp.p_poly_set(B);
906 
907  while (Bpp.is_zero()==0)
908  {
909  R.p_poly_div_rem(App,Bpp);
910  App.p_poly_set(Bpp);
911  Bpp.p_poly_set(R);
912  }
913  p_poly_set(App);
914  }
915 
916 }
917 
918 //Nach nach Fieker 2.12 Symbolisches Rechnen (2012)
919 // gibt g=s*A+t*B aus
920 void p_poly::p_poly_extgcd(p_poly &s,p_poly &t,p_poly &g, p_poly A, p_poly B)
921 {
922 
923  //Reduktion mod p
924  A.p_poly_reduce(A,A.mod);
925  B.p_poly_reduce(B,B.mod);
926 
927 
928  if (A.deg<B.deg)
929  p_poly_extgcd(t,s,g,B,A);
930  else if (B.is_zero()==1)
931  {
932  g.p_poly_set(A);
933  t.p_poly_set_zero();
934 
935  mpz_t temp;
936  mpz_init_set_ui(temp,1);
937  s.p_poly_set(temp,A.mod);
938  }
939 
940  else
941  {
942  mpz_t temp;
943  mpz_init_set_ui(temp,1);
944 
945  p_poly R1;
946  R1.p_poly_set(A);
947  p_poly R2;
948  R2.p_poly_set(B);
949  p_poly R3;
950  R3.mod=A.mod;
951 
952 
953  p_poly S1;
954  S1.p_poly_set(temp,A.mod);
955  p_poly S2;
956  S2.p_poly_set_zero();
957  S2.mod=A.mod;
958  p_poly S3;
959  S3.mod=A.mod;
960 
961  p_poly T1;
962  T1.p_poly_set_zero();
963  T1.mod=A.mod;
964  p_poly T2;
965  T2.p_poly_set(temp,A.mod);
966  p_poly T3;
967  T3.mod=A.mod;
968 
969  p_poly Q;
970 
971  while (R2.is_zero()!=1)
972  {
973  p_poly_div(Q,R3,R1,R2);
974 
975  S3.p_poly_mult_n(Q,S2);
976  S3.p_poly_neg();
977  S3.p_poly_add_to(S1);
978 
979  T3.p_poly_mult_n(Q,T2);
980  T3.p_poly_neg();
981  T3.p_poly_add_to(T1);
982 
983  R1.p_poly_set(R2);
984  R2.p_poly_set(R3);
985 
986  S1.p_poly_set(S2);
987  S2.p_poly_set(S3);
988 
989  T1.p_poly_set(T2);
990  T2.p_poly_set(T3);
991  }
992  t.p_poly_set(T1);
993  s.p_poly_set(S1);
994  g.p_poly_set(R1);
995  }
996 }
997 
998 
999 //Ein & Ausgabe
1000 
1001 //Eingabe
1002 
1003 void p_poly::p_poly_insert()
1004 {
1005 #if 0
1006  cout << "Bitte geben Sie ein p_polynom ein! Zunächst den Grad: " << endl;
1007  cin >> deg;
1008  cout << "Jetzt den modul: " << endl;
1009  cin >> mod;
1010 
1011  for ( int i=0; i<=deg;i++)
1012  {
1013  mpz_init_set_ui(coef[i],0);
1014  printf("Geben Sie nun f[%i] ein:",i);
1015  mpz_inp_str(coef[i],stdin, 10);
1016  }
1017  //Reduktion
1018  this->p_poly_reduce(*this,mod);
1019 #endif
1020 }
1021 
1022 
1023 //Ausgabe
1024 void p_poly::p_poly_print()
1025 {
1026 #if 0
1027 
1028  //Reduktion
1029  this->p_poly_reduce(*this,mod);
1030 
1031 
1032  if (is_zero()==1)
1033  cout << "0" << "\n" <<endl;
1034  else
1035  {
1036  for (int i=deg;i>=1;i--)
1037  {
1038  mpz_out_str(stdout,10, coef[i]);
1039  printf("X%i+",i);
1040  }
1041  mpz_out_str(stdout,10, coef[0]);
1042  printf("\n");
1043  }
1044 #endif
1045 }
1046 
1047 #endif
Rational pow(const Rational &a, int e)
Definition: GMPrat.cc:411
All the auxiliary stuff.
CF_NO_INLINE CanonicalForm mod(const CanonicalForm &, const CanonicalForm &)
Definition: cf_inline.cc:564
int i
Definition: cfEzgcd.cc:125
int k
Definition: cfEzgcd.cc:92
int p
Definition: cfModGcd.cc:4019
g
Definition: cfModGcd.cc:4031
CanonicalForm b
Definition: cfModGcd.cc:4044
FILE * f
Definition: checklibs.c:9
CanonicalForm & mod(const CanonicalForm &)
void T1(ideal h)
Definition: cohomo.cc:2838
void T2(ideal h)
Definition: cohomo.cc:3097
const CanonicalForm int s
Definition: facAbsFact.cc:55
b *CanonicalForm B
Definition: facBivar.cc:52
int j
Definition: facHensel.cc:105
STATIC_VAR jList * Q
Definition: janet.cc:30
gmp_float log(const gmp_float &a)
Definition: mpr_complex.cc:343
const signed long ceil(const ampf< Precision > &x)
Definition: amp.h:788
#define R
Definition: sirandom.c:27
#define A
Definition: sirandom.c:24