source: git/kernel/p_polys.cc @ ab4778

spielwiese
Last change on this file since ab4778 was ab4778, checked in by Hans Schönemann <hannes@…>, 17 years ago
*hannes: optim. git-svn-id: file:///usr/local/Singular/svn/trunk@10175 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 20.0 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/***************************************************************
5 *  File:    p_polys.cc
6 *  Purpose: implementation of currRing independent poly procedures
7 *  Author:  obachman (Olaf Bachmann)
8 *  Created: 8/00
9 *  Version: $Id: p_polys.cc,v 1.8 2007-07-05 16:14:37 Singular Exp $
10 *******************************************************************/
11
12#include "mod2.h"
13#include "structs.h"
14#include "structs.h"
15#include "p_polys.h"
16#include "ring.h"
17#include "int64vec.h"
18#ifndef NDEBUG
19#include "febase.h"
20#endif
21
22/***************************************************************
23 *
24 * Completing what needs to be set for the monomial
25 *
26 ***************************************************************/
27// this is special for the syz stuff
28static int* _Components = NULL;
29static long* _ShiftedComponents = NULL;
30static int _ExternalComponents = 0;
31
32BOOLEAN pSetm_error=0;
33
34void p_Setm_General(poly p, const ring r)
35{
36  p_LmCheckPolyRing(p, r);
37  int pos=0;
38  if (r->typ!=NULL)
39  {
40    loop
41    {
42      long ord=0;
43      sro_ord* o=&(r->typ[pos]);
44      switch(o->ord_typ)
45      {
46        case ro_dp:
47        {
48          int a,e;
49          a=o->data.dp.start;
50          e=o->data.dp.end;
51          for(int i=a;i<=e;i++) ord+=p_GetExp(p,i,r);
52          p->exp[o->data.dp.place]=ord;
53          break;
54        }
55        case ro_wp_neg:
56          ord=POLY_NEGWEIGHT_OFFSET;
57          // no break;
58        case ro_wp:
59        {
60          int a,e;
61          a=o->data.wp.start;
62          e=o->data.wp.end;
63          int *w=o->data.wp.weights;
64#if 1
65          for(int i=a;i<=e;i++) ord+=p_GetExp(p,i,r)*w[i-a];
66#else
67          long ai;
68          int ei,wi;
69          for(int i=a;i<=e;i++)
70          {
71             ei=p_GetExp(p,i,r);
72             wi=w[i-a];
73             ai=ei*wi;
74             if (ai/ei!=wi) pSetm_error=TRUE;
75             ord+=ai;
76             if (ord<ai) pSetm_error=TRUE;
77          }
78#endif
79          p->exp[o->data.wp.place]=ord;
80          break;
81        }
82      case ro_wp64:
83        {
84          int64 ord=0;
85          int a,e;
86          a=o->data.wp64.start;
87          e=o->data.wp64.end;
88          int64 *w=o->data.wp64.weights64;
89          int64 ei,wi,ai;
90          for(int i=a;i<=e;i++) {
91            //Print("exp %d w %d \n",p_GetExp(p,i,r),(int)w[i-a]);
92            //ord+=((int64)p_GetExp(p,i,r))*w[i-a];
93            ei=(int64)p_GetExp(p,i,r);
94            wi=w[i-a];
95            ai=ei*wi;
96            if(ei!=0 && ai/ei!=wi){
97              pSetm_error=TRUE;
98              Print("ai %lld, wi %lld\n",ai,wi);
99            }
100            ord+=ai;
101            if (ord<ai){
102               pSetm_error=TRUE;
103               Print("ai %lld, ord %lld\n",ai,ord);
104            }
105          }
106          int64 mask=(int64)0x7fffffff;
107          long a_0=(long)(ord&mask); //2^31
108          long a_1=(long)(ord >>31 ); /*(ord/(mask+1));*/
109
110          //Print("mask: %x,  ord: %d,  a_0: %d,  a_1: %d\n"
111          //,(int)mask,(int)ord,(int)a_0,(int)a_1);
112                    //Print("mask: %d",mask);
113
114          p->exp[o->data.wp64.place]=a_1;
115          p->exp[o->data.wp64.place+1]=a_0;
116//            if(p_Setm_error) Print("***************************\n
117//                                    ***************************\n
118//                                    **WARNING: overflow error**\n
119//                                    ***************************\n
120//                                    ***************************\n");
121          break;
122        }
123        case ro_cp:
124        {
125          int a,e;
126          a=o->data.cp.start;
127          e=o->data.cp.end;
128          int pl=o->data.cp.place;
129          for(int i=a;i<=e;i++) { p->exp[pl]=p_GetExp(p,i,r); pl++; }
130          break;
131        }
132        case ro_syzcomp:
133        {
134          int c=p_GetComp(p,r);
135          long sc = c;
136          int* Components = (_ExternalComponents ? _Components :
137                             o->data.syzcomp.Components);
138          long* ShiftedComponents = (_ExternalComponents ? _ShiftedComponents:
139                                     o->data.syzcomp.ShiftedComponents);
140          if (ShiftedComponents != NULL)
141          {
142            assume(Components != NULL);
143            assume(c == 0 || Components[c] != 0);
144            sc = ShiftedComponents[Components[c]];
145            assume(c == 0 || sc != 0);
146          }
147          p->exp[o->data.syzcomp.place]=sc;
148          break;
149        }
150        case ro_syz:
151        {
152          int c=p_GetComp(p, r);
153          if (c > o->data.syz.limit)
154            p->exp[o->data.syz.place] = o->data.syz.curr_index;
155          else if (c > 0)
156            p->exp[o->data.syz.place]= o->data.syz.syz_index[c];
157          else
158          {
159            assume(c == 0);
160            p->exp[o->data.syz.place]= 0;
161          }
162          break;
163        }
164        default:
165          dReportError("wrong ord in rSetm:%d\n",o->ord_typ);
166          return;
167      }
168      pos++;
169      if (pos == r->OrdSize) return;
170    }
171  }
172}
173
174void p_Setm_Syz(poly p, ring r, int* Components, long* ShiftedComponents)
175{
176  _Components = Components;
177  _ShiftedComponents = ShiftedComponents;
178  _ExternalComponents = 1;
179  p_Setm_General(p, r);
180  _ExternalComponents = 0;
181}
182
183// dummy for lp, ls, etc
184void p_Setm_Dummy(poly p, const ring r)
185{
186  p_LmCheckPolyRing(p, r);
187}
188
189// for dp, Dp, ds, etc
190void p_Setm_TotalDegree(poly p, const ring r)
191{
192  p_LmCheckPolyRing(p, r);
193  p->exp[r->pOrdIndex] = p_ExpVectorQuerSum(p, r);
194}
195
196// for wp, Wp, ws, etc
197void p_Setm_WFirstTotalDegree(poly p, const ring r)
198{
199  p_LmCheckPolyRing(p, r);
200  p->exp[r->pOrdIndex] = pWFirstTotalDegree(p, r);
201}
202
203p_SetmProc p_GetSetmProc(ring r)
204{
205  // covers lp, rp, ls,
206  if (r->typ == NULL) return p_Setm_Dummy;
207
208  if (r->OrdSize == 1)
209  {
210    if (r->typ[0].ord_typ == ro_dp &&
211        r->typ[0].data.dp.start == 1 &&
212        r->typ[0].data.dp.end == r->N &&
213        r->typ[0].data.dp.place == r->pOrdIndex)
214      return p_Setm_TotalDegree;
215    if (r->typ[0].ord_typ == ro_wp &&
216        r->typ[0].data.wp.start == 1 &&
217        r->typ[0].data.wp.end == r->N &&
218        r->typ[0].data.wp.place == r->pOrdIndex &&
219        r->typ[0].data.wp.weights == r->firstwv)
220      return p_Setm_WFirstTotalDegree;
221  }
222  return p_Setm_General;
223}
224
225
226/* -------------------------------------------------------------------*/
227/* several possibilities for pFDeg: the degree of the head term       */
228/*2
229* compute the degree of the leading monomial of p
230* the ordering is compatible with degree, use a->order
231*/
232inline long _pDeg(poly a, ring r)
233{
234  p_LmCheckPolyRing(a, r);
235  assume(p_GetOrder(a, r) == pWTotaldegree(a, r));
236  return p_GetOrder(a, r);
237}
238
239long pDeg(poly a, ring r)
240{
241  return _pDeg(a, r);
242}
243
244/*2
245* compute the degree of the leading monomial of p
246* with respect to weigths 1
247* (all are 1 so save multiplications or they are of different signs)
248* the ordering is not compatible with degree so do not use p->Order
249*/
250inline long _pTotaldegree(poly p, ring r)
251{
252  p_LmCheckPolyRing(p, r);
253  return (long) p_ExpVectorQuerSum(p, r);
254}
255long pTotaldegree(poly p, ring r)
256{
257  return (long) _pTotaldegree(p, r);
258}
259
260// pWTotalDegree for weighted orderings
261// whose first block covers all variables
262inline long _pWFirstTotalDegree(poly p, ring r)
263{
264  int i;
265  long sum = 0;
266
267  for (i=1; i<= r->firstBlockEnds; i++)
268  {
269    sum += p_GetExp(p, i, r)*r->firstwv[i-1];
270  }
271  return sum;
272}
273
274long pWFirstTotalDegree(poly p, ring r)
275{
276  return (long) _pWFirstTotalDegree(p, r);
277}
278
279/*2
280* compute the degree of the leading monomial of p
281* with respect to weigths from the ordering
282* the ordering is not compatible with degree so do not use p->Order
283*/
284long pWTotaldegree(poly p, ring r)
285{
286  p_LmCheckPolyRing(p, r);
287  int i, k;
288  long j =0;
289
290  // iterate through each block:
291  for (i=0;r->order[i]!=0;i++)
292  {
293    int b0=r->block0[i];
294    int b1=r->block1[i];
295    switch(r->order[i])
296    {
297      case ringorder_M:
298        for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
299        { // in jedem block:
300          j+= p_GetExp(p,k,r)*r->wvhdl[i][k - b0 /*r->block0[i]*/]*r->OrdSgn;
301        }
302        break;
303      case ringorder_wp:
304      case ringorder_ws:
305      case ringorder_Wp:
306      case ringorder_Ws:
307        for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
308        { // in jedem block:
309          j+= p_GetExp(p,k,r)*r->wvhdl[i][k - b0 /*r->block0[i]*/];
310        }
311        break;
312      case ringorder_lp:
313      case ringorder_ls:
314      case ringorder_dp:
315      case ringorder_ds:
316      case ringorder_Dp:
317      case ringorder_Ds:
318      case ringorder_rp:
319        for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
320        {
321          j+= p_GetExp(p,k,r);
322        }
323        break;
324      case ringorder_a64:
325        {
326          int64* w=(int64*)r->wvhdl[i];
327          for (k=0;k<=(b1 /*r->block1[i]*/ - b0 /*r->block0[i]*/);k++)
328          {
329            //there should be added a line which checks if w[k]>2^31
330            j+= p_GetExp(p,k+1, r)*(long)w[k];
331          }
332          //break;
333          return j;
334        }
335      case ringorder_c:
336      case ringorder_C:
337      case ringorder_S:
338      case ringorder_s:
339      case ringorder_aa:
340        break;
341      case ringorder_a:
342        for (k=b0 /*r->block0[i]*/;k<=b1 /*r->block1[i]*/;k++)
343        { // only one line
344          j+= p_GetExp(p,k, r)*r->wvhdl[i][ k- b0 /*r->block0[i]*/];
345        }
346        //break;
347        return j;
348
349#ifndef NDEBUG
350      default:
351        Print("missing order %d in pWTotaldegree\n",r->order[i]);
352        break;
353#endif
354    }
355  }
356  return  j;
357}
358
359int pWeight(int i, const ring r)
360{
361  if ((r->firstwv==NULL) || (i>r->firstBlockEnds))
362  {
363    return 1;
364  }
365  return r->firstwv[i-1];
366}
367
368long pWDegree(poly p, ring r)
369{
370  if (r->firstwv==NULL) return pTotaldegree(p, r);
371  p_LmCheckPolyRing(p, r);
372  int i, k;
373  long j =0;
374
375  for(i=1;i<=r->firstBlockEnds;i++)
376    j+=p_GetExp(p, i, r)*r->firstwv[i-1];
377
378  for (;i<=r->N;i++)
379    j+=p_GetExp(p,i, r)*pWeight(i, r);
380
381  return j;
382}
383
384
385/* ---------------------------------------------------------------------*/
386/* several possibilities for pLDeg: the maximal degree of a monomial in p*/
387/*  compute in l also the pLength of p                                   */
388
389/*2
390* compute the length of a polynomial (in l)
391* and the degree of the monomial with maximal degree: the last one
392*/
393long pLDeg0(poly p,int *l, ring r)
394{
395  p_CheckPolyRing(p, r);
396  Exponent_t k= p_GetComp(p, r);
397  int ll=1;
398
399  if (k > 0)
400  {
401    while ((pNext(p)!=NULL) && (p_GetComp(pNext(p), r)==k))
402    {
403      pIter(p);
404      ll++;
405    }
406  }
407  else
408  {
409     while (pNext(p)!=NULL)
410     {
411       pIter(p);
412       ll++;
413     }
414  }
415  *l=ll;
416  return r->pFDeg(p, r);
417}
418
419/*2
420* compute the length of a polynomial (in l)
421* and the degree of the monomial with maximal degree: the last one
422* but search in all components before syzcomp
423*/
424long pLDeg0c(poly p,int *l, ring r)
425{
426  assume(p!=NULL);
427#ifdef PDEBUG
428  _p_Test(p,r,PDEBUG);
429#endif
430  p_CheckPolyRing(p, r);
431  long o;
432  int ll=1;
433
434  if (! rIsSyzIndexRing(r))
435  {
436    while (pNext(p) != NULL)
437    {
438      pIter(p);
439      ll++;
440    }
441    o = r->pFDeg(p, r);
442  }
443  else
444  {
445    int curr_limit = rGetCurrSyzLimit(r);
446    poly pp = p;
447    while ((p=pNext(p))!=NULL)
448    {
449      if (p_GetComp(p, r)<=curr_limit/*syzComp*/)
450        ll++;
451      else break;
452      pp = p;
453    }
454#ifdef PDEBUG
455    _p_Test(pp,r,PDEBUG);
456#endif
457    o = r->pFDeg(pp, r);
458  }
459  *l=ll;
460  return o;
461}
462
463/*2
464* compute the length of a polynomial (in l)
465* and the degree of the monomial with maximal degree: the first one
466* this works for the polynomial case with degree orderings
467* (both c,dp and dp,c)
468*/
469long pLDegb(poly p,int *l, ring r)
470{
471  p_CheckPolyRing(p, r);
472  Exponent_t k= p_GetComp(p, r);
473  long o = r->pFDeg(p, r);
474  int ll=1;
475
476  if (k != 0)
477  {
478    while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
479    {
480      ll++;
481    }
482  }
483  else
484  {
485    while ((p=pNext(p)) !=NULL)
486    {
487      ll++;
488    }
489  }
490  *l=ll;
491  return o;
492}
493
494/*2
495* compute the length of a polynomial (in l)
496* and the degree of the monomial with maximal degree:
497* this is NOT the last one, we have to look for it
498*/
499long pLDeg1(poly p,int *l, ring r)
500{
501  p_CheckPolyRing(p, r);
502  Exponent_t k= p_GetComp(p, r);
503  int ll=1;
504  long  t,max;
505
506  max=r->pFDeg(p, r);
507  if (k > 0)
508  {
509    while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
510    {
511      t=r->pFDeg(p, r);
512      if (t>max) max=t;
513      ll++;
514    }
515  }
516  else
517  {
518    while ((p=pNext(p))!=NULL)
519    {
520      t=r->pFDeg(p, r);
521      if (t>max) max=t;
522      ll++;
523    }
524  }
525  *l=ll;
526  return max;
527}
528
529/*2
530* compute the length of a polynomial (in l)
531* and the degree of the monomial with maximal degree:
532* this is NOT the last one, we have to look for it
533* in all components
534*/
535long pLDeg1c(poly p,int *l, ring r)
536{
537  p_CheckPolyRing(p, r);
538  int ll=1;
539  long  t,max;
540
541  max=r->pFDeg(p, r);
542  if (rIsSyzIndexRing(r))
543  {
544    long limit = rGetCurrSyzLimit(r);
545    while ((p=pNext(p))!=NULL)
546    {
547      if (p_GetComp(p, r)<=limit)
548      {
549        if ((t=r->pFDeg(p, r))>max) max=t;
550        ll++;
551      }
552      else break;
553    }
554  }
555  else
556  {
557    while ((p=pNext(p))!=NULL)
558    {
559      if ((t=r->pFDeg(p, r))>max) max=t;
560      ll++;
561    }
562  }
563  *l=ll;
564  return max;
565}
566
567// like pLDeg1, only pFDeg == pDeg
568long pLDeg1_Deg(poly p,int *l, ring r)
569{
570  assume(r->pFDeg == pDeg);
571  p_CheckPolyRing(p, r);
572  Exponent_t k= p_GetComp(p, r);
573  int ll=1;
574  long  t,max;
575
576  max=_pDeg(p, r);
577  if (k > 0)
578  {
579    while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
580    {
581      t=_pDeg(p, r);
582      if (t>max) max=t;
583      ll++;
584    }
585  }
586  else
587  {
588    while ((p=pNext(p))!=NULL)
589    {
590      t=_pDeg(p, r);
591      if (t>max) max=t;
592      ll++;
593    }
594  }
595  *l=ll;
596  return max;
597}
598
599long pLDeg1c_Deg(poly p,int *l, ring r)
600{
601  assume(r->pFDeg == pDeg);
602  p_CheckPolyRing(p, r);
603  int ll=1;
604  long  t,max;
605
606  max=_pDeg(p, r);
607  if (rIsSyzIndexRing(r))
608  {
609    long limit = rGetCurrSyzLimit(r);
610    while ((p=pNext(p))!=NULL)
611    {
612      if (p_GetComp(p, r)<=limit)
613      {
614        if ((t=_pDeg(p, r))>max) max=t;
615        ll++;
616      }
617      else break;
618    }
619  }
620  else
621  {
622    while ((p=pNext(p))!=NULL)
623    {
624      if ((t=_pDeg(p, r))>max) max=t;
625      ll++;
626    }
627  }
628  *l=ll;
629  return max;
630}
631
632// like pLDeg1, only pFDeg == pTotoalDegree
633long pLDeg1_Totaldegree(poly p,int *l, ring r)
634{
635  p_CheckPolyRing(p, r);
636  Exponent_t k= p_GetComp(p, r);
637  int ll=1;
638  long  t,max;
639
640  max=_pTotaldegree(p, r);
641  if (k > 0)
642  {
643    while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
644    {
645      t=_pTotaldegree(p, r);
646      if (t>max) max=t;
647      ll++;
648    }
649  }
650  else
651  {
652    while ((p=pNext(p))!=NULL)
653    {
654      t=_pTotaldegree(p, r);
655      if (t>max) max=t;
656      ll++;
657    }
658  }
659  *l=ll;
660  return max;
661}
662
663long pLDeg1c_Totaldegree(poly p,int *l, ring r)
664{
665  p_CheckPolyRing(p, r);
666  int ll=1;
667  long  t,max;
668
669  max=_pTotaldegree(p, r);
670  if (rIsSyzIndexRing(r))
671  {
672    long limit = rGetCurrSyzLimit(r);
673    while ((p=pNext(p))!=NULL)
674    {
675      if (p_GetComp(p, r)<=limit)
676      {
677        if ((t=_pTotaldegree(p, r))>max) max=t;
678        ll++;
679      }
680      else break;
681    }
682  }
683  else
684  {
685    while ((p=pNext(p))!=NULL)
686    {
687      if ((t=_pTotaldegree(p, r))>max) max=t;
688      ll++;
689    }
690  }
691  *l=ll;
692  return max;
693}
694
695// like pLDeg1, only pFDeg == pWFirstTotalDegree
696long pLDeg1_WFirstTotalDegree(poly p,int *l, ring r)
697{
698  p_CheckPolyRing(p, r);
699  Exponent_t k= p_GetComp(p, r);
700  int ll=1;
701  long  t,max;
702
703  max=_pWFirstTotalDegree(p, r);
704  if (k > 0)
705  {
706    while (((p=pNext(p))!=NULL) && (p_GetComp(p, r)==k))
707    {
708      t=_pWFirstTotalDegree(p, r);
709      if (t>max) max=t;
710      ll++;
711    }
712  }
713  else
714  {
715    while ((p=pNext(p))!=NULL)
716    {
717      t=_pWFirstTotalDegree(p, r);
718      if (t>max) max=t;
719      ll++;
720    }
721  }
722  *l=ll;
723  return max;
724}
725
726long pLDeg1c_WFirstTotalDegree(poly p,int *l, ring r)
727{
728  p_CheckPolyRing(p, r);
729  int ll=1;
730  long  t,max;
731
732  max=_pWFirstTotalDegree(p, r);
733  if (rIsSyzIndexRing(r))
734  {
735    long limit = rGetCurrSyzLimit(r);
736    while ((p=pNext(p))!=NULL)
737    {
738      if (p_GetComp(p, r)<=limit)
739      {
740        if ((t=_pTotaldegree(p, r))>max) max=t;
741        ll++;
742      }
743      else break;
744    }
745  }
746  else
747  {
748    while ((p=pNext(p))!=NULL)
749    {
750      if ((t=_pTotaldegree(p, r))>max) max=t;
751      ll++;
752    }
753  }
754  *l=ll;
755  return max;
756}
757
758/***************************************************************
759 *
760 * Maximal Exponent business
761 *
762 ***************************************************************/
763
764static inline unsigned long
765p_GetMaxExpL2(unsigned long l1, unsigned long l2, ring r,
766              unsigned long number_of_exp)
767{
768  const unsigned long bitmask = r->bitmask;
769  unsigned long ml1 = l1 & bitmask;
770  unsigned long ml2 = l2 & bitmask;
771  unsigned long max = (ml1 > ml2 ? ml1 : ml2);
772  unsigned long j = number_of_exp - 1;
773
774  if (j > 0)
775  {
776    unsigned long mask = bitmask << r->BitsPerExp;
777    while (1)
778    {
779      ml1 = l1 & mask;
780      ml2 = l2 & mask;
781      max |= ((ml1 > ml2 ? ml1 : ml2) & mask);
782      j--;
783      if (j == 0) break;
784      mask = mask << r->BitsPerExp;
785    }
786  }
787  return max;
788}
789
790static inline unsigned long
791p_GetMaxExpL2(unsigned long l1, unsigned long l2, ring r)
792{
793  return p_GetMaxExpL2(l1, l2, r, r->ExpPerLong);
794}
795
796poly p_GetMaxExpP(poly p, ring r)
797{
798  p_CheckPolyRing(p, r);
799  if (p == NULL) return p_Init(r);
800  poly max = p_LmInit(p, r);
801  pIter(p);
802  if (p == NULL) return max;
803  int i, offset;
804  unsigned long l_p, l_max;
805  unsigned long divmask = r->divmask;
806
807  do
808  {
809    offset = r->VarL_Offset[0];
810    l_p = p->exp[offset];
811    l_max = max->exp[offset];
812    // do the divisibility trick to find out whether l has an exponent
813    if (l_p > l_max ||
814        (((l_max & divmask) ^ (l_p & divmask)) != ((l_max-l_p) & divmask)))
815      max->exp[offset] = p_GetMaxExpL2(l_max, l_p, r);
816
817    for (i=1; i<r->VarL_Size; i++)
818    {
819      offset = r->VarL_Offset[i];
820      l_p = p->exp[offset];
821      l_max = max->exp[offset];
822      // do the divisibility trick to find out whether l has an exponent
823      if (l_p > l_max ||
824          (((l_max & divmask) ^ (l_p & divmask)) != ((l_max-l_p) & divmask)))
825        max->exp[offset] = p_GetMaxExpL2(l_max, l_p, r);
826    }
827    pIter(p);
828  }
829  while (p != NULL);
830  return max;
831}
832
833unsigned long p_GetMaxExpL(poly p, ring r, unsigned long l_max)
834{
835  unsigned long l_p, divmask = r->divmask;
836  int i;
837
838  while (p != NULL)
839  {
840    l_p = p->exp[r->VarL_Offset[0]];
841    if (l_p > l_max ||
842        (((l_max & divmask) ^ (l_p & divmask)) != ((l_max-l_p) & divmask)))
843      l_max = p_GetMaxExpL2(l_max, l_p, r);
844    for (i=1; i<r->VarL_Size; i++)
845    {
846      l_p = p->exp[r->VarL_Offset[i]];
847      // do the divisibility trick to find out whether l has an exponent
848      if (l_p > l_max ||
849          (((l_max & divmask) ^ (l_p & divmask)) != ((l_max-l_p) & divmask)))
850        l_max = p_GetMaxExpL2(l_max, l_p, r);
851    }
852    pIter(p);
853  }
854  return l_max;
855}
856
857
858
859
860/***************************************************************
861 *
862 * Misc things
863 *
864 ***************************************************************/
865// returns TRUE, if all monoms have the same component
866BOOLEAN p_OneComp(poly p, ring r)
867{
868  if(p!=NULL)
869  {
870    long i = p_GetComp(p, r);
871    while (pNext(p)!=NULL)
872    {
873      pIter(p);
874      if(i != p_GetComp(p, r)) return FALSE;
875    }
876  }
877  return TRUE;
878}
879
880/*2
881*test if a monomial /head term is a pure power
882*/
883int p_IsPurePower(const poly p, const ring r)
884{
885  int i,k=0;
886
887  for (i=r->N;i;i--)
888  {
889    if (p_GetExp(p,i, r)!=0)
890    {
891      if(k!=0) return 0;
892      k=i;
893    }
894  }
895  return k;
896}
897
898/*2
899* returns a polynomial representing the integer i
900*/
901poly p_ISet(int i, ring r)
902{
903  poly rc = NULL;
904  if (i!=0)
905  {
906    rc = p_Init(r);
907    pSetCoeff0(rc,r->cf->nInit(i));
908    if (r->cf->nIsZero(p_GetCoeff(rc,r)))
909      p_DeleteLm(&rc,r);
910  }
911  return rc;
912}
913
914/*2
915* returns a polynomial representing the number n
916* destroys n
917*/
918poly p_NSet(number n, ring r)
919{
920  if (r->cf->nIsZero(n))
921  {
922    r->cf->cfDelete(&n, r);
923    return NULL;
924  }
925  else
926  {
927    poly rc = p_Init(r);
928    pSetCoeff0(rc,n);
929    return rc;
930  }
931}
932
933/***************************************************************
934 *
935 * p_ShallowDelete
936 *
937 ***************************************************************/
938#undef LINKAGE
939#define LINKAGE
940#undef p_Delete
941#define p_Delete p_ShallowDelete
942#undef n_Delete
943#define n_Delete(n, r) ((void)0)
944
945#include "p_Delete__T.cc"
946
Note: See TracBrowser for help on using the repository browser.