source: git/kernel/polys.cc @ a2466f

spielwiese
Last change on this file since a2466f was a2466f, checked in by Oliver Wienand <wienand@…>, 18 years ago
[oliver] @hannes: Please check changes to kutil.cc(pDivComp). kutil.cc: * pDivComp now distinguishes between equal and incomparable monoms * changes to chainCritRing, in progress ringgb.*: * changes to testGB, now tests also inclusion of I polys.cc: * just doc git-svn-id: file:///usr/local/Singular/svn/trunk@9195 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 18.4 KB
RevLine 
[35aab3]1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
[a2466f]4/* $Id: polys.cc,v 1.9 2006-06-12 00:07:11 wienand Exp $ */
[35aab3]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
[68349d]22#ifdef HAVE_PLURAL
23#include "gring.h"
24#endif
25
[35aab3]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
[8c5988]52void pSetGlobals(const ring r, BOOLEAN complete)
[35aab3]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      s = eati(s,&i);
198      p_AddExp(rc,1+j, (Exponent_t)i, r);
199    }
200    else
201    {
202      s--;
203      return s;
204    }
205  }
206done:
207  if (r->cf->nIsZero(pGetCoeff(rc))) p_DeleteLm(&rc,r);
208  else
209  {
210    p_Setm(rc,r);
211  }
212  return s;
213}
214
215poly pmInit(char *st, BOOLEAN &ok)
216{
217  poly p;
218  char *s=p_Read(st,p,currRing);
219  if (*s!='\0')
220  {
221    if ((s!=st)&&isdigit(st[0]))
222    {
223      errorreported=TRUE;
224    }
225    ok=FALSE;
226    pDelete(&p);
227    return NULL;
228  }
229  ok=!errorreported;
230  return p;
231}
232
233/*2
234*make p homogeneous by multiplying the monomials by powers of x_varnum
235*/
236poly pHomogen (poly p, int varnum)
237{
238  poly q=NULL, qn;
239  int  o,ii;
240  sBucket_pt bp;
241
242  if (p!=NULL)
243  {
244    if ((varnum < 1) || (varnum > pVariables))
245    {
246      return NULL;
247    }
248    o=pWTotaldegree(p);
249    q=pNext(p);
250    while (q != NULL)
251    {
252      ii=pWTotaldegree(q);
253      if (ii>o) o=ii;
254      pIter(q);
255    }
256    q = pCopy(p);
257    bp = sBucketCreate(currRing);
258    while (q != NULL)
259    {
260      ii = o-pWTotaldegree(q);
261      if (ii!=0)
262      {
263        pAddExp(q,varnum, (Exponent_t)ii);
264        pSetm(q);
265      }
266      qn = pNext(q);
267      pNext(q) = NULL;
268      sBucket_Add_p(bp, q, 1);
269      q = qn;
270    }
271    sBucketDestroyAdd(bp, &q, &ii);
272  }
273  return q;
274}
275
276/*2
277*replaces the maximal powers of the leading monomial of p2 in p1 by
278*the same powers of n, utility for dehomogenization
279*/
280poly pDehomogen (poly p1,poly p2,number n)
281{
282  polyset P;
283  int     SizeOfSet=5;
284  int     i;
285  poly    p;
286  number  nn;
287
288  P = (polyset)omAlloc0(5*sizeof(poly));
289  //for (i=0; i<5; i++)
290  //{
291  //  P[i] = NULL;
292  //}
293  pCancelPolyByMonom(p1,p2,&P,&SizeOfSet);
294  p = P[0];
295  //P[0] = NULL ;// for safety, may be removed later
296  for (i=1; i<SizeOfSet; i++)
297  {
298    if (P[i] != NULL)
299    {
300      nPower(n,i,&nn);
301      pMult_nn(P[i],nn);
302      p = pAdd(p,P[i]);
303      //P[i] =NULL; // for safety, may be removed later
304      nDelete(&nn);
305    }
306  }
307  omFreeSize((ADDRESS)P,SizeOfSet*sizeof(poly));
308  return p;
309}
310
311/*4
312*Returns the exponent of the maximal power of the leading monomial of
313*p2 in that of p1
314*/
315static int pGetMaxPower (poly p1,poly p2)
316{
317  int     i,k,res = 32000; /*a very large integer*/
318
319  if (p1 == NULL) return 0;
320  for (i=1; i<=pVariables; i++)
321  {
322    if ( pGetExp(p2,i) != 0)
323    {
324      k =  pGetExp(p1,i) /  pGetExp(p2,i);
325      if (k < res) res = k;
326    }
327  }
328  return res;
329}
330
331/*2
332*Returns as i-th entry of P the coefficient of the (i-1) power of
333*the leading monomial of p2 in p1
334*/
335void pCancelPolyByMonom (poly p1,poly p2,polyset * P,int * SizeOfSet)
336{
337  int   maxPow;
338  poly  p,qp,Coeff;
339
340  if (*P == NULL)
341  {
342    *P = (polyset) omAlloc(5*sizeof(poly));
343    *SizeOfSet = 5;
344  }
345  p = pCopy(p1);
346  while (p != NULL)
347  {
348    qp = p->next;
349    p->next = NULL;
350    maxPow = pGetMaxPower(p,p2);
351    Coeff = pDivByMonom(p,p2);
352    if (maxPow > *SizeOfSet)
353    {
354      pEnlargeSet(P,*SizeOfSet,maxPow+1-*SizeOfSet);
355      *SizeOfSet = maxPow+1;
356    }
357    (*P)[maxPow] = pAdd((*P)[maxPow],Coeff);
358    pDelete(&p);
359    p = qp;
360  }
361}
362
363/*2
364*returns the leading monomial of p1 divided by the maximal power of that
365*of p2
366*/
367poly pDivByMonom (poly p1,poly p2)
368{
369  int     k, i;
370
371  if (p1 == NULL) return NULL;
372  k = pGetMaxPower(p1,p2);
373  if (k == 0)
374    return pHead(p1);
375  else
376  {
377    number n;
378    poly p = pInit();
379
380    p->next = NULL;
381    for (i=1; i<=pVariables; i++)
382    {
383       pSetExp(p,i, pGetExp(p1,i)-k* pGetExp(p2,i));
384    }
385    nPower(p2->coef,k,&n);
386    pSetCoeff0(p,nDiv(p1->coef,n));
387    nDelete(&n);
388    pSetm(p);
389    return p;
390  }
391}
392/*----------utilities for syzygies--------------*/
393poly pTakeOutComp(poly * p, int k)
394{
395  poly q = *p,qq=NULL,result = NULL;
396
397  if (q==NULL) return NULL;
398  if (pGetComp(q)==k)
399  {
400    result = q;
401    while ((q!=NULL) && (pGetComp(q)==k))
402    {
403      pSetComp(q,0);
404      pSetmComp(q);
405      qq = q;
406      pIter(q);
407    }
408    *p = q;
409    pNext(qq) = NULL;
410  }
411  if (q==NULL) return result;
412  if (pGetComp(q) > k)
413  {
414    pDecrComp(q);
415    pSetmComp(q);
416  }
417  poly pNext_q;
418  while ((pNext_q=pNext(q))!=NULL)
419  {
420    if (pGetComp(pNext_q)==k)
421    {
422      if (result==NULL)
423      {
424        result = pNext_q;
425        qq = result;
426      }
427      else
428      {
429        pNext(qq) = pNext_q;
430        pIter(qq);
431      }
432      pNext(q) = pNext(pNext_q);
433      pNext(qq) =NULL;
434      pSetComp(qq,0);
435      pSetmComp(qq);
436    }
437    else
438    {
439      /*pIter(q);*/ q=pNext_q;
440      if (pGetComp(q) > k)
441      {
442        pDecrComp(q);
443        pSetmComp(q);
444      }
445    }
446  }
447  return result;
448}
449
450// Splits *p into two polys: *q which consists of all monoms with
451// component == comp and *p of all other monoms *lq == pLength(*q)
452void pTakeOutComp(poly *r_p, Exponent_t comp, poly *r_q, int *lq)
453{
454  spolyrec pp, qq;
455  poly p, q, p_prev;
456  int l = 0;
457
458#ifdef HAVE_ASSUME
459  int lp = pLength(*r_p);
460#endif
461
462  pNext(&pp) = *r_p;
463  p = *r_p;
464  p_prev = &pp;
465  q = &qq;
466
467  while(p != NULL)
468  {
469    while (pGetComp(p) == comp)
470    {
471      pNext(q) = p;
472      pIter(q);
473      pSetComp(p, 0);
474      pSetmComp(p);
475      pIter(p);
476      l++;
477      if (p == NULL)
478      {
479        pNext(p_prev) = NULL;
480        goto Finish;
481      }
482    }
483    pNext(p_prev) = p;
484    p_prev = p;
485    pIter(p);
486  }
487
488  Finish:
489  pNext(q) = NULL;
490  *r_p = pNext(&pp);
491  *r_q = pNext(&qq);
492  *lq = l;
493#ifdef HAVE_ASSUME
494  assume(pLength(*r_p) + pLength(*r_q) == lp);
495#endif
496  pTest(*r_p);
497  pTest(*r_q);
498}
499
500void pDecrOrdTakeOutComp(poly *r_p, Exponent_t comp, Order_t order,
501                         poly *r_q, int *lq)
502{
503  spolyrec pp, qq;
504  poly p, q, p_prev;
505  int l = 0;
506
507  pNext(&pp) = *r_p;
508  p = *r_p;
509  p_prev = &pp;
510  q = &qq;
511
512#ifdef HAVE_ASSUME
513  if (p != NULL)
514  {
515    while (pNext(p) != NULL)
516    {
517      assume(pGetOrder(p) >= pGetOrder(pNext(p)));
518      pIter(p);
519    }
520  }
521  p = *r_p;
522#endif
523
524  while (p != NULL && pGetOrder(p) > order) pIter(p);
525
526  while(p != NULL && pGetOrder(p) == order)
527  {
528    while (pGetComp(p) == comp)
529    {
530      pNext(q) = p;
531      pIter(q);
532      pIter(p);
533      pSetComp(p, 0);
534      pSetmComp(p);
535      l++;
536      if (p == NULL || pGetOrder(p) != order)
537      {
538        pNext(p_prev) = p;
539        goto Finish;
540      }
541    }
542    pNext(p_prev) = p;
543    p_prev = p;
544    pIter(p);
545  }
546
547  Finish:
548  pNext(q) = NULL;
549  *r_p = pNext(&pp);
550  *r_q = pNext(&qq);
551  *lq = l;
552}
553
554#if 1
555poly pTakeOutComp1(poly * p, int k)
556{
557  poly q = *p;
558
559  if (q==NULL) return NULL;
560
561  poly qq=NULL,result = NULL;
562
563  if (pGetComp(q)==k)
564  {
565    result = q; /* *p */
566    while ((q!=NULL) && (pGetComp(q)==k))
567    {
568      pSetComp(q,0);
569      pSetmComp(q);
570      qq = q;
571      pIter(q);
572    }
573    *p = q;
574    pNext(qq) = NULL;
575  }
576  if (q==NULL) return result;
577//  if (pGetComp(q) > k) pGetComp(q)--;
578  while (pNext(q)!=NULL)
579  {
580    if (pGetComp(pNext(q))==k)
581    {
582      if (result==NULL)
583      {
584        result = pNext(q);
585        qq = result;
586      }
587      else
588      {
589        pNext(qq) = pNext(q);
590        pIter(qq);
591      }
592      pNext(q) = pNext(pNext(q));
593      pNext(qq) =NULL;
594      pSetComp(qq,0);
595      pSetmComp(qq);
596    }
597    else
598    {
599      pIter(q);
600//      if (pGetComp(q) > k) pGetComp(q)--;
601    }
602  }
603  return result;
604}
605#endif
606
607void pDeleteComp(poly * p,int k)
608{
609  poly q;
610
611  while ((*p!=NULL) && (pGetComp(*p)==k)) pDeleteLm(p);
612  if (*p==NULL) return;
613  q = *p;
614  if (pGetComp(q)>k)
615  {
616    pDecrComp(q);
617    pSetmComp(q);
618  }
619  while (pNext(q)!=NULL)
620  {
621    if (pGetComp(pNext(q))==k)
622      pDeleteLm(&(pNext(q)));
623    else
624    {
625      pIter(q);
626      if (pGetComp(q)>k)
627      {
628        pDecrComp(q);
629        pSetmComp(q);
630      }
631    }
632  }
633}
634/*----------end of utilities for syzygies--------------*/
635
636/*2
637* pair has no common factor ? or is no polynomial
638*/
639BOOLEAN pHasNotCF(poly p1, poly p2)
640{
641
642  if (pGetComp(p1) > 0 || pGetComp(p2) > 0)
643    return FALSE;
644  int i = 1;
645  loop
646  {
647    if ((pGetExp(p1, i) > 0) && (pGetExp(p2, i) > 0))   return FALSE;
648    if (i == pVariables)                                return TRUE;
649    i++;
650  }
651}
652
[585bbcb]653#ifdef HAVE_RING2TOM
654number nGetUnit(number k) {
655  long test = (long) k;
656  while (test%2 == 0) {
657    test = test / 2;
658  }
659  return (number) test;
660}
661#endif
[35aab3]662
663/*2
664*divides p1 by its leading coefficient
665*/
666void pNorm(poly p1)
667{
668  poly h;
669  number k, c;
[585bbcb]670#ifdef HAVE_RING2TOM
671  if (currRing->cring != 0)
672  {
673    if (p1!=NULL)
674    {
675      k = nGetUnit(pGetCoeff(p1));
676      if (!nIsOne(k))
677      {
678        k = nGetUnit(pGetCoeff(p1));
679        c = nDiv(pGetCoeff(p1), k);
680        pSetCoeff0(p1, c);
681        h = pNext(p1);
682        while (h != NULL)
683        {
684          c = nDiv(pGetCoeff(h), k);
685          pSetCoeff(h, c);
686          pIter(h);
687        }
688        nDelete(&k);
689      }
690     return;
691    }
692  }
[f92547]693#endif
[35aab3]694  if (p1!=NULL)
695  {
[5e8fe91]696    if (pNext(p1)==NULL)
697    {
[98938c]698      pSetCoeff(p1,nInit(1));
[5e8fe91]699      return;
700    }
[35aab3]701    if (!nIsOne(pGetCoeff(p1)))
702    {
703      nNormalize(pGetCoeff(p1));
[585bbcb]704      k = pGetCoeff(p1);
[35aab3]705      c = nInit(1);
706      pSetCoeff0(p1,c);
707      h = pNext(p1);
708      while (h!=NULL)
709      {
710        c=nDiv(pGetCoeff(h),k);
711        if (!nIsOne(c)) nNormalize(c);
712        pSetCoeff(h,c);
713        pIter(h);
714      }
715      nDelete(&k);
716    }
717    else
718    {
719      if (nNormalize != nDummy2)
720      {
721        h = pNext(p1);
722        while (h!=NULL)
723        {
724          nNormalize(pGetCoeff(h));
725          pIter(h);
726        }
727      }
728    }
729  }
730}
731
732/*2
733*normalize all coefficients
734*/
[585bbcb]735void p_Normalize(poly p, ring r)
[35aab3]736{
737  if (rField_has_simple_inverse(r)) return; /* Z/p, GF(p,n), R, long R/C */
738  while (p!=NULL)
739  {
740    if (currRing==r) {nTest(pGetCoeff(p));}
741    n_Normalize(pGetCoeff(p),r);
742    pIter(p);
743  }
744}
745
746// splits p into polys with Exp(n) == 0 and Exp(n) != 0
747// Poly with Exp(n) != 0 is reversed
748static void pSplitAndReversePoly(poly p, int n, poly *non_zero, poly *zero)
749{
750  if (p == NULL)
751  {
752    *non_zero = NULL;
753    *zero = NULL;
754    return;
755  }
756  spolyrec sz;
757  poly z, n_z, next;
758  z = &sz;
759  n_z = NULL;
760
761  while(p != NULL)
762  {
763    next = pNext(p);
764    if (pGetExp(p, n) == 0)
765    {
766      pNext(z) = p;
767      pIter(z);
768    }
769    else
770    {
771      pNext(p) = n_z;
772      n_z = p;
773    }
774    p = next;
775  }
776  pNext(z) = NULL;
777  *zero = pNext(&sz);
778  *non_zero = n_z;
779  return;
780}
781
782/*3
783* substitute the n-th variable by 1 in p
784* destroy p
785*/
786static poly pSubst1 (poly p,int n)
787{
[788529d]788  poly qq=NULL, result = NULL;
789  poly zero=NULL, non_zero=NULL;
[35aab3]790
791  // reverse, so that add is likely to be linear
792  pSplitAndReversePoly(p, n, &non_zero, &zero);
793
794  while (non_zero != NULL)
795  {
796    assume(pGetExp(non_zero, n) != 0);
797    qq = non_zero;
798    pIter(non_zero);
799    qq->next = NULL;
800    pSetExp(qq,n,0);
801    pSetm(qq);
802    result = pAdd(result,qq);
803  }
804  p = pAdd(result, zero);
805  pTest(p);
806  return p;
807}
808
809/*3
810* substitute the n-th variable by number e in p
811* destroy p
812*/
813static poly pSubst2 (poly p,int n, number e)
814{
815  assume( ! nIsZero(e) );
816  poly qq,result = NULL;
817  number nn, nm;
818  poly zero, non_zero;
819
820  // reverse, so that add is likely to be linear
821  pSplitAndReversePoly(p, n, &non_zero, &zero);
822
823  while (non_zero != NULL)
824  {
825    assume(pGetExp(non_zero, n) != 0);
826    qq = non_zero;
827    pIter(non_zero);
828    qq->next = NULL;
829    nPower(e, pGetExp(qq, n), &nn);
830    nm = nMult(nn, pGetCoeff(qq));
831    pSetCoeff(qq, nm);
832    nDelete(&nn);
833    pSetExp(qq, n, 0);
834    pSetm(qq);
835    result = pAdd(result,qq);
836  }
837  p = pAdd(result, zero);
838  pTest(p);
839  return p;
840}
841
842
843/* delete monoms whose n-th exponent is different from zero */
844poly pSubst0(poly p, int n)
845{
846  spolyrec res;
847  poly h = &res;
848  pNext(h) = p;
849
850  while (pNext(h)!=NULL)
851  {
852    if (pGetExp(pNext(h),n)!=0)
853    {
854      pDeleteLm(&pNext(h));
855    }
856    else
857    {
858      pIter(h);
859    }
860  }
861  pTest(pNext(&res));
862  return pNext(&res);
863}
864
865/*2
866* substitute the n-th variable by e in p
867* destroy p
868*/
869poly pSubst(poly p, int n, poly e)
870{
871  if (e == NULL) return pSubst0(p, n);
872
873  if (pIsConstant(e))
874  {
875    if (nIsOne(pGetCoeff(e))) return pSubst1(p,n);
876    else return pSubst2(p, n, pGetCoeff(e));
877  }
878
[68349d]879#ifdef HAVE_PLURAL
880  if (rIsPluralRing(currRing))
881  {
882    return nc_pSubst(p,n,e);
883  }
884#endif
885
[35aab3]886  int exponent,i;
887  poly h, res, m;
888  int *me,*ee;
889  number nu,nu1;
890
891  me=(int *)omAlloc((pVariables+1)*sizeof(int));
892  ee=(int *)omAlloc((pVariables+1)*sizeof(int));
893  if (e!=NULL) pGetExpV(e,ee);
894  res=NULL;
895  h=p;
896  while (h!=NULL)
897  {
898    if ((e!=NULL) || (pGetExp(h,n)==0))
899    {
900      m=pHead(h);
901      pGetExpV(m,me);
902      exponent=me[n];
903      me[n]=0;
904      for(i=pVariables;i>0;i--)
905        me[i]+=exponent*ee[i];
906      pSetExpV(m,me);
907      if (e!=NULL)
908      {
909        nPower(pGetCoeff(e),exponent,&nu);
910        nu1=nMult(pGetCoeff(m),nu);
911        nDelete(&nu);
912        pSetCoeff(m,nu1);
913      }
914      res=pAdd(res,m);
915    }
916    pDeleteLm(&h);
917  }
918  omFreeSize((ADDRESS)me,(pVariables+1)*sizeof(int));
919  omFreeSize((ADDRESS)ee,(pVariables+1)*sizeof(int));
920  return res;
921}
922
[a2466f]923/* Returns TRUE if
924     * LM(p) | LM(lcm)
925     * LC(p) | LC(lcm) only if ring
926     * Exists i, j:
927         * LE(p, i)  != LE(lcm, i)
928         * LE(p1, i) != LE(lcm, i)   ==> LCM(p1, p) != lcm
929         * LE(p, j)  != LE(lcm, j)
930         * LE(p2, j) != LE(lcm, j)   ==> LCM(p2, p) != lcm
931*/
[35aab3]932BOOLEAN pCompareChain (poly p,poly p1,poly p2,poly lcm)
933{
934  int k, j;
935
936  if (lcm==NULL) return FALSE;
[a2466f]937  #ifdef HAVE_RING2TOM
[f92547]938  // In coefficient rings, the coefficient plays a role in chain crit TODO
939  if (currRing->cring == 1 && !pLmDivisibleByNoComp(p, lcm)) return FALSE;
940  #endif
[35aab3]941
942  for (j=pVariables; j; j--)
943    if ( pGetExp(p,j) >  pGetExp(lcm,j)) return FALSE;
944  if ( pGetComp(p) !=  pGetComp(lcm)) return FALSE;
945  for (j=pVariables; j; j--)
946  {
947    if (pGetExp(p1,j)!=pGetExp(lcm,j))
948    {
949      if (pGetExp(p,j)!=pGetExp(lcm,j))
950      {
951        for (k=pVariables; k>j; k--)
952        {
953          if ((pGetExp(p,k)!=pGetExp(lcm,k))
954          && (pGetExp(p2,k)!=pGetExp(lcm,k)))
955            return TRUE;
956        }
957        for (k=j-1; k; k--)
958        {
959          if ((pGetExp(p,k)!=pGetExp(lcm,k))
960          && (pGetExp(p2,k)!=pGetExp(lcm,k)))
961            return TRUE;
962        }
963        return FALSE;
964      }
965    }
966    else if (pGetExp(p2,j)!=pGetExp(lcm,j))
967    {
968      if (pGetExp(p,j)!=pGetExp(lcm,j))
969      {
970        for (k=pVariables; k>j; k--)
971        {
972          if ((pGetExp(p,k)!=pGetExp(lcm,k))
973          && (pGetExp(p1,k)!=pGetExp(lcm,k)))
974            return TRUE;
975        }
976        for (k=j-1; k!=0 ; k--)
977        {
978          if ((pGetExp(p,k)!=pGetExp(lcm,k))
979          && (pGetExp(p1,k)!=pGetExp(lcm,k)))
980            return TRUE;
981        }
982        return FALSE;
983      }
984    }
985  }
986  return FALSE;
987}
Note: See TracBrowser for help on using the repository browser.