source: git/Singular/mpsr_GetPoly.cc @ 0df59c8

spielwiese
Last change on this file since 0df59c8 was 762407, checked in by Oleksandr Motsak <motsak@…>, 12 years ago
config.h is for sources files only FIX: config.h should only be used by source (not from inside kernel/mod2.h!) NOTE: each source file should better include mod2.h right after config.h, while headers should better not include mod2.h.
  • Property mode set to 100644
File size: 25.7 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4
5/* $Id$ */
6
7/***************************************************************
8 *
9 * File:       mpsr_GetPoly.cc
10 * Purpose:    rotines which get polys and polynomails (i.e. ring) annotations
11 * Author:     Olaf Bachmann (2/97)
12 *
13 *
14 ***************************************************************/
15#include "config.h"
16#include <kernel/mod2.h>
17
18#ifdef HAVE_MPSR
19#include <misc/mylimits.h>
20
21#include <Singular/mpsr_Get.h>
22
23#include <coeffs/si_gmp.h>
24#include <omalloc/omalloc.h>
25#include <Singular/tok.h>
26#include <Singular/ipid.h>
27#include <polys/monomials/ring.h>
28//#include "kernel/longalg.h"
29#include <polys/monomials/maps.h>
30#include <kernel/ideals.h>
31#include <Singular/grammar.h>
32#include <kernel/febase.h>
33#include "kernel/modulop.h"
34
35#include <Singular/mpsr_Tok.h>
36
37#ifdef PARI_BIGINT_TEST
38#include <Singular/MP_PariBigInt.h>
39
40MP_Status_t IMP_MyGetApInt(MP_Link_pt link, MP_ApInt_t *apint)
41{
42  GEN pnum;
43  mp_failr(IMP_GetApInt(link, (MP_ApInt_t *) &pnum));
44  _pari_to_gmp(pnum, (mpz_ptr *) apint);
45
46  return MP_Success;
47}
48
49#else
50
51#define IMP_MyGetApInt IMP_GetApInt
52
53#endif
54
55
56/***************************************************************
57 *
58 * global variable definitions
59 *
60 ***************************************************************/
61
62static mpsr_Status_t (*GetCoeff)(MP_Link_pt link, number *x);
63static mpsr_Status_t (*GetAlgNumberNumber)(MP_Link_pt link, number *x);
64static MP_Sint32_t gNalgvars = 0;
65static MP_Sint32_t gNvars = 0;
66static ring        currGetRing = NULL;
67
68
69/***************************************************************
70 *
71 * prototype declarations
72 *
73 ***************************************************************/
74static void        SetGetFuncs(ring r);
75static mpsr_Status_t GetModuloNumber(MP_Link_pt link, number *a);
76static mpsr_Status_t GetGaloisNumber(MP_Link_pt link, number *a);
77static mpsr_Status_t GetFloatNumber(MP_Link_pt link, number *a);
78static mpsr_Status_t GetApInt(MP_Link_pt link, mpz_ptr ap);
79static mpsr_Status_t GetRationalNumber(MP_Link_pt link, number *a);
80static mpsr_Status_t GetAlgNumber(MP_Link_pt link, number *a);
81static mpsr_Status_t GetVarNumberAnnot(MPT_Node_pt node, ring r, BOOLEAN mv);
82static mpsr_Status_t GetProtoTypeAnnot(MPT_Node_pt node, ring r, BOOLEAN mv,
83                                     ring &subring);
84static mpsr_Status_t GetOrderingAnnot(MPT_Node_pt node, ring r, BOOLEAN mv,
85                                      BOOLEAN &IsUnOrdered);
86static mpsr_Status_t GetSimpleOrdering(MPT_Node_pt node, ring r, short i);
87static mpsr_Status_t GetVarNamesAnnot(MPT_Node_pt node, ring r);
88static mpsr_Status_t GetDefRelsAnnot(MPT_Node_pt node, ring r);
89
90/***************************************************************
91 *
92 * Setting the global Put Functions
93 *
94 ***************************************************************/
95
96static inline BOOLEAN IsCurrGetRing(ring r)
97{
98  if (r == currGetRing && r == currRing) return 1;
99  else return 0;
100}
101
102static void SetGetFuncs(ring r)
103{
104  currGetRing = r;
105  // first, we set the PutNumber function
106  gNvars = r->N;
107  mpsr_InitTempArray(gNvars + 1);
108
109  if (rField_is_Q(r))
110    // rational numbers
111    GetCoeff = GetRationalNumber;
112  else if (rField_is_Zp(r))
113    GetCoeff = GetModuloNumber;
114  else if (rField_is_GF(r))
115      GetCoeff = GetGaloisNumber;
116  else if (rField_is_R(r))
117    GetCoeff = GetFloatNumber;
118  else
119  {
120    // now we come to algebraic numbers
121    gNalgvars = rPar(r);
122    mpsr_InitTempArray(gNalgvars);
123    GetCoeff = GetAlgNumber;
124    if (rField_is_Zp_a(r))
125      // first, Z/p(a)
126      GetAlgNumberNumber = GetModuloNumber;
127    else
128      GetAlgNumberNumber = GetRationalNumber;
129  }
130
131  // still need to set the global ring
132  mpsr_SetCurrRing(r, TRUE);
133}
134
135
136
137/***************************************************************
138 *
139 * Routines for Getting coeffs
140 *
141 ***************************************************************/
142// we always Get modulo numbers without a node, since
143// we have type-spec this before
144static mpsr_Status_t GetModuloNumber(MP_Link_pt link, number *a)
145{
146  MP_Uint32_t x;
147  mp_failr(IMP_GetUint32(link, &x));
148  *a=npInit((int)x, currRing);
149  return mpsr_Success;
150}
151
152static mpsr_Status_t GetGaloisNumber(MP_Link_pt link, number *a)
153{
154  mp_return(IMP_GetUint32(link, (MP_Uint32_t *) a));
155}
156
157static mpsr_Status_t GetFloatNumber(MP_Link_pt link, number *a)
158{
159   mp_return( IMP_GetReal32(link , (MP_Real32_t *) a));
160}
161
162static mpsr_Status_t GetApInt(MP_Link_pt link, mpz_ptr ap)
163{
164  MP_NodeType_t     node;
165  MP_DictTag_t      dict;
166  MP_NumChild_t     num_child;
167  MP_NumAnnot_t     num_annots;
168  MP_Common_t       cvalue;
169  MP_Boolean_t      req = 0;
170
171  mp_failr(IMP_GetNodeHeader(link,&node,&dict, &cvalue, &num_annots,
172                             &num_child));
173
174  if (node == MP_ApIntType)
175  {
176    mpz_init(ap);
177    mp_failr(IMP_MyGetApInt(link, (MP_ApInt_t *) &ap));
178  }
179  else if (node == MP_Sint32Type || node == MP_Sint8Type)
180  {
181    MP_Sint32_t i;
182    if (node == MP_Sint8Type)
183      i = (int) ((signed char) cvalue);
184    else
185      mp_failr(IMP_GetSint32(link, &i));
186    mpz_init_set_si((mpz_ptr) ap, i);
187  }
188  else if (node == MP_Uint32Type || node == MP_Uint8Type)
189  {
190    MP_Uint32_t i;
191    if (node == MP_Uint8Type)
192      i = cvalue;
193    else
194      mp_failr(IMP_GetUint32(link, &i));
195    mpz_init_set_ui((mpz_ptr) ap, i);
196  }
197  else
198    return mpsr_SetError(mpsr_WrongNodeType);
199
200  if (num_annots > 0)
201  {
202    mpt_failr(MPT_SkipAnnots(link, num_annots, &req));
203    if (req) return mpsr_SetError(mpsr_ReqAnnotSkip);
204  }
205
206  return mpsr_Success;
207}
208
209// This supposes that number is of char 0, i.e. a rational number
210static mpsr_Status_t GetRationalNumber(MP_Link_pt link, number *x)
211{
212  MP_NodeType_t     node;
213  MP_DictTag_t      dict;
214  MP_NumChild_t     num_child;
215  MP_NumAnnot_t     num_annots;
216  MP_Sint32_t       i;
217  MP_Common_t       cvalue;
218  number            y;
219  MP_Boolean_t      req;
220
221  mp_failr(IMP_GetNodeHeader(link,&node,&dict, &cvalue, &num_annots,
222                             &num_child));
223
224  // start with the most frequent cases
225  if (node == MP_Sint32Type)
226  {
227    mp_failr(IMP_GetSint32(link, &i));
228    *x = nlInit(i, currRing);
229  }
230  else if (node == MP_ApIntType)
231  {
232    mpz_ptr gnum;
233    y =  ALLOC0_RNUMBER();
234#if defined(LDEBUG)
235    y->debug = 123456;
236#endif
237    y->s = 3;
238    gnum = y->z;
239    mpz_init(gnum);
240    mp_failr(IMP_MyGetApInt(link, (MP_ApInt_t *) &gnum));
241    nlNormalize(y);
242    *x = y;
243  }
244  // fraction of numbers
245  else if (node == MP_CommonOperatorType &&
246           dict== MP_BasicDict &&
247           cvalue == MP_CopBasicDiv)
248  {
249    if (num_annots > 0)
250    {
251      mpt_failr(MPT_SkipAnnots(link, num_annots, &req));
252      if (req) return mpsr_SetError(mpsr_ReqAnnotSkip);
253    }
254    *x =  (number) ALLOC0_RNUMBER();
255    y = (number) *x;
256#if defined(LDEBUG)
257    y->debug = 123456;
258#endif
259    y->s = 1;
260    failr(GetApInt(link, y->z));
261    return GetApInt(link, y->n);
262  }
263  // check for some more esoteric cases
264  else if (node == MP_Uint8Type)
265    *x = nlInit(cvalue, currRing);
266  else if (node == MP_Sint8Type)
267    // be careful -- need to handle the value "-2", for example
268    *x = nlInit((int) ((MP_Sint8_t) cvalue), currRing);
269  else if (node == MP_Uint32Type)
270  {
271    MP_Uint32_t ui;
272    mp_failr(IMP_GetUint32(link, &ui));
273    // check whether u_int can be casted safely to int
274    if (ui < MAX_INT_VAL)
275      *x = nlInit(ui, currRing);
276    else
277    {
278      // otherwise, make an apint out of it
279      *x =  (number) ALLOC0_RNUMBER();
280      y = (number) *x;
281#if defined(LDEBUG)
282      y->debug = 123456;
283#endif
284      mpz_init_set_ui(y->z, ui);
285      y->s = 3;
286    }
287  }
288  else
289    return mpsr_SetError(mpsr_WrongNodeType);
290
291  if (num_annots > 0)
292  {
293    mpt_failr(MPT_SkipAnnots(link, num_annots, &req));
294    if (req) return mpsr_SetError(mpsr_ReqAnnotSkip);
295  }
296
297  return mpsr_Success;
298}
299
300/***************************************************************
301 *
302 * Algebraic Numbers (a la Singular)
303 *
304 ***************************************************************/
305static inline mpsr_Status_t GetAlgPoly(MP_Link_pt link, napoly *p)
306{
307  MP_Uint32_t j, nm;
308  int i;
309  napoly a;
310  int *exp;
311
312  IMP_GetUint32(link, &nm);
313
314  if (nm == 0)
315  {
316    *p = NULL;
317    return mpsr_Success;
318  }
319  a = napNew();
320  *p = a;
321
322  failr(GetAlgNumberNumber(link, &(napGetCoeff(a))));
323  mp_failr(IMP_GetSint32Vector(link, (MP_Sint32_t **) &gTa, rPar(currRing)));
324  for (i=0; i<rPar(currRing); i++)
325    napSetExp(a,i+1,gTa[i]);
326
327  for (j=1; j<nm; j++)
328  {
329    pNext(a) = napNew();
330    pIter(a);
331    failr(GetAlgNumberNumber(link, &(napGetCoeff(a))));
332    mp_failr(IMP_GetSint32Vector(link, (MP_Sint32_t **) &gTa, rPar(currRing)));
333    for (i=0; i<rPar(currRing); i++)
334      napSetExp(a,i+1,gTa[i]);
335  }
336  pNext(a) = NULL;
337
338  return mpsr_Success;
339}
340
341static mpsr_Status_t GetAlgNumber(MP_Link_pt link, number *a)
342{
343  lnumber b;
344  MP_Uint32_t ut;
345
346  // Get the union tag
347  mp_failr(IMP_GetUint32(link, &ut));
348  if (ut == 0)
349  {
350    *a = NULL;
351    return mpsr_Success;
352  }
353  else if (ut == 1 || ut == 2)
354  {
355    // single number
356    b = (lnumber) ALLOC0_RNUMBER();
357    *a = (number) b;
358    failr(GetAlgPoly(link, &(b->z)));
359    if (ut == 2)
360      return GetAlgPoly(link, &(b->n));
361    else
362      return mpsr_Success;
363  }
364  else
365    return mpsr_SetError(mpsr_WrongUnionDiscriminator);
366}
367
368/***************************************************************
369 *
370 *  Getting polys
371 *
372 ***************************************************************/
373mpsr_Status_t mpsr_GetPoly(MP_Link_pt link, poly &p, MP_Uint32_t nmon,
374                         ring cring)
375{
376  poly pp;
377  MP_Sint32_t i;
378  MP_Uint32_t j;
379
380  if (! IsCurrGetRing(cring))
381    SetGetFuncs(cring);
382
383  if (nmon == 0)
384  {
385    p = NULL;
386    return mpsr_Success;
387  }
388
389  pp = pInit();
390  p = pp;
391  failr(GetCoeff(link, &(pp->coef)));
392  if (gNvars > 1)
393  {
394    MP_Sint32_t* Ta = &gTa[1];
395    mp_failr(IMP_GetSint32Vector(link, &Ta, gNvars));
396    for (i=1; i<=gNvars; i++)
397      pSetExp(pp,i , gTa[i]);
398    pSetm(pp);
399
400    for (j=1; j<nmon; j++)
401    {
402      pp->next = pInit();
403      pp = pp->next;
404      failr(GetCoeff(link, &(pp->coef)));
405      mp_failr(IMP_GetSint32Vector(link, &Ta, gNvars));
406      for (i=1; i<=gNvars; i++)
407        pSetExp(pp, i, gTa[i]);
408      pSetm(pp);
409    }
410  }
411  else
412  {
413    mp_failr(IMP_GetSint32(link, &i));
414    pSetExp(pp,1, i);
415    pSetm(pp);
416
417    for (j=1; j<nmon; j++)
418    {
419      pp->next = pInit();
420      pp = pp->next;
421      failr(GetCoeff(link, &(pp->coef)));
422      mp_failr(IMP_GetSint32(link, &i));
423      pSetExp(pp,1, i);
424      pSetm(pp);
425    }
426  }
427
428  pp->next = NULL;
429
430  return mpsr_Success;
431}
432
433mpsr_Status_t mpsr_GetPolyVector(MP_Link_pt link, poly &p, MP_Uint32_t nmon,
434                          ring cring)
435{
436  poly pp;
437  MP_Sint32_t i, n1;
438  MP_Uint32_t j;
439
440  if (!IsCurrGetRing(cring))
441    SetGetFuncs(cring);
442
443  n1 = gNvars + 1;
444  if (nmon == 0)
445  {
446    p = NULL;
447    return mpsr_Success;
448  }
449
450  pp = pInit();
451  p = pp;
452  failr(GetCoeff(link, &(pp->coef)));
453  if (gNvars > 1)
454  {
455    mp_failr(IMP_GetSint32Vector(link, &gTa, n1));
456    pSetComp(pp, gTa[0]);
457    for (i=1; i<n1; i++)
458      pSetExp(pp, i, gTa[i]);
459    pSetm(pp);
460
461    for (j=1; j<nmon; j++)
462    {
463      pp->next = pInit();
464      pp = pp->next;
465      failr(GetCoeff(link, &(pp->coef)));
466      mp_failr(IMP_GetSint32Vector(link, &gTa, n1));
467      pSetComp(pp, gTa[0]);
468      for (i=1; i<n1; i++)
469        pSetExp(pp,i, gTa[i]);
470      pSetm(pp);
471    }
472  }
473  else
474  {
475    mp_failr(IMP_GetSint32(link, &i));
476    pSetComp(pp, i);
477    mp_failr(IMP_GetSint32(link, &i));
478    pSetExp(pp,1, i);
479    pSetm(pp);
480
481    for (j=1; j<nmon; j++)
482    {
483      pp->next = pInit();
484      pp = pp->next;
485      failr(GetCoeff(link, &(pp->coef)));
486      mp_failr(IMP_GetSint32(link, &i));
487      pSetComp(pp, i);
488      mp_failr(IMP_GetSint32(link, &i));
489      pSetExp(pp,1, i);
490      pSetm(pp);
491    }
492  }
493  pp->next = NULL;
494
495  return mpsr_Success;
496}
497
498/***************************************************************
499 *
500 *  The Getting annotation buisness
501 *
502 ***************************************************************/
503#define falser(x)                               \
504do                                              \
505{                                               \
506  if (!(x)) return mpsr_Failure;                   \
507}                                               \
508while (0)
509
510// We assume that the node is that of a DDP: This returns
511// MP_Succcess, if annots of node can be used to construct a
512// Singular ring (in which case r is the respective ring) or,
513// MP_Failure, if not
514mpsr_Status_t mpsr_GetRingAnnots(MPT_Node_pt node, ring &r,
515                                 BOOLEAN &mv, BOOLEAN &IsUnOrdered)
516{
517  sip_sring r1, *subring;
518  poly minpoly = NULL;
519
520  memset(&r1, 0, sizeof(sip_sring));
521
522  r = NULL;
523  if (MPT_Annot(node, MP_PolyDict, MP_AnnotPolyModuleVector) != NULL)
524    mv = 1;
525  else
526    mv = 0;
527
528  // sets r->N
529  if (GetVarNumberAnnot(node, &r1, mv) != mpsr_Success)
530    Warn("GetVarNumberAnnot: using the one found in the prototype");
531
532  // sets r->char and r->minpoly, r->parameter; if necessary
533  failr(GetProtoTypeAnnot(node, &r1, mv, subring));
534
535  // if we are still here, then we are successful in constructing the ring
536  r = (ring) omAllocBin(sip_sring_bin);
537  memcpy(r, &r1, sizeof(sip_sring));
538
539  if (GetVarNamesAnnot(node, r) != mpsr_Success)
540    Warn("GetVarNamesAnnot: using default variable names");
541
542  if (GetOrderingAnnot(node,r, mv, IsUnOrdered) != mpsr_Success)
543    Warn("GetOrderingAnnot: using unspec ordering");
544
545  rComplete(r);
546
547  if (GetDefRelsAnnot(node, r) != mpsr_Success)
548    Warn("GetDefRelsAnnot: using no defining relations");
549
550  // check on whether or not I have to set a minpoly
551  if (subring != NULL)
552  {
553    if ((subring->qideal != NULL) &&
554        ((minpoly = subring->qideal->m[0]) != NULL))
555    {
556      mpsr_SetCurrRing(r, TRUE);
557      minpoly = maIMap(subring, minpoly);
558      r->minpoly = minpoly->coef;
559      pLmFree(minpoly);
560    }
561    rKill(subring);
562  }
563
564  // complete ring constructions
565  return mpsr_Success;
566}
567
568
569static mpsr_Status_t GetVarNumberAnnot(MPT_Node_pt node, ring r, BOOLEAN mv)
570{
571  MPT_Annot_pt annot = MPT_Annot(node, MP_PolyDict, MP_AnnotPolyVarNumber);
572
573  if (annot != NULL)
574  {
575    if (annot->value != NULL && annot->value->node->type == MP_Uint32Type)
576    {
577      // Hm.. should check that r->N is not too big for Singular
578      r->N = (short) MP_UINT32_T(annot->value->node->nvalue);
579      if (mv) (r->N)--;
580      return mpsr_Success;
581    }
582  }
583  return mpsr_Failure;
584}
585
586
587static mpsr_Status_t GetProtoTypeAnnot(MPT_Node_pt node, ring r, BOOLEAN mv,
588                                       ring &subring)
589{
590  MPT_Annot_pt annot = NULL;
591  MPT_Tree_pt  val;
592  MPT_Tree_pt  *ta;
593
594  subring = NULL;
595
596  // look for prototype annot
597  if ((val = MPT_ProtoAnnotValue(node)) == NULL)
598    return mpsr_Failure;
599
600  // check value of annot
601  node = val->node;
602  if (! (NodeCheck(node, MP_CommonOperatorType, MP_ProtoDict,
603                   MP_CopProtoStruct) && node->numchild == 2))
604    return mpsr_Failure;
605  // get the two args of the value
606  ta = (MPT_Tree_pt *) val->args;
607
608
609  // We get the exponent vector specification first
610  node = ta[1]->node;
611  if (! (NodeCheck(node, MP_CommonMetaOperatorType, MP_ProtoDict,
612                   MP_CopProtoArray) && node->numchild > 0))
613    return mpsr_Failure;
614  // check r->N and reset, if necessary
615  if (mv)
616  {
617    if (r->N != (int) (node->numchild - 1))
618    {
619      Warn("GetProtoAnnot: Inconsistent NumVars specification");
620      r->N = (node->numchild -1);
621    }
622  }
623  else
624  {
625    if (r->N != (int) node->numchild)
626    {
627      Warn("GetProtoAnnot: Inconsistent NumVars specification");
628      r->N = (node->numchild);
629    }
630  }
631  // check for type of exponent
632  if ((val = MPT_ProtoAnnotValue(node)) == NULL)
633    return mpsr_Failure;
634
635  node = val->node;
636  falser(NodeCheck(node, MP_CommonMetaType, MP_ProtoDict, MP_CmtProtoIMP_Sint32));
637
638  // consider the first arg -- which specify the coeffs
639  val = ta[0];
640  node = val->node;
641  if (node->type == MP_CommonMetaType)
642  {
643    // char 0
644    if (MP_COMMON_T(node->nvalue) == MP_CmtNumberRational &&
645        node->dict == MP_NumberDict)
646    {
647      r->ch = 0;
648      // Hmm ... we should check for the normalized annot
649    }
650    else if (MP_COMMON_T(node->nvalue) == MP_CmtProtoIMP_Uint32 &&
651             node->dict == MP_ProtoDict &&
652             (annot = MPT_Annot(node,MP_NumberDict,MP_AnnotNumberModulos))
653              != NULL)
654    {
655      // char p || GF(p,n)
656      falser(annot->value != NULL &&
657             annot->value->node->type == MP_Uint32Type);
658      r->ch = MP_UINT32_T(annot->value->node->nvalue);
659
660      if (MPT_Annot(annot->value->node,
661                        MP_NumberDict, MP_AnnotNumberIsPrime) == NULL)
662      {
663        // GF(p,n)
664        falser((annot = MPT_Annot(annot->value->node, 129,
665                                  MP_AnnotSingularGalois)) != NULL &&
666           (annot->value != NULL) &&
667           (annot->value->node->type == MP_StringType));
668        r->parameter = (char **)omAllocBin(char_ptr_bin);
669        r->parameter[0] = omStrDup(MP_STRING_T(annot->value->node->nvalue));
670        r->P = 1;
671      }
672    }
673    else if (MP_COMMON_T(node->nvalue) == MP_CmtProtoIMP_Real32 &&
674             node->dict == MP_ProtoDict)
675    {
676      // floats
677      r->ch = -1;
678    }
679    else
680      return mpsr_SetError(mpsr_UnknownCoeffDomain);
681
682    return mpsr_Success;
683  }
684  else
685  {
686    // alg numbers
687    BOOLEAN mv2, IsUnOrdered;
688    int i;
689
690    // DDP Frac Node check
691    falser(NodeCheck(node, MP_CommonMetaOperatorType, MP_BasicDict,
692                    MP_CopBasicDiv) &&
693           node->numchild == 0);
694    falser((val = MPT_ProtoAnnotValue(node)) != NULL);
695    node = val->node;
696    mpsr_assume(node != NULL);
697    falser(NodeCheck(node, MP_CommonMetaOperatorType, MP_PolyDict,
698                     MP_CopPolyDenseDistPoly) &&
699           node->numchild == 0);
700    // GetRingAnnots
701    failr(mpsr_GetRingAnnots(node, subring, mv2, IsUnOrdered));
702    // Check whether the ring can be "coerced" to an algebraic number
703    falser( (rField_is_Zp(subring)||rField_is_Q(subring)) &&
704           // orig: subring->ch >= 0 &&a ???
705           subring->order[0] == ringorder_lp &&
706           subring->order[2] == 0 &&
707           mv2 == FALSE &&
708           IsUnOrdered == FALSE);
709
710    // Now do the coercion
711    r->ch = (rField_is_Q(subring) ? 1 : - rChar(subring));
712    r->parameter = (char **) omAlloc((subring->N)*sizeof(char*));
713    r->P = subring->N;
714    for (i=0; i < subring->N; i++)
715      r->parameter[i] = omStrDup(subring->names[i]);
716
717    // everything is ok
718    return mpsr_Success;
719  }
720}
721
722static mpsr_Status_t GetVarNamesAnnot(MPT_Node_pt node, ring r)
723{
724  MPT_Annot_pt annot = MPT_Annot(node, MP_PolyDict, MP_AnnotPolyVarNames);
725  short num_vars = 0, N, lb, offset, nc;
726
727  mpsr_assume(r != NULL);
728  N = r->N;
729  r->names = (char **) omAlloc0(N * sizeof(char *));
730
731  // fill in varnames from the back
732  if (annot != NULL && annot->value != NULL)
733  {
734    node = annot->value->node;
735    nc = (short) node->numchild;
736    if (NodeCheck(node, MP_CommonOperatorType, MP_ProtoDict, MP_CopProtoArray))
737    {
738      MPT_Tree_pt val = MPT_ProtoAnnotValue(node);
739      if (val != NULL &&
740          NodeCheck(val->node, MP_CommonMetaType, MP_ProtoDict,
741                    MP_CmtProtoIMP_Identifier))
742      {
743        MPT_Arg_pt arg_pt = annot->value->args;
744        lb = si_min(nc, N);
745        offset = N - (short) nc;
746        if (offset < 0) offset = 0;
747        for (; num_vars < lb; num_vars++)
748          r->names[offset + num_vars] =
749            omStrDup(MP_STRING_T(arg_pt[num_vars]));
750      }
751    }
752    else if (node->type == MP_IdentifierType)
753    {
754      r->names[N-1] = omStrDup(MP_STRING_T(annot->value->node->nvalue));
755      num_vars = 1;
756    }
757  }
758
759  // fill in all remaining varnames
760  if (num_vars < N)
761  {
762    char vn[10];
763    offset = N - num_vars;
764    for (nc = 0; nc < offset; nc++)
765    {
766      sprintf(vn, "x(%d)", nc);
767      r->names[nc] = omStrDup(vn);
768    }
769  }
770
771  if (num_vars < N) return mpsr_Failure;
772  else return mpsr_Success;
773}
774
775static mpsr_Status_t GetOrderingAnnot(MPT_Node_pt node, ring r,
776                                      BOOLEAN mv, BOOLEAN &IsUnOrdered)
777{
778  MPT_Annot_pt annot = MPT_Annot(node, MP_PolyDict,
779                                 MP_AnnotShouldHavePolyOrdering);
780  IsUnOrdered = FALSE;
781  mpsr_Status_t status = mpsr_Success;
782  if (annot == NULL)
783  {
784    annot = MPT_Annot(node, MP_PolyDict,MP_AnnotPolyOrdering);
785    if (annot == NULL) status = mpsr_Failure;
786  }
787  else
788  {
789    IsUnOrdered = TRUE;
790  }
791
792
793  if (status == mpsr_Success) node =  annot->value->node;
794
795  // Check for BlockOrdering
796  if (status == mpsr_Success &&
797      NodeCheck(annot->value->node, MP_CommonOperatorType,
798               MP_BasicDict, MP_CopBasicList))
799  {
800    MP_NumChild_t nc = node->numchild, i;
801    MPT_Tree_pt *tarray = (MPT_Tree_pt *) annot->value->args, *tarray2, tree;
802
803    if (! mv) nc += 2; else nc++;
804    r->block0 = (int *) omAlloc0(nc*sizeof(int *));
805    r->block1 = (int *) omAlloc0(nc*sizeof(int *));
806    r->wvhdl  = (int **) omAlloc0(nc*sizeof(int *));
807    r->order  = (int *) omAlloc0(nc*sizeof(int *));
808
809    if (! mv)
810    {
811      r->order[nc-2] = ringorder_C;
812      nc = nc - 2;
813    }
814    else
815      nc--;
816
817    for (i=0; i<nc; i++)
818    {
819      tree = tarray[i];
820      if (NodeCheck(tree->node, MP_CommonOperatorType,
821                   MP_BasicDict, MP_CopBasicList) &&
822          tree->node->numchild == 3)
823      {
824        tarray2 = (MPT_Tree_pt *) tree->args;
825        if (GetSimpleOrdering(tarray2[0]->node, r, i) != mpsr_Success ||
826            tarray2[1]->node->type != MP_Uint32Type ||
827            tarray2[2]->node->type != MP_Uint32Type)
828        {
829          status = mpsr_Failure;
830          break;
831        }
832        else
833        {
834          r->block0[i] = MP_SINT32_T(tarray2[1]->node->nvalue);
835          r->block1[i] = MP_SINT32_T(tarray2[2]->node->nvalue);
836        }
837      }
838      else
839      {
840          status = mpsr_Failure;
841          break;
842      }
843    }
844
845    if (status == mpsr_Success) status = mpsr_rSetOrdSgn(r);
846
847    // Clean up if sth went wrong
848    if (status == mpsr_Failure)
849    {
850      if (mv) nc++;
851      else nc += 2;
852      omFreeSize(r->block0, nc*sizeof(int *));
853      omFreeSize(r->block1, nc*sizeof(int *));
854      omFreeSize(r->order, nc*sizeof(int *));
855      omFreeSize(r->wvhdl, nc*sizeof(short *));
856    }
857    else
858      return mpsr_Success;
859  }
860
861  // Either Simple Ordering, or sth failed from before
862  r->wvhdl = (int **)omAlloc0(3 * sizeof(int *));
863  r->order = (int *) omAlloc0(3 * sizeof(int *));
864  r->block0 = (int *)omAlloc0(3 * sizeof(int *));
865  r->block1 = (int *)omAlloc0(3 * sizeof(int *));
866  r->order[1] = ringorder_C;
867  r->block0[0] = 1;
868  r->block1[0] = r->N;
869
870  // Check for simple Ordering
871  if (status == mpsr_Success)
872    status = GetSimpleOrdering(node, r, 0);
873
874  if (status != mpsr_Success)
875  {
876    r->order[0] = ringorder_unspec;
877    IsUnOrdered = FALSE;
878  }
879
880  return mpsr_rSetOrdSgn(r);
881}
882
883static mpsr_Status_t GetSimpleOrdering(MPT_Node_pt node, ring r, short i)
884{
885  if (node->type != MP_CommonConstantType)
886    return mpsr_Failure;
887
888  int sr_ord =  mpsr_mp2ord(MP_COMMON_T(node->nvalue));
889
890  r->order[i] = sr_ord;
891  if (r->order[i] == ringorder_unspec) return mpsr_Failure;
892
893  MPT_Annot_pt annot = MPT_Annot(node, MP_PolyDict, MP_AnnotPolyWeights);
894
895  if (annot == NULL) return mpsr_Success;
896  if (annot->value == NULL) return mpsr_Failure;
897
898  node = annot->value->node;
899  if (r->order[i] == ringorder_M)
900  {
901    if (! NodeCheck(node, MP_CommonOperatorType, MP_MatrixDict,
902                   MP_CopMatrixDenseMatrix))
903      return mpsr_Failure;
904  }
905  else
906  {
907    if (! NodeCheck(node, MP_CommonOperatorType, MP_MatrixDict,
908                    MP_CopMatrixDenseVector))
909      return mpsr_Failure;
910    if (sr_ord == ringorder_lp) r->order[i] = ringorder_Wp;
911    else if (sr_ord == ringorder_ls) r->order[i] = ringorder_Ws;
912    else if (sr_ord != ringorder_wp && sr_ord != ringorder_ws &&
913             sr_ord != ringorder_a)
914      return mpsr_Failure;
915  }
916
917  MPT_Annot_pt
918    annot2 = MPT_Annot(node, MP_ProtoDict, MP_AnnotProtoPrototype);
919
920  if (annot2 == NULL ||
921      ! NodeCheck(annot2->value->node, MP_CommonMetaType, MP_ProtoDict,
922                 MP_CmtProtoIMP_Sint32))
923    return mpsr_Failure;
924
925  MP_Uint32_t nc = node->numchild, j;
926  MP_Sint32_t *w = (MP_Sint32_t *) annot->value->args;
927  int *w2 = (int *) omAlloc(nc*sizeof(int));
928
929  r->wvhdl[i] = w2;
930  for (j = 0; j < nc ; j++)
931    w2[j] = w[j];
932
933  return mpsr_Success;
934}
935
936static mpsr_Status_t GetDefRelsAnnot(MPT_Node_pt node, ring r)
937{
938  MPT_Annot_pt annot = MPT_Annot(node, MP_PolyDict, MP_AnnotPolyDefRel);
939  mpsr_leftv mlv;
940  leftv lv;
941  ring r1;
942
943  if (annot == NULL) return mpsr_Success;
944
945  node = annot->value->node;
946  if (node->type != MPT_ExternalDataType) return mpsr_Failure;
947
948  mlv = (mpsr_leftv) annot->value->args;
949  r1 = mlv->r;
950  lv = mlv->lv;
951
952  if (! rEqual(r1, r)) return mpsr_Failure;
953
954  if (lv->rtyp == POLY_CMD)
955  {
956    r->qideal = idInit(1,1);
957    r->qideal->m[0] = (poly) lv->data;
958    lv->data = NULL;
959  }
960  else if (lv->rtyp == IDEAL_CMD)
961  {
962    r->qideal = (ideal) lv->data;
963    lv->data = NULL;
964  }
965  else return mpsr_Failure;
966
967  return mpsr_Success;
968}
969
970extern mpsr_Status_t mpsr_rSetOrdSgn(ring r)
971{
972  short i = 0, order;
973  r->OrdSgn = 1;
974
975  while ((order = r->order[i]) != ringorder_no)
976  {
977    if (order == ringorder_ls ||
978        order == ringorder_Ws ||
979        order == ringorder_ws ||
980        order == ringorder_Ds ||
981        order == ringorder_ds)
982    {
983      r->OrdSgn = -1;
984      return mpsr_Success;
985    }
986    if (order == ringorder_M)
987    {
988      int sz = r->block1[i] - r->block0[i] + 1, j, k=0;
989      int *matrix = r->wvhdl[i];
990
991      while (k < sz)
992      {
993        j = 0;
994        while ((j < sz) && matrix[j*sz+k]==0) j++;
995        if (j>=sz)
996        {
997          Warn("Matrix order not complete");
998          r->OrdSgn = 0;
999          return mpsr_Failure;
1000        }
1001        else if (matrix[j*sz+k]<0)
1002        {
1003          r->OrdSgn = -1;
1004          return mpsr_Success;
1005        }
1006        else
1007          k++;
1008      }
1009    }
1010    i++;
1011  }
1012  return mpsr_Success;
1013}
1014#endif
Note: See TracBrowser for help on using the repository browser.