source: git/factory/libfac/factor/Truefactor.cc @ 741bb9

spielwiese
Last change on this file since 741bb9 was 16055bd, checked in by Martin Lee <martinlee84@…>, 12 years ago
chg: move libfac into factory
  • Property mode set to 100644
File size: 12.7 KB
Line 
1///////////////////////////////////////////////////////////////////////////////
2// emacs edit mode for this file is -*- C++ -*-
3///////////////////////////////////////////////////////////////////////////////
4// Factory - Includes
5#include <factory.h>
6#ifndef NOSTREAMIO
7#ifdef HAVE_IOSTREAM
8#include <iostream>
9#define OSTREAM std::ostream
10#define ISTREAM std::istream
11#define CERR std::cerr
12#define CIN std::cin
13#elif defined(HAVE_IOSTREAM_H)
14#include <iostream.h>
15#define OSTREAM ostream
16#define ISTREAM istream
17#define CERR cerr
18#define CIN cin
19#endif
20#endif
21// Factor - Includes
22#include "tmpl_inst.h"
23#include "helpstuff.h"
24// some CC's need this:
25#include "Truefactor.h"
26
27#ifdef TRUEFACTORDEBUG
28#  define DEBUGOUTPUT
29#else
30#  undef DEBUGOUTPUT
31#endif
32
33#include <libfac/factor/debug.h>
34#include "timing.h"
35
36int hasAlgVar(const CanonicalForm &f, const Variable &v)
37{
38  if (f.inBaseDomain()) return 0;
39  if (f.inCoeffDomain())
40  {
41    if (f.mvar()==v) return 1;
42    return hasAlgVar(f.LC(),v);
43  }
44  if (f.inPolyDomain())
45  {
46    if (hasAlgVar(f.LC(),v)) return 1;
47    for( CFIterator i=f; i.hasTerms(); i++)
48    {
49      if (hasAlgVar(i.coeff(),v)) return 1;
50    }
51  }
52  return 0;
53}
54
55int hasVar(const CanonicalForm &f, const Variable &v)
56{
57  if (f.inBaseDomain()) return 0;
58  if (f.inCoeffDomain())
59  {
60    if (f.mvar()==v) return 1;
61    return hasAlgVar(f.LC(),v);
62  }
63  if (f.inPolyDomain())
64  {
65    if (f.mvar()==v) return 1;
66    if (hasVar(f.LC(),v)) return 1;
67    for( CFIterator i=f; i.hasTerms(); i++)
68    {
69      if (hasVar(i.coeff(),v)) return 1;
70    }
71  }
72  return 0;
73}
74
75int hasAlgVar(const CanonicalForm &f)
76{
77  if (f.inBaseDomain()) return 0;
78  if (f.inCoeffDomain())
79  {
80    if (f.level()!=0)
81    {
82      //CERR << "hasAlgVar:" << f.mvar() <<"\n";
83      return 1;
84    }
85    return hasAlgVar(f.LC());
86  }
87  if (f.inPolyDomain())
88  {
89    if (hasAlgVar(f.LC())) return 1;
90    for( CFIterator i=f; i.hasTerms(); i++)
91    {
92      if (hasAlgVar(i.coeff())) return 1;
93    }
94  }
95  return 0;
96}
97
98///////////////////////////////////////////////////////////////
99// generate all different k-subsets of the set with n        //
100// elements and return them in returnlist.                   //
101///////////////////////////////////////////////////////////////
102static void
103combinat( int k, int n, List<IntList> & returnlist ){
104  ListIntList ListofLists;
105  IntList intermediate,intermediate2;
106  int value,j;
107
108  if ( k == 1 ){                     // k=1
109    for ( j=1; j<=n; j++)
110      returnlist.append( IntList(j) );
111  }
112  else{                              // k-1 --> k
113    combinat(k-1,n,returnlist);  // generate (k-1,n)
114    for ( ListIntListIterator l=returnlist; l.hasItem(); l++ ){
115      intermediate = l.getItem();
116      value = intermediate.getLast();
117      if ( value != n )
118        for ( j=value+1; j<=n; j++ ){
119          intermediate2 = intermediate; intermediate2.append(j);
120          ListofLists.append( intermediate2 );
121        }
122    }
123    returnlist = ListofLists;
124  }
125}
126
127///////////////////////////////////////////////////////////////
128// Return the CanonicalForm number nr in  Factorlist.        //
129///////////////////////////////////////////////////////////////
130static CanonicalForm
131getItemNr(int nr, const CFFList & Factorlist ){
132  ListIterator<CFFactor> i=Factorlist;
133  int Nr=nr;
134
135  for ( Nr=1; Nr<nr; Nr++ ) i++;
136  return i.getItem().factor();
137}
138
139///////////////////////////////////////////////////////////////
140// Generate all sets of m factors out of LiftedFactors list. //
141///////////////////////////////////////////////////////////////
142static CFFList
143combine( int m, const CFFList & LiftedFactors ){
144  CFFList result;
145  ListIntList CombinatList;
146  CanonicalForm intermediate;
147
148  combinat(m, LiftedFactors.length(), CombinatList);
149  for ( ListIntListIterator j=CombinatList ; j.hasItem(); j++ ){
150    intermediate=1;
151    for ( IntListIterator k=j.getItem(); k.hasItem(); k++ )
152      intermediate *= getItemNr(k.getItem(), LiftedFactors);
153    if (!hasAlgVar(intermediate))
154    result.append(CFFactor(intermediate,1));
155  }
156  return result;
157}
158
159///////////////////////////////////////////////////////////////
160// Remove element elem from the list L.                      //
161///////////////////////////////////////////////////////////////
162static CFFList
163Remove_from_List( const CFFList & L, const CanonicalForm & elem ){
164  CFFList Returnlist;
165
166  DEBOUTLN(CERR, "Remove_from_List called with L= ",L);
167  DEBOUTLN(CERR, "                     and  elem= ",elem);
168  for ( ListIterator<CFFactor> i = L ; i.hasItem(); i++)
169    if ( i.getItem().factor() != elem )
170      Returnlist.append( i.getItem() );
171
172  return Returnlist;
173}
174
175///////////////////////////////////////////////////////////////
176// Here we solve:          G= F mod ( P, S^h )               //
177///////////////////////////////////////////////////////////////
178static CanonicalForm
179Multmod_power( const CanonicalForm & F, const SFormList & Substituionlist, int h, int levelF){
180  CanonicalForm G;
181
182  G= change_poly(F, Substituionlist, 0);
183  G= mod_power(G, h, levelF);
184  G= change_poly(G, Substituionlist, 1);
185
186  return G;
187}
188
189///////////////////////////////////////////////////////////////
190// Determine the right degree for the list of combinations   //
191// of factors, i.e. delete any element from list CombL which //
192// degree in the main variable levelU exceeeds degU.         //
193///////////////////////////////////////////////////////////////
194static CFFList
195Rightdegree( const CFFList & CombL, int degU, int levelU ){
196  CFFList Returnlist;
197  CFFactor factor;
198
199  for ( ListIterator<CFFactor> i= CombL; i.hasItem(); i++ ){
200    factor= i.getItem();
201    if ( degree(factor.factor(), levelU) <= degU )
202      Returnlist.append(factor);
203  }
204
205  return Returnlist;
206}
207
208///////////////////////////////////////////////////////////////
209// We have properly lifted all the specialized factors. See  //
210// which one works.                                          //
211// We use the (modified) algorithm TRUEFACTORS given by      //
212// Paul S. Wang and Linda Preiss Rothschild:                 //
213// Factoring Multivariate Polynomials Over the Integers      //
214// Math. Comp. V29 Nr131 (July 1975) p. 935-950              //
215///////////////////////////////////////////////////////////////
216CFFList
217Truefactors( const CanonicalForm Ua, int levelU, const SFormList & SubstitutionList, const CFFList & PiList){
218  CanonicalForm U=Ua,a,b,Y;
219  CFFactor factor;
220  CFFList L,FAC,E_all;
221  int c,r = PiList.length(),degU, onemore,M, h = subvardegree(Ua,levelU) + 1;
222  ListIterator<CFFactor> i;
223
224  //CERR << "SubstitutionList="<< SubstitutionList<<"\n";
225// step 1: simply test the generated factors alone
226  for ( i= PiList; i.hasItem();i++){
227    factor = i.getItem();
228    //CanonicalForm test_f=change_poly(factor.factor(),SubstitutionList,0);
229    CanonicalForm test_f=factor.factor();
230    //CERR <<"f:" << factor.factor() << " -> test_f:"<<test_f <<"\n";
231    //CERR << "           1:" << change_poly(factor.factor(),SubstitutionList,1) <<"\n";
232    c= mydivremt(U,test_f,a,b);
233    if (  c  && b == U.genZero() && !hasAlgVar(test_f))
234    // factor.getFactor() divides U
235    {
236      //CERR << " teilt:" << test_f <<"\n";
237      U= a;
238      FAC.append(factor);
239    }
240    else{
241      //CERR << " teilt nicht:" << test_f <<"\n";
242      L.append(factor);
243    }
244  }
245  DEBOUTLN(CERR,"Truefactors: (step1) L  = ", L);
246  DEBOUTLN(CERR,"                     FAC= ", FAC);
247
248// step 2: Do we have to check combinations?
249  degU = L.length();
250  if ( degU == 0 ) // No elements: Return
251    return FAC;
252  else
253    if ( degU < 4 ){ // Less then four elements: no combinations possible
254      FAC.append( CFFactor(U,1) );
255      return FAC;
256    }
257    else {
258      M = 1; r = r - FAC.length(); degU = degree(U, levelU)/2;
259    }
260
261  DEBOUTLN(CERR,"Truefactors: (step2) M   = ", M);
262  DEBOUTLN(CERR,"                     r   = ", r);
263  DEBOUTLN(CERR,"                     degU= ", degU);
264
265// Now do the real work!
266// Test all the combinations of possible factors.
267
268  onemore=1;
269// steps 3 to 6
270  while (1)
271  {
272    // step 3 iff onemore == 1
273    if ( onemore ) M+= 1;  else onemore = 1;
274    // step 4
275    if ( U.isOne() ) break; // Return FAC
276    if ( ( r == 1 ) || ( M >= ( r-1 ) ) || ( M > degU ) )
277    {
278      FAC.append( CFFactor(U,1) );
279      break; // Return FAC union {U}
280    }
281    // step 5
282    E_all = combine(M, L); // generate all combinations of M elements from L
283    DEBOUTLN(CERR,"Truefactors: (step5) E_all= ", E_all);
284    // select combinations with the degree not to exceed degU:
285    E_all = Rightdegree( E_all, degU, levelU );
286    DEBOUTLN(CERR,"Truefactors: (step5) E_all(Rightdegree)= ", E_all);
287    if ( E_all.length() == 0 )
288    {
289      FAC.append( CFFactor(U,1) );
290      break; // Return FAC union {U}
291    }
292    for ( i=E_all; i.hasItem(); i++)
293    {
294      factor = i.getItem();
295      Y = Multmod_power( factor.factor(), SubstitutionList, h, levelU);
296      DEBOUTLN(CERR, "Truefactors: (step6) Testing Y  = ", Y);
297      c = mydivremt(U,Y,a,b);
298      //      if (  c  && b == U.genZero()) { // Y divides U
299      if ( c && b.isZero() )
300      {
301        DEBOUT(CERR,"Truefactors: (step6): ",Y );
302        DEBOUTLN(CERR, "  divides  ",U);
303        U = a;
304        FAC.append(Y); // Y is a real factor
305        onemore = 0;
306        degU = degree(U, levelU)/2; // new degU
307        // L = L \ {factor}
308        // Hier ist noch etwas faul; wir muessen (f=prod(f_i)) die f_i
309        // entfernen und nicht f!
310        L = Remove_from_List( L, factor.factor() );
311        r -= 1;
312        // delete from L any element with degree greater than degU
313        L = Rightdegree( L, degU, levelU );
314      }
315    }
316  }
317  return FAC;
318}
319
320///////////////////////////////////////////////////////////////
321// Check if poly f is in Fp (returns true) or in Fp(a)       //
322///////////////////////////////////////////////////////////////
323static bool is_in_Fp( const CanonicalForm & f )
324{
325  if ( f.inCoeffDomain() )
326    return f.inBaseDomain() ;
327  else
328  {
329    CFIterator i=f;
330    bool ok=true;
331    while ( ok && i.hasTerms() )
332    {
333      ok = is_in_Fp( i.coeff() );
334      i++ ;
335    }
336    return ok;
337  }
338}
339
340///////////////////////////////////////////////////////////////
341// We have factors which possibly lie in an extension of the //
342// base field. If one of these is not over the base field,   //
343// find its norm by (the theoretically handicapped method    //
344// of) multiplying by other elements.                        //
345///////////////////////////////////////////////////////////////
346CFFList TakeNorms(const CFFList & PiList)
347{
348  CFFList CopyPossibleFactors, PossibleFactors, TrueFactors;
349  CFFListIterator i;
350  CFFactor Factor;
351  CanonicalForm intermediate;
352  ListIntList CombinatList;
353  ListIntListIterator j;
354  IntListIterator k;
355
356  // First check if the factors in PiList already lie in Fp-Domain
357  for ( i=PiList; i.hasItem(); i++ )
358  {
359    Factor = i.getItem();
360    if ( is_in_Fp( Factor.factor() ) )
361      TrueFactors.append(Factor);
362    else
363      PossibleFactors.append(Factor);
364  }
365  // Now we have to check if combinations of the remaining factors
366  // (now in PossibleFactors) do lie in Fp-Domain
367  if ( PossibleFactors.length() > 0 ) // there are (at least two) items
368  {
369    int n=2;
370    if ( PossibleFactors.length() < n )  // a little check
371    {
372      factoryError("libfac: ERROR: TakeNorms less then two items remaining!");
373    }
374    while ( n < PossibleFactors.length() )
375    {
376      // generate all combinations of n elements
377      combinat(n, PossibleFactors.length(), CombinatList);
378      for ( j=CombinatList ; j.hasItem(); j++ )
379      {
380        intermediate=1;
381        for ( k=j.getItem(); k.hasItem(); k++ )
382          intermediate *= getItemNr( k.getItem(), PossibleFactors );
383        if ( is_in_Fp( intermediate ) )
384        {
385          TrueFactors.append(intermediate); // found a true factor
386          CopyPossibleFactors=PossibleFactors; // save list
387          for ( k=j.getItem(); k.hasItem(); k++ )
388            //remove combined factors from PossibleFactors
389            PossibleFactors=Remove_from_List(PossibleFactors,
390                                getItemNr( k.getItem(), CopyPossibleFactors ));
391          n-=1; // look for the same number of combined factors:
392          break;
393        }
394        else
395        {
396          //CERR << "Schade!" << "\n";
397        }
398        DEBOUT(CERR, "Truefactor: Combined ", n);
399        DEBOUTLN(CERR, " factors to: ", intermediate);
400      }
401      n += 1;
402    }
403  // All remaining factors in PossibleFactors multiplied
404  // should lie in Fp domain
405    if ( PossibleFactors.length() >=1 )
406    {
407      for ( i=PossibleFactors; i.hasItem(); i++ )
408        intermediate *= i.getItem().factor();
409      // a last check:
410      if ( is_in_Fp(intermediate) )
411      {
412        TrueFactors.append(CFFactor(intermediate,1));
413      }
414      else
415      {
416        factoryError("libfac: TakeNorms: somethings wrong with remaining factors!");
417      }
418    }
419  }
420  return TrueFactors;
421}
422
423////////////////////////////////////////////////////////////
Note: See TracBrowser for help on using the repository browser.