source: git/libfac/factor/SqrFree.cc @ 4a81ec

spielwiese
Last change on this file since 4a81ec was 4a81ec, checked in by Hans Schönemann <hannes@…>, 26 years ago
* hannes/michael: libfac-0.3.0 git-svn-id: file:///usr/local/Singular/svn/trunk@708 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 11.2 KB
Line 
1/* Copyright 1996 Michael Messollen. All rights reserved. */
2///////////////////////////////////////////////////////////////////////////////
3// emacs edit mode for this file is -*- C++ -*-
4static char * rcsid = "$Id: SqrFree.cc,v 1.3 1997-09-12 07:19:50 Singular Exp $";
5static char * errmsg = "\nYou found a bug!\nPlease inform (Michael Messollen) michael@math.uni-sb.de .\n Please include above information and your input (the ideal/polynomial and characteristic) in your bug-report.\nThank you.";
6///////////////////////////////////////////////////////////////////////////////
7// FACTORY - Includes
8#include<factory.h>
9// Factor - Includes
10#include "tmpl_inst.h"
11#include "helpstuff.h"
12// some CC's need this:
13#include "SqrFree.h"
14
15#ifdef SINGULAR
16#  define HAVE_SINGULAR
17#endif
18
19#ifdef SQRFREEDEBUG
20# define DEBUGOUTPUT
21#else
22# undef DEBUGOUTPUT
23#endif
24
25#include "debug.h"
26#include "timing.h"
27TIMING_DEFINE_PRINT(squarefree_time);
28TIMING_DEFINE_PRINT(gcd_time);
29
30
31static inline CFFactor
32Powerup( const CFFactor & F , int exp=1){ 
33  return CFFactor(F.factor(), exp*F.exp()) ; 
34}
35
36static CFFList
37Powerup( const CFFList & Inputlist , int exp=1 ){
38  CFFList Outputlist;
39
40  for ( CFFListIterator i=Inputlist; i.hasItem(); i++ )
41    Outputlist.append(Powerup(i.getItem(), exp));
42  return Outputlist ;
43}
44
45///////////////////////////////////////////////////////////////
46// Compute the Pth root of a polynomial in characteristic p  //
47// f must be a polynomial which we can take the Pth root of. //
48// Domain is q=p^m , f a uni/multivariate polynomial         //
49///////////////////////////////////////////////////////////////
50static CanonicalForm
51PthRoot( const CanonicalForm & f ){
52  CanonicalForm RES, R = f;
53  int n= max(level(R),getNumVars(R)), p= getCharacteristic();
54 
55  if (n==0){ // constant
56    if (R.inExtension()) // not in prime field; f over |F(q=p^k)
57      R = power(R,Powerup(p,getGFDegree() - 1)) ;
58    // if f in prime field, do nothing
59    return R;
60  }
61  // we assume R is a Pth power here
62  RES = R.genZero();
63  Variable x(n);
64  for (int i=0; i<= (int) (degree(R,level(R))/p) ; i++)
65    RES += PthRoot( R[i*p] ) * power(x,i);
66  return RES;
67}
68
69///////////////////////////////////////////////////////////////
70// A uni/multivariate SqrFreeTest routine.                   //
71// Cheaper to run if all you want is a test.                 //
72// Works for charcteristic 0 and q=p^m                       //
73// Returns 1 if poly r is SqrFree, 0 if SqrFree will do some //
74// kind of factorization.                                    //
75// Would be much more effcient iff we had *good*             //
76//  uni/multivariate gcd's and/or gcdtest's                  //
77///////////////////////////////////////////////////////////////
78int
79SqrFreeTest( const CanonicalForm & r, int opt){
80  CanonicalForm f=r, g;
81  int n=level(f);
82
83  if (getNumVars(f)==0) return 1 ; // a constant is SqrFree
84  if ( f.isUnivariate() ) {
85    g= f.deriv();
86    if ( getCharacteristic() > 0 && g.isZero() ) return 0 ;
87    // Next: it would be best to have a *univariate* gcd-test which returns
88    // 0 iff gcdtest(f,g) == 1 or a constant ( for real Polynomials )
89    g = mygcd(f,g);
90    if ( g.isOne() || (-g).isOne() ) return 1;
91    else 
92      if ( getNumVars(g) == 0 ) return 1;// <- totaldegree!!!
93      else return 0 ;
94  }
95  else { // multivariate case
96    for ( int k=1; k<=n; k++ ) {
97      g = swapvar(f,k,n); g = content(g);
98      // g = 1 || -1 : sqr-free, g poly : not sqr-free, g number : opt helps
99      if ( ! (g.isOne() || (-g).isOne() || getNumVars(g)==0 ) ) {
100        if ( opt==0 ) return 0;
101        else {
102          if ( SqrFreeTest(g,1) == 0 ) return 0;
103          g = swapvar(g,k,n);
104          f /=g ;
105        }
106      }
107    }
108    // Now f is primitive
109    n = level(f); // maybe less indeterminants
110    //    if ( totaldegree(f) <= 1 ) return 1;
111   
112    // Let`s look if it is a Pth root
113    if ( getCharacteristic() > 0 )
114      for (int k=1; k<=n; k++ ) {
115        g=swapvar(f,k,n); g=g.deriv();
116        if ( ! g.isZero() ) break ;
117        else if ( k==n) return 0 ; // really is Pth root
118      }
119    g = f.deriv() ;
120    // Next: it would be best to have a *multivariate* gcd-test which returns
121    // 0 iff gcdtest(f,g) == 1 or a constant ( for real Polynomials )
122    g= mygcd(f,g);
123    if ( g.isOne() || (-g).isOne() || (g==f) || (getNumVars(g)==0) ) return 1 ;
124    else return 0 ;
125  }
126#ifdef HAVE_SINGULAR
127  extern void WerrorS(char *);
128  WerrorS("libfac: ERROR: SqrFreeTest: we should never fall trough here!");
129#else
130  cerr << "\nlibfac: ERROR: SqrFreeTest: we should never fall trough here!\n" 
131       << rcsid << errmsg << endl;
132#endif
133  return 0;
134}
135
136///////////////////////////////////////////////////////////////
137// A uni/multivariate SqrFree routine.Works for polynomials  //
138// which don\'t have a constant as the content.              //
139// Works for charcteristic 0 and q=p^m                       //
140// returns a list of polys each of sqrfree, but gcd(f_i,f_j) //
141// needs not to be 1 !!!!!                                   //
142///////////////////////////////////////////////////////////////
143static CFFList
144SqrFreed( const CanonicalForm & r ){
145  CanonicalForm h, g, f = r;
146  CFFList Outputlist;
147  int n = level(f);
148
149  DEBINCLEVEL(cout, "SqrFreed");
150  DEBOUTLN(cout, "Called with r= ", r);
151  if (getNumVars(f)==0 ) { // just a constant; return it
152    Outputlist= myappend(Outputlist,CFFactor(f,1));
153    return Outputlist ;
154  }
155
156// We look if we do have a content; if so, SqrFreed the content
157// and continue computations with pp(f)
158  for (int k=1; k<=n; k++) {
159    g = swapvar(f,k,n); g = content(g);
160    if ( ! (g.isOne() || (-g).isOne() || degree(g)==0 )) {
161      g = swapvar(g,k,n);
162      DEBOUTLN(cout, "We have a content: ", g);
163      Outputlist = myUnion(InternalSqrFree(g),Outputlist); // should we add a
164                                                // SqrFreeTest(g) first ?
165      DEBOUTLN(cout, "Outputlist is now: ", Outputlist);
166      f /=g;
167      DEBOUTLN(cout, "f is now: ", f);
168    }
169  }
170
171// Now f is primitive; Let`s look if f is univariate
172  if ( f.isUnivariate() ) {
173    DEBOUTLN(cout, "f is univariate: ", f);
174    g = content(g);
175    if ( ! (g.isOne() || (-g).isOne() ) ){
176      Outputlist= myappend(Outputlist,CFFactor(g,1)) ;
177      f /= g;
178    }
179    Outputlist = Union(sqrFree(f),Outputlist) ; 
180    DEBOUTLN(cout, "Outputlist after univ. sqrFree(f) = ", Outputlist);
181    DEBDECLEVEL(cout, "SqrFreed");
182    return Outputlist ;
183  }
184
185// Linear?
186  if ( totaldegree(f) <= 1 ) {
187    Outputlist= myappend(Outputlist,CFFactor(f,1)) ;
188    DEBDECLEVEL(cout, "SqrFreed");
189    return Outputlist ;
190  }
191
192// is it Pth root?
193  n=level(f); // maybe less indeterminants
194  g= f.deriv();
195  if ( getCharacteristic() > 0 && g.isZero() ){  // Pth roots only apply to char > 0
196    for (int k=1; k<=n; k++) {
197      g=swapvar(f,k,n) ; g = g.deriv();
198      if ( ! g.isZero() ){ // can`t be Pth root
199        CFFList Outputlist2= SqrFreed(swapvar(f,k,n)); 
200        for (CFFListIterator inter=Outputlist2; inter.hasItem(); inter++){
201          Outputlist= myappend(Outputlist, CFFactor(swapvar(inter.getItem().factor(),k,n), inter.getItem().exp()));
202        }
203        return Outputlist;
204      }
205      else 
206        if ( k==n ) { // really is Pth power
207          DEBOUTLN(cout, "f is a p'th root: ", f);
208          CFMap m;
209          g = compress(f,m);
210          f = m(PthRoot(g));
211          DEBOUTLN(cout, "  that is       : ", f);
212          // now : Outputlist union ( SqrFreed(f) )^getCharacteristic()
213          Outputlist=myUnion(Powerup(InternalSqrFree(f),getCharacteristic()),Outputlist);
214          DEBDECLEVEL(cout, "SqrFreed");
215          return Outputlist ;
216        }
217    }
218  }
219  g = f.deriv();
220  DEBOUTLN(cout, "calculating mygcd of ", f);
221  DEBOUTLN(cout, "               and ", g);
222  h = mygcd(f,pp(g));  h /= lc(h);
223  DEBOUTLN(cout,"mygcd(f,g)= ",h);
224  if ( (h.isOne()) || ( h==f) || ((-h).isOne()) || getNumVars(h)==0 ) { // no common factor
225    Outputlist= myappend(Outputlist,CFFactor(f,1)) ;
226    DEBOUTLN(cout, "Outputlist= ", Outputlist);
227    DEBDECLEVEL(cout, "SqrFreed");
228    return Outputlist ;
229  }
230  else { // we can split into two nontrivial pieces
231    f /= h; // Now we have split the poly into f and h
232    g = lc(f);
233    if ( g != f.genOne() && getNumVars(g) == 0 ){
234       Outputlist= myappend(Outputlist,CFFactor(g,1)) ;
235       f /= g;
236    }
237    DEBOUTLN(cout, "Split into f= ", f);
238    DEBOUTLN(cout, "       and h= ", h);
239    // For char > 0 the polys f and h can have Pth roots; so we need a test
240    // Test is disabled because timing is the same
241//    if ( SqrFreeTest(f,0) )
242//      Outputlist= myappend(Outputlist,CFFactor(f,1)) ;
243//    else
244    Outputlist=myUnion(Outputlist, InternalSqrFree(f));
245//    if ( SqrFreeTest(h,0) )
246//      Outputlist= myappend(Outputlist,CFFactor(h,1)) ;
247//    else
248    Outputlist=myUnion(Outputlist,InternalSqrFree(h));
249    DEBOUTLN(cout, "Returning list ", Outputlist);
250    DEBDECLEVEL(cout, "SqrFreed");
251    return Outputlist ;
252  }
253#ifdef HAVE_SINGULAR
254  extern void WerrorS(char *);
255  WerrorS("libfac: ERROR: SqrFreed: we should never fall trough here!");
256#else
257  cerr << "\nlibfac: ERROR: SqrFreed: we should never fall trough here!\n" 
258       << rcsid << errmsg << endl;
259#endif
260  DEBDECLEVEL(cout, "SqrFreed");
261  return Outputlist; // for safety purpose
262}
263
264///////////////////////////////////////////////////////////////
265// The user front-end for the SqrFreed routine.              //
266// Input can have a constant as content                      //
267///////////////////////////////////////////////////////////////
268CFFList
269InternalSqrFree( const CanonicalForm & r ){
270  CanonicalForm g=icontent(r), f = r;
271  CFFList Outputlist, Outputlist2;
272
273  DEBINCLEVEL(cout, "InternalSqrFree");
274  DEBOUTMSG(cout, rcsid);
275  DEBOUTLN(cout,"Called with f= ", f);
276
277  // Take care of stupid users giving us constants
278  if ( getNumVars(f) == 0 ) { // a constant ; Exp==1 even if f==0
279      Outputlist= myappend(Outputlist,CFFactor(f,1));
280  }
281  else{
282      // Now we are sure: we have a nonconstant polynomial
283      g = lc(f);
284      while ( getNumVars(g) != 0 ) g=content(g);
285      if ( ! g.isOne() ) Outputlist= myappend(Outputlist,CFFactor(g,1)) ;
286      f /= g;
287      if ( getNumVars(f) != 0 ) // a real polynomial
288        Outputlist=myUnion(SqrFreed(f),Outputlist) ;
289  }
290  DEBOUTLN(cout,"Outputlist = ", Outputlist);
291  for ( CFFListIterator i=Outputlist; i.hasItem(); i++ )
292    if ( getNumVars(i.getItem().factor()) > 0 )
293      Outputlist2.append(i.getItem());
294
295  DEBOUTLN(cout,"Outputlist2 = ", Outputlist2);
296  DEBDECLEVEL(cout, "InternalSqrFree");
297  return Outputlist2 ;
298}
299
300CFFList
301SqrFree(const CanonicalForm & r ){
302  CFFList outputlist, sqrfreelist = InternalSqrFree(r);
303  CFFListIterator i;
304  CanonicalForm elem;
305  int n=totaldegree(r);
306
307  DEBINCLEVEL(cout, "SqrFree");
308
309  if ( sqrfreelist.length() < 2 ){
310    DEBDECLEVEL(cout, "SqrFree");
311    return sqrfreelist;
312  }
313  for ( int j=1; j<=n; j++ ){
314    elem =1;
315    for ( i = sqrfreelist; i.hasItem() ; i++ ){
316      if ( i.getItem().exp() == j ) elem *= i.getItem().factor();
317    }
318    if ( elem != 1 ) outputlist.append(CFFactor(elem,j));
319  }
320  elem=1;
321  for ( i=outputlist; i.hasItem(); i++ )
322    if ( getNumVars(i.getItem().factor()) > 0 )
323      elem*= power(i.getItem().factor(),i.getItem().exp());
324  elem= r/elem;
325  outputlist.insert(CFFactor(elem,1));
326
327  DEBOUTLN(cout, "SqrFree returns list:", outputlist);
328  DEBDECLEVEL(cout, "SqrFree");
329  return outputlist;
330}
331
332/*
333$Log: not supported by cvs2svn $
334Revision 1.4  1997/04/25 22:19:46  michael
335changed cerr and cout messages for use with Singular
336Version for libfac-0.2.1
337
338*/
Note: See TracBrowser for help on using the repository browser.