ncSAMult.h
Go to the documentation of this file.
1 #ifndef GRING_SA_MULT_H
2 #define GRING_SA_MULT_H
3 /*****************************************
4  * Computer Algebra System SINGULAR *
5  *****************************************/
6 
7 #ifdef HAVE_PLURAL
8 
9 // #include <ncSAMult.h> // for CMultiplier etc classes
10 
11 #include <misc/options.h>
12 #include <polys/monomials/ring.h>
13 #include <polys/nc/summator.h>// for CPolynomialSummator class
14 #include <reporter/reporter.h> // for Print!
17 
18 #include <polys/coeffrings.h>
19 
20 //#include <polys/nc/ncSACache.h> // for CCacheHash etc classes
21 #include <polys/nc/ncSAFormula.h> // for CFormulaPowerMultiplier and enum Enum_ncSAType
22 
23 // //////////////////////////////////////////////////////////////////////// //
24 //
25 
27 
28 
29 template <typename CExponent>
31 {
32  protected:
33  const ring m_basering;
34  const int m_NVars; // N = number of variables
35 
36  public:
37  CMultiplier(ring rBaseRing): m_basering(rBaseRing), m_NVars(rBaseRing->N) {};
38  virtual ~CMultiplier() {};
39 
40  inline ring GetBasering() const { return m_basering; };
41  inline int NVars() const { return m_NVars; }
42 
43 
44  inline poly LM(const poly pTerm, const ring r, int i = 1) const
45  {
46  poly pMonom = p_LmInit(pTerm, r);
47  pSetCoeff0(pMonom, n_Init(i, r));
48  return pMonom;
49  }
50 
51  // Term * Exponent -> Monom * Exponent
52  inline poly MultiplyTE(const poly pTerm, const CExponent expRight)
53  {
54  const ring r = GetBasering();
55  poly pMonom = LM(pTerm, r);
56 
57  poly result = p_Mult_nn(MultiplyME(pMonom, expRight), p_GetCoeff(pTerm, r), r);
58 
59  p_Delete(&pMonom, r);
60 
61  return result;
62  }
63 
64 
65  // Exponent * Term -> Exponent * Monom
66  inline poly MultiplyET(const CExponent expLeft, const poly pTerm)
67  {
68  const ring r = GetBasering();
69  poly pMonom = LM(pTerm, r);
70 
71  poly result = p_Mult_nn(MultiplyEM(expLeft, pMonom), p_GetCoeff(pTerm, r), r);
72 
73  p_Delete(&pMonom, r);
74  return result;
75 
76 
77  }
78 
79 // protected:
80 
81  // Exponent * Exponent
82  virtual poly MultiplyEE(const CExponent expLeft, const CExponent expRight) = 0;
83 
84  // Monom * Exponent
85  virtual poly MultiplyME(const poly pMonom, const CExponent expRight) = 0;
86 
87  // Exponent * Monom
88  virtual poly MultiplyEM(const CExponent expLeft, const poly pMonom) = 0;
89 
90  private: // no copy constuctors!
91  CMultiplier();
92  CMultiplier(const CMultiplier&);
94 
95 };
96 
97 
99 {
100  private:
101  int m_i; // 2-gen subalgebra in these variables...
102  int m_j;
103 
104 // poly m_c_ij;
105 // poly m_d_ij;
106 
107 
108  public:
109  // 1 <= i < j <= NVars()
110  CSpecialPairMultiplier(ring r, int i, int j);
111  virtual ~CSpecialPairMultiplier();
112 
113  inline int GetI() const { return m_i; } // X
114  inline int GetJ() const { return m_j; } // Y > X!
115 
116 // protected:
117  typedef int CExponent;
118 
119  // Exponent * Exponent
120  // Computes: var(j)^{expLeft} * var(i)^{expRight}
121  virtual poly MultiplyEE(const CExponent expLeft, const CExponent expRight) = 0;
122 
123  // Monom * Exponent
124  // pMonom must be of the form: var(j)^{n}
125  virtual poly MultiplyME(const poly pMonom, const CExponent expRight);
126 
127  // Exponent * Monom
128  // pMonom must be of the form: var(i)^{m}
129  virtual poly MultiplyEM(const CExponent expLeft, const poly pMonom);
130 
131 };
132 
133 
134 
135 
136 
137 struct CPower // represents var(iVar)^{iPower}
138 {
139  int Var;
140  int Power;
141 
142  CPower(int i, int n): Var(i), Power(n) {};
143 
144 /*
145  inline poly GetPoly(const ring r) const // TODO: search for GetPoly(r, 1) and remove "1"!
146  {
147  poly p = p_One(r);
148  p_SetExp(p, Var, Power, r);
149  p_Setm(p, r);
150  return p;
151  };
152  inline poly GetPoly(const ring r, int c) const
153  {
154  poly p = p_ISet(c, r);
155  p_SetExp(p, Var, Power, r);
156  p_Setm(p, r);
157  return p;
158  };
159 */
160 
161 };
162 
163 
164 
165 
166 
167 
168 class CPowerMultiplier: public CMultiplier<CPower>
169 {
170  private:
171  CSpecialPairMultiplier** m_specialpairs; // upper triangular submatrix of pairs 1 <= i < j <= N of a N x N matrix.
172 
173 
174  public:
175  CPowerMultiplier(ring r);
176  virtual ~CPowerMultiplier();
177 
178  inline CSpecialPairMultiplier* GetPair(int i, int j) const
179  {
180  assume( m_specialpairs != NULL );
181  assume( i > 0 );
182  assume( i < j );
183  assume( j <= NVars() );
184 
185  return m_specialpairs[( (NVars() * ((i)-1) - ((i) * ((i)-1))/2 + (j)-1) - (i) )];
186  }
187 
188  inline CSpecialPairMultiplier*& GetPair(int i, int j)
189  {
190  assume( m_specialpairs != NULL );
191  assume( i > 0 );
192  assume( i < j );
193  assume( j <= NVars() );
194 
195  return m_specialpairs[( (NVars() * ((i)-1) - ((i) * ((i)-1))/2 + (j)-1) - (i) )];
196  }
197 
198 // protected:
199  typedef CPower CExponent;
200 
201  // Exponent * Exponent
202  // Computes: var(j)^{expLeft} * var(i)^{expRight}
203  virtual poly MultiplyEE(const CExponent expLeft, const CExponent expRight);
204 
205  // Monom * Exponent
206  // pMonom may NOT be of the form: var(j)^{n}!
207  virtual poly MultiplyME(const poly pMonom, const CExponent expRight);
208 
209  // Exponent * Monom
210  // pMonom may NOT be of the form: var(i)^{m}!
211  virtual poly MultiplyEM(const CExponent expLeft, const poly pMonom);
212 
213  // Main templates:
214 
215  // Poly * Exponent
216  inline poly MultiplyPE(const poly pPoly, const CExponent expRight)
217  {
218  bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);
219  CPolynomialSummator sum(GetBasering(), bUsePolynomial);
220 
221  for( poly q = pPoly; q !=NULL; q = pNext(q) )
222  sum += MultiplyTE(q, expRight);
223 
224  return sum;
225  }
226 
227  // Exponent * Poly
228  inline poly MultiplyEP(const CExponent expLeft, const poly pPoly)
229  {
230  bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);
231  CPolynomialSummator sum(GetBasering(), bUsePolynomial);
232 
233  for( poly q = pPoly; q !=NULL; q = pNext(q) )
234  sum += MultiplyET(expLeft, q);
235 
236  return sum;
237  }
238 
239  // Poly * Exponent
240  inline poly MultiplyPEDestroy(poly pPoly, const CExponent expRight)
241  {
242  bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);
243  CPolynomialSummator sum(GetBasering(), bUsePolynomial);
244 
245  for( ; pPoly!=NULL; pPoly = p_LmDeleteAndNext(pPoly, GetBasering()) )
246  sum += MultiplyTE(pPoly, expRight);
247 
248  return sum;
249  }
250 
251  // Exponent * Poly
252  inline poly MultiplyEPDestroy(const CExponent expLeft, poly pPoly)
253  {
254  bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);
255  CPolynomialSummator sum(GetBasering(), bUsePolynomial);
256 
257  for( ; pPoly!=NULL; pPoly = p_LmDeleteAndNext(pPoly, GetBasering()) )
258  sum += MultiplyET(expLeft, pPoly);
259 
260  return sum;
261  }
262 
263 
264 };
265 
266 
267 
268 class CGlobalMultiplier: public CMultiplier<poly>
269 {
270  private:
271 // CGlobalCacheHash* m_cache;
274 
275  public:
277 
278  CGlobalMultiplier(ring r);
279  virtual ~CGlobalMultiplier();
280 
281 
282 // protected:
283  typedef poly CExponent;
284 
285  // the following methods are literally equal!
286 
287  // Exponent * Exponent
288  // TODO: handle components!!!
289  virtual poly MultiplyEE(const CExponent expLeft, const CExponent expRight);
290 
291  // Monom * Exponent
292  virtual poly MultiplyME(const poly pMonom, const CExponent expRight);
293 
294  // Exponent * Monom
295  virtual poly MultiplyEM(const CExponent expLeft, const poly pMonom);
296 
297 
298  // Main templates:
299 
300  // Poly * Exponent
301  inline poly MultiplyPE(const poly pPoly, const CExponent expRight)
302  {
303  assume( pPoly != NULL ); assume( expRight != NULL );
304  const int iComponentMonom = p_GetComp(expRight, GetBasering());
305 
306  bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);
307  CPolynomialSummator sum(GetBasering(), bUsePolynomial);
308 
309 
310  if( iComponentMonom!=0 )
311  {
312  for( poly q = pPoly; q !=NULL; q = pNext(q) )
313  {
314 #ifdef PDEBUG
315  {
316  const int iComponent = p_GetComp(q, GetBasering());
317  assume(iComponent == 0);
318  if( iComponent!=0 )
319  {
320  Werror("MultiplyPE: both sides have non-zero components: %d and %d!\n", iComponent, iComponentMonom);
321  // what should we do further?!?
322  return NULL;
323  }
324 
325  }
326 #endif
327  sum += MultiplyTE(q, expRight); // NO Component!!!
328  }
329  poly t = sum; p_SetCompP(t, iComponentMonom, GetBasering());
330  return t;
331  } // iComponentMonom != 0!
332  else
333  { // iComponentMonom == 0!
334  for( poly q = pPoly; q !=NULL; q = pNext(q) )
335  {
336  const int iComponent = p_GetComp(q, GetBasering());
337 
338 #ifdef PDEBUG
339  if( iComponent!=0 )
340  {
341  Warn("MultiplyPE: Multiplication in the left module from the right by component %d!\n", iComponent);
342  // what should we do further?!?
343  }
344 #endif
345  poly t = MultiplyTE(q, expRight); // NO Component!!!
346  p_SetCompP(t, iComponent, GetBasering());
347  sum += t;
348  }
349  return sum;
350  } // iComponentMonom == 0!
351  }
352 
353  // Exponent * Poly
354  inline poly MultiplyEP(const CExponent expLeft, const poly pPoly)
355  {
356  assume( pPoly != NULL ); assume( expLeft != NULL );
357  const int iComponentMonom = p_GetComp(expLeft, GetBasering());
358 
359  bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);
360  CPolynomialSummator sum(GetBasering(), bUsePolynomial);
361 
362  if( iComponentMonom!=0 )
363  {
364  for( poly q = pPoly; q !=NULL; q = pNext(q) )
365  {
366 #ifdef PDEBUG
367  {
368  const int iComponent = p_GetComp(q, GetBasering());
369  assume(iComponent == 0);
370  if( iComponent!=0 )
371  {
372  Werror("MultiplyEP: both sides have non-zero components: %d and %d!\n", iComponent, iComponentMonom);
373  // what should we do further?!?
374  return NULL;
375  }
376  }
377 #endif
378  sum += MultiplyET(expLeft, q);
379  }
380  poly t = sum; p_SetCompP(t, iComponentMonom, GetBasering());
381  return t;
382  } // iComponentMonom != 0!
383  else
384  { // iComponentMonom == 0!
385  for( poly q = pPoly; q !=NULL; q = pNext(q) )
386  {
387  const int iComponent = p_GetComp(q, GetBasering());
388 
389  poly t = MultiplyET(expLeft, q); // NO Component!!!
390  p_SetCompP(t, iComponent, GetBasering());
391  sum += t;
392  }
393  return sum;
394  } // iComponentMonom == 0!
395  }
396 
397  // Poly * Exponent
398  inline poly MultiplyPEDestroy(poly pPoly, const CExponent expRight)
399  {
400  assume( pPoly != NULL ); assume( expRight != NULL );
401  const int iComponentMonom = p_GetComp(expRight, GetBasering());
402 
403  bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);
404  CPolynomialSummator sum(GetBasering(), bUsePolynomial);
405 
406 
407  if( iComponentMonom!=0 )
408  {
409  for(poly q = pPoly ; q!=NULL; q = p_LmDeleteAndNext(q, GetBasering()) )
410  {
411 #ifdef PDEBUG
412  {
413  const int iComponent = p_GetComp(q, GetBasering());
414  assume(iComponent == 0);
415  if( iComponent!=0 )
416  {
417  Werror("MultiplyPEDestroy: both sides have non-zero components: %d and %d!\n", iComponent, iComponentMonom);
418  // what should we do further?!?
419  return NULL;
420  }
421 
422  }
423 #endif
424  sum += MultiplyTE(q, expRight); // NO Component!!!
425  }
426  poly t = sum; p_SetCompP(t, iComponentMonom, GetBasering());
427  return t;
428  } // iComponentMonom != 0!
429  else
430  { // iComponentMonom == 0!
431  for(poly q = pPoly ; q!=NULL; q = p_LmDeleteAndNext(q, GetBasering()) )
432  {
433  const int iComponent = p_GetComp(q, GetBasering());
434 
435 #ifdef PDEBUG
436  if( iComponent!=0 )
437  {
438  Warn("MultiplyPEDestroy: Multiplication in the left module from the right by component %d!\n", iComponent);
439  // what should we do further?!?
440  }
441 #endif
442  poly t = MultiplyTE(q, expRight); // NO Component!!!
443  p_SetCompP(t, iComponent, GetBasering());
444  sum += t;
445  }
446  return sum;
447  } // iComponentMonom == 0!
448 
449  }
450 
451  // Exponent * Poly
452  inline poly MultiplyEPDestroy(const CExponent expLeft, poly pPoly)
453  {
454 
455  assume( pPoly != NULL ); assume( expLeft != NULL );
456  const int iComponentMonom = p_GetComp(expLeft, GetBasering());
457 
458  bool bUsePolynomial = TEST_OPT_NOT_BUCKETS || (pLength(pPoly) < MIN_LENGTH_BUCKET);
459  CPolynomialSummator sum(GetBasering(), bUsePolynomial);
460 
461  if( iComponentMonom!=0 )
462  {
463  for(poly q = pPoly ; q!=NULL; q = p_LmDeleteAndNext(q, GetBasering()) )
464  {
465 #ifdef PDEBUG
466  {
467  const int iComponent = p_GetComp(q, GetBasering());
468  assume(iComponent == 0);
469  if( iComponent!=0 )
470  {
471  Werror("MultiplyEPDestroy: both sides have non-zero components: %d and %d!\n", iComponent, iComponentMonom);
472  // what should we do further?!?
473  return NULL;
474  }
475  }
476 #endif
477  sum += MultiplyET(expLeft, q);
478  }
479  poly t = sum; p_SetCompP(t, iComponentMonom, GetBasering());
480  return t;
481  } // iComponentMonom != 0!
482  else
483  { // iComponentMonom == 0!
484  for(poly q = pPoly ; q!=NULL; q = p_LmDeleteAndNext(q, GetBasering()) )
485  {
486  const int iComponent = p_GetComp(q, GetBasering());
487 
488  poly t = MultiplyET(expLeft, q); // NO Component!!!
489  p_SetCompP(t, iComponent, GetBasering());
490  sum += t;
491  }
492  return sum;
493  } // iComponentMonom == 0!
494 
495  }
496 
497 
498 
499 
500 };
501 
502 
503 
504 //////////////////////////////////////////////////////////////////////////
506 {
507  public:
508  CCommutativeSpecialPairMultiplier(ring r, int i, int j);
510 
511  // Exponent * Exponent
512  virtual poly MultiplyEE(const int expLeft, const int expRight);
513 };
514 
515 //////////////////////////////////////////////////////////////////////////
517 {
518  public:
519  CAntiCommutativeSpecialPairMultiplier(ring r, int i, int j);
521 
522  // Exponent * Exponent
523  virtual poly MultiplyEE(const int expLeft, const int expRight);
524 };
525 
526 
527 //////////////////////////////////////////////////////////////////////////
529 {
530  private:
531  const number m_q;
532  // TODO: make cache for some 'good' powers!?
533 
534  public:
535  CQuasiCommutativeSpecialPairMultiplier(ring r, int i, int j, number q);
537 
538  // Exponent * Exponent
539  virtual poly MultiplyEE(const int expLeft, const int expRight);
540 };
541 
542 
543 //////////////////////////////////////////////////////////////////////////
545 {
546  private:
547  const number m_g;
548  // TODO: make cache for some 'good' powers!?
549 
550  public:
551  CWeylSpecialPairMultiplier(ring r, int i, int j, number g);
552  virtual ~CWeylSpecialPairMultiplier();
553 
554  // Exponent * Exponent
555  virtual poly MultiplyEE(const int expLeft, const int expRight);
556 };
557 
558 //////////////////////////////////////////////////////////////////////////
560 {
561  private:
562  const int m_k;
563  // TODO: make cache for some 'good' powers!?
564 
565  public:
566  CHWeylSpecialPairMultiplier(ring r, int i, int j, int k);
567  virtual ~CHWeylSpecialPairMultiplier();
568 
569  // Exponent * Exponent
570  virtual poly MultiplyEE(const int expLeft, const int expRight);
571 };
572 
573 
574 //////////////////////////////////////////////////////////////////////////
576 {
577  private:
578  const number m_shiftCoef;
579  const int m_shiftVar;
580  // TODO: make cache for some 'good' powers!?
581 
582  public:
583  CShiftSpecialPairMultiplier(ring r, int i, int j, int s, number c);
584  virtual ~CShiftSpecialPairMultiplier();
585 
586  // Exponent * Exponent
587  virtual poly MultiplyEE(const int expLeft, const int expRight);
588 };
589 
590 
591 
592 // need: enum Enum_ncSAType;
593 
594 //////////////////////////////////////////////////////////////////////////
595 // Using external 'formula' routins
597 {
598  private:
600  public:
601  CExternalSpecialPairMultiplier(ring r, int i, int j, Enum_ncSAType type);
603 
604  // Exponent * Exponent
605  virtual poly MultiplyEE(const int expLeft, const int expRight);
606 };
607 
608 
609 #endif // HAVE_PLURAL :(
610 #endif //
const CanonicalForm int s
Definition: facAbsFact.cc:55
virtual poly MultiplyME(const poly pMonom, const CExponent expRight)=0
int GetJ() const
Definition: ncSAMult.h:114
static poly p_LmDeleteAndNext(poly p, const ring r)
Definition: p_polys.h:720
int NVars() const
Definition: ncSAMult.h:41
virtual ~CMultiplier()
Definition: ncSAMult.h:38
poly MultiplyTE(const poly pTerm, const CExponent expRight)
Definition: ncSAMult.h:52
CSpecialPairMultiplier ** m_specialpairs
Definition: ncSAMult.h:171
#define p_GetComp(p, r)
Definition: monomials.h:72
virtual poly MultiplyEM(const CExponent expLeft, const poly pMonom)=0
static FORCE_INLINE number n_Init(long i, const coeffs r)
a number representing i in the given coeff field/ring r
Definition: coeffs.h:539
poly MultiplyEPDestroy(const CExponent expLeft, poly pPoly)
Definition: ncSAMult.h:452
#define MIN_LENGTH_BUCKET
Definition: p_Mult_q.h:22
const CFormulaPowerMultiplier * m_RingFormulaMultiplier
Definition: ncSAMult.h:273
CPower(int i, int n)
Definition: ncSAMult.h:142
int Var
Definition: ncSAMult.h:139
poly MultiplyPEDestroy(poly pPoly, const CExponent expRight)
Definition: ncSAMult.h:240
virtual poly MultiplyEE(const CExponent expLeft, const CExponent expRight)=0
g
Definition: cfModGcd.cc:4031
int k
Definition: cfEzgcd.cc:93
int Power
Definition: ncSAMult.h:140
CPower CExponent
Definition: ncSAMult.h:199
const int m_NVars
Definition: ncSAMult.h:34
CSpecialPairMultiplier * GetPair(int i, int j) const
Definition: ncSAMult.h:178
CSpecialPairMultiplier *& GetPair(int i, int j)
Definition: ncSAMult.h:188
int GetI() const
Definition: ncSAMult.h:113
static int pLength(poly a)
Definition: p_polys.h:189
CMultiplier< poly > CBaseType
Definition: ncSAMult.h:276
#define TEST_OPT_NOT_BUCKETS
Definition: options.h:100
CPowerMultiplier * m_powers
Definition: ncSAMult.h:272
poly MultiplyPE(const poly pPoly, const CExponent expRight)
Definition: ncSAMult.h:216
const ring r
Definition: syzextra.cc:208
poly MultiplyEP(const CExponent expLeft, const poly pPoly)
Definition: ncSAMult.h:228
poly MultiplyPEDestroy(poly pPoly, const CExponent expRight)
Definition: ncSAMult.h:398
static void p_SetCompP(poly p, int i, ring r)
Definition: p_polys.h:249
const CanonicalForm CFMap CFMap & N
Definition: cfEzgcd.cc:49
int j
Definition: myNF.cc:70
#define assume(x)
Definition: mod2.h:405
CMultiplier & operator=(const CMultiplier &)
poly MultiplyEP(const CExponent expLeft, const poly pPoly)
Definition: ncSAMult.h:354
int i
Definition: cfEzgcd.cc:123
static poly p_Mult_nn(poly p, number n, const ring r)
Definition: p_polys.h:901
ring GetBasering() const
Definition: ncSAMult.h:40
static void p_Delete(poly *p, const ring r)
Definition: p_polys.h:849
const number m_shiftCoef
Definition: ncSAMult.h:578
#define NULL
Definition: omList.c:10
const ring m_basering
Definition: ncSAMult.h:33
#define pNext(p)
Definition: monomials.h:43
static poly p_LmInit(poly p, const ring r)
Definition: p_polys.h:1264
#define pSetCoeff0(p, n)
Definition: monomials.h:67
#define p_GetCoeff(p, r)
Definition: monomials.h:57
poly MultiplyET(const CExponent expLeft, const poly pTerm)
Definition: ncSAMult.h:66
poly LM(const poly pTerm, const ring r, int i=1) const
Definition: ncSAMult.h:44
poly MultiplyEPDestroy(const CExponent expLeft, poly pPoly)
Definition: ncSAMult.h:252
BOOLEAN ncInitSpecialPairMultiplication(ring r)
Definition: ncSAMult.cc:267
poly MultiplyPE(const poly pPoly, const CExponent expRight)
Definition: ncSAMult.h:301
CMultiplier(ring rBaseRing)
Definition: ncSAMult.h:37
CPolynomialSummator: unifies bucket and polynomial summation as the later is brocken in buckets :(...
Definition: summator.h:29
polyrec * poly
Definition: hilb.h:10
int BOOLEAN
Definition: auxiliary.h:131
void Werror(const char *fmt,...)
Definition: reporter.cc:199
return result
Definition: facAbsBiFact.cc:76
Enum_ncSAType
Definition: ncSAFormula.h:17
#define Warn
Definition: emacs.cc:80