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