Changeset 3dbee61 in git


Ignore:
Timestamp:
Sep 29, 2011, 5:34:55 PM (12 years ago)
Author:
Hans Schoenemann <hannes@…>
Branches:
(u'spielwiese', '91fdef05f09f54b8d58d92a472e9c4a43aa4656f')
Children:
61eb53a68205301a8f2e44993354fbe213d648eb
Parents:
229530c0e416b111939f7e2475c01b2a091b0e55
Message:
merging SINGULAR-algnumbers: longalg rev 13634

git-svn-id: file:///usr/local/Singular/svn/trunk@14399 2c84dea3-7e68-4137-9b89-c4e89433aadc
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • Singular/iparith.cc

    r229530 r3dbee61  
    20522052        if ((r->minpoly != NULL) || (r->minideal != NULL))
    20532053          naSetChar(rInternalChar(r),r);
    2054         else ntSetChar(rInternalChar(r),r);
     2054        else naSetChar(rInternalChar(r),r);
    20552055        nSetChar(currRing);
    20562056        test=save_test;
  • Singular/ipassign.cc

    r229530 r3dbee61  
    191191    /* redefine function pointers due to switch from
    192192       transcendental to algebraic field extension */
    193     redefineFunctionPointers();
     193    //redefineFunctionPointers();
    194194    // and now, normalize all already defined objects in this ring
    195195    idhdl h=currRing->idroot;
  • Singular/ipshell.cc

    r229530 r3dbee61  
    23362336            R->algring->qideal->m[0]=NULL;
    23372337            idDelete(&(R->algring->qideal));
    2338             redefineFunctionPointers();
     2338            //redefineFunctionPointers();
    23392339          }
    23402340          else
     
    24092409            if ((orig_ring->minpoly != NULL) || (orig_ring->minideal != NULL))
    24102410              naSetChar(rInternalChar(orig_ring),orig_ring);
    2411             else ntSetChar(rInternalChar(orig_ring),orig_ring);
     2411            else naSetChar(rInternalChar(orig_ring),orig_ring);
    24122412            nSetChar(currRing);
    24132413            test=save_test;
  • kernel/longalg.cc

    r229530 r3dbee61  
    2626#endif
    2727#include <kernel/longalg.h>
    28 #include <kernel/longtrans.h>
     28
     29struct snaIdeal
     30{
     31  int anz;
     32  napoly *liste;
     33};
     34typedef struct snaIdeal * naIdeal;
     35
     36naIdeal naI=NULL;
     37
     38omBin snaIdeal_bin = omGetSpecBin(sizeof(snaIdeal));
     39
     40int naNumbOfPar;
     41napoly naMinimalPoly;
     42#define naParNames (currRing->parameter)
     43static int naIsChar0;
     44static ring naMapRing;
    2945
    3046#ifdef LDEBUG
     
    3551#endif
    3652
    37 naIdeal naI = NULL;
    38 napoly  naMinimalPoly;
    39 omBin   snaIdeal_bin = omGetSpecBin(sizeof(snaIdeal));
    40 number  (*naMap)(number from);
    41 //omBin lnumber_bin = omGetSpecBin(sizeof(slnumber));
    42 //omBin rnumber_bin = omGetSpecBin(sizeof(snumber));
    43 
    44 void redefineFunctionPointers()
    45 {
    46   n_Procs_s* n = currRing->cf;
    47   /* re-defining function pointers */
    48   n->cfDelete       = naDelete;
    49   n->nNormalize     = naNormalize;
    50   n->cfInit         = naInit;
    51   n->nPar           = naPar;
    52   n->nParDeg        = naParDeg;
    53   n->n_Int          = naInt;
    54   n->nAdd           = naAdd;
    55   n->nSub           = naSub;
    56   n->nMult          = naMult;
    57   n->nDiv           = naDiv;
    58   n->nExactDiv      = naDiv;
    59   n->nIntDiv        = naIntDiv;
    60   n->nNeg           = naNeg;
    61   n->nInvers        = naInvers;
    62   n->nCopy          = naCopy;
    63   n->cfCopy         = na_Copy;
    64   n->nGreater       = naGreater;
    65   n->nEqual         = naEqual;
    66   n->nIsZero        = naIsZero;
    67   n->nIsOne         = naIsOne;
    68   n->nIsMOne        = naIsMOne;
    69   n->nGreaterZero   = naGreaterZero;
    70   n->nGreater       = naGreater;
    71   n->cfWrite        = naWrite;
    72   n->nRead          = naRead;
    73   n->nPower         = naPower;
    74   n->nGcd           = naGcd;
    75   n->nLcm           = naLcm;
    76   n->cfSetMap       = naSetMap;
    77   n->nName          = naName;
    78   n->nSize          = naSize;
    79   n->cfGetDenom     = napGetDenom;
    80   n->cfGetNumerator = napGetNumerator;
    81 #ifdef LDEBUG
    82   n->nDBTest        = naDBTest;
    83 #endif
    84   /* re-defining global function pointers */
    85   nNormalize=naNormalize;
    86   nPar   = naPar;
    87   nParDeg= nParDeg;
    88   n_Int  = naInt;
    89   nAdd   = naAdd;
    90   nSub   = naSub;
    91   nMult  = naMult;
    92   nDiv   = naDiv;
    93   nExactDiv= naDiv;
    94   nIntDiv= naIntDiv;
    95   nNeg   = naNeg;
    96   nInvers= naInvers;
    97   nCopy  = naCopy;
    98   nGreater = naGreater;
    99   nEqual = naEqual;
    100   nIsZero = naIsZero;
    101   nIsOne = naIsOne;
    102   nIsMOne = naIsMOne;
    103   nGreaterZero = naGreaterZero;
    104   nGreater = naGreater;
    105   nRead = naRead;
    106   nPower = naPower;
    107   nGcd  = naGcd;
    108   nLcm  = naLcm;
    109   nName= naName;
    110   nSize  = naSize;
    111 }
     53number (*naMap)(number from);
     54/* procedure variables */
     55static numberfunc
     56                nacMult, nacSub, nacAdd, nacDiv, nacIntDiv;
     57static number   (*nacGcd)(number a, number b, const ring r);
     58static number   (*nacLcm)(number a, number b, const ring r);
     59static number   (*nacInit)(int i, const ring r);
     60static int      (*nacInt)(number &n, const ring r);
     61static void     (*nacDelete)(number *a, const ring r);
     62#undef n_Delete
     63#define n_Delete(A,R) nacDelete(A,R)
     64       void     (*nacNormalize)(number &a);
     65static number   (*nacNeg)(number a);
     66       number   (*nacCopy)(number a);
     67static number   (*nacInvers)(number a);
     68       BOOLEAN  (*nacIsZero)(number a);
     69static BOOLEAN  (*nacIsOne)(number a);
     70static BOOLEAN  (*nacIsMOne)(number a);
     71static BOOLEAN  (*nacGreaterZero)(number a);
     72static const char   * (*nacRead) (const char *s, number *a);
     73static napoly napRedp(napoly q);
     74static napoly napTailred(napoly q);
     75static BOOLEAN napDivPoly(napoly p, napoly q);
     76static int napExpi(int i, napoly a, napoly b);
     77       ring nacRing;
     78
     79void naCoefNormalize(number pp);
     80
     81#define napCopy(p)       p_Copy(p,nacRing)
    11282
    11383static number nadGcd( number a, number b, const ring r) { return nacInit(1,r); }
     
    11787void naSetChar(int i, ring r)
    11888{
    119   assume((r->minpoly  != NULL) ||
    120          (r->minideal != NULL)    );
    121  
    12289  if (naI!=NULL)
    12390  {
     
    140107  else
    141108    naMinimalPoly = NULL;
    142    
    143109  if (r->minideal!=NULL)
    144110  {
     
    154120  }
    155121
    156   ntNumbOfPar=rPar(r);
     122  naNumbOfPar=rPar(r);
    157123  if (i == 1)
    158     ntIsChar0 = 1;
     124  {
     125    naIsChar0 = 1;
     126  }
    159127  else if (i < 0)
    160128  {
    161     ntIsChar0 = 0;
     129    naIsChar0 = 0;
    162130    npSetChar(-i, r->algring); // to be changed HS
    163131  }
     
    177145  nacNeg         = nacRing->cf->nNeg;
    178146  nacIsZero      = nacRing->cf->nIsZero;
     147  nacRead        = nacRing->cf->nRead;
    179148  nacGreaterZero = nacRing->cf->nGreaterZero;
    180   nacGreater     = nacRing->cf->nGreater;
    181149  nacIsOne       = nacRing->cf->nIsOne;
     150  nacIsMOne      = nacRing->cf->nIsMOne;
    182151  nacGcd         = nacRing->cf->nGcd;
    183152  nacLcm         = nacRing->cf->nLcm;
     
    186155  nacIntDiv      = nacRing->cf->nIntDiv;
    187156  nacInvers      = nacRing->cf->nInvers;
    188 }
     157  nacDelete      = nacRing->cf->cfDelete;
     158}
     159
     160/*============= procedure for polynomials: napXXXX =======================*/
     161
     162
     163
     164#ifdef LDEBUG
     165static void napTest(napoly p)
     166{
     167  while (p != NULL)
     168  {
     169    if (naIsChar0)
     170      nlDBTest(pGetCoeff(p), "", 0);
     171    pIter(p);
     172  }
     173}
     174#else
     175#define napTest(p) ((void) 0)
     176#endif
     177
     178#define napSetCoeff(p,n) {n_Delete(&pGetCoeff(p),nacRing);pGetCoeff(p)=n;}
     179#define napComp(p,q)     p_LmCmp((poly)p,(poly)q, nacRing)
     180#define napMultT(A,E)    A=(napoly)p_Mult_mm((poly)A,(poly)E,nacRing)
     181#define napDeg(p)        (int)p_Totaldegree(p, nacRing)
     182
     183/*3
     184* creates  an napoly
     185*/
     186napoly napInitz(number z)
     187{
     188  napoly a = (napoly)p_Init(nacRing);
     189  pGetCoeff(a) = z;
     190  return a;
     191}
     192
     193/*3
     194* copy a napoly. poly
     195*/
     196static napoly napCopyNeg(napoly p)
     197{
     198  napoly r=napCopy(p);
     199  r=(napoly)p_Neg((poly)r, nacRing);
     200  return r;
     201}
     202
     203/*3
     204* returns ph * z
     205*/
     206static void napMultN(napoly p, number z)
     207{
     208  number t;
     209
     210  while (p!=NULL)
     211  {
     212    t = nacMult(pGetCoeff(p), z);
     213    nacNormalize(t);
     214    n_Delete(&pGetCoeff(p),nacRing);
     215    pGetCoeff(p) = t;
     216    pIter(p);
     217  }
     218}
     219
     220/*3
     221*  division with rest; f = g*q + r, returns r and destroy f
     222*/
     223napoly napRemainder(napoly f, const napoly  g)
     224{
     225  napoly a, h, qq;
     226
     227  qq = (napoly)p_Init(nacRing);
     228  pNext(qq) = NULL;
     229  p_Normalize(g, nacRing);
     230  p_Normalize(f, nacRing);
     231  a = f;
     232  do
     233  {
     234    napSetExp(qq,1, p_GetExp(a,1,nacRing) - p_GetExp(g,1,nacRing));
     235    napSetm(qq);
     236    pGetCoeff(qq) = nacDiv(pGetCoeff(a), pGetCoeff(g));
     237    pGetCoeff(qq) = nacNeg(pGetCoeff(qq));
     238    nacNormalize(pGetCoeff(qq));
     239    h = napCopy(g);
     240    napMultT(h, qq);
     241    p_Normalize(h,nacRing);
     242    n_Delete(&pGetCoeff(qq),nacRing);
     243    a = napAdd(a, h);
     244  }
     245  while ((a!=NULL) && (p_GetExp(a,1,nacRing) >= p_GetExp(g,1,nacRing)));
     246  omFreeBinAddr(qq);
     247  return a;
     248}
     249
     250/*3
     251*  division with rest; f = g*q + r,  destroy f
     252*/
     253static void napDivMod(napoly f, napoly  g, napoly *q, napoly *r)
     254{
     255  napoly a, h, b, qq;
     256
     257  qq = (napoly)p_Init(nacRing);
     258  pNext(qq) = b = NULL;
     259  p_Normalize(g,nacRing);
     260  p_Normalize(f,nacRing);
     261  a = f;
     262  do
     263  {
     264    napSetExp(qq,1, p_GetExp(a,1,nacRing) - p_GetExp(g,1,nacRing));
     265    p_Setm(qq,nacRing);
     266    pGetCoeff(qq) = nacDiv(pGetCoeff(a), pGetCoeff(g));
     267    nacNormalize(pGetCoeff(qq));
     268    b = napAdd(b, napCopy(qq));
     269    pGetCoeff(qq) = nacNeg(pGetCoeff(qq));
     270    h = napCopy(g);
     271    napMultT(h, qq);
     272    p_Normalize(h,nacRing);
     273    n_Delete(&pGetCoeff(qq),nacRing);
     274    a = napAdd(a, h);
     275  }
     276  while ((a!=NULL) && (p_GetExp(a,1,nacRing) >= p_GetExp(g,1,nacRing)));
     277  omFreeBinAddr(qq);
     278  *q = b;
     279  *r = a;
     280}
     281
     282/*3
     283*  returns z with z*x mod c = 1
     284*/
     285static napoly napInvers(napoly x, const napoly c)
     286{
     287  napoly y, r, qa, qn, q;
     288  number t, h;
     289
     290  if (p_GetExp(x,1,nacRing) >= p_GetExp(c,1,nacRing))
     291    x = napRemainder(x, c);
     292  if (x==NULL)
     293  {
     294    goto zero_divisor;
     295  }
     296  if (p_GetExp(x,1,nacRing)==0)
     297  {
     298    if (!nacIsOne(pGetCoeff(x)))
     299    {
     300      nacNormalize(pGetCoeff(x));
     301      t = nacInvers(pGetCoeff(x));
     302      nacNormalize(t);
     303      n_Delete(&pGetCoeff(x),nacRing);
     304      pGetCoeff(x) = t;
     305    }
     306    return x;
     307  }
     308  y = napCopy(c);
     309  napDivMod(y, x, &qa, &r);
     310  if (r==NULL)
     311  {
     312    goto zero_divisor;
     313  }
     314  if (p_GetExp(r,1,nacRing)==0)
     315  {
     316    nacNormalize(pGetCoeff(r));
     317    t = nacInvers(pGetCoeff(r));
     318    nacNormalize(t);
     319    t = nacNeg(t);
     320    napMultN(qa, t);
     321    n_Delete(&t,nacRing);
     322    p_Normalize(qa,nacRing);
     323    p_Delete(&x,nacRing);
     324    p_Delete(&r,nacRing);
     325    return qa;
     326  }
     327  y = x;
     328  x = r;
     329  napDivMod(y, x, &q, &r);
     330  if (r==NULL)
     331  {
     332    goto zero_divisor;
     333  }
     334  if (p_GetExp(r,1,nacRing)==0)
     335  {
     336    q = p_Mult_q(q, qa,nacRing);
     337    q = napAdd(q, p_ISet(1,nacRing));
     338    nacNormalize(pGetCoeff(r));
     339    t = nacInvers(pGetCoeff(r));
     340    napMultN(q, t);
     341    p_Normalize(q,nacRing);
     342    n_Delete(&t,nacRing);
     343    p_Delete(&x,nacRing);
     344    p_Delete(&r,nacRing);
     345    if (p_GetExp(q,1,nacRing) >= p_GetExp(c,1,nacRing))
     346      q = napRemainder(q, c);
     347    return q;
     348  }
     349  q = p_Mult_q(q, napCopy(qa),nacRing);
     350  q = napAdd(q, p_ISet(1,nacRing));
     351  qa = napNeg(qa);
     352  loop
     353  {
     354    y = x;
     355    x = r;
     356    napDivMod(y, x, &qn, &r);
     357    if (r==NULL)
     358    {
     359      break;
     360    }
     361    if (p_GetExp(r,1,nacRing)==0)
     362    {
     363      q = p_Mult_q(q, qn,nacRing);
     364      q = napNeg(q);
     365      q = napAdd(q, qa);
     366      nacNormalize(pGetCoeff(r));
     367      t = nacInvers(pGetCoeff(r));
     368      //nacNormalize(t);
     369      napMultN(q, t);
     370      p_Normalize(q,nacRing);
     371      n_Delete(&t,nacRing);
     372      p_Delete(&x,nacRing);
     373      p_Delete(&r,nacRing);
     374      if (p_GetExp(q,1,nacRing) >= p_GetExp(c,1,nacRing))
     375        q = napRemainder(q, c);
     376      return q;
     377    }
     378    y = q;
     379    q = p_Mult_q(napCopy(q), qn, nacRing);
     380    q = napNeg(q);
     381    q = napAdd(q, qa);
     382    qa = y;
     383  }
     384// zero divisor found:
     385zero_divisor:
     386  Werror("zero divisor found - your minpoly is not irreducible");
     387  return x;
     388}
     389
     390/*3
     391* the max degree of an napoly poly (used for test of "simple" et al.)
     392*/
     393static int  napMaxDeg(napoly p)
     394{
     395  int  d = 0;
     396  while(p!=NULL)
     397  {
     398    d=si_max(d,napDeg(p));
     399    pIter(p);
     400  }
     401  return d;
     402}
     403
     404/*3
     405* the max degree of an napoly poly (used for test of "simple" et al.)
     406*/
     407static int  napMaxDegLen(napoly p, int &l)
     408{
     409  int  d = 0;
     410  int ll=0;
     411  while(p!=NULL)
     412  {
     413    d=si_max(d,napDeg(p));
     414    pIter(p);
     415    ll++;
     416  }
     417  l=ll;
     418  return d;
     419}
     420
     421
     422/*3
     423*writes a polynomial number
     424*/
     425void napWrite(napoly p,const BOOLEAN has_denom, const ring r)
     426{
     427  ring nacring=r->algring;
     428  if (p==NULL)
     429    StringAppendS("0");
     430  else if (p_LmIsConstant(p,nacring))
     431  {
     432    BOOLEAN kl=FALSE;
     433    if (has_denom)
     434    {
     435      number den=nacring->cf->cfGetDenom(pGetCoeff(p),nacring );
     436      kl=!n_IsOne(den,nacring);
     437      n_Delete(&den, nacring);
     438    }
     439    if (kl) StringAppendS("(");
     440    //StringAppendS("-1");
     441    n_Write(pGetCoeff(p),nacring);
     442    if (kl) StringAppendS(")");
     443  }
     444  else
     445  {
     446    StringAppendS("(");
     447    loop
     448    {
     449      BOOLEAN wroteCoeff=FALSE;
     450      if ((p_LmIsConstant(p,nacring))
     451      || ((!n_IsOne(pGetCoeff(p),nacring))
     452        && (!n_IsMOne(pGetCoeff(p),nacring))))
     453      {
     454        n_Write(pGetCoeff(p),nacring);
     455        wroteCoeff=(r->ShortOut==0);
     456      }
     457      else if (n_IsMOne(pGetCoeff(p),nacring))
     458      {
     459        StringAppendS("-");
     460      }
     461      int  i;
     462      for (i = 0; i < r->P; i++)
     463      {
     464        int e=p_GetExp(p,i+1,nacring);
     465        if (e > 0)
     466        {
     467          if (wroteCoeff)
     468            StringAppendS("*");
     469          else
     470            wroteCoeff=(r->ShortOut==0);
     471          StringAppendS(r->parameter[i]);
     472          if (e > 1)
     473          {
     474            if (r->ShortOut == 0)
     475              StringAppendS("^");
     476            StringAppend("%d", e);
     477          }
     478        }
     479      } /*for*/
     480      pIter(p);
     481      if (p==NULL)
     482        break;
     483      if (n_GreaterZero(pGetCoeff(p),nacring))
     484        StringAppendS("+");
     485    }
     486    StringAppendS(")");
     487  }
     488}
     489
     490
     491static const char *napHandleMons(const char *s, int i, napoly ex)
     492{
     493  int  j;
     494  if (strncmp(s,naParNames[i],strlen(naParNames[i]))==0)
     495  {
     496    s+=strlen(naParNames[i]);
     497    if ((*s >= '0') && (*s <= '9'))
     498    {
     499      s = eati(s, &j);
     500      napAddExp(ex,i+1,j);
     501    }
     502    else
     503      napAddExp(ex,i+1,1);
     504  }
     505  return s;
     506}
     507static const char *napHandlePars(const char *s, int i, napoly ex)
     508{
     509  int  j;
     510  if (strcmp(s,naParNames[i])==0)
     511  {
     512    s+=strlen(naParNames[i]);
     513    napSetExp(ex,i+1,1);
     514  }
     515  return s;
     516}
     517
     518/*3  reads a monomial  */
     519static const char  *napRead(const char *s, napoly *b)
     520{
     521  napoly a;
     522  int  i;
     523  a = (napoly)p_Init(nacRing);
     524  if ((*s >= '0') && (*s <= '9'))
     525  {
     526    s = nacRead(s, &pGetCoeff(a));
     527    if (nacIsZero(pGetCoeff(a)))
     528    {
     529      p_LmDelete(&a,nacRing);
     530      *b = NULL;
     531      return s;
     532    }
     533  }
     534  else
     535    pGetCoeff(a) = nacInit(1,nacRing);
     536  i = 0;
     537  const char  *olds=s;
     538  loop
     539  {
     540    s = napHandlePars(s, i, a);
     541    if (olds == s)
     542      i++;
     543    else if (*s == '\0')
     544    {
     545      *b = a;
     546      return s;
     547    }
     548    if (i >= naNumbOfPar)
     549      break;
     550  }
     551  i=0;
     552  loop
     553  {
     554    olds = s;
     555    s = napHandleMons(s, i, a);
     556    if (olds == s)
     557      i++;
     558    else
     559      i = 0;
     560    if ((*s == '\0') || (i >= naNumbOfPar))
     561      break;
     562  }
     563  *b = a;
     564  return s;
     565}
     566
     567static int napExp(napoly a, napoly b)
     568{
     569  while (pNext(a)!=NULL) pIter(a);
     570  int m = p_GetExp(a,1,nacRing);
     571  if (m==0) return 0;
     572  while (pNext(b)!=NULL) pIter(b);
     573  int mm=p_GetExp(b,1,nacRing);
     574  if (m > mm) m = mm;
     575  return m;
     576}
     577
     578/*
     579* finds the smallest i-th exponent in a and b
     580* used to find it in a fraction
     581*/
     582static int napExpi(int i, napoly a, napoly b)
     583{
     584  if (a==NULL || b==NULL) return 0;
     585  int m = p_GetExp(a,i+1,nacRing);
     586  if (m==0) return 0;
     587  while (pNext(a) != NULL)
     588  {
     589    pIter(a);
     590    if (m > p_GetExp(a,i+1,nacRing))
     591    {
     592      m = p_GetExp(a,i+1,nacRing);
     593      if (m==0) return 0;
     594    }
     595  }
     596  do
     597  {
     598    if (m > p_GetExp(b,i+1,nacRing))
     599    {
     600      m = p_GetExp(b,i+1,nacRing);
     601      if (m==0) return 0;
     602    }
     603    pIter(b);
     604  }
     605  while (b != NULL);
     606  return m;
     607}
     608
     609static void napContent(napoly ph)
     610{
     611  number h,d;
     612  napoly p;
     613
     614  p = ph;
     615  if (nacIsOne(pGetCoeff(p)))
     616    return;
     617  h = nacCopy(pGetCoeff(p));
     618  pIter(p);
     619  do
     620  {
     621    d=nacGcd(pGetCoeff(p), h, nacRing);
     622    if(nacIsOne(d))
     623    {
     624      n_Delete(&h,nacRing);
     625      n_Delete(&d,nacRing);
     626      return;
     627    }
     628    n_Delete(&h,nacRing);
     629    h = d;
     630    pIter(p);
     631  }
     632  while (p!=NULL);
     633  h = nacInvers(d);
     634  n_Delete(&d,nacRing);
     635  p = ph;
     636  while (p!=NULL)
     637  {
     638    d = nacMult(pGetCoeff(p), h);
     639    n_Delete(&pGetCoeff(p),nacRing);
     640    pGetCoeff(p) = d;
     641    pIter(p);
     642  }
     643  n_Delete(&h,nacRing);
     644}
     645
     646static void napCleardenom(napoly ph)
     647{
     648  number d, h;
     649  napoly p;
     650
     651  if (!naIsChar0)
     652    return;
     653  p = ph;
     654  h = nacInit(1,nacRing);
     655  while (p!=NULL)
     656  {
     657    d = nacLcm(h, pGetCoeff(p), nacRing);
     658    n_Delete(&h,nacRing);
     659    h = d;
     660    pIter(p);
     661  }
     662  if(!nacIsOne(h))
     663  {
     664    p = ph;
     665    while (p!=NULL)
     666    {
     667      d=nacMult(h, pGetCoeff(p));
     668      n_Delete(&pGetCoeff(p),nacRing);
     669      nacNormalize(d);
     670      pGetCoeff(p) = d;
     671      pIter(p);
     672    }
     673    n_Delete(&h,nacRing);
     674  }
     675  napContent(ph);
     676}
     677
     678static napoly napGcd0(napoly a, napoly b)
     679{
     680  number x, y;
     681  if (!naIsChar0)
     682    return p_ISet(1,nacRing);
     683  x = nacCopy(pGetCoeff(a));
     684  if (nacIsOne(x))
     685    return napInitz(x);
     686  while (pNext(a)!=NULL)
     687  {
     688    pIter(a);
     689    y = nacGcd(x, pGetCoeff(a), nacRing);
     690    n_Delete(&x,nacRing);
     691    x = y;
     692    if (nacIsOne(x))
     693      return napInitz(x);
     694  }
     695  do
     696  {
     697    y = nacGcd(x, pGetCoeff(b), nacRing);
     698    n_Delete(&x,nacRing);
     699    x = y;
     700    if (nacIsOne(x))
     701      return napInitz(x);
     702    pIter(b);
     703  }
     704  while (b!=NULL);
     705  return napInitz(x);
     706}
     707
     708/*3
     709* result =gcd(a,b)
     710*/
     711static napoly napGcd(napoly a, napoly b)
     712{
     713  int i;
     714  napoly g, x, y, h;
     715  if ((a==NULL)
     716  || ((pNext(a)==NULL)&&(nacIsZero(pGetCoeff(a)))))
     717  {
     718    if ((b==NULL)
     719    || ((pNext(b)==NULL)&&(nacIsZero(pGetCoeff(b)))))
     720    {
     721      return p_ISet(1,nacRing);
     722    }
     723    return napCopy(b);
     724  }
     725  else
     726  if ((b==NULL)
     727  || ((pNext(b)==NULL)&&(nacIsZero(pGetCoeff(b)))))
     728  {
     729    return napCopy(a);
     730  }
     731  if (naMinimalPoly != NULL)
     732  {
     733    if (p_GetExp(a,1,nacRing) >= p_GetExp(b,1,nacRing))
     734    {
     735      x = a;
     736      y = b;
     737    }
     738    else
     739    {
     740      x = b;
     741      y = a;
     742    }
     743    if (!naIsChar0) g = p_ISet(1,nacRing);
     744    else            g = napGcd0(x, y);
     745    if (pNext(y)==NULL)
     746    {
     747      napSetExp(g,1, napExp(x, y));
     748      p_Setm(g,nacRing);
     749      return g;
     750    }
     751    x = napCopy(x);
     752    y = napCopy(y);
     753    loop
     754    {
     755      h = napRemainder(x, y);
     756      if (h==NULL)
     757      {
     758        napCleardenom(y);
     759        if (!nacIsOne(pGetCoeff(g)))
     760          napMultN(y, pGetCoeff(g));
     761        p_LmDelete(&g,nacRing);
     762        return y;
     763      }
     764      else if (pNext(h)==NULL)
     765        break;
     766      x = y;
     767      y = h;
     768    }
     769    p_Delete(&y,nacRing);
     770    p_LmDelete(&h,nacRing);
     771    napSetExp(g,1, napExp(a, b));
     772    p_Setm(g,nacRing);
     773    return g;
     774  }
     775  // Hmm ... this is a memory leak
     776  // x = (napoly)p_Init(nacRing);
     777  g=a;
     778  h=b;
     779  if (!naIsChar0) x = p_ISet(1,nacRing);
     780  else            x = napGcd0(g,h);
     781  for (i=(naNumbOfPar-1); i>=0; i--)
     782  {
     783    napSetExp(x,i+1, napExpi(i,a,b));
     784    p_Setm(x,nacRing);
     785  }
     786  return x;
     787}
     788
     789
     790number napLcm(napoly a)
     791{
     792  number h = nacInit(1,nacRing);
     793
     794  if (naIsChar0)
     795  {
     796    number d;
     797    napoly b = a;
     798
     799    while (b!=NULL)
     800    {
     801      d = nacLcm(h, pGetCoeff(b), nacRing);
     802      n_Delete(&h,nacRing);
     803      h = d;
     804      pIter(b);
     805    }
     806  }
     807  return h;
     808}
     809
     810
     811/*2
     812* meins  (for reduction in algebraic extension)
     813* checks if head of p divides head of q
     814* doesn't delete p and q
     815*/
     816static BOOLEAN napDivPoly (napoly p, napoly q)
     817{
     818  int j=1;   /* evtl. von naNumber.. -1 abwaerts zaehlen */
     819
     820  while (p_GetExp(p,j,nacRing) <= p_GetExp(q,j,nacRing))
     821  {
     822    j++;
     823    if (j > naNumbOfPar)
     824      return 1;
     825  }
     826  return 0;
     827}
     828
     829
     830/*2
     831* meins  (for reduction in algebraic extension)
     832* Normalform of poly with naI
     833* changes q and returns it
     834*/
     835napoly napRedp (napoly q)
     836{
     837  napoly h = (napoly)p_Init(nacRing);
     838  int i=0,j;
     839
     840  loop
     841  {
     842    if (napDivPoly (naI->liste[i], q))
     843    {
     844      /* h = lt(q)/lt(naI->liste[i])*/
     845      pGetCoeff(h) = nacCopy(pGetCoeff(q));
     846      for (j=naNumbOfPar; j>0; j--)
     847        napSetExp(h,j, p_GetExp(q,j,nacRing) - p_GetExp(naI->liste[i],j,nacRing));
     848      p_Setm(h,nacRing);
     849      h = p_Mult_q(h, napCopy(naI->liste[i]),nacRing);
     850      h = napNeg (h);
     851      q = napAdd (q, napCopy(h));
     852      p_Delete (&pNext(h),nacRing);
     853      if (q == NULL)
     854      {
     855        p_Delete(&h,nacRing);
     856        return q;
     857      }
     858      /* try to reduce further */
     859      i = 0;
     860    }
     861    else
     862    {
     863      i++;
     864      if (i >= naI->anz)
     865      {
     866        p_Delete(&h,nacRing);
     867        return q;
     868      }
     869    }
     870  }
     871}
     872
     873
     874/*2
     875* meins  (for reduction in algebraic extension)
     876* reduces the tail of Poly q
     877* needs q != NULL
     878* changes q and returns it
     879*/
     880napoly napTailred (napoly q)
     881{
     882  napoly h;
     883
     884  h = pNext(q);
     885  while (h != NULL)
     886  {
     887     h = napRedp (h);
     888     if (h == NULL)
     889        return q;
     890     pIter(h);
     891  }
     892  return q;
     893}
     894
    189895
    190896/*================ procedure for rational functions: naXXXX =================*/
     
    202908      poly z=p_Init(r->algring);
    203909      pSetCoeff0(z,c);
    204       lnumber l = ALLOC_LNUMBER();
     910      lnumber l = (lnumber)omAllocBin(rnumber_bin);
    205911      l->z = z;
    206912      l->s = 2;
     
    215921number  naPar(int i)
    216922{
    217   lnumber l = ALLOC_LNUMBER();
     923  lnumber l = (lnumber)omAllocBin(rnumber_bin);
    218924  l->s = 2;
    219925  l->z = p_ISet(1,nacRing);
     
    272978    p_Delete(&(l->z),r->algring);
    273979    p_Delete(&(l->n),r->algring);
    274     FREE_LNUMBER(l);
     980    omFreeBin((ADDRESS)l,  rnumber_bin);
    275981  }
    276982  *p = NULL;
     
    286992  lnumber erg;
    287993  lnumber src = (lnumber)p;
    288   erg = ALLOC_LNUMBER();
     994  erg = (lnumber)omAlloc0Bin(rnumber_bin);
    289995  erg->z = p_Copy(src->z, nacRing);
    290996  erg->n = p_Copy(src->n, nacRing);
     
    2971003  lnumber erg;
    2981004  lnumber src = (lnumber)p;
    299   erg = (lnumber)ALLOC_LNUMBER();
     1005  erg = (lnumber)omAlloc0Bin(rnumber_bin);
    3001006  erg->z = p_Copy(src->z,r->algring);
    3011007  erg->n = p_Copy(src->n,r->algring);
     
    3171023  lnumber b = (lnumber)lb;
    3181024  #ifdef LDEBUG
    319   omCheckAddrSize(a,sizeof(slnumber));
    320   omCheckAddrSize(b,sizeof(slnumber));
     1025  omCheckAddrSize(a,sizeof(snumber));
     1026  omCheckAddrSize(b,sizeof(snumber));
    3211027  #endif
    3221028  if (b->n!=NULL) x = pp_Mult_qq(a->z, b->n,nacRing);
     
    3291035    return (number)NULL;
    3301036  }
    331   lu = ALLOC_LNUMBER();
     1037  lu = (lnumber)omAllocBin(rnumber_bin);
    3321038  lu->z=res;
    3331039  if (a->n!=NULL)
     
    3851091
    3861092  #ifdef LDEBUG
    387   omCheckAddrSize(a,sizeof(slnumber));
    388   omCheckAddrSize(b,sizeof(slnumber));
     1093  omCheckAddrSize(a,sizeof(snumber));
     1094  omCheckAddrSize(b,sizeof(snumber));
    3891095  #endif
    3901096
     
    3991105    return (number)NULL;
    4001106  }
    401   lu = ALLOC_LNUMBER();
     1107  lu = (lnumber)omAllocBin(rnumber_bin);
    4021108  lu->z=res;
    4031109  if (a->n!=NULL)
     
    4311137number naMult(number la, number lb)
    4321138{
    433   if ((la==NULL) || (lb==NULL))   /* never occurs even when la or lb
    434                                      represents zero??? */
     1139  if ((la==NULL) || (lb==NULL))
    4351140    return NULL;
    4361141
     
    4411146
    4421147  #ifdef LDEBUG
    443   omCheckAddrSize(a,sizeof(slnumber));
    444   omCheckAddrSize(b,sizeof(slnumber));
     1148  omCheckAddrSize(a,sizeof(snumber));
     1149  omCheckAddrSize(b,sizeof(snumber));
    4451150  #endif
    4461151  naTest(la);
    4471152  naTest(lb);
    4481153
    449   lo = ALLOC_LNUMBER();
     1154  lo = (lnumber)omAllocBin(rnumber_bin);
    4501155  lo->z = pp_Mult_qq(a->z, b->z,nacRing);
    4511156
     
    4701175  if (naMinimalPoly!=NULL)
    4711176  {
    472     if ((lo->z != NULL) &&
    473         (p_GetExp(lo->z,1,nacRing) >= p_GetExp(naMinimalPoly,1,nacRing)))
     1177    if (p_GetExp(lo->z,1,nacRing) >= p_GetExp(naMinimalPoly,1,nacRing))
    4741178      lo->z = napRemainder(lo->z, naMinimalPoly);
    475     if ((x!=NULL) &&
    476         (p_GetExp(x,1,nacRing) >= p_GetExp(naMinimalPoly,1,nacRing)))
     1179    if ((x!=NULL) && (p_GetExp(x,1,nacRing) >= p_GetExp(naMinimalPoly,1,nacRing)))
    4771180      x = napRemainder(x, naMinimalPoly);
    4781181  }
     
    4951198  if(lo->z==NULL)
    4961199  {
    497     FREE_LNUMBER(lo);
     1200    omFreeBin((ADDRESS)lo, rnumber_bin);
    4981201    lo=NULL;
    4991202  }
     
    5271230  assume(a->z!=NULL && b->z!=NULL);
    5281231  assume(a->n==NULL && b->n==NULL);
    529   res = ALLOC_LNUMBER();
     1232  res = (lnumber)omAllocBin(rnumber_bin);
    5301233  res->z = napCopy(a->z);
    5311234  res->n = napCopy(b->z);
     
    5581261  }
    5591262  #ifdef LDEBUG
    560   omCheckAddrSize(a,sizeof(slnumber));
    561   omCheckAddrSize(b,sizeof(slnumber));
     1263  omCheckAddrSize(a,sizeof(snumber));
     1264  omCheckAddrSize(b,sizeof(snumber));
    5621265  #endif
    563   lo = ALLOC_LNUMBER();
     1266  lo = (lnumber)omAllocBin(rnumber_bin);
    5641267  if (b->n!=NULL)
    5651268    lo->z = pp_Mult_qq(a->z, b->n,nacRing);
     
    6351338  }
    6361339  #ifdef LDEBUG
    637   omCheckAddrSize(b,sizeof(slnumber));
     1340  omCheckAddrSize(b,sizeof(snumber));
    6381341  #endif
    639   lo = ALLOC0_LNUMBER();
     1342  lo = (lnumber)omAlloc0Bin(rnumber_bin);
    6401343  lo->s = b->s;
    6411344  if (b->n!=NULL)
     
    7441447
    7451448
    746 /* This method will only consider the numerators of a and b.
    747    Moreover it may return TRUE only if one or both numerators
    748    are zero or if their degrees are equal. Then TRUE is returned iff
    749    coeff(numerator(a)) > coeff(numerator(b));
    750    In all other cases, FALSE will be returned. */
    7511449BOOLEAN naGreater (number a, number b)
    7521450{
    753   int az = 0; int ad = 0;
    754   if (naIsZero(a)) az = 1;
    755   else ad = napDeg(((lnumber)a)->z);
    756   int bz = 0; int bd = 0;
    757   if (naIsZero(b)) bz = 1;
    758   else bd = napDeg(((lnumber)b)->z);
    759  
    760   if ((az == 1) && (bz == 1)) /* a = b = 0 */ return FALSE;
    761   if (az == 1) /* a = 0, b != 0 */
    762   {
    763     return (!nacGreaterZero(pGetCoeff(((lnumber)b)->z)));
    764   }
    765   if (bz == 1) /* a != 0, b = 0 */
    766   {
    767     return (nacGreaterZero(pGetCoeff(((lnumber)a)->z)));
    768   }
    769   if (ad == bd) 
    770     return nacGreater(pGetCoeff(((lnumber)a)->z),
    771                       pGetCoeff(((lnumber)b)->z));
    772   return FALSE;
     1451  if (naIsZero(a))
     1452    return FALSE;
     1453  if (naIsZero(b))
     1454    return TRUE; /* a!= 0)*/
     1455  return napDeg(((lnumber)a)->z)>napDeg(((lnumber)b)->z);
    7731456}
    7741457
     
    7861469    return s;
    7871470  }
    788   *p = (number)ALLOC0_LNUMBER();
     1471  *p = (number)omAlloc0Bin(rnumber_bin);
    7891472  a = (lnumber)*p;
    7901473  if ((naMinimalPoly!=NULL)
     
    8011484  if(a->z==NULL)
    8021485  {
    803     FREE_LNUMBER(a);
     1486    omFreeBin((ADDRESS)*p, rnumber_bin);
    8041487    *p=NULL;
    8051488  }
     
    8221505    return NULL;
    8231506  int i;
    824   char *s=(char *)omAlloc(4* ntNumbOfPar);
     1507  char *s=(char *)omAlloc(4* naNumbOfPar);
    8251508  char *t=(char *)omAlloc(8);
    8261509  s[0]='\0';
    827   for (i = 0; i <= ntNumbOfPar - 1; i++)
     1510  for (i = 0; i <= naNumbOfPar - 1; i++)
    8281511  {
    8291512    int e=p_GetExp(ph->z,i+1,nacRing);
     
    8321515      if (e >1)
    8331516      {
    834         sprintf(t,"%s%d",ntParNames[i],e);
     1517        sprintf(t,"%s%d",naParNames[i],e);
    8351518        strcat(s,t);
    8361519      }
    8371520      else
    8381521      {
    839         strcat(s,ntParNames[i]);
     1522        strcat(s,naParNames[i]);
    8401523      }
    8411524    }
     
    8811564  if (a==NULL) return FALSE;
    8821565#ifdef LDEBUG
    883   omCheckAddrSize(a,sizeof(slnumber));
     1566  omCheckAddrSize(a,sizeof(snumber));
    8841567  if (a->z==NULL)
    8851568  {
     
    9391622  if (a==NULL) return FALSE;
    9401623#ifdef LDEBUG
    941   omCheckAddrSize(a,sizeof(slnumber));
     1624  omCheckAddrSize(a,sizeof(snumber));
    9421625  if (a->z==NULL)
    9431626  {
     
    9481631  if (a->n==NULL)
    9491632  {
    950     if (p_LmIsConstant(a->z,nacRing)) return n_IsMOne(pGetCoeff(a->z),nacRing);
     1633    if (p_LmIsConstant(a->z,nacRing)) return nacIsMOne(pGetCoeff(a->z));
    9511634    /*else                   return FALSE;*/
    9521635  }
     
    9781661
    9791662  lnumber x, y;
    980   lnumber result = ALLOC0_LNUMBER();
     1663  lnumber result = (lnumber)omAlloc0Bin(rnumber_bin);
    9811664
    9821665  x = (lnumber)a;
    9831666  y = (lnumber)b;
    984   if ((ntNumbOfPar == 1) && (naMinimalPoly!=NULL))
     1667  if ((naNumbOfPar == 1) && (naMinimalPoly!=NULL))
    9851668  {
    9861669    if (pNext(x->z)!=NULL)
     
    10221705
    10231706/*2
    1024 * ntNumbOfPar = 1:
     1707* naNumbOfPar = 1:
    10251708* clears denominator         algebraic case;
    10261709* tries to simplify ratio    transcendental case;
     
    10281711* cancels monomials
    10291712* occuring in denominator
    1030 * and enumerator  ?          ntNumbOfPar != 1;
     1713* and enumerator  ?          naNumbOfPar != 1;
    10311714*
    10321715* #defines for Factory:
     
    11651848    p->z = x;
    11661849    p->n = y = NULL;
    1167     norm=ntIsChar0;
     1850    norm=naIsChar0;
    11681851  }
    11691852
     
    11751858    x = napRemainder(x, naMinimalPoly);
    11761859    p->z = x;
    1177     norm=ntIsChar0;
     1860    norm=naIsChar0;
    11781861  }
    11791862  /* normalize all coefficients in n and z (if in Q) */
     
    11891872  {
    11901873    int i;
    1191     for (i=ntNumbOfPar-1; i>=0; i--)
     1874    for (i=naNumbOfPar-1; i>=0; i--)
    11921875    {
    11931876      napoly xx=x;
     
    12281911  }
    12291912#ifndef FACTORY_GCD_TEST
    1230   if (ntNumbOfPar == 1) /* apply built-in gcd */
     1913  if (naNumbOfPar == 1) /* apply built-in gcd */
    12311914  {
    12321915    napoly x1,y1;
     
    12631946    y = p->n;
    12641947    /* collect all denoms from y and multiply x and y by it */
    1265     if (ntIsChar0)
     1948    if (naIsChar0)
    12661949    {
    12671950      number n=napLcm(y);
     
    13782061  lnumber a = (lnumber)la;
    13792062  lnumber b = (lnumber)lb;
    1380   result = ALLOC0_LNUMBER();
     2063  result = (lnumber)omAlloc0Bin(rnumber_bin);
     2064  //if (((naMinimalPoly==NULL) && (naI==NULL)) || !naIsChar0)
     2065  //{
     2066  //  result->z = p_ISet(1,nacRing);
     2067  //  return (number)result;
     2068  //}
     2069  //naNormalize(lb);
    13812070  naTest(la);
    13822071  naTest(lb);
     
    14642153{
    14652154  if (npIsZero(c)) return NULL;
    1466   lnumber l=ALLOC_LNUMBER();
     2155  lnumber l=(lnumber)omAllocBin(rnumber_bin);
    14672156  l->s=2;
    14682157  l->z=(napoly)p_Init(nacRing);
    14692158  int i=(int)((long)c);
    1470   if (i>((long)ntMapRing->ch>>2)) i-=(long)ntMapRing->ch;
     2159  if (i>((long)naMapRing->ch>>2)) i-=(long)naMapRing->ch;
    14712160  pGetCoeff(l->z)=nlInit(i, nacRing);
    14722161  l->n=NULL;
     
    14802169{
    14812170  if (nlIsZero(c)) return NULL;
    1482   lnumber l=ALLOC_LNUMBER();
     2171  lnumber l=(lnumber)omAllocBin(rnumber_bin);
    14832172  l->s=0;
    14842173  l->z=(napoly)p_Init(nacRing);
     
    14942183{
    14952184  if (npIsZero(c)) return NULL;
    1496   lnumber l=ALLOC_LNUMBER();
     2185  lnumber l=(lnumber)omAllocBin(rnumber_bin);
    14972186  l->s=2;
    14982187  l->z=(napoly)p_Init(nacRing);
     
    15092198  if (npIsZero(c)) return NULL;
    15102199  int i=(int)((long)c);
    1511   if (i>(long)ntMapRing->ch) i-=(long)ntMapRing->ch;
    1512   number n=npInit(i,ntMapRing);
     2200  if (i>(long)naMapRing->ch) i-=(long)naMapRing->ch;
     2201  number n=npInit(i,naMapRing);
    15132202  if (npIsZero(n)) return NULL;
    1514   lnumber l=ALLOC_LNUMBER();
     2203  lnumber l=(lnumber)omAllocBin(rnumber_bin);
    15152204  l->s=2;
    15162205  l->z=(napoly)p_Init(nacRing);
     
    15292218  if (npIsZero(n)) return NULL;
    15302219  npTest(n);
    1531   lnumber l=ALLOC_LNUMBER();
     2220  lnumber l=(lnumber)omAllocBin(rnumber_bin);
    15322221  l->s=2;
    15332222  l->z=(napoly)p_Init(nacRing);
     
    15372226}
    15382227
     2228static number (*nacMap)(number);
     2229static int naParsToCopy;
     2230static napoly napMap(napoly p)
     2231{
     2232  napoly w, a;
     2233
     2234  if (p==NULL) return NULL;
     2235  a = w = (napoly)p_Init(nacRing);
     2236  int i;
     2237  for(i=1;i<=naParsToCopy;i++)
     2238    napSetExp(a,i,napGetExpFrom(p,i,naMapRing));
     2239  p_Setm(a,nacRing);
     2240  pGetCoeff(w) = nacMap(pGetCoeff(p));
     2241  loop
     2242  {
     2243    pIter(p);
     2244    if (p==NULL) break;
     2245    pNext(a) = (napoly)p_Init(nacRing);
     2246    pIter(a);
     2247    for(i=1;i<=naParsToCopy;i++)
     2248      napSetExp(a,i,napGetExpFrom(p,i,naMapRing));
     2249    p_Setm(a,nacRing);
     2250    pGetCoeff(a) = nacMap(pGetCoeff(p));
     2251  }
     2252  pNext(a) = NULL;
     2253  return w;
     2254}
     2255
     2256static napoly napPerm(napoly p,const int *par_perm,const ring src_ring,const nMapFunc nMap)
     2257{
     2258  napoly w, a;
     2259
     2260  if (p==NULL) return NULL;
     2261  w = (napoly)p_Init(nacRing);
     2262  int i;
     2263  BOOLEAN not_null=TRUE;
     2264  loop
     2265  {
     2266    for(i=1;i<=rPar(src_ring);i++)
     2267    {
     2268      int e;
     2269      if (par_perm!=NULL) e=par_perm[i-1];
     2270      else                e=-i;
     2271      int ee=napGetExpFrom(p,i,src_ring);
     2272      if (e<0)
     2273        napSetExp(w,-e,ee);
     2274      else if (ee>0)
     2275        not_null=FALSE;
     2276    }
     2277    pGetCoeff(w) = nMap(pGetCoeff(p));
     2278    p_Setm(w,nacRing);
     2279    pIter(p);
     2280    if (!not_null)
     2281    {
     2282      if (p==NULL)
     2283      {
     2284        p_Delete(&w,nacRing);
     2285        return NULL;
     2286      }
     2287      /* else continue*/
     2288      n_Delete(&(pGetCoeff(w)),nacRing);
     2289    }
     2290    else
     2291    {
     2292      if (p==NULL) return w;
     2293      else
     2294      {
     2295        pNext(w)=napPerm(p,par_perm,src_ring,nMap);
     2296        return w;
     2297      }
     2298    }
     2299  }
     2300}
     2301
    15392302/*2
    15402303* map _(a) -> _(b)
     
    15432306{
    15442307  if (c==NULL) return NULL;
    1545   lnumber erg= ALLOC_LNUMBER();
     2308  lnumber erg= (lnumber)omAlloc0Bin(rnumber_bin);
    15462309  lnumber src =(lnumber)c;
    15472310  erg->s=src->s;
     
    15732336nMapFunc naSetMap(const ring src, const ring dst)
    15742337{
    1575   ntMapRing=src;
     2338  naMapRing=src;
    15762339  if (rField_is_Q_a(dst)) /* -> Q(a) */
    15772340  {
     
    15872350    {
    15882351      int i;
    1589       ntParsToCopy=0;
     2352      naParsToCopy=0;
    15902353      for(i=0;i<rPar(src);i++)
    15912354      {
     
    15932356        ||(strcmp(src->parameter[i],dst->parameter[i])!=0))
    15942357           return NULL;
    1595         ntParsToCopy++;
     2358        naParsToCopy++;
    15962359      }
    15972360      nacMap=nacCopy;
    1598       if ((ntParsToCopy==rPar(dst))&&(ntParsToCopy==rPar(src)))
     2361      if ((naParsToCopy==rPar(dst))&&(naParsToCopy==rPar(src)))
    15992362        return naCopy;    /* Q(a) -> Q(a) */
    16002363      return naMapQaQb;   /* Q(a..) -> Q(a..) */
     
    16302393      }
    16312394      int i;
    1632       ntParsToCopy=0;
     2395      naParsToCopy=0;
    16332396      for(i=0;i<rPar(src);i++)
    16342397      {
     
    16362399        ||(strcmp(src->parameter[i],dst->parameter[i])!=0))
    16372400           return NULL;
    1638         ntParsToCopy++;
    1639       }
    1640       if ((ntParsToCopy==rPar(dst))&&(ntParsToCopy==rPar(src))
     2401        naParsToCopy++;
     2402      }
     2403      if ((naParsToCopy==rPar(dst))&&(naParsToCopy==rPar(src))
    16412404      && (nacMap==nacCopy))
    16422405        return naCopy;    /* Z/p(a) -> Z/p(a) */
     
    16452408  }
    16462409  return NULL;      /* default */
     2410}
     2411
     2412/*2
     2413* convert a napoly number into a poly
     2414*/
     2415poly naPermNumber(number z, int * par_perm, int P, ring oldRing)
     2416{
     2417  if (z==NULL) return NULL;
     2418  poly res=NULL;
     2419  poly p;
     2420  napoly za=((lnumber)z)->z;
     2421  napoly zb=((lnumber)z)->n;
     2422  nMapFunc nMap=naSetMap(oldRing,currRing);
     2423  if (currRing->parameter!=NULL)
     2424    nMap=currRing->algring->cf->cfSetMap(oldRing->algring, nacRing);
     2425  else
     2426    nMap=currRing->cf->cfSetMap(oldRing->algring, currRing);
     2427  if (nMap==NULL) return NULL; /* emergency exit only */
     2428  do
     2429  {
     2430    p = pInit();
     2431    pNext(p)=NULL;
     2432    nNew(&pGetCoeff(p));
     2433    int i;
     2434    for(i=pVariables;i;i--)
     2435       pSetExp(p,i, 0);
     2436    if (rRing_has_Comp(currRing)) pSetComp(p, 0);
     2437    napoly pa=NULL;
     2438    lnumber pan;
     2439    if (currRing->parameter!=NULL)
     2440    {
     2441      assume(oldRing->algring!=NULL);
     2442      pGetCoeff(p)=(number)omAlloc0Bin(rnumber_bin);
     2443      pan=(lnumber)pGetCoeff(p);
     2444      pan->s=2;
     2445      pan->z=napInitz(nMap(pGetCoeff(za)));
     2446      pa=pan->z;
     2447    }
     2448    else
     2449    {
     2450      pGetCoeff(p)=nMap(pGetCoeff(za));
     2451    }
     2452    for(i=0;i<P;i++)
     2453    {
     2454      if(napGetExpFrom(za,i+1,oldRing)!=0)
     2455      {
     2456        if(par_perm==NULL)
     2457        {
     2458          if ((rPar(currRing)>=i) && (pa!=NULL))
     2459          {
     2460            napSetExp(pa,i+1,napGetExpFrom(za,i+1,oldRing));
     2461            p_Setm(pa,nacRing);
     2462          }
     2463          else
     2464          {
     2465            pDelete(&p);
     2466            break;
     2467          }
     2468        }
     2469        else if(par_perm[i]>0)
     2470          pSetExp(p,par_perm[i],napGetExpFrom(za,i+1,oldRing));
     2471        else if((par_perm[i]<0)&&(pa!=NULL))
     2472        {
     2473          napSetExp(pa,-par_perm[i], napGetExpFrom(za,i+1,oldRing));
     2474          p_Setm(pa,nacRing);
     2475        }
     2476        else
     2477        {
     2478          pDelete(&p);
     2479          break;
     2480        }
     2481      }
     2482    }
     2483    if (p!=NULL)
     2484    {
     2485      pSetm(p);
     2486      if (zb!=NULL)
     2487      {
     2488        if  (currRing->P>0)
     2489        {
     2490          pan->n=napPerm(zb,par_perm,oldRing,nMap);
     2491          if(pan->n==NULL) /* error in mapping or mapping to variable */
     2492            pDelete(&p);
     2493        }
     2494        else
     2495          pDelete(&p);
     2496      }
     2497      pTest(p);
     2498      res=pAdd(res,p);
     2499    }
     2500    pIter(za);
     2501  }
     2502  while (za!=NULL);
     2503  pTest(res);
     2504  return res;
     2505}
     2506
     2507number   naGetDenom(number &n, const ring r)
     2508{
     2509  lnumber x=(lnumber)n;
     2510  if (x->n!=NULL)
     2511  {
     2512    lnumber rr=(lnumber)omAlloc0Bin(rnumber_bin);
     2513    rr->z=p_Copy(x->n,r->algring);
     2514    rr->s = 2;
     2515    return (number)rr;
     2516  }
     2517  return n_Init(1,r);
     2518}
     2519
     2520number   naGetNumerator(number &n, const ring r)
     2521{
     2522  lnumber x=(lnumber)n;
     2523  lnumber rr=(lnumber)omAlloc0Bin(rnumber_bin);
     2524  rr->z=p_Copy(x->z,r->algring);
     2525  rr->s = 2;
     2526  return (number)rr;
    16472527}
    16482528
     
    16542534    return TRUE;
    16552535  #ifdef LDEBUG
    1656   omCheckAddrSize(a, sizeof(slnumber));
     2536  omCheckAddrSize(a, sizeof(snumber));
    16572537  #endif
    16582538  napoly p = x->z;
     
    16642544  while(p!=NULL)
    16652545  {
    1666     if (( ntIsChar0 && nlIsZero(pGetCoeff(p)))
    1667     || ((!ntIsChar0) && npIsZero(pGetCoeff(p))))
     2546    if ((naIsChar0 && nlIsZero(pGetCoeff(p)))
     2547    || ((!naIsChar0) && npIsZero(pGetCoeff(p))))
    16682548    {
    16692549      Print("coeff 0 in %s:%d\n",f,l);
     
    16772557      return FALSE;
    16782558    }
    1679     //if (ntIsChar0 && (((int)p->ko &3) == 0) && (p->ko->s==0) && (x->s==2))
     2559    //if (naIsChar0 && (((int)p->ko &3) == 0) && (p->ko->s==0) && (x->s==2))
    16802560    //{
    16812561    //  Print("normalized with non-normal coeffs in %s:%d\n",f,l);
    16822562    //  return FALSE;
    16832563    //}
    1684     if (ntIsChar0 && !(nlDBTest(pGetCoeff(p),f,l)))
     2564    if (naIsChar0 && !(nlDBTest(pGetCoeff(p),f,l)))
    16852565      return FALSE;
    16862566    pIter(p);
     
    16892569  while(p!=NULL)
    16902570  {
    1691     if (ntIsChar0 && !(nlDBTest(pGetCoeff(p),f,l)))
     2571    if (naIsChar0 && !(nlDBTest(pGetCoeff(p),f,l)))
    16922572      return FALSE;
    16932573    pIter(p);
     
    16962576}
    16972577#endif
     2578
     2579/*2
     2580* convert a napoly number into a poly
     2581*/
     2582poly napPermNumber(number z, int * par_perm, int P, ring oldRing)
     2583{
     2584  if (z==NULL) return NULL;
     2585  poly res=NULL;
     2586  poly p;
     2587  napoly za=((lnumber)z)->z;
     2588  napoly zb=((lnumber)z)->n;
     2589  nMapFunc nMap=naSetMap(oldRing,currRing); /* todo: check naSetMap
     2590                                                     vs. ntSetMap */
     2591  if (currRing->parameter!=NULL)
     2592    nMap=currRing->algring->cf->cfSetMap(oldRing->algring, nacRing);
     2593  else
     2594    nMap=currRing->cf->cfSetMap(oldRing->algring, currRing);
     2595  if (nMap==NULL) return NULL; /* emergency exit only */
     2596  while(za!=NULL)
     2597  {
     2598    p = pInit();
     2599    pNext(p)=NULL;
     2600    //nNew(&pGetCoeff(p));
     2601    int i;
     2602    //for(i=pVariables;i;i--) pSetExp(p,i, 0); // done by pInit
     2603    //if (rRing_has_Comp(currRing)) pSetComp(p, 0); // done by pInit
     2604    napoly pa=NULL;
     2605    lnumber pan;
     2606    if (currRing->parameter!=NULL)
     2607    {
     2608      assume(oldRing->algring!=NULL);
     2609      pGetCoeff(p)=(number)ALLOC0_LNUMBER();
     2610      pan=(lnumber)pGetCoeff(p);
     2611      pan->s=2;
     2612      pan->z=napInitz(nMap(pGetCoeff(za)));
     2613      pa=pan->z;
     2614    }
     2615    else
     2616    {
     2617      pGetCoeff(p)=nMap(pGetCoeff(za));
     2618    }
     2619    for(i=0;i<P;i++)
     2620    {
     2621      if(napGetExpFrom(za,i+1,oldRing)!=0)
     2622      {
     2623        if(par_perm==NULL)
     2624        {
     2625          if ((rPar(currRing)>=i) && (pa!=NULL))
     2626          {
     2627            napSetExp(pa,i+1,napGetExpFrom(za,i+1,oldRing));
     2628            p_Setm(pa,nacRing);
     2629          }
     2630          else
     2631          {
     2632            pDelete(&p);
     2633            break;
     2634          }
     2635        }
     2636        else if(par_perm[i]>0)
     2637          pSetExp(p,par_perm[i],napGetExpFrom(za,i+1,oldRing));
     2638        else if((par_perm[i]<0)&&(pa!=NULL))
     2639        {
     2640          napSetExp(pa,-par_perm[i], napGetExpFrom(za,i+1,oldRing));
     2641          p_Setm(pa,nacRing);
     2642        }
     2643        else
     2644        {
     2645          pDelete(&p);
     2646          break;
     2647        }
     2648      }
     2649    }
     2650    if (p!=NULL)
     2651    {
     2652      pSetm(p);
     2653      if (zb!=NULL)
     2654      {
     2655        if  (currRing->P>0)
     2656        {
     2657          pan->n=napPerm(zb,par_perm,oldRing,nMap);
     2658          if(pan->n==NULL) /* error in mapping or mapping to variable */
     2659            pDelete(&p);
     2660        }
     2661        else
     2662          pDelete(&p);
     2663      }
     2664      nNormalize(pGetCoeff(p));
     2665      if (nIsZero(pGetCoeff(p)))
     2666        pDelete(&p);
     2667      else
     2668      {
     2669        pTest(p);
     2670        res=pAdd(res,p);
     2671      }
     2672    }
     2673    pIter(za);
     2674  }
     2675  pTest(res);
     2676  return res;
     2677}
     2678number   napGetDenom(number &n, const ring r)
     2679{
     2680  lnumber x=(lnumber)n;
     2681  if (x->n!=NULL)
     2682  {
     2683    lnumber rr=ALLOC0_LNUMBER();
     2684    rr->z=p_Copy(x->n,r->algring);
     2685    rr->s = 2;
     2686    return (number)rr;
     2687  }
     2688  return n_Init(1,r);
     2689}
     2690
     2691number   napGetNumerator(number &n, const ring r)
     2692{
     2693  lnumber x=(lnumber)n;
     2694  lnumber rr=ALLOC0_LNUMBER();
     2695  rr->z=p_Copy(x->z,r->algring);
     2696  rr->s = 2;
     2697  return (number)rr;
     2698}
  • kernel/longalg.h

    r229530 r3dbee61  
    1111#include <kernel/longrat.h>
    1212#include <kernel/polys-impl.h>
    13 #include <kernel/longtrans.h>
    1413
    15 /*
    16    IMPORTANT INFORMATION:
    17    Instantiation of an algebraic field extension in SINGULAR
    18    works by first creating a transcendental field extension
    19    and then providing a minimal polynomial / minimal ideal.
    20    Consequently, first the code for transcendental field
    21    extensions will be activated; see longtrans.*.
    22    When providing a minimal polynomial / minimal ideal, all
    23    function pointers will be re-defined, using the below method
    24    redefineFunctionPointers(). After that, the code for algebraic
    25    field extensions is active; see longalg.*.
    26    
    27    NOTE:
    28    Elements of algebraic and transcendental field extensions
    29    are polynomials or quotients of two polynomials, respectively.
    30    All these polynomials, as well as the minimal polynomial (in
    31    the case of an algebraic extension) live in the globally
    32    accessible ring 'nacRing', defined in longtrans.*.
    33    
    34    METHOD NAMING CONVENTIONS
    35    (not true for types / structs / global variables):
    36    nap*   macros and methods operating on polynomials living in
    37           nacRing (defined in longtrans.*),
    38    na*    (but not nap*) methods and macros for algebraic field
    39           extensions (defined in longalg.*),
    40    nt*    methods and macros for transcendental field extensions,
    41           (defined in longtrans.*)
    42    nac*   function pointers for computing with the coefficients of
    43           polynomials living in nacRing (defined in longtrans.*)
    44 */
    45 
    46 struct snaIdeal
    47 {
    48   int anz;
    49   napoly *liste;
    50 };
    51 typedef struct snaIdeal * naIdeal;
    52 extern omBin snaIdeal_bin;
    53 extern naIdeal naI;
    54 extern napoly naMinimalPoly;
    55 
    56 /* for re-defining function pointers when switching from
    57    transcendental to algebraic extension (by providing
    58    a minpoly) */
    59 void          redefineFunctionPointers();
    60 
    61 /* specific methods / macros for algebraic field extensions */
    62 void          naSetChar(int p, ring r);
    63 void          naDelete (number *p, const ring r);
    64 number        naInit(int i, const ring r);                /* z := i      */
    65 number        naPar(int i);                               /* z := par(i) */
    66 int           naParDeg(number n);                         /* i := deg(n) */
    67 int           naSize(number n);                           /* size desc.  */
    68 int           naInt(number &n, const ring r);
    69 BOOLEAN       naIsZero(number za);                        /* za = 0 ?    */
    70 BOOLEAN       naIsOne(number  za);                        /* za = 1 ?    */
    71 BOOLEAN       naIsMOne(number  za);                       /* za = -1 ?   */
    72 BOOLEAN       naEqual(number a, number b);                /* a = b ?     */
    73 BOOLEAN       naGreater(number a, number b);              /* dummy       */
    74 number        naNeg(number za);                           /* za := - za  */
    75 number        naInvers(number a);
    76 void          naPower(number x, int exp, number *lo);
    77 BOOLEAN       naGreaterZero(number a);
    78 number        naCopy(number p);                           /* erg:= p     */
    79 number        na_Copy(number p, const ring r);            /* erg:= p     */
    80 number        naAdd(number la, number li);                /* lu := la+li */
    81 number        naMult(number la, number li);               /* lo := la*li */
    82 number        naDiv(number la, number li);                /* lo := la/li */
    83 number        naIntDiv(number la, number li);             /* lo := la/li */
    84 //number        naIntMod(number la, number li);           /* lo := la/li */
    85 number        naSub(number la, number li);                /* lu := la-li */
    86 void          naNormalize(number &p);
    87 number        naGcd(number a, number b, const ring r);
    88 number        naLcm(number a, number b, const ring r);
    89 const char *  naRead(const char * s, number * p);
    90 void          naWrite(number &p, const ring r);
    91 char *        naName(number n);
    92 nMapFunc      naSetMap(const ring src, const ring dst);
    93 number        naMap0P(number c);
    94 number        naMap00(number c);
    95 #ifdef LDEBUG
    96 BOOLEAN       naDBTest(number a, const char *f,const int l);
    97 #endif
    98 void          naSetIdeal(ideal I);
    99 void          naCoefNormalize(number pp);
    100 extern number (*naMap)(number from);
    101 
    102 //extern omBin lnumber_bin;
    103 //#define ALLOC_LNUMBER() (lnumber)omAllocBin(lnumber_bin)
    104 //#define ALLOC0_LNUMBER() (lnumber)omAlloc0Bin(lnumber_bin)
    105 //#define FREE_LNUMBER(x) omFreeBin((ADDRESS)x, lnumber_bin)
    106 extern omBin rnumber_bin;
    10714#define ALLOC_LNUMBER() (lnumber)omAllocBin(rnumber_bin)
    10815#define ALLOC0_LNUMBER() (lnumber)omAlloc0Bin(rnumber_bin)
    10916#define FREE_LNUMBER(x) omFreeBin((ADDRESS)x, rnumber_bin)
     17
     18
     19typedef polyrec *   napoly;
     20
     21struct slnumber;
     22typedef struct slnumber * lnumber;
     23
     24struct slnumber
     25{
     26  napoly z;
     27  napoly n;
     28  BOOLEAN s;
     29};
     30
     31extern napoly naMinimalPoly;
     32extern ring nacRing;
     33
     34
     35void naSetChar(int p, ring r);
     36void    naDelete (number *p, const ring r);
     37number  naInit(int i, const ring r);                /* z := i */
     38number  naPar(int i);                               /* z := par(i) */
     39int     naParDeg(number n);                         /* i := deg(n) */
     40int     naSize(number n);                           /* size desc. */
     41int     naInt(number &n, const ring r);
     42BOOLEAN naIsZero(number za);                         /* za = 0 ? */
     43BOOLEAN naIsOne(number  za);                         /* za = 1 ? */
     44BOOLEAN naIsMOne(number  za);                        /* za = -1 ? */
     45BOOLEAN naEqual(number a, number b);                  /* a = b ? */
     46BOOLEAN naGreater(number a, number b);                  /* dummy */
     47number  naNeg(number za);                           /* za := - za */
     48number  naInvers(number a);
     49void    naPower(number x, int exp, number *lo);
     50BOOLEAN naGreaterZero(number a);
     51number  naCopy(number p);                              /* erg:= p */
     52number  na_Copy(number p, const ring r);               /* erg:= p */
     53number  naAdd(number la, number li);               /* lu := la+li */
     54number  naMult(number la, number li);              /* lo := la*li */
     55number  naDiv(number la, number li);               /* lo := la/li */
     56number  naIntDiv(number la, number li);            /* lo := la/li */
     57//number  naIntMod(number la, number li);            /* lo := la/li */
     58number  naSub(number la, number li);               /* lu := la-li */
     59void    naNormalize(number &p);
     60number  naGcd(number a, number b, const ring r);
     61number  naLcm(number a, number b, const ring r);
     62const char *  naRead(const char * s, number * p);
     63void    naWrite(number &p, const ring r);
     64char *  naName(number n);
     65nMapFunc naSetMap(const ring src, const ring dst);
     66number naMap0P(number c);
     67number naMap00(number c);
     68#ifdef LDEBUG
     69BOOLEAN naDBTest(number a, const char *f,const int l);
    11070#endif
    11171
     72void    naSetIdeal(ideal I);
     73
     74// external access to the interna
     75poly naPermNumber(number z, int * par_perm, int P, ring r);
     76#define napAddExp(p,i,e)       (p_AddExp(p,i,e,currRing->algring))
     77#define napLength(p)           pLength(p)
     78#define napNeg(p)              (p_Neg(p,currRing->algring))
     79#define napVariables           naNumbOfPar
     80#define napGetCoeff(p)         pGetCoeff(p)
     81#define napGetExpFrom(p,i,r)   (p_GetExp(p,i,r->algring))
     82#define napSetExp(p,i,e)       (p_SetExp(p,i,e,currRing->algring))
     83#define napNew()               (p_Init(currRing->algring))
     84#define napAdd(p1,p2)          (p_Add_q(p1,p2,currRing->algring))
     85#define napSetm(p)             p_Setm(p,currRing->algring)
     86#define nanumber               lnumber
     87napoly napRemainder(napoly f, const napoly  g);
     88number   naGetDenom(number &n, const ring r);
     89number   naGetNumerator(number &n, const ring r);
     90number  napGetDenom(number &n, const ring r);
     91number  napGetNumerator(number &n, const ring r);
     92poly    napPermNumber(number z, int * par_perm, int P, ring r);
     93#endif
     94
  • kernel/longtrans.cc

    r229530 r3dbee61  
    88*/
    99
     10#if 0
    1011#include <stdio.h>
    1112#include <string.h>
     
    25262527}
    25272528#endif
     2529#endif
  • kernel/longtrans.h

    r229530 r3dbee61  
    1010*/
    1111#include <kernel/structs.h>
     12#include <kernel/longalg.h>
     13#if 0
    1214#include <kernel/longrat.h>
    1315#include <kernel/polys-impl.h>
     
    168170
    169171#endif
     172#endif
    170173
  • kernel/numbers.cc

    r229530 r3dbee61  
    140140  else if (rField_is_Extension(r))
    141141  {
    142     if (r->minpoly != NULL)
    143     {
    144142      naSetChar(c,r);
    145143      if (rField_is_Q_a()) nInit_bigint=naMap00;
    146144      if (rField_is_Zp_a()) nInit_bigint=naMap0P;
    147     }
    148     else
    149     {
    150       ntSetChar(c,r);
    151       if (rField_is_Q_a()) nInit_bigint=ntMap00;
    152       if (rField_is_Zp_a()) nInit_bigint=ntMap0P;
    153     }
    154145  }
    155146#ifdef HAVE_RINGS
     
    324315  if (rField_is_Extension(r))
    325316  {
    326     //ntInitChar(c,TRUE,r);
    327     n->cfDelete       = ntDelete;
    328     n->nNormalize     = ntNormalize;
    329     n->cfInit         = ntInit;
    330     n->nPar           = ntPar;
    331     n->nParDeg        = ntParDeg;
    332     n->n_Int          = ntInt;
    333     n->nAdd           = ntAdd;
    334     n->nSub           = ntSub;
    335     n->nMult          = ntMult;
    336     n->nDiv           = ntDiv;
    337     n->nExactDiv      = ntDiv;
    338     n->nIntDiv        = ntIntDiv;
    339     n->nNeg           = ntNeg;
    340     n->nInvers        = ntInvers;
    341     n->nCopy          = ntCopy;
    342     n->cfCopy         = nt_Copy;
    343     n->nGreater       = ntGreater;
    344     n->nEqual         = ntEqual;
    345     n->nIsZero        = ntIsZero;
    346     n->nIsOne         = ntIsOne;
    347     n->nIsMOne        = ntIsMOne;
    348     n->nGreaterZero   = ntGreaterZero;
    349     n->cfWrite        = ntWrite;
    350     n->nRead          = ntRead;
    351     n->nPower         = ntPower;
    352     n->nGcd           = ntGcd;
    353     n->nLcm           = ntLcm;
    354     n->cfSetMap       = ntSetMap;
    355     n->nName          = ntName;
    356     n->nSize          = ntSize;
     317    //naInitChar(c,TRUE,r);
     318    n->cfDelete       = naDelete;
     319    n->nNormalize     = naNormalize;
     320    n->cfInit         = naInit;
     321    n->nPar           = naPar;
     322    n->nParDeg        = naParDeg;
     323    n->n_Int          = naInt;
     324    n->nAdd           = naAdd;
     325    n->nSub           = naSub;
     326    n->nMult          = naMult;
     327    n->nDiv           = naDiv;
     328    n->nExactDiv      = naDiv;
     329    n->nIntDiv        = naIntDiv;
     330    n->nNeg           = naNeg;
     331    n->nInvers        = naInvers;
     332    n->nCopy          = naCopy;
     333    n->cfCopy         = na_Copy;
     334    n->nGreater       = naGreater;
     335    n->nEqual         = naEqual;
     336    n->nIsZero        = naIsZero;
     337    n->nIsOne         = naIsOne;
     338    n->nIsMOne        = naIsMOne;
     339    n->nGreaterZero   = naGreaterZero;
     340    n->cfWrite        = naWrite;
     341    n->nRead          = naRead;
     342    n->nPower         = naPower;
     343    n->nGcd           = naGcd;
     344    n->nLcm           = naLcm;
     345    n->cfSetMap       = naSetMap;
     346    n->nName          = naName;
     347    n->nSize          = naSize;
    357348    n->cfGetDenom     = napGetDenom;
    358349    n->cfGetNumerator = napGetNumerator;
    359350#ifdef LDEBUG
    360     n->nDBTest        = ntDBTest;
     351    n->nDBTest        = naDBTest;
    361352#endif
    362353  }
Note: See TracChangeset for help on using the changeset viewer.