Changeset 4d9dc5 in git


Ignore:
Timestamp:
Jun 8, 1999, 11:15:20 AM (25 years ago)
Author:
Mathias Schulze <mschulze@…>
Branches:
(u'spielwiese', 'fe61d9c35bf7c61f2b6cbf1b56e25e2f08d536cc')
Children:
559d1ee5651bb2aab60891dce7d68b2348256ca7
Parents:
bb21e33b4dde3d1a2778e90aee11ade770bea780
Message:
*** empty log message ***


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

Legend:

Unmodified
Added
Removed

  • Singular/LIB/jordan.lib

    rbb21e3 r4d9dc5  
    11///////////////////////////////////////////////////////////////////////////////
    22
    3 version="$Id: jordan.lib,v 1.9 1999-02-18 19:29:01 mschulze Exp $";
     3version="$Id: jordan.lib,v 1.10 1999-06-08 09:15:20 mschulze Exp $";
    44info="
    55LIBRARY: jordan.lib  PROCEDURES TO COMPUTE THE JORDAN NORMAL FORM
    66AUTHOR:  Mathias Schulze, email: mschulze@mathematik.uni-kl.de
    77
    8  jordan(M[,opt]);  eigenvalues, Jordan block sizes, transformation matrix of M
     8 jordan(M[,opt]);  eigenvalues, Jordan block sizes, Jordan transformation of M
    99 jordanmatrix(l);  Jordan matrix with eigenvalues l[1], Jordan block sizes l[2]
    1010 jordanform(M);    Jordan normal form of constant square matrix M
    11  inversemat(M);    inverse matrix of invertible constant matrix M
     11 invmat(M);        inverse matrix of invertible constant matrix M
    1212";
    1313
    … …  
    1616
    1717proc jordan(matrix M,list #)
    18 "USAGE:   jordan(M[,opt]); M constant square matrix, opt integer
    19 ASSUME:  the eigenvalues of M are in the coefficient field
    20 RETURN:  a list l with 3 entries of type ideal, list of intvecs, matrix with
    21            l[1]       : eigenvalues of M
    22            l[2][i][j] : size of j-th Jordan block with eigenvalue l[1][i]
    23            l[3] : transformation matrix U with U^(-1)*M*U in Jordan normal form
    24          which entries of l are computed depends on opt:
    25            opt=-1 : compute l[1]
    26            opt= 0 : compute l[1] and l[2] (default)
    27            opt= 1 : compute l[1],l[2] and l[3]
    28            opt= 2 : compute l[1] and l[3]
    29 NOTE:    a non constant polynomial matrix M is replaced by its constant part
    30 EXAMPLE: example jordan; shows an example
     18"USAGE:   jordan(M[,opt]); with M constant square matrix and opt integer.
     19ASSUME:  The eigenvalues of M are in the coefficient field.
     20RETURN:  The procedure returns a list l with 3 entries of type
     21         ideal, list of intvecs, matrix with
     22         l[1] eigenvalues of M,
     23         l[2][i][j] size of j-th Jordan block with eigenvalue l[1][i], and
     24         l[3]^(-1)*M*l[3] in Jordan normal form.
     25         Depending on opt, only certain entries of l are computed.
     26         If opt=-1, l[1] is computed,
     27         if opt= 0, l[1] and l[2] are computed,
     28         if opt= 1, l[1],l[2], and l[3] are computed, and,
     29         if opt= 2, l[1] and l[3] are computed.
     30         By default, opt=0.
     31NOTE:    A non constant polynomial matrix M is replaced by its constant part.
     32EXAMPLE: example jordan; shows an example.
    3133"
    3234{
    33   // test if square matrix
     35  // test if M is a square matrix
    3436  int n=nrows(M);
    3537  if(n!=ncols(M))
    3638  {
    37     "//no square matrix";
    38     return();
    39   }
    40 
    41   // get constant part
     39    print("//no square matrix");
     40    return(list());
     41  }
     42
     43  // replace M by its constant part
     44  M=jet(M,0);
     45
     46  // change to polynomial ring for factorization
    4247  def br=basering;
    43   map zero=br,0;
    44   M=zero(M);
    45   kill zero;
    46 
    47   // change to polynomial ring for factorization
    4848  changeord("pr","dp");
    4949  matrix M=imap(br,M);
    5050
    51   // factorize characteristic polynomial
     51  // factorize characteristic polynomial of M
    5252  list l=factorize(det(M-var(1)*freemodule(n)),2);
    5353
    54   // get multiplicities mM
     54  // get multiplicities mM of M
    5555  def eM,mM=l[1..2];
    5656  kill l;
    5757
    58   // test if eigenvalues in the coefficient field
     58  // test if eigenvalues of M are in the coefficient field
    5959  int i;
    6060  for(i=size(eM);i>=1;i--)
    … …  
    6262    if(deg(eM[i])>1)
    6363    {
     64      print("//eigenvalues not in the coefficient field");
     65      setring br;
    6466      kill pr;
    65       "//eigenvalues not in the coefficient field";
    66       return();
    67     }
    68   }
    69 
    70   // get eigenvalues eM
     67      return(list());
     68    }
     69  }
     70
     71  // get eigenvalues eM of M
    7172  map inv=pr,-var(1);
    72   eM=simplify(inv(eM),1);
     73  eM=jet(simplify(inv(eM),1),0);
    7374  setring br;
    74   map zero=pr,0;
    75   ideal eM=zero(eM);
     75  ideal eM=fetch(pr,eM);
    7676  kill pr;
    7777
    78   // sort eigenvalues
     78  // sort eigenvalues eM
    7979  int j;
    8080  poly e;
    … …  
    9797  kill e,m;
    9898
    99   // get option parameter
     99  // get optional parameter opt
    100100  int opt=0;
    101101  if(size(#)>0)
    … …  
    111111  }
    112112
    113   // define needed variables
     113  // define variables
    114114  int k,l;
    115115  matrix E=freemodule(n);
    … …  
    128128  }
    129129
    130   // do the following for all eigenvalues eM[i]
     130  // do for all eigenvalues eM[i]
    131131  for(i=ncols(eM);i>=1;i--)
    132132  {
    … …  
    143143    if(opt<=1)
    144144    {
    145       // compute Jordan block size vector corresponding to eigenvalue eM[i]
     145      // compute Jordan block size vector bMi corresponding to eigenvalue eM[i]
    146146      bMi=0;
    147147      bMi[size(K[2])]=0;
    … …  
    158158    if(opt>=1)
    159159    {
    160       // compute generating vectors for Jordan basis vectors corresponding to
    161       // eigenvalue eM[i]
     160      // compute Jordan basis vectors K corresponding to eigenvalue eM[i]
    162161      if(size(K)>1)
    163162      {
    … …  
    170169        K[j]=interred(reduce(K[j],std(K1)));
    171170      }
    172 
    173       // compute Jordan basis vectors corresponding to eigenvalue eM[i] from
    174       // generating vectors
    175171      for(j=size(K);j>=2;j--)
    176172      {
    … …  
    189185  kill l;
    190186
    191   // create return list
     187  // create return list l
    192188  list l=eM;
    193189  if(opt<=1)
    … …  
    212208
    213209proc jordanmatrix(list l)
    214 "USAGE:   jordanmatrix(l); l list of ideal and list of intvecs
    215 RETURN:  the Jordan matrix J defined by l:
    216            l[1]       : eigenvalues of J
    217            l[2][i][j] : size of j-th Jordan block with eigenvalue l[1][i]
    218 EXAMPLE: example jordanmatrix; shows an example
     210"USAGE:   jordanmatrix(l); with l list of ideal and list of intvecs.
     211RETURN:  The procedure returns the Jordan matrix J with eigenvalues l[1] and
     212         size l[2][i][j] of j-th Jordan block with eigenvalue l[1][i].
     213EXAMPLE: example jordanmatrix; shows an example.
    219214"
    220215{
    … …  
    222217  if(size(l)<2)
    223218  {
    224     "//not enough entries in argument list";
    225     return();
     219    print("//not enough entries in argument list");
     220    matrix J[1][0];
     221    return(J);
    226222  }
    227223  def eJ,bJ=l[1..2];
    … …  
    229225  if(typeof(eJ)!="ideal")
    230226  {
    231     "//first entry in argument list not an ideal";
    232     return();
     227    print("//first entry in argument list not an ideal");
     228    matrix J[1][0];
     229    return(J);
    233230  }
    234231  if(typeof(bJ)!="list")
    235232  {
    236     "//second entry in argument list not a list";
    237     return();
     233    print("//second entry in argument list not a list");
     234    matrix J[1][0];
     235    return(J);
    238236  }
    239237  if(size(eJ)<size(bJ))
    … …  
    252250    if(typeof(bJ[i])!="intvec")
    253251    {
    254       "//second entry in argument list not a list of integer vectors";
    255       return();
     252      print("//second entry in argument list not a list of intvecs");
     253      matrix J[1][0];
     254      return(J);
    256255    }
    257256    else
    … …  
    267266    }
    268267  }
    269   matrix J[n][n];
    270268
    271269  // create Jordan matrix
    272270  int l;
     271  matrix J[n][n];
    273272  for(i,l=1,1;i<=s;i++)
    274273  {
    … …  
    303302
    304303proc jordanform(matrix M)
    305 "USAGE:   jordanform(M); M constant square matrix
    306 ASSUME:  the eigenvalues of M are in the coefficient field
    307 RETURN:  the Jordan normal form of M
    308 NOTE:    a non constant polynomial matrix M is replaced by its constant part
    309 EXAMPLE: example jordanform; shows an example
     304"USAGE:   jordanform(M); with M constant square matrix.
     305ASSUME:  The eigenvalues of M are in the coefficient field.
     306RETURN:  The procedure returns the Jordan normal form of M.
     307NOTE:    A non constant polynomial matrix M is replaced by its constant part.
     308EXAMPLE: example jordanform; shows an example.
    310309"
    311310{
    … …  
    321320///////////////////////////////////////////////////////////////////////////////
    322321
    323 proc inversemat(matrix M)
    324 "USAGE:   inversemat(M); M constant square matrix
    325 ASSUME:  M is invertible
    326 RETURN:  the inverse matrix of M
    327 NOTE:    a non constant polynomial matrix M is replaced by its constant part
    328 EXAMPLE: example inversemat; shows an example
     322proc invmat(matrix M)
     323"USAGE:   invmat(M); with M constant square matrix.
     324ASSUME:  M is invertible.
     325RETURN:  The procedure returns the inverse matrix of M.
     326NOTE:    A non constant polynomial matrix M is replaced by its constant part.
     327EXAMPLE: example invmat; shows an example.
    329328"
    330329{
    331   // test if square matrix
    332   int n=nrows(M);
    333   if(n!=ncols(M))
    334   {
    335     "//no square matrix";
    336     return();
    337   }
    338 
    339   // get constant part
    340   def br=basering;
    341   map zero=br,0;
    342   M=zero(M);
    343 
    344   // compute inverse
    345   return(lift(M,freemodule(n)));
     330  if(nrows(M)==ncols(M))
     331  {
     332    matrix invM=lift(jet(M,0),freemodule(nrows(M)));
     333  }
     334  else
     335  {
     336    print("//no square matrix");
     337    matrix[1][0]=invM;
     338  }
     339  return(invM);
    346340}
    347341example
    … …  
    350344  matrix M[3][3]=3,2,1,0,2,1,0,0,3;
    351345  print(M);
    352   print(inversemat(M));
    353 }
    354 ///////////////////////////////////////////////////////////////////////////////
    355 
     346  print(invmat(M));
     347}
     348///////////////////////////////////////////////////////////////////////////////

  • Singular/pcv.cc

    rbb21e3 r4d9dc5  
    22*  Computer Algebra System SINGULAR      *
    33*****************************************/
    4 /* $Id: pcv.cc,v 1.19 1999-04-29 15:09:32 Singular Exp $ */
     4/* $Id: pcv.cc,v 1.20 1999-06-08 09:13:54 mschulze Exp $ */
    55/*
    66* ABSTRACT: conversion between polys and coef vectors
    … …  
    88
    99#include "mod2.h"
     10
    1011#ifdef HAVE_PCV
    1112#if !defined(HAVE_DYNAMIC_LOADING) || defined(BUILD_MODULE)
    … …  
    1617#include "polys.h"
    1718#include "lists.h"
     19#include "matpol.h"
    1820#include "febase.h"
    1921#include "pcv.h"
    2022
    21 static int pcvDegBound;
     23static int pcvMaxDegree;
    2224static int pcvTableSize;
    2325static int pcvIndexSize;
    … …  
    2527static unsigned** pcvIndex=NULL;
    2628
    27 #ifndef HAVE_DYNAMIC_LOADING
    28 /* Without dynamic-loading we need to provides following functions */
     29int pcvDeg(poly p)
     30{
     31  int d=0;
     32  for(int i=pVariables;i>=1;i--) d+=pGetExp(p,i);
     33  return d;
     34}
     35
     36int pcvMinDeg(poly p)
     37{
     38  if(!p) return -1;
     39  int md=pcvDeg(p);
     40  pIter(p);
     41  while(p)
     42  {
     43    int d=pcvDeg(p);
     44    if(d<md) md=d;
     45    pIter(p);
     46  }
     47  return md;
     48}
     49
     50int pcvMinDeg(matrix m)
     51{
     52  int i,j,d;
     53  int md=-1;
     54  for(i=1;i<=MATROWS(m);i++)
     55  {
     56    for(j=1;j<=MATCOLS(m);j++)
     57    {
     58      d=pcvMinDeg(MATELEM(m,i,j));
     59      if((d>=0&&md>d)||md==-1) md=d;
     60    }
     61  }
     62  return(md);
     63}
     64
     65int pcvMaxDeg(poly p)
     66{
     67  if(!p) return -1;
     68  int md=pcvDeg(p);
     69  pIter(p);
     70  while(p)
     71  {
     72    int d=pcvDeg(p);
     73    if(d>md) md=d;
     74    pIter(p);
     75  }
     76  return md;
     77}
     78
     79int pcvMaxDeg(matrix m)
     80{
     81  int i,j,d;
     82  int md=-1;
     83  for(i=1;i<=MATROWS(m);i++)
     84  {
     85    for(j=1;j<=MATCOLS(m);j++)
     86    {
     87      d=pcvMinDeg(MATELEM(m,i,j));
     88      if((d>=0&&md<d)||md==-1) md=d;
     89    }
     90  }
     91  return(md);
     92}
    2993
    3094BOOLEAN pcvMinDeg(leftv res,leftv h)
    … …  
    38102      return FALSE;
    39103    }
    40   }
    41   WerrorS("<poly> expected");
     104    else
     105    if(h->Typ()==MATRIX_CMD)
     106    {
     107      res->rtyp=INT_CMD;
     108      res->data=(void*)pcvMinDeg((matrix)h->Data());     
     109      return FALSE;
     110    }
     111  }
     112  WerrorS("<poly> or <matrix> expected");
    42113  return TRUE;
    43114}
    … …  
    53124      return FALSE;
    54125    }
    55   }
    56   WerrorS("<poly> expected");
    57   return TRUE;
     126    else
     127    if(h->Typ()==MATRIX_CMD)
     128    {
     129      res->rtyp=INT_CMD;
     130      res->data=(void*)pcvMaxDeg((matrix)h->Data());     
     131      return FALSE;
     132    }
     133  }
     134  WerrorS("<poly> or <matrix> expected");
     135  return TRUE;
     136}
     137
     138void pcvInit(int d)
     139{
     140  if(d<0) d=0;
     141  pcvMaxDegree=d;
     142  pcvTableSize=pVariables*pcvMaxDegree*sizeof(unsigned);
     143  pcvTable=(unsigned*)Alloc0(pcvTableSize);
     144  pcvIndexSize=pVariables*sizeof(unsigned*);
     145  pcvIndex=(unsigned**)Alloc(pcvIndexSize);
     146  for(int i=0;i<pVariables;i++)
     147    pcvIndex[i]=pcvTable+i*pcvMaxDegree;
     148  for(int i=0;i<pcvMaxDegree;i++)
     149    pcvIndex[0][i]=i;
     150  for(int i=1;i<pVariables;i++)
     151  {
     152    unsigned x=0;
     153    for(int j=0;j<pcvMaxDegree;j++)
     154    {
     155      x+=pcvIndex[i-1][j];
     156      pcvIndex[i][j]=x;
     157    }
     158  }
     159}
     160
     161void pcvClean()
     162{
     163  if(pcvTable)
     164  {
     165    Free(pcvTable,pcvTableSize);
     166    pcvTable=NULL;
     167  }
     168  if(pcvIndex)
     169  {
     170    Free(pcvIndex,pcvIndexSize);
     171    pcvIndex=NULL;
     172  }
     173}
     174
     175int pcvM2N(poly m)
     176{
     177  unsigned n=0,d=0;
     178  for(int i=0;i<pVariables;i++)
     179  {
     180    d+=pGetExp(m,i+1);
     181    n+=pcvIndex[i][d];
     182  }
     183  return n+1;
     184}
     185
     186poly pcvN2M(int n)
     187{
     188  n--;
     189  poly m=pOne();
     190  int i,j,k;
     191  for(i=pVariables-1;i>=0;i--)
     192  {
     193    k=j;
     194    for(j=0;j<pcvMaxDegree&&pcvIndex[i][j]<=n;j++);
     195    j--;
     196    n-=pcvIndex[i][j];
     197    if(i<pVariables-1) pSetExp(m,i+2,k-j);
     198  }
     199  if(n==0)
     200  {
     201    pSetExp(m,1,j);
     202    pSetm(m);
     203    return m;
     204  }
     205  else
     206  {
     207    pDelete1(&m);
     208    return NULL;
     209  }
     210}
     211
     212poly pcvP2CV(poly p,int d0,int d1)
     213{
     214  poly cv=NULL;
     215  while(p)
     216  {
     217    int d=pcvDeg(p);
     218    if(d0<=d&&d<d1)
     219    {
     220      poly c=pOne();
     221      pSetComp(c,pcvM2N(p));
     222      pSetCoeff(c,nCopy(pGetCoeff(p)));
     223      cv=pAdd(cv,c);
     224    }
     225    pIter(p);
     226  }
     227  return cv;
     228}
     229
     230poly pcvCV2P(poly cv,int d0,int d1)
     231{
     232  poly p=NULL;
     233  while(cv)
     234  {
     235    poly m=pcvN2M(pGetComp(cv));
     236    if(m)
     237    {
     238      int d=pcvDeg(m);
     239      if(d0<=d&&d<d1)
     240      {
     241        pSetCoeff(m,nCopy(pGetCoeff(cv)));
     242        p=pAdd(p,m);
     243      }
     244    }
     245    pIter(cv);
     246  }
     247  return p;
     248}
     249
     250lists pcvP2CV(lists pl,int d0,int d1)
     251{
     252  lists cvl=(lists)Alloc(sizeof(slists));
     253  cvl->Init(pl->nr+1);
     254  pcvInit(d1);
     255  for(int i=pl->nr;i>=0;i--)
     256  {
     257    if(pl->m[i].rtyp==POLY_CMD)
     258    {
     259      cvl->m[i].rtyp=VECTOR_CMD;
     260      cvl->m[i].data=pcvP2CV((poly)pl->m[i].data,d0,d1);
     261    }
     262  }
     263  pcvClean();
     264  return cvl;
     265}
     266
     267lists pcvCV2P(lists cvl,int d0,int d1)
     268{
     269  lists pl=(lists)Alloc(sizeof(slists));
     270  pl->Init(cvl->nr+1);
     271  pcvInit(d1);
     272  for(int i=cvl->nr;i>=0;i--)
     273  {
     274    if(cvl->m[i].rtyp==VECTOR_CMD)
     275    {
     276      pl->m[i].rtyp=POLY_CMD;
     277      pl->m[i].data=pcvCV2P((poly)cvl->m[i].data,d0,d1);
     278    }
     279  }
     280  pcvClean();
     281  return pl;
    58282}
    59283
    … …  
    114338}
    115339
     340int pcvDim(int d0,int d1)
     341{
     342  if(d0<0) d0=0;
     343  if(d1<0) d1=0;
     344  pcvInit(d1+1);
     345  int d=pcvIndex[pVariables-1][d1]-pcvIndex[pVariables-1][d0];
     346  pcvClean();
     347  return d;
     348}
     349
    116350BOOLEAN pcvDim(leftv res,leftv h)
    117351{
    … …  
    137371}
    138372
     373int pcvBasis(lists b,int i,poly m,int d,int n)
     374{
     375  if(n<pVariables)
     376  {
     377    for(int k=0,l=d;k<=l;k++,d--)
     378    {
     379      pSetExp(m,n,k);
     380      i=pcvBasis(b,i,m,d,n+1);
     381    }
     382  }
     383  else
     384  {
     385    pSetExp(m,n,d);
     386    pSetm(m);
     387    b->m[i].rtyp=POLY_CMD;
     388    b->m[i++].data=pCopy(m);
     389  }
     390  return i;
     391}
     392
     393lists pcvBasis(int d0,int d1)
     394{
     395  if(d0<0) d0=0;
     396  if(d1<0) d1=0;
     397  lists b=(lists)Alloc(sizeof(slists));
     398  b->Init(pcvDim(d0,d1));
     399  poly m=pOne();
     400  for(int d=d0,i=0;d<d1;d++)
     401    i=pcvBasis(b,i,m,d,1);
     402  pDelete1(&m);
     403  return b;
     404}
     405
    139406BOOLEAN pcvBasis(leftv res,leftv h)
    140407{
    … …  
    160427}
    161428
    162 #endif /* HAVE_DYNAMIC_LOADING */
    163 
    164 int pcvDeg(poly p)
    165 {
    166   int dp=0;
    167   for(int i=1;i<=pVariables;i++) dp+=pGetExp(p,i);
    168   return dp;
    169 }
    170 
    171 int pcvMinDeg(poly p)
    172 {
    173   if(!p) return 0;
    174   int md=pcvDeg(p);
    175   pIter(p);
    176   while(p)
    177   {
    178     int d=pcvDeg(p);
    179     if(d<md) md=d;
    180     pIter(p);
    181   }
    182   return md;
    183 }
    184 
    185 int pcvMaxDeg(poly p)
    186 {
    187   if(!p) return 0;
    188   int md=pcvDeg(p);
    189   pIter(p);
    190   while(p)
    191   {
    192     int d=pcvDeg(p);
    193     if(d>md) md=d;
    194     pIter(p);
    195   }
    196   return md;
    197 }
    198 
    199 void pcvInit(int d)
    200 {
    201   if(d<0) d=0;
    202   pcvDegBound=d;
    203   pcvTableSize=pVariables*pcvDegBound*sizeof(unsigned);
    204   pcvTable=(unsigned*)Alloc0(pcvTableSize);
    205   pcvIndexSize=pVariables*sizeof(unsigned*);
    206   pcvIndex=(unsigned**)Alloc(pcvIndexSize);
    207   int i;
    208   for(i=0;i<pVariables;i++)
    209     pcvIndex[i]=pcvTable+i*pcvDegBound;
    210   for(i=0;i<pcvDegBound;i++)
    211     pcvIndex[0][i]=i;
    212   for(i=1;i<pVariables;i++)
    213   {
    214     unsigned x=0;
    215     for(int j=0;j<pcvDegBound;j++)
    216     {
    217       x+=pcvIndex[i-1][j];
    218       pcvIndex[i][j]=x;
    219     }
    220   }
    221 }
    222 
    223 void pcvClean()
    224 {
    225   if(pcvTable)
    226   {
    227     Free(pcvTable,pcvTableSize);
    228     pcvTable=NULL;
    229   }
    230   if(pcvIndex)
    231   {
    232     Free(pcvIndex,pcvIndexSize);
    233     pcvIndex=NULL;
    234   }
    235 }
    236 
    237 int pcvM2N(poly m)
    238 {
    239   unsigned n=0,d=0;
    240   for(int i=0;i<pVariables;i++)
    241   {
    242     d+=pGetExp(m,i+1);
    243     n+=pcvIndex[i][d];
    244   }
    245   return n+1;
    246 }
    247 
    248 poly pcvN2M(int n)
    249 {
    250   n--;
    251   poly m=pOne();
    252   int i,j,k;
    253   for(i=pVariables-1;i>=0;i--)
    254   {
    255     k=j;
    256     for(j=0;j<pcvDegBound&& ((int)pcvIndex[i][j])<=n;j++);
    257     j--;
    258     n-=pcvIndex[i][j];
    259     if(i<pVariables-1) pSetExp(m,i+2,k-j);
    260   }
    261   if(n==0)
    262   {
    263     pSetExp(m,1,j);
    264     pSetm(m);
    265     return m;
    266   }
    267   else
    268   {
    269     pDelete1(&m);
    270     return NULL;
    271   }
    272 }
    273 
    274 poly pcvP2CV(poly p,int d0,int d1)
    275 {
    276   poly cv=NULL;
    277   while(p)
    278   {
    279     int d=pcvDeg(p);
    280     if(d0<=d&&d<d1)
    281     {
    282       poly c=pOne();
    283       pSetComp(c,pcvM2N(p));
    284       pSetCoeff(c,nCopy(pGetCoeff(p)));
    285       cv=pAdd(cv,c);
    286     }
    287     pIter(p);
    288   }
    289   return cv;
    290 }
    291 
    292 poly pcvCV2P(poly cv,int d0,int d1)
    293 {
    294   poly p=NULL;
    295   while(cv)
    296   {
    297     poly m=pcvN2M(pGetComp(cv));
    298     if(m)
    299     {
    300       int d=pcvDeg(m);
    301       if(d0<=d&&d<d1)
    302       {
    303         pSetCoeff(m,nCopy(pGetCoeff(cv)));
    304         p=pAdd(p,m);
    305       }
    306     }
    307     pIter(cv);
    308   }
    309   return p;
    310 }
    311 
    312 lists pcvP2CV(lists pl,int d0,int d1)
    313 {
    314   lists cvl=(lists)Alloc(sizeof(slists));
    315   cvl->Init(pl->nr+1);
    316   pcvInit(d1);
    317   for(int i=pl->nr;i>=0;i--)
    318   {
    319     if(pl->m[i].rtyp==POLY_CMD)
    320     {
    321       cvl->m[i].rtyp=VECTOR_CMD;
    322       cvl->m[i].data=pcvP2CV((poly)pl->m[i].data,d0,d1);
    323     }
    324   }
    325   pcvClean();
    326   return cvl;
    327 }
    328 
    329 lists pcvCV2P(lists cvl,int d0,int d1)
    330 {
    331   lists pl=(lists)Alloc(sizeof(slists));
    332   pl->Init(cvl->nr+1);
    333   pcvInit(d1);
    334   for(int i=cvl->nr;i>=0;i--)
    335   {
    336     if(cvl->m[i].rtyp==VECTOR_CMD)
    337     {
    338       pl->m[i].rtyp=POLY_CMD;
    339       pl->m[i].data=pcvCV2P((poly)cvl->m[i].data,d0,d1);
    340     }
    341   }
    342   pcvClean();
    343   return pl;
    344 }
    345 
    346 int pcvDim(int d0,int d1)
    347 {
    348   if(d0<0) d0=0;
    349   if(d1<0) d1=0;
    350   pcvInit(d1+1);
    351   int d=pcvIndex[pVariables-1][d1]-pcvIndex[pVariables-1][d0];
    352   pcvClean();
    353   return d;
    354 }
    355 
    356 int pcvBasis(lists b,int i,poly m,int d,int n)
    357 {
    358   if(n<pVariables)
    359   {
    360     for(int k=0,l=d;k<=l;k++,d--)
    361     {
    362       pSetExp(m,n,k);
    363       i=pcvBasis(b,i,m,d,n+1);
    364     }
    365   }
    366   else
    367   {
    368     pSetExp(m,n,d);
    369     pSetm(m);
    370     b->m[i].rtyp=POLY_CMD;
    371     b->m[i++].data=pCopy(m);
    372   }
    373   return i;
    374 }
    375 
    376 lists pcvBasis(int d0,int d1)
    377 {
    378   if(d0<0) d0=0;
    379   if(d1<0) d1=0;
    380   lists b=(lists)Alloc(sizeof(slists));
    381   b->Init(pcvDim(d0,d1));
    382   poly m=pOne();
    383   for(int d=d0,i=0;d<d1;d++)
    384     i=pcvBasis(b,i,m,d,1);
    385   pDelete1(&m);
    386   return b;
    387 }
    388 
    389429#endif /* !defined(HAVE_DYNAMIC_LOADING) || defined(BUILD_MODULE) */
    390430#endif /* HAVE_PCV */

  • Singular/pcv.h

    rbb21e3 r4d9dc5  
    22*  Computer Algebra System SINGULAR     *
    33****************************************/
    4 /* $Id: pcv.h,v 1.10 1998-12-21 10:54:59 mschulze Exp $ */
     4/* $Id: pcv.h,v 1.11 1999-06-08 09:14:01 mschulze Exp $ */
    55/*
    66* ABSTRACT: conversion between polys and coef vectors
    … …  
    1212int pcvDeg(poly p);
    1313int pcvMinDeg(poly p);
     14int pcvMinDeg(matrix m);
    1415int pcvMaxDeg(poly p);
     16int pcvMaxDeg(matrix m);
    1517BOOLEAN pcvMinDeg(leftv res,leftv h);
    1618BOOLEAN pcvMaxDeg(leftv res,leftv h);
Note: See TracChangeset for help on using the changeset viewer.