# source:git/factory/facFqSquarefree.cc@8c5f4b

spielwiese
Last change on this file since 8c5f4b was 8c5f4b, checked in by Martin Lee <martinlee84@…>, 9 years ago
deleted some unnecessary const
• Property mode set to 100644
File size: 7.8 KB
Line
1/*****************************************************************************\
2 * Computer Algebra System SINGULAR
3\*****************************************************************************/
4/** @file facFqSquarefree.cc
5 *
6 * This file provides functions for squarefrees factorizing over
7 * \f$F_{p} \f$ , \f$F_{p}(\alpha ) \f$ or GF, as decribed in "Factoring
8 * multivariate polynomials over a finite field" by L. Bernardin
9 *
10 * @author Martin Lee
11 *
12 **/
13/*****************************************************************************/
14
15
16#include "config.h"
17
18
19#include "canonicalform.h"
20
21#include "cf_iter.h"
22#include "cf_map.h"
23#include "cf_util.h"
24#include "cf_factory.h"
25#include "facFqSquarefree.h"
26
27#ifdef HAVE_NTL
28#include "NTLconvert.h"
29#endif
30
31#ifdef HAVE_FLINT
32#include "FLINTconvert.h"
33#endif
34
35static inline
36CanonicalForm
37pthRoot (const CanonicalForm & F, int q)
38{
39  CanonicalForm A= F;
40  int p= getCharacteristic ();
41  if (A.inCoeffDomain())
42  {
43    A= power (A, q/p);
44    return A;
45  }
46  else
47  {
48    CanonicalForm buf= 0;
49    for (CFIterator i= A; i.hasTerms(); i++)
50      buf= buf + power(A.mvar(), i.exp()/p)*pthRoot (i.coeff(), q);
51    return buf;
52  }
53}
54
55#ifdef HAVE_NTL
56CanonicalForm
57pthRoot (const CanonicalForm & F, const ZZ& q, const Variable& alpha)
58{
59  CanonicalForm A= F;
60  int p= getCharacteristic ();
61  if (A.inCoeffDomain())
62  {
63    zz_p::init (p);
64    zz_pX NTLMipo= convertFacCF2NTLzzpX (getMipo (alpha));
65    zz_pE::init (NTLMipo);
66    zz_pX NTLA= convertFacCF2NTLzzpX (A);
67    zz_pE NTLA2= to_zz_pE (NTLA);
68    power (NTLA2, NTLA2, q/p);
69    A= convertNTLzzpE2CF (NTLA2, alpha);
70    return A;
71  }
72  else
73  {
74    CanonicalForm buf= 0;
75    for (CFIterator i= A; i.hasTerms(); i++)
76      buf= buf + power(A.mvar(), i.exp()/p)*pthRoot (i.coeff(), q, alpha);
77    return buf;
78  }
79}
80#endif
81
82#if (HAVE_FLINT && __FLINT_VERSION_MINOR >= 4)
83CanonicalForm
84pthRoot (const CanonicalForm & F, const fmpz_t& q, const Variable& alpha)
85{
86  CanonicalForm A= F;
87  int p= getCharacteristic ();
88  if (A.inCoeffDomain())
89  {
90    nmod_poly_t FLINTmipo;
91    fq_nmod_ctx_t fq_con;
92    fmpz_t qp;
93    fq_nmod_t FLINTA;
94
95    nmod_poly_init (FLINTmipo, p);
96    convertFacCF2nmod_poly_t (FLINTmipo, getMipo (alpha));
97
98    fq_nmod_ctx_init_modulus (fq_con, FLINTmipo, "Z");
99
100    fq_nmod_init2 (FLINTA, fq_con);
101
102    convertFacCF2Fq_nmod_t (FLINTA, A, fq_con);
103
104    fmpz_init_set (qp, q);
105    fmpz_divexact_si (qp, qp, p);
106
107    fq_nmod_pow (FLINTA, FLINTA, qp, fq_con);
108    A= convertFq_nmod_t2FacCF (FLINTA, alpha);
109
110    fmpz_clear (qp);
111    nmod_poly_clear (FLINTmipo);
112    fq_nmod_clear (FLINTA, fq_con);
113    fq_nmod_ctx_clear (fq_con);
114    return A;
115  }
116  else
117  {
118    CanonicalForm buf= 0;
119    for (CFIterator i= A; i.hasTerms(); i++)
120      buf= buf + power(A.mvar(), i.exp()/p)*pthRoot (i.coeff(), q, alpha);
121    return buf;
122  }
123}
124#endif
125
126CanonicalForm
127maxpthRoot (const CanonicalForm & F, int q, int& l)
128{
129  CanonicalForm result= F;
130  bool derivZero= true;
131  l= 0;
132  while (derivZero)
133  {
134    for (int i= 1; i <= result.level(); i++)
135    {
136      if (!deriv (result, Variable (i)).isZero())
137      {
138        derivZero= false;
139        break;
140      }
141    }
142    if (!derivZero)
143      break;
144    result= pthRoot (result, q);
145    l++;
146  }
147  return result;
148}
149
150static inline
151CFFList
152sqrfPosDer (const CanonicalForm & F, const Variable & x,
153            CanonicalForm & c)
154{
155  CanonicalForm b= deriv (F, x);
156  c= gcd (F, b);
157  CanonicalForm w= F/c;
158  CanonicalForm v= b/c;
159  CanonicalForm u= v - deriv (w, x);
160  int j= 1;
161  int p= getCharacteristic();
162  CanonicalForm g;
163  CFFList result;
164  while (j < p - 1 && degree(u) >= 0)
165  {
166    g= gcd (w, u);
167    if (!g.inCoeffDomain())
168      result.append (CFFactor (g, j));
169    w= w/g;
170    c= c/w;
171    v= u/g;
172    u= v - deriv (w, x);
173    j++;
174  }
175  if (!w.inCoeffDomain())
176    result.append (CFFactor (w, j));
177  return result;
178}
179
180CFFList
181squarefreeFactorization (const CanonicalForm & F, const Variable & alpha)
182{
183  int p= getCharacteristic();
184  CanonicalForm A= F;
185  CFMap M;
186  A= compress (A, M);
187  Variable x= A.mvar();
188  int l= x.level();
189  int k;
190  if (CFFactory::gettype() == GaloisFieldDomain)
191    k= getGFDegree();
192  else if (alpha.level() != 1)
193    k= degree (getMipo (alpha));
194  else
195    k= 1;
196  Variable buf, buf2;
197  CanonicalForm tmp;
198
199  CFFList tmp1, tmp2;
200  bool found;
201  for (int i= l; i > 0; i--)
202  {
203    buf= Variable (i);
204    if (degree (deriv (A, buf)) >= 0)
205    {
206      tmp1= sqrfPosDer (A, buf, tmp);
207      A= tmp;
208      for (CFFListIterator j= tmp1; j.hasItem(); j++)
209      {
210        found= false;
211        CFFListIterator k= tmp2;
212        if (!k.hasItem() && !j.getItem().factor().inCoeffDomain()) tmp2.append (j.getItem());
213        else
214        {
215          for (; k.hasItem(); k++)
216          {
217            if (k.getItem().exp() == j.getItem().exp())
218            {
219              k.getItem()= CFFactor (k.getItem().factor()*j.getItem().factor(),
220                                     j.getItem().exp());
221              found= true;
222            }
223          }
224          if (found == false && !j.getItem().factor().inCoeffDomain())
225            tmp2.append(j.getItem());
226        }
227      }
228    }
229  }
230
231  bool degcheck= false;;
232  for (int i= l; i > 0; i--)
233  if (degree (A, Variable(i)) >= p)
234    degcheck= true;
235
236  if (degcheck == false && tmp1.isEmpty() && tmp2.isEmpty())
237    return CFFList (CFFactor (F/Lc(F), 1));
238
239  CanonicalForm buffer;
240#if defined(HAVE_NTL) || (HAVE_FLINT && __FLINT_VERSION_MINOR >= 4)
241  if (alpha.level() == 1)
242#endif
243    buffer= pthRoot (A, ipower (p, k));
244#if (HAVE_FLINT && __FLINT_VERSION_MINOR >= 4)
245  else
246  {
247    fmpz_t qq;
248    fmpz_init_set_ui (qq, p);
249    fmpz_pow_ui (qq, qq, k);
250    buffer= pthRoot (A, qq, alpha);
251    fmpz_clear (qq);
252  }
253#elif defined(HAVE_NTL)
254  else
255  {
256    ZZ q;
257    power (q, p, k);
258    buffer= pthRoot (A, q, alpha);
259  }
260#endif
261
262  tmp1= squarefreeFactorization (buffer, alpha);
263
264  CFFList result;
265  buf= alpha;
266  for (CFFListIterator i= tmp2; i.hasItem(); i++)
267  {
268    for (CFFListIterator j= tmp1; j.hasItem(); j++)
269    {
270      tmp= gcd (i.getItem().factor(), j.getItem().factor());
271      i.getItem()= CFFactor (i.getItem().factor()/tmp, i.getItem().exp());
272      j.getItem()= CFFactor (j.getItem().factor()/tmp, j.getItem().exp());
273      if (!tmp.inCoeffDomain())
274      {
275        tmp= M (tmp);
276        result.append (CFFactor (tmp/Lc(tmp),
277                       j.getItem().exp()*p + i.getItem().exp()));
278      }
279    }
280  }
281  for (CFFListIterator i= tmp2; i.hasItem(); i++)
282  {
283    if (!i.getItem().factor().inCoeffDomain())
284    {
285      tmp= M (i.getItem().factor());
286      result.append (CFFactor (tmp/Lc(tmp), i.getItem().exp()));
287    }
288  }
289  for (CFFListIterator j= tmp1; j.hasItem(); j++)
290  {
291    if (!j.getItem().factor().inCoeffDomain())
292    {
293      tmp= M (j.getItem().factor());
294      result.append (CFFactor (tmp/Lc(tmp), j.getItem().exp()*p));
295    }
296  }
297  return result;
298}
299
300CanonicalForm
301sqrfPart (const CanonicalForm& F, CanonicalForm& pthPower,
302          const Variable& alpha)
303{
304  if (F.inCoeffDomain())
305  {
306    pthPower= 1;
307    return F;
308  }
309  CFMap M;
310  CanonicalForm A= compress (F, M);
311  Variable vBuf= alpha;
312  CanonicalForm w, v, b;
313  pthPower= 1;
314  CanonicalForm result;
315  int i= 1;
316  bool allZero= true;
317  for (; i <= A.level(); i++)
318  {
319    if (!deriv (A, Variable (i)).isZero())
320    {
321      allZero= false;
322      break;
323    }
324  }
325  if (allZero)
326  {
327    pthPower= F;
328    return 1;
329  }
330  w= gcd (A, deriv (A, Variable (i)));
331
332  b= A/w;
333  result= b;
334  if (degree (w) < 1)
335    return M (result);
336  i++;
337  for (; i <= A.level(); i++)
338  {
339    if (!deriv (w, Variable (i)).isZero())
340    {
341      b= w;
342      w= gcd (w, deriv (w, Variable (i)));
343      b /= w;
344      if (degree (b) < 1)
345        break;
346      CanonicalForm g;
347      g= gcd (b, result);
348      if (degree (g) > 0)
349        result *= b/g;
350      if (degree (g) <= 0)
351        result *= b;
352    }
353  }
354  result= M (result);
355  return result;
356}
357
Note: See TracBrowser for help on using the repository browser.