source: git/kernel/clapsing.cc @ 64a88e

spielwiese
Last change on this file since 64a88e was 64a88e, checked in by Hans Schönemann <hannes@…>, 14 years ago
*hannes: remove unused parameter git-svn-id: file:///usr/local/Singular/svn/trunk@12060 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 39.8 KB
Line 
1// emacs edit mode for this file is -*- C++ -*-
2/****************************************
3*  Computer Algebra System SINGULAR     *
4****************************************/
5// $Id: clapsing.cc,v 1.43 2009-08-13 15:31:27 Singular Exp $
6/*
7* ABSTRACT: interface between Singular and factory
8*/
9
10//#define FACTORIZE2_DEBUG
11#include "mod2.h"
12#include "omalloc.h"
13#ifdef HAVE_FACTORY
14#define SI_DONT_HAVE_GLOBAL_VARS
15#include "structs.h"
16#include "clapsing.h"
17#include "numbers.h"
18#include "ring.h"
19#include "ffields.h"
20#include <factory.h>
21#include "clapconv.h"
22#include <factor.h>
23//CanonicalForm algcd(const CanonicalForm & F, const CanonicalForm & g, const CFList & as, const Varlist & order);
24CanonicalForm alg_gcd(const CanonicalForm &, const CanonicalForm &, const CFList &);
25#include "ring.h"
26
27void out_cf(char *s1,const CanonicalForm &f,char *s2);
28
29poly singclap_gcd_r ( poly f, poly g, const ring r )
30{
31  // assume pCleardenom is done
32  // assume f!=0, g!=0
33  poly res=NULL;
34
35  if (p_IsConstantPoly(f,r) || p_IsConstantPoly(g,r))
36  {
37    return p_One(r);
38  }
39
40  // for now there is only the possibility to handle polynomials over
41  // Q and Fp ...
42  Off(SW_RATIONAL);
43  if (rField_is_Q(r) || (rField_is_Zp(r)))
44  {
45    CanonicalForm newGCD(const CanonicalForm & A, const CanonicalForm & B);
46    setCharacteristic( n_GetChar(r) );
47    CanonicalForm F( convSingPFactoryP( f,r ) ), G( convSingPFactoryP( g, r ) );
48    //if (nGetChar() > 1 )
49    //{
50    //  res=convFactoryPSingP( newGCD( F,G ));
51    //  if (!nGreaterZero(pGetCoeff(res))) res=pNeg(res);
52    //}
53    //else
54      res=convFactoryPSingP( gcd( F, G ) , r);
55  }
56  // and over Q(a) / Fp(a)
57  else if (( n_GetChar(r)==1 ) /* Q(a) */
58  || (n_GetChar(r) <-1))       /* Fp(a) */
59  {
60    if (n_GetChar(r)==1) setCharacteristic( 0 );
61    else                 setCharacteristic( -n_GetChar(r) );
62    if (r->minpoly!=NULL)
63    {
64    #if 0
65      if (( n_GetChar(r)==1 ) /* Q(a) */ && (!isOn(SW_USE_QGCD)))
66      {
67      //  WerrorS( feNotImplemented );
68        CanonicalForm mipo=convSingPFactoryP(((lnumber)r->minpoly)->z,r->algring);
69        //Varlist ord;
70        //ord.append(mipo.mvar());
71        CFList as(mipo);
72        Variable a=rootOf(mipo);
73        //CanonicalForm F( convSingAPFactoryAP( f,a,r ) ), G( convSingAPFactoryAP( g,a,r ) );
74        CanonicalForm F( convSingTrPFactoryP(f,r) ), G( convSingTrPFactoryP(g,r) );
75        //res= convFactoryAPSingAP( algcd( F, G, as, ord) );
76        //res= convFactoryAPSingAP( alg_gcd( F, G, as) );
77        res= convFactoryAPSingAP( alg_gcd( F, G, as),r );
78      }
79      else
80      #endif
81      {
82        bool b=isOn(SW_USE_QGCD);
83        if ( nGetChar()==1 ) On(SW_USE_QGCD);
84        CanonicalForm mipo=convSingPFactoryP(((lnumber)r->minpoly)->z,
85                                           r->algring);
86        Variable a=rootOf(mipo);
87        CanonicalForm F( convSingAPFactoryAP( f,a ) ),
88                      G( convSingAPFactoryAP( g,a ) );
89        res= convFactoryAPSingAP( gcd( F, G ) );
90        if (!b) Off(SW_USE_QGCD);
91      }
92    }
93    else
94    {
95      CanonicalForm F( convSingTrPFactoryP( f,r ) ), G( convSingTrPFactoryP( g,r ) );
96      res= convFactoryPSingTrP( gcd( F, G ),r );
97    }
98  }
99  #if 0
100  else if (( n_GetChar(r)>1 )&&(r->parameter!=NULL)) /* GF(q) */
101  {
102    int p=rChar(r);
103    int n=2;
104    int t=p*p;
105    while (t!=n_Char(r)) { t*=p;n++; }
106    setCharacteristic(p,n,'a');
107    CanonicalForm F( convSingGFFactoryGF( f,r ) ), G( convSingGFFactoryGF( g,r ) );
108    res= convFactoryGFSingGF( gcd( F, G ),r );
109  }
110  #endif
111  else
112    WerrorS( feNotImplemented );
113
114  Off(SW_RATIONAL);
115  return res;
116}
117
118poly singclap_gcd ( poly f, poly g )
119{
120  poly res=NULL;
121
122  if (f!=NULL) pCleardenom(f);
123  if (g!=NULL) pCleardenom(g);
124  else         return f; // g==0 => gcd=f (but do a pCleardenom)
125  if (f==NULL) return g; // f==0 => gcd=g (but do a pCleardenom)
126
127  res=singclap_gcd_r(f,g,currRing);
128  pDelete(&f);
129  pDelete(&g);
130  return res;
131}
132
133/*2 find the maximal exponent of var(i) in poly p*/
134int pGetExp_Var(poly p, int i)
135{
136  int m=0;
137  int mm;
138  while (p!=NULL)
139  {
140    mm=pGetExp(p,i);
141    if (mm>m) m=mm;
142    pIter(p);
143  }
144  return m;
145}
146
147poly singclap_resultant ( poly f, poly g , poly x)
148{
149  int i=pIsPurePower(x);
150  if (i==0)
151  {
152    WerrorS("3rd argument must be a ring variable");
153    return NULL;
154  }
155  if ((f==NULL) || (g==NULL))
156    return NULL;
157  // for now there is only the possibility to handle polynomials over
158  // Q and Fp ...
159  if (( nGetChar() == 0 || nGetChar() > 1 )
160  && (currRing->parameter==NULL))
161  {
162    Variable X(i);
163    setCharacteristic( nGetChar() );
164    CanonicalForm F( convSingPFactoryP( f ) ), G( convSingPFactoryP( g ) );
165    poly res=convFactoryPSingP( resultant( F, G, X ) );
166    Off(SW_RATIONAL);
167    return res;
168  }
169  // and over Q(a) / Fp(a)
170  else if (( nGetChar()==1 ) /* Q(a) */
171  || (nGetChar() <-1))       /* Fp(a) */
172  {
173    if (nGetChar()==1) setCharacteristic( 0 );
174    else               setCharacteristic( -nGetChar() );
175    poly res;
176    Variable X(i+rPar(currRing));
177    if (currRing->minpoly!=NULL)
178    {
179      //Variable X(i);
180      CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
181                                                 currRing->algring);
182      Variable a=rootOf(mipo);
183      CanonicalForm F( convSingAPFactoryAP( f,a ) ), G( convSingAPFactoryAP( g,a ) );
184      res= convFactoryAPSingAP( resultant( F, G, X ) );
185    }
186    else
187    {
188      //Variable X(i+rPar(currRing));
189      number nf,ng;
190      pCleardenom_n(f,nf);pCleardenom_n(g,ng);
191      int ef,eg;
192      ef=pGetExp_Var(f,i);
193      eg=pGetExp_Var(g,i);
194      CanonicalForm F( convSingTrPFactoryP( f ) ), G( convSingTrPFactoryP( g ) );
195      res= convFactoryPSingTrP( resultant( F, G, X ) );
196      if ((nf!=NULL)&&(!nIsOne(nf))&&(!nIsZero(nf)))
197      {
198        number n=nInvers(nf);
199        while(eg>0)
200        {
201          res=pMult_nn(res,n);
202          eg--;
203        }
204        nDelete(&n);
205      }
206      nDelete(&nf);
207      if ((ng!=NULL)&&(!nIsOne(ng))&&(!nIsZero(ng)))
208      {
209        number n=nInvers(ng);
210        while(ef>0)
211        {
212          res=pMult_nn(res,n);
213          ef--;
214        }
215        nDelete(&n);
216      }
217      nDelete(&ng);
218    }
219    Off(SW_RATIONAL);
220    return res;
221  }
222  else
223    WerrorS( feNotImplemented );
224  return NULL;
225}
226//poly singclap_resultant ( poly f, poly g , poly x)
227//{
228//  int i=pVar(x);
229//  if (i==0)
230//  {
231//    WerrorS("ringvar expected");
232//    return NULL;
233//  }
234//  ideal I=idInit(1,1);
235//
236//  // get the coeffs von f wrt. x:
237//  I->m[0]=pCopy(f);
238//  matrix ffi=mpCoeffs(I,i);
239//  ffi->rank=1;
240//  ffi->ncols=ffi->nrows;
241//  ffi->nrows=1;
242//  ideal fi=(ideal)ffi;
243//
244//  // get the coeffs von g wrt. x:
245//  I->m[0]=pCopy(g);
246//  matrix ggi=mpCoeffs(I,i);
247//  ggi->rank=1;
248//  ggi->ncols=ggi->nrows;
249//  ggi->nrows=1;
250//  ideal gi=(ideal)ggi;
251//
252//  // contruct the matrix:
253//  int fn=IDELEMS(fi); //= deg(f,x)+1
254//  int gn=IDELEMS(gi); //= deg(g,x)+1
255//  matrix m=mpNew(fn+gn-2,fn+gn-2);
256//  if(m==NULL)
257//  {
258//    return NULL;
259//  }
260//
261//  // enter the coeffs into m:
262//  int j;
263//  for(i=0;i<gn-1;i++)
264//  {
265//    for(j=0;j<fn;j++)
266//    {
267//      MATELEM(m,i+1,fn-j+i)=pCopy(fi->m[j]);
268//    }
269//  }
270//  for(i=0;i<fn-1;i++)
271//  {
272//    for(j=0;j<gn;j++)
273//    {
274//      MATELEM(m,gn+i,gn-j+i)=pCopy(gi->m[j]);
275//    }
276//  }
277//
278//  poly r=mpDet(m);
279//
280//  idDelete(&fi);
281//  idDelete(&gi);
282//  idDelete((ideal *)&m);
283//  return r;
284//}
285
286BOOLEAN singclap_extgcd ( poly f, poly g, poly &res, poly &pa, poly &pb )
287{
288  // for now there is only the possibility to handle univariate
289  // polynomials over
290  // Q and Fp ...
291  res=NULL;pa=NULL;pb=NULL;
292  On(SW_SYMMETRIC_FF);
293  if (( nGetChar() == 0 || nGetChar() > 1 )
294  && (currRing->parameter==NULL))
295  {
296    setCharacteristic( nGetChar() );
297    CanonicalForm F( convSingPFactoryP( f ) ), G( convSingPFactoryP( g ) );
298    CanonicalForm FpG=F+G;
299    if (!(FpG.isUnivariate()|| FpG.inCoeffDomain()))
300    //if (!F.isUnivariate() || !G.isUnivariate() || F.mvar()!=G.mvar())
301    {
302      Off(SW_RATIONAL);
303      WerrorS("not univariate");
304      return TRUE;
305    }
306    CanonicalForm Fa,Gb;
307    On(SW_RATIONAL);
308    res=convFactoryPSingP( extgcd( F, G, Fa, Gb ) );
309    pa=convFactoryPSingP(Fa);
310    pb=convFactoryPSingP(Gb);
311    Off(SW_RATIONAL);
312  }
313  // and over Q(a) / Fp(a)
314  else if (( nGetChar()==1 ) /* Q(a) */
315  || (nGetChar() <-1))       /* Fp(a) */
316  {
317    if (nGetChar()==1) setCharacteristic( 0 );
318    else               setCharacteristic( -nGetChar() );
319    CanonicalForm Fa,Gb;
320    if (currRing->minpoly!=NULL)
321    {
322      CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
323                                                 currRing->algring);
324      Variable a=rootOf(mipo);
325      CanonicalForm F( convSingAPFactoryAP( f,a ) ), G( convSingAPFactoryAP( g,a ) );
326      CanonicalForm FpG=F+G;
327      if (!(FpG.isUnivariate()|| FpG.inCoeffDomain()))
328      //if (!F.isUnivariate() || !G.isUnivariate() || F.mvar()!=G.mvar())
329      {
330        WerrorS("not univariate");
331        return TRUE;
332      }
333      res= convFactoryAPSingAP( extgcd( F, G, Fa, Gb ) );
334      pa=convFactoryAPSingAP(Fa);
335      pb=convFactoryAPSingAP(Gb);
336    }
337    else
338    {
339      CanonicalForm F( convSingTrPFactoryP( f ) ), G( convSingTrPFactoryP( g ) );
340      CanonicalForm FpG=F+G;
341      if (!(FpG.isUnivariate()|| FpG.inCoeffDomain()))
342      //if (!F.isUnivariate() || !G.isUnivariate() || F.mvar()!=G.mvar())
343      {
344        Off(SW_RATIONAL);
345        WerrorS("not univariate");
346        return TRUE;
347      }
348      res= convFactoryPSingTrP( extgcd( F, G, Fa, Gb ) );
349      pa=convFactoryPSingTrP(Fa);
350      pb=convFactoryPSingTrP(Gb);
351    }
352    Off(SW_RATIONAL);
353  }
354  else
355  {
356    WerrorS( feNotImplemented );
357    return TRUE;
358  }
359  return FALSE;
360}
361
362BOOLEAN singclap_extgcd_r ( poly f, poly g, poly &res, poly &pa, poly &pb, const ring r )
363{
364  // for now there is only the possibility to handle univariate
365  // polynomials over
366  // Q and Fp ...
367  res=NULL;pa=NULL;pb=NULL;
368  On(SW_SYMMETRIC_FF);
369  if (( n_GetChar(r) == 0 || n_GetChar(r) > 1 )
370  && (r->parameter==NULL))
371  {
372    setCharacteristic( n_GetChar(r) );
373    CanonicalForm F( convSingPFactoryP( f,r ) ), G( convSingPFactoryP( g,r) );
374    CanonicalForm FpG=F+G;
375    if (!(FpG.isUnivariate()|| FpG.inCoeffDomain()))
376    //if (!F.isUnivariate() || !G.isUnivariate() || F.mvar()!=G.mvar())
377    {
378      Off(SW_RATIONAL);
379      WerrorS("not univariate");
380      return TRUE;
381    }
382    CanonicalForm Fa,Gb;
383    On(SW_RATIONAL);
384    res=convFactoryPSingP( extgcd( F, G, Fa, Gb ),r );
385    pa=convFactoryPSingP(Fa,r);
386    pb=convFactoryPSingP(Gb,r);
387    Off(SW_RATIONAL);
388  }
389  // and over Q(a) / Fp(a)
390  else if (( n_GetChar(r)==1 ) /* Q(a) */
391  || (n_GetChar(r) <-1))       /* Fp(a) */
392  {
393    if (n_GetChar(r)==1) setCharacteristic( 0 );
394    else                 setCharacteristic( -n_GetChar(r) );
395    CanonicalForm Fa,Gb;
396    if (r->minpoly!=NULL)
397    {
398      CanonicalForm mipo=convSingPFactoryP(((lnumber)r->minpoly)->z,
399                                            r->algring);
400      Variable a=rootOf(mipo);
401      CanonicalForm F( convSingAPFactoryAP( f,a ) ),
402                    G( convSingAPFactoryAP( g,a ) );
403      CanonicalForm FpG=F+G;
404      if (!(FpG.isUnivariate()|| FpG.inCoeffDomain()))
405      //if (!F.isUnivariate() || !G.isUnivariate() || F.mvar()!=G.mvar())
406      {
407        WerrorS("not univariate");
408        return TRUE;
409      }
410      res= convFactoryAPSingAP( extgcd( F, G, Fa, Gb ) );
411      pa=convFactoryAPSingAP(Fa);
412      pb=convFactoryAPSingAP(Gb);
413    }
414    else
415    {
416      CanonicalForm F( convSingTrPFactoryP( f, r ) ), G( convSingTrPFactoryP( g, r ) );
417      CanonicalForm FpG=F+G;
418      if (!(FpG.isUnivariate()|| FpG.inCoeffDomain()))
419      //if (!F.isUnivariate() || !G.isUnivariate() || F.mvar()!=G.mvar())
420      {
421        Off(SW_RATIONAL);
422        WerrorS("not univariate");
423        return TRUE;
424      }
425      res= convFactoryPSingTrP( extgcd( F, G, Fa, Gb ), r );
426      pa=convFactoryPSingTrP(Fa, r);
427      pb=convFactoryPSingTrP(Gb, r);
428    }
429    Off(SW_RATIONAL);
430  }
431  else
432  {
433    WerrorS( feNotImplemented );
434    return TRUE;
435  }
436  return FALSE;
437}
438
439poly singclap_pdivide ( poly f, poly g )
440{
441  // for now there is only the possibility to handle polynomials over
442  // Q and Fp ...
443  poly res=NULL;
444  On(SW_RATIONAL);
445  if (( nGetChar() == 0 || nGetChar() > 1 )
446  && (currRing->parameter==NULL))
447  {
448    setCharacteristic( nGetChar() );
449    CanonicalForm F( convSingPFactoryP( f ) ), G( convSingPFactoryP( g ) );
450    res = convFactoryPSingP( F / G );
451  }
452  // and over Q(a) / Fp(a)
453  else if (( nGetChar()==1 ) /* Q(a) */
454  || (nGetChar() <-1))       /* Fp(a) */
455  {
456    if (nGetChar()==1) setCharacteristic( 0 );
457    else               setCharacteristic( -nGetChar() );
458    if (currRing->minpoly!=NULL)
459    {
460      CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
461                                                 currRing->algring);
462      Variable a=rootOf(mipo);
463      CanonicalForm F( convSingAPFactoryAP( f,a ) ), G( convSingAPFactoryAP( g,a ) );
464      res= convFactoryAPSingAP(  F / G  );
465    }
466    else
467    {
468      CanonicalForm F( convSingTrPFactoryP( f ) ), G( convSingTrPFactoryP( g ) );
469      res= convFactoryPSingTrP(  F / G  );
470    }
471  }
472  #if 0 // not yet working
473  else if (rField_is_GF())
474  {
475    //Print("GF(%d^%d)\n",nfCharP,nfMinPoly[0]);
476    setCharacteristic( nfCharP,nfMinPoly[0], currRing->parameter[0][0] );
477    CanonicalForm F( convSingGFFactoryGF( f ) ), G( convSingGFFactoryGF( g ) );
478    res = convFactoryGFSingGF( F / G );
479  }
480  #endif
481  else
482    WerrorS( feNotImplemented );
483  Off(SW_RATIONAL);
484  return res;
485}
486
487void singclap_divide_content ( poly f )
488{
489  if ( f==NULL )
490  {
491    return;
492  }
493  else  if ( pNext( f ) == NULL )
494  {
495    pSetCoeff( f, nInit( 1 ) );
496    return;
497  }
498  else
499  {
500    if ( nGetChar() == 1 )
501      setCharacteristic( 0 );
502    else  if ( nGetChar() == -1 )
503      return; /* not implemented for R */
504    else  if ( nGetChar() < 0 )
505      setCharacteristic( -nGetChar() );
506    else if (currRing->parameter==NULL) /* not GF(q) */
507      setCharacteristic( nGetChar() );
508    else
509      return; /* not implemented*/
510
511    CFList L;
512    CanonicalForm g, h;
513    poly p = pNext(f);
514
515    // first attemp: find 2 smallest g:
516
517    number g1=pGetCoeff(f);
518    number g2=pGetCoeff(p); // p==pNext(f);
519    pIter(p);
520    int sz1=nSize(g1);
521    int sz2=nSize(g2);
522    if (sz1>sz2)
523    {
524      number gg=g1;
525      g1=g2; g2=gg;
526      int sz=sz1;
527      sz1=sz2; sz2=sz;
528    }
529    while (p!=NULL)
530    {
531      int n_sz=nSize(pGetCoeff(p));
532      if (n_sz<sz1)
533      {
534        sz2=sz1;
535        g2=g1;
536        g1=pGetCoeff(p);
537        sz1=n_sz;
538        if (sz1<=3) break;
539      }
540      else if(n_sz<sz2)
541      {
542        sz2=n_sz;
543        g2=pGetCoeff(p);
544        sz2=n_sz;
545      }
546      pIter(p);
547    }
548    g = convSingPFactoryP( ((lnumber)g1)->z, currRing->algring );
549    g = gcd( g, convSingPFactoryP( ((lnumber)g2)->z , currRing->algring));
550
551    // second run: gcd's
552
553    p = f;
554    while ( (p != NULL) && (g != 1)  && ( g != 0))
555    {
556      h = convSingPFactoryP( ((lnumber)pGetCoeff(p))->z, currRing->algring );
557      pIter( p );
558
559      g = gcd( g, h );
560
561      L.append( h );
562    }
563    if (( g == 1 ) || (g == 0))
564    {
565      // pTest(f);
566      return;
567    }
568    else
569    {
570      CFListIterator i;
571      for ( i = L, p = f; i.hasItem(); i++, p=pNext(p) )
572      {
573        lnumber c=(lnumber)pGetCoeff(p);
574        napDelete(&c->z);
575        c->z=convFactoryPSingP( i.getItem() / g, currRing->algring );
576        //nTest((number)c);
577        //#ifdef LDEBUG
578        //number cn=(number)c;
579        //StringSetS(""); nWrite(nt); StringAppend(" ==> ");
580        //nWrite(cn);PrintS(StringAppend("\n"));
581        //#endif
582      }
583    }
584    // pTest(f);
585  }
586}
587
588static int primepower(int c)
589{
590  int p=1;
591  int cc=c;
592  while(cc!= rInternalChar(currRing)) { cc*=c; p++; }
593  return p;
594}
595
596BOOLEAN count_Factors(ideal I, intvec *v,int j, poly &f, poly fac)
597{
598  pTest(f);
599  pTest(fac);
600  int e=0;
601  if (!pIsConstantPoly(fac))
602  {
603#ifdef FACTORIZE2_DEBUG
604    printf("start count_Factors(%d), Fdeg=%d, factor deg=%d\n",j,pTotaldegree(f),pTotaldegree(fac));
605    p_wrp(fac,currRing);PrintLn();
606#endif
607    On(SW_RATIONAL);
608    CanonicalForm F, FAC,Q,R;
609    Variable a;
610    if (( nGetChar() == 0 || nGetChar() > 1 )
611    && (currRing->parameter==NULL))
612    {
613      F=convSingPFactoryP( f );
614      FAC=convSingPFactoryP( fac );
615    }
616    // and over Q(a) / Fp(a)
617    else if (( nGetChar()==1 ) /* Q(a) */
618    || (nGetChar() <-1))       /* Fp(a) */
619    {
620      if (currRing->minpoly!=NULL)
621      {
622        CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
623                                    currRing->algring);
624        a=rootOf(mipo);
625        F=convSingAPFactoryAP( f,a );
626        FAC=convSingAPFactoryAP( fac,a );
627      }
628      else
629      {
630        F=convSingTrPFactoryP( f );
631        FAC=convSingTrPFactoryP( fac );
632      }
633    }
634    else
635      WerrorS( feNotImplemented );
636
637    poly q;
638    loop
639    {
640      Q=F;
641      Q/=FAC;
642      R=Q;
643      R*=FAC;
644      R-=F;
645      if (R.isZero())
646      {
647        if (( nGetChar() == 0 || nGetChar() > 1 )
648        && (currRing->parameter==NULL))
649        {
650          q = convFactoryPSingP( Q );
651        }
652        else if (( nGetChar()==1 ) /* Q(a) */
653        || (nGetChar() <-1))       /* Fp(a) */
654        {
655          if (currRing->minpoly!=NULL)
656          {
657            q= convFactoryAPSingAP(  Q  );
658          }
659          else
660          {
661            q= convFactoryPSingTrP(  Q  );
662          }
663        }
664        e++; pDelete(&f); f=q; q=NULL; F=Q;
665      }
666      else
667      {
668        break;
669      }
670    }
671    if (e==0)
672    {
673      Off(SW_RATIONAL);
674      return FALSE;
675    }
676  }
677  else e=1;
678  I->m[j]=fac;
679  if (v!=NULL) (*v)[j]=e;
680  Off(SW_RATIONAL);
681  return TRUE;
682}
683
684int singclap_factorize_retry;
685extern int libfac_interruptflag;
686
687ideal singclap_factorize ( poly f, intvec ** v , int with_exps)
688{
689  pTest(f);
690#ifdef FACTORIZE2_DEBUG
691  printf("singclap_factorize, degree %d\n",pTotaldegree(f));
692#endif
693  // with_exps: 3,1 return only true factors, no exponents
694  //            2 return true factors and exponents
695  //            0 return coeff, factors and exponents
696  BOOLEAN save_errorreported=errorreported;
697
698  ideal res=NULL;
699
700  // handle factorize(0) =========================================
701  if (f==NULL)
702  {
703    res=idInit(1,1);
704    if (with_exps!=1)
705    {
706      (*v)=new intvec(1);
707      (**v)[0]=1;
708    }
709    return res;
710  }
711  // handle factorize(mon) =========================================
712  if (pNext(f)==NULL)
713  {
714    int i=0;
715    int n=0;
716    int e;
717    for(i=pVariables;i>0;i--) if(pGetExp(f,i)!=0) n++;
718    if (with_exps==0) n++; // with coeff
719    res=idInit(si_max(n,1),1);
720    switch(with_exps)
721    {
722      case 0: // with coef & exp.
723        res->m[0]=pOne();
724        pSetCoeff(res->m[0],nCopy(pGetCoeff(f)));
725        // no break
726      case 2: // with exp.
727        (*v)=new intvec(si_max(1,n));
728        (**v)[0]=1;
729        // no break
730      case 1: ;
731      #ifdef TEST
732      default: ;
733      #endif
734    }
735    if (n==0)
736    {
737      res->m[0]=pOne();
738      // (**v)[0]=1; is already done
739      return res;
740    }
741    for(i=pVariables;i>0;i--)
742    {
743      e=pGetExp(f,i);
744      if(e!=0)
745      {
746        n--;
747        poly p=pOne();
748        pSetExp(p,i,1);
749        pSetm(p);
750        res->m[n]=p;
751        if (with_exps!=1) (**v)[n]=e;
752      }
753    }
754    return res;
755  }
756  //PrintS("S:");pWrite(f);PrintLn();
757  // use factory/libfac in general ==============================
758  Off(SW_RATIONAL);
759  On(SW_SYMMETRIC_FF);
760  #ifdef HAVE_NTL
761  extern int prime_number;
762  if(rField_is_Q()) prime_number=0;
763  #endif
764  CFFList L;
765  number N=NULL;
766  number NN=NULL;
767  number old_lead_coeff=nCopy(pGetCoeff(f));
768
769  if (!rField_is_Zp() && !rField_is_Zp_a()) /* Q, Q(a) */
770  {
771    //if (f!=NULL) // already tested at start of routine
772    {
773      number n0=nCopy(pGetCoeff(f));
774      if (with_exps==0)
775        N=nCopy(n0);
776      pCleardenom(f);
777      NN=nDiv(n0,pGetCoeff(f));
778      nDelete(&n0);
779      if (with_exps==0)
780      {
781        nDelete(&N);
782        N=nCopy(NN);
783      }
784    }
785  }
786  else if (rField_is_Zp_a())
787  {
788    //if (f!=NULL) // already tested at start of routine
789    if (singclap_factorize_retry==0)
790    {
791      number n0=nCopy(pGetCoeff(f));
792      if (with_exps==0)
793        N=nCopy(n0);
794      pNorm(f);
795      pCleardenom(f);
796      NN=nDiv(n0,pGetCoeff(f));
797      nDelete(&n0);
798      if (with_exps==0)
799      {
800        nDelete(&N);
801        N=nCopy(NN);
802      }
803    }
804  }
805  if (rField_is_Q() || rField_is_Zp())
806  {
807    setCharacteristic( nGetChar() );
808    CanonicalForm F( convSingPFactoryP( f ) );
809    if (nGetChar()==0) /* Q */
810    {
811      L = factorize( F );
812    }
813    else /* Fp */
814    {
815      do
816      {
817        libfac_interruptflag=0;
818        L = Factorize( F );
819      }
820      while ((libfac_interruptflag!=0) ||(L.isEmpty()));
821    }
822  }
823  #if 0
824  else if (rField_is_GF())
825  {
826    int c=rChar(currRing);
827    setCharacteristic( c, primepower(c) );
828    CanonicalForm F( convSingGFFactoryGF( f ) );
829    if (F.isUnivariate())
830    {
831      L = factorize( F );
832    }
833    else
834    {
835      goto notImpl;
836    }
837  }
838  #endif
839  // and over Q(a) / Fp(a)
840  else if (rField_is_Extension())
841  {
842    if (rField_is_Q_a()) setCharacteristic( 0 );
843    else                 setCharacteristic( -nGetChar() );
844    if (currRing->minpoly!=NULL)
845    {
846      CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
847                   currRing->algring);
848      Variable a=rootOf(mipo);
849      CanonicalForm F( convSingAPFactoryAP( f,a ) );
850      if (rField_is_Zp_a() && F.isUnivariate())
851      {
852        L = factorize( F, a );
853      }
854      else
855      {
856        CanonicalForm G( convSingTrPFactoryP( f ) );
857        //  over Q(a) / multivariate over Fp(a)
858        do
859        {
860          libfac_interruptflag=0;
861          L=Factorize2(G, mipo);
862        }
863        while ((libfac_interruptflag!=0) ||(L.isEmpty()));
864        #ifdef FACTORIZE2_DEBUG
865        printf("while okay\n");
866        #endif
867        libfac_interruptflag=0;
868      }
869    }
870    else
871    {
872      CanonicalForm F( convSingTrPFactoryP( f ) );
873      if (rField_is_Q_a())
874      {
875        L = factorize( F );
876      }
877      else /* Fp(a) */
878      {
879        L = Factorize( F );
880      }
881    }
882  }
883  else
884  {
885    goto notImpl;
886  }
887  {
888    poly ff=pCopy(f); // a copy for the retry stuff
889    // the first factor should be a constant
890    if ( ! L.getFirst().factor().inCoeffDomain() )
891      L.insert(CFFactor(1,1));
892    // convert into ideal
893    int n = L.length();
894    if (n==0) n=1;
895    CFFListIterator J=L;
896    int j=0;
897    if (with_exps!=1)
898    {
899      if ((with_exps==2)&&(n>1))
900      {
901        n--;
902        J++;
903      }
904      *v = new intvec( n );
905    }
906    res = idInit( n ,1);
907    for ( ; J.hasItem(); J++, j++ )
908    {
909      poly p;
910      if (with_exps!=1) (**v)[j] = J.getItem().exp();
911      if (rField_is_Zp() || rField_is_Q())           /* Q, Fp */
912      {
913        //count_Factors(res,*v,f, j, convFactoryPSingP( J.getItem().factor() );
914        res->m[j] = convFactoryPSingP( J.getItem().factor() );
915      }
916      #if 0
917      else if (rField_is_GF())
918        res->m[j] = convFactoryGFSingGF( J.getItem().factor() );
919      #endif
920      else if (rField_is_Extension())     /* Q(a), Fp(a) */
921      {
922        intvec *w=NULL;
923        if (v!=NULL) w=*v;
924        if (currRing->minpoly==NULL)
925        {
926          if(!count_Factors(res,w,j,ff,convFactoryPSingTrP( J.getItem().factor() )))
927          {
928            if (w!=NULL)
929              (*w)[j]=1;
930            res->m[j]=pOne();
931          }
932        }
933        else
934        {
935          if (!count_Factors(res,w,j,ff,convFactoryAPSingAP( J.getItem().factor() )))
936          {
937            if (w!=NULL)
938              (*w)[j]=1;
939            res->m[j]=pOne();
940          }
941        }
942      }
943    }
944    if (rField_is_Extension() && (!pIsConstantPoly(ff)))
945    {
946      singclap_factorize_retry++;
947      if (singclap_factorize_retry<3)
948      {
949        int jj;
950        #ifdef FACTORIZE2_DEBUG
951        printf("factorize_retry\n");
952        #endif
953        intvec *ww=NULL;
954        idTest(res);
955        ideal h=singclap_factorize ( ff, &ww , with_exps);
956        idTest(h);
957        int l=(*v)->length();
958        (*v)->resize(l+ww->length());
959        for(jj=0;jj<ww->length();jj++)
960          (**v)[jj+l]=(*ww)[jj];
961        delete ww;
962        ideal hh=idInit(IDELEMS(res)+IDELEMS(h),1);
963        for(jj=IDELEMS(res)-1;jj>=0;jj--)
964        {
965          hh->m[jj]=res->m[jj];
966          res->m[jj]=NULL;
967        }
968        for(jj=IDELEMS(h)-1;jj>=0;jj--)
969        {
970          hh->m[jj+IDELEMS(res)]=h->m[jj];
971          h->m[jj]=NULL;
972        }
973        idDelete(&res);
974        idDelete(&h);
975        res=hh;
976        idTest(res);
977        ff=NULL;
978      }
979      else
980      {
981        WarnS("problem with factorize");
982        #if 0
983        pWrite(ff);
984        idShow(res);
985        #endif
986        idDelete(&res);
987        res=idInit(2,1);
988        res->m[0]=pOne();
989        res->m[1]=ff; ff=NULL;
990      }
991    }
992    pDelete(&ff);
993    if (N!=NULL)
994    {
995      pMult_nn(res->m[0],N);
996      nDelete(&N);
997      N=NULL;
998    }
999    // delete constants
1000    if (res!=NULL)
1001    {
1002      int i=IDELEMS(res)-1;
1003      int j=0;
1004      for(;i>=0;i--)
1005      {
1006        if ((res->m[i]!=NULL)
1007        && (pNext(res->m[i])==NULL)
1008        && (pIsConstant(res->m[i])))
1009        {
1010          if (with_exps!=0)
1011          {
1012            pDelete(&(res->m[i]));
1013            if ((v!=NULL) && ((*v)!=NULL))
1014              (**v)[i]=0;
1015            j++;
1016          }
1017          else if (i!=0)
1018          {
1019            while ((v!=NULL) && ((*v)!=NULL) && ((**v)[i]>1))
1020            {
1021              res->m[0]=pMult(res->m[0],pCopy(res->m[i]));
1022              (**v)[i]--;
1023            }
1024            res->m[0]=pMult(res->m[0],res->m[i]);
1025            res->m[i]=NULL;
1026            if ((v!=NULL) && ((*v)!=NULL))
1027              (**v)[i]=1;
1028            j++;
1029          }
1030        }
1031      }
1032      if (j>0)
1033      {
1034        idSkipZeroes(res);
1035        if ((v!=NULL) && ((*v)!=NULL))
1036        {
1037          intvec *w=*v;
1038          int len=IDELEMS(res);
1039          *v = new intvec( len );
1040          for (i=0,j=0;i<si_min(w->length(),len);i++)
1041          {
1042            if((*w)[i]!=0)
1043            {
1044              (**v)[j]=(*w)[i]; j++;
1045            }
1046          }
1047          delete w;
1048        }
1049      }
1050      if (res->m[0]==NULL)
1051      {
1052        res->m[0]=pOne();
1053      }
1054    }
1055  }
1056  if (rField_is_Q_a() && (currRing->minpoly!=NULL))
1057  {
1058    int i=IDELEMS(res)-1;
1059    int stop=1;
1060    if (with_exps!=0) stop=0;
1061    for(;i>=stop;i--)
1062    {
1063      pNorm(res->m[i]);
1064    }
1065    if (with_exps==0) pSetCoeff(res->m[0],old_lead_coeff);
1066    else nDelete(&old_lead_coeff);
1067  }
1068  else
1069    nDelete(&old_lead_coeff);
1070  errorreported=save_errorreported;
1071notImpl:
1072  if (res==NULL)
1073    WerrorS( feNotImplemented );
1074  if (NN!=NULL)
1075  {
1076    nDelete(&NN);
1077  }
1078  if (N!=NULL)
1079  {
1080    nDelete(&N);
1081  }
1082  //if (f!=NULL) pDelete(&f);
1083  //PrintS("......S\n");
1084  return res;
1085}
1086ideal singclap_sqrfree ( poly f)
1087{
1088  pTest(f);
1089#ifdef FACTORIZE2_DEBUG
1090  printf("singclap_sqrfree, degree %d\n",pTotaldegree(f));
1091#endif
1092  // with_exps: 3,1 return only true factors, no exponents
1093  //            2 return true factors and exponents
1094  //            0 return coeff, factors and exponents
1095  BOOLEAN save_errorreported=errorreported;
1096
1097  ideal res=NULL;
1098
1099  // handle factorize(0) =========================================
1100  if (f==NULL)
1101  {
1102    res=idInit(1,1);
1103    return res;
1104  }
1105  // handle factorize(mon) =========================================
1106  if (pNext(f)==NULL)
1107  {
1108    int i=0;
1109    int n=0;
1110    int e;
1111    for(i=pVariables;i>0;i--) if(pGetExp(f,i)!=0) n++;
1112    n++; // with coeff
1113    res=idInit(si_max(n,1),1);
1114    res->m[0]=pOne();
1115    pSetCoeff(res->m[0],nCopy(pGetCoeff(f)));
1116    if (n==0)
1117    {
1118      res->m[0]=pOne();
1119      // (**v)[0]=1; is already done
1120      return res;
1121    }
1122    for(i=pVariables;i>0;i--)
1123    {
1124      e=pGetExp(f,i);
1125      if(e!=0)
1126      {
1127        n--;
1128        poly p=pOne();
1129        pSetExp(p,i,1);
1130        pSetm(p);
1131        res->m[n]=p;
1132      }
1133    }
1134    return res;
1135  }
1136  //PrintS("S:");pWrite(f);PrintLn();
1137  // use factory/libfac in general ==============================
1138  Off(SW_RATIONAL);
1139  On(SW_SYMMETRIC_FF);
1140  #ifdef HAVE_NTL
1141  extern int prime_number;
1142  if(rField_is_Q()) prime_number=0;
1143  #endif
1144  CFFList L;
1145
1146  if (!rField_is_Zp() && !rField_is_Zp_a()) /* Q, Q(a) */
1147  {
1148    //if (f!=NULL) // already tested at start of routine
1149    {
1150      pCleardenom(f);
1151    }
1152  }
1153  else if (rField_is_Zp_a())
1154  {
1155    //if (f!=NULL) // already tested at start of routine
1156    if (singclap_factorize_retry==0)
1157    {
1158      pNorm(f);
1159      pCleardenom(f);
1160    }
1161  }
1162  if (rField_is_Q() || rField_is_Zp())
1163  {
1164    setCharacteristic( nGetChar() );
1165    CanonicalForm F( convSingPFactoryP( f ) );
1166    L = sqrFree( F );
1167  }
1168  #if 0
1169  else if (rField_is_GF())
1170  {
1171    int c=rChar(currRing);
1172    setCharacteristic( c, primepower(c) );
1173    CanonicalForm F( convSingGFFactoryGF( f ) );
1174    if (F.isUnivariate())
1175    {
1176      L = factorize( F );
1177    }
1178    else
1179    {
1180      goto notImpl;
1181    }
1182  }
1183  #endif
1184  // and over Q(a) / Fp(a)
1185  else if (rField_is_Extension())
1186  {
1187    if (rField_is_Q_a()) setCharacteristic( 0 );
1188    else                 setCharacteristic( -nGetChar() );
1189    if (currRing->minpoly!=NULL)
1190    {
1191      CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
1192                    currRing->algring);
1193      Variable a=rootOf(mipo);
1194      CanonicalForm F( convSingAPFactoryAP( f,a ) );
1195      CFFList SqrFreeMV( const CanonicalForm & f , const CanonicalForm & mipo=0) ;
1196
1197      L = SqrFreeMV( F,mipo );
1198      //WarnS("L = sqrFree( F,mipo );");
1199      //L = sqrFree( F );
1200    }
1201    else
1202    {
1203      CanonicalForm F( convSingTrPFactoryP( f ) );
1204      L = sqrFree( F );
1205    }
1206  }
1207  else
1208  {
1209    goto notImpl;
1210  }
1211  {
1212    poly ff=pCopy(f); // a copy for the retry stuff
1213    // convert into ideal
1214    int n = L.length();
1215    if (n==0) n=1;
1216    CFFListIterator J=L;
1217    int j=0;
1218    res = idInit( n ,1);
1219    for ( ; J.hasItem(); J++, j++ )
1220    {
1221      poly p;
1222      if (rField_is_Zp() || rField_is_Q())           /* Q, Fp */
1223        //count_Factors(res,*v,f, j, convFactoryPSingP( J.getItem().factor() );
1224        res->m[j] = convFactoryPSingP( J.getItem().factor() );
1225      #if 0
1226      else if (rField_is_GF())
1227        res->m[j] = convFactoryGFSingGF( J.getItem().factor() );
1228      #endif
1229      else if (rField_is_Extension())     /* Q(a), Fp(a) */
1230      {
1231        if (currRing->minpoly==NULL)
1232          res->m[j]=convFactoryPSingTrP( J.getItem().factor() );
1233        else
1234          res->m[j]=convFactoryAPSingAP( J.getItem().factor() );
1235      }
1236    }
1237    if (res->m[0]==NULL)
1238    {
1239      res->m[0]=pOne();
1240    }
1241  }
1242  errorreported=save_errorreported;
1243notImpl:
1244  if (res==NULL)
1245    WerrorS( feNotImplemented );
1246  return res;
1247}
1248matrix singclap_irrCharSeries ( ideal I)
1249{
1250  if (idIs0(I)) return mpNew(1,1);
1251
1252  // for now there is only the possibility to handle polynomials over
1253  // Q and Fp ...
1254  matrix res=NULL;
1255  int i;
1256  Off(SW_RATIONAL);
1257  On(SW_SYMMETRIC_FF);
1258  CFList L;
1259  ListCFList LL;
1260  if (((nGetChar() == 0) || (nGetChar() > 1) )
1261  && (currRing->parameter==NULL))
1262  {
1263    setCharacteristic( nGetChar() );
1264    for(i=0;i<IDELEMS(I);i++)
1265    {
1266      poly p=I->m[i];
1267      if (p!=NULL)
1268      {
1269        p=pCopy(p);
1270        pCleardenom(p);
1271        L.append(convSingPFactoryP(p));
1272      }
1273    }
1274  }
1275  // and over Q(a) / Fp(a)
1276  else if (( nGetChar()==1 ) /* Q(a) */
1277  || (nGetChar() <-1))       /* Fp(a) */
1278  {
1279    if (nGetChar()==1) setCharacteristic( 0 );
1280    else               setCharacteristic( -nGetChar() );
1281    for(i=0;i<IDELEMS(I);i++)
1282    {
1283      poly p=I->m[i];
1284      if (p!=NULL)
1285      {
1286        p=pCopy(p);
1287        pCleardenom(p);
1288        L.append(convSingTrPFactoryP(p));
1289      }
1290    }
1291  }
1292  else
1293  {
1294    WerrorS( feNotImplemented );
1295    return res;
1296  }
1297
1298  // a very bad work-around --- FIX IT in libfac
1299  // should be fixed as of 2001/6/27
1300  int tries=0;
1301  int m,n;
1302  ListIterator<CFList> LLi;
1303  loop
1304  {
1305    LL=IrrCharSeries(L);
1306    m= LL.length(); // Anzahl Zeilen
1307    n=0;
1308    for ( LLi = LL; LLi.hasItem(); LLi++ )
1309    {
1310      n = si_max(LLi.getItem().length(),n);
1311    }
1312    if ((m!=0) && (n!=0)) break;
1313    tries++;
1314    if (tries>=5) break;
1315  }
1316  if ((m==0) || (n==0))
1317  {
1318    Warn("char_series returns %d x %d matrix from %d input polys (%d)",
1319      m,n,IDELEMS(I)+1,LL.length());
1320    iiWriteMatrix((matrix)I,"I",2,0);
1321    m=si_max(m,1);
1322    n=si_max(n,1);
1323  }
1324  res=mpNew(m,n);
1325  CFListIterator Li;
1326  for ( m=1, LLi = LL; LLi.hasItem(); LLi++, m++ )
1327  {
1328    for (n=1, Li = LLi.getItem(); Li.hasItem(); Li++, n++)
1329    {
1330      if ( (nGetChar() == 0) || (nGetChar() > 1) )
1331        MATELEM(res,m,n)=convFactoryPSingP(Li.getItem());
1332      else
1333        MATELEM(res,m,n)=convFactoryPSingTrP(Li.getItem());
1334    }
1335  }
1336  Off(SW_RATIONAL);
1337  return res;
1338}
1339
1340char* singclap_neworder ( ideal I)
1341{
1342  int i;
1343  Off(SW_RATIONAL);
1344  On(SW_SYMMETRIC_FF);
1345  CFList L;
1346  if (((nGetChar() == 0) || (nGetChar() > 1) )
1347  && (currRing->parameter==NULL))
1348  {
1349    setCharacteristic( nGetChar() );
1350    for(i=0;i<IDELEMS(I);i++)
1351    {
1352      L.append(convSingPFactoryP(I->m[i]));
1353    }
1354  }
1355  // and over Q(a) / Fp(a)
1356  else if (( nGetChar()==1 ) /* Q(a) */
1357  || (nGetChar() <-1))       /* Fp(a) */
1358  {
1359    if (nGetChar()==1) setCharacteristic( 0 );
1360    else               setCharacteristic( -nGetChar() );
1361    for(i=0;i<IDELEMS(I);i++)
1362    {
1363      L.append(convSingTrPFactoryP(I->m[i]));
1364    }
1365  }
1366  else
1367  {
1368    WerrorS( feNotImplemented );
1369    return NULL;
1370  }
1371
1372  List<int> IL=neworderint(L);
1373  ListIterator<int> Li;
1374  StringSetS("");
1375  Li = IL;
1376  int offs=rPar(currRing);
1377  int* mark=(int*)omAlloc0((pVariables+offs)*sizeof(int));
1378  int cnt=pVariables+offs;
1379  loop
1380  {
1381    if(! Li.hasItem()) break;
1382    BOOLEAN done=TRUE;
1383    i=Li.getItem()-1;
1384    mark[i]=1;
1385    if (i<offs)
1386    {
1387      done=FALSE;
1388      //StringAppendS(currRing->parameter[i]);
1389    }
1390    else
1391    {
1392      StringAppendS(currRing->names[i-offs]);
1393    }
1394    Li++;
1395    cnt--;
1396    if(cnt==0) break;
1397    if (done) StringAppendS(",");
1398  }
1399  for(i=0;i<pVariables+offs;i++)
1400  {
1401    BOOLEAN done=TRUE;
1402    if(mark[i]==0)
1403    {
1404      if (i<offs)
1405      {
1406        done=FALSE;
1407        //StringAppendS(currRing->parameter[i]);
1408      }
1409      else
1410      {
1411        StringAppendS(currRing->names[i-offs]);
1412      }
1413      cnt--;
1414      if(cnt==0) break;
1415      if (done) StringAppendS(",");
1416    }
1417  }
1418  char * s=omStrDup(StringAppendS(""));
1419  if (s[strlen(s)-1]==',') s[strlen(s)-1]='\0';
1420  return s;
1421}
1422
1423BOOLEAN singclap_isSqrFree(poly f)
1424{
1425  BOOLEAN b=FALSE;
1426  Off(SW_RATIONAL);
1427  //  Q / Fp
1428  if (((nGetChar() == 0) || (nGetChar() > 1) )
1429  &&(currRing->parameter==NULL))
1430  {
1431    setCharacteristic( nGetChar() );
1432    CanonicalForm F( convSingPFactoryP( f ) );
1433    if((nGetChar()>1)&&(!F.isUnivariate()))
1434      goto err;
1435    b=(BOOLEAN)isSqrFree(F);
1436  }
1437  // and over Q(a) / Fp(a)
1438  else if (( nGetChar()==1 ) /* Q(a) */
1439  || (nGetChar() <-1))       /* Fp(a) */
1440  {
1441    if (nGetChar()==1) setCharacteristic( 0 );
1442    else               setCharacteristic( -nGetChar() );
1443    //if (currRing->minpoly!=NULL)
1444    //{
1445    //  CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
1446    //                                             currRing->algring);
1447    //  Variable a=rootOf(mipo);
1448    //  CanonicalForm F( convSingAPFactoryAP( f,a ) );
1449    //  ...
1450    //}
1451    //else
1452    {
1453      CanonicalForm F( convSingTrPFactoryP( f ) );
1454      b=(BOOLEAN)isSqrFree(F);
1455    }
1456    Off(SW_RATIONAL);
1457  }
1458  else
1459  {
1460err:
1461    WerrorS( feNotImplemented );
1462  }
1463  return b;
1464}
1465
1466poly singclap_det( const matrix m )
1467{
1468  int r=m->rows();
1469  if (r!=m->cols())
1470  {
1471    Werror("det of %d x %d matrix",r,m->cols());
1472    return NULL;
1473  }
1474  poly res=NULL;
1475  if (( nGetChar() == 0 || nGetChar() > 1 )
1476  && (currRing->parameter==NULL))
1477  {
1478    setCharacteristic( nGetChar() );
1479    CFMatrix M(r,r);
1480    int i,j;
1481    for(i=r;i>0;i--)
1482    {
1483      for(j=r;j>0;j--)
1484      {
1485        M(i,j)=convSingPFactoryP(MATELEM(m,i,j));
1486      }
1487    }
1488    res= convFactoryPSingP( determinant(M,r) ) ;
1489  }
1490  // and over Q(a) / Fp(a)
1491  else if (( nGetChar()==1 ) /* Q(a) */
1492  || (nGetChar() <-1))       /* Fp(a) */
1493  {
1494    if (nGetChar()==1) setCharacteristic( 0 );
1495    else               setCharacteristic( -nGetChar() );
1496    CFMatrix M(r,r);
1497    poly res;
1498    if (currRing->minpoly!=NULL)
1499    {
1500      CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
1501                                            currRing->algring);
1502      Variable a=rootOf(mipo);
1503      int i,j;
1504      for(i=r;i>0;i--)
1505      {
1506        for(j=r;j>0;j--)
1507        {
1508          M(i,j)=convSingAPFactoryAP(MATELEM(m,i,j),a);
1509        }
1510      }
1511      res= convFactoryAPSingAP( determinant(M,r) ) ;
1512    }
1513    else
1514    {
1515      int i,j;
1516      for(i=r;i>0;i--)
1517      {
1518        for(j=r;j>0;j--)
1519        {
1520          M(i,j)=convSingTrPFactoryP(MATELEM(m,i,j));
1521        }
1522      }
1523      res= convFactoryPSingTrP( determinant(M,r) );
1524    }
1525  }
1526  else
1527    WerrorS( feNotImplemented );
1528  Off(SW_RATIONAL);
1529  return res;
1530}
1531
1532int singclap_det_i( intvec * m )
1533{
1534  setCharacteristic( 0 );
1535  CFMatrix M(m->rows(),m->cols());
1536  int i,j;
1537  for(i=m->rows();i>0;i--)
1538  {
1539    for(j=m->cols();j>0;j--)
1540    {
1541      M(i,j)=IMATELEM(*m,i,j);
1542    }
1543  }
1544  int res= convFactoryISingI( determinant(M,m->rows())) ;
1545  Off(SW_RATIONAL);
1546  return res;
1547}
1548napoly singclap_alglcm ( napoly f, napoly g )
1549{
1550
1551 // over Q(a) / Fp(a)
1552 if (nGetChar()==1) setCharacteristic( 0 );
1553 else               setCharacteristic( -nGetChar() );
1554 napoly res;
1555
1556 if (currRing->minpoly!=NULL)
1557 {
1558   CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
1559                                         currRing->algring);
1560   Variable a=rootOf(mipo);
1561   CanonicalForm F( convSingAFactoryA( f,a ) ), G( convSingAFactoryA( g,a ) );
1562   CanonicalForm GCD;
1563
1564   // calculate gcd
1565   GCD = gcd( F, G );
1566
1567   // calculate lcm
1568   res= convFactoryASingA( (F/GCD)*G );
1569 }
1570 else
1571 {
1572   CanonicalForm F( convSingPFactoryP( f,currRing->algring ) ),
1573                 G( convSingPFactoryP( g,currRing->algring ) );
1574   CanonicalForm GCD;
1575   // calculate gcd
1576   GCD = gcd( F, G );
1577
1578   // calculate lcm
1579   res= convFactoryPSingP( (F/GCD)*G, currRing->algring );
1580 }
1581
1582 Off(SW_RATIONAL);
1583 return res;
1584}
1585
1586void singclap_algdividecontent ( napoly f, napoly g, napoly &ff, napoly &gg )
1587{
1588 // over Q(a) / Fp(a)
1589 if (nGetChar()==1) setCharacteristic( 0 );
1590 else               setCharacteristic( -nGetChar() );
1591 ff=gg=NULL;
1592 On(SW_RATIONAL);
1593
1594 if (currRing->minpoly!=NULL)
1595 {
1596   CanonicalForm mipo=convSingPFactoryP(((lnumber)currRing->minpoly)->z,
1597                         currRing->algring);
1598   Variable a=rootOf(mipo);
1599   CanonicalForm F( convSingAFactoryA( f,a ) ), G( convSingAFactoryA( g,a ) );
1600   CanonicalForm GCD;
1601
1602   GCD=gcd( F, G );
1603
1604   if ((GCD!=1) && (GCD!=0))
1605   {
1606     ff= convFactoryASingA( F/ GCD );
1607     gg= convFactoryASingA( G/ GCD );
1608   }
1609 }
1610 else
1611 {
1612   CanonicalForm F( convSingPFactoryP( f,currRing->algring ) ),
1613                 G( convSingPFactoryP( g,currRing->algring ) );
1614   CanonicalForm GCD;
1615
1616   GCD=gcd( F, G );
1617
1618   if ((GCD!=1) && (GCD!=0))
1619   {
1620     ff= convFactoryPSingP( F/ GCD, currRing->algring );
1621     gg= convFactoryPSingP( G/ GCD, currRing->algring );
1622   }
1623 }
1624
1625 Off(SW_RATIONAL);
1626}
1627
1628#if 0
1629lists singclap_chineseRemainder(lists x, lists q)
1630{
1631  //assume(x->nr == q->nr);
1632  //assume(x->nr >= 0);
1633  int n=x->nr+1;
1634  if ((x->nr<0) || (x->nr!=q->nr))
1635  {
1636    WerrorS("list are empty or not of equal length");
1637    return NULL;
1638  }
1639  lists res=(lists)omAlloc0Bin(slists_bin);
1640  CFArray X(1,n), Q(1,n);
1641  int i;
1642  for(i=0; i<n; i++)
1643  {
1644    if (x->m[i-1].Typ()==INT_CMD)
1645    {
1646      X[i]=(int)x->m[i-1].Data();
1647    }
1648    else if (x->m[i-1].Typ()==NUMBER_CMD)
1649    {
1650      number N=(number)x->m[i-1].Data();
1651      X[i]=convSingNFactoryN(N);
1652    }
1653    else
1654    {
1655      WerrorS("illegal type in chineseRemainder");
1656      omFreeBin(res,slists_bin);
1657      return NULL;
1658    }
1659    if (q->m[i-1].Typ()==INT_CMD)
1660    {
1661      Q[i]=(int)q->m[i-1].Data();
1662    }
1663    else if (q->m[i-1].Typ()==NUMBER_CMD)
1664    {
1665      number N=(number)x->m[i-1].Data();
1666      Q[i]=convSingNFactoryN(N);
1667    }
1668    else
1669    {
1670      WerrorS("illegal type in chineseRemainder");
1671      omFreeBin(res,slists_bin);
1672      return NULL;
1673    }
1674  }
1675  CanonicalForm r, prod;
1676  chineseRemainder( X, Q, r, prod );
1677  res->Init(2);
1678  res->m[0].rtyp=NUMBER_CMD;
1679  res->m[1].rtyp=NUMBER_CMD;
1680  res->m[0].data=(char *)convFactoryNSingN( r );
1681  res->m[1].data=(char *)convFactoryNSingN( prod );
1682  return res;
1683}
1684#endif
1685#endif
Note: See TracBrowser for help on using the repository browser.