Changeset c6d1377 in git


Ignore:
Timestamp:
Jan 29, 2003, 5:04:18 PM (21 years ago)
Author:
Viktor Levandovskyy <levandov@…>
Branches:
(u'spielwiese', '8e0ad00ce244dfd0756200662572aef8402f13d5')
Children:
18e9f2f7714c20737d239d02c81a2a64fcc2e65b
Parents:
1c3523c71debb3d63cefdfa2e1548fd8078e8416
Message:
new NC multiplication and analysis tools


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

Legend:

Unmodified
Added
Removed

  • Singular/extra.cc

    r1c3523c rc6d1377  
    22*  Computer Algebra System SINGULAR      *
    33*****************************************/
    4 /* $Id: extra.cc,v 1.191 2002-11-21 13:37:20 Singular Exp $ */
     4/* $Id: extra.cc,v 1.192 2003-01-29 16:04:18 levandov Exp $ */
    55/*
    66* ABSTRACT: general interface to internals of Singular ("system" command)
    … …  
    13151315    else
    13161316#ifdef HAVE_PLURAL
     1317/*==================== PrintMat  =================*/
     1318      if (strcmp(sys_cmd, "PrintMat") == 0)
     1319      {
     1320        int a;
     1321        int b;
     1322        ring r;
     1323        int metric;
     1324        if ((h!=NULL) && (h->Typ()==INT_CMD))
     1325        {
     1326          a=(int)h->CopyD();
     1327          h=h->next;
     1328        }
     1329        else return TRUE;
     1330        if ((h!=NULL) && (h->Typ()==INT_CMD))
     1331        {
     1332          b=(int)h->CopyD();
     1333          h=h->next;
     1334        }
     1335        else return TRUE;
     1336        if ((h!=NULL) && (h->Typ()==RING_CMD))
     1337        {
     1338          r=(ring)h->Data();
     1339          h=h->next;
     1340        }
     1341        else return TRUE;
     1342        if ((h!=NULL) && (h->Typ()==INT_CMD))
     1343        {
     1344          metric=(int)h->CopyD();
     1345        }
     1346        res->rtyp=MATRIX_CMD;
     1347        if (rIsPluralRing(r)) res->data=nc_PrintMat(a,b,r,metric);
     1348        else res->data=NULL;
     1349        return FALSE;
     1350      }
    13171351/*==================== twostd  =================*/
    13181352      if (strcmp(sys_cmd, "twostd") == 0)

  • Singular/gring.cc

    r1c3523c rc6d1377  
    77 *  Author:  levandov (Viktor Levandovsky)
    88 *  Created: 8/00 - 11/00
    9  *  Version: $Id: gring.cc,v 1.18 2002-12-20 17:16:58 levandov Exp $
     9 *  Version: $Id: gring.cc,v 1.19 2003-01-29 16:04:17 levandov Exp $
    1010 *******************************************************************/
    1111#include "mod2.h"
    … …  
    593593}
    594594
    595 poly nc_uu_Mult_ww (int i, int a, int j, int b, const ring r)
     595poly nc_uu_Mult_ww_vert (int i, int a, int j, int b, const ring r)
    596596{
    597597  poly out=pOne();
    … …  
    707707  t=MATELEM(cMT,a,b);
    708708  return(p_Copy(t,r));  /* as the last computed element was cMT[a,b] */
     709}
     710
     711poly nc_uu_Mult_ww_horvert (int i, int a, int j, int b, const ring r)
     712  /* (x_i)^a times (x_j)^b */
     713  /* x_i = y,  x_j = x ! */
     714{
     715  /* Check zero exeptions, (q-)commutativity and should we do something at all?  */
     716  assume(a!=0);
     717  assume(b!=0);
     718  poly out=pOne();
     719  if (i<=j)
     720  {
     721    p_SetExp(out,i,a,r);
     722    p_AddExp(out,j,b,r);
     723    p_Setm(out,r);
     724    return(out);
     725  }/* zero exeptions and usual case */
     726  /*  if ((a==0)||(b==0)||(i<=j)) return(out); */
     727 
     728  if (MATELEM(r->nc->COM,j,i)!=NULL)
     729    /* commutative or quasicommutative case */
     730  {
     731    if (r->cf->nIsOne(p_GetCoeff(MATELEM(r->nc->COM,j,i),r))) /* commutative case */
     732    {
     733      return(out);
     734    }     
     735    else
     736    {
     737      number tmp_number=p_GetCoeff(MATELEM(r->nc->COM,j,i),r); /* quasicommutative case */
     738      nPower(tmp_number,a*b,&tmp_number);
     739      p_SetCoeff(out,tmp_number,r);
     740      return(out);
     741    }
     742  }/* end_of commutative or quasicommutative case */
     743  p_Delete(&out,r);
     744
     745  /* we are here if  i>j and variables do not commute or quasicommute */
     746  /* in fact, now a>=1 and b>=1; and j<i */
     747  /* now check whether the polynomial is already computed */
     748
     749  int vik = UPMATELEM(j,i,r->N);
     750  int cMTsize=r->nc->MTsize[vik];
     751  int newcMTsize=0;
     752  if (a>b) {newcMTsize=a;} else {newcMTsize=b;}
     753
     754  if (newcMTsize<=cMTsize)
     755  {
     756    if ( MATELEM(r->nc->MT[vik],a,b)!=NULL)
     757    {
     758      out=p_Copy(MATELEM(r->nc->MT[vik],a,b),r);
     759      return (out);
     760    }
     761  }
     762  int k,m;   
     763  if (newcMTsize > cMTsize)
     764  {
     765    number nM = nInit(newcMTsize);
     766    number cM = nInit(7);
     767    int inM;
     768    number nM2 = nDiv(nM,cM);
     769    nDelete(&nM);
     770    nDelete(&cM);
     771    if (nIsZero(nM2))
     772    {
     773      /* Indicates |nM/cM| > MAXINT */
     774      /* Should never get there */
     775      WarnS("too big matrix in nc_uu_Mult_ww");
     776    }
     777    else
     778    {
     779      inM = nInt(nM2);
     780      inM=inM*7;/* 7=DefMTSize in extra.cc */
     781      if (inM<newcMTsize) inM=inM+7;
     782    }
     783    nDelete(&nM2);
     784    newcMTsize = inM;
     785    //    matrix tmp = (matrix)omAlloc0(inM*inM*sizeof(poly));
     786    matrix tmp = mpNew(newcMTsize,newcMTsize);
     787     
     788    for (k=1;k<=cMTsize;k++)
     789    {
     790      for (m=1;m<=cMTsize;m++)
     791      {
     792        if ( MATELEM(r->nc->MT[UPMATELEM(j,i,r->N)],k,m) != NULL )
     793        {
     794          MATELEM(tmp,k,m) = MATELEM(r->nc->MT[UPMATELEM(j,i,r->N)],k,m);
     795          //       omCheckAddr(tmp->m);
     796          MATELEM(r->nc->MT[UPMATELEM(j,i,r->N)],k,m)=NULL;
     797          //       omCheckAddr(r->nc->MT[UPMATELEM(j,i,r->N)]->m);
     798        }
     799      }
     800    }
     801    id_Delete((ideal *)&(r->nc->MT[UPMATELEM(j,i,r->N)]),r);
     802    r->nc->MT[UPMATELEM(j,i,r->N)] = tmp;
     803    tmp=NULL;
     804    r->nc->MTsize[UPMATELEM(j,i,r->N)] = newcMTsize;
     805  }  /* The update of multiplication matrix is finished */
     806
     807  matrix cMT=r->nc->MT[UPMATELEM(j,i,r->N)];         /* cMT=current MT */
     808
     809  poly x=pOne();p_SetExp(x,j,1,r);p_Setm(x,r);p_Test(x,r); /* var(j); */
     810  poly y=pOne();p_SetExp(y,i,1,r);p_Setm(y,r);p_Test(y,r); /*var(i);  for convenience */
     811 
     812  poly t=NULL;
     813
     814  int toXY;
     815  int toYX;
     816
     817  if (a==1) /* y*x^b, b>=2 */
     818  {
     819    toXY=b-1;
     820    while ( (MATELEM(cMT,1,toXY)==NULL) && (toXY>=2)) toXY--;
     821    for (m=toXY+1;m<=b;m++)
     822    {
     823      t=MATELEM(cMT,1,m);
     824      if (t==NULL)   /* remove after debug */
     825      {
     826        t = p_Copy(MATELEM(cMT,1,m-1),r);
     827        t = nc_p_Mult_mm(t,x,r);
     828        MATELEM(cMT,1,m) = t;
     829        /*      omCheckAddr(cMT->m); */
     830      }
     831      else
     832      {   
     833        /* Error, should never get there */
     834        WarnS("Error: a=1; MATELEM!=0");
     835      }
     836      t=NULL;
     837    }
     838    return(p_Copy(MATELEM(cMT,1,b),r));
     839  }
     840
     841  if (b==1) /* y^a*x, a>=2 */
     842  {
     843    toYX=a-1;
     844    while ( (MATELEM(cMT,toYX,1)==NULL) && (toYX>=2)) toYX--;
     845    for (m=toYX+1;m<=a;m++)
     846    {
     847      t=MATELEM(cMT,m,1);
     848      if (t==NULL)   /* remove after debug */
     849      {
     850        t = p_Copy(MATELEM(cMT,m-1,1),r);
     851        t = nc_mm_Mult_p(y,t,r);
     852        MATELEM(cMT,m,1) = t;
     853        /*      omCheckAddr(cMT->m); */
     854      }
     855      else
     856      {   
     857        /* Error, should never get there */
     858        WarnS("Error: b=1, MATELEM!=0");
     859      }
     860      t=NULL;
     861    }
     862    return(p_Copy(MATELEM(cMT,a,1),r));
     863  }
     864
     865/* ------------ Main Cycles ----------------------------*/
     866  /*            a>1, b>1              */
     867
     868  int dXY=0; int dYX=0;
     869  /* dXY = distance for computing x-mult, then y-mult */
     870  /* dYX = distance for computing y-mult, then x-mult */
     871  int toX=a-1; int toY=b-1; /* toX = to axe X, toY = to axe Y */
     872  toXY=b-1; toYX=a-1;
     873  /* if toX==0, toXY = dist. to computed y * x^toXY */
     874  /* if toY==0, toYX = dist. to computed y^toYX * x */
     875  while ( (MATELEM(cMT,toX,b)==NULL) && (toX>=1)) toX--;
     876  if (toX==0) /* the whole column is not computed yet */
     877  {
     878    while ( (MATELEM(cMT,1,toXY)==NULL) && (toXY>=1)) toXY--;
     879    /* toXY >=1 */
     880    dXY=b-1-toXY;
     881  }
     882  dXY=dXY+a-toX; /* the distance to nearest computed y^toX x^b */
     883
     884  while ( (MATELEM(cMT,a,toY)==NULL) && (toY>=1)) toY--;
     885  if (toY==0) /* the whole row is not computed yet */
     886  {
     887    while ( (MATELEM(cMT,toYX,1)==NULL) && (toYX>=1)) toYX--;
     888    /* toYX >=1 */
     889    dYX=a-1-toYX;
     890  }
     891  dYX=dYX+b-toY; /* the distance to nearest computed y^a x^toY */
     892 
     893  if (dYX>=dXY)
     894  {
     895    /* first x, then y */
     896    if (toX==0) /* start with the row*/
     897    {
     898      for (m=toXY+1;m<=b;m++)
     899      {
     900        t=MATELEM(cMT,1,m);
     901        if (t==NULL)   /* remove after debug */
     902        {
     903          t = p_Copy(MATELEM(cMT,1,m-1),r);
     904          t = nc_p_Mult_mm(t,x,r);
     905          MATELEM(cMT,1,m) = t;
     906          /*    omCheckAddr(cMT->m); */
     907        }
     908        else
     909        {         
     910          /* Error, should never get there */
     911          WarnS("dYX>=dXY,toXY; MATELEM==0");
     912        }
     913        t=NULL;
     914      }
     915      toX=1; /* y*x^b is computed */
     916    }
     917    /* Now toX>=1 */
     918    for (k=toX+1;k<=a;k++)
     919    {
     920      t=MATELEM(cMT,k,b);
     921      if (t==NULL)   /* remove after debug */
     922      {
     923        t = p_Copy(MATELEM(cMT,k-1,b),r);
     924        t = nc_mm_Mult_p(y,t,r);
     925        MATELEM(cMT,k,b) = t;
     926        /*      omCheckAddr(cMT->m); */
     927      }
     928      else
     929      {
     930        /* Error, should never get there */
     931        WarnS("dYX>=dXY,toX; MATELEM==0");
     932      }
     933      t=NULL;
     934    }
     935  } /* endif (dYX>=dXY) */
     936
     937
     938  if (dYX<dXY)
     939  {
     940    /* first y, then x */
     941    if (toY==0) /* start with the column*/
     942    {
     943      for (m=toYX+1;m<=a;m++)
     944      {
     945        t=MATELEM(cMT,m,1);
     946        if (t==NULL)   /* remove after debug */
     947        {
     948          t = p_Copy(MATELEM(cMT,m-1,1),r);
     949          t = nc_mm_Mult_p(y,t,r);
     950          MATELEM(cMT,m,1) = t;
     951          /*    omCheckAddr(cMT->m); */
     952        }
     953        else
     954        {         
     955          /* Error, should never get there */
     956          WarnS("dYX<dXY,toYX; MATELEM==0");
     957        }
     958        t=NULL;
     959      }
     960      toY=1; /* y^a*x is computed */
     961    }
     962    /* Now toY>=1 */
     963    for (k=toY+1;k<=b;k++)
     964    {
     965      t=MATELEM(cMT,a,k);
     966      if (t==NULL)   /* remove after debug */
     967      {
     968        t = p_Copy(MATELEM(cMT,a,k-1),r);
     969        t = nc_p_Mult_mm(t,x,r);
     970        MATELEM(cMT,a,k) = t;
     971        /*      omCheckAddr(cMT->m); */
     972      }
     973      else
     974      {
     975        /* Error, should never get there */
     976        WarnS("dYX<dXY,toY; MATELEM==0");
     977      }
     978      t=NULL;
     979    }
     980  } /* endif (dYX<dXY) */
     981       
     982  p_Delete(&x,r);
     983  p_Delete(&y,r);
     984  t=p_Copy(MATELEM(cMT,a,b),r);
     985  return(t);  /* as the last computed element was cMT[a,b] */
    709986}
    710987
    … …  
    11411418}
    11421419
     1420matrix nc_PrintMat(int a, int b, ring r, int metric)
     1421  /* returns matrix with the info on noncomm multiplication */
     1422{
     1423 
     1424  if ( (a==b) || !rIsPluralRing(r) ) return(NULL);
     1425  int i;
     1426  int j;
     1427  if (a>b) {j=b; i=a;}
     1428  else {j=a; i=b;}
     1429  /* i<j */
     1430  int size=r->nc->MTsize[UPMATELEM(i,j,r->N)];
     1431  matrix M = r->nc->MT[UPMATELEM(i,j,r->N)];
     1432  /*  return(M); */
     1433  int sizeofres;
     1434  if (metric==0)
     1435  {
     1436    sizeofres=sizeof(int);
     1437  }
     1438  if (metric==1)
     1439  {
     1440    sizeofres=sizeof(number);
     1441  }
     1442  matrix res=mpNew(size,size);
     1443  int s;
     1444  int t;
     1445  int length;
     1446  long totdeg;
     1447  poly p;
     1448  for(s=1;s<=size;s++)
     1449  {
     1450    for(t=1;t<=size;t++)
     1451    {
     1452      p=MATELEM(M,s,t);
     1453      if (p==NULL)
     1454      {
     1455        MATELEM(res,s,t)=0;
     1456      }
     1457      else
     1458      {
     1459        length = pLength(p);
     1460        if (metric==0) /* length */
     1461        {       
     1462          MATELEM(res,s,t)= p_ISet(length,r);
     1463        }
     1464        else if (metric==1) /* sum of deg divided by the length */
     1465        {
     1466          totdeg=0;
     1467          while (p!=NULL)
     1468          {
     1469            totdeg=totdeg+pDeg(p,r);
     1470            pIter(p);
     1471          }
     1472          number ntd = nInit(totdeg);
     1473          number nln = nInit(length);
     1474          number nres=nDiv(ntd,nln);
     1475          nDelete(&ntd);
     1476          nDelete(&nln);
     1477          MATELEM(res,s,t)=p_NSet(nres,r);       
     1478        }
     1479      }
     1480    }
     1481  }
     1482  return(res);
     1483}
     1484
    11431485#endif

  • Singular/gring.h

    r1c3523c rc6d1377  
    44*  Computer Algebra System SINGULAR     *
    55****************************************/
    6 /* $Id: gring.h,v 1.14 2002-12-06 20:51:42 levandov Exp $ */
     6/* $Id: gring.h,v 1.15 2003-01-29 16:04:18 levandov Exp $ */
    77/*
    88* ABSTRACT additional defines etc for --with-plural
    … …  
    2121poly nc_mm_Mult_nn (Exponent_t *F, Exponent_t *G, const ring r);
    2222poly nc_mm_Mult_uu (Exponent_t *F,int jG,int bG, const ring r);
    23 poly nc_uu_Mult_ww (int i, int a, int j, int b, const ring r);
     23
     24#define nc_uu_Mult_ww nc_uu_Mult_ww_horvert
     25poly nc_uu_Mult_ww_vert (int i, int a, int j, int b, const ring r);
     26poly nc_uu_Mult_ww_horvert (int i, int a, int j, int b, const ring r);
     27poly nc_uu_Mult_ww_hvdiag (int i, int a, int j, int b, const ring r);
     28
     29
    2430poly _nc_p_Mult_q(poly p, poly q, const int copy, const ring r);
     31
    2532//syzygies :
    2633poly nc_spGSpolyCreate(poly p1, poly p2,poly spNoether, const ring r);
    … …  
    4653void nc_PolyPolyRed(poly &b, poly p, number *c);
    4754
     55matrix nc_PrintMat(int a, int b, ring r, int metric);
     56
    4857#else
    4958// dummy definition to make gcc happy
Note: See TracChangeset for help on using the changeset viewer.