source: git/factory/cfCharSets.cc @ 8d1432e

spielwiese
Last change on this file since 8d1432e was c65b73, checked in by Martin Lee <martinlee84@…>, 9 years ago
format: fix indentation
  • Property mode set to 100644
File size: 15.4 KB
Line 
1/*****************************************************************************\
2 * Computer Algebra System SINGULAR
3\*****************************************************************************/
4/** @file cfCharSets.cc
5 *
6 * This file provides functions to compute characteristic sets
7 *
8 * @note some of the code is code from libfac or derived from code from libfac.
9 * Libfac is written by M. Messollen. See also COPYING for license information
10 * and README for general information on characteristic sets.
11 *
12 * ABSTRACT: Descriptions can be found in Wang "On the Parallelization of
13 * characteristic-set based algorithms" or Greuel/Pfister "A Singular
14 * Introduction to Commutative Algebra".
15 *
16 * @authors Martin Lee
17 *
18 **/
19/*****************************************************************************/
20
21#include "config.h"
22
23#include "timing.h"
24
25#include "canonicalform.h"
26#include "cfCharSets.h"
27#include "cfCharSetsUtil.h"
28#include "cf_algorithm.h"
29#include "facAlgFunc.h"
30
31TIMING_DEFINE_PRINT(neworder_time)
32
33// set up a new orderd list of Variables.
34// we try to reorder the variables heuristically optimal.
35Varlist
36neworder (const CFList & PolyList)
37{
38  CFList PS= PolyList, PS1=PolyList;
39  Varlist oldorder, reorder, difference;
40
41  TIMING_START (neworder_time);
42
43  int highest_level= level (get_max_var (PS));
44
45  // set up oldorder and first criterion: only_in_one
46  for (int i= highest_level; i>=1; i--)
47  {
48    oldorder.insert (Variable (i));
49    CFList is_one= only_in_one (PS1, Variable (i));
50    if (is_one.length() == 1)
51    {
52      reorder.insert (Variable (i));
53      PS1= Difference (PS1, is_one);
54    }
55    else if (is_one.length() == 0)
56    {
57      reorder.append (Variable (i)); // assigne it the highest level
58      PS1= Difference (PS1, is_one);
59    }
60  }
61  difference= Difference (oldorder, reorder);
62
63  // rearrange the ordering of the variables!
64  difference= reorderb (difference, PS, highest_level);
65  reorder= Union (reorder, difference);
66  TIMING_END(neworder_time);
67
68  TIMING_PRINT(neworder_time, "\ntime used for neworder   : ");
69
70  return Union (reorder, Difference (oldorder, reorder));
71}
72
73// the same as above, only returning a list of CanonicalForms
74CFList
75newordercf (const CFList & PolyList)
76{
77  Varlist reorder= neworder (PolyList);
78  CFList output;
79
80  for (VarlistIterator i=reorder; i.hasItem(); i++)
81    output.append (CanonicalForm (i.getItem()));
82
83  return output;
84}
85
86// the same as above, only returning a list of ints
87IntList
88neworderint (const CFList & PolyList)
89{
90  Varlist reorder= neworder (PolyList);
91  IntList output;
92
93  for (VarlistIterator i= reorder; i.hasItem(); i++)
94    output.append (level (i.getItem()));
95
96  return output;
97}
98
99// a library function: we reorganize the global variable ordering
100CFList
101reorder (const Varlist & betterorder, const CFList & PS)
102{
103  int i= 1, n= betterorder.length();
104  Intarray v (1, n);
105  CFList ps= PS;
106
107  //initalize:
108  for (VarlistIterator j= betterorder; j.hasItem(); j++)
109  {
110    v[i]= level (j.getItem());
111    i++;
112  }
113  // reorder:
114  for (i= 1; i <= n; i++)
115    ps= swapvar (ps, Variable (v[i]), Variable (n + i));
116  return ps;
117}
118
119CFFList
120reorder (const Varlist & betterorder, const CFFList & PS)
121{
122  int i= 1, n= betterorder.length();
123  Intarray v (1, n);
124  CFFList ps= PS;
125
126  //initalize:
127  for (VarlistIterator j= betterorder; j.hasItem(); j++)
128  {
129    v[i]= level (j.getItem());
130    i++;
131  }
132
133  // reorder:
134  for (i= 1; i <= n; i++)
135    ps= swapvar (ps, Variable (v[i]), Variable (n + i));
136  return ps;
137}
138
139ListCFList
140reorder (const Varlist & betterorder, const ListCFList & Q)
141{
142  ListCFList Q1;
143
144  for (ListCFListIterator i= Q; i.hasItem(); i++)
145    Q1.append (reorder (betterorder, i.getItem()));
146  return Q1;
147}
148
149CFList
150basicSet (const CFList &PS)
151{
152  CFList QS= PS, BS, RS;
153  CanonicalForm b;
154  int cb, degb;
155
156  if (PS.length() < 2)
157    return PS;
158
159  CFListIterator i;
160
161  while (!QS.isEmpty())
162  {
163    b= lowestRank (QS);
164    cb= b.level();
165
166    BS= Union(CFList (b), BS);
167
168    if (cb <= 0)
169      return CFList();
170    else
171    {
172      degb= degree (b);
173      RS= CFList();
174      for (i= QS; i.hasItem(); i++)
175      {
176        if (degree (i.getItem(), cb) < degb)
177          RS= Union (CFList (i.getItem()), RS);
178      }
179      QS= RS;
180    }
181  }
182
183  return BS;
184}
185
186CFList
187charSet (const CFList &PS)
188{
189  CFList QS= PS, RS= PS, CSet, tmp;
190  CFListIterator i;
191  CanonicalForm r;
192
193  while (!RS.isEmpty())
194  {
195    CSet= basicSet (QS);
196
197    RS= CFList();
198    if (CSet.length() > 0 && CSet.getFirst().level() > 0)
199    {
200      tmp= Difference (QS, CSet);
201      for (i= tmp; i.hasItem(); i++)
202      {
203        r= Prem (i.getItem(), CSet);
204        if (r != 0)
205          RS= Union (RS, CFList (r));
206      }
207      QS= Union (QS, RS);
208    }
209  }
210
211  return CSet;
212}
213
214/// medial set
215CFList
216charSetN (const CFList &PS)
217{
218  CFList QS= PS, RS= PS, CSet, tmp;
219  CFListIterator i;
220  CanonicalForm r;
221
222  while (!RS.isEmpty())
223  {
224    QS= uniGcd (QS);
225    CSet= basicSet (QS);
226
227    RS= CFList();
228    if (CSet.length() > 0 && CSet.getFirst().level() > 0)
229    {
230      tmp= Difference (QS, CSet);
231      for (i= tmp; i.hasItem(); i++)
232      {
233        r= Prem (i.getItem(), CSet);
234        if (!r.isZero())
235          RS= Union (RS, CFList (r));
236      }
237      QS= Union (CSet, RS);
238    }
239  }
240
241  return CSet;
242}
243
244/// compute a characteristic set via medial set
245CFList
246charSetViaCharSetN (const CFList& PS)
247{
248  CFList L;
249  CFFList sqrfFactors;
250  CanonicalForm sqrf;
251  CFFListIterator iter2;
252  for (CFListIterator iter= PS; iter.hasItem(); iter++)
253  {
254    sqrf= 1;
255    sqrfFactors= sqrFree (iter.getItem());
256    for (iter2= sqrfFactors; iter2.hasItem(); iter2++)
257      sqrf *= iter2.getItem().factor();
258    L= Union (L, CFList (normalize (sqrf)));
259  }
260
261  CFList result= charSetN (L);
262
263  if (result.isEmpty() || result.getFirst().inCoeffDomain())
264    return CFList(1);
265
266  CanonicalForm r;
267  CFList RS;
268  CFList tmp= Difference (L, result);
269
270  for (CFListIterator i= tmp; i.hasItem(); i++)
271  {
272    r= Premb (i.getItem(), result);
273    if (!r.isZero())
274      RS= Union (RS, CFList (r));
275  }
276  if (RS.isEmpty())
277    return result;
278
279  return charSetViaCharSetN (Union (L, Union (RS, result)));
280}
281
282/// modified medial set
283CFList
284modCharSet (const CFList& L, StoreFactors& StoredFactors, bool removeContents)
285{
286  CFList QS, RS= L, CSet, tmp, contents, initial, removedFactors;
287  CFListIterator i;
288  CanonicalForm r, cF;
289  bool noRemainder= true;
290  StoreFactors StoredFactors2;
291
292  QS= uniGcd (L);
293
294  while (!RS.isEmpty())
295  {
296    noRemainder= true;
297    CSet= basicSet (QS);
298
299    initial= factorsOfInitials (CSet);
300
301    StoredFactors2.FS1= StoredFactors.FS1;
302    StoredFactors2.FS2= Union (StoredFactors2.FS2, initial);
303
304    RS= CFList();
305
306    if (CSet.length() > 0 && CSet.getFirst().level() > 0)
307    {
308      tmp= Difference (QS, CSet);
309
310      for (i= tmp; i.hasItem(); i++)
311      {
312        r= Prem (i.getItem(), CSet);
313        if (!r.isZero())
314        {
315          noRemainder= false;
316          if (removeContents)
317          {
318            removeContent (r, cF);
319
320            if (!cF.isZero())
321              contents= Union (contents, factorPSet (CFList(cF))); //factorPSet maybe too much it should suffice to do a squarefree factorization instead
322          }
323
324          removeFactors (r, StoredFactors2, removedFactors);
325          StoredFactors2.FS1= Union (StoredFactors2.FS1, removedFactors);
326          StoredFactors2.FS2= Difference (StoredFactors2.FS2, removedFactors);
327
328          removedFactors= CFList();
329
330          RS= Union (RS, CFList (r));
331        }
332      }
333
334      if (removeContents && !noRemainder)
335      {
336        StoredFactors.FS1= Union (StoredFactors2.FS1, contents);
337        StoredFactors.FS2= StoredFactors2.FS2;
338      }
339      else
340        StoredFactors= StoredFactors2;
341
342      QS= Union (CSet, RS);
343
344      contents= CFList();
345      removedFactors= CFList();
346    }
347    else
348      StoredFactors= StoredFactors2;
349  }
350
351  return CSet;
352}
353
354/// characteristic set via modified medial set
355CFList
356charSetViaModCharSet (const CFList& PS, StoreFactors& StoredFactors,
357                      bool removeContents)
358{
359  CFList L;
360  CFFList sqrfFactors;
361  CanonicalForm sqrf;
362  CFFListIterator iter2;
363  for (CFListIterator iter= PS; iter.hasItem(); iter++)
364  {
365    sqrf= 1;
366    sqrfFactors= sqrFree (iter.getItem());
367    for (iter2= sqrfFactors; iter2.hasItem(); iter2++)
368      sqrf *= iter2.getItem().factor();
369    L= Union (L, CFList (normalize (sqrf)));
370  }
371
372  L= uniGcd (L);
373
374  CFList result= modCharSet (L, StoredFactors, removeContents);
375
376  if (result.isEmpty() || result.getFirst().inCoeffDomain())
377    return CFList(1);
378
379  CanonicalForm r;
380  CFList RS;
381  CFList tmp= Difference (L, result);
382
383  for (CFListIterator i= tmp; i.hasItem(); i++)
384  {
385    r= Premb (i.getItem(), result);
386    if (!r.isZero())
387      RS= Union (RS, CFList (r));
388  }
389  if (RS.isEmpty())
390    return result;
391
392  return charSetViaModCharSet (Union (L, Union (RS, result)), StoredFactors,
393                               removeContents);
394}
395
396CFList
397charSetViaModCharSet (const CFList& PS, bool removeContents)
398{
399  StoreFactors tmp;
400  return charSetViaModCharSet (PS, tmp, removeContents);
401}
402
403CFList
404modCharSet (const CFList& PS, bool removeContents)
405{
406  StoreFactors tmp;
407  return modCharSet (PS, tmp, removeContents);
408}
409
410ListCFList
411charSeries (const CFList& L)
412{
413  ListCFList tmp, result, tmp2, ppi1, ppi2, qqi, ppi, alreadyConsidered;
414  CFList l, charset, ini;
415
416  int count= 0;
417  int highestLevel= 1;
418  CFListIterator iter;
419
420  StoreFactors StoredFactors;
421
422  l= L;
423
424  for (iter= l; iter.hasItem(); iter++)
425  {
426    iter.getItem()= normalize (iter.getItem());
427    if (highestLevel < iter.getItem().level())
428      highestLevel= iter.getItem().level();
429  }
430
431  tmp= ListCFList (l);
432
433  while (!tmp.isEmpty())
434  {
435    sortListCFList (tmp);
436
437    l= tmp.getFirst();
438
439    tmp= Difference (tmp, l);
440
441    select (ppi, l.length(), ppi1, ppi2);
442
443    inplaceUnion (ppi2, qqi);
444
445    if (count > 0)
446      ppi= Union (ppi1, ListCFList (l));
447    else
448      ppi= ListCFList();
449
450    if (l.length() - 3 < highestLevel)
451      charset= charSetViaModCharSet (l, StoredFactors);
452    else
453      charset= charSetViaCharSetN (l);
454
455    if (charset.length() > 0 && charset.getFirst().level() > 0)
456    {
457      result= Union (ListCFList (charset), result);
458      ini= factorsOfInitials (charset);
459
460      ini= Union (ini, factorPSet (StoredFactors.FS1));
461      sortCFListByLevel (ini);
462    }
463    else
464    {
465      ini= factorPSet (StoredFactors.FS1);
466      sortCFListByLevel (ini);
467    }
468
469    tmp2= adjoin (ini, l, qqi);
470    tmp= Union (tmp2, tmp);
471
472    StoredFactors.FS1= CFList();
473    StoredFactors.FS2= CFList();
474
475    ppi1= ListCFList();
476    ppi2= ListCFList();
477
478    count++;
479  }
480
481  //TODO need to remove superflous components
482
483  return result;
484}
485
486static bool
487irreducible (const CFList & AS)
488{
489// AS is given by AS = { A1, A2, .. Ar }, d_i = degree(Ai)
490
491// 1) we test: if d_i > 1, d_j =1 for all j<>i, then AS is irreducible.
492  bool deg1= true;
493  for (CFListIterator i= AS ; i.hasItem(); i++)
494  {
495    if (degree (i.getItem()) > 1)
496    {
497      if (deg1)
498        deg1= false;
499      else
500        return false; // found 2nd poly with deg > 1
501    }
502  }
503  return true;
504}
505
506static CFList
507irredAS (CFList & AS, int & indexRed, CanonicalForm & reducible)
508{
509  CFFList qs;
510  CFList ts, as;
511  CanonicalForm elem;
512  bool ind= true;
513  int nr= 0;
514  CFListIterator i;
515
516  indexRed= 0;
517  for (i= AS; i.hasItem(); i++ )
518  {
519    nr += 1;
520    qs= factorize (i.getItem());
521    if (qs.getFirst().factor().inCoeffDomain())
522      qs.removeFirst();
523
524    if ((qs.length() >= 2 ) || (qs.getFirst().exp() > 1))
525    {
526      indexRed= nr;
527      ind= false;
528      reducible= i.getItem();
529      break;
530    }
531  }
532
533  if (ind)
534  {
535    if (irreducible (AS))  // as quasilinear? => irreducible!
536      indexRed= 0;
537    else
538    {
539      i= AS;
540      for (nr= 1; nr< AS.length(); nr++)
541      {
542        as.append (i.getItem());
543        i++;
544        if (degree (i.getItem()) > 1)
545        {  // search for a non linear elem
546          qs= facAlgFunc2 (i.getItem(), as);
547          if (qs.length() > 0)
548          {
549            if (qs.getFirst().factor().inCoeffDomain())
550              qs.removeFirst();
551            if (qs.length() > 1 || qs.getFirst().exp() > 1)
552            { //found elem is reducible
553              reducible= i.getItem();
554              indexRed= nr + 1;
555              break;
556            }
557          }
558        }
559      }
560    }
561  }
562  for (CFFListIterator k= qs; k.hasItem(); k++)
563    ts.append (normalize (k.getItem().factor()));
564  return ts;
565}
566
567ListCFList
568irrCharSeries (const CFList & PS)
569{
570  CanonicalForm reducible, reducible2;
571  CFList qs, cs, factorset, is, ts, L;
572  CanonicalForm sqrf;
573  CFFList sqrfFactors;
574  CFFListIterator iter2;
575  for (CFListIterator iter= PS; iter.hasItem(); iter++)
576  {
577    sqrf= 1;
578    sqrfFactors= sqrFree (iter.getItem());
579    if (sqrfFactors.getFirst().factor().inCoeffDomain())
580      sqrfFactors.removeFirst();
581    for (iter2= sqrfFactors; iter2.hasItem(); iter2++)
582      sqrf *= iter2.getItem().factor();
583    sqrf= normalize (sqrf);
584    L= Union (CFList (sqrf), L);
585  }
586
587  ListCFList pi, ppi, qqi, qsi, iss, qhi= ListCFList(L);
588
589  int nr_of_iteration= 0, indexRed, highestlevel= 0;
590
591  for (CFListIterator iter= PS; iter.hasItem(); iter++)
592  {
593    if (level (iter.getItem()) > highestlevel)
594      highestlevel= level(iter.getItem());
595  }
596
597  while (!qhi.isEmpty())
598  {
599    sortListCFList (qhi);
600
601    qs= qhi.getFirst();
602
603    ListCFList ppi1,ppi2;
604    select (ppi, qs.length(), ppi1, ppi2);
605
606    inplaceUnion (ppi2, qqi);
607
608    if (nr_of_iteration == 0)
609    {
610      nr_of_iteration += 1;
611      ppi= ListCFList();
612    }
613    else
614    {
615      nr_of_iteration += 1;
616      ppi= Union (ppi1, ListCFList (qs));
617    }
618
619    StoreFactors StoredFactors;
620    if (qs.length() - 3 < highestlevel)
621      cs= modCharSet (qs, StoredFactors, false);
622    else
623      cs= charSetN (qs);
624    cs= removeContent (cs, StoredFactors);
625
626    factorset= StoredFactors.FS1;
627
628    if (!cs.isEmpty() && cs.getFirst().level() > 0)
629    {
630      ts= irredAS (cs, indexRed, reducible);
631
632      if (indexRed <= 0) // irreducible
633      {
634        if (!isSubset (cs,qs))
635          cs= charSetViaCharSetN (Union (qs,cs));
636        if (!find (pi, cs))
637        {
638          pi= Union (ListCFList (cs), pi);
639          if (cs.getFirst().level() > 0)
640          {
641            ts= irredAS (cs, indexRed, reducible);
642
643            if (indexRed <= 0) //irreducible
644            {
645              qsi= Union (ListCFList(cs), qsi);
646              if (cs.length() == highestlevel)
647                is= factorPSet (factorset);
648              else
649                is= Union (factorsOfInitials (cs), factorPSet (factorset));
650              iss= adjoin (is, qs, qqi);
651            }
652          }
653          else
654            iss= adjoin (factorPSet (factorset), qs, qqi);
655        }
656        else
657          iss= adjoin (factorPSet (factorset), qs, qqi);
658      }
659
660      if (indexRed > 0)
661      {
662        is= factorPSet (factorset);
663        if (indexRed > 1)
664        {
665          CFList cst;
666          for (CFListIterator i= cs ; i.hasItem(); i++)
667          {
668            if (i.getItem() == reducible)
669              break;
670            else
671              cst.append (i.getItem());
672          }
673          is= Union (factorsOfInitials (Union (cst,  CFList (reducible))), is);
674          iss= Union (adjoinb (ts, qs, qqi, cst), adjoin (is, qs, qqi));
675        }
676        else
677          iss= adjoin (Union (is, ts), qs, qqi);
678      }
679    }
680    else
681      iss= adjoin (factorPSet (factorset), qs, qqi);
682    if (qhi.length() > 1)
683    {
684      qhi.removeFirst();
685      qhi= Union (iss, qhi);
686    }
687    else
688      qhi= iss;
689  }
690  if (!qsi.isEmpty())
691    return contract (qsi);
692  return ListCFList(CFList (1)) ;
693}
694
Note: See TracBrowser for help on using the repository browser.