source: git/Singular/mpsr_GetMisc.cc @ 8a150b

spielwiese
Last change on this file since 8a150b was 8a150b, checked in by Hans Schönemann <hannes@…>, 25 years ago
* hannes: added long reals git-svn-id: file:///usr/local/Singular/svn/trunk@3012 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 12.6 KB
Line 
1/****************************************
2*  Computer Algebra System SINGULAR     *
3****************************************/
4/* $Id: mpsr_GetMisc.cc,v 1.14 1999-04-29 11:38:52 Singular Exp $ */
5
6/***************************************************************
7 *
8 * File:       mpsr_GetMisc.cc
9 * Purpose:    Miscellanous routines which are needed by mpsr_Get
10 * Author:     Olaf Bachmann (10/95)
11 *
12 * Change History (most recent first):
13 *
14 ***************************************************************/
15
16#include "mod2.h"
17
18#ifdef HAVE_MPSR
19
20#include "mpsr_Get.h"
21#include "longalg.h"
22#include "tok.h"
23#include "maps.h"
24#include "lists.h"
25
26BOOLEAN currComplete = FALSE;
27
28static ring rDefault(short ch, char *name);
29static char* GenerateRingName();
30
31
32// use the varname so that they are compatible with the ones generated
33// by GetRingAnnots
34#define MPSR_DEFAULT_VARNAME "x(1)"
35
36// returns some default ring
37int mpsr_rDefault(short ch, char *name, ring &r)
38{
39  // check for currRing
40  if (currRing != NULL && rInternalChar(currRing) == ch)
41  // orig: currRing->ch==ch ???
42  {
43    int i, n = currRing->N;
44    char **names = currRing->names;
45
46    for (i=0; i<n; i++)
47    {
48      if (strcmp(names[i], name) == 0)
49      {
50        (currRing->ref)++;
51        r = currRing;
52        return i;
53      }
54    }
55  }
56  r = rDefault(ch, name);
57  return 0;
58}
59
60ring mpsr_rDefault(short ch)
61{
62  if (currRing != NULL && rInternalChar(currRing) == ch)
63  // orig: currRing->ch==ch ???
64  {
65    (currRing->ref)++;
66    return currRing;
67  }
68  else
69    return rDefault(ch, MPSR_DEFAULT_VARNAME);
70}
71
72static ring rDefault(short ch, char *name)
73{
74  ring r = (ring) Alloc0(sizeof(sip_sring));
75  r->ch = ch;
76  r->N = 1;
77  r->names = (char **) Alloc(sizeof(char *));
78  r->names[0] = mstrdup(name);
79
80  r->wvhdl = (short **)Alloc0(3 * sizeof(short *));
81  /*order: dp,C,0*/
82  r->order = (int *) Alloc(3 * sizeof(int *));
83  r->block0 = (int *)Alloc(3 * sizeof(int *));
84  r->block1 = (int *)Alloc(3 * sizeof(int *));
85  /* ringorder dp for the first block: var 1..3 */
86  r->order[0]  = ringorder_unspec;
87  r->block0[0] = 1;
88  r->block1[0] = 1;
89  /* ringorder C for the second block: no vars */
90  r->order[1]  = ringorder_C;
91  r->block0[1] = 0;
92  r->block1[1] = 0;
93  /* the last block: everything is 0 */
94  r->order[2]  = 0;
95  r->block0[2] = 0;
96  r->block1[2] = 0;
97  rComplete(r);
98  return r;
99}
100
101// returns TRUE, if r1 equals r2
102// FALSE, otherwise
103// Equality is determined componentwise
104BOOLEAN mpsr_RingEqual(ring r1, ring r2)
105{
106  int i, j;
107
108  if (r1 == r2) return 1;
109
110  if (r1 == NULL || r2 == NULL) return 0;
111
112  if ((rInternalChar(r1) != rInternalChar(r2))
113  // orig: r1->ch == r2->ch ???
114  || (r1->N != r2->N) || (r1->OrdSgn != r2->OrdSgn)
115      || (rPar(r1) != rPar(r2)))
116    return 0;
117
118  for (i=0; i<r1->N; i++)
119    if (strcmp(r1->names[i], r2->names[i])) return 0;
120
121  i=0;
122  while (r1->order[i] != 0)
123  {
124    if (r2->order[i] == 0) return 0;
125    if ((r1->order[i] != r2->order[i]) ||
126        (r1->block0[i] != r2->block0[i]) || (r2->block0[i] != r1->block0[i]))
127      return 0;
128    if (r1->wvhdl[i] != NULL)
129    {
130      if (r2->wvhdl[i] == NULL)
131        return 0;
132      for (j=0; j<r1->block1[i]-r1->block0[i]+1; j++)
133        if (r2->wvhdl[i][j] != r1->wvhdl[i][j])
134          return 0;
135    }
136    else if (r2->wvhdl[i] != NULL) return 0;
137    i++;
138  }
139
140  for (i=0; i<rPar(r1);i++)
141  {
142      if (strcmp(r1->parameter[i], r2->parameter[i])!=0)
143        return 0;
144  }
145
146  if (r1->minpoly != NULL)
147  {
148    if (r2->minpoly == NULL) return 0;
149    mpsr_SetCurrRing(r1, FALSE);
150    if (! naEqual(r1->minpoly, r2->minpoly)) return 0;
151  }
152  else if (r2->minpoly != NULL) return 0;
153
154  if (r1->qideal != NULL)
155  {
156    ideal id1 = r1->qideal, id2 = r2->qideal;
157    int i, n;
158    poly *m1, *m2;
159
160    if (id2 == NULL) return 0;
161    if ((n = IDELEMS(id1)) != IDELEMS(id2)) return 0;
162
163    mpsr_SetCurrRing(r1, FALSE);
164    m1 = id1->m;
165    m2 = id2->m;
166    for (i=0; i<n; i++)
167      if (! pEqualPolys(m1[i],m2[i])) return 0;
168  }
169  else if (r2->qideal != NULL) return 0;
170
171  return 1;
172}
173
174
175// returns TRUE, if r1 less or equals r2
176// FALSE, otherwise
177// Less or equal means that r1 is a strong subring of r2
178inline BOOLEAN RingLessEqual(ring r1, ring r2)
179{
180  int i, j;
181
182  if (r1 == r2) return 1;
183
184  if (r1 == NULL) return 1;
185
186  if (r2 == NULL) return 0;
187
188  if ((r1->N > r2->N) || (r1->OrdSgn != r2->OrdSgn) || (rPar(r1) > rPar(r2)))
189    return 0;
190
191  if (!rField_is_Q(r1) && rInternalChar(r1) != rInternalChar(r2)) return 0;
192  // orig: if (r1->ch != 0 && r1->ch != r2->ch) return 0;
193
194  for (i=0, j=0; j<r1->N && i<r2->N; i++)
195    if (strcmp(r1->names[j], r2->names[i]) == 0) j++;
196  if (j < r1->N) return 0;
197
198  // for ordering, suppose that they are only simple orderings
199  if (r1->order[2] != 0 || r2->order[2] != 0 ||
200      (r1->order[0] != r2->order[0] && r1->order[0] != ringorder_unspec) ||
201      r1->order[1] != r2->order[1])
202    return 0;
203
204  for (i=0; i<rPar(r1);i++)
205  {
206      if (strcmp(r1->parameter[i], r2->parameter[i])!=0)
207        return 0;
208  }
209  // r1->parameter == NULL && r2->parameter != NULL  is ok
210
211  if (r1->minpoly != NULL)
212  {
213    if (r2->minpoly == NULL) return 0;
214    mpsr_SetCurrRing(r1, FALSE);
215    if (! naEqual(r1->minpoly, r2->minpoly)) return 0;
216  }
217  return 1;
218}
219
220// returns MP_Success and lv2 appended to lv1, both over the same ring,
221// or MP_Failure
222mpsr_Status_t mpsr_MergeLeftv(mpsr_leftv mlv1, mpsr_leftv mlv2)
223{
224  ring r, r1 = mlv1->r, r2 = mlv2->r;
225  leftv lv;
226
227  if (mpsr_RingEqual(r1,r2))
228  {
229    if (r2 != NULL) rKill(r2);
230    r = r1;
231  }
232  else if (RingLessEqual(r1, r2))
233  {
234    r = r2;
235    if (r1 != NULL)
236    {
237      mpsr_MapLeftv(mlv1->lv, r1, r);
238      rKill(r1);
239    }
240  }
241  else if (RingLessEqual(r2, r1))
242  {
243    r = r1;
244    if (r2 != NULL)
245    {
246      mpsr_MapLeftv(mlv2->lv, r2, r);
247      rKill(r2);
248    }
249  }
250  else if (rSum(r1, r2, r) >= 0)
251  {
252    mpsr_MapLeftv(mlv1->lv, r1, r);
253    mpsr_MapLeftv(mlv2->lv, r2, r);
254    rKill(r1);
255    rKill(r2);
256  }
257  else return mpsr_Failure;
258
259  lv = mlv1->lv;
260  while (lv->next != NULL) lv = lv->next;
261
262  lv->next = mlv2->lv;
263  mlv1->r = r;
264
265#ifdef RDEBUG
266  if (r!= NULL) rTest(r);
267#endif
268  // this is an optimization for the mpsr_rDefault routines
269  currRing = r;
270  return mpsr_Success;
271}
272
273void mpsr_MapLeftv(leftv l, ring from_ring, ring to_ring)
274{
275  int i, n;
276
277  while (l != NULL)
278  {
279    short typ = l->Typ();
280
281    switch(typ)
282    {
283        case POLY_CMD:
284        case VECTOR_CMD:
285        {
286          poly p = (poly) l->data;
287          mpsr_SetCurrRing(to_ring, TRUE);
288          l->data = (void *) maIMap(from_ring, (poly) l->data);
289          ppDelete(&p, from_ring);
290          break;
291        }
292
293        case MODUL_CMD:
294        case IDEAL_CMD:
295        case MATRIX_CMD:
296        case MAP_CMD:
297        {
298          ideal id = (ideal) l->Data();
299          n = IDELEMS(id);
300          poly *m = id->m, *m1 = (poly *) Alloc(n*sizeof(poly));
301          mpsr_SetCurrRing(to_ring, TRUE);
302          for (i=0; i<n; i++)
303          {
304            m1[i] = m[i];
305            m[i] = maIMap(from_ring, m[i]);
306          }
307          mpsr_SetCurrRing(from_ring, FALSE);
308          for (i=0; i<n; i++) pDelete(&(m1[i]));
309          Free(m1, n*sizeof(poly));
310          break;
311        }
312
313        case LIST_CMD:
314        {
315          lists ll = (lists) l->Data();
316          n = ll->nr + 1;
317          for (i=0; i<n; i++) mpsr_MapLeftv(&(ll->m[i]), from_ring, to_ring);
318          break;
319        }
320
321        case COMMAND:
322        {
323          command cmd = (command) l->Data();
324          if (cmd->op == PROC_CMD && cmd->argc == 2)
325            mpsr_MapLeftv(&(cmd->arg2), from_ring, to_ring);
326          else if (cmd->argc > 0)
327          {
328            mpsr_MapLeftv(&(cmd->arg1), from_ring, to_ring);
329            if (cmd->argc > 1)
330            {
331              mpsr_MapLeftv(&(cmd->arg2), from_ring, to_ring);
332              if (cmd->argc > 2)
333                mpsr_MapLeftv(&(cmd->arg3), from_ring, to_ring);
334            }
335          }
336          break;
337        }
338
339        case NUMBER_CMD:
340        {
341          number nn = (number) l->data;
342          mpsr_SetCurrRing(to_ring, TRUE);
343          nSetMap(rInternalChar(from_ring), from_ring->parameter,
344            rPar(from_ring), from_ring->minpoly);
345          l->data = (void *) nMap(nn);
346          mpsr_SetCurrRing(from_ring, FALSE);
347          nDelete(&nn);
348        }
349    }
350    l = l->next;
351  }
352}
353
354
355// searches for a ring handle which has a ring which is equal to r
356// if one is found, then this one is set to the new global ring
357// otherwise, a ring name is generated, and a new idhdl is created
358void mpsr_SetCurrRingHdl(ring r)
359{
360  idhdl h = IDROOT, rh = NULL;
361
362  if (r == NULL)
363  {
364    if (currRingHdl != NULL && currRing != IDRING(currRingHdl))
365      mpsr_SetCurrRing(IDRING(currRingHdl), TRUE);
366    return;
367  }
368
369  rTest(r);
370  // try to find an idhdl which is an equal ring
371  while (h != NULL)
372  {
373    if ((IDTYP(h) == RING_CMD || IDTYP(h) == QRING_CMD) &&
374        (mpsr_RingEqual(IDRING(h), r)))
375    {
376      // found one
377      rh = h;
378      break;
379    }
380    h = h->next;
381  }
382
383  if (rh != NULL)
384  {
385    // found an idhdl to an equal ring
386    // we better reset currRing, so that rSetHdl does not choke (see
387    // sLastPrinted)
388    if (currRingHdl != NULL && IDRING(currRingHdl) != currRing)
389      mpsr_SetCurrRing(IDRING(currRingHdl), TRUE);
390
391    rSetHdl(rh, TRUE);
392
393    if (currRing != r)
394    {
395      mpsr_assume(r->ref <= 0);
396      rKill(r);
397    }
398  }
399  else
400  {
401    rh = mpsr_InitIdhdl((r->qideal == NULL ? (short) RING_CMD
402                         : (short) QRING_CMD),
403                        (void *) r, GenerateRingName());
404    // reset currRing for reasons explained above
405    if (currRingHdl != NULL) mpsr_SetCurrRing(IDRING(currRingHdl), TRUE);
406    rSetHdl(rh, TRUE);
407    rh->next = IDROOT;
408    IDROOT = rh;
409    r->ref = 0;
410  }
411}
412
413
414int gringcounter = 0;
415char grname[14];
416
417static char* GenerateRingName()
418{
419  sprintf(grname, "mpsr_r%d", gringcounter++);
420  return grname;
421}
422
423// searches through the Singular namespace for a matching name:
424// the first found is returned together witht the respective ring
425idhdl mpsr_FindIdhdl(char *name, ring &r)
426{
427#ifdef HAVE_NAMESPACES
428  idhdl h = (NSROOT(namespaceroot) != NULL ? namespaceroot->get(name, 0): (idhdl) NULL), h2;
429#else /* HAVE_NAMESPACES */
430  idhdl h = (idroot != NULL ? idroot->get(name, 0): (idhdl) NULL), h2;
431#endif /* HAVE_NAMESPACES */
432  r = NULL;
433
434  if (h != NULL)
435  {
436    r = NULL;
437    return h;
438  }
439
440  h = IDROOT;
441  while ( h != NULL)
442  {
443    if (IDTYP(h) == RING_CMD || IDTYP(h) == QRING_CMD)
444    {
445      h2 = (IDRING(h)->idroot!=NULL ?IDRING(h)->idroot->get(name, 0) :
446            (idhdl) NULL);
447
448      if (h2 != NULL)
449      {
450        r = IDRING(h);
451        r->ref++;
452        return h2;
453      }
454    }
455    h = h->next;
456  }
457  return NULL;
458}
459
460
461/***************************************************************
462 *
463 * Stuff which deals with External Data
464 *
465 ***************************************************************/
466
467void mpsr_DeleteExternalData(MPT_ExternalData_t edata)
468{
469  mpsr_leftv mlv = (mpsr_leftv) edata;
470
471  if (edata != NULL)
472  {
473    if (mlv->r != NULL) mpsr_SetCurrRing(mlv->r, FALSE);
474    if (mlv->lv != NULL)
475    {
476      mlv->lv->CleanUp();
477      Free(mlv->lv, sizeof(sleftv));
478    }
479    if (mlv->r != NULL) rKill(mlv->r);
480  }
481  Free(mlv, sizeof(mpsr_sleftv));
482}
483
484void mpsr_CopyExternalData(MPT_ExternalData_t *dest,
485                           MPT_ExternalData_t src)
486{
487  mpsr_leftv slv = (mpsr_leftv) src, dlv;
488
489  if (slv != NULL)
490  {
491    dlv = (mpsr_leftv) Alloc0(sizeof(mpsr_sleftv));
492    dlv->r = rCopy(slv->r);
493    dlv->lv = (leftv) Alloc0(sizeof(sleftv));
494    if (slv->lv != NULL) dlv->lv->Copy(slv->lv);
495    else dlv->lv = NULL;
496
497    *dest = (MPT_ExternalData_t) dlv;
498  }
499  else
500    *dest = NULL;
501}
502
503/***************************************************************
504 *
505 * mpsr initialization
506 *
507 ***************************************************************/
508
509#ifdef MDEBUG
510void * mpAllocBlock( size_t t)
511{
512  return mmDBAllocBlock(t,"mp",0);
513}
514void mpFreeBlock( void* a, size_t t)
515{
516  mmDBFreeBlock(a,t,"mp",0);
517}
518
519void * mpAlloc( size_t t)
520{
521  return mmDBAlloc(t,"mp",0);
522}
523void mpFree(void* a)
524{
525  mmDBFree(a,"mp",0);
526}
527#endif
528
529void mpsr_Init()
530{
531#ifndef EXTERNAL_MALLOC_H
532  // memory management functions of MP (and MPT)
533#ifndef MDEBUG
534  IMP_RawMemAllocFnc = mmAlloc;
535  IMP_RawMemFreeFnc = mmFree;
536  IMP_MemAllocFnc = mmAllocBlock;
537  IMP_MemFreeFnc = mmFreeBlock;
538#else
539  IMP_RawMemAllocFnc = mpAlloc;
540  IMP_RawMemFreeFnc = mpFree;
541  IMP_MemAllocFnc = mpAllocBlock;
542  IMP_MemFreeFnc = mpFreeBlock;
543#endif
544#endif
545
546  // Init of the MPT External Data functions
547  MPT_GetExternalData = mpsr_GetExternalData;
548  MPT_DeleteExternalData = mpsr_DeleteExternalData;
549
550#ifdef  PARI_BIGINT_TEST
551  init(4000000, 2);
552#endif
553}
554
555#ifdef MPSR_DEBUG
556// this is just a dummy function, where we can set a debugger breakpoint
557void mpsr_Break()
558{
559  Werror("mpsr_Error");
560}
561#endif
562
563#endif // HAVE_MPSR
Note: See TracBrowser for help on using the repository browser.