source: git/factory/facFqBivarUtil.cc @ 9614bb

spielwiese
Last change on this file since 9614bb was 9614bb, checked in by Martin Lee <martinlee84@…>, 11 years ago
chg: some optimizations to computation of logarithmic derivatives
  • Property mode set to 100644
File size: 21.4 KB
RevLine 
[7bf145]1/*****************************************************************************\
[806c18]2 * Computer Algebra System SINGULAR
[7bf145]3\*****************************************************************************/
4/** @file facFqBivarUtil.cc
[806c18]5 *
6 * This file provides utility functions for bivariate factorization
[7bf145]7 *
8 * @author Martin Lee
9 *
10 **/
11/*****************************************************************************/
12
[e4fe2b]13#include "config.h"
[7bf145]14
15#include "cf_map.h"
[81d96c]16#include "algext.h"
[7bf145]17#include "cf_map_ext.h"
[6db552]18#include "templates/ftmpl_functions.h"
[7bf145]19#include "ExtensionInfo.h"
[fecc08]20#include "cf_algorithm.h"
21#include "cf_factory.h"
[07718c3]22#include "cf_util.h"
[fecc08]23#include "imm.h"
24#include "cf_iter.h"
[7bf145]25#include "facFqBivarUtil.h"
[f876a66]26#include "cfNewtonPolygon.h"
27#include "facHensel.h"
[0e2e23]28#include "facMul.h"
[7bf145]29
30
[806c18]31void append (CFList& factors1, const CFList& factors2)
[7bf145]32{
[806c18]33  for (CFListIterator i= factors2; i.hasItem(); i++)
[7bf145]34  {
35    if (!i.getItem().inCoeffDomain())
36      factors1.append (i.getItem());
37  }
38  return;
39}
40
[806c18]41void decompress (CFList& factors, const CFMap& N)
[7bf145]42{
[806c18]43  for (CFListIterator i= factors; i.hasItem(); i++)
[7bf145]44    i.getItem()= N (i.getItem());
[f876a66]45}
46
47void decompress (CFFList& factors, const CFMap& N)
48{
49  for (CFFListIterator i= factors; i.hasItem(); i++)
50    i.getItem()= CFFactor (N (i.getItem().factor()), i.getItem().exp());
[7bf145]51}
52
[806c18]53void appendSwapDecompress (CFList& factors1, const CFList& factors2,
54                           const CFList& factors3, const bool swap1,
55                           const bool swap2, const CFMap& N)
[7bf145]56{
57  Variable x= Variable (1);
58  Variable y= Variable (2);
[806c18]59  for (CFListIterator i= factors1; i.hasItem(); i++)
[7bf145]60  {
[806c18]61    if (swap1)
[7bf145]62    {
[806c18]63      if (!swap2)
[7bf145]64        i.getItem()= swapvar (i.getItem(), x, y);
65    }
[806c18]66    else
[7bf145]67    {
[806c18]68      if (swap2)
[7bf145]69        i.getItem()= swapvar (i.getItem(), y, x);
70    }
71    i.getItem()= N (i.getItem());
72  }
[806c18]73  for (CFListIterator i= factors2; i.hasItem(); i++)
[7bf145]74    factors1.append (N (i.getItem()));
[806c18]75  for (CFListIterator i= factors3; i.hasItem(); i++)
[7bf145]76    factors1.append (N (i.getItem()));
77  return;
78}
79
[806c18]80void swapDecompress (CFList& factors, const bool swap, const CFMap& N)
[7bf145]81{
82  Variable x= Variable (1);
83  Variable y= Variable (2);
[806c18]84  for (CFListIterator i= factors; i.hasItem(); i++)
[7bf145]85  {
[806c18]86    if (swap)
[7bf145]87      i.getItem()= swapvar (i.getItem(), x, y);
88    i.getItem()= N (i.getItem());
89  }
90  return;
91}
92
[806c18]93static inline
94bool GFInExtensionHelper (const CanonicalForm& F, const int number)
[7bf145]95{
96  if (F.isOne()) return false;
97  InternalCF* buf;
98  int exp;
99  bool result= false;
[806c18]100  if (F.inBaseDomain())
[7bf145]101  {
102    buf= F.getval();
103    exp= imm2int (buf);
[806c18]104    if (exp%number != 0)
[7bf145]105      return true;
[806c18]106    else
[7bf145]107      return result;
108  }
[806c18]109  else
[7bf145]110  {
[806c18]111    for (CFIterator i= F; i.hasTerms(); i++)
[7bf145]112    {
113      result= GFInExtensionHelper (i.coeff(), number);
114      if (result == true)
115        return result;
116    }
117  }
118  return result;
119}
120
[806c18]121static inline
[f876a66]122bool FqInExtensionHelper (const CanonicalForm& F, const CanonicalForm& gamma,
123                          const CanonicalForm& delta, CFList& source,
124                          CFList& dest)
[7bf145]125{
126  bool result= false;
[806c18]127  if (F.inBaseDomain())
[7bf145]128    return result;
[806c18]129  else if (F.inCoeffDomain())
[7bf145]130  {
[806c18]131    if (!fdivides (gamma, F))
[7bf145]132      return true;
[806c18]133    else
[f876a66]134    {
135      int pos= findItem (source, F);
136      if (pos > 0)
137        return false;
138      Variable a;
139      hasFirstAlgVar (F, a);
140      int bound= ipower (getCharacteristic(), degree (getMipo (a)));
141      CanonicalForm buf= 1;
142      for (int i= 1; i < bound; i++)
143      {
144        buf *= gamma;
145        if (buf == F)
146        {
147          source.append (buf);
148          dest.append (power (delta, i));
149          return false;
150        }
151      }
152      return true;
153    }
[7bf145]154  }
[806c18]155  else
[7bf145]156  {
[806c18]157    for (CFIterator i= F; i.hasTerms(); i++)
[7bf145]158    {
[f876a66]159      result= FqInExtensionHelper (i.coeff(), gamma, delta, source, dest);
[7bf145]160      if (result == true)
[806c18]161        return result;
[7bf145]162    }
163  }
164  return result;
165}
166
167bool isInExtension (const CanonicalForm& F, const CanonicalForm& gamma,
[f876a66]168                    const int k, const CanonicalForm& delta,
169                    CFList& source, CFList& dest)
[7bf145]170{
171  bool result;
[806c18]172  if (CFFactory::gettype() == GaloisFieldDomain)
[7bf145]173  {
174    int p= getCharacteristic();
175    int orderFieldExtension= ipower (p, getGFDegree()) - 1;
176    int order= ipower (p, k) - 1;
177    int number= orderFieldExtension/order;
178    result= GFInExtensionHelper (F, number);
179    return result;
180  }
[806c18]181  else
[7bf145]182  {
[f876a66]183    result= FqInExtensionHelper (F, gamma, delta, source, dest);
[7bf145]184    return result;
185  }
186}
187
[806c18]188CanonicalForm
189mapDown (const CanonicalForm& F, const ExtensionInfo& info, CFList& source,
190         CFList& dest)
[7bf145]191{
192  int k= info.getGFDegree();
193  Variable beta= info.getAlpha();
194  CanonicalForm primElem= info.getGamma();
195  CanonicalForm imPrimElem= info.getDelta();
196  if (k > 1)
197    return GFMapDown (F, k);
198  else if (k == 1)
199    return F;
[c1b9927]200  if (/*k==0 &&*/ beta == Variable (1))
[7bf145]201    return F;
[c1b9927]202  else /*if (k==0 && beta != Variable (1))*/
[618da5]203    return mapDown (F, imPrimElem, primElem, beta, source, dest);
[7bf145]204}
205
[806c18]206void appendTestMapDown (CFList& factors, const CanonicalForm& f,
207                        const ExtensionInfo& info, CFList& source, CFList& dest)
[7bf145]208{
209  int k= info.getGFDegree();
210  Variable beta= info.getBeta();
211  Variable alpha= info.getAlpha();
212  CanonicalForm delta= info.getDelta();
213  CanonicalForm gamma= info.getGamma();
214  CanonicalForm g= f;
215  int degMipoBeta;
[1372ae]216  if (!k && beta.level() == 1)
[7bf145]217    degMipoBeta= 1;
[1372ae]218  else if (!k && beta.level() != 1)
[7bf145]219    degMipoBeta= degree (getMipo (beta));
[806c18]220  if (k > 1)
[7bf145]221  {
[f876a66]222    if (!isInExtension (g, gamma, k, delta, source, dest))
[7bf145]223    {
224      g= GFMapDown (g, k);
225      factors.append (g);
226    }
227  }
[806c18]228  else if (k == 1)
229  {
[1709e1]230    if (!isInExtension (g, gamma, k, delta, source, dest))
[7bf145]231      factors.append (g);
232  }
[806c18]233  else if (!k && beta == Variable (1))
[7bf145]234  {
235    if (degree (g, alpha) < degMipoBeta)
236      factors.append (g);
237  }
[806c18]238  else if (!k && beta != Variable (1))
[7bf145]239  {
[f876a66]240    if (!isInExtension (g, gamma, k, delta, source, dest))
[7bf145]241    {
[618da5]242      g= mapDown (g, delta, gamma, alpha, source, dest);
[7bf145]243      factors.append (g);
[806c18]244    }
[7bf145]245  }
246  return;
247}
248
[806c18]249void
250appendMapDown (CFList& factors, const CanonicalForm& g,
251               const ExtensionInfo& info, CFList& source, CFList& dest)
[7bf145]252{
253  int k= info.getGFDegree();
254  Variable beta= info.getBeta();
255  Variable alpha= info.getAlpha();
256  CanonicalForm gamma= info.getGamma();
257  CanonicalForm delta= info.getDelta();
[806c18]258  if (k > 1)
[7bf145]259    factors.append (GFMapDown (g, k));
[806c18]260  else if (k == 1)
[7bf145]261    factors.append (g);
[806c18]262  else if (!k && beta == Variable (1))
[7bf145]263    factors.append (g);
[806c18]264  else if (!k && beta != Variable (1))
[618da5]265    factors.append (mapDown (g, delta, gamma, alpha, source, dest));
[7bf145]266  return;
267}
268
[806c18]269void normalize (CFList& factors)
[7bf145]270{
[1682691]271  CanonicalForm lcinv;
[806c18]272  for (CFListIterator i= factors; i.hasItem(); i++)
[1682691]273  {
274    lcinv= 1/Lc (i.getItem());
275    i.getItem() *= lcinv;
276  }
[7bf145]277  return;
278}
279
[f876a66]280void normalize (CFFList& factors)
281{
[1682691]282  CanonicalForm lcinv;
[c1b9927]283  for (CFFListIterator i= factors; i.hasItem(); i++)
[1682691]284  {
285    lcinv= 1/ Lc (i.getItem().factor());
286    i.getItem()= CFFactor (i.getItem().factor()*lcinv,
[f876a66]287                           i.getItem().exp());
[1682691]288  }
[f876a66]289  return;
290}
291
[806c18]292CFList subset (int index [], const int& s, const CFArray& elements,
293               bool& noSubset)
[7bf145]294{
295  int r= elements.size();
296  int i= 0;
297  CFList result;
298  noSubset= false;
[806c18]299  if (index[s - 1] == 0)
[7bf145]300  {
[806c18]301    while (i < s)
[7bf145]302    {
303      index[i]= i + 1;
304      result.append (elements[i]);
305      i++;
[806c18]306    }
[7bf145]307    return result;
308  }
309  int buf;
310  int k;
311  bool found= false;
[806c18]312  if (index[s - 1] == r)
[7bf145]313  {
[806c18]314    if (index[0] == r - s + 1)
[7bf145]315    {
316      noSubset= true;
317      return result;
318    }
319    else {
[806c18]320      while (found == false)
[7bf145]321      {
322        if (index[s - 2 - i] < r - i - 1)
323          found= true;
[806c18]324        i++;
[7bf145]325      }
326      buf= index[s - i - 1];
327      k= 0;
[806c18]328      while (s - i - 1 + k < s)
[7bf145]329      {
330        index[s - i - 1 + k]= buf + k + 1;
331        k++;
[806c18]332      }
[7bf145]333    }
[806c18]334    for (int j= 0; j < s; j++)
[7bf145]335      result.append (elements[index[j] - 1]);
336    return result;
337  }
[806c18]338  else
[7bf145]339  {
340    index[s - 1] += 1;
341    for (int j= 0; j < s; j++)
342      result.append (elements[index[j] - 1]);
343    return result;
344  }
[806c18]345}
[7bf145]346
[806c18]347CFArray copy (const CFList& list)
[7bf145]348{
349  CFArray array= CFArray (list.length());
350  int j= 0;
[806c18]351  for (CFListIterator i= list; i.hasItem(); i++, j++)
[7bf145]352    array[j]= i.getItem();
353  return array;
354}
355
[806c18]356void indexUpdate (int index [], const int& subsetSize, const int& setSize,
357                   bool& noSubset)
[7bf145]358{
359  noSubset= false;
[806c18]360  if (subsetSize > setSize)
[7bf145]361  {
362    noSubset= true;
363    return;
364  }
[1673386]365  int * v= new int [setSize];
[806c18]366  for (int i= 0; i < setSize; i++)
[7bf145]367    v[i]= index[i];
[806c18]368  if (subsetSize == 1)
[7bf145]369  {
370    v[0]= v[0] - 1;
[806c18]371    if (v[0] >= setSize)
[7bf145]372    {
373      noSubset= true;
[1673386]374      delete [] v;
[7bf145]375      return;
376    }
377  }
[806c18]378  else
[7bf145]379  {
[806c18]380    if (v[subsetSize - 1] - v[0] + 1 == subsetSize && v[0] > 1)
[7bf145]381    {
[806c18]382      if (v[0] + subsetSize - 1 > setSize)
383      {
[7bf145]384        noSubset= true;
[1673386]385        delete [] v;
[7bf145]386        return;
387      }
388      v[0]= v[0] - 1;
[806c18]389      for (int i= 1; i < subsetSize - 1; i++)
[7bf145]390        v[i]= v[i - 1] + 1;
391      v[subsetSize - 1]= v[subsetSize - 2];
392    }
[806c18]393    else
[7bf145]394    {
[806c18]395      if (v[0] + subsetSize - 1 > setSize)
396      {
[7bf145]397        noSubset= true;
[1673386]398        delete [] v;
[7bf145]399        return;
400      }
401      for (int i= 1; i < subsetSize - 1; i++)
402        v[i]= v[i - 1] + 1;
403      v[subsetSize - 1]= v[subsetSize - 2];
404    }
405  }
406  for (int i= 0; i < setSize; i++)
407    index[i]= v[i];
[1673386]408  delete [] v;
[7bf145]409}
410
[806c18]411int subsetDegree (const CFList& S)
[7bf145]412{
413  int result= 0;
414  for (CFListIterator i= S; i.hasItem(); i++)
415    result += degree (i.getItem(), Variable (1));
416  return result;
417}
418
419CFFList multiplicity (CanonicalForm& F, const CFList& factors)
420{
421  if (F.inCoeffDomain())
422    return CFFList (CFFactor (F, 1));
423  CFFList result;
424  int multi= 0;
[21b8f4]425  CanonicalForm quot;
[7bf145]426  for (CFListIterator i= factors; i.hasItem(); i++)
427  {
[21b8f4]428    while (fdivides (i.getItem(), F, quot))
[7bf145]429    {
430      multi++;
[21b8f4]431      F= quot;
[7bf145]432    }
[575e1d]433    if (multi > 0)
434      result.append (CFFactor (i.getItem(), multi));
[7bf145]435    multi= 0;
436  }
437  return result;
438}
439
[615ca8]440#ifdef HAVE_NTL
[f876a66]441CFArray
442logarithmicDerivative (const CanonicalForm& F, const CanonicalForm& G, int l,
443                       CanonicalForm& Q
444                      )
445{
446  Variable x= Variable (2);
447  Variable y= Variable (1);
448  CanonicalForm xToL= power (x, l);
449  CanonicalForm q,r;
450  CanonicalForm logDeriv;
451
452  q= newtonDiv (F, G, xToL);
453
454  logDeriv= mulMod2 (q, deriv (G, y), xToL);
455
[9614bb]456  int j= degree (logDeriv, y) + 1;
[f876a66]457  CFArray result= CFArray (j);
[9614bb]458  CFIterator ii;
[f876a66]459  for (CFIterator i= logDeriv; i.hasTerms() && !logDeriv.isZero(); i++)
460  {
[9614bb]461    for (ii= i.coeff(); ii.hasTerms(); ii++)
462      result[ii.exp()] += ii.coeff()*power (x,i.exp());
[f876a66]463  }
464  Q= q;
465  return result;
466}
467
468CFArray
469logarithmicDerivative (const CanonicalForm& F, const CanonicalForm& G, int l,
470                       int oldL, const CanonicalForm& oldQ, CanonicalForm& Q
471                      )
472{
473  Variable x= Variable (2);
474  Variable y= Variable (1);
475  CanonicalForm xToL= power (x, l);
476  CanonicalForm xToOldL= power (x, oldL);
477  CanonicalForm xToLOldL= power (x, l-oldL);
478  CanonicalForm q,r;
479  CanonicalForm logDeriv;
480
[ac3fcca]481  CanonicalForm bufF;
482  if ((oldL > 100 && l - oldL < 50) || (oldL < 100 && l - oldL < 30))
483  {
484    bufF= F;
485    CanonicalForm oldF= mulMod2 (G, oldQ, xToL);
486    bufF -= oldF;
487    bufF= div (bufF, xToOldL);
488  }
489  else
490  {
491    //middle product style computation of [G*oldQ]^{l}_{oldL}
492    CanonicalForm G3= div (G, xToOldL);
493    CanonicalForm Up= mulMod2 (G3, oldQ, xToLOldL);
[fd80670]494    CanonicalForm xToOldL2= power (x, (oldL+1)/2);
[ac3fcca]495    CanonicalForm G2= mod (G, xToOldL);
496    CanonicalForm G1= div (G2, xToOldL2);
497    CanonicalForm G0= mod (G2, xToOldL2);
498    CanonicalForm oldQ1= div (oldQ, xToOldL2);
499    CanonicalForm oldQ0= mod (oldQ, xToOldL2);
[fd80670]500    CanonicalForm Mid;
501    if (oldL % 2 == 1)
502      Mid= mulMod2 (G1, oldQ1*x, xToLOldL);
503    else
504      Mid= mulMod2 (G1, oldQ1, xToLOldL);
[ac3fcca]505    //computation of Low might be faster using a real middle product?
506    CanonicalForm Low= mulMod2 (G0, oldQ1, xToOldL)+mulMod2 (G1, oldQ0, xToOldL);
[fd80670]507    Low= div (Low, power (x, oldL/2));
508    Low= mod (Low, xToLOldL);
[ac3fcca]509    Up += Mid + Low;
510    bufF= div (F, xToOldL);
511    bufF -= Up;
512  }
[f876a66]513
[72f1e4b]514  if (l-oldL > 0)
515    q= newtonDiv (bufF, G, xToLOldL);
516  else
517    q= 0;
[f876a66]518  q *= xToOldL;
519  q += oldQ;
520
521  logDeriv= mulMod2 (q, deriv (G, y), xToL);
522
[9614bb]523  int j= degree (logDeriv,y) + 1;
[f876a66]524  CFArray result= CFArray (j);
[9614bb]525  CFIterator ii;
[f876a66]526  for (CFIterator i= logDeriv; i.hasTerms() && !logDeriv.isZero(); i++)
527  {
[9614bb]528    for (ii= i.coeff(); ii.hasTerms(); ii++)
529      result[ii.exp()] += ii.coeff()*power (x,i.exp());
[f876a66]530  }
531  Q= q;
532  return result;
533}
[615ca8]534#endif
[f876a66]535
536void
537writeInMatrix (CFMatrix& M, const CFArray& A, const int column,
538               const int startIndex
539              )
540{
[050d1b]541  ASSERT (A.size () - startIndex >= 0, "wrong starting index");
542  ASSERT (A.size () - startIndex <= M.rows(), "wrong starting index");
[f876a66]543  ASSERT (column > 0 && column <= M.columns(), "wrong column");
544  if (A.size() - startIndex <= 0) return;
545  int j= 1;
546  for (int i= startIndex; i < A.size(); i++, j++)
547    M (j, column)= A [i];
548}
549
550CFArray getCoeffs (const CanonicalForm& F, const int k)
551{
[050d1b]552  ASSERT (F.isUnivariate() || F.inCoeffDomain(), "univariate input expected");
[f876a66]553  if (degree (F, 2) < k)
554    return CFArray();
555
556  CFArray result= CFArray (degree (F) - k + 1);
557  CFIterator j= F;
558  for (int i= degree (F); i >= k; i--)
559  {
560    if (j.exp() == i)
561    {
562      result [i - k]= j.coeff();
[c1b9927]563      j++;
[f876a66]564      if (!j.hasTerms())
565        return result;
566    }
567    else
568      result[i - k]= 0;
569  }
570  return result;
571}
572
573CFArray getCoeffs (const CanonicalForm& F, const int k, const Variable& alpha)
574{
[050d1b]575  ASSERT (F.isUnivariate() || F.inCoeffDomain(), "univariate input expected");
[f876a66]576  if (degree (F, 2) < k)
577    return CFArray ();
578
579  int d= degree (getMipo (alpha));
580  CFArray result= CFArray ((degree (F) - k + 1)*d);
581  CFIterator j= F;
582  CanonicalForm buf;
583  CFIterator iter;
584  for (int i= degree (F); i >= k; i--)
585  {
586    if (j.exp() == i)
587    {
588      iter= j.coeff();
589      for (int l= degree (j.coeff(), alpha); l >= 0; l--)
590      {
591        if (iter.exp() == l)
592        {
593          result [(i - k)*d + l]= iter.coeff();
594          iter++;
595          if (!iter.hasTerms())
596            break;
597        }
598      }
599      j++;
600      if (!j.hasTerms())
601        return result;
602    }
603    else
604    {
605      for (int l= 0; l < d; l++)
606        result[(i - k)*d + l]= 0;
607    }
608  }
609  return result;
610}
611
612#ifdef HAVE_NTL
613CFArray
614getCoeffs (const CanonicalForm& G, const int k, const int l, const int degMipo,
615           const Variable& alpha, const CanonicalForm& evaluation,
616           const mat_zz_p& M)
617{
[050d1b]618  ASSERT (G.isUnivariate() || G.inCoeffDomain(), "univariate input expected");
[f876a66]619  CanonicalForm F= G (G.mvar() - evaluation, G.mvar());
620  if (F.isZero())
621    return CFArray ();
622
623  Variable y= Variable (2);
624  F= F (power (y, degMipo), y);
625  F= F (y, alpha);
626  zz_pX NTLF= convertFacCF2NTLzzpX (F);
627  NTLF.rep.SetLength (l*degMipo);
628  NTLF.rep= M*NTLF.rep;
629  NTLF.normalize();
630  F= convertNTLzzpX2CF (NTLF, y);
631
632  if (degree (F, 2) < k)
633    return CFArray();
634
635  CFArray result= CFArray (degree (F) - k + 1);
636
637  CFIterator j= F;
638  for (int i= degree (F); i >= k; i--)
639  {
640    if (j.exp() == i)
641    {
642      result [i - k]= j.coeff();
[c1b9927]643      j++;
[f876a66]644      if (!j.hasTerms())
645        return result;
646    }
647    else
648      result[i - k]= 0;
649  }
650  return result;
651}
652#endif
653
[542864]654int * computeBounds (const CanonicalForm& F, int& n, bool& isIrreducible)
[f876a66]655{
656  n= degree (F, 1);
657  int* result= new int [n];
658  int sizeOfNewtonPolygon;
659  int** newtonPolyg= newtonPolygon (F, sizeOfNewtonPolygon);
660
[542864]661  isIrreducible= false;
662  if (sizeOfNewtonPolygon == 3)
663  {
664    bool check1=
665        (newtonPolyg[0][0]==0 || newtonPolyg[1][0]==0 || newtonPolyg[2][0]==0);
666    if (check1)
667    {
668      bool check2=
669        (newtonPolyg[0][1]==0 || newtonPolyg[1][1]==0 || newtonPolyg[2][0]==0);
670      if (check2)
671      {
672        int p=getCharacteristic();
673        int d=1;
674        char bufGFName='Z';
675        bool GF= (CFFactory::gettype()==GaloisFieldDomain);
676        if (GF)
677        {
678          d= getGFDegree();
679          bufGFName=gf_name;
680        }
681        setCharacteristic(0);
682        CanonicalForm tmp= gcd (newtonPolyg[0][0],newtonPolyg[0][1]);
683        tmp= gcd (tmp, newtonPolyg[1][0]);
684        tmp= gcd (tmp, newtonPolyg[1][1]);
685        tmp= gcd (tmp, newtonPolyg[2][0]);
686        tmp= gcd (tmp, newtonPolyg[2][1]);
687        isIrreducible= (tmp==1);
688        if (GF)
689          setCharacteristic (p, d, bufGFName);
690        else
691          setCharacteristic(p);
692      }
693    }
694  }
695
[f876a66]696  int minX, minY, maxX, maxY;
697  minX= newtonPolyg [0] [0];
698  minY= newtonPolyg [0] [1];
699  maxX= minX;
700  maxY= minY;
701  for (int i= 1; i < sizeOfNewtonPolygon; i++)
702  {
703    if (minX > newtonPolyg [i] [0])
704      minX= newtonPolyg [i] [0];
705    if (maxX < newtonPolyg [i] [0])
706      maxX= newtonPolyg [i] [0];
707    if (minY > newtonPolyg [i] [1])
708      minY= newtonPolyg [i] [1];
709    if (maxY < newtonPolyg [i] [1])
710      maxY= newtonPolyg [i] [1];
711  }
712
713  int k= maxX;
714  for (int i= 0; i < n; i++)
715  {
716    if (i + 1 > maxY || i + 1 < minY)
717    {
718      result [i]= 0;
719      continue;
720    }
721    int* point= new int [2];
722    point [0]= k;
723    point [1]= i + 1;
724    while (!isInPolygon (newtonPolyg, sizeOfNewtonPolygon, point) && k > 0)
725    {
726      k--;
727      point [0]= k;
728    }
729    result [i]= k;
730    k= maxX;
731    delete [] point;
732  }
733
734  return result;
735}
736
[dd8047]737int
738substituteCheck (const CanonicalForm& F, const Variable& x)
739{
740  if (F.inCoeffDomain())
741    return 0;
742  if (degree (F, x) < 0)
743    return 0;
744  CanonicalForm f= swapvar (F, F.mvar(), x);
745  int sizef= 0;
746  for (CFIterator i= f; i.hasTerms(); i++, sizef++)
747  {
748    if (i.exp() == 1)
749      return 0;
750  }
751  int * expf= new int [sizef];
752  int j= 0;
753  for (CFIterator i= f; i.hasTerms(); i++, j++)
754    expf [j]= i.exp();
755
756  int indf= sizef - 1;
757  if (expf[indf] == 0)
758    indf--;
759
760  int result= expf[indf];
761  for (int i= indf - 1; i >= 0; i--)
762  {
763    if (expf [i]%result != 0)
764    {
765      delete [] expf;
766      return 0;
767    }
768  }
769
770  delete [] expf;
771  return result;
772}
773
[686ce3]774static int
775substituteCheck (const CanonicalForm& F, const CanonicalForm& G)
776{
777  if (F.inCoeffDomain() || G.inCoeffDomain())
778    return 0;
779  Variable x= Variable (1);
780  if (degree (F, x) <= 1 || degree (G, x) <= 1)
781    return 0;
782  CanonicalForm f= swapvar (F, F.mvar(), x);
783  CanonicalForm g= swapvar (G, G.mvar(), x);
784  int sizef= 0;
785  int sizeg= 0;
786  for (CFIterator i= f; i.hasTerms(); i++, sizef++)
787  {
788    if (i.exp() == 1)
789      return 0;
790  }
791  for (CFIterator i= g; i.hasTerms(); i++, sizeg++)
792  {
793    if (i.exp() == 1)
794      return 0;
795  }
796  int * expf= new int [sizef];
797  int * expg= new int [sizeg];
798  int j= 0;
799  for (CFIterator i= f; i.hasTerms(); i++, j++)
800  {
801    expf [j]= i.exp();
802  }
803  j= 0;
804  for (CFIterator i= g; i.hasTerms(); i++, j++)
805  {
806    expg [j]= i.exp();
807  }
808
809  int indf= sizef - 1;
810  int indg= sizeg - 1;
811  if (expf[indf] == 0)
812    indf--;
813  if (expg[indg] == 0)
814    indg--;
815
816  if ((expg[indg]%expf [indf] != 0 && expf[indf]%expg[indg] != 0) ||
817      (expg[indg] == 1 && expf[indf] == 1))
818  {
819    delete [] expg;
820    delete [] expf;
821    return 0;
822  }
823
824  int result;
825  if (expg [indg]%expf [indf] == 0)
826    result= expf[indf];
827  else
828    result= expg[indg];
829  for (int i= indf - 1; i >= 0; i--)
830  {
831    if (expf [i]%result != 0)
832    {
833      delete [] expf;
834      delete [] expg;
835      return 0;
836    }
837  }
838
839  for (int i= indg - 1; i >= 0; i--)
840  {
841    if (expg [i]%result != 0)
842    {
843      delete [] expf;
844      delete [] expg;
845      return 0;
846    }
847  }
848
849  delete [] expg;
850  delete [] expf;
851  return result;
852}
853
854int recSubstituteCheck (const CanonicalForm& F, const int d)
855{
856  if (F.inCoeffDomain())
857    return 0;
858  Variable x= Variable (1);
859  if (degree (F, x) <= 1)
860    return 0;
861  CanonicalForm f= swapvar (F, F.mvar(), x);
862  int sizef= 0;
863  for (CFIterator i= f; i.hasTerms(); i++, sizef++)
864  {
865    if (i.exp() == 1)
866      return 0;
867  }
868  int * expf= new int [sizef];
869  int j= 0;
870  for (CFIterator i= f; i.hasTerms(); i++, j++)
871  {
872    expf [j]= i.exp();
873  }
874
875  int indf= sizef - 1;
876  if (expf[indf] == 0)
877    indf--;
878
879  if ((d%expf [indf] != 0 && expf[indf]%d != 0) || (expf[indf] == 1))
880  {
881    delete [] expf;
882    return 0;
883  }
884
885  int result;
886  if (d%expf [indf] == 0)
887    result= expf[indf];
888  else
889    result= d;
890  for (int i= indf - 1; i >= 0; i--)
891  {
892    if (expf [i]%result != 0)
893    {
894      delete [] expf;
895      return 0;
896    }
897  }
898
899  delete [] expf;
900  return result;
901}
902
903int substituteCheck (const CFList& L)
904{
905  ASSERT (L.length() > 1, "expected a list of at least two elements");
906  if (L.length() < 2)
907    return 0;
908  CFListIterator i= L;
909  i++;
910  int result= substituteCheck (L.getFirst(), i.getItem());
911  if (result <= 1)
912    return result;
913  i++;
914  for (;i.hasItem(); i++)
915  {
916    result= recSubstituteCheck (i.getItem(), result);
917    if (result <= 1)
918      return result;
919  }
920  return result;
921}
922
923void
924subst (const CanonicalForm& F, CanonicalForm& A, const int d, const Variable& x)
925{
926  if (d <= 1)
927  {
928    A= F;
929    return;
930  }
931  if (degree (F, x) <= 0)
932  {
933    A= F;
934    return;
935  }
936  CanonicalForm C= 0;
937  CanonicalForm f= swapvar (F, x, F.mvar());
938  for (CFIterator i= f; i.hasTerms(); i++)
939    C += i.coeff()*power (f.mvar(), i.exp()/ d);
940  A= swapvar (C, x, F.mvar());
941}
942
943CanonicalForm
944reverseSubst (const CanonicalForm& F, const int d, const Variable& x)
945{
946  if (d <= 1)
947    return F;
948  if (degree (F, x) <= 0)
949    return F;
950  CanonicalForm f= swapvar (F, x, F.mvar());
951  CanonicalForm result= 0;
952  for (CFIterator i= f; i.hasTerms(); i++)
953    result += i.coeff()*power (f.mvar(), d*i.exp());
954  return swapvar (result, x, F.mvar());
955}
956
957void
958reverseSubst (CFList& L, const int d, const Variable& x)
959{
960  for (CFListIterator i= L; i.hasItem(); i++)
961    i.getItem()= reverseSubst (i.getItem(), d, x);
962}
963
Note: See TracBrowser for help on using the repository browser.