source: git/kernel/Poly.h @ fd283f

spielwiese
Last change on this file since fd283f was fd283f, checked in by Hans Schönemann <hannes@…>, 19 years ago
*bricken: invalid use of forward declarated Poly fixed git-svn-id: file:///usr/local/Singular/svn/trunk@8588 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 14.4 KB
Line 
1//$Id: Poly.h,v 1.23 2005-09-02 12:46:36 Singular Exp $
2
3
4
5#ifndef POLYCPP_HEADER
6#define POLYCPP_HEADER
7#include "mod2.h"
8
9#include "numbers.h"
10#include "Number.h"
11#include "febase.h"
12#include "polys.h"
13#include "ring.h"
14
15
16#include <boost/shared_ptr.hpp>
17
18#include <vector>
19
20#define BOOST_DISABLE_THREADS
21
22
23//PolyImpl is a 08/15 poly wrapper
24//Poly wraps around PolyImpl with reference counting using boost
25
26class PolyImpl{
27  friend class PolyBase<POLY_VARIANT_RING,Poly>;
28  friend class PolyBase<POLY_VARIANT_MODUL,Vector>;
29  //friend class PolyBase<POLY_VARIANT_MODUL>;
30  friend class Poly;
31  friend class Vector;
32  //friend class Number;
33 protected:
34  poly getInternalReference() const{
35    return p;
36  }
37 public:
38  ring getRing() const{
39    return r;
40  }
41  friend PolyImpl operator+(const PolyImpl& p1, const PolyImpl& n2);
42  friend PolyImpl operator-(const PolyImpl& p1, const PolyImpl& n2);
43  friend PolyImpl operator/(const PolyImpl& p1, const PolyImpl& n2);
44  friend PolyImpl operator*(const PolyImpl& p1, const PolyImpl& n2);
45  friend bool operator==(const PolyImpl& p1, const PolyImpl& n2);
46  friend PolyImpl operator+(const PolyImpl& p1, int n2);
47  friend PolyImpl operator-(const PolyImpl& p1, int n2);
48  friend PolyImpl operator/(const PolyImpl& p1, int n2);
49  friend PolyImpl operator*(const PolyImpl& p1, int n2);
50  friend bool operator==(const PolyImpl& p1, int n2);
51  Number leadCoef(){
52    return Number(p->coef,r);
53  }
54  PolyImpl& operator=(const PolyImpl& p2){
55    //durch Reihenfolge Selbstzuweisungen berücksichtigt
56    if (this==&p2) return *this;
57    poly pc=p_Copy(p2.p,p2.r);
58    p_Delete(&p,r);
59    r=p2.r;
60    p=pc;
61    return *this;
62  }
63  PolyImpl operator-(){
64    PolyImpl t(*this);
65    t.p=p_Copy(p,r);
66    t.p=p_Neg(t.p,r);
67    return t;
68  }
69  PolyImpl& operator+=(const PolyImpl & p2){
70    if (r!=p2.r){
71      Werror("not the same ring");
72      return *this;
73    }
74    if (this==&p2){
75      number two=n_Init(2,r);
76      p_Mult_nn(p,two,r);
77      return *this;
78    }
79    p=p_Add_q(p,p_Copy(p2.p,p2.r),r);
80   
81    return *this;
82  }
83  PolyImpl& operator*=(const PolyImpl & p2){
84    if (r!=p2.r){
85      Werror("not the same ring");
86      return *this;
87    }
88    if (this==&p2){
89      poly pc=p_Copy(p,r);
90      p=p_Mult_q(p,p2.p,r);
91      return *this;
92    }
93    p=p_Mult_q(p,p_Copy(p2.p,p2.r),r);
94    return *this;
95  }
96  PolyImpl& operator*=(const Number & n){
97    if (r!=n.r){
98      Werror("not the same ring");
99      return *this;
100    }
101   
102    p=p_Mult_nn(p,n.n,r);
103    return *this;
104  }
105  PolyImpl& operator-=(const PolyImpl & p2){
106    if (r!=p2.r){
107      Werror("not the same ring");
108      return *this;
109    }
110    if (this==&p2){
111      p_Delete(&p,r);
112      p=NULL;
113      return *this;
114    }
115
116    poly pc=p_Copy(p2.p,p2.r);
117    pc=p_Neg(pc,r);
118    p=p_Add_q(p,pc,r);
119
120   
121    return *this;
122  }
123  //Div not available for rings other than currRing
124/*   PolyImpl& operator/=(const PolyImpl & p2){ */
125/*     if (r!=p2.r){ */
126/*       Werror("not the same ring"); */
127/*       return *this; */
128/*     } */
129/*     if (this==&p2){ */
130/*       poly one=p_ISet(1,r); */
131/*       p_Delete(&p,r); */
132/*       p=one; */
133/*       return *this; */
134/*     } */
135/*     number nv=n_Div(n,p2.n,r); */
136/*     n_Delete(&n,r); */
137/*     n=nv; */
138/*     return *this; */
139/*   } */
140
141
142
143
144
145
146
147
148
149
150  PolyImpl& operator=(int n){
151 
152    p_Delete(&p,r);
153    p=p_ISet(n,r);
154    return *this;
155 
156  }
157 
158
159  PolyImpl(){
160    r=currRing;
161    p=NULL;
162  }
163  PolyImpl(const PolyImpl & p){
164    r=p.r;
165    this->p=p_Copy(p.p,r);
166  }
167  PolyImpl(poly p, ring r){
168    this->p=p_Copy(p,r);
169    this->r=r;
170  }
171  PolyImpl(poly p, ring r,int){
172    this->p=p;
173    this->r=r;
174  }
175  PolyImpl(int n, ring r){
176    this->p=p_ISet(n,r);
177    this->r=r;
178  }
179  PolyImpl(const Number & n){
180   
181    r=n.r;
182    this->p=p_NSet(n_Copy(n.n,r),r);
183   
184  }
185  explicit PolyImpl(int n){
186    r=currRing;
187    this->p=p_ISet(n,r);
188  }
189  void print(){
190    p_Write(p,r,r);
191  }
192
193  virtual ~PolyImpl(){
194    p_Delete(&p,r);
195  }
196
197 protected:
198  poly p;
199  ring r;
200
201};
202
203inline PolyImpl operator+(const PolyImpl &p1, const PolyImpl& p2){
204  PolyImpl erg(p1);
205  erg+=p2;
206  return erg;
207}
208inline PolyImpl operator*(const PolyImpl &p1, const PolyImpl& p2){
209  PolyImpl erg(p1);
210  erg*=p2;
211  return erg;
212}
213inline PolyImpl operator-(const PolyImpl &p1, const PolyImpl& p2){
214  PolyImpl erg(p1);
215  erg-=p2;
216  return erg;
217}
218/*PolyImpl operator/(const PolyImpl &p1, const PolyImpl& p2){
219  PolyImpl erg(p1);
220  erg/=p2;
221  return erg;
222  }*/
223/*
224bool operator==(const PolyImpl &p1, const PolyImpl& p2){
225  if(p1.r!=p2.r)
226    return false;
227  return n_Equal(p1.n,p2.n,p1.r);
228  }*/
229//Equal Polys not available for oth. rings than currRing
230
231
232inline PolyImpl operator+(const PolyImpl &p1, int p2){
233  PolyImpl erg(p1);
234  erg+=PolyImpl(p2,p1.r);
235  return erg;
236}
237inline PolyImpl operator*(const PolyImpl &p1, int p2){
238  PolyImpl erg(p1);
239  erg*=PolyImpl(p2,p1.r);
240  return erg;
241}
242inline PolyImpl operator-(const PolyImpl &p1, int p2){
243  PolyImpl erg(p1);
244  erg-=PolyImpl(p2,p1.r);
245  return erg;
246}
247/*PolyImpl operator/(const PolyImpl &p1, int p2){
248  PolyImpl erg(p1);
249  erg/=PolyImpl(p2,p1.r);
250  return erg;
251  }*/
252
253/*bool operator==(const PolyImpl &p1, int p2){
254  return n_Equal(p1.n,PolyImpl(p2,p1.r).n,p1.r);
255  }*/
256inline PolyImpl operator+(int p1, const PolyImpl& p2){
257  PolyImpl erg(p2);
258  return erg+=PolyImpl(p1,p2.getRing());
259}
260/*PolyImpl operator-(int p1, const PolyImpl& p2){
261
262  PolyImpl erg(p1,p2.r);
263  return erg-=p2;
264  }*/
265/*PolyImpl operator/(int p1, const PolyImpl& p2){
266  PolyImpl erg(p1,p2.r);
267  return erg/=p2;
268  }*/
269
270inline PolyImpl operator*(int p1, const PolyImpl& p2){
271  PolyImpl erg(p2);
272  return erg*=PolyImpl(p1,p2.getRing());
273}
274/*bool operator==(int p1, const PolyImpl& p2){
275  return p2==PolyImpl(p1,p2.r);
276  }*/
277using namespace boost;
278
279
280template<class T> class ConstTermReference{
281 private:
282  ring r;
283  poly t;
284 public:
285  operator T() const {
286    return T(p_Head(t,r),r);
287  }
288  ConstTermReference(poly p, ring r){
289    this->t=p;
290    this->r=r;
291  }
292  bool isConstant() const{
293    return p_LmIsConstant(t,r);
294  }
295 
296};
297
298template<class T> class PolyInputIterator:
299public std::iterator<std::input_iterator_tag,T,int, shared_ptr<const T>,ConstTermReference<T> >
300{
301
302 
303 private:
304  poly t;
305  ring r;
306  public:
307  bool operator==(const PolyInputIterator& t2){
308    return t2.t==t;
309  }
310  bool operator!=(const PolyInputIterator& t2){
311    return t2.t!=t;
312  }
313  PolyInputIterator(poly p, ring r){
314    t=p;
315    this->r=r;
316  }
317  PolyInputIterator(const PolyInputIterator& it){
318    t=it.t;
319    r=it.r;
320  }
321  PolyInputIterator& operator++(){
322    t=t->next;
323    return *this;
324  }
325  PolyInputIterator operator++(int){
326    PolyInputIterator it(*this);
327    ++(*this);
328    return it;
329  }
330  const ConstTermReference<T> operator*(){
331    return ConstTermReference<T> (t,r);
332  }
333  shared_ptr<const T> operator->(){
334    return shared_ptr<const T>(new T(p_Head(t,r),r,0));
335  }
336
337};
338
339
340//template <poly_variant v> inline PolyBase<v> operator+(const PolyBase<v>& p1, const PolyBase<v>& p2);
341//template <poly_variant v> inline PolyBase<v> operator*(const PolyBase<v>& p1, const PolyBase<v>& p2);
342//template <poly_variant v>inline PolyBase<v> operator*(const PolyBase<v>& p1, const Number& n);
343template<poly_variant variant, class create_type_input> class PolyBase{
344 
345 public:
346  poly as_poly() const{
347    return p_Copy(ptr->p,ptr->r);
348  }
349  typedef create_type_input create_type;
350  typedef PolyInputIterator<create_type> iterator;
351  void copy_on_write(){
352    if (!ptr.unique()){
353      ptr.reset(new PolyImpl(*ptr));
354    }
355  }
356  void print() const {
357    ptr->print();
358  }
359  //* ressource managed by Singular
360  char* c_string(){
361
362    return p_String(ptr->p,ptr->r,ptr->r);
363  }
364
365  PolyBase(ring r=currRing):ptr(new PolyImpl((poly) NULL,r)){
366  }
367
368  PolyBase(const char* c, ring r=currRing):ptr(new PolyImpl((poly)NULL,r)){
369    //p_Read takes no const so do
370    char* cp=(char*) omalloc((strlen(c)+1)*sizeof(char));
371    strcpy(cp,c);
372    p_Read(cp,ptr->p,r);
373    omfree(cp);
374  }
375  PolyBase(const PolyBase&p):ptr(p.ptr){
376  }
377
378
379
380  PolyBase& operator+=(const PolyBase& p2){
381    copy_on_write();
382    *ptr += *p2.ptr;
383   
384    return *this;
385  }
386  PolyBase& operator*=(const Poly & p2);
387  PolyBase& operator*=(Number n){
388    copy_on_write();
389    *ptr *=n;
390   
391    return *this;
392  }
393  /*  void print(){
394     StringSetS("");
395     write();
396     Print(StringAppendS(""));
397     }*/
398  virtual ~PolyBase(){}
399  PolyBase(poly p, ring r):ptr(new PolyImpl(p_Copy(p,r),r)){
400  }
401  PolyBase(poly p, ring r,int):ptr(new PolyImpl(p,r,0)){
402  }
403  /*Poly(Poly& p){
404    ptr=p.ptr;
405    }*/
406
407  PolyInputIterator<create_type> begin(){
408    return PolyInputIterator<create_type>(ptr->p,ptr->r);
409  }
410  PolyInputIterator<create_type> end(){
411    return PolyInputIterator<create_type>(NULL, ptr->r);
412  }
413  ring getRing() const{
414    return ptr->getRing();
415  }
416  int lmTotalDegree() const{
417    return pTotaldegree(ptr->p);
418  }
419  Number leadCoef(){
420    return ptr->leadCoef();
421  }
422  create_type operator-(){
423    create_type erg(*this);
424    erg*=Number(-1,ptr->r);
425    return erg;
426  }
427 protected:
428 
429  PolyBase(PolyImpl& impl):ptr(&impl){
430   
431  }
432  poly getInternalReference(){
433    return ptr->getInternalReference();
434  }
435 protected:
436
437  shared_ptr<PolyImpl> ptr;
438  //friend inline inline Poly operator+(const Poly& p1, const Poly& p2);
439  ///friend inline PolyBase operator*(const Poly& p1, const Poly& p2);
440  //friend inline PolyBase operator*(const Poly& p1, const Number& n);
441  // friend inline inline Poly operator*(const Poly& p1, const Number& n);
442  //  friend inline template PolyBase<poly_variant variant> operator+(const PolyBase<v>& p1, const PolyBase<v>& p2);
443  //friend PolyBase<variant> operator+<>(const PolyBase<variant>& p1, const PolyBase<variant>& p2);
444  //friend PolyBase<variant> operator*<>(const PolyBase<variant>& p1, const PolyBase<variant>& p2);
445  //friend PolyBase<variant> operator*<>(const PolyBase<variant>& p1, const Number& p2);
446};
447
448class Poly: public PolyBase<POLY_VARIANT_RING, Poly>{
449  friend class Vector;
450  friend class PolyBase<POLY_VARIANT_MODUL,Vector>;
451 public:
452
453  Poly(ring r=currRing):PolyBase<POLY_VARIANT_RING, Poly> ((poly)NULL,r,0){
454  }
455  Poly(int n, ring r=currRing):PolyBase<POLY_VARIANT_RING, Poly>(*(new PolyImpl(n,r))){
456   
457  }
458  Poly(const char* c, ring r=currRing):PolyBase<POLY_VARIANT_RING, Poly>(c,r){
459
460  }
461  Poly(const PolyBase<POLY_VARIANT_RING, Poly>& p):PolyBase<POLY_VARIANT_RING, Poly>(p){
462  }
463 
464  Poly(const Number& n):PolyBase<POLY_VARIANT_RING, Poly>(*(new PolyImpl(n))){
465   
466  }
467  Poly(poly p, ring r):PolyBase<POLY_VARIANT_RING, Poly>(p,r){
468   
469  }
470  Poly(poly p, ring r, int):PolyBase<POLY_VARIANT_RING, Poly>(p,r,0){
471  }
472  Poly(std::vector<int> v, ring r=currRing):PolyBase<POLY_VARIANT_RING, Poly>(*(new PolyImpl((poly) NULL,r))){
473    unsigned int i;
474    int s=v.size();
475    poly p=p_ISet(1,r);
476    for(i=0;i<v.size();i++){
477      pSetExp(p,i+1,v[i]);
478    }
479    pSetm(p);
480    ptr.reset(new PolyImpl(p,r));
481  }
482  /*  Poly& operator+=(const Number& n){
483  Poly p2(n);
484  ((PolyBase<POLY_VARIANT_RING, Poly>&) (*this))+=p2;
485  return *this;
486  }*/
487  Poly& operator+=(const Poly& p ){
488
489    ((PolyBase<POLY_VARIANT_RING, Poly>&)*this)+=p;
490    return *this;
491  }
492  Poly& operator+=(const PolyBase<POLY_VARIANT_RING, Poly>& p ){
493
494    ((PolyBase<POLY_VARIANT_RING, Poly>&)*this)+=p;
495    return *this;
496  }
497
498
499};
500class Vector: public PolyBase<POLY_VARIANT_MODUL, Vector>{
501 
502 public:
503
504  Vector(ring r=currRing):PolyBase<POLY_VARIANT_MODUL, Vector> ((poly)NULL,r,0){
505  }
506  Vector(int n, ring r=currRing):PolyBase<POLY_VARIANT_MODUL, Vector>(*(new PolyImpl(n,r))){
507   
508  }
509  Vector(const char* c, ring r=currRing):PolyBase<POLY_VARIANT_MODUL, Vector>(c,r){
510
511  }
512  Vector(const PolyBase<POLY_VARIANT_MODUL, Vector>& p):PolyBase<POLY_VARIANT_MODUL, Vector>(p){
513  }
514 
515
516  Vector(poly p, ring r):PolyBase<POLY_VARIANT_MODUL, Vector>(p,r){
517   
518  }
519  Vector(poly p, ring r, int):PolyBase<POLY_VARIANT_MODUL, Vector>(p,r,0){
520  }
521  Vector(std::vector<int> v, ring r=currRing):PolyBase<POLY_VARIANT_MODUL, Vector>(*(new PolyImpl((poly) NULL,r))){
522    unsigned int i;
523    int s=v.size();
524    poly p=p_ISet(1,r);
525    for(i=0;i<v.size();i++){
526      pSetExp(p,i+1,v[i]);
527    }
528    pSetm(p);
529    ptr.reset(new PolyImpl(p,r));
530  }
531  /*  Poly& operator+=(const Number& n){
532  Poly p2(n);
533  ((PolyBase<POLY_VARIANT_MODUL, Poly>&) (*this))+=p2;
534  return *this;
535  }*/
536  Vector& operator+=(const Vector& p ){
537
538    ((PolyBase<POLY_VARIANT_MODUL, Vector>&)*this)+=p;
539    return *this;
540  }
541  Vector& operator+=(const PolyBase<POLY_VARIANT_MODUL, Vector>& p ){
542
543    ((PolyBase<POLY_VARIANT_MODUL, Vector>&)*this)+=p;
544    return *this;
545  }
546
547};
548
549//typedef Poly PolyBase<POLY_VARIANT_RING>::create_type;
550/*template <poly_variant v, class c> inline PolyBase<v> operator+(const PolyBase<v,c>& p1, const PolyBase<v,c>& p2){
551    PolyImpl* res=new PolyImpl(*p1.ptr);
552    *res+=*p2.ptr;
553    return(PolyBase<v,c>(*res));
554    }*/
555/*template <poly_variant v> inline PolyBase<v> operator*(const PolyBase<v>& p1, const PolyBase<v>& p2){
556    PolyImpl* res=new PolyImpl(*p1.ptr);
557    *res *= *p2.ptr;
558    return(PolyBase<v> (*res));
559    }*/
560/*template <class c> inline PolyBase<POLY_VARIANT_MODUL> operator*(const PolyBase<POLY_VARIANT_MODUL>& p1, const Number& n){
561  PolyBase<POLY_VARIANT_MODUL> erg(p1);
562  erg*=n;
563  return erg;
564  }*/
565inline Poly operator*(const Poly& p, const Poly& p2){
566  Poly erg=p;
567  erg*=p2;
568  return erg;
569}
570inline Vector operator*(const Number& n, const Vector& v){
571  Vector res=v;
572  res*=n;
573  return res;
574}
575
576//assumes monomials commute with numbers
577template <poly_variant variant, class create_type> 
578  inline typename PolyBase<variant,create_type>::create_type
579  operator*
580  (const Number& n, 
581   const PolyBase<variant,create_type>& p)
582{
583  typename PolyBase<variant, create_type>::create_type erg(p);
584  erg*=n;
585  return erg;
586}
587
588inline Vector operator*(const Poly& p, const Vector& v){
589  Vector res(v);
590  res*=p;
591  return res;
592}
593inline Poly operator+(const Poly& p1, const Number& n){
594 Poly f(p1);
595  f+=n;
596  return f;
597  }
598template <poly_variant variant, class create_type> 
599  inline typename PolyBase<variant,create_type>::create_type
600  operator+
601  (const PolyBase<variant,create_type>& b1, 
602   const PolyBase<variant,create_type>& b2)
603{
604  typename PolyBase<variant, create_type>::create_type erg(b1);
605  erg+=b2;
606  return erg;
607}
608inline Vector unitVector(int i,ring r=currRing){
609  poly p=p_ISet(1,r);
610  p_SetComp(p,i,r);
611  return Vector(p,r,0);
612}
613
614template <poly_variant variant, class create_type> 
615   
616inline PolyBase<variant, create_type>& 
617PolyBase<variant, create_type>::operator*=(const Poly & p2){
618    copy_on_write();
619    *ptr *= *p2.ptr;
620   
621    return *this;
622  }
623#endif
Note: See TracBrowser for help on using the repository browser.