source: git/dyn_modules/syzextra/mod_main.cc @ 7b7c2c

spielwiese
Last change on this file since 7b7c2c was 7b7c2c, checked in by Oleksandr Motsak <motsak@…>, 12 years ago
further Syz/Resolutions: chg: sorting: "ds" ->> "c,ds" (since reversed order!) chg: SSComputeLeadingSyzygyTerms moved into the dynamic module (ComputeLeadingSyzygyTerms)! TODO: add id_Sort wrt "c,ds" to ComputeLeadingSyzygyTerms???
  • Property mode set to 100644
File size: 28.2 KB
Line 
1
2
3
4
5#include <kernel/mod2.h>
6
7#include <omalloc/omalloc.h>
8
9#include <misc/intvec.h>
10
11#include <coeffs/coeffs.h>
12
13#include <polys/PolyEnumerator.h>
14
15#include <polys/monomials/p_polys.h>
16#include <polys/monomials/ring.h>
17// #include <kernel/longrat.h>
18#include <kernel/GBEngine/kstd1.h>
19
20#include <kernel/polys.h>
21
22#include <kernel/GBEngine/syz.h>
23
24#include <Singular/tok.h>
25#include <Singular/ipid.h>
26#include <Singular/lists.h>
27#include <Singular/attrib.h>
28
29#include <Singular/ipid.h> 
30#include <Singular/ipshell.h> // For iiAddCproc
31
32// extern coeffs coeffs_BIGINT
33
34#include "singularxx_defs.h"
35#include "DebugPrint.h"
36#include "myNF.h"
37
38
39#include <Singular/mod_lib.h>
40
41
42#if GOOGLE_PROFILE_ENABLED
43#include <google/profiler.h>
44#endif // #if GOOGLE_PROFILE_ENABLED
45
46
47#include <stdio.h>
48#include <stdlib.h>
49#include <string.h>
50
51
52extern void pISUpdateComponents(ideal F, const intvec *const V, const int MIN, const ring r);
53// extern ring rCurrRingAssure_SyzComp();
54extern ring rAssure_InducedSchreyerOrdering(const ring r, BOOLEAN complete, int sign);
55extern int rGetISPos(const int p, const ring r);
56
57
58USING_NAMESPACE( SINGULARXXNAME :: DEBUG )
59USING_NAMESPACE( SINGULARXXNAME :: NF )
60
61BEGIN_NAMESPACE_NONAME
62static inline void NoReturn(leftv& res)
63{
64  res->rtyp = NONE;
65  res->data = NULL;
66}
67
68/// wrapper around n_ClearContent
69static BOOLEAN _ClearContent(leftv res, leftv h)
70{
71  NoReturn(res);
72
73  const char *usage = "'ClearContent' needs a (non-zero!) poly or vector argument...";
74 
75  if( h == NULL )
76  {
77    WarnS(usage);
78    return TRUE;
79  }
80
81  assume( h != NULL );
82
83  if( !( h->Typ() == POLY_CMD || h->Typ() == VECTOR_CMD) )
84  {
85    WarnS(usage);
86    return TRUE;
87  }
88
89  assume (h->Next() == NULL);
90 
91  poly ph = reinterpret_cast<poly>(h->Data());
92 
93  if( ph == NULL )
94  {
95    WarnS(usage);
96    return TRUE;
97  }
98 
99  const ring r =  currRing;
100  assume( r != NULL ); assume( r->cf != NULL ); const coeffs C = r->cf;
101
102  number n;
103
104  // experimentall (recursive enumerator treatment) of alg. ext
105  CPolyCoeffsEnumerator itr(ph);
106  n_ClearContent(itr, n, C);
107
108  res->data = n;
109  res->rtyp = NUMBER_CMD;
110
111  return FALSE;
112}
113
114/// wrapper around n_ClearDenominators
115static BOOLEAN _ClearDenominators(leftv res, leftv h)
116{
117  NoReturn(res);
118
119  const char *usage = "'ClearDenominators' needs a (non-zero!) poly or vector argument...";
120
121  if( h == NULL )
122  {
123    WarnS(usage);
124    return TRUE;
125  }
126
127  assume( h != NULL );
128
129  if( !( h->Typ() == POLY_CMD || h->Typ() == VECTOR_CMD) )
130  {
131    WarnS(usage);
132    return TRUE;
133  }
134
135  assume (h->Next() == NULL);
136
137  poly ph = reinterpret_cast<poly>(h->Data());
138
139  if( ph == NULL )
140  {
141    WarnS(usage);
142    return TRUE;
143  }
144
145  const ring r =  currRing;
146  assume( r != NULL ); assume( r->cf != NULL ); const coeffs C = r->cf;
147
148  number n;
149
150  // experimentall (recursive enumerator treatment) of alg. ext.
151  CPolyCoeffsEnumerator itr(ph);
152  n_ClearDenominators(itr, n, C);
153
154  res->data = n;
155  res->rtyp = NUMBER_CMD;
156
157  return FALSE;
158}
159
160
161/// try to get an optional (simple) integer argument out of h
162/// or return the default value
163static int getOptionalInteger(const leftv& h, const int _n)
164{
165  if( h!= NULL && h->Typ() == INT_CMD )
166  {
167    int n = (int)(long)(h->Data());
168
169    if( n < 0 )
170      Warn("Negative (%d) optional integer argument", n);
171
172    return (n);
173  }
174
175  return (_n); 
176}
177
178static BOOLEAN noop(leftv __res, leftv /*__v*/)
179{
180  NoReturn(__res);
181  return FALSE;
182}
183
184static BOOLEAN _ProfilerStart(leftv __res, leftv h)
185{
186  NoReturn(__res);
187#if GOOGLE_PROFILE_ENABLED
188  if( h!= NULL && h->Typ() == STRING_CMD )
189  {
190    const char* name = (char*)(h->Data());
191    assume( name != NULL );   
192    ProfilerStart(name);
193  } else
194    WerrorS("ProfilerStart requires a string [name] argument"); 
195#else
196  WarnS("Sorry no google profiler support (GOOGLE_PROFILE_ENABLE!=1)...");
197//  return TRUE; // ?
198#endif // #if GOOGLE_PROFILE_ENABLED
199  return FALSE;
200  (void)h;
201}
202static BOOLEAN _ProfilerStop(leftv __res, leftv /*__v*/)
203{
204  NoReturn(__res);
205#if GOOGLE_PROFILE_ENABLED
206  ProfilerStop();
207#else
208  WarnS("Sorry no google profiler support (GOOGLE_PROFILE_ENABLED!=1)...");
209//  return TRUE; // ?
210#endif // #if GOOGLE_PROFILE_ENABLED
211  return FALSE;
212}
213
214static inline number jjLONG2N(long d)
215{
216  return n_Init(d, coeffs_BIGINT);
217}
218
219static inline void view(const intvec* v)
220{
221#ifndef SING_NDEBUG
222  v->view();
223#else
224  // This code duplication is only due to Hannes's #ifndef SING_NDEBUG!
225  Print ("intvec: {rows: %d, cols: %d, length: %d, Values: \n", v->rows(), v->cols(), v->length());
226
227  for (int i = 0; i < v->rows(); i++)
228  {
229    Print ("Row[%3d]:", i);
230    for (int j = 0; j < v->cols(); j++)
231      Print (" %5d", (*v)[j + i * (v->cols())] );
232    PrintLn ();
233  }
234  PrintS ("}\n");
235#endif
236
237}
238
239                   
240
241static BOOLEAN DetailedPrint(leftv __res, leftv h)
242{
243  NoReturn(__res);
244
245  if( h == NULL )
246  {
247    WarnS("DetailedPrint needs an argument...");
248    return TRUE;
249  }
250
251  if( h->Typ() == NUMBER_CMD)
252  {
253    number n = (number)h->Data(); 
254
255    const ring r = currRing;
256
257#ifdef LDEBUG
258    r->cf->cfDBTest(n,__FILE__,__LINE__,r->cf);
259#endif
260
261    StringSetS("");
262    n_Write(n, r->cf);
263    PrintS(StringEndS());
264    PrintLn();
265
266    return FALSE;
267  }
268 
269  if( h->Typ() == RING_CMD)
270  {
271    const ring r = (const ring)h->Data();
272    rWrite(r, TRUE);
273    PrintLn();
274#ifdef RDEBUG
275    rDebugPrint(r);
276#endif
277    return FALSE;
278  }
279
280  if( h->Typ() == POLY_CMD || h->Typ() == VECTOR_CMD)
281  {
282    const poly p = (const poly)h->Data(); h = h->Next();
283
284    dPrint(p, currRing, currRing, getOptionalInteger(h, 3));
285
286    return FALSE;
287  }
288
289  if( h->Typ() == IDEAL_CMD || h->Typ() == MODUL_CMD)
290  {
291    const ideal id = (const ideal)h->Data(); h = h->Next(); 
292
293    dPrint(id, currRing, currRing, getOptionalInteger(h, 3));
294   
295    return FALSE;           
296  }
297
298  if( h->Typ() == RESOLUTION_CMD )
299  {
300    const syStrategy syzstr = reinterpret_cast<const syStrategy>(h->Data());
301
302    h = h->Next();
303
304    int nTerms = getOptionalInteger(h, 1);
305
306
307    Print("RESOLUTION_CMD(%p): ", reinterpret_cast<const void*>(syzstr)); PrintLn();
308
309    const ring save = currRing;
310    const ring r = syzstr->syRing;
311    const ring rr = (r != NULL) ? r: save;
312
313
314    const int iLength = syzstr->length;
315
316    Print("int 'length': %d", iLength); PrintLn();
317    Print("int 'regularity': %d", syzstr->regularity); PrintLn();
318    Print("short 'list_length': %hd", syzstr->list_length); PrintLn();
319    Print("short 'references': %hd", syzstr->references); PrintLn();
320
321
322#define PRINT_pINTVECTOR(s, v) Print("intvec '%10s'(%p)", #v, reinterpret_cast<const void*>((s)->v)); \
323if( (s)->v != NULL ){ PrintS(": "); view((s)->v); }; \
324PrintLn();
325
326    PRINT_pINTVECTOR(syzstr, resolution);
327    PRINT_pINTVECTOR(syzstr, betti);
328    PRINT_pINTVECTOR(syzstr, Tl);
329    PRINT_pINTVECTOR(syzstr, cw);
330#undef PRINT_pINTVECTOR
331
332    if (r == NULL)
333      Print("ring '%10s': NULL", "syRing");
334    else 
335      if (r == currRing)
336        Print("ring '%10s': currRing", "syRing");
337      else
338        if (r != NULL && r != save)
339        {
340          Print("ring '%10s': ", "syRing");
341          rWrite(r);
342#ifdef RDEBUG
343          //              rDebugPrint(r);
344#endif
345          // rChangeCurrRing(r);
346        }           
347    PrintLn();
348
349    const SRes rP = syzstr->resPairs;
350    Print("SRes 'resPairs': %p", reinterpret_cast<const void*>(rP)); PrintLn();
351
352    if (rP != NULL)
353      for (int iLevel = 0; (iLevel < iLength) && (rP[iLevel] != NULL) && ((*syzstr->Tl)[iLevel] >= 0); iLevel++)
354      {
355        int n = 0;
356        const int iTl = (*syzstr->Tl)[iLevel];
357        for (int j = 0; (j < iTl) && ((rP[iLevel][j].lcm!=NULL) || (rP[iLevel][j].syz!=NULL)); j++)
358        {
359          if (rP[iLevel][j].isNotMinimal==NULL)
360            n++;
361        }
362        Print("minimal-resPairs-Size[1+%d]: %d", iLevel, n); PrintLn();
363      }
364
365
366    //  const ring rrr = (iLevel > 0) ? rr : save; ?
367#define PRINT_RESOLUTION(s, v) Print("resolution '%12s': %p", #v, reinterpret_cast<const void*>((s)->v)); PrintLn(); \
368if ((s)->v != NULL) \
369  for (int iLevel = 0; (iLevel < iLength) && ( ((s)->v)[iLevel] != NULL ); iLevel++) \
370  { \
371    /* const ring rrr = (iLevel > 0) ? save : save; */ \
372    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)*/ ); \
373    PrintLn(); \
374  } \
375  PrintLn();
376
377    // resolvente:
378    PRINT_RESOLUTION(syzstr, minres);
379    PRINT_RESOLUTION(syzstr, fullres);
380
381    assume (id_RankFreeModule (syzstr->res[1], rr) == syzstr->res[1]->rank);
382
383    PRINT_RESOLUTION(syzstr, res);
384    PRINT_RESOLUTION(syzstr, orderedRes);
385#undef PRINT_RESOLUTION
386
387#define PRINT_POINTER(s, v) Print("pointer '%17s': %p", #v, reinterpret_cast<const void*>((s)->v)); PrintLn();
388    // 2d arrays:
389    PRINT_POINTER(syzstr, truecomponents);
390    PRINT_POINTER(syzstr, ShiftedComponents);
391    PRINT_POINTER(syzstr, backcomponents);
392    PRINT_POINTER(syzstr, Howmuch);
393    PRINT_POINTER(syzstr, Firstelem);
394    PRINT_POINTER(syzstr, elemLength);
395    PRINT_POINTER(syzstr, sev);
396
397    // arrays of intvects:
398    PRINT_POINTER(syzstr, weights);
399    PRINT_POINTER(syzstr, hilb_coeffs);
400#undef PRINT_POINTER
401
402
403    if (syzstr->fullres==NULL)
404    {
405      PrintS("resolution 'fullres': (NULL) => resolution not computed yet");
406      PrintLn();
407    } else
408    {
409      Print("resolution 'fullres': (%p) => resolution seems to be computed already", reinterpret_cast<const void*>(syzstr->fullres));
410      PrintLn();
411      dPrint(*syzstr->fullres, save, save, nTerms);
412    }
413
414
415
416
417    if (syzstr->minres==NULL)
418    {
419      PrintS("resolution 'minres': (NULL) => resolution not minimized yet");
420      PrintLn();
421    } else
422    {
423      Print("resolution 'minres': (%p) => resolution seems to be minimized already", reinterpret_cast<const void*>(syzstr->minres));
424      PrintLn();
425      dPrint(*syzstr->minres, save, save, nTerms);
426    }
427
428
429
430
431    /*
432    int ** truecomponents;
433    long** ShiftedComponents;
434    int ** backcomponents;
435    int ** Howmuch;
436    int ** Firstelem;
437    int ** elemLength;
438    unsigned long ** sev;
439
440    intvec ** weights;
441    intvec ** hilb_coeffs;
442
443    SRes resPairs;               //polynomial data for internal use only
444
445    resolvente fullres;
446    resolvente minres;
447    resolvente res;              //polynomial data for internal use only
448    resolvente orderedRes;       //polynomial data for internal use only
449*/
450
451    //            if( currRing != save ) rChangeCurrRing(save);
452  }
453
454
455  return FALSE;
456}
457
458
459
460
461
462
463static BOOLEAN Tail(leftv res, leftv h)
464{
465  NoReturn(res);
466
467  if( h == NULL )
468  {
469    WarnS("Tail needs an argument...");
470    return TRUE;
471  }
472
473  assume( h != NULL );
474
475  if( h->Typ() == POLY_CMD || h->Typ() == VECTOR_CMD)
476  {
477    const poly p = (const poly)h->Data();
478
479    res->rtyp = h->Typ();
480
481    if( p == NULL)
482      res->data = NULL;
483    else
484      res->data = p_Copy( pNext(p), currRing );
485
486    h = h->Next(); assume (h == NULL);
487   
488    return FALSE;
489  }
490
491  if( h->Typ() == IDEAL_CMD || h->Typ() == MODUL_CMD)
492  {
493    const ideal id = (const ideal)h->Data();
494
495    res->rtyp = h->Typ();
496
497    if( id == NULL)
498      res->data = NULL;
499    else
500    {
501      const ideal newid = idInit(IDELEMS(id),id->rank);
502      for (int i=IDELEMS(id) - 1; i >= 0; i--)
503        if( id->m[i] != NULL )
504          newid->m[i] = p_Copy(pNext(id->m[i]), currRing);
505        else
506          newid->m[i] = NULL;
507     
508      res->data = newid; 
509    }
510   
511    h = h->Next(); assume (h == NULL);
512   
513    return FALSE;
514  }
515
516  WarnS("Tail needs a single poly/vector/ideal/module argument...");
517  return TRUE;
518}
519
520
521static BOOLEAN ComputeLeadingSyzygyTerms(leftv res, leftv h)
522{
523  const ring r = currRing;
524  NoReturn(res);
525
526  if( h == NULL )
527  {
528    WarnS("ComputeLeadingSyzygyTerms needs an argument...");
529    return TRUE;
530  }
531
532  assume( h != NULL ); 
533
534#ifndef NDEBUG
535  const BOOLEAN __DEBUG__ = (BOOLEAN)((long)(atGet(currRingHdl,"DEBUG",INT_CMD, (void*)TRUE)));
536#else
537  const BOOLEAN __DEBUG__ = (BOOLEAN)((long)(atGet(currRingHdl,"DEBUG",INT_CMD, (void*)FALSE)));
538#endif
539
540  if( h->Typ() == IDEAL_CMD || h->Typ() == MODUL_CMD)
541  {
542    const ideal id = (const ideal)h->Data();
543
544    assume(id != NULL);
545
546    if( __DEBUG__ )
547    {
548      PrintS("ComputeLeadingSyzygyTerms::Input: \n");
549     
550      const BOOLEAN __LEAD2SYZ__ = (BOOLEAN)((long)(atGet(currRingHdl,"LEAD2SYZ",INT_CMD, (void*)0)));
551      const BOOLEAN __TAILREDSYZ__ = (BOOLEAN)((long)(atGet(currRingHdl,"TAILREDSYZ",INT_CMD, (void*)0)));
552      const BOOLEAN __SYZCHECK__ = (BOOLEAN)((long)(atGet(currRingHdl,"SYZCHECK",INT_CMD, (void*)0)));   
553
554      Print("\nSYZCHECK: \t%d", __SYZCHECK__);
555      Print(", DEBUG: \t%d", __DEBUG__);
556      Print(", LEAD2SYZ: \t%d", __LEAD2SYZ__);
557      Print(", TAILREDSYZ: \t%d\n", __TAILREDSYZ__);
558
559      dPrint(id, r, r, 1);
560    }
561
562    h = h->Next(); assume (h == NULL);
563
564    // 1. set of components S?
565    // 2. for each component c from S: set of indices of leading terms
566    // with this component?
567    // 3. short exp. vectors for each leading term?
568
569    const int size = IDELEMS(id);
570
571    if( size < 2 )
572    {
573      const ideal newid = idInit(1, 0); // zero ideal...
574
575      newid->m[0] = NULL;
576     
577      res->data = newid;
578      res->rtyp = MODUL_CMD;
579
580      return FALSE;
581    }
582   
583
584    // TODO/NOTE: input is supposed to be sorted wrt "C,ds"!??
585
586    // components should come in groups: count elements in each group
587    // && estimate the real size!!!
588
589
590    // use just a vector instead???
591    const ideal newid = idInit( (size * (size-1))/2, size); // maximal size: ideal case!
592   
593    int k = 0;
594
595    for (int j = 0; j < size; j++)
596    {
597      const poly p = id->m[j];
598      assume( p != NULL );
599      const int  c = p_GetComp(p, r);
600
601      for (int i = j - 1; i >= 0; i--)
602      {
603        const poly pp = id->m[i];
604        assume( pp != NULL );
605        const int  cc = p_GetComp(pp, r);
606
607        if( c != cc )
608          continue;
609
610        const poly m = p_Init(r); // p_New???
611
612        // m = LCM(p, pp) / p! // TODO: optimize: knowing the ring structure: (C/lp)!
613        for (int v = rVar(r); v > 0; v--)
614        {
615          assume( v > 0 );
616          assume( v <= rVar(r) );
617         
618          const short e1 = p_GetExp(p , v, r);
619          const short e2 = p_GetExp(pp, v, r);
620
621          if( e1 >= e2 )
622            p_SetExp(m, v, 0, r);
623          else
624            p_SetExp(m, v, e2 - e1, r);
625           
626        }
627
628        assume( (j > i) && (i >= 0) );
629
630        p_SetComp(m, j + 1, r);
631        pNext(m) = NULL;
632        p_SetCoeff0(m, n_Init(1, r->cf), r); // for later...
633
634        p_Setm(m, r); // should not do anything!!!
635
636        newid->m[k++] = m;
637      }
638    }
639
640    if( __DEBUG__ && FALSE )
641    {
642      PrintS("ComputeLeadingSyzygyTerms::Temp0: \n");
643      dPrint(newid, r, r, 1);
644    }
645
646    // the rest of newid is assumed to be zeroes...
647
648    // simplify(newid, 2 + 32)??
649    // sort(newid, "C,ds")[1]???
650    id_DelDiv(newid, r); // #define SIMPL_LMDIV 32
651
652    if( __DEBUG__ && FALSE )
653    {
654      PrintS("ComputeLeadingSyzygyTerms::Temp1: \n");
655      dPrint(newid, r, r, 1);
656    }
657
658    idSkipZeroes(newid); // #define SIMPL_NULL 2
659
660    if( __DEBUG__ )
661    {
662      PrintS("ComputeLeadingSyzygyTerms::Output: \n");
663      dPrint(newid, r, r, 1);
664    }
665
666
667    // TODO: add sorting wrt <c,ds> & reversing...
668
669    res->data = newid;
670    res->rtyp = MODUL_CMD;
671   
672    return FALSE;
673  }
674
675  WarnS("ComputeLeadingSyzygyTerms needs a single ideal/module argument...");
676  return TRUE;
677}
678
679
680
681
682
683/// Get leading term without a module component
684static BOOLEAN leadmonom(leftv res, leftv h)
685{
686  NoReturn(res);
687
688  if ((h!=NULL) && (h->Typ()==VECTOR_CMD || h->Typ()==POLY_CMD) && (h->Data() != NULL))
689  {
690    const ring r = currRing;
691    const poly p = (poly)(h->Data());
692
693    poly m = NULL;
694
695    if( p != NULL )
696    {
697      assume( p != NULL );
698      assume( p_LmTest(p, r) );
699
700      m = p_LmInit(p, r);
701      p_SetCoeff0(m, n_Copy(p_GetCoeff(p, r), r), r);
702
703      //            if( p_GetComp(m, r) != 0 )
704      //            {
705      p_SetComp(m, 0, r);
706      p_Setm(m, r);
707      //            }
708
709      assume( p_GetComp(m, r) == 0 );
710      assume( m != NULL );
711      assume( p_LmTest(m, r) );
712    }
713
714    res->rtyp = POLY_CMD;
715    res->data = reinterpret_cast<void *>(m);
716
717    return FALSE;
718  }
719
720  WerrorS("`leadmonom(<poly/vector>)` expected");
721  return TRUE;
722}
723
724/// Get leading component
725static BOOLEAN leadcomp(leftv res, leftv h)
726{
727  NoReturn(res);
728
729  if ((h!=NULL) && (h->Typ()==VECTOR_CMD || h->Typ()==POLY_CMD))
730  {
731    const ring r = currRing;
732
733    const poly p = (poly)(h->Data());
734
735    if (p != NULL )
736    {
737      assume( p != NULL );
738      assume( p_LmTest(p, r) );
739
740      const unsigned long iComp = p_GetComp(p, r);
741
742  //    assume( iComp > 0 ); // p is a vector
743
744      res->data = reinterpret_cast<void *>(jjLONG2N(iComp));
745    } else
746      res->data = reinterpret_cast<void *>(jjLONG2N(0));
747     
748
749    res->rtyp = BIGINT_CMD;
750    return FALSE;
751  }
752
753  WerrorS("`leadcomp(<poly/vector>)` expected");
754  return TRUE;
755}
756
757
758
759
760/// Get raw leading exponent vector
761static BOOLEAN leadrawexp(leftv res, leftv h)
762{
763  NoReturn(res);
764
765  if ((h!=NULL) && (h->Typ()==VECTOR_CMD || h->Typ()==POLY_CMD) && (h->Data() != NULL))
766  {
767    const ring r = currRing;
768    const poly p = (poly)(h->Data());
769
770    assume( p != NULL );
771    assume( p_LmTest(p, r) );
772
773    const int iExpSize = r->ExpL_Size;
774
775//    intvec *iv = new intvec(iExpSize);
776
777    lists l=(lists)omAllocBin(slists_bin);
778    l->Init(iExpSize);
779
780    for(int i = iExpSize-1; i >= 0; i--)
781    {
782      l->m[i].rtyp = BIGINT_CMD;
783      l->m[i].data = reinterpret_cast<void *>(jjLONG2N(p->exp[i])); // longs...
784    }
785
786    res->rtyp = LIST_CMD; // list of bigints
787    res->data = reinterpret_cast<void *>(l);
788    return FALSE;
789  }
790
791  WerrorS("`leadrawexp(<poly/vector>)` expected");
792  return TRUE;
793}
794
795
796/// Endowe the current ring with additional (leading) Syz-component ordering
797static BOOLEAN MakeSyzCompOrdering(leftv res, leftv /*h*/)
798{
799
800  NoReturn(res);
801
802  //    res->data = rCurrRingAssure_SyzComp(); // changes current ring! :(
803  res->data = reinterpret_cast<void *>(rAssure_SyzComp(currRing, TRUE));
804  res->rtyp = RING_CMD; // return new ring!
805  // QRING_CMD?
806
807  return FALSE;
808}
809
810
811/// Same for Induced Schreyer ordering (ordering on components is defined by sign!)
812static BOOLEAN MakeInducedSchreyerOrdering(leftv res, leftv h)
813{
814
815  NoReturn(res);
816
817  int sign = 1;
818  if ((h!=NULL) && (h->Typ()==INT_CMD))
819  {
820    const int s = (int)((long)(h->Data()));
821
822    if( s != -1 && s != 1 )
823    {
824      WerrorS("`MakeInducedSchreyerOrdering(<int>)` called with wrong integer argument (must be +-1)!");
825      return TRUE;
826    }
827
828    sign = s;           
829  }
830
831  assume( sign == 1 || sign == -1 );
832  res->data = reinterpret_cast<void *>(rAssure_InducedSchreyerOrdering(currRing, TRUE, sign));
833  res->rtyp = RING_CMD; // return new ring!
834  // QRING_CMD?
835  return FALSE;
836}
837
838
839/// Returns old SyzCompLimit, can set new limit
840static BOOLEAN SetSyzComp(leftv res, leftv h)
841{
842  NoReturn(res);
843
844  const ring r = currRing;
845
846  if( !rIsSyzIndexRing(r) )
847  {
848    WerrorS("`SetSyzComp(<int>)` called on incompatible ring (not created by 'MakeSyzCompOrdering'!)");
849    return TRUE;
850  }
851
852  res->rtyp = INT_CMD;
853  res->data = reinterpret_cast<void *>(rGetCurrSyzLimit(r)); // return old syz limit
854
855  if ((h!=NULL) && (h->Typ()==INT_CMD))
856  {
857    const int iSyzComp = (int)reinterpret_cast<long>(h->Data());
858    assume( iSyzComp > 0 );
859    rSetSyzComp(iSyzComp, currRing);
860  }
861
862  return FALSE;
863}
864
865/// ?
866static BOOLEAN GetInducedData(leftv res, leftv h)
867{
868  NoReturn(res);
869
870  const ring r = currRing;
871
872  int p = 0; // which IS-block? p^th!
873
874  if ((h!=NULL) && (h->Typ()==INT_CMD))
875  {
876    p = (int)((long)(h->Data())); h=h->next;
877    assume(p >= 0);
878  }
879
880  const int pos = rGetISPos(p, r);
881
882  if(  /*(*/ -1 == pos /*)*/  )
883  {
884    WerrorS("`GetInducedData([int])` called on incompatible ring (not created by 'MakeInducedSchreyerOrdering'!)");
885    return TRUE;
886  }
887
888
889  const int iLimit = r->typ[pos].data.is.limit;
890  const ideal F = r->typ[pos].data.is.F;
891  ideal FF = id_Copy(F, r);
892
893  lists l=(lists)omAllocBin(slists_bin);
894  l->Init(2);
895
896  l->m[0].rtyp = INT_CMD;
897  l->m[0].data = reinterpret_cast<void *>(iLimit);
898
899
900  //        l->m[1].rtyp = MODUL_CMD;
901
902  if( idIsModule(FF, r) )
903  {
904    l->m[1].rtyp = MODUL_CMD;
905
906    //          Print("before: %d\n", FF->nrows);
907    //          FF->nrows = id_RankFreeModule(FF, r); // ???
908    //          Print("after: %d\n", FF->nrows);
909  }
910  else
911    l->m[1].rtyp = IDEAL_CMD;
912
913  l->m[1].data = reinterpret_cast<void *>(FF);
914
915  res->rtyp = LIST_CMD; // list of int/module
916  res->data = reinterpret_cast<void *>(l);
917
918  return FALSE;
919
920}
921
922
923/// Returns old SyzCompLimit, can set new limit
924static BOOLEAN SetInducedReferrence(leftv res, leftv h)
925{
926  NoReturn(res);
927
928  const ring r = currRing;
929
930  if( !( (h!=NULL) && ( (h->Typ()==IDEAL_CMD) || (h->Typ()==MODUL_CMD))) )
931  {
932    WerrorS("`SetInducedReferrence(<ideal/module>, [int[, int]])` expected");
933    return TRUE;
934  }
935
936  const ideal F = (ideal)h->Data(); ; // No copy!
937  h=h->next;
938
939  int rank = 0;
940
941  if ((h!=NULL) && (h->Typ()==INT_CMD))
942  {
943    rank = (int)((long)(h->Data())); h=h->next;
944    assume(rank >= 0);
945  } else
946    rank = id_RankFreeModule(F, r); // Starting syz-comp (1st: i+1)
947
948  int p = 0; // which IS-block? p^th!
949
950  if ((h!=NULL) && (h->Typ()==INT_CMD))
951  {
952    p = (int)((long)(h->Data())); h=h->next;
953    assume(p >= 0);
954  }
955
956  const int posIS = rGetISPos(p, r);
957
958  if(  /*(*/ -1 == posIS /*)*/  )
959  {
960    WerrorS("`SetInducedReferrence(<ideal/module>, [int[, int]])` called on incompatible ring (not created by 'MakeInducedSchreyerOrdering'!)");
961    return TRUE;
962  }
963
964
965
966  // F & componentWeights belong to that ordering block of currRing now:
967  rSetISReference(r, F, rank, p); // F will be copied!
968  return FALSE;
969}
970
971
972//    F = ISUpdateComponents( F, V, MIN );
973//    // replace gen(i) -> gen(MIN + V[i-MIN]) for all i > MIN in all terms from F!
974static BOOLEAN ISUpdateComponents(leftv res, leftv h)
975{
976  NoReturn(res);
977
978  PrintS("ISUpdateComponents:.... \n");
979
980  if ((h!=NULL) && (h->Typ()==MODUL_CMD))
981  {
982    ideal F = (ideal)h->Data(); ; // No copy!
983    h=h->next;
984
985    if ((h!=NULL) && (h->Typ()==INTVEC_CMD))
986    {
987      const intvec* const V = (const intvec* const) h->Data();
988      h=h->next;
989
990      if ((h!=NULL) && (h->Typ()==INT_CMD))
991      {
992        const int MIN = (int)((long)(h->Data()));
993
994        pISUpdateComponents(F, V, MIN, currRing);
995        return FALSE;
996      }
997    }
998  }
999
1000  WerrorS("`ISUpdateComponents(<module>, intvec, int)` expected");
1001  return TRUE;
1002}
1003
1004
1005/// NF using length
1006static BOOLEAN reduce_syz(leftv res, leftv h)
1007{
1008  // const ring r = currRing;
1009
1010  if ( !( (h!=NULL) && (h->Typ()==VECTOR_CMD || h->Typ()==POLY_CMD) ) )
1011  {
1012    WerrorS("`reduce_syz(<poly/vector>!, <ideal/module>, <int>, [int])` expected");
1013    return TRUE;
1014  }
1015
1016  res->rtyp = h->Typ();
1017  const poly v = reinterpret_cast<poly>(h->Data());
1018  h=h->next;
1019
1020  if ( !( (h!=NULL) && (h->Typ()==MODUL_CMD || h->Typ()==IDEAL_CMD ) ) )
1021  {
1022    WerrorS("`reduce_syz(<poly/vector>, <ideal/module>!, <int>, [int])` expected");
1023    return TRUE;
1024  }
1025
1026  assumeStdFlag(h);
1027  const ideal M = reinterpret_cast<ideal>(h->Data()); h=h->next;
1028
1029
1030  if ( !( (h!=NULL) && (h->Typ()== INT_CMD)  ) )
1031  {
1032    WerrorS("`reduce_syz(<poly/vector>, <ideal/module>, <int>!, [int])` expected");
1033    return TRUE;
1034  }
1035
1036  const int iSyzComp = (int)((long)(h->Data())); h=h->next;
1037
1038  int iLazyReduce = 0;
1039
1040  if ( ( (h!=NULL) && (h->Typ()== INT_CMD)  ) )
1041    iLazyReduce = (int)((long)(h->Data())); 
1042
1043  res->data = (void *)kNFLength(M, currQuotient, v, iSyzComp, iLazyReduce);
1044  return FALSE;
1045}
1046
1047
1048/// Get raw syzygies (idPrepare)
1049static BOOLEAN idPrepare(leftv res, leftv h)
1050{
1051  //        extern int rGetISPos(const int p, const ring r);
1052
1053  const ring r = currRing;
1054
1055  const bool isSyz = rIsSyzIndexRing(r);
1056  const int posIS = rGetISPos(0, r);
1057
1058
1059  if ( !( (h!=NULL) && (h->Typ()==MODUL_CMD) && (h->Data() != NULL) ) )
1060  {
1061    WerrorS("`idPrepare(<module>)` expected");
1062    return TRUE;
1063  }
1064
1065  const ideal I = reinterpret_cast<ideal>(h->Data());
1066
1067  assume( I != NULL );
1068  idTest(I);
1069
1070  int iComp = -1;
1071
1072  h=h->next;
1073  if ( (h!=NULL) && (h->Typ()==INT_CMD) )
1074  {
1075    iComp = (int)((long)(h->Data()));
1076  } else
1077  {
1078      if( (!isSyz) && (-1 == posIS) )
1079      {
1080        WerrorS("`idPrepare(<...>)` called on incompatible ring (not created by 'MakeSyzCompOrdering' or 'MakeInducedSchreyerOrdering'!)");
1081        return TRUE;
1082      }
1083
1084    if( isSyz )
1085      iComp = rGetCurrSyzLimit(r);
1086    else
1087      iComp = id_RankFreeModule(r->typ[posIS].data.is.F, r); // ;
1088  }
1089 
1090  assume(iComp >= 0);
1091
1092
1093  intvec* w = reinterpret_cast<intvec *>(atGet(h, "isHomog", INTVEC_CMD));
1094  tHomog hom = testHomog;
1095
1096  //           int add_row_shift = 0;
1097  //
1098  if (w!=NULL)
1099  {
1100    w = ivCopy(w);
1101  //             add_row_shift = ww->min_in();
1102  //
1103  //             (*ww) -= add_row_shift;
1104  //             
1105  //             if (idTestHomModule(I, currQuotient, ww))
1106  //             {
1107    hom = isHomog;
1108  //               w = ww;
1109  //             }
1110  //             else
1111  //             {
1112  //               //WarnS("wrong weights");
1113  //               delete ww;
1114  //               w = NULL;
1115  //               hom=testHomog;
1116  //             }
1117  }
1118
1119
1120  // computes syzygies of h1,
1121  // works always in a ring with ringorder_s
1122  // NOTE: rSetSyzComp(syzcomp) should better be called beforehand
1123  //        ideal idPrepare (ideal  h1, tHomog hom, int syzcomp, intvec **w);
1124
1125  ideal J = // idPrepare( I, hom, iComp, &w);
1126           kStd(I, currQuotient, hom, &w, NULL, iComp);
1127
1128  idTest(J);
1129
1130  if (w!=NULL)
1131    atSet(res, omStrDup("isHomog"), w, INTVEC_CMD);
1132  //             if (w!=NULL) delete w;
1133
1134  res->rtyp = MODUL_CMD;
1135  res->data = reinterpret_cast<void *>(J);
1136  return FALSE;
1137}
1138
1139/// Get raw syzygies (idPrepare)
1140static BOOLEAN _p_Content(leftv res, leftv h)
1141{
1142  if ( !( (h!=NULL) && (h->Typ()==POLY_CMD) && (h->Data() != NULL) ) )
1143  {
1144    WerrorS("`p_Content(<poly-var>)` expected");
1145    return TRUE;
1146  }
1147
1148
1149  const poly p = reinterpret_cast<poly>(h->Data());
1150
1151 
1152  pTest(p);  pWrite(p); PrintLn();
1153
1154 
1155  p_Content( p, currRing);     
1156
1157  pTest(p);
1158  pWrite(p); PrintLn();
1159 
1160  NoReturn(res);
1161  return false;
1162}
1163
1164static BOOLEAN _m2_end(leftv res, leftv h)
1165{
1166  int ret = 0;
1167 
1168  if ( (h!=NULL) && (h->Typ()!=INT_CMD) )
1169  {
1170    WerrorS("`m2_end([<int>])` expected");
1171    return TRUE;
1172  }
1173  ret = (int)(long)(h->Data());
1174
1175  m2_end( ret );
1176
1177  NoReturn(res);
1178  return FALSE;
1179}
1180
1181   
1182
1183END_NAMESPACE
1184
1185
1186int SI_MOD_INIT(syzextra)(SModulFunctions* psModulFunctions) 
1187{
1188#define ADD0(A,B,C,D,E) A(B, (char*)C, D, E)
1189// #define ADD(A,B,C,D,E) ADD0(iiAddCproc, "", C, D, E)
1190  #define ADD(A,B,C,D,E) ADD0(A->iiAddCproc, B, C, D, E)
1191  ADD(psModulFunctions, currPack->libname, "ClearContent", FALSE, _ClearContent);
1192  ADD(psModulFunctions, currPack->libname, "ClearDenominators", FALSE, _ClearDenominators);
1193
1194  ADD(psModulFunctions, currPack->libname, "DetailedPrint", FALSE, DetailedPrint);
1195  ADD(psModulFunctions, currPack->libname, "leadmonomial", FALSE, leadmonom);
1196  ADD(psModulFunctions, currPack->libname, "leadcomp", FALSE, leadcomp);
1197  ADD(psModulFunctions, currPack->libname, "leadrawexp", FALSE, leadrawexp);
1198
1199  ADD(psModulFunctions, currPack->libname, "ISUpdateComponents", FALSE, ISUpdateComponents);
1200  ADD(psModulFunctions, currPack->libname, "SetInducedReferrence", FALSE, SetInducedReferrence);
1201  ADD(psModulFunctions, currPack->libname, "GetInducedData", FALSE, GetInducedData);
1202  ADD(psModulFunctions, currPack->libname, "SetSyzComp", FALSE, SetSyzComp);
1203  ADD(psModulFunctions, currPack->libname, "MakeInducedSchreyerOrdering", FALSE, MakeInducedSchreyerOrdering);
1204  ADD(psModulFunctions, currPack->libname, "MakeSyzCompOrdering", FALSE, MakeSyzCompOrdering);
1205
1206  ADD(psModulFunctions, currPack->libname, "ProfilerStart", FALSE, _ProfilerStart); ADD(psModulFunctions, currPack->libname, "ProfilerStop",  FALSE, _ProfilerStop );
1207 
1208  ADD(psModulFunctions, currPack->libname, "noop", FALSE, noop);
1209  ADD(psModulFunctions, currPack->libname, "idPrepare", FALSE, idPrepare);
1210  ADD(psModulFunctions, currPack->libname, "reduce_syz", FALSE, reduce_syz);
1211
1212  ADD(psModulFunctions, currPack->libname, "p_Content", FALSE, _p_Content);
1213
1214  ADD(psModulFunctions, currPack->libname, "Tail", FALSE, Tail);
1215  ADD(psModulFunctions, currPack->libname, "ComputeLeadingSyzygyTerms", FALSE, ComputeLeadingSyzygyTerms);
1216 
1217  //  ADD(psModulFunctions, currPack->libname, "", FALSE, );
1218
1219  ADD(psModulFunctions, currPack->libname, "m2_end", FALSE, _m2_end);
1220
1221#undef ADD 
1222  return 0;
1223}
1224
1225#ifndef EMBED_PYTHON
1226extern "C" { 
1227int mod_init(SModulFunctions* psModulFunctions)
1228{ 
1229  return SI_MOD_INIT(syzextra)(psModulFunctions); 
1230}
1231}
1232#endif
Note: See TracBrowser for help on using the repository browser.