source: git/factory/cfGcdUtil.cc @ a723558

spielwiese
Last change on this file since a723558 was a723558, checked in by Hans Schoenemann <hannes@…>, 4 years ago
fix: typo
  • Property mode set to 100644
File size: 8.0 KB
Line 
1#include "config.h"
2
3#include "canonicalform.h"
4#include "cf_factory.h"
5#include "cf_reval.h"
6#include "cf_util.h"
7#include "cf_iter.h"
8#include "gfops.h"
9#include "cf_map_ext.h"
10#include "templates/ftmpl_functions.h"
11
12#ifdef HAVE_NTL
13#include "NTLconvert.h"
14#endif
15
16#ifdef HAVE_FLINT
17#include "FLINTconvert.h"
18#endif
19
20/// Coprimality Check. f and g are assumed to have the same level. If swap is
21/// true, the main variables of f and g are swapped with Variable(1). If the
22/// result is false, d is set to the degree of the gcd of f and g evaluated at a
23/// random point in K^n-1. This gcd is a gcd of univariate polynomials.
24#ifdef HAVE_NTL // mapPrimElem
25bool
26gcd_test_one ( const CanonicalForm & f, const CanonicalForm & g, bool swap, int & d )
27{
28    d= 0;
29    int count = 0;
30    // assume polys have same level;
31
32    Variable v= Variable (1);
33    bool algExtension= (hasFirstAlgVar (f, v) || hasFirstAlgVar (g, v));
34    CanonicalForm lcf, lcg;
35    if ( swap )
36    {
37        lcf = swapvar( LC( f ), Variable(1), f.mvar() );
38        lcg = swapvar( LC( g ), Variable(1), f.mvar() );
39    }
40    else
41    {
42        lcf = LC( f, Variable(1) );
43        lcg = LC( g, Variable(1) );
44    }
45
46    CanonicalForm F, G;
47    if ( swap )
48    {
49        F=swapvar( f, Variable(1), f.mvar() );
50        G=swapvar( g, Variable(1), g.mvar() );
51    }
52    else
53    {
54        F = f;
55        G = g;
56    }
57
58    #define TEST_ONE_MAX 50
59    int p= getCharacteristic();
60    bool passToGF= false;
61    int k= 1;
62    bool extOfExt= false;
63    Variable v3;
64    if (p > 0 && p < TEST_ONE_MAX && CFFactory::gettype() != GaloisFieldDomain && !algExtension)
65    {
66      if (p == 2)
67        setCharacteristic (2, 6, 'Z');
68      else if (p == 3)
69        setCharacteristic (3, 4, 'Z');
70      else if (p == 5 || p == 7)
71        setCharacteristic (p, 3, 'Z');
72      else
73        setCharacteristic (p, 2, 'Z');
74      passToGF= true;
75    }
76    else if (p > 0 && CFFactory::gettype() == GaloisFieldDomain && ipower (p , getGFDegree()) < TEST_ONE_MAX)
77    {
78      k= getGFDegree();
79      if (ipower (p, 2*k) > TEST_ONE_MAX)
80        setCharacteristic (p, 2*k, gf_name);
81      else
82        setCharacteristic (p, 3*k, gf_name);
83      F= GFMapUp (F, k);
84      G= GFMapUp (G, k);
85      lcf= GFMapUp (lcf, k);
86      lcg= GFMapUp (lcg, k);
87    }
88    else if (p > 0 && p < TEST_ONE_MAX && algExtension)
89    {
90#if defined(HAVE_NTL) || defined(HAVE_FLINT)
91      int d= degree (getMipo (v));
92      CFList source, dest;
93      Variable v2;
94      CanonicalForm primElem, imPrimElem;
95      #if defined(HAVE_NTL) && !defined(HAVE_FLINT)
96      if (fac_NTL_char != p)
97      {
98        fac_NTL_char= p;
99        zz_p::init (p);
100      }
101      #endif
102      if (p == 2 && d < 6)
103      {
104        bool primFail= false;
105        Variable vBuf;
106        primElem= primitiveElement (v, vBuf, primFail);
107        ASSERT (!primFail, "failure in integer factorizer");
108        if (d < 3)
109        {
110          #ifdef HAVE_FLINT
111          nmod_poly_t Irredpoly;
112          nmod_poly_init(Irredpoly,p);
113          nmod_poly_randtest_monic_irreducible(Irredpoly, FLINTrandom, 3*d+1);
114          CanonicalForm newMipo=convertnmod_poly_t2FacCF(Irredpoly,Variable(1));
115          nmod_poly_clear(Irredpoly);
116          #elif defined(HAVE_NTL)
117          zz_pX NTLIrredpoly;
118          BuildIrred (NTLIrredpoly, d*3);
119          CanonicalForm newMipo= convertNTLzzpX2CF (NTLIrredpoly, Variable (1));
120          #endif
121          v2= rootOf (newMipo);
122        }
123        else
124        {
125          #ifdef HAVE_FLINT
126          nmod_poly_t Irredpoly;
127          nmod_poly_init(Irredpoly,p);
128          nmod_poly_randtest_monic_irreducible(Irredpoly, FLINTrandom, 3*d+1);
129          CanonicalForm newMipo=convertnmod_poly_t2FacCF(Irredpoly,Variable(1));
130          nmod_poly_clear(Irredpoly);
131          #elif defined(HAVE_NTL)
132          zz_pX NTLIrredpoly;
133          BuildIrred (NTLIrredpoly, d*2);
134          CanonicalForm newMipo= convertNTLzzpX2CF (NTLIrredpoly, Variable (1));
135          #endif
136          v2= rootOf (newMipo);
137        }
138        imPrimElem= mapPrimElem (primElem, v, v2);
139        extOfExt= true;
140      }
141      else if ((p == 3 && d < 4) || ((p == 5 || p == 7) && d < 3))
142      {
143        bool primFail= false;
144        Variable vBuf;
145        primElem= primitiveElement (v, vBuf, primFail);
146        ASSERT (!primFail, "failure in integer factorizer");
147        #ifdef HAVE_FLINT
148        nmod_poly_t Irredpoly;
149        nmod_poly_init(Irredpoly,p);
150        nmod_poly_randtest_monic_irreducible(Irredpoly, FLINTrandom, 2*d+1);
151        CanonicalForm newMipo=convertnmod_poly_t2FacCF(Irredpoly,Variable(1));
152        nmod_poly_clear(Irredpoly);
153        #elif defined(HAVE_NTL)
154        zz_pX NTLIrredpoly;
155        BuildIrred (NTLIrredpoly, d*2);
156        CanonicalForm newMipo= convertNTLzzpX2CF (NTLIrredpoly, Variable (1));
157        #endif
158        v2= rootOf (newMipo);
159        imPrimElem= mapPrimElem (primElem, v, v2);
160        extOfExt= true;
161      }
162      if (extOfExt)
163      {
164        v3= v;
165        F= mapUp (F, v, v2, primElem, imPrimElem, source, dest);
166        G= mapUp (G, v, v2, primElem, imPrimElem, source, dest);
167        lcf= mapUp (lcf, v, v2, primElem, imPrimElem, source, dest);
168        lcg= mapUp (lcg, v, v2, primElem, imPrimElem, source, dest);
169        v= v2;
170      }
171#endif
172    }
173
174    CFRandom * sample;
175    if ((!algExtension && p > 0) || p == 0)
176      sample = CFRandomFactory::generate();
177    else
178      sample = AlgExtRandomF (v).clone();
179
180    REvaluation e( 2, tmax( f.level(), g.level() ), *sample );
181    delete sample;
182
183    if (passToGF)
184    {
185      lcf= lcf.mapinto();
186      lcg= lcg.mapinto();
187    }
188
189    CanonicalForm eval1, eval2;
190    if (passToGF)
191    {
192      eval1= e (lcf);
193      eval2= e (lcg);
194    }
195    else
196    {
197      eval1= e (lcf);
198      eval2= e (lcg);
199    }
200
201    while ( ( eval1.isZero() || eval2.isZero() ) && count < TEST_ONE_MAX )
202    {
203        e.nextpoint();
204        count++;
205        eval1= e (lcf);
206        eval2= e (lcg);
207    }
208    if ( count >= TEST_ONE_MAX )
209    {
210        if (passToGF)
211          setCharacteristic (p);
212        if (k > 1)
213          setCharacteristic (p, k, gf_name);
214        if (extOfExt)
215          prune1 (v3);
216        return false;
217    }
218
219
220    if (passToGF)
221    {
222      F= F.mapinto();
223      G= G.mapinto();
224      eval1= e (F);
225      eval2= e (G);
226    }
227    else
228    {
229      eval1= e (F);
230      eval2= e (G);
231    }
232
233    CanonicalForm c= gcd (eval1, eval2);
234    d= c.degree();
235    bool result= d < 1;
236    if (d < 0)
237      d= 0;
238
239    if (passToGF)
240      setCharacteristic (p);
241    if (k > 1)
242      setCharacteristic (p, k, gf_name);
243    if (extOfExt)
244      prune1 (v3);
245    return result;
246}
247#endif
248
249/**
250 * same as balance_p ( const CanonicalForm & f, const CanonicalForm & q )
251 * but qh= q/2 is provided, too.
252**/
253CanonicalForm
254balance_p ( const CanonicalForm & f, const CanonicalForm & q, const CanonicalForm & qh )
255{
256    Variable x = f.mvar();
257    CanonicalForm result = 0;
258    CanonicalForm c;
259    CFIterator i;
260    for ( i = f; i.hasTerms(); i++ )
261    {
262        c = i.coeff();
263        if ( c.inCoeffDomain())
264        {
265          if ( c > qh )
266            result += power( x, i.exp() ) * (c - q);
267          else
268            result += power( x, i.exp() ) * c;
269        }
270        else
271          result += power( x, i.exp() ) * balance_p(c,q,qh);
272    }
273    return result;
274}
275
276/** static CanonicalForm balance_p ( const CanonicalForm & f, const CanonicalForm & q )
277 *
278 * balance_p() - map f from positive to symmetric representation
279 *   mod q.
280 *
281 * This makes sense for polynomials over Z only.
282 * q should be an integer.
283 *
284**/
285CanonicalForm
286balance_p ( const CanonicalForm & f, const CanonicalForm & q )
287{
288    CanonicalForm qh = q / 2;
289    return balance_p (f, q, qh);
290}
291
292
293/*static CanonicalForm
294balance ( const CanonicalForm & f, const CanonicalForm & q )
295{
296    Variable x = f.mvar();
297    CanonicalForm result = 0, qh = q / 2;
298    CanonicalForm c;
299    CFIterator i;
300    for ( i = f; i.hasTerms(); i++ ) {
301        c = mod( i.coeff(), q );
302        if ( c > qh )
303            result += power( x, i.exp() ) * (c - q);
304        else
305            result += power( x, i.exp() ) * c;
306    }
307    return result;
308}*/
Note: See TracBrowser for help on using the repository browser.