source: git/Singular/mpsr_GetPoly.cc @ 1629ab

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