source: git/kernel/polys.cc @ 225d94

spielwiese
Last change on this file since 225d94 was 225d94, checked in by Hans Schönemann <hannes@…>, 17 years ago
*hannes: V_IDLIFT, idLift redLazy: select shortest poly git-svn-id: file:///usr/local/Singular/svn/trunk@9598 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 18.4 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/* $Id: polys.cc,v 1.13 2006-12-15 17:16:07 Singular Exp $ */
5
6/*
7* ABSTRACT - all basic methods to manipulate polynomials
8*/
9
10/* includes */
11#include <stdio.h>
12#include <string.h>
13#include <ctype.h>
14#include "mod2.h"
15#include "structs.h"
16#include "omalloc.h"
17#include "febase.h"
18#include "numbers.h"
19#include "polys.h"
20#include "ring.h"
21
22#ifdef HAVE_PLURAL
23#include "gring.h"
24#endif
25
26/* ----------- global variables, set by pSetGlobals --------------------- */
27/* computes length and maximal degree of a POLYnomial */
28pLDegProc pLDeg;
29/* computes the degree of the initial term, used for std */
30pFDegProc pFDeg;
31/* the monomial ordering of the head monomials a and b */
32/* returns -1 if a comes before b, 0 if a=b, 1 otherwise */
33
34int pVariables;     // number of variables
35
36/* 1 for polynomial ring, -1 otherwise */
37int     pOrdSgn;
38// it is of type int, not BOOLEAN because it is also in ip
39/* TRUE if the monomial ordering is not compatible with pFDeg */
40BOOLEAN pLexOrder;
41
42/* ----------- global variables, set by procedures from hecke/kstd1 ----- */
43/* the highest monomial below pHEdge */
44poly      ppNoether = NULL;
45
46/* -------------------------------------------------------- */
47/*2
48* change all global variables to fit the description of the new ring
49*/
50
51
52void pSetGlobals(const ring r, BOOLEAN complete)
53{
54  int i;
55  if (ppNoether!=NULL) pDelete(&ppNoether);
56  pVariables = r->N;
57  pOrdSgn = r->OrdSgn;
58  pFDeg=r->pFDeg;
59  pLDeg=r->pLDeg;
60  pLexOrder=r->LexOrder;
61
62  if (complete)
63  {
64    test &= ~ TEST_RINGDEP_OPTS;
65    test |= r->options;
66  }
67}
68
69// resets the pFDeg and pLDeg: if pLDeg is not given, it is
70// set to currRing->pLDegOrig, i.e. to the respective LDegProc which
71// only uses pFDeg (and not pDeg, or pTotalDegree, etc)
72void pSetDegProcs(pFDegProc new_FDeg, pLDegProc new_lDeg)
73{
74  assume(new_FDeg != NULL);
75  pFDeg = new_FDeg;
76  currRing->pFDeg = new_FDeg;
77
78  if (new_lDeg == NULL)
79    new_lDeg = currRing->pLDegOrig;
80
81  pLDeg = new_lDeg;
82  currRing->pLDeg = new_lDeg;
83}
84
85
86// restores pFDeg and pLDeg:
87extern void pRestoreDegProcs(pFDegProc old_FDeg, pLDegProc old_lDeg)
88{
89  assume(old_FDeg != NULL && old_lDeg != NULL);
90  pFDeg = old_FDeg;
91  currRing->pFDeg = old_FDeg;
92  pLDeg = old_lDeg;
93  currRing->pLDeg = old_lDeg;
94}
95
96/*2
97* assumes that the head term of b is a multiple of the head term of a
98* and return the multiplicant *m
99*/
100poly pDivide(poly a, poly b)
101{
102  int i;
103  poly result = pInit();
104
105  for(i=(int)pVariables; i; i--)
106    pSetExp(result,i, pGetExp(a,i)- pGetExp(b,i));
107  pSetComp(result, pGetComp(a) - pGetComp(b));
108  pSetm(result);
109  return result;
110}
111
112/*2
113* divides a by the monomial b, ignores monomials which are not divisible
114* assumes that b is not NULL
115*/
116poly pDivideM(poly a, poly b)
117{
118  if (a==NULL) return NULL;
119  poly result=a;
120  poly prev=NULL;
121  int i;
122  number inv=nInvers(pGetCoeff(b));
123
124  while (a!=NULL)
125  {
126    if (pDivisibleBy(b,a))
127    {
128      for(i=(int)pVariables; i; i--)
129         pSubExp(a,i, pGetExp(b,i));
130      pSubComp(a, pGetComp(b));
131      pSetm(a);
132      prev=a;
133      pIter(a);
134    }
135    else
136    {
137      if (prev==NULL)
138      {
139        pDeleteLm(&result);
140        a=result;
141      }
142      else
143      {
144        pDeleteLm(&pNext(prev));
145        a=pNext(prev);
146      }
147    }
148  }
149  pMult_nn(result,inv);
150  nDelete(&inv);
151  pDelete(&b);
152  return result;
153}
154
155/*2
156* returns the LCM of the head terms of a and b in *m
157*/
158void pLcm(poly a, poly b, poly m)
159{
160  int i;
161  for (i=pVariables; i; i--)
162  {
163    pSetExp(m,i, si_max( pGetExp(a,i), pGetExp(b,i)));
164  }
165  pSetComp(m, si_max(pGetComp(a), pGetComp(b)));
166  /* Don't do a pSetm here, otherwise hres/lres chockes */
167}
168
169/*2
170* convert monomial given as string to poly, e.g. 1x3y5z
171*/
172char * p_Read(char *st, poly &rc, ring r)
173{
174  if (r==NULL) { rc=NULL;return st;}
175  int i,j;
176  rc = p_Init(r);
177  char *s = r->cf->nRead(st,&(rc->coef));
178  if (s==st)
179  /* i.e. it does not start with a coeff: test if it is a ringvar*/
180  {
181    j = r_IsRingVar(s,r);
182    if (j >= 0)
183    {
184      p_IncrExp(rc,1+j,r);
185      while (*s!='\0') s++;
186      goto done;
187    }
188  }
189  while (*s!='\0')
190  {
191    char ss[2];
192    ss[0] = *s++;
193    ss[1] = '\0';
194    j = r_IsRingVar(ss,r);
195    if (j >= 0)
196    {
197      char *s_save=s;
198      s = eati(s,&i);
199      if (((unsigned long)i) >  r->bitmask)
200      {
201        return s_save;
202      }
203      p_AddExp(rc,1+j, (Exponent_t)i, r);
204    }
205    else
206    {
207      s--;
208      return s;
209    }
210  }
211done:
212  if (r->cf->nIsZero(pGetCoeff(rc))) p_DeleteLm(&rc,r);
213  else
214  {
215    p_Setm(rc,r);
216  }
217  return s;
218}
219
220poly pmInit(char *st, BOOLEAN &ok)
221{
222  poly p;
223  char *s=p_Read(st,p,currRing);
224  if (*s!='\0')
225  {
226    if ((s!=st)&&isdigit(st[0]))
227    {
228      errorreported=TRUE;
229    }
230    ok=FALSE;
231    pDelete(&p);
232    return NULL;
233  }
234  ok=!errorreported;
235  return p;
236}
237
238/*2
239*make p homogeneous by multiplying the monomials by powers of x_varnum
240*/
241poly pHomogen (poly p, int varnum)
242{
243  poly q=NULL, qn;
244  int  o,ii;
245  sBucket_pt bp;
246
247  if (p!=NULL)
248  {
249    if ((varnum < 1) || (varnum > pVariables))
250    {
251      return NULL;
252    }
253    o=pWTotaldegree(p);
254    q=pNext(p);
255    while (q != NULL)
256    {
257      ii=pWTotaldegree(q);
258      if (ii>o) o=ii;
259      pIter(q);
260    }
261    q = pCopy(p);
262    bp = sBucketCreate(currRing);
263    while (q != NULL)
264    {
265      ii = o-pWTotaldegree(q);
266      if (ii!=0)
267      {
268        pAddExp(q,varnum, (Exponent_t)ii);
269        pSetm(q);
270      }
271      qn = pNext(q);
272      pNext(q) = NULL;
273      sBucket_Add_p(bp, q, 1);
274      q = qn;
275    }
276    sBucketDestroyAdd(bp, &q, &ii);
277  }
278  return q;
279}
280
281/*2
282*replaces the maximal powers of the leading monomial of p2 in p1 by
283*the same powers of n, utility for dehomogenization
284*/
285poly pDehomogen (poly p1,poly p2,number n)
286{
287  polyset P;
288  int     SizeOfSet=5;
289  int     i;
290  poly    p;
291  number  nn;
292
293  P = (polyset)omAlloc0(5*sizeof(poly));
294  //for (i=0; i<5; i++)
295  //{
296  //  P[i] = NULL;
297  //}
298  pCancelPolyByMonom(p1,p2,&P,&SizeOfSet);
299  p = P[0];
300  //P[0] = NULL ;// for safety, may be removed later
301  for (i=1; i<SizeOfSet; i++)
302  {
303    if (P[i] != NULL)
304    {
305      nPower(n,i,&nn);
306      pMult_nn(P[i],nn);
307      p = pAdd(p,P[i]);
308      //P[i] =NULL; // for safety, may be removed later
309      nDelete(&nn);
310    }
311  }
312  omFreeSize((ADDRESS)P,SizeOfSet*sizeof(poly));
313  return p;
314}
315
316/*4
317*Returns the exponent of the maximal power of the leading monomial of
318*p2 in that of p1
319*/
320static int pGetMaxPower (poly p1,poly p2)
321{
322  int     i,k,res = 32000; /*a very large integer*/
323
324  if (p1 == NULL) return 0;
325  for (i=1; i<=pVariables; i++)
326  {
327    if ( pGetExp(p2,i) != 0)
328    {
329      k =  pGetExp(p1,i) /  pGetExp(p2,i);
330      if (k < res) res = k;
331    }
332  }
333  return res;
334}
335
336/*2
337*Returns as i-th entry of P the coefficient of the (i-1) power of
338*the leading monomial of p2 in p1
339*/
340void pCancelPolyByMonom (poly p1,poly p2,polyset * P,int * SizeOfSet)
341{
342  int   maxPow;
343  poly  p,qp,Coeff;
344
345  if (*P == NULL)
346  {
347    *P = (polyset) omAlloc(5*sizeof(poly));
348    *SizeOfSet = 5;
349  }
350  p = pCopy(p1);
351  while (p != NULL)
352  {
353    qp = p->next;
354    p->next = NULL;
355    maxPow = pGetMaxPower(p,p2);
356    Coeff = pDivByMonom(p,p2);
357    if (maxPow > *SizeOfSet)
358    {
359      pEnlargeSet(P,*SizeOfSet,maxPow+1-*SizeOfSet);
360      *SizeOfSet = maxPow+1;
361    }
362    (*P)[maxPow] = pAdd((*P)[maxPow],Coeff);
363    pDelete(&p);
364    p = qp;
365  }
366}
367
368/*2
369*returns the leading monomial of p1 divided by the maximal power of that
370*of p2
371*/
372poly pDivByMonom (poly p1,poly p2)
373{
374  int     k, i;
375
376  if (p1 == NULL) return NULL;
377  k = pGetMaxPower(p1,p2);
378  if (k == 0)
379    return pHead(p1);
380  else
381  {
382    number n;
383    poly p = pInit();
384
385    p->next = NULL;
386    for (i=1; i<=pVariables; i++)
387    {
388       pSetExp(p,i, pGetExp(p1,i)-k* pGetExp(p2,i));
389    }
390    nPower(p2->coef,k,&n);
391    pSetCoeff0(p,nDiv(p1->coef,n));
392    nDelete(&n);
393    pSetm(p);
394    return p;
395  }
396}
397/*----------utilities for syzygies--------------*/
398poly pTakeOutComp(poly * p, int k)
399{
400  poly q = *p,qq=NULL,result = NULL;
401
402  if (q==NULL) return NULL;
403  if (pGetComp(q)==k)
404  {
405    result = q;
406    do
407    {
408      pSetComp(q,0);
409      pSetmComp(q);
410      qq = q;
411      pIter(q);
412    }
413    while ((q!=NULL) && (pGetComp(q)==k));
414    *p = q;
415    pNext(qq) = NULL;
416  }
417  if (q==NULL) return result;
418  if (pGetComp(q) > k)
419  {
420    pDecrComp(q);
421    pSetmComp(q);
422  }
423  poly pNext_q;
424  while ((pNext_q=pNext(q))!=NULL)
425  {
426    if (pGetComp(pNext_q)==k)
427    {
428      if (result==NULL)
429      {
430        result = pNext_q;
431        qq = result;
432      }
433      else
434      {
435        pNext(qq) = pNext_q;
436        pIter(qq);
437      }
438      pNext(q) = pNext(pNext_q);
439      pNext(qq) =NULL;
440      pSetComp(qq,0);
441      pSetmComp(qq);
442    }
443    else
444    {
445      /*pIter(q);*/ q=pNext_q;
446      if (pGetComp(q) > k)
447      {
448        pDecrComp(q);
449        pSetmComp(q);
450      }
451    }
452  }
453  return result;
454}
455
456// Splits *p into two polys: *q which consists of all monoms with
457// component == comp and *p of all other monoms *lq == pLength(*q)
458void pTakeOutComp(poly *r_p, Exponent_t comp, poly *r_q, int *lq)
459{
460  spolyrec pp, qq;
461  poly p, q, p_prev;
462  int l = 0;
463
464#ifdef HAVE_ASSUME
465  int lp = pLength(*r_p);
466#endif
467
468  pNext(&pp) = *r_p;
469  p = *r_p;
470  p_prev = &pp;
471  q = &qq;
472
473  while(p != NULL)
474  {
475    while (pGetComp(p) == comp)
476    {
477      pNext(q) = p;
478      pIter(q);
479      pSetComp(p, 0);
480      pSetmComp(p);
481      pIter(p);
482      l++;
483      if (p == NULL)
484      {
485        pNext(p_prev) = NULL;
486        goto Finish;
487      }
488    }
489    pNext(p_prev) = p;
490    p_prev = p;
491    pIter(p);
492  }
493
494  Finish:
495  pNext(q) = NULL;
496  *r_p = pNext(&pp);
497  *r_q = pNext(&qq);
498  *lq = l;
499#ifdef HAVE_ASSUME
500  assume(pLength(*r_p) + pLength(*r_q) == lp);
501#endif
502  pTest(*r_p);
503  pTest(*r_q);
504}
505
506void pDecrOrdTakeOutComp(poly *r_p, Exponent_t comp, Order_t order,
507                         poly *r_q, int *lq)
508{
509  spolyrec pp, qq;
510  poly p, q, p_prev;
511  int l = 0;
512
513  pNext(&pp) = *r_p;
514  p = *r_p;
515  p_prev = &pp;
516  q = &qq;
517
518#ifdef HAVE_ASSUME
519  if (p != NULL)
520  {
521    while (pNext(p) != NULL)
522    {
523      assume(pGetOrder(p) >= pGetOrder(pNext(p)));
524      pIter(p);
525    }
526  }
527  p = *r_p;
528#endif
529
530  while (p != NULL && pGetOrder(p) > order) pIter(p);
531
532  while(p != NULL && pGetOrder(p) == order)
533  {
534    while (pGetComp(p) == comp)
535    {
536      pNext(q) = p;
537      pIter(q);
538      pIter(p);
539      pSetComp(p, 0);
540      pSetmComp(p);
541      l++;
542      if (p == NULL || pGetOrder(p) != order)
543      {
544        pNext(p_prev) = p;
545        goto Finish;
546      }
547    }
548    pNext(p_prev) = p;
549    p_prev = p;
550    pIter(p);
551  }
552
553  Finish:
554  pNext(q) = NULL;
555  *r_p = pNext(&pp);
556  *r_q = pNext(&qq);
557  *lq = l;
558}
559
560#if 1
561poly pTakeOutComp1(poly * p, int k)
562{
563  poly q = *p;
564
565  if (q==NULL) return NULL;
566
567  poly qq=NULL,result = NULL;
568
569  if (pGetComp(q)==k)
570  {
571    result = q; /* *p */
572    while ((q!=NULL) && (pGetComp(q)==k))
573    {
574      pSetComp(q,0);
575      pSetmComp(q);
576      qq = q;
577      pIter(q);
578    }
579    *p = q;
580    pNext(qq) = NULL;
581  }
582  if (q==NULL) return result;
583//  if (pGetComp(q) > k) pGetComp(q)--;
584  while (pNext(q)!=NULL)
585  {
586    if (pGetComp(pNext(q))==k)
587    {
588      if (result==NULL)
589      {
590        result = pNext(q);
591        qq = result;
592      }
593      else
594      {
595        pNext(qq) = pNext(q);
596        pIter(qq);
597      }
598      pNext(q) = pNext(pNext(q));
599      pNext(qq) =NULL;
600      pSetComp(qq,0);
601      pSetmComp(qq);
602    }
603    else
604    {
605      pIter(q);
606//      if (pGetComp(q) > k) pGetComp(q)--;
607    }
608  }
609  return result;
610}
611#endif
612
613void pDeleteComp(poly * p,int k)
614{
615  poly q;
616
617  while ((*p!=NULL) && (pGetComp(*p)==k)) pDeleteLm(p);
618  if (*p==NULL) return;
619  q = *p;
620  if (pGetComp(q)>k)
621  {
622    pDecrComp(q);
623    pSetmComp(q);
624  }
625  while (pNext(q)!=NULL)
626  {
627    if (pGetComp(pNext(q))==k)
628      pDeleteLm(&(pNext(q)));
629    else
630    {
631      pIter(q);
632      if (pGetComp(q)>k)
633      {
634        pDecrComp(q);
635        pSetmComp(q);
636      }
637    }
638  }
639}
640/*----------end of utilities for syzygies--------------*/
641
642/*2
643* pair has no common factor ? or is no polynomial
644*/
645BOOLEAN pHasNotCF(poly p1, poly p2)
646{
647
648  if (!TEST_OPT_IDLIFT)
649  {
650    if (pGetComp(p1) > 0 || pGetComp(p2) > 0)
651      return FALSE;
652  }
653  int i = 1;
654  loop
655  {
656    if ((pGetExp(p1, i) > 0) && (pGetExp(p2, i) > 0))   return FALSE;
657    if (i == pVariables)                                return TRUE;
658    i++;
659  }
660}
661
662#ifdef HAVE_RING2TOM
663number nGetUnit(number k) {
664  long test = (long) k;
665  while (test%2 == 0) {
666    test = test / 2;
667  }
668  return (number) test;
669}
670#endif
671
672/*2
673*divides p1 by its leading coefficient
674*/
675void pNorm(poly p1)
676{
677  poly h;
678  number k, c;
679#ifdef HAVE_RING2TOM
680  if (currRing->cring != 0)
681  {
682    if (p1!=NULL)
683    {
684      k = nGetUnit(pGetCoeff(p1));
685      if (!nIsOne(k))
686      {
687        k = nGetUnit(pGetCoeff(p1));
688        c = nDiv(pGetCoeff(p1), k);
689        pSetCoeff0(p1, c);
690        h = pNext(p1);
691        while (h != NULL)
692        {
693          c = nDiv(pGetCoeff(h), k);
694          pSetCoeff(h, c);
695          pIter(h);
696        }
697        nDelete(&k);
698      }
699     return;
700    }
701  }
702#endif
703  if (p1!=NULL)
704  {
705    if (pNext(p1)==NULL)
706    {
707      pSetCoeff(p1,nInit(1));
708      return;
709    }
710    if (!nIsOne(pGetCoeff(p1)))
711    {
712      nNormalize(pGetCoeff(p1));
713      k = pGetCoeff(p1);
714      c = nInit(1);
715      pSetCoeff0(p1,c);
716      h = pNext(p1);
717      while (h!=NULL)
718      {
719        c=nDiv(pGetCoeff(h),k);
720        if (!nIsOne(c)) nNormalize(c);
721        pSetCoeff(h,c);
722        pIter(h);
723      }
724      nDelete(&k);
725    }
726    else
727    {
728      if (nNormalize != nDummy2)
729      {
730        h = pNext(p1);
731        while (h!=NULL)
732        {
733          nNormalize(pGetCoeff(h));
734          pIter(h);
735        }
736      }
737    }
738  }
739}
740
741/*2
742*normalize all coefficients
743*/
744void p_Normalize(poly p, ring r)
745{
746  if (rField_has_simple_inverse(r)) return; /* Z/p, GF(p,n), R, long R/C */
747  while (p!=NULL)
748  {
749    if (currRing==r) {nTest(pGetCoeff(p));}
750    n_Normalize(pGetCoeff(p),r);
751    pIter(p);
752  }
753}
754
755// splits p into polys with Exp(n) == 0 and Exp(n) != 0
756// Poly with Exp(n) != 0 is reversed
757static void pSplitAndReversePoly(poly p, int n, poly *non_zero, poly *zero)
758{
759  if (p == NULL)
760  {
761    *non_zero = NULL;
762    *zero = NULL;
763    return;
764  }
765  spolyrec sz;
766  poly z, n_z, next;
767  z = &sz;
768  n_z = NULL;
769
770  while(p != NULL)
771  {
772    next = pNext(p);
773    if (pGetExp(p, n) == 0)
774    {
775      pNext(z) = p;
776      pIter(z);
777    }
778    else
779    {
780      pNext(p) = n_z;
781      n_z = p;
782    }
783    p = next;
784  }
785  pNext(z) = NULL;
786  *zero = pNext(&sz);
787  *non_zero = n_z;
788  return;
789}
790
791/*3
792* substitute the n-th variable by 1 in p
793* destroy p
794*/
795static poly pSubst1 (poly p,int n)
796{
797  poly qq=NULL, result = NULL;
798  poly zero=NULL, non_zero=NULL;
799
800  // reverse, so that add is likely to be linear
801  pSplitAndReversePoly(p, n, &non_zero, &zero);
802
803  while (non_zero != NULL)
804  {
805    assume(pGetExp(non_zero, n) != 0);
806    qq = non_zero;
807    pIter(non_zero);
808    qq->next = NULL;
809    pSetExp(qq,n,0);
810    pSetm(qq);
811    result = pAdd(result,qq);
812  }
813  p = pAdd(result, zero);
814  pTest(p);
815  return p;
816}
817
818/*3
819* substitute the n-th variable by number e in p
820* destroy p
821*/
822static poly pSubst2 (poly p,int n, number e)
823{
824  assume( ! nIsZero(e) );
825  poly qq,result = NULL;
826  number nn, nm;
827  poly zero, non_zero;
828
829  // reverse, so that add is likely to be linear
830  pSplitAndReversePoly(p, n, &non_zero, &zero);
831
832  while (non_zero != NULL)
833  {
834    assume(pGetExp(non_zero, n) != 0);
835    qq = non_zero;
836    pIter(non_zero);
837    qq->next = NULL;
838    nPower(e, pGetExp(qq, n), &nn);
839    nm = nMult(nn, pGetCoeff(qq));
840    pSetCoeff(qq, nm);
841    nDelete(&nn);
842    pSetExp(qq, n, 0);
843    pSetm(qq);
844    result = pAdd(result,qq);
845  }
846  p = pAdd(result, zero);
847  pTest(p);
848  return p;
849}
850
851
852/* delete monoms whose n-th exponent is different from zero */
853poly pSubst0(poly p, int n)
854{
855  spolyrec res;
856  poly h = &res;
857  pNext(h) = p;
858
859  while (pNext(h)!=NULL)
860  {
861    if (pGetExp(pNext(h),n)!=0)
862    {
863      pDeleteLm(&pNext(h));
864    }
865    else
866    {
867      pIter(h);
868    }
869  }
870  pTest(pNext(&res));
871  return pNext(&res);
872}
873
874/*2
875* substitute the n-th variable by e in p
876* destroy p
877*/
878poly pSubst(poly p, int n, poly e)
879{
880  if (e == NULL) return pSubst0(p, n);
881
882  if (pIsConstant(e))
883  {
884    if (nIsOne(pGetCoeff(e))) return pSubst1(p,n);
885    else return pSubst2(p, n, pGetCoeff(e));
886  }
887
888#ifdef HAVE_PLURAL
889  if (rIsPluralRing(currRing))
890  {
891    return nc_pSubst(p,n,e);
892  }
893#endif
894
895  int exponent,i;
896  poly h, res, m;
897  int *me,*ee;
898  number nu,nu1;
899
900  me=(int *)omAlloc((pVariables+1)*sizeof(int));
901  ee=(int *)omAlloc((pVariables+1)*sizeof(int));
902  if (e!=NULL) pGetExpV(e,ee);
903  res=NULL;
904  h=p;
905  while (h!=NULL)
906  {
907    if ((e!=NULL) || (pGetExp(h,n)==0))
908    {
909      m=pHead(h);
910      pGetExpV(m,me);
911      exponent=me[n];
912      me[n]=0;
913      for(i=pVariables;i>0;i--)
914        me[i]+=exponent*ee[i];
915      pSetExpV(m,me);
916      if (e!=NULL)
917      {
918        nPower(pGetCoeff(e),exponent,&nu);
919        nu1=nMult(pGetCoeff(m),nu);
920        nDelete(&nu);
921        pSetCoeff(m,nu1);
922      }
923      res=pAdd(res,m);
924    }
925    pDeleteLm(&h);
926  }
927  omFreeSize((ADDRESS)me,(pVariables+1)*sizeof(int));
928  omFreeSize((ADDRESS)ee,(pVariables+1)*sizeof(int));
929  return res;
930}
931
932/* Returns TRUE if
933     * LM(p) | LM(lcm)
934     * LC(p) | LC(lcm) only if ring
935     * Exists i, j:
936         * LE(p, i)  != LE(lcm, i)
937         * LE(p1, i) != LE(lcm, i)   ==> LCM(p1, p) != lcm
938         * LE(p, j)  != LE(lcm, j)
939         * LE(p2, j) != LE(lcm, j)   ==> LCM(p2, p) != lcm
940*/
941BOOLEAN pCompareChain (poly p,poly p1,poly p2,poly lcm)
942{
943  int k, j;
944
945  if (lcm==NULL) return FALSE;
946
947  for (j=pVariables; j; j--)
948    if ( pGetExp(p,j) >  pGetExp(lcm,j)) return FALSE;
949  if ( pGetComp(p) !=  pGetComp(lcm)) return FALSE;
950  for (j=pVariables; j; j--)
951  {
952    if (pGetExp(p1,j)!=pGetExp(lcm,j))
953    {
954      if (pGetExp(p,j)!=pGetExp(lcm,j))
955      {
956        for (k=pVariables; k>j; k--)
957        {
958          if ((pGetExp(p,k)!=pGetExp(lcm,k))
959          && (pGetExp(p2,k)!=pGetExp(lcm,k)))
960            return TRUE;
961        }
962        for (k=j-1; k; k--)
963        {
964          if ((pGetExp(p,k)!=pGetExp(lcm,k))
965          && (pGetExp(p2,k)!=pGetExp(lcm,k)))
966            return TRUE;
967        }
968        return FALSE;
969      }
970    }
971    else if (pGetExp(p2,j)!=pGetExp(lcm,j))
972    {
973      if (pGetExp(p,j)!=pGetExp(lcm,j))
974      {
975        for (k=pVariables; k>j; k--)
976        {
977          if ((pGetExp(p,k)!=pGetExp(lcm,k))
978          && (pGetExp(p1,k)!=pGetExp(lcm,k)))
979            return TRUE;
980        }
981        for (k=j-1; k!=0 ; k--)
982        {
983          if ((pGetExp(p,k)!=pGetExp(lcm,k))
984          && (pGetExp(p1,k)!=pGetExp(lcm,k)))
985            return TRUE;
986        }
987        return FALSE;
988      }
989    }
990  }
991  return FALSE;
992}
Note: See TracBrowser for help on using the repository browser.