source: git/dyn_modules/syzextra/mod_main.cc @ 74afe1f

spielwiese
Last change on this file since 74afe1f was 74afe1f, checked in by Oleksandr Motsak <motsak@…>, 12 years ago
moved 'Compute2LeadingSyzygyTerms' into the dyn.mod. TODO/TEST: interred & tailred
  • Property mode set to 100644
File size: 39.7 KB
Line 
1
2
3
4
5#include <kernel/mod2.h>
6
7#include <omalloc/omalloc.h>
8
9#include <misc/intvec.h>
10#include <misc/options.h>
11
12#include <coeffs/coeffs.h>
13
14#include <polys/PolyEnumerator.h>
15
16#include <polys/monomials/p_polys.h>
17#include <polys/monomials/ring.h>
18// #include <kernel/longrat.h>
19#include <kernel/GBEngine/kstd1.h>
20
21#include <kernel/polys.h>
22
23#include <kernel/GBEngine/syz.h>
24
25#include <Singular/tok.h>
26#include <Singular/ipid.h>
27#include <Singular/lists.h>
28#include <Singular/attrib.h>
29
30#include <Singular/ipid.h> 
31#include <Singular/ipshell.h> // For iiAddCproc
32
33// extern coeffs coeffs_BIGINT
34
35#include "singularxx_defs.h"
36#include "DebugPrint.h"
37#include "myNF.h"
38
39
40#include <Singular/mod_lib.h>
41
42
43#if GOOGLE_PROFILE_ENABLED
44#include <google/profiler.h>
45#endif // #if GOOGLE_PROFILE_ENABLED
46
47
48#include <stdio.h>
49#include <stdlib.h>
50#include <string.h>
51
52
53extern void pISUpdateComponents(ideal F, const intvec *const V, const int MIN, const ring r);
54// extern ring rCurrRingAssure_SyzComp();
55extern ring rAssure_InducedSchreyerOrdering(const ring r, BOOLEAN complete, int sign);
56extern int rGetISPos(const int p, const ring r);
57
58USING_NAMESPACE( SINGULARXXNAME :: DEBUG )
59USING_NAMESPACE( SINGULARXXNAME :: NF )
60
61
62BEGIN_NAMESPACE_NONAME
63
64
65static inline void NoReturn(leftv& res)
66{
67  res->rtyp = NONE;
68  res->data = NULL;
69}
70
71/// wrapper around n_ClearContent
72static BOOLEAN _ClearContent(leftv res, leftv h)
73{
74  NoReturn(res);
75
76  const char *usage = "'ClearContent' needs a (non-zero!) poly or vector argument...";
77 
78  if( h == NULL )
79  {
80    WarnS(usage);
81    return TRUE;
82  }
83
84  assume( h != NULL );
85
86  if( !( h->Typ() == POLY_CMD || h->Typ() == VECTOR_CMD) )
87  {
88    WarnS(usage);
89    return TRUE;
90  }
91
92  assume (h->Next() == NULL);
93 
94  poly ph = reinterpret_cast<poly>(h->Data());
95 
96  if( ph == NULL )
97  {
98    WarnS(usage);
99    return TRUE;
100  }
101 
102  const ring r =  currRing;
103  assume( r != NULL ); assume( r->cf != NULL ); const coeffs C = r->cf;
104
105  number n;
106
107  // experimentall (recursive enumerator treatment) of alg. ext
108  CPolyCoeffsEnumerator itr(ph);
109  n_ClearContent(itr, n, C);
110
111  res->data = n;
112  res->rtyp = NUMBER_CMD;
113
114  return FALSE;
115}
116
117/// wrapper around n_ClearDenominators
118static BOOLEAN _ClearDenominators(leftv res, leftv h)
119{
120  NoReturn(res);
121
122  const char *usage = "'ClearDenominators' needs a (non-zero!) poly or vector argument...";
123
124  if( h == NULL )
125  {
126    WarnS(usage);
127    return TRUE;
128  }
129
130  assume( h != NULL );
131
132  if( !( h->Typ() == POLY_CMD || h->Typ() == VECTOR_CMD) )
133  {
134    WarnS(usage);
135    return TRUE;
136  }
137
138  assume (h->Next() == NULL);
139
140  poly ph = reinterpret_cast<poly>(h->Data());
141
142  if( ph == NULL )
143  {
144    WarnS(usage);
145    return TRUE;
146  }
147
148  const ring r =  currRing;
149  assume( r != NULL ); assume( r->cf != NULL ); const coeffs C = r->cf;
150
151  number n;
152
153  // experimentall (recursive enumerator treatment) of alg. ext.
154  CPolyCoeffsEnumerator itr(ph);
155  n_ClearDenominators(itr, n, C);
156
157  res->data = n;
158  res->rtyp = NUMBER_CMD;
159
160  return FALSE;
161}
162
163
164/// try to get an optional (simple) integer argument out of h
165/// or return the default value
166static int getOptionalInteger(const leftv& h, const int _n)
167{
168  if( h!= NULL && h->Typ() == INT_CMD )
169  {
170    int n = (int)(long)(h->Data());
171
172    if( n < 0 )
173      Warn("Negative (%d) optional integer argument", n);
174
175    return (n);
176  }
177
178  return (_n); 
179}
180
181static BOOLEAN noop(leftv __res, leftv /*__v*/)
182{
183  NoReturn(__res);
184  return FALSE;
185}
186
187static BOOLEAN _ProfilerStart(leftv __res, leftv h)
188{
189  NoReturn(__res);
190#if GOOGLE_PROFILE_ENABLED
191  if( h!= NULL && h->Typ() == STRING_CMD )
192  {
193    const char* name = (char*)(h->Data());
194    assume( name != NULL );   
195    ProfilerStart(name);
196  } else
197    WerrorS("ProfilerStart requires a string [name] argument"); 
198#else
199  WarnS("Sorry no google profiler support (GOOGLE_PROFILE_ENABLE!=1)...");
200//  return TRUE; // ?
201#endif // #if GOOGLE_PROFILE_ENABLED
202  return FALSE;
203  (void)h;
204}
205static BOOLEAN _ProfilerStop(leftv __res, leftv /*__v*/)
206{
207  NoReturn(__res);
208#if GOOGLE_PROFILE_ENABLED
209  ProfilerStop();
210#else
211  WarnS("Sorry no google profiler support (GOOGLE_PROFILE_ENABLED!=1)...");
212//  return TRUE; // ?
213#endif // #if GOOGLE_PROFILE_ENABLED
214  return FALSE;
215}
216
217static inline number jjLONG2N(long d)
218{
219  return n_Init(d, coeffs_BIGINT);
220}
221
222static inline void view(const intvec* v)
223{
224#ifndef SING_NDEBUG
225  v->view();
226#else
227  // This code duplication is only due to Hannes's #ifndef SING_NDEBUG!
228  Print ("intvec: {rows: %d, cols: %d, length: %d, Values: \n", v->rows(), v->cols(), v->length());
229
230  for (int i = 0; i < v->rows(); i++)
231  {
232    Print ("Row[%3d]:", i);
233    for (int j = 0; j < v->cols(); j++)
234      Print (" %5d", (*v)[j + i * (v->cols())] );
235    PrintLn ();
236  }
237  PrintS ("}\n");
238#endif
239
240}
241
242                   
243
244static BOOLEAN DetailedPrint(leftv __res, leftv h)
245{
246  NoReturn(__res);
247
248  if( h == NULL )
249  {
250    WarnS("DetailedPrint needs an argument...");
251    return TRUE;
252  }
253
254  if( h->Typ() == NUMBER_CMD)
255  {
256    number n = (number)h->Data(); 
257
258    const ring r = currRing;
259
260#ifdef LDEBUG
261    r->cf->cfDBTest(n,__FILE__,__LINE__,r->cf);
262#endif
263
264    StringSetS("");
265    n_Write(n, r->cf);
266    PrintS(StringEndS());
267    PrintLn();
268
269    return FALSE;
270  }
271 
272  if( h->Typ() == RING_CMD)
273  {
274    const ring r = (const ring)h->Data();
275    rWrite(r, TRUE);
276    PrintLn();
277#ifdef RDEBUG
278    rDebugPrint(r);
279#endif
280    return FALSE;
281  }
282
283  if( h->Typ() == POLY_CMD || h->Typ() == VECTOR_CMD)
284  {
285    const poly p = (const poly)h->Data(); h = h->Next();
286
287    dPrint(p, currRing, currRing, getOptionalInteger(h, 3));
288
289    return FALSE;
290  }
291
292  if( h->Typ() == IDEAL_CMD || h->Typ() == MODUL_CMD)
293  {
294    const ideal id = (const ideal)h->Data(); h = h->Next(); 
295
296    dPrint(id, currRing, currRing, getOptionalInteger(h, 3));
297   
298    return FALSE;           
299  }
300
301  if( h->Typ() == RESOLUTION_CMD )
302  {
303    const syStrategy syzstr = reinterpret_cast<const syStrategy>(h->Data());
304
305    h = h->Next();
306
307    int nTerms = getOptionalInteger(h, 1);
308
309
310    Print("RESOLUTION_CMD(%p): ", reinterpret_cast<const void*>(syzstr)); PrintLn();
311
312    const ring save = currRing;
313    const ring r = syzstr->syRing;
314    const ring rr = (r != NULL) ? r: save;
315
316
317    const int iLength = syzstr->length;
318
319    Print("int 'length': %d", iLength); PrintLn();
320    Print("int 'regularity': %d", syzstr->regularity); PrintLn();
321    Print("short 'list_length': %hd", syzstr->list_length); PrintLn();
322    Print("short 'references': %hd", syzstr->references); PrintLn();
323
324
325#define PRINT_pINTVECTOR(s, v) Print("intvec '%10s'(%p)", #v, reinterpret_cast<const void*>((s)->v)); \
326if( (s)->v != NULL ){ PrintS(": "); view((s)->v); }; \
327PrintLn();
328
329    PRINT_pINTVECTOR(syzstr, resolution);
330    PRINT_pINTVECTOR(syzstr, betti);
331    PRINT_pINTVECTOR(syzstr, Tl);
332    PRINT_pINTVECTOR(syzstr, cw);
333#undef PRINT_pINTVECTOR
334
335    if (r == NULL)
336      Print("ring '%10s': NULL", "syRing");
337    else 
338      if (r == currRing)
339        Print("ring '%10s': currRing", "syRing");
340      else
341        if (r != NULL && r != save)
342        {
343          Print("ring '%10s': ", "syRing");
344          rWrite(r);
345#ifdef RDEBUG
346          //              rDebugPrint(r);
347#endif
348          // rChangeCurrRing(r);
349        }           
350    PrintLn();
351
352    const SRes rP = syzstr->resPairs;
353    Print("SRes 'resPairs': %p", reinterpret_cast<const void*>(rP)); PrintLn();
354
355    if (rP != NULL)
356      for (int iLevel = 0; (iLevel < iLength) && (rP[iLevel] != NULL) && ((*syzstr->Tl)[iLevel] >= 0); iLevel++)
357      {
358        int n = 0;
359        const int iTl = (*syzstr->Tl)[iLevel];
360        for (int j = 0; (j < iTl) && ((rP[iLevel][j].lcm!=NULL) || (rP[iLevel][j].syz!=NULL)); j++)
361        {
362          if (rP[iLevel][j].isNotMinimal==NULL)
363            n++;
364        }
365        Print("minimal-resPairs-Size[1+%d]: %d", iLevel, n); PrintLn();
366      }
367
368
369    //  const ring rrr = (iLevel > 0) ? rr : save; ?
370#define PRINT_RESOLUTION(s, v) Print("resolution '%12s': %p", #v, reinterpret_cast<const void*>((s)->v)); PrintLn(); \
371if ((s)->v != NULL) \
372  for (int iLevel = 0; (iLevel < iLength) && ( ((s)->v)[iLevel] != NULL ); iLevel++) \
373  { \
374    /* const ring rrr = (iLevel > 0) ? save : save; */ \
375    Print("id '%10s'[%d]: (%p) ncols = %d / size: %d; nrows = %d, rank = %ld / rk: %ld", #v, iLevel, reinterpret_cast<const void*>(((s)->v)[iLevel]), ((s)->v)[iLevel]->ncols, idSize(((s)->v)[iLevel]), ((s)->v)[iLevel]->nrows, ((s)->v)[iLevel]->rank, -1L/*id_RankFreeModule(((s)->v)[iLevel], rrr)*/ ); \
376    PrintLn(); \
377  } \
378  PrintLn();
379
380    // resolvente:
381    PRINT_RESOLUTION(syzstr, minres);
382    PRINT_RESOLUTION(syzstr, fullres);
383
384    assume (id_RankFreeModule (syzstr->res[1], rr) == syzstr->res[1]->rank);
385
386    PRINT_RESOLUTION(syzstr, res);
387    PRINT_RESOLUTION(syzstr, orderedRes);
388#undef PRINT_RESOLUTION
389
390#define PRINT_POINTER(s, v) Print("pointer '%17s': %p", #v, reinterpret_cast<const void*>((s)->v)); PrintLn();
391    // 2d arrays:
392    PRINT_POINTER(syzstr, truecomponents);
393    PRINT_POINTER(syzstr, ShiftedComponents);
394    PRINT_POINTER(syzstr, backcomponents);
395    PRINT_POINTER(syzstr, Howmuch);
396    PRINT_POINTER(syzstr, Firstelem);
397    PRINT_POINTER(syzstr, elemLength);
398    PRINT_POINTER(syzstr, sev);
399
400    // arrays of intvects:
401    PRINT_POINTER(syzstr, weights);
402    PRINT_POINTER(syzstr, hilb_coeffs);
403#undef PRINT_POINTER
404
405
406    if (syzstr->fullres==NULL)
407    {
408      PrintS("resolution 'fullres': (NULL) => resolution not computed yet");
409      PrintLn();
410    } else
411    {
412      Print("resolution 'fullres': (%p) => resolution seems to be computed already", reinterpret_cast<const void*>(syzstr->fullres));
413      PrintLn();
414      dPrint(*syzstr->fullres, save, save, nTerms);
415    }
416
417
418
419
420    if (syzstr->minres==NULL)
421    {
422      PrintS("resolution 'minres': (NULL) => resolution not minimized yet");
423      PrintLn();
424    } else
425    {
426      Print("resolution 'minres': (%p) => resolution seems to be minimized already", reinterpret_cast<const void*>(syzstr->minres));
427      PrintLn();
428      dPrint(*syzstr->minres, save, save, nTerms);
429    }
430
431
432
433
434    /*
435    int ** truecomponents;
436    long** ShiftedComponents;
437    int ** backcomponents;
438    int ** Howmuch;
439    int ** Firstelem;
440    int ** elemLength;
441    unsigned long ** sev;
442
443    intvec ** weights;
444    intvec ** hilb_coeffs;
445
446    SRes resPairs;               //polynomial data for internal use only
447
448    resolvente fullres;
449    resolvente minres;
450    resolvente res;              //polynomial data for internal use only
451    resolvente orderedRes;       //polynomial data for internal use only
452*/
453
454    //            if( currRing != save ) rChangeCurrRing(save);
455  }
456
457
458  return FALSE;
459}
460
461
462
463
464
465
466static BOOLEAN Tail(leftv res, leftv h)
467{
468  NoReturn(res);
469
470  if( h == NULL )
471  {
472    WarnS("Tail needs an argument...");
473    return TRUE;
474  }
475
476  assume( h != NULL );
477
478  if( h->Typ() == POLY_CMD || h->Typ() == VECTOR_CMD)
479  {
480    const poly p = (const poly)h->Data();
481
482    res->rtyp = h->Typ();
483
484    if( p == NULL)
485      res->data = NULL;
486    else
487      res->data = p_Copy( pNext(p), currRing );
488
489    h = h->Next(); assume (h == NULL);
490   
491    return FALSE;
492  }
493
494  if( h->Typ() == IDEAL_CMD || h->Typ() == MODUL_CMD)
495  {
496    const ideal id = (const ideal)h->Data();
497
498    res->rtyp = h->Typ();
499
500    if( id == NULL)
501      res->data = NULL;
502    else
503    {
504      const ideal newid = idInit(IDELEMS(id),id->rank);
505      for (int i=IDELEMS(id) - 1; i >= 0; i--)
506        if( id->m[i] != NULL )
507          newid->m[i] = p_Copy(pNext(id->m[i]), currRing);
508        else
509          newid->m[i] = NULL;
510     
511      res->data = newid; 
512    }
513   
514    h = h->Next(); assume (h == NULL);
515   
516    return FALSE;
517  }
518
519  WarnS("Tail needs a single poly/vector/ideal/module argument...");
520  return TRUE;
521}
522
523
524static int cmp_c_ds(const void *p1, const void *p2, void *R)
525{
526  const int YES = 1;
527  const int NO = -1;
528
529  const ring r =  (const ring) R; // TODO/NOTE: the structure is known: C, lp!!!
530
531  const poly a = *(const poly*)p1;
532  const poly b = *(const poly*)p2;
533
534  assume( a != NULL );
535  assume( b != NULL );
536 
537  assume( p_LmTest(a, r) );
538  assume( p_LmTest(b, r) );
539
540
541  const signed long iCompDiff = p_GetComp(a, r) - p_GetComp(b, r);
542
543  // TODO: test this!!!!!!!!!!!!!!!!
544
545  //return -( compare (c, ds) )
546
547  const int __DEBUG__ = 0;
548
549#ifndef NDEBUG
550  if( __DEBUG__ )
551  {
552    PrintS("cmp_c_ds: a, b: \np1: "); dPrint(a, r, r, 2);
553    PrintS("b: "); dPrint(b, r, r, 2);
554    PrintLn();
555  }
556#endif
557
558
559  if( iCompDiff > 0 )
560    return YES;
561
562  if( iCompDiff < 0 )
563    return  NO;
564
565  assume( iCompDiff == 0 );
566
567  const signed long iDegDiff = p_Totaldegree(a, r) - p_Totaldegree(b, r);
568
569  if( iDegDiff > 0 )
570    return YES;
571
572  if( iDegDiff < 0 )
573    return  NO;
574
575  assume( iDegDiff == 0 );
576
577#ifndef NDEBUG
578  if( __DEBUG__ )
579  {
580    PrintS("   a & b have the same comp & deg! "); PrintLn();
581  }
582#endif
583
584  for (int v = rVar(r); v > 0; v--)
585  {
586    assume( v > 0 );
587    assume( v <= rVar(r) );
588
589    const signed int d = p_GetExp(a, v, r) - p_GetExp(b, v, r);
590
591    if( d > 0 )
592      return YES;
593
594    if( d < 0 )
595      return NO;
596
597    assume( d == 0 );
598  }
599
600  return 0; 
601}
602
603
604static BOOLEAN ComputeLeadingSyzygyTerms(leftv res, leftv h)
605{
606  const ring r = currRing;
607  NoReturn(res);
608
609  if( h == NULL )
610  {
611    WarnS("ComputeLeadingSyzygyTerms needs an argument...");
612    return TRUE;
613  }
614
615  assume( h != NULL ); 
616
617#ifndef NDEBUG
618  const BOOLEAN __DEBUG__ = (BOOLEAN)((long)(atGet(currRingHdl,"DEBUG",INT_CMD, (void*)TRUE)));
619#else
620  const BOOLEAN __DEBUG__ = (BOOLEAN)((long)(atGet(currRingHdl,"DEBUG",INT_CMD, (void*)FALSE)));
621#endif
622
623  if( h->Typ() == IDEAL_CMD || h->Typ() == MODUL_CMD)
624  {
625    const ideal id = (const ideal)h->Data();
626
627    assume(id != NULL);
628
629    if( __DEBUG__ )
630    {
631      PrintS("ComputeLeadingSyzygyTerms::Input: \n");
632     
633      const BOOLEAN __LEAD2SYZ__ = (BOOLEAN)((long)(atGet(currRingHdl,"LEAD2SYZ",INT_CMD, (void*)0)));
634      const BOOLEAN __TAILREDSYZ__ = (BOOLEAN)((long)(atGet(currRingHdl,"TAILREDSYZ",INT_CMD, (void*)0)));
635      const BOOLEAN __SYZCHECK__ = (BOOLEAN)((long)(atGet(currRingHdl,"SYZCHECK",INT_CMD, (void*)0)));
636
637
638      Print("\nSYZCHECK: \t%d", __SYZCHECK__);
639      Print(", DEBUG: \t%d", __DEBUG__);
640      Print(", LEAD2SYZ: \t%d", __LEAD2SYZ__);
641      Print(", TAILREDSYZ: \t%d\n", __TAILREDSYZ__);
642
643      dPrint(id, r, r, 1);
644
645      assume( !__LEAD2SYZ__ );
646    }
647
648    h = h->Next(); assume (h == NULL);
649
650    // 1. set of components S?
651    // 2. for each component c from S: set of indices of leading terms
652    // with this component?
653    // 3. short exp. vectors for each leading term?
654
655    const int size = IDELEMS(id);
656
657    if( size < 2 )
658    {
659      const ideal newid = idInit(1, 0); // zero ideal...
660
661      newid->m[0] = NULL;
662     
663      res->data = newid;
664      res->rtyp = MODUL_CMD;
665
666      return FALSE;
667    }
668   
669
670    // TODO/NOTE: input is supposed to be (reverse-) sorted wrt "(c,ds)"!??
671
672    // components should come in groups: count elements in each group
673    // && estimate the real size!!!
674
675
676    // use just a vector instead???
677    const ideal newid = idInit( (size * (size-1))/2, size); // maximal size: ideal case!
678   
679    int k = 0;
680
681    for (int j = 0; j < size; j++)
682    {
683      const poly p = id->m[j];
684      assume( p != NULL );
685      const int  c = p_GetComp(p, r);
686
687      for (int i = j - 1; i >= 0; i--)
688      {
689        const poly pp = id->m[i];
690        assume( pp != NULL );
691        const int  cc = p_GetComp(pp, r);
692
693        if( c != cc )
694          continue;
695
696        const poly m = p_Init(r); // p_New???
697
698        // m = LCM(p, pp) / p! // TODO: optimize: knowing the ring structure: (C/lp)!
699        for (int v = rVar(r); v > 0; v--)
700        {
701          assume( v > 0 );
702          assume( v <= rVar(r) );
703         
704          const short e1 = p_GetExp(p , v, r);
705          const short e2 = p_GetExp(pp, v, r);
706
707          if( e1 >= e2 )
708            p_SetExp(m, v, 0, r);
709          else
710            p_SetExp(m, v, e2 - e1, r);
711           
712        }
713
714        assume( (j > i) && (i >= 0) );
715
716        p_SetComp(m, j + 1, r);
717        pNext(m) = NULL;
718        p_SetCoeff0(m, n_Init(1, r->cf), r); // for later...
719
720        p_Setm(m, r); // should not do anything!!!
721
722        newid->m[k++] = m;
723      }
724    }
725
726    if( __DEBUG__ && FALSE )
727    {
728      PrintS("ComputeLeadingSyzygyTerms::Temp0: \n");
729      dPrint(newid, r, r, 1);
730    }
731
732    // the rest of newid is assumed to be zeroes...
733
734    // simplify(newid, 2 + 32)??
735    // sort(newid, "C,ds")[1]???
736    id_DelDiv(newid, r); // #define SIMPL_LMDIV 32
737
738    if( __DEBUG__ && FALSE )
739    {
740      PrintS("ComputeLeadingSyzygyTerms::Temp1: \n");
741      dPrint(newid, r, r, 1);
742    }
743
744    idSkipZeroes(newid); // #define SIMPL_NULL 2
745
746    if( __DEBUG__ )
747    {
748      PrintS("ComputeLeadingSyzygyTerms::Output: \n");
749      dPrint(newid, r, r, 1);
750    }
751
752    const int sizeNew = IDELEMS(newid);
753
754    if( sizeNew >= 2 )
755      qsort_r(newid->m, sizeNew, sizeof(poly), cmp_c_ds, r);
756
757    if (IDELEMS(newid) == 0 || (IDELEMS(newid) == 1 && newid->m[0] == NULL) )
758      newid->rank = 1;   
759   
760    // TODO: add sorting wrt <c,ds> & reversing...
761
762    res->data = newid;
763    res->rtyp = MODUL_CMD;
764   
765    return FALSE;
766  }
767
768  WarnS("ComputeLeadingSyzygyTerms needs a single ideal/module argument...");
769  return TRUE;
770}
771
772///  sorting wrt <c,ds> & reversing...
773/// change the input inplace!!!
774// TODO: use a ring with >_{c, ds}!???
775static BOOLEAN Sort_c_ds(leftv res, leftv h)
776{
777  NoReturn(res);
778
779  const ring r = currRing;
780  NoReturn(res);
781
782  if( h == NULL )
783  {
784    WarnS("Sort_c_ds needs an argument...");
785    return TRUE;
786  }
787
788  assume( h != NULL ); 
789
790#ifndef NDEBUG
791  const BOOLEAN __DEBUG__ = (BOOLEAN)((long)(atGet(currRingHdl,"DEBUG",INT_CMD, (void*)TRUE)));
792#else
793  const BOOLEAN __DEBUG__ = (BOOLEAN)((long)(atGet(currRingHdl,"DEBUG",INT_CMD, (void*)FALSE)));
794#endif
795
796  if(    (h->Typ() == IDEAL_CMD || h->Typ() == MODUL_CMD)
797      && (h->rtyp  == IDHDL) // must be a variable!
798      && (h->e == NULL) // not a list element
799      ) 
800  {
801    const ideal id = (const ideal)h->Data();
802
803    assume(id != NULL);
804
805    if( __DEBUG__ )
806    {
807      PrintS("Sort_c_ds::Input: \n");
808
809      const BOOLEAN __LEAD2SYZ__ = (BOOLEAN)((long)(atGet(currRingHdl,"LEAD2SYZ",INT_CMD, (void*)0)));
810      const BOOLEAN __TAILREDSYZ__ = (BOOLEAN)((long)(atGet(currRingHdl,"TAILREDSYZ",INT_CMD, (void*)0)));
811      const BOOLEAN __SYZCHECK__ = (BOOLEAN)((long)(atGet(currRingHdl,"SYZCHECK",INT_CMD, (void*)0)));   
812
813      Print("\nSYZCHECK: \t%d", __SYZCHECK__);
814      Print(", DEBUG: \t%d", __DEBUG__);
815      Print(", LEAD2SYZ: \t%d", __LEAD2SYZ__);
816      Print(", TAILREDSYZ: \t%d\n", __TAILREDSYZ__);
817
818      dPrint(id, r, r, 1);     
819    }
820
821    assume (h->Next() == NULL);
822
823    id_Test(id, r);
824
825    const int size = IDELEMS(id);
826
827    const ideal newid = id; // id_Copy(id, r); // copy???
828
829    if( size >= 2 )
830      qsort_r(newid->m, size, sizeof(poly), cmp_c_ds, r);
831   
832//    res->data = newid;
833//    res->rtyp = h->Typ();
834   
835    if( __DEBUG__ )
836    {
837      PrintS("Sort_c_ds::Output: \n");
838      dPrint(newid, r, r, 1);
839    }
840
841    return FALSE;
842  }
843
844  WarnS("ComputeLeadingSyzygyTerms needs a single ideal/module argument (must be a variable!)...");
845  return TRUE; 
846}
847
848
849
850static BOOLEAN Compute2LeadingSyzygyTerms(leftv res, leftv h)
851{
852  const ring r = currRing;
853  NoReturn(res);
854
855  if( h == NULL )
856  {
857    WarnS("Compute2LeadingSyzygyTerms needs an argument...");
858    return TRUE;
859  }
860
861  assume( h != NULL ); 
862
863#ifndef NDEBUG
864  const BOOLEAN __DEBUG__ = (BOOLEAN)((long)(atGet(currRingHdl,"DEBUG",INT_CMD, (void*)TRUE)));
865#else
866  const BOOLEAN __DEBUG__ = (BOOLEAN)((long)(atGet(currRingHdl,"DEBUG",INT_CMD, (void*)FALSE)));
867#endif
868
869  const BOOLEAN __TAILREDSYZ__ = (BOOLEAN)((long)(atGet(currRingHdl,"TAILREDSYZ",INT_CMD, (void*)0)));
870  const BOOLEAN __LEAD2SYZ__ = (BOOLEAN)((long)(atGet(currRingHdl,"LEAD2SYZ",INT_CMD, (void*)0)));
871  const BOOLEAN __SYZCHECK__ = (BOOLEAN)((long)(atGet(currRingHdl,"SYZCHECK",INT_CMD, (void*)0)));   
872
873
874  assume( __LEAD2SYZ__ );
875 
876  if( h->Typ() == IDEAL_CMD || h->Typ() == MODUL_CMD)
877  {
878    const ideal id = (const ideal)h->Data();
879
880    assume(id != NULL);
881
882    if( __DEBUG__ )
883    {
884      PrintS("Compute2LeadingSyzygyTerms::Input: \n");
885
886      Print("\nSYZCHECK: \t%d", __SYZCHECK__);
887      Print(", DEBUG: \t%d", __DEBUG__);
888      Print(", LEAD2SYZ: \t%d", __LEAD2SYZ__);
889      Print(", TAILREDSYZ: \t%d\n", __TAILREDSYZ__);
890
891      dPrint(id, r, r, 1);
892    }
893
894    h = h->Next(); assume (h == NULL);
895
896    // 1. set of components S?
897    // 2. for each component c from S: set of indices of leading terms
898    // with this component?
899    // 3. short exp. vectors for each leading term?
900
901    const int size = IDELEMS(id);
902
903    if( size < 2 )
904    {
905      const ideal newid = idInit(1, 1); // zero module...
906
907      newid->m[0] = NULL;
908
909      res->data = newid;
910      res->rtyp = MODUL_CMD;
911
912      return FALSE;
913    }
914
915
916    // TODO/NOTE: input is supposed to be sorted wrt "C,ds"!??
917
918    // components should come in groups: count elements in each group
919    // && estimate the real size!!!
920
921
922    // use just a vector instead???
923    ideal newid = idInit( (size * (size-1))/2, size); // maximal size: ideal case!
924
925    int k = 0;
926
927    for (int j = 0; j < size; j++)
928    {
929      const poly p = id->m[j];
930      assume( p != NULL );
931      const int  c = p_GetComp(p, r);
932
933      for (int i = j - 1; i >= 0; i--)
934      {
935        const poly pp = id->m[i];
936        assume( pp != NULL );
937        const int  cc = p_GetComp(pp, r);
938
939        if( c != cc )
940          continue;
941
942        // allocate memory & zero it out!
943        const poly m = p_Init(r); const poly mm = p_Init(r);
944       
945
946        // m = LCM(p, pp) / p! mm = LCM(p, pp) / pp!
947        // TODO: optimize: knowing the ring structure: (C/lp)!
948       
949        for (int v = rVar(r); v > 0; v--)
950        {
951          assume( v > 0 );
952          assume( v <= rVar(r) );
953
954          const short e1 = p_GetExp(p , v, r);
955          const short e2 = p_GetExp(pp, v, r);
956
957          if( e1 >= e2 )
958            p_SetExp(mm, v, e1 - e2, r); //            p_SetExp(m, v, 0, r);
959          else
960            p_SetExp(m, v, e2 - e1, r); //            p_SetExp(mm, v, 0, r);
961
962        }
963
964        assume( (j > i) && (i >= 0) );
965
966        p_SetComp(m, j + 1, r);
967        p_SetComp(mm, i + 1, r);
968       
969        const number& lc1 = p_GetCoeff(p , r);
970        const number& lc2 = p_GetCoeff(pp, r);
971
972        number g = n_Lcm( lc1, lc2, r );
973
974        p_SetCoeff0(m ,       n_Div(g, lc1, r), r);
975        p_SetCoeff0(mm, n_Neg(n_Div(g, lc2, r), r), r);
976
977        n_Delete(&g, r);
978
979        p_Setm(m, r); // should not do anything!!!
980        p_Setm(mm, r); // should not do anything!!!
981
982        pNext(m) = mm; //        pNext(mm) = NULL;
983       
984        newid->m[k++] = m;
985      }
986    }
987
988    if( __DEBUG__ && FALSE )
989    {
990      PrintS("Compute2LeadingSyzygyTerms::Temp0: \n");
991      dPrint(newid, r, r, 1);
992    }
993
994    if( !__TAILREDSYZ__ )
995    {
996      // simplify(newid, 2 + 32)??
997      // sort(newid, "C,ds")[1]???
998      id_DelDiv(newid, r); // #define SIMPL_LMDIV 32
999
1000      if( __DEBUG__ && FALSE )
1001      {
1002        PrintS("Compute2LeadingSyzygyTerms::Temp1 (deldiv): \n");
1003        dPrint(newid, r, r, 1);
1004      }
1005    }
1006    else
1007    {
1008      //      option(redSB); option(redTail);
1009      //      TEST_OPT_REDSB
1010      //      TEST_OPT_REDTAIL
1011      assume( r == currRing );
1012      BITSET _save_test = test;
1013      test |= (Sy_bit(OPT_REDTAIL) | Sy_bit(OPT_REDSB));
1014
1015      intvec* w=new intvec(IDELEMS(newid));
1016      ideal tmp = kStd(newid, currQuotient, isHomog, &w);
1017      delete w;
1018     
1019      test = _save_test;
1020     
1021      id_Delete(&newid, r);
1022      newid = tmp;
1023
1024      if( __DEBUG__ && FALSE )
1025      {
1026        PrintS("Compute2LeadingSyzygyTerms::Temp1 (std): \n");
1027        dPrint(newid, r, r, 1);
1028      }
1029     
1030    }
1031
1032    idSkipZeroes(newid);
1033
1034    const int sizeNew = IDELEMS(newid);
1035
1036    if( sizeNew >= 2 )
1037      qsort_r(newid->m, sizeNew, sizeof(poly), cmp_c_ds, r);
1038
1039    if (IDELEMS(newid) == 0 || (IDELEMS(newid) == 1 && newid->m[0] == NULL) )
1040      newid->rank = 1;
1041
1042    // TODO: add sorting wrt <c,ds> & reversing...
1043
1044    res->data = newid;
1045    res->rtyp = MODUL_CMD;
1046
1047    return FALSE;
1048  }
1049
1050  WarnS("Compute2LeadingSyzygyTerms needs a single ideal/module argument...");
1051  return TRUE;
1052}
1053
1054
1055/// Get leading term without a module component
1056static BOOLEAN leadmonom(leftv res, leftv h)
1057{
1058  NoReturn(res);
1059
1060  if ((h!=NULL) && (h->Typ()==VECTOR_CMD || h->Typ()==POLY_CMD) && (h->Data() != NULL))
1061  {
1062    const ring r = currRing;
1063    const poly p = (poly)(h->Data());
1064
1065    poly m = NULL;
1066
1067    if( p != NULL )
1068    {
1069      assume( p != NULL );
1070      assume( p_LmTest(p, r) );
1071
1072      m = p_LmInit(p, r);
1073      p_SetCoeff0(m, n_Copy(p_GetCoeff(p, r), r), r);
1074
1075      //            if( p_GetComp(m, r) != 0 )
1076      //            {
1077      p_SetComp(m, 0, r);
1078      p_Setm(m, r);
1079      //            }
1080
1081      assume( p_GetComp(m, r) == 0 );
1082      assume( m != NULL );
1083      assume( p_LmTest(m, r) );
1084    }
1085
1086    res->rtyp = POLY_CMD;
1087    res->data = reinterpret_cast<void *>(m);
1088
1089    return FALSE;
1090  }
1091
1092  WerrorS("`leadmonom(<poly/vector>)` expected");
1093  return TRUE;
1094}
1095
1096/// Get leading component
1097static BOOLEAN leadcomp(leftv res, leftv h)
1098{
1099  NoReturn(res);
1100
1101  if ((h!=NULL) && (h->Typ()==VECTOR_CMD || h->Typ()==POLY_CMD))
1102  {
1103    const ring r = currRing;
1104
1105    const poly p = (poly)(h->Data());
1106
1107    if (p != NULL )
1108    {
1109      assume( p != NULL );
1110      assume( p_LmTest(p, r) );
1111
1112      const unsigned long iComp = p_GetComp(p, r);
1113
1114  //    assume( iComp > 0 ); // p is a vector
1115
1116      res->data = reinterpret_cast<void *>(jjLONG2N(iComp));
1117    } else
1118      res->data = reinterpret_cast<void *>(jjLONG2N(0));
1119     
1120
1121    res->rtyp = BIGINT_CMD;
1122    return FALSE;
1123  }
1124
1125  WerrorS("`leadcomp(<poly/vector>)` expected");
1126  return TRUE;
1127}
1128
1129
1130
1131
1132/// Get raw leading exponent vector
1133static BOOLEAN leadrawexp(leftv res, leftv h)
1134{
1135  NoReturn(res);
1136
1137  if ((h!=NULL) && (h->Typ()==VECTOR_CMD || h->Typ()==POLY_CMD) && (h->Data() != NULL))
1138  {
1139    const ring r = currRing;
1140    const poly p = (poly)(h->Data());
1141
1142    assume( p != NULL );
1143    assume( p_LmTest(p, r) );
1144
1145    const int iExpSize = r->ExpL_Size;
1146
1147//    intvec *iv = new intvec(iExpSize);
1148
1149    lists l=(lists)omAllocBin(slists_bin);
1150    l->Init(iExpSize);
1151
1152    for(int i = iExpSize-1; i >= 0; i--)
1153    {
1154      l->m[i].rtyp = BIGINT_CMD;
1155      l->m[i].data = reinterpret_cast<void *>(jjLONG2N(p->exp[i])); // longs...
1156    }
1157
1158    res->rtyp = LIST_CMD; // list of bigints
1159    res->data = reinterpret_cast<void *>(l);
1160    return FALSE;
1161  }
1162
1163  WerrorS("`leadrawexp(<poly/vector>)` expected");
1164  return TRUE;
1165}
1166
1167
1168/// Endowe the current ring with additional (leading) Syz-component ordering
1169static BOOLEAN MakeSyzCompOrdering(leftv res, leftv /*h*/)
1170{
1171
1172  NoReturn(res);
1173
1174  //    res->data = rCurrRingAssure_SyzComp(); // changes current ring! :(
1175  res->data = reinterpret_cast<void *>(rAssure_SyzComp(currRing, TRUE));
1176  res->rtyp = RING_CMD; // return new ring!
1177  // QRING_CMD?
1178
1179  return FALSE;
1180}
1181
1182
1183/// Same for Induced Schreyer ordering (ordering on components is defined by sign!)
1184static BOOLEAN MakeInducedSchreyerOrdering(leftv res, leftv h)
1185{
1186
1187  NoReturn(res);
1188
1189  int sign = 1;
1190  if ((h!=NULL) && (h->Typ()==INT_CMD))
1191  {
1192    const int s = (int)((long)(h->Data()));
1193
1194    if( s != -1 && s != 1 )
1195    {
1196      WerrorS("`MakeInducedSchreyerOrdering(<int>)` called with wrong integer argument (must be +-1)!");
1197      return TRUE;
1198    }
1199
1200    sign = s;           
1201  }
1202
1203  assume( sign == 1 || sign == -1 );
1204  res->data = reinterpret_cast<void *>(rAssure_InducedSchreyerOrdering(currRing, TRUE, sign));
1205  res->rtyp = RING_CMD; // return new ring!
1206  // QRING_CMD?
1207  return FALSE;
1208}
1209
1210
1211/// Returns old SyzCompLimit, can set new limit
1212static BOOLEAN SetSyzComp(leftv res, leftv h)
1213{
1214  NoReturn(res);
1215
1216  const ring r = currRing;
1217
1218  if( !rIsSyzIndexRing(r) )
1219  {
1220    WerrorS("`SetSyzComp(<int>)` called on incompatible ring (not created by 'MakeSyzCompOrdering'!)");
1221    return TRUE;
1222  }
1223
1224  res->rtyp = INT_CMD;
1225  res->data = reinterpret_cast<void *>(rGetCurrSyzLimit(r)); // return old syz limit
1226
1227  if ((h!=NULL) && (h->Typ()==INT_CMD))
1228  {
1229    const int iSyzComp = (int)reinterpret_cast<long>(h->Data());
1230    assume( iSyzComp > 0 );
1231    rSetSyzComp(iSyzComp, currRing);
1232  }
1233
1234  return FALSE;
1235}
1236
1237/// ?
1238static BOOLEAN GetInducedData(leftv res, leftv h)
1239{
1240  NoReturn(res);
1241
1242  const ring r = currRing;
1243
1244  int p = 0; // which IS-block? p^th!
1245
1246  if ((h!=NULL) && (h->Typ()==INT_CMD))
1247  {
1248    p = (int)((long)(h->Data())); h=h->next;
1249    assume(p >= 0);
1250  }
1251
1252  const int pos = rGetISPos(p, r);
1253
1254  if(  /*(*/ -1 == pos /*)*/  )
1255  {
1256    WerrorS("`GetInducedData([int])` called on incompatible ring (not created by 'MakeInducedSchreyerOrdering'!)");
1257    return TRUE;
1258  }
1259
1260
1261  const int iLimit = r->typ[pos].data.is.limit;
1262  const ideal F = r->typ[pos].data.is.F;
1263  ideal FF = id_Copy(F, r);
1264
1265  lists l=(lists)omAllocBin(slists_bin);
1266  l->Init(2);
1267
1268  l->m[0].rtyp = INT_CMD;
1269  l->m[0].data = reinterpret_cast<void *>(iLimit);
1270
1271
1272  //        l->m[1].rtyp = MODUL_CMD;
1273
1274  if( idIsModule(FF, r) )
1275  {
1276    l->m[1].rtyp = MODUL_CMD;
1277
1278    //          Print("before: %d\n", FF->nrows);
1279    //          FF->nrows = id_RankFreeModule(FF, r); // ???
1280    //          Print("after: %d\n", FF->nrows);
1281  }
1282  else
1283    l->m[1].rtyp = IDEAL_CMD;
1284
1285  l->m[1].data = reinterpret_cast<void *>(FF);
1286
1287  res->rtyp = LIST_CMD; // list of int/module
1288  res->data = reinterpret_cast<void *>(l);
1289
1290  return FALSE;
1291
1292}
1293
1294
1295/* // the following turned out to be unnecessary...   
1296/// Finds p^th AM ordering, and returns its position in r->typ[] AND
1297/// corresponding &r->wvhdl[]
1298/// returns FALSE if something went wrong!
1299/// p - starts with 0!
1300BOOLEAN rGetAMPos(const ring r, const int p, int &typ_pos, int &wvhdl_pos, const BOOLEAN bSearchWvhdl = FALSE)
1301{
1302#if MYTEST
1303  Print("rGetAMPos(p: %d...)\nF:", p);
1304  PrintLn();
1305#endif
1306  typ_pos = -1;
1307  wvhdl_pos = -1;
1308
1309  if (r->typ==NULL)
1310    return FALSE;
1311
1312
1313  int j = p; // Which IS record to use...
1314  for( int pos = 0; pos < r->OrdSize; pos++ )
1315    if( r->typ[pos].ord_typ == ro_am)
1316      if( j-- == 0 )
1317      {
1318        typ_pos = pos;
1319
1320        if( bSearchWvhdl )
1321        {
1322          const int nblocks = rBlocks(r) - 1;
1323          const int* w = r->typ[pos].data.am.weights; // ?
1324
1325          for( pos = 0; pos <= nblocks; pos ++ )
1326            if (r->order[pos] == ringorder_am)
1327              if( r->wvhdl[pos] == w )
1328              {
1329                wvhdl_pos = pos;
1330                break;
1331              }
1332          if (wvhdl_pos < 0)
1333            return FALSE;
1334
1335          assume(wvhdl_pos >= 0);
1336        }
1337        assume(typ_pos >= 0);
1338        return TRUE;
1339      }
1340
1341  return FALSE;
1342}
1343
1344// // ?
1345// static BOOLEAN GetAMData(leftv res, leftv h)
1346// {
1347//   NoReturn(res);
1348//
1349//   const ring r = currRing;
1350//
1351//   int p = 0; // which IS-block? p^th!
1352//
1353//   if ((h!=NULL) && (h->Typ()==INT_CMD))
1354//     p = (int)((long)(h->Data())); h=h->next;
1355//
1356//   assume(p >= 0);
1357//
1358//   int d, w;
1359//   
1360//   if( !rGetAMPos(r, p, d, w, TRUE) )
1361//   {
1362//     Werror("`GetAMData([int])`: no %d^th _am block-ordering!", p);
1363//     return TRUE;
1364//   }
1365//
1366//   assume( r->typ[d].ord_typ == ro_am );
1367//   assume( r->order[w] == ringorder_am );
1368//
1369//
1370//   const short start = r->typ[d].data.am.start;  // bounds of ordering (in E)
1371//   const short end = r->typ[d].data.am.end;
1372//   const short len_gen = r->typ[d].data.am.len_gen; // i>len_gen: weight(gen(i)):=0
1373//   const int *weights = r->typ[d].data.am.weights; // pointers into wvhdl field of length (end-start+1) + len_gen
1374//   // contents w_1,... w_n, len, mod_w_1, .. mod_w_len, 0
1375//
1376//   assume( weights == r->wvhdl[w] );
1377//
1378//   
1379//   lists l=(lists)omAllocBin(slists_bin);
1380//   l->Init(2);
1381//
1382//   const short V = end-start+1;
1383//   intvec* ww_vars = new intvec(V);
1384//   intvec* ww_gens = new intvec(len_gen);
1385//
1386//   for (int i = 0; i < V; i++ )
1387//     (*ww_vars)[i] = weights[i];
1388//
1389//   assume( weights[V] == len_gen );
1390//
1391//   for (int i = 0; i < len_gen; i++ )
1392//     (*ww_gens)[i] = weights[i - V - 1];
1393//   
1394//
1395//   l->m[0].rtyp = INTVEC_CMD;
1396//   l->m[0].data = reinterpret_cast<void *>(ww_vars);
1397//
1398//   l->m[1].rtyp = INTVEC_CMD;
1399//   l->m[1].data = reinterpret_cast<void *>(ww_gens);
1400//
1401//
1402//   return FALSE;
1403//
1404// }
1405*/
1406
1407/// Returns old SyzCompLimit, can set new limit
1408static BOOLEAN SetInducedReferrence(leftv res, leftv h)
1409{
1410  NoReturn(res);
1411
1412  const ring r = currRing;
1413
1414  if( !( (h!=NULL) && ( (h->Typ()==IDEAL_CMD) || (h->Typ()==MODUL_CMD))) )
1415  {
1416    WerrorS("`SetInducedReferrence(<ideal/module>, [int[, int]])` expected");
1417    return TRUE;
1418  }
1419
1420  const ideal F = (ideal)h->Data(); ; // No copy!
1421  h=h->next;
1422
1423  int rank = 0;
1424
1425  if ((h!=NULL) && (h->Typ()==INT_CMD))
1426  {
1427    rank = (int)((long)(h->Data())); h=h->next;
1428    assume(rank >= 0);
1429  } else
1430    rank = id_RankFreeModule(F, r); // Starting syz-comp (1st: i+1)
1431
1432  int p = 0; // which IS-block? p^th!
1433
1434  if ((h!=NULL) && (h->Typ()==INT_CMD))
1435  {
1436    p = (int)((long)(h->Data())); h=h->next;
1437    assume(p >= 0);
1438  }
1439
1440  const int posIS = rGetISPos(p, r);
1441
1442  if(  /*(*/ -1 == posIS /*)*/  )
1443  {
1444    WerrorS("`SetInducedReferrence(<ideal/module>, [int[, int]])` called on incompatible ring (not created by 'MakeInducedSchreyerOrdering'!)");
1445    return TRUE;
1446  }
1447
1448
1449
1450  // F & componentWeights belong to that ordering block of currRing now:
1451  rSetISReference(r, F, rank, p); // F will be copied!
1452  return FALSE;
1453}
1454
1455
1456//    F = ISUpdateComponents( F, V, MIN );
1457//    // replace gen(i) -> gen(MIN + V[i-MIN]) for all i > MIN in all terms from F!
1458static BOOLEAN ISUpdateComponents(leftv res, leftv h)
1459{
1460  NoReturn(res);
1461
1462  PrintS("ISUpdateComponents:.... \n");
1463
1464  if ((h!=NULL) && (h->Typ()==MODUL_CMD))
1465  {
1466    ideal F = (ideal)h->Data(); ; // No copy!
1467    h=h->next;
1468
1469    if ((h!=NULL) && (h->Typ()==INTVEC_CMD))
1470    {
1471      const intvec* const V = (const intvec* const) h->Data();
1472      h=h->next;
1473
1474      if ((h!=NULL) && (h->Typ()==INT_CMD))
1475      {
1476        const int MIN = (int)((long)(h->Data()));
1477
1478        pISUpdateComponents(F, V, MIN, currRing);
1479        return FALSE;
1480      }
1481    }
1482  }
1483
1484  WerrorS("`ISUpdateComponents(<module>, intvec, int)` expected");
1485  return TRUE;
1486}
1487
1488
1489/// NF using length
1490static BOOLEAN reduce_syz(leftv res, leftv h)
1491{
1492  // const ring r = currRing;
1493
1494  if ( !( (h!=NULL) && (h->Typ()==VECTOR_CMD || h->Typ()==POLY_CMD) ) )
1495  {
1496    WerrorS("`reduce_syz(<poly/vector>!, <ideal/module>, <int>, [int])` expected");
1497    return TRUE;
1498  }
1499
1500  res->rtyp = h->Typ();
1501  const poly v = reinterpret_cast<poly>(h->Data());
1502  h=h->next;
1503
1504  if ( !( (h!=NULL) && (h->Typ()==MODUL_CMD || h->Typ()==IDEAL_CMD ) ) )
1505  {
1506    WerrorS("`reduce_syz(<poly/vector>, <ideal/module>!, <int>, [int])` expected");
1507    return TRUE;
1508  }
1509
1510  assumeStdFlag(h);
1511  const ideal M = reinterpret_cast<ideal>(h->Data()); h=h->next;
1512
1513
1514  if ( !( (h!=NULL) && (h->Typ()== INT_CMD)  ) )
1515  {
1516    WerrorS("`reduce_syz(<poly/vector>, <ideal/module>, <int>!, [int])` expected");
1517    return TRUE;
1518  }
1519
1520  const int iSyzComp = (int)((long)(h->Data())); h=h->next;
1521
1522  int iLazyReduce = 0;
1523
1524  if ( ( (h!=NULL) && (h->Typ()== INT_CMD)  ) )
1525    iLazyReduce = (int)((long)(h->Data())); 
1526
1527  res->data = (void *)kNFLength(M, currQuotient, v, iSyzComp, iLazyReduce);
1528  return FALSE;
1529}
1530
1531
1532/// Get raw syzygies (idPrepare)
1533static BOOLEAN idPrepare(leftv res, leftv h)
1534{
1535  //        extern int rGetISPos(const int p, const ring r);
1536
1537  const ring r = currRing;
1538
1539  const bool isSyz = rIsSyzIndexRing(r);
1540  const int posIS = rGetISPos(0, r);
1541
1542
1543  if ( !( (h!=NULL) && (h->Typ()==MODUL_CMD) && (h->Data() != NULL) ) )
1544  {
1545    WerrorS("`idPrepare(<module>)` expected");
1546    return TRUE;
1547  }
1548
1549  const ideal I = reinterpret_cast<ideal>(h->Data());
1550
1551  assume( I != NULL );
1552  idTest(I);
1553
1554  int iComp = -1;
1555
1556  h=h->next;
1557  if ( (h!=NULL) && (h->Typ()==INT_CMD) )
1558  {
1559    iComp = (int)((long)(h->Data()));
1560  } else
1561  {
1562      if( (!isSyz) && (-1 == posIS) )
1563      {
1564        WerrorS("`idPrepare(<...>)` called on incompatible ring (not created by 'MakeSyzCompOrdering' or 'MakeInducedSchreyerOrdering'!)");
1565        return TRUE;
1566      }
1567
1568    if( isSyz )
1569      iComp = rGetCurrSyzLimit(r);
1570    else
1571      iComp = id_RankFreeModule(r->typ[posIS].data.is.F, r); // ;
1572  }
1573 
1574  assume(iComp >= 0);
1575
1576
1577  intvec* w = reinterpret_cast<intvec *>(atGet(h, "isHomog", INTVEC_CMD));
1578  tHomog hom = testHomog;
1579
1580  //           int add_row_shift = 0;
1581  //
1582  if (w!=NULL)
1583  {
1584    w = ivCopy(w);
1585  //             add_row_shift = ww->min_in();
1586  //
1587  //             (*ww) -= add_row_shift;
1588  //             
1589  //             if (idTestHomModule(I, currQuotient, ww))
1590  //             {
1591    hom = isHomog;
1592  //               w = ww;
1593  //             }
1594  //             else
1595  //             {
1596  //               //WarnS("wrong weights");
1597  //               delete ww;
1598  //               w = NULL;
1599  //               hom=testHomog;
1600  //             }
1601  }
1602
1603
1604  // computes syzygies of h1,
1605  // works always in a ring with ringorder_s
1606  // NOTE: rSetSyzComp(syzcomp) should better be called beforehand
1607  //        ideal idPrepare (ideal  h1, tHomog hom, int syzcomp, intvec **w);
1608
1609  ideal J = // idPrepare( I, hom, iComp, &w);
1610           kStd(I, currQuotient, hom, &w, NULL, iComp);
1611
1612  idTest(J);
1613
1614  if (w!=NULL)
1615    atSet(res, omStrDup("isHomog"), w, INTVEC_CMD);
1616  //             if (w!=NULL) delete w;
1617
1618  res->rtyp = MODUL_CMD;
1619  res->data = reinterpret_cast<void *>(J);
1620  return FALSE;
1621}
1622
1623/// Get raw syzygies (idPrepare)
1624static BOOLEAN _p_Content(leftv res, leftv h)
1625{
1626  if ( !( (h!=NULL) && (h->Typ()==POLY_CMD) && (h->Data() != NULL) ) )
1627  {
1628    WerrorS("`p_Content(<poly-var>)` expected");
1629    return TRUE;
1630  }
1631
1632
1633  const poly p = reinterpret_cast<poly>(h->Data());
1634
1635 
1636  pTest(p);  pWrite(p); PrintLn();
1637
1638 
1639  p_Content( p, currRing);     
1640
1641  pTest(p);
1642  pWrite(p); PrintLn();
1643 
1644  NoReturn(res);
1645  return false;
1646}
1647
1648static BOOLEAN _m2_end(leftv res, leftv h)
1649{
1650  int ret = 0;
1651 
1652  if ( (h!=NULL) && (h->Typ()!=INT_CMD) )
1653  {
1654    WerrorS("`m2_end([<int>])` expected");
1655    return TRUE;
1656  }
1657  ret = (int)(long)(h->Data());
1658
1659  m2_end( ret );
1660
1661  NoReturn(res);
1662  return FALSE;
1663}
1664
1665   
1666
1667END_NAMESPACE
1668
1669
1670int SI_MOD_INIT(syzextra)(SModulFunctions* psModulFunctions) 
1671{
1672#define ADD0(A,B,C,D,E) A(B, (char*)C, D, E)
1673// #define ADD(A,B,C,D,E) ADD0(iiAddCproc, "", C, D, E)
1674  #define ADD(A,B,C,D,E) ADD0(A->iiAddCproc, B, C, D, E)
1675  ADD(psModulFunctions, currPack->libname, "ClearContent", FALSE, _ClearContent);
1676  ADD(psModulFunctions, currPack->libname, "ClearDenominators", FALSE, _ClearDenominators);
1677
1678  ADD(psModulFunctions, currPack->libname, "m2_end", FALSE, _m2_end);
1679
1680  ADD(psModulFunctions, currPack->libname, "DetailedPrint", FALSE, DetailedPrint);
1681  ADD(psModulFunctions, currPack->libname, "leadmonomial", FALSE, leadmonom);
1682  ADD(psModulFunctions, currPack->libname, "leadcomp", FALSE, leadcomp);
1683  ADD(psModulFunctions, currPack->libname, "leadrawexp", FALSE, leadrawexp);
1684
1685  ADD(psModulFunctions, currPack->libname, "ISUpdateComponents", FALSE, ISUpdateComponents);
1686  ADD(psModulFunctions, currPack->libname, "SetInducedReferrence", FALSE, SetInducedReferrence);
1687  ADD(psModulFunctions, currPack->libname, "GetInducedData", FALSE, GetInducedData);
1688  ADD(psModulFunctions, currPack->libname, "SetSyzComp", FALSE, SetSyzComp);
1689  ADD(psModulFunctions, currPack->libname, "MakeInducedSchreyerOrdering", FALSE, MakeInducedSchreyerOrdering);
1690  ADD(psModulFunctions, currPack->libname, "MakeSyzCompOrdering", FALSE, MakeSyzCompOrdering);
1691
1692  ADD(psModulFunctions, currPack->libname, "ProfilerStart", FALSE, _ProfilerStart); ADD(psModulFunctions, currPack->libname, "ProfilerStop",  FALSE, _ProfilerStop );
1693 
1694  ADD(psModulFunctions, currPack->libname, "noop", FALSE, noop);
1695  ADD(psModulFunctions, currPack->libname, "idPrepare", FALSE, idPrepare);
1696  ADD(psModulFunctions, currPack->libname, "reduce_syz", FALSE, reduce_syz);
1697
1698  ADD(psModulFunctions, currPack->libname, "p_Content", FALSE, _p_Content);
1699
1700  ADD(psModulFunctions, currPack->libname, "Tail", FALSE, Tail);
1701 
1702  ADD(psModulFunctions, currPack->libname, "ComputeLeadingSyzygyTerms", FALSE, ComputeLeadingSyzygyTerms);
1703  ADD(psModulFunctions, currPack->libname, "Compute2LeadingSyzygyTerms", FALSE, Compute2LeadingSyzygyTerms);
1704 
1705  ADD(psModulFunctions, currPack->libname, "Sort_c_ds", FALSE, Sort_c_ds);
1706 
1707  //  ADD(psModulFunctions, currPack->libname, "GetAMData", FALSE, GetAMData);
1708
1709  //  ADD(psModulFunctions, currPack->libname, "", FALSE, );
1710
1711#undef ADD 
1712  return 0;
1713}
1714
1715#ifndef EMBED_PYTHON
1716extern "C" { 
1717int mod_init(SModulFunctions* psModulFunctions)
1718{ 
1719  return SI_MOD_INIT(syzextra)(psModulFunctions); 
1720}
1721}
1722#endif
Note: See TracBrowser for help on using the repository browser.