source: git/libfac/charset/charset.cc @ a13956

fieker-DuValspielwiese
Last change on this file since a13956 was a13956, checked in by Hans Schoenemann <hannes@…>, 14 years ago
factoryError, removed dep. on SINGULAR git-svn-id: file:///usr/local/Singular/svn/trunk@13616 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 16.6 KB
Line 
1////////////////////////////////////////////////////////////
2// emacs edit mode for this file is -*- C++ -*-
3/* $Id$ */
4/////////////////////////////////////////////////////////////
5// FACTORY - Includes
6#include <factory.h>
7#ifdef HAVE_IOSTREAM
8#define CIN std::cin
9#else
10#define CIN cin
11#endif
12// Factor - Includes
13#include <SqrFree.h>
14#include <Factor.h>
15#include <interrupt.h>
16// Charset - Includes
17#include "csutil.h"
18#include "algfactor.h"
19#include "alg_factor.h"
20// Some CC's need this:
21#include "charset.h"
22
23#ifdef BASICSETDEBUG
24#  define DEBUGOUTPUT
25#else
26#  undef DEBUGOUTPUT
27#endif
28
29#include "debug.h"
30#include "timing.h"
31TIMING_DEFINE_PRINT(subfactorize_time);
32
33
34// forward declarations:
35static CFList     irras(CFList & AS, int &ja, CanonicalForm & reducible);
36
37#ifdef BASICSETDEBUG
38#  define DEBUGOUTPUT
39#else
40#  undef DEBUGOUTPUT
41#endif
42#include "debug.h"
43
44// the next computes a characteristic set (a basic set in Wang's sense)
45CFList
46BasicSet( const CFList &PS )
47{
48    CFList QS = PS, BS, RS;
49    CanonicalForm b;
50    int cb;
51
52    DEBOUTLN(CERR, "BasicSet: called with ps= ", PS);
53    if ( PS.length() < 2 ) return PS;
54    while ( ! QS.isEmpty() ) {
55        b = lowestRank( QS );
56        cb = rank( b );
57        DEBOUTLN(CERR, "BasicSet: choose b  = ", b);
58        DEBOUTLN(CERR, "BasicSet: it's rank = ", cb);
59        BS=Union(CFList(b),BS);//BS.append( b );
60        if ( rank( b ) == 0 )
61            return Union(PS, CFList(b)) ; // b should be the first elem!
62        else {
63            RS = CFList();
64            // QS:= {q \in QS -{B} | q is reduced wrt b}
65            // We can process whole QS, because b is not reduced wrt. b
66            for ( CFListIterator i = QS; i.hasItem(); ++i )
67                if ( degree( i.getItem(), cb ) < degree( b ) )
68                    //RS.append( i.getItem() );
69                    RS = Union(CFList(i.getItem()),RS);
70            QS = RS;
71        }
72    }
73    DEBOUTLN(CERR, "BasicSet: returning bs= ", BS);
74    return BS;
75}
76
77int
78checkok( const CFList & PS, CFList & FS2){
79  CanonicalForm elem;
80
81  for ( CFListIterator i=PS; i.hasItem(); i++){
82    elem= i.getItem();
83    for (CFListIterator j=FS2; j.hasItem(); j++){
84      if (elem == j.getItem()){
85         //        FS2= Difference(FS2,CFList(elem));
86         return 0;
87      }
88    }
89  }
90  return 1;
91}
92
93#ifdef MCHARSETNDEBUG
94#  define DEBUGOUTPUT
95#else
96#  undef DEBUGOUTPUT
97#endif
98#include "debug.h"
99
100// The modified CharSet (an extended characteristic set with certain factors
101// canceled; this is a characteristic set in Wang's sense)
102CFList
103MCharSetN( const CFList &PS, PremForm & Remembern ){
104  CFList QS = PS, RS = PS, CSet, OLDCS;
105
106  DEBOUTLN(CERR, "MCharSetN: called with ps= ", PS);
107  while ( ! RS.isEmpty() ) {
108    CSet = BasicSet( QS );
109    OLDCS=CSet;
110    DEBOUTLN(CERR, "MCharSetN: CS= ", CSet);
111//     if ( getNumVars(CSet.getFirst()) > 1 ){
112//       //CSet = removecontent(CSet, Remembern);
113// #ifdef MCHARSETNDEBUG
114//       CERR << "MCharSetN: CSet= " << CSet << "\n";
115// #endif
116//     }
117    Remembern.FS1 = Union(Remembern.FS1, initalset1(CSet));
118    DEBOUTLN(CERR, "MCharSetN: Remembern.FS1= ", Remembern.FS1);
119    DEBOUTLN(CERR, "MCharSetN: Remembern.FS2= ", Remembern.FS2);
120    RS = CFList();
121    if ( rank( CSet.getFirst() ) != 0 ) {
122      CFList D = Difference( QS, CSet );
123      DEBOUT(CERR, "MCharSetN: Difference( ", QS);
124      DEBOUT(CERR, " , ", CSet);
125      DEBOUTLN(CERR, " ) = ", D);
126//CERR << "MCharSetN: Difference( " << QS << " , " << CSet << " ) = " << D << "\n";
127      //PremForm Oldremember=Remembern;
128      //PremForm Newremember=Remembern;
129      for ( CFListIterator i = D; i.hasItem(); ++i ) {
130        CanonicalForm r = Prem( i.getItem(), CSet );
131        DEBOUT(CERR,"MCharSetN: Prem(", i.getItem()  );
132        DEBOUT(CERR, ",", CSet);
133        DEBOUTLN(CERR,") = ", r);
134//CERR << "MCharSetN: Prem("<< i.getItem() << "," << CSet << ") = " << r << "\n";
135        if ( r != 0 ){
136          //removefactor( r, Newremember );
137          removefactor( r, Remembern );
138          //Remembern.FS2 = Union(Remembern.FS2, Newremember.FS2);
139          //Newremember = Oldremember;
140          //if ( cls(r) > 0 )
141            //RS=Union(CFList(r),RS);//RS.append( r );
142            RS=Union(RS,CFList(r));
143        }
144      }
145      if ( ! checkok(RS,Remembern.FS2)) return CFList(CanonicalForm(1));
146      DEBOUTLN(CERR, "MCharSetN: RS= ", RS);
147      //QS = Union( QS, RS );
148      QS = Union(OLDCS,RS);
149      DEBOUTLN(CERR, "MCharSetN: QS= Union(QS,RS)= ", QS);
150    }
151    else{ return CFList(CanonicalForm(1)); }
152  }
153  DEBOUTLN(CERR, "MCharSetN: Removed factors: ", Remembern.FS2);
154  DEBOUTLN(CERR, "MCharSetN: Remembern.FS1: ", Remembern.FS1);
155
156  return CSet;
157}
158
159
160CFList
161mcharset( const CFList &PS, PremForm & Remembern ){
162  CFList cs= MCharSetN(PS, Remembern );
163  CFList rs= remsetb(Difference(PS,cs),cs);
164
165  DEBOUTLN(CERR, "mcharset: cs= ", cs);
166  DEBOUTLN(CERR, "mcharset: rs= ", rs);
167  if ( rs.length() > 0 )
168    cs= mcharset(Union(PS,Union(cs,rs)), Remembern);
169  return cs;
170}
171
172// the "original" extended characteristic set
173CFList
174CharSet( const CFList &PS ){
175  CFList QS = PS, RS = PS, CSet;
176
177  while ( ! RS.isEmpty() ) {
178    CSet = BasicSet( QS );
179    DEBOUTLN(CERR, "CharSet: CSet= ", CSet);
180    RS = CFList();
181    if ( rank( CSet.getFirst() ) != 0 ) {
182      CFList D = Difference( QS, CSet );
183      for ( CFListIterator i = D; i.hasItem(); ++i ) {
184        CanonicalForm r = Prem( i.getItem(), CSet );
185        if ( r != 0 )  RS=Union(CFList(r),RS);//RS.append( r );
186      }
187      QS = Union( QS, RS );
188    }
189  }
190  return CSet;
191}
192
193static CFList
194charseta( const CFList & PS ){
195  CFList QS = PS, RS = PS, CSet;
196
197  while ( ! RS.isEmpty() ) {
198    CSet = CharSet( QS );
199    RS = CFList();
200    if ( rank( CSet.getFirst() ) != 0 ) {
201      CFList D = Difference( QS, CSet );
202      for ( CFListIterator i = D; i.hasItem(); ++i ) {
203        CanonicalForm r = Prem( i.getItem(), CSet );
204        if ( r != 0 )  RS=Union(CFList(r),RS);//RS.append( r );
205      }
206      QS = Union(CSet,Union( QS, RS ));
207    }
208    else return CFList(CanonicalForm(1));
209  }
210  return CSet;
211}
212
213static bool
214contractsub( const CFList & cs1, const CFList & cs2 ){
215  CFListIterator i;
216
217  for ( i=cs1; i.hasItem(); i++ )
218    if ( Prem(i.getItem(),cs2 ) != 0 ) return false;
219  CFList is=initalset1(cs1);
220  for ( i=is; i.hasItem(); i++ )
221    if ( Prem(i.getItem(),cs2 ) == 0 ) return false;
222  return true;
223}
224
225static ListCFList
226contract( const ListCFList & cs){
227  ListCFList mem,ts;
228  CFList iitem,jitem;
229
230  if ( cs.length() < 2 ) return cs;
231
232  for ( ListCFListIterator i=cs; i.hasItem(); i++ ){
233    iitem=i.getItem();
234    if ( ! member(iitem, mem))
235      for ( ListCFListIterator j=i; j.hasItem(); j++){
236        jitem=j.getItem();
237        if ( ! same( iitem, jitem ) )
238          if ( ! member(jitem, mem))
239            if ( contractsub(iitem, jitem) ){
240              ts.append(jitem); mem.append(jitem);
241            }
242            else
243              if ( contractsub(jitem, iitem) ){
244                ts.append(iitem);
245              }
246      }
247  }
248  return Minus(cs,ts);
249}
250
251static ListCFList
252adjoin(const CFList & is, const CFList & qs, const ListCFList & qh ){
253  ListCFList iss,qhi;
254  ListCFListIterator j;
255  CFList iscopy,itt;
256  CFListIterator i;
257  CanonicalForm elem;
258  int ind, length;
259
260  for ( i=is ; i.hasItem(); i++ ){
261    elem=i.getItem();
262    if ( cls(elem) > 0 ) iscopy=Union(CFList(elem),iscopy);
263  }
264  if ( iscopy.isEmpty() ) return iss;
265  qhi = MyDifference(qh,qs);
266  length = qhi.length();
267  for ( i =iscopy; i.hasItem(); i++){
268    itt = Union(qs,CFList(i.getItem()));
269    ind = 0;
270    if ( length > 0 )
271      for ( j=qhi; j.hasItem(); j++ )
272        if ( subset(j.getItem(),itt )) ind=1;
273    if ( ind == 0 ) iss.append(itt);
274  }
275  return iss;
276}
277
278static ListCFList
279adjoinb(const CFList & is, const CFList & qs, const ListCFList & qh ,const CFList & cs){
280  ListCFList iss,qhi;
281  ListCFListIterator j;
282  CFList iscopy,itt;
283  CFListIterator i;
284  CanonicalForm elem;
285  int ind, length;
286
287  for ( i=is ; i.hasItem(); i++ ){
288    elem=i.getItem();
289    if ( cls(elem) > 0 ) iscopy=Union(CFList(elem),iscopy);
290  }
291  if ( iscopy.isEmpty() ) return iss;
292  qhi = MyDifference(qh,qs);
293  length = qhi.length();
294  for ( i =iscopy; i.hasItem(); i++){
295    itt = Union(Union(qs,CFList(i.getItem())),cs);
296    ind = 0;
297    if ( length > 0 )
298      for ( j=qhi; j.hasItem(); j++ )
299        if ( subset(j.getItem(),itt )) ind=1;
300    if ( ind == 0 ) {iss.append(itt);}
301  }
302  return iss;
303}
304
305static ListCFList
306sort( const ListCFList & list_to_sort ){
307  ListCFList output,copy=list_to_sort;
308  CFList l,qs1,elem;
309
310  l = copy.getLast(); copy.removeLast();
311  if ( copy.length() == 0 ){ return ListCFList(l); }
312  for ( ListCFListIterator i=copy ; i.hasItem(); i++ ){
313    elem = i.getItem();
314    if ( elem.length() > l.length() ) {
315      output = MyUnion( ListCFList(l), output);
316      l= elem;
317    }
318    else{ output = MyUnion(output, ListCFList(elem) ); }
319  }
320  //output = MyUnion(ListCFList(l),sort(output));
321  output = MyUnion(ListCFList(l),output);
322  return output;
323}
324
325#ifdef EXPERIMENTAL
326static CFList
327getItemNr( int nr, const ListCFList & copy){
328  int i =1;
329  CFList elem;
330
331  for ( ListCFListIterator j=copy; j.hasItem(); j++ )
332    if ( i == nr ) { elem=j.getItem(); break; }
333    else { i+= 1; }
334  return elem;
335}
336
337static int
338choosefrom(){
339int choice;
340    CERR << "choose from qhi! ->";
341    CIN >> choice;
342return choice;
343}
344
345static ListCFList
346msort( const ListCFList & list_to_sort ){
347  int nr, number = list_to_sort.length();
348  ListCFList output;
349
350  CERR << "Sort: list to sort is: " <<  list_to_sort << "\n";
351  for (int i=1; i<= number; i++){
352    CERR << " Next elem = "; CIN >> nr;
353    output.append(getItemNr(nr,list_to_sort));
354  }
355  return output;
356}
357#endif
358
359#ifdef IRRCHARSERIESDEBUG
360#  define DEBUGOUTPUT
361#else
362#  undef DEBUGOUTPUT
363#endif
364#include "debug.h"
365
366ListCFList
367IrrCharSeries( const CFList &PS, int opt ){
368  CanonicalForm reducible,reducible2;
369  CFList qs,cs,factorset,is,ts;
370  ListCFList pi,ppi,qqi,qsi,iss,qhi= ListCFList(PS);
371  int nr_of_iteration=0,ts2,highestlevel=0;
372#ifdef EXPERIMENTAL
373  int choice=1;;
374#endif
375
376//  CERR << getCharacteristic() << "\n";
377  for ( CFListIterator Ps=PS; Ps.hasItem(); Ps++ )
378    if ( level(Ps.getItem() ) > highestlevel ) highestlevel = level(Ps.getItem()) ;
379//  for ( int xx=1; xx <= highestlevel; xx++)
380//   CERR << Variable(xx) ;
381//  CERR << "\n";
382//  for ( CFListIterator Ps=PS; Ps.hasItem(); Ps++ )
383//    CERR << Ps.getItem() << ", " ;//<< "\n";
384//  CERR <<  "\n";
385  while ( ! qhi.isEmpty() ) {
386    qhi=sort(qhi);
387    DEBOUTLN(CERR, "qhi is: ", qhi);
388#ifdef EXPERIMENTAL
389    choice=choosefrom();
390    CERR <<"/n Choose " << choice << "\n";
391    qs= getItemNr(choice, qhi);
392#else
393    qs=qhi.getFirst();
394#endif
395    DEBOUTLN(CERR, "qs  is: ", qs);
396    DEBOUTLN(CERR, "ppi is: ", ppi);
397    ListCFList ppi1,ppi2;
398    select(ppi,qs.length(),ppi1,ppi2);
399    DEBOUTLN(CERR, "ppi1 is: ", ppi1);
400    DEBOUTLN(CERR, "ppi2 is: ", ppi2);
401    qqi = MyUnion(ppi2,qqi);
402    DEBOUTLN(CERR, "qqi is: ", qqi);
403    if ( nr_of_iteration == 0 ){ nr_of_iteration += 1; ppi = ListCFList(); }
404    else{ nr_of_iteration += 1; ppi = MyUnion(ListCFList(qs),ppi1); }
405    DEBOUTLN(CERR,"ppi is: ", ppi);
406    PremForm Remembern;
407    cs = MCharSetN(qs,Remembern);
408    DEBOUTLN(CERR, "cs is: ", cs);
409    DEBOUTLN(CERR, "factorset is: ", Remembern.FS2);
410    cs = removecontent(cs,Remembern);
411    factorset=Remembern.FS2;
412    DEBOUTLN(CERR, "cs (after removecontent) is: ", cs);
413    DEBOUTLN(CERR, "factorset is: ", factorset);
414    // Hier: removecontent einfuegen!!!!
415    if ( cls(cs.getFirst()) > 0 ){
416      ts = irras(cs,ts2, reducible);
417
418      // INTERRUPTHANDLER
419      if ( interrupt_handle() ) return ListCFList() ;
420      // INTERRUPTHANDLER
421
422      DEBOUTLN(CERR, "ts is: ", ts);
423      DEBOUTLN(CERR, "ts2 is: ", ts2);
424      // next is preliminary: should be ==0
425      if ( ts2 <= 0 ){ //irreducible
426        if ( ! subset(cs,qs) ){
427          DEBOUTMSG(CERR, "cs is not a subset of qs");
428          cs = charseta(Union(qs,cs));
429          DEBOUTLN(CERR, "new cs is: ", cs);
430        }
431        if ( ! member(cs,pi) ){
432          pi = MyUnion(pi, ListCFList(cs));
433          DEBOUTMSG(CERR, "cs is not a member of pi");
434          DEBOUTLN(CERR, "pi is: ", pi);
435          if ( cls(cs.getFirst()) > 0 ){
436            ts = irras(cs,ts2,reducible);
437
438            // INTERRUPTHANDLER
439            if ( interrupt_handle() ) return ListCFList() ;
440            // INTERRUPTHANDLER
441
442            DEBOUTLN(CERR, "ts is: ", ts);
443            DEBOUTLN(CERR, "ts2 is: ", ts2);
444            // next is preliminary: should be ==0
445            if ( ts2 <= 0 ){ //irreducible
446              qsi = MyUnion(qsi,ListCFList(cs));
447              DEBOUTLN(CERR, "qsi is: ", qsi);
448              if ( cs.length() == highestlevel ){
449                DEBOUTLN(CERR, "cs.length() == nops(ord) :", cs.length());
450                is = factorps(factorset);
451              }
452              else{
453                DEBOUT(CERR,"cs.length() != nops(ord) :", cs.length());
454                DEBOUTLN(CERR, "  nops(ord)= ", highestlevel);
455                is = Union(initalset1(cs),factorps(factorset));
456              }
457              DEBOUTLN(CERR, "is is: ", is);
458              iss = adjoin(is,qs,qqi);
459              DEBOUTLN(CERR, "iss is: ", iss);
460            }
461          }
462          else{ iss = adjoin(factorps(factorset),qs,qqi); }
463        }
464        else{
465          DEBOUTMSG(CERR, "cs is a member of pi");
466          iss = adjoin(factorps(factorset),qs,qqi); }
467        DEBOUTLN(CERR, "iss is: ", iss);
468        DEBOUTLN(CERR, "   factorps(factorset)= ", factorps(factorset));
469        DEBOUTLN(CERR, "   qs= ", qs);
470        DEBOUTLN(CERR, "   qqi= ", qqi);
471      }
472      // next is preliminary: should be !=0
473      if ( ts2 > 0 ){
474        is = factorps(factorset);
475        DEBOUTLN(CERR, "is is: ", is);
476        if ( ts2 > 1 ){
477          // setup cst: need it later for adjoinb
478          CFList cst;
479          for ( CFListIterator i=cs ; i.hasItem(); i++){
480            if ( i.getItem() == reducible ) { break; }
481            else { cst.append(i.getItem()); }
482          }
483          is = Union(initalset1(cst), is);
484          iss = MyUnion(adjoin(is,qs,qqi), adjoinb(ts,qs,qqi,cst));
485        }
486        else{ iss = adjoin(Union(is,ts),qs,qqi); }
487        DEBOUTLN(CERR, "iss is: ", iss);
488      }
489    }
490    else{
491      iss = adjoin(factorps(factorset),qs,qqi);
492      DEBOUTMSG(CERR, "case: cs is a constant.");
493      DEBOUTLN(CERR, "  qs = ", qs);
494      DEBOUTLN(CERR, "  qqi = ", qqi);
495      DEBOUTLN(CERR, "  iss = adjoin(factorps(factorset),qs,qqi) = ",iss);
496    }
497    if ( qhi.length() > 1 ){ qhi.removeFirst(); qhi = MyUnion(qhi,iss); }
498    else{ qhi = iss; }
499    DEBOUTLN(CERR, "iss is: ", iss);
500  }
501  if ( ! qsi.isEmpty() ){
502    DEBOUTLN(CERR, "qsi before contract= ", qsi);
503    if ( opt == 0 ){
504       return contract( qsi );
505    }
506    else { return qsi; }
507  }
508  else{ return ListCFList() ; }
509}
510
511// tests for characteristic sets
512//////////////////////////////////
513#ifdef IRRASDEBUG
514#  define DEBUGOUTPUT
515#else
516#  undef DEBUGOUTPUT
517#endif
518
519#include "debug.h"
520
521static CFList
522irras( CFList & AS, int & ja, CanonicalForm & reducible)
523{
524  CFFList qs;
525  CFList ts,as;
526  CanonicalForm elem;
527  int ind=1,nr=0, success=-1;
528  CFListIterator i;
529
530  ja = 0;
531  DEBOUTLN(CERR, "irras: called with: AS= ", AS);
532  for ( i=AS; i.hasItem(); i++ )
533  {
534    elem = i.getItem();
535    nr += 1;
536    DEBOUT(CERR, "irras: factoring: ", elem);
537    if ( degree(elem) > 1 ) // linear poly's are irreduzible
538    {
539      qs = Factorize(elem);
540      // remove a constant
541      if (qs.getFirst().factor().degree()==0) qs.removeFirst();
542    }
543    else
544    {
545      qs=(CFFactor(elem,1));
546    }
547    DEBOUTLN(CERR, "  = ", qs);
548    // INTERRUPTHANDLER
549    if ( interrupt_handle() ) return CFList() ;
550    // INTERRUPTHANDLER
551
552    if ( (qs.length() >= 2 ) || (qs.getFirst().exp() > 1))
553    {
554      DEBOUTLN(CERR, "irras: Setting ind=0, ja= ", nr);
555      ja=nr; ind=0; reducible= elem;
556      break;
557    }
558    //    else{ as.append(elem) ; }
559  }
560  //  CERR << "ind= " << ind << "\n";
561  if ( (ind == 1) ) //&& ( as.length() > 1) )
562  {
563    if ( irreducible(AS) )
564    { // as quasilinear? => irreducible!
565      ja = 0;
566      DEBOUTLN(CERR, "as is irreducible. as= ", AS);
567    }
568    else
569    {
570      i=AS;
571      for ( nr=1; nr< AS.length(); nr++)
572      {
573        as.append(i.getItem());
574        i++;
575        if ( degree(i.getItem()) > 1 )
576        {  // search for a non linear elem
577          elem=i.getItem();
578          //if (as.length()==1)
579          //  qs = Factorize2(elem,as.getFirst());
580          //else
581            qs= newfactoras(elem,as,success);
582          if ( qs.length() > 1 || qs.getFirst().exp() > 1 )
583          { //found elem is reducible
584            reducible=elem;
585            ja=nr+1;
586            break;
587          }
588        }
589      }
590    }
591  }
592  for ( CFFListIterator k=qs; k.hasItem();k++)
593    ts.append(k.getItem().factor());
594  return ts;
595}
596
597///////////////////////////////////////////////////////////////////////////////
Note: See TracBrowser for help on using the repository browser.