Changeset efaa52 in git


Ignore:
Timestamp:
Aug 13, 1999, 7:12:22 PM (25 years ago)
Author:
Hans Schönemann <hannes@…>
Branches:
(u'spielwiese', '17f1d200f27c5bd38f5dfc6e8a0879242279d1d8')
Children:
a848f8e428c6a9689038df1f3558643b6f843484
Parents:
578da371b36ce1fa90f41579b61226e98cc72a90
Message:
 * hannes: changes to spectrum stuff ( apdapt it to Singular convetions,
                                       use existing procedures)
          part 1 (to be continued)


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

Legend:

Unmodified
Added
Removed
  • Singular/iparith.cc

    r578da3 refaa52  
    22*  Computer Algebra System SINGULAR     *
    33****************************************/
    4 /* $Id: iparith.cc,v 1.168 1999-08-13 16:31:12 Singular Exp $ */
     4/* $Id: iparith.cc,v 1.169 1999-08-13 17:12:20 Singular Exp $ */
    55
    66/*
     
    281281  { "sortvec",     0, SORTVEC_CMD ,       CMD_1},
    282282  #ifdef HAVE_SPECTRUM
    283   { "spectrum",    0, SPECTRUM_CMD ,      CMD_1},
    284   { "spectrumnd",  0, SPECTRUMF_CMD ,     CMD_1},
    285   { "spadd",       0, SPADD_CMD ,         CMD_2},
    286   { "spmul",       0, SPMUL_CMD ,         CMD_2},
     283  { "spectrum",    0, SPECTRUM_CMD ,      CMD_123},
    287284  { "semic",       0, SEMIC_CMD ,         CMD_2},
    288285  { "semicsqh",    0, SEMICH_CMD ,        CMD_2},
     
    22862283,{jjSIMPL_ID,  SIMPLIFY_CMD,   MODUL_CMD,      MODUL_CMD,  INT_CMD PROFILER}
    22872284#ifdef HAVE_SPECTRUM
    2288 ,{spaddProc,   SPADD_CMD,      LIST_CMD,       LIST_CMD,   LIST_CMD PROFILER}
    2289 ,{spmulProc,   SPMUL_CMD,      LIST_CMD,       LIST_CMD,   INT_CMD  PROFILER}
     2285,{spectrumProc2,SPECTRUM_CMD,  LIST_CMD,       POLY_CMD,   INT_CMD PROFILER}
    22902286,{semicProc,   SEMIC_CMD,      INT_CMD,        LIST_CMD,   LIST_CMD PROFILER}
    22912287,{semichProc,  SEMICH_CMD,     INT_CMD,        LIST_CMD,   LIST_CMD PROFILER}
     
    36053601#ifdef HAVE_SPECTRUM
    36063602,{spectrumProc, SPECTRUM_CMD,    LIST_CMD,       POLY_CMD }
    3607 ,{spectrumfProc,SPECTRUMF_CMD,   LIST_CMD,       POLY_CMD }
    36083603#endif
    36093604,{jjSTD,        STD_CMD,         IDEAL_CMD,      IDEAL_CMD }
     
    43264321#else
    43274322,{jjWRONG3,         RESULTANT_CMD, POLY_CMD,POLY_CMD,   POLY_CMD,   POLY_CMD }
     4323#endif
     4324#ifdef HAVE_SPECTRUM
     4325,{spectrumOp3,      SPECTRUM_CMD, LIST_CMD, LIST_CMD,   STRING_CMD, LIST_CMD }
    43284326#endif
    43294327#ifdef OLD_RES
  • Singular/mmisc.c

    r578da3 refaa52  
    22*  Computer Algebra System SINGULAR     *
    33****************************************/
    4 /* $Id: mmisc.c,v 1.7 1999-08-06 14:06:41 obachman Exp $ */
     4/* $Id: mmisc.c,v 1.8 1999-08-13 17:12:21 Singular Exp $ */
    55
    66/*
     
    1313#include "mmpage.h"
    1414#include "febase.h"
    15 #include "version.h"
    1615#ifdef MTRACK
    1716#include "mmbt.h"
  • Singular/polys1.cc

    r578da3 refaa52  
    22*  Computer Algebra System SINGULAR     *
    33****************************************/
    4 /* $Id: polys1.cc,v 1.23 1999-07-23 13:45:17 Singular Exp $ */
     4/* $Id: polys1.cc,v 1.24 1999-08-13 17:12:21 Singular Exp $ */
    55
    66/*
     
    6161
    6262/*2
    63 *test if a monomial is a pure power
     63*test if a monomial /head term is a pure power
    6464*/
    6565int pIsPurePower(poly p)
     
    11161116  }
    11171117  result=pOrdPolyMerge(result);
    1118 #else 
     1118#else
    11191119  //  if (qq!=NULL)
    11201120  //  {
  • Singular/spectrum.cc

    r578da3 refaa52  
    3939
    4040// ----------------------------------------------------------------------------
     41//  test if the polynomial  h  has a term of total degree d
     42// ----------------------------------------------------------------------------
     43
     44BOOLEAN hasTermOfDegree( poly h, int d )
     45{
     46  do
     47  {
     48    if( pTotaldegree( h )== d )
     49      return  TRUE;
     50    pIter(h);
     51  }
     52  while( h!=(poly)NULL );
     53
     54  return  FALSE;
     55}
     56
     57// ----------------------------------------------------------------------------
    4158//  test if the polynomial  h  has a constant term
    4259// ----------------------------------------------------------------------------
    4360
    44 BOOLEAN hasConstTerm( poly h )
    45 {
    46     do
    47     {
    48         if( pTotaldegree( h )==0 )
     61static BOOLEAN inline hasConstTerm( poly h )
     62{
     63  return  hasTermOfDegree(h,0);
     64}
     65
     66// ----------------------------------------------------------------------------
     67//  test if the polynomial  h  has a linear term
     68// ----------------------------------------------------------------------------
     69
     70static BOOLEAN inline hasLinearTerm( poly h )
     71{
     72  return  hasTermOfDegree(h,1);
     73}
     74
     75// ----------------------------------------------------------------------------
     76//  test if the ideal  J  has a lead monomial on the axis number  k
     77// ----------------------------------------------------------------------------
     78
     79BOOLEAN hasAxis( ideal J,int k )
     80{
     81  int i;
     82
     83  for( i=0; i<IDELEMS(J); i++ )
     84  {
     85    if (pIsPurePower( J->m[i] ) == k) return TRUE;
     86  }
     87  return  FALSE;
     88}
     89
     90// ----------------------------------------------------------------------------
     91//  test if the ideal  J  has  1  as a lead monomial
     92// ----------------------------------------------------------------------------
     93
     94int     hasOne( ideal J )
     95{
     96  int i;
     97
     98  for( i=0; i<IDELEMS(J); i++ )
     99  {
     100    if (pIsConstant(J->m[i])) return TRUE;
     101  }
     102  return  FALSE;
     103}
     104// ----------------------------------------------------------------------------
     105//  test if  m  is a multiple of one of the monomials of  f
     106// ----------------------------------------------------------------------------
     107
     108int     isMultiple( poly f,poly m )
     109{
     110  while( f!=(poly)NULL )
     111  {
     112    // ---------------------------------------------------
     113    //  for a local order  f|m  is only possible if  f>=m
     114    // ---------------------------------------------------
     115
     116    if( pComp0( f,m )>=0 )
     117    {
     118      if( pDivisibleBy2( f,m ) )
     119      {
     120        return  TRUE;
     121      }
     122      else
     123      {
     124        pIter( f );
     125      }
     126    }
     127    else
     128    {
     129      return FALSE;
     130    }
     131  }
     132
     133  return  FALSE;
     134}
     135
     136// ----------------------------------------------------------------------------
     137//  compute the minimal monomial of minimmal  weight>=max_weight
     138// ----------------------------------------------------------------------------
     139
     140poly    computeWC( const newtonPolygon &np,Rational max_weight )
     141{
     142  poly m  = pOne();
     143  poly wc = (poly)NULL;
     144  int  mdegree;
     145
     146  for( int i=1; i<=pVariables; i++ )
     147  {
     148    mdegree = 1;
     149    pSetExp( m,i,mdegree );
     150    // pSetm( m );
     151    // np.weight_shift does not need pSetm( m ), postpone it
     152
     153    while(  np.weight_shift( m )<max_weight )
     154    {
     155      mdegree++;
     156      pSetExp( m,i,mdegree );
     157      // pSetm( m );
     158      // np.weight_shift does not need pSetm( m ), postpone it
     159    }
     160    pSetm( m );
     161
     162    if( i==1 || pComp( m,wc )<0 )
     163    {
     164      pDelete( &wc );
     165      wc = pHead( m );
     166    }
     167
     168    pSetExp( m,i,0 );
     169  }
     170
     171  pDelete( &m );
     172
     173  return  wc;
     174}
     175
     176// ----------------------------------------------------------------------------
     177//  deletes all monomials of  f  which are less than  hc
     178// ----------------------------------------------------------------------------
     179
     180static inline  poly    normalFormHC( poly f,poly hc )
     181{
     182  poly    *ptr = &f;
     183
     184  while( (*ptr)!=(poly)NULL )
     185  {
     186    if( pComp0( *ptr,hc )>=0 )
     187    {
     188      ptr = &(pNext( *ptr ));
     189    }
     190    else
     191    {
     192      pDelete( ptr );
     193      return  f;
     194    }
     195  }
     196
     197  return f;
     198}
     199
     200// ----------------------------------------------------------------------------
     201//  deletes all monomials of  f  which are multiples of monomials of  Z
     202// ----------------------------------------------------------------------------
     203
     204static inline  poly    normalFormZ( poly f,poly Z )
     205{
     206  poly    *ptr = &f;
     207
     208  while( (*ptr)!=(poly)NULL )
     209  {
     210    if( !isMultiple( Z,*ptr ) )
     211    {
     212      ptr = &(pNext( *ptr ));
     213    }
     214    else
     215    {
     216      pDelete1( ptr );
     217    }
     218  }
     219
     220  return f;
     221}
     222
     223// ?? HS:
     224// Looks for the shortest polynomial f in stdJ which is divisible
     225// by the monimial m, return its index in stdJ
     226// ----------------------------------------------------------------------------
     227//  Looks for the first polynomial f in stdJ which satisfies m=LT(f)*eta
     228//  for a monomial eta. The return eta*f-m and cancel all monomials
     229//  which are smaller than the highest corner hc
     230// ----------------------------------------------------------------------------
     231
     232static inline  int     isLeadMonomial( poly m,ideal stdJ )
     233{
     234  int     length = INT_MAX;
     235  int     result = -1;
     236
     237  for( int i=0; i<IDELEMS(stdJ); i++ )
     238  {
     239    if( pComp( stdJ->m[i],m )>=0 && pDivisibleBy( stdJ->m[i],m ) )
     240    {
     241      int     tmp = pLength( stdJ->m[i] );
     242
     243      if( tmp < length )
     244      {
     245        length = tmp;
     246        result = i;
     247      }
     248    }
     249  }
     250
     251  return  result;
     252}
     253
     254// ----------------------------------------------------------------------------
     255//  set the exponent of a monomial t an integer array
     256// ----------------------------------------------------------------------------
     257
     258static void    setExp( poly m,int *r )
     259{
     260  for( int i=pVariables; i>0; i-- )
     261  {
     262    pSetExp( m,i,r[i-1] );
     263  }
     264  pSetm( m );
     265}
     266
     267// ----------------------------------------------------------------------------
     268//  check if the ordering is a reverse wellordering, i.e. every monomial
     269//  is smaller than only finitely monomials
     270// ----------------------------------------------------------------------------
     271
     272static BOOLEAN isWell( void )
     273{
     274  int b = rBlocks( currRing );
     275
     276  if( b==3 &&
     277      ( currRing->order[0] == ringorder_ds ||
     278        currRing->order[0] == ringorder_Ds ||
     279        currRing->order[0] == ringorder_ws ||
     280        currRing->order[0] == ringorder_Ws ) )
     281  {
     282    return  TRUE;
     283  }
     284  else if( b>=3
     285  && (( currRing->order[0] ==ringorder_a
     286        && currRing->block1[0]==pVariables  )
     287    || (currRing->order[0]==ringorder_M
     288        && currRing->block1[0]==pVariables*pVariables )))
     289  {
     290    for( int i=pVariables-1; i>=0; i-- )
     291    {
     292      if( currRing->wvhdl[0][i]>=0 )
     293      {
     294        return  FALSE;
     295      }
     296    }
     297    return  TRUE;
     298  }
     299
     300  return  FALSE;
     301}
     302
     303// ----------------------------------------------------------------------------
     304//  compute all monomials not in  stdJ  and their normal forms
     305// ----------------------------------------------------------------------------
     306
     307static void    computeNF( ideal stdJ,poly hc,poly wc,spectrumPolyList *NF )
     308{
     309  int         carry,k;
     310  multiCnt    C( pVariables,0 );
     311  poly        Z = (poly)NULL;
     312
     313  int         well = isWell( );
     314
     315  do
     316  {
     317    poly    m = pOne();
     318    setExp( m,C.cnt );
     319
     320    carry = FALSE;
     321
     322    k = isLeadMonomial( m,stdJ );
     323
     324    if( k < 0 )
     325    {
     326      // ---------------------------
     327      //  m  is not a lead monomial
     328      // ---------------------------
     329
     330      NF->insert_node( m,(poly)NULL );
     331    }
     332    else if( isMultiple( Z,m ) )
     333    {
     334      // ------------------------------------
     335      //  m  is trivially in the ideal  stdJ
     336      // ------------------------------------
     337
     338      pDelete( &m );
     339      carry = TRUE;
     340    }
     341    else if( pComp( m,hc ) < 0 || pComp( m,wc ) < 0 )
     342    {
     343      // -------------------
     344      //  we do not need  m
     345      // -------------------
     346
     347      pDelete( &m );
     348      carry = TRUE;
     349    }
     350    else
     351    {
     352      // --------------------------
     353      //  compute lazy normal form
     354      // --------------------------
     355
     356      poly    multiplicant = pDivide( m,stdJ->m[k] );
     357      pGetCoeff( multiplicant ) = nInit(1);
     358
     359      poly    nf = pMultT( pCopy( stdJ->m[k] ), multiplicant );
     360
     361      pDelete( &multiplicant );
     362
     363      nf = normalFormHC( nf,hc );
     364
     365      if( pNext( nf )==(poly)NULL )
     366      {
     367        // ----------------------------------
     368        //  normal form of  m  is  m  itself
     369        // ----------------------------------
     370
     371        pDelete( &nf );
     372        NF->delete_monomial( m );
     373        Z = pAdd( Z,m );
     374        carry = TRUE;
     375      }
     376      else
     377      {
     378        nf = normalFormZ( nf,Z );
     379
     380        if( pNext( nf )==(poly)NULL )
    49381        {
    50             return  TRUE;
     382          // ----------------------------------
     383          //  normal form of  m  is  m  itself
     384          // ----------------------------------
     385
     386          pDelete( &nf );
     387          NF->delete_monomial( m );
     388          Z = pAdd( Z,m );
     389          carry = TRUE;
    51390        }
    52 
    53         h = pNext( h );
    54     }
    55     while( h!=(poly)NULL );
    56 
    57     return  FALSE;
    58 }
    59 
    60 // ----------------------------------------------------------------------------
    61 //  test if the polynomial  h  has a linear term
    62 // ----------------------------------------------------------------------------
    63 
    64 BOOLEAN hasLinearTerm( poly h )
    65 {
    66     do
    67     {
    68         if( pTotaldegree( h )==1 )
     391        else
    69392        {
    70             return  TRUE;
    71         }
    72 
    73         h = pNext( h );
    74     }
    75     while( h!=(poly)NULL );
    76 
    77     return  FALSE;
    78 }
    79 
    80 // ----------------------------------------------------------------------------
    81 //  test if the ideal  J  has a lead monomial on the axis number  k
    82 // ----------------------------------------------------------------------------
    83 
    84 BOOLEAN hasAxis( ideal J,int k )
    85 {
    86     int i,j,found;
    87 
    88     for( i=0; i<IDELEMS(J); i++ )
    89     {
    90         found = TRUE;
    91 
    92         for( j=1; j<k && found==TRUE; j++ )
    93         {
    94             if( pGetExp( J->m[i],j ) > 0 )
     393          // ------------------------------------
     394          //  normal form of  m  is a polynomial
     395          // ------------------------------------
     396
     397          pNorm( nf );
     398          if( well==TRUE )
     399          {
     400            NF->insert_node( m,nf );
     401          }
     402          else
     403          {
     404            poly    nfhard = kNF( stdJ,(ideal)NULL,pNext( nf ),0,0 );
     405            nfhard = normalFormHC( nfhard,hc );
     406            nfhard = normalFormZ ( nfhard,Z );
     407
     408            if( nfhard==(poly)NULL )
    95409            {
    96                 found = FALSE;
    97             }
    98         }
    99 
    100         if( pGetExp( J->m[i],k )==0 )
    101         {
    102             found = FALSE;
    103         }
    104 
    105         for( j=k+1; j<=pVariables && found==TRUE; j++ )
    106         {
    107             if( pGetExp( J->m[i],j ) > 0 )
    108             {
    109                 found = FALSE;
    110             }
    111         }
    112 
    113         if( found==TRUE )
    114         {
    115             return  TRUE;
    116         }
    117     }
    118     return  FALSE;
    119 }
    120 
    121 // ----------------------------------------------------------------------------
    122 //  test if the ideal  J  has  1  as a lead monomial
    123 // ----------------------------------------------------------------------------
    124 
    125 int     hasOne( ideal J )
    126 {
    127     int i,j,found;
    128 
    129     for( i=0; i<IDELEMS(J); i++ )
    130     {
    131         found = TRUE;
    132 
    133         for( j=1; j<=pVariables && found==TRUE; j++ )
    134         {
    135             if( pGetExp( J->m[i],j ) > 0 )
    136             {
    137                 found = FALSE;
    138             }
    139         }
    140 
    141         if( found==TRUE )
    142         {
    143             return  TRUE;
    144         }
    145     }
    146     return  FALSE;
    147 }
    148 // ----------------------------------------------------------------------------
    149 //  test if  m  is a multiple of one of the monomials of  f
    150 // ----------------------------------------------------------------------------
    151 
    152 int     isMultiple( poly f,poly m )
    153 {
    154     while( f!=(poly)NULL )
    155     {
    156         // ---------------------------------------------------
    157         //  for a local order  f|m  is only possible if  f>=m
    158         // ---------------------------------------------------
    159 
    160         if( pComp( f,m )>=0 )
    161         {
    162             if( pDivisibleBy2( f,m ) )
    163             {
    164                 return  TRUE;
     410              NF->delete_monomial( m );
     411              Z = pAdd( Z,m );
     412              carry = TRUE;
    165413            }
    166414            else
    167415            {
    168                 f = pNext( f );
     416              pDelete( &pNext( nf ) );
     417              pNext( nf ) = nfhard;
     418              NF->insert_node( m,nf );
    169419            }
     420          }
    170421        }
    171         else
    172         {
    173             return FALSE;
    174         }
    175     }
    176 
    177     return  FALSE;
    178 }
    179 
    180 // ----------------------------------------------------------------------------
    181 //  compute the minimal monomial of minimmal  weight>=max_weight
    182 // ----------------------------------------------------------------------------
    183 
    184 poly    computeWC( const newtonPolygon &np,Rational max_weight )
    185 {
    186     poly m  = pISet( 1 );
    187     poly wc = (poly)NULL;
    188 
    189     int  mdegree;
    190 
    191     for( int i=1; i<=pVariables; i++ )
    192     {
    193         mdegree = 1;
    194         pSetExp( m,i,mdegree );
    195         pSetm( m );
    196 
    197         while(  np.weight_shift( m )<max_weight )
    198         {
    199             mdegree++;
    200             pSetExp( m,i,mdegree );
    201             pSetm( m );
    202         }
    203 
    204         if( i==1 || pComp( m,wc )<0 )
    205         {
    206             pDelete( &wc );
    207             wc = pCopy( m );
    208         }
    209 
    210         pSetExp( m,i,0 );
    211     }
    212 
    213     pDelete( &m );
    214 
    215     return  wc;
    216 }
    217 
    218 // ----------------------------------------------------------------------------
    219 //  deletes all monomials of  f  which are less than  hc
    220 // ----------------------------------------------------------------------------
    221 
    222 inline  poly    normalFormHC( poly f,poly hc )
    223 {
    224     poly    *ptr = &f;
    225 
    226     while( (*ptr)!=(poly)NULL )
    227     {
    228         if( pComp( *ptr,hc )>=0 )
    229         {
    230             ptr = &(pNext( *ptr ));
    231         }
    232         else
    233         {
    234             pDelete( ptr );
    235             return  f;
    236         }
    237     }
    238 
    239     return f;
    240 }
    241 
    242 // ----------------------------------------------------------------------------
    243 //  deletes all monomials of  f  which are multiples of monomials of  Z
    244 // ----------------------------------------------------------------------------
    245 
    246 inline  poly    normalFormZ( poly f,poly Z )
    247 {
    248     poly    *ptr = &f;
    249 
    250     while( (*ptr)!=(poly)NULL )
    251     {
    252         if( !isMultiple( Z,*ptr ) )
    253         {
    254             ptr = &(pNext( *ptr ));
    255         }
    256         else
    257         {
    258             pDelete1( ptr );
    259         }
    260     }
    261 
    262     return f;
    263 }
    264 
    265 // ----------------------------------------------------------------------------
    266 //  Looks for the first polynomial f in stdJ which satisfies m=LT(f)*eta
    267 //  for a monomial eta. The return eta*f-m and cancel all monomials
    268 //  which are smaller than the highest corner hc
    269 // ----------------------------------------------------------------------------
    270 
    271 inline  int     isLeadMonomial( poly m,ideal stdJ )
    272 {
    273     int     length = INT_MAX;
    274     int     result = -1;
    275 
    276     for( int i=0; i<IDELEMS(stdJ); i++ )
    277     {
    278         if( pComp( stdJ->m[i],m )>=0 && pDivisibleBy( stdJ->m[i],m ) )
    279         {
    280             int     tmp = pLength( stdJ->m[i] );
    281 
    282             if( tmp < length )
    283             {
    284                 length = tmp;
    285                 result = i;
    286             }
    287         }
    288     }
    289 
    290     return  result;
    291 }
    292 
    293 // ----------------------------------------------------------------------------
    294 //  set the exponent of a monomial t an integer array
    295 // ----------------------------------------------------------------------------
    296 
    297 void    setExp( poly m,int *r )
    298 {
    299     for( int i=0; i<pVariables; i++ )
    300     {
    301         pSetExp( m,i+1,r[i] );
    302     }
    303 }
    304 
    305 // ----------------------------------------------------------------------------
    306 //  check if the ordering is a reverse wellordering, i.e. every monomial
    307 //  is smaller than only finitely monomials
    308 // ----------------------------------------------------------------------------
    309 
    310 BOOLEAN isWell( void )
    311 {
    312     int b = rBlocks( currRing );
    313 
    314     if( b==3 &&
    315         ( currRing->order[0] == ringorder_ds ||
    316           currRing->order[0] == ringorder_Ds ||
    317           currRing->order[0] == ringorder_ws ||
    318           currRing->order[0] == ringorder_Ws ) )
    319     {
    320         return  TRUE;
    321     }
    322     else if( b>=3 && currRing->order[0] ==ringorder_a
    323                   && currRing->block1[0]==pVariables  )
    324     {
    325         for( int i=0; i<pVariables; i++ )
    326         {
    327             if( currRing->wvhdl[0][i]>=0 )
    328             {
    329                 return  FALSE;
    330                 }
    331         }
    332         return  TRUE;
    333     }
    334     else if( b>=3 && currRing->order[0]==ringorder_M
    335                   && currRing->block1[0]==pVariables*pVariables )
    336     {
    337         for( int i=0; i<pVariables; i++ )
    338         {
    339             if( currRing->wvhdl[0][i]>=0 )
    340             {
    341                 return  FALSE;
    342                 }
    343         }
    344         return  TRUE;
    345     }
    346 
    347     return  FALSE;
    348 }
    349 
    350 // ----------------------------------------------------------------------------
    351 //  compute all monomials not in  stdJ  and their normal forms
    352 // ----------------------------------------------------------------------------
    353 
    354 void    computeNF( ideal stdJ,poly hc,poly wc,spectrumPolyList *NF )
    355 {
    356     int         carry,k;
    357     multiCnt    C( pVariables,0 );
    358     poly        Z = (poly)NULL;
    359 
    360     int         well = isWell( );
    361 
    362     do
    363     {
    364         poly    m = pISet( 1 );
    365         setExp( m,C.cnt );
    366         pSetm( m );
    367 
    368         carry = FALSE;
    369 
    370         k = isLeadMonomial( m,stdJ );
    371 
    372         if( k < 0 )
    373         {
    374             // ---------------------------
    375             //  m  is not a lead monomial
    376             // ---------------------------
    377 
    378             NF->insert_node( m,(poly)NULL );
    379         }
    380         else if( isMultiple( Z,m ) )
    381         {
    382             // ------------------------------------
    383             //  m  is trivially in the ideal  stdJ
    384             // ------------------------------------
    385 
    386             pDelete( &m );
    387             carry = TRUE;
    388         }
    389         else if( pComp( m,hc ) < 0 || pComp( m,wc ) < 0 )
    390         {
    391             // -------------------
    392             //  we do not need  m
    393             // -------------------
    394 
    395             pDelete( &m );
    396             carry = TRUE;
    397         }
    398         else
    399         {
    400             // --------------------------
    401             //  compute lazy normal form
    402             // --------------------------
    403 
    404             poly    multiplicant = pDivide( m,stdJ->m[k] );
    405             pGetCoeff( multiplicant ) = nInit(1);
    406 
    407             poly    nf = pMultT( pCopy( stdJ->m[k] ), multiplicant );
    408 
    409             pDelete( &multiplicant );
    410 
    411             nf = normalFormHC( nf,hc );
    412 
    413             if( pNext( nf )==(poly)NULL )
    414             {
    415                 // ----------------------------------
    416                 //  normal form of  m  is  m  itself
    417                 // ----------------------------------
    418 
    419                 pDelete( &nf );
    420                 NF->delete_monomial( m );
    421                 Z = pAdd( Z,m );
    422                 carry = TRUE;
    423             }
    424             else
    425             {
    426                 nf = normalFormZ( nf,Z );
    427 
    428                 if( pNext( nf )==(poly)NULL )
    429                 {
    430                     // ----------------------------------
    431                     //  normal form of  m  is  m  itself
    432                     // ----------------------------------
    433 
    434                     pDelete( &nf );
    435                     NF->delete_monomial( m );
    436                     Z = pAdd( Z,m );
    437                     carry = TRUE;
    438                 }
    439                 else
    440                 {
    441                       // ------------------------------------
    442                     //  normal form of  m  is a polynomial
    443                     // ------------------------------------
    444 
    445                     if( well==TRUE )
    446                     {
    447                         pNorm( nf );
    448                          NF->insert_node( m,nf );
    449                     }
    450                     else
    451                     {
    452                         pNorm( nf );
    453                         poly    nfhard = kNF( stdJ,(ideal)NULL,
    454                                               pNext( nf ),0,0 );
    455                         nfhard = normalFormHC( nfhard,hc );
    456                         nfhard = normalFormZ ( nfhard,Z );
    457 
    458                         if( nfhard==(poly)NULL )
    459                         {
    460                             NF->delete_monomial( m );
    461                             Z = pAdd( Z,m );
    462                             carry = TRUE;
    463                         }
    464                         else
    465                         {
    466                             pDelete( &pNext( nf ) );
    467                             pNext( nf ) = nfhard;
    468                             NF->insert_node( m,nf );
    469                         }
    470                     }
    471                 }
    472             }
    473         }
    474     }
    475     while( C.inc( carry ) );
    476 
    477     // delete single monomials
    478 
    479     int  finished = FALSE;
    480 
    481     while( finished==FALSE )
    482     {
    483         finished = TRUE;
    484 
    485         spectrumPolyNode  *node = NF->root;
    486 
    487         while( node!=(spectrumPolyNode*)NULL )
    488         {
    489             if( node->nf!=(poly)NULL && pNext( node->nf )==(poly)NULL )
    490             {
    491                 NF->delete_monomial( node->nf );
    492                 finished = FALSE;
    493                 node = (spectrumPolyNode*)NULL;
    494             }
    495             else
    496             {
    497                 node = node->next;
    498             }
    499         }
    500     }
    501 
    502     pDelete( &Z );
     422      }
     423    }
     424  }
     425  while( C.inc( carry ) );
     426
     427  // delete single monomials
     428
     429  BOOLEAN  not_finished;
     430
     431  do
     432  {
     433    not_finished = FALSE;
     434
     435    spectrumPolyNode  *node = NF->root;
     436
     437    while( node!=(spectrumPolyNode*)NULL )
     438    {
     439      if( node->nf!=(poly)NULL && pNext( node->nf )==(poly)NULL )
     440      {
     441        NF->delete_monomial( node->nf );
     442        not_finished = TRUE;
     443        node = (spectrumPolyNode*)NULL;
     444      }
     445      else
     446      {
     447        node = node->next;
     448      }
     449    }
     450  } while( not_finished );
     451
     452  pDelete( &Z );
    503453}
    504454
     
    507457// ----------------------------------------------------------------------------
    508458
    509 spectrumState   spectrumCompute( poly h,lists *L,int fast )
    510 {
    511     int i,j;
     459static spectrumState   spectrumCompute( poly h,lists *L,int fast )
     460{
     461  int i,j;
     462
     463  #ifdef SPECTRUM_DEBUG
     464  #ifdef SPECTRUM_PRINT
     465  #ifdef SPECTRUM_IOSTREAM
     466    cout << "spectrumCompute\n";
     467    if( fast==0 ) cout << "    no optimization" << endl;
     468    if( fast==1 ) cout << "    weight optimization" << endl;
     469    if( fast==2 ) cout << "    symmetry optimization" << endl;
     470  #else
     471    fprintf( stdout,"spectrumCompute\n" );
     472    if( fast==0 ) fprintf( stdout,"    no optimization\n" );
     473    if( fast==1 ) fprintf( stdout,"    weight optimization\n" );
     474    if( fast==2 ) fprintf( stdout,"    symmetry optimization\n" );
     475  #endif
     476  #endif
     477  #endif
     478
     479  // ----------------------
     480  //  check if  h  is zero
     481  // ----------------------
     482
     483  if( h==(poly)NULL )
     484  {
     485    return  spectrumZero;
     486  }
     487
     488  // ----------------------------------
     489  //  check if  h  has a constant term
     490  // ----------------------------------
     491
     492  if( hasConstTerm( h ) )
     493  {
     494    return  spectrumBadPoly;
     495  }
     496
     497  // --------------------------------
     498  //  check if  h  has a linear term
     499  // --------------------------------
     500
     501  if( hasLinearTerm( h ) )
     502  {
     503    *L = (lists)Alloc( sizeof( slists ) );
     504    (*L)->Init( 1 );
     505    (*L)->m[0].rtyp = INT_CMD;    //  milnor number
     506    /* (*L)->m[0].data = (void*)0;a  -- done by Init */
     507
     508    return  spectrumNoSingularity;
     509  }
     510
     511  // ----------------------------------
     512  //  compute the jacobi ideal of  (h)
     513  // ----------------------------------
     514
     515  ideal J = NULL;
     516  J = idInit( pVariables,1 );
     517
     518  #ifdef SPECTRUM_DEBUG
     519  #ifdef SPECTRUM_PRINT
     520  #ifdef SPECTRUM_IOSTREAM
     521    cout << "\n   computing the Jacobi ideal...\n";
     522  #else
     523    fprintf( stdout,"\n   computing the Jacobi ideal...\n" );
     524  #endif
     525  #endif
     526  #endif
     527
     528  for( i=0; i<pVariables; i++ )
     529  {
     530    J->m[i] = pDiff( h,i+1); //j );
    512531
    513532    #ifdef SPECTRUM_DEBUG
    514533    #ifdef SPECTRUM_PRINT
    515534    #ifdef SPECTRUM_IOSTREAM
    516         cout << "spectrumCompute" << endl;
    517         if( fast==0 ) cout << "    no optimization" << endl;
    518         if( fast==1 ) cout << "    weight optimization" << endl;
    519         if( fast==2 ) cout << "    symmetry optimization" << endl;
     535      cout << "        ";
    520536    #else
    521         fprintf( stdout,"spectrumCompute\n" );
    522         if( fast==0 ) fprintf( stdout,"    no optimization\n" );
    523         if( fast==1 ) fprintf( stdout,"    weight optimization\n" );
    524         if( fast==2 ) fprintf( stdout,"    symmetry optimization\n" );
     537      fprintf( stdout,"        " );
     538    #endif
     539      pWrite( J->m[i] );
    525540    #endif
    526541    #endif
    527     #endif
    528 
    529     // ----------------------
    530     //  check if  h  is zero
    531     // ----------------------
    532 
    533     if( h==(poly)NULL )
    534     {
    535         return  spectrumZero;
    536     }
    537 
    538     // ----------------------------------
    539     //  check if  h  has a constant term
    540     // ----------------------------------
    541 
    542     if( hasConstTerm( h )==TRUE )
    543     {
    544         return  spectrumBadPoly;
    545     }
    546 
    547     // --------------------------------
    548     //  check if  h  has a linear term
    549     // --------------------------------
    550 
    551     if( hasLinearTerm( h )==TRUE )
    552     {
    553         *L = (lists)Alloc( sizeof( slists ) );
    554         (*L)->Init( 1 );
    555         (*L)->m[0].rtyp = INT_CMD;    //  milnor number
    556         (*L)->m[0].data = (void*)0;
    557 
    558         return  spectrumNoSingularity;
    559     }
    560 
    561     // ----------------------------------
    562     //  compute the jacobi ideal of  (h)
    563     // ----------------------------------
    564 
    565     ideal J = NULL;
    566     J = idInit( pVariables,1 );
    567 
    568     #ifdef SPECTRUM_DEBUG
    569     #ifdef SPECTRUM_PRINT
    570     #ifdef SPECTRUM_IOSTREAM
    571         cout << endl;
    572         cout << "   computing the Jacobi ideal..." << endl;
    573     #else
    574         fprintf( stdout,"\n" );
    575         fprintf( stdout,"   computing the Jacobi ideal...\n" );
    576     #endif
    577     #endif
    578     #endif
    579 
    580     for( i=0,j=1; i<pVariables; i++,j++ )
    581     {
    582         J->m[i] = pDiff( h,j );
    583 
    584         #ifdef SPECTRUM_DEBUG
    585         #ifdef SPECTRUM_PRINT
    586         #ifdef SPECTRUM_IOSTREAM
    587             cout << "        ";
    588         #else
    589             fprintf( stdout,"        " );
    590         #endif
    591 
    592             pWrite( J->m[i] );
    593         #endif
    594         #endif
    595     }
    596 
    597     // --------------------------------------------
    598     //  compute a standard basis  stdJ  of  jac(h)
    599     // --------------------------------------------
    600 
    601     #ifdef SPECTRUM_DEBUG
    602     #ifdef SPECTRUM_PRINT
    603     #ifdef SPECTRUM_IOSTREAM
    604         cout << endl;
    605         cout << "    computing a standard basis..." << endl;
    606     #else
    607         fprintf( stdout,"\n" );
    608         fprintf( stdout,"    computing a standard basis...\n" );
    609     #endif
    610     #endif
    611     #endif
    612 
    613     ideal stdJ = NULL;
    614     stdJ = kStd(J,currQuotient,isNotHomog,NULL);
    615     idSkipZeroes( stdJ );
    616 
    617     #ifdef SPECTRUM_DEBUG
    618     #ifdef SPECTRUM_PRINT
    619         for( i=0; i<IDELEMS(stdJ); i++ )
    620         {
    621             #ifdef SPECTRUM_IOSTREAM
    622                 cout << "        ";
    623             #else
    624                 fprintf( stdout,"        " );
    625             #endif
    626 
    627             pWrite( stdJ->m[i] );
    628         }
    629     #endif
    630     #endif
    631 
    632     idDelete( &J );
    633 
    634     // ------------------------------------------
    635     //  check if the  h  has a singularity
    636     // ------------------------------------------
    637 
    638     if( hasOne( stdJ )==TRUE )
    639     {
    640         // -------------------------------
    641         //  h is smooth in the origin
    642         //  return only the Milnor number
    643         // -------------------------------
    644 
    645         *L = (lists)Alloc( sizeof( slists ) );
    646         (*L)->Init( 1 );
    647         (*L)->m[0].rtyp = INT_CMD;    //  milnor number
    648         (*L)->m[0].data = (void*)0;
    649 
    650         return  spectrumNoSingularity;
    651     }
    652 
    653     // ------------------------------------------
    654     //  check if the singularity  h  is isolated
    655     // ------------------------------------------
     542  }
     543
     544  // --------------------------------------------
     545  //  compute a standard basis  stdJ  of  jac(h)
     546  // --------------------------------------------
     547
     548  #ifdef SPECTRUM_DEBUG
     549  #ifdef SPECTRUM_PRINT
     550  #ifdef SPECTRUM_IOSTREAM
     551    cout << endl;
     552    cout << "    computing a standard basis..." << endl;
     553  #else
     554    fprintf( stdout,"\n" );
     555    fprintf( stdout,"    computing a standard basis...\n" );
     556  #endif
     557  #endif
     558  #endif
     559
     560  ideal stdJ = kStd(J,currQuotient,isNotHomog,NULL);
     561  idSkipZeroes( stdJ );
     562
     563  #ifdef SPECTRUM_DEBUG
     564  #ifdef SPECTRUM_PRINT
     565    for( i=0; i<IDELEMS(stdJ); i++ )
     566    {
     567      #ifdef SPECTRUM_IOSTREAM
     568        cout << "        ";
     569      #else
     570        fprintf( stdout,"        " );
     571      #endif
     572
     573      pWrite( stdJ->m[i] );
     574    }
     575  #endif
     576  #endif
     577
     578  idDelete( &J );
     579
     580  // ------------------------------------------
     581  //  check if the  h  has a singularity
     582  // ------------------------------------------
     583
     584  if( hasOne( stdJ ) )
     585  {
     586    // -------------------------------
     587    //  h is smooth in the origin
     588    //  return only the Milnor number
     589    // -------------------------------
     590
     591    *L = (lists)Alloc( sizeof( slists ) );
     592    (*L)->Init( 1 );
     593    (*L)->m[0].rtyp = INT_CMD;    //  milnor number
     594    /* (*L)->m[0].data = (void*)0;a  -- done by Init */
     595
     596    return  spectrumNoSingularity;
     597  }
     598
     599  // ------------------------------------------
     600  //  check if the singularity  h  is isolated
     601  // ------------------------------------------
     602
     603  for( i=pVariables; i>0; i-- )
     604  {
     605    if( hasAxis( stdJ,i )==FALSE )
     606    {
     607      return  spectrumNotIsolated;
     608    }
     609  }
     610
     611  // ------------------------------------------
     612  //  compute the highest corner  hc  of  stdJ
     613  // ------------------------------------------
     614
     615  #ifdef SPECTRUM_DEBUG
     616  #ifdef SPECTRUM_PRINT
     617  #ifdef SPECTRUM_IOSTREAM
     618    cout << "\n    computing the highest corner...\n";
     619  #else
     620    fprintf( stdout,"\n    computing the highest corner...\n" );
     621  #endif
     622  #endif
     623  #endif
     624
     625  poly hc = (poly)NULL;
     626
     627  scComputeHC( stdJ,0,hc );
     628
     629  if( hc!=(poly)NULL )
     630  {
     631    pGetCoeff(hc) = nInit(1);
    656632
    657633    for( i=pVariables; i>0; i-- )
    658634    {
    659         if( hasAxis( stdJ,i )==FALSE )
    660         {
    661             return  spectrumNotIsolated;
    662         }
    663     }
    664 
    665     // ------------------------------------------
    666     //  compute the highest corner  hc  of  stdJ
    667     // ------------------------------------------
    668 
    669     #ifdef SPECTRUM_DEBUG
    670     #ifdef SPECTRUM_PRINT
    671     #ifdef SPECTRUM_IOSTREAM
    672         cout << endl;
    673         cout << "    computing the highest corner..."  << endl;
    674     #else
    675         fprintf( stdout,"\n" );
    676         fprintf( stdout,"    computing the highest corner...\n" );
    677     #endif
    678     #endif
    679     #endif
    680 
    681     poly hc = (poly)NULL;
    682 
    683     scComputeHC( stdJ,0,hc );
    684 
    685     if( hc!=(poly)NULL )
    686     {
    687         pGetCoeff(hc) = nInit(1);
    688 
    689         for( i=pVariables; i>0; i-- )
    690         {
    691             if( pGetExp( hc,i )>0 ) pDecrExp( hc,i );
    692         }
    693         pSetm( hc );
    694     }
    695     else
    696     {
    697         return  spectrumNoHC;
    698     }
    699 
    700     #ifdef SPECTRUM_DEBUG
    701     #ifdef SPECTRUM_PRINT
    702     #ifdef SPECTRUM_IOSTREAM
    703         cout << "       ";
    704     #else
    705         fprintf( stdout,"       " );
    706     #endif
    707         pWrite( hc );
    708     #endif
    709     #endif
    710 
    711     // ----------------------------------------
    712     //  compute the Newton polygon  nph  of  h
    713     // ----------------------------------------
    714 
    715     #ifdef SPECTRUM_DEBUG
    716     #ifdef SPECTRUM_PRINT
    717     #ifdef SPECTRUM_IOSTREAM
    718         cout << endl;
    719         cout << "    computing the newton polygon..." << endl;
    720     #else
    721         fprintf( stdout,"\n" );
    722         fprintf( stdout,"    computing the newton polygon...\n" );
    723     #endif
    724     #endif
    725     #endif
    726 
    727     newtonPolygon nph( h );
    728 
    729     #ifdef SPECTRUM_DEBUG
    730     #ifdef SPECTRUM_PRINT
    731         cout << nph;
    732     #endif
    733     #endif
    734 
    735     // -----------------------------------------------
    736     //  compute the weight corner  wc  of  (stdj,nph)
    737     // -----------------------------------------------
    738 
    739     #ifdef SPECTRUM_DEBUG
    740     #ifdef SPECTRUM_PRINT
    741     #ifdef SPECTRUM_IOSTREAM
    742         cout << endl;
    743         cout << "    computing the weight corner..." << endl;
    744     #else
    745         fprintf( stdout,"\n" );
    746         fprintf( stdout,"    computing the weight corner...\n" );
    747     #endif
    748     #endif
    749     #endif
    750 
    751     poly    wc = ( fast==0 ? pCopy( hc ) :
    752                  ( fast==1 ? computeWC( nph,(Rational)pVariables ) :
    753                 /* fast==2 */computeWC( nph,((Rational)pVariables)/(Rational)2 ) ) );
    754 
    755     #ifdef SPECTRUM_DEBUG
    756     #ifdef SPECTRUM_PRINT
    757     #ifdef SPECTRUM_IOSTREAM
    758         cout << "        ";
    759     #else
    760         fprintf( stdout,"        " );
    761     #endif
    762         pWrite( wc );
    763     #endif
    764     #endif
    765 
    766     // -------------
    767     //  compute  NF
    768     // -------------
    769 
    770     #ifdef SPECTRUM_DEBUG
    771     #ifdef SPECTRUM_PRINT
    772     #ifdef SPECTRUM_IOSTREAM
    773         cout << endl;
    774         cout << "    computing NF..." << endl;
    775     #else
    776         fprintf( stdout,"\n" );
    777         fprintf( stdout,"    computing NF...\n" );
    778     #endif
    779     #endif
    780     #endif
    781 
    782     spectrumPolyList NF( &nph );
    783 
    784     computeNF( stdJ,hc,wc,&NF );
    785 
    786     #ifdef SPECTRUM_DEBUG
    787     #ifdef SPECTRUM_PRINT
    788         cout << NF;
    789     #ifdef SPECTRUM_IOSTREAM
    790         cout << endl;
    791     #else
    792         fprintf( stdout,"n" );
    793     #endif
    794     #endif
    795     #endif
    796 
    797     // ----------------------------
    798     //  compute the spectrum of  h
    799     // ----------------------------
    800 
    801     return  NF.spectrum( L,fast );
     635      if( pGetExp( hc,i )>0 ) pDecrExp( hc,i );
     636    }
     637    pSetm( hc );
     638  }
     639  else
     640  {
     641    return  spectrumNoHC;
     642  }
     643
     644  #ifdef SPECTRUM_DEBUG
     645  #ifdef SPECTRUM_PRINT
     646  #ifdef SPECTRUM_IOSTREAM
     647    cout << "       ";
     648  #else
     649    fprintf( stdout,"       " );
     650  #endif
     651    pWrite( hc );
     652  #endif
     653  #endif
     654
     655  // ----------------------------------------
     656  //  compute the Newton polygon  nph  of  h
     657  // ----------------------------------------
     658
     659  #ifdef SPECTRUM_DEBUG
     660  #ifdef SPECTRUM_PRINT
     661  #ifdef SPECTRUM_IOSTREAM
     662    cout << "\n    computing the newton polygon...\n";
     663  #else
     664    fprintf( stdout,"\n    computing the newton polygon...\n" );
     665  #endif
     666  #endif
     667  #endif
     668
     669  newtonPolygon nph( h );
     670
     671  #ifdef SPECTRUM_DEBUG
     672  #ifdef SPECTRUM_PRINT
     673    cout << nph;
     674  #endif
     675  #endif
     676
     677  // -----------------------------------------------
     678  //  compute the weight corner  wc  of  (stdj,nph)
     679  // -----------------------------------------------
     680
     681  #ifdef SPECTRUM_DEBUG
     682  #ifdef SPECTRUM_PRINT
     683  #ifdef SPECTRUM_IOSTREAM
     684    cout << "\n    computing the weight corner...\n";
     685  #else
     686    fprintf( stdout,"\n    computing the weight corner...\n" );
     687  #endif
     688  #endif
     689  #endif
     690
     691  poly    wc = ( fast==0 ? pCopy( hc ) :
     692               ( fast==1 ? computeWC( nph,(Rational)pVariables ) :
     693              /* fast==2 */computeWC( nph,((Rational)pVariables)/(Rational)2 ) ) );
     694
     695  #ifdef SPECTRUM_DEBUG
     696  #ifdef SPECTRUM_PRINT
     697  #ifdef SPECTRUM_IOSTREAM
     698    cout << "        ";
     699  #else
     700    fprintf( stdout,"        " );
     701  #endif
     702    pWrite( wc );
     703  #endif
     704  #endif
     705
     706  // -------------
     707  //  compute  NF
     708  // -------------
     709
     710  #ifdef SPECTRUM_DEBUG
     711  #ifdef SPECTRUM_PRINT
     712  #ifdef SPECTRUM_IOSTREAM
     713    cout << "\n    computing NF...\n" << endl;
     714  #else
     715    fprintf( stdout,"\n    computing NF...\n" );
     716  #endif
     717  #endif
     718  #endif
     719
     720  spectrumPolyList NF( &nph );
     721
     722  computeNF( stdJ,hc,wc,&NF );
     723
     724  #ifdef SPECTRUM_DEBUG
     725  #ifdef SPECTRUM_PRINT
     726    cout << NF;
     727  #ifdef SPECTRUM_IOSTREAM
     728    cout << endl;
     729  #else
     730    fprintf( stdout,"\n" );
     731  #endif
     732  #endif
     733  #endif
     734
     735  // ----------------------------
     736  //  compute the spectrum of  h
     737  // ----------------------------
     738
     739  return  NF.spectrum( L,fast );
    802740}
    803741
     
    806744// ----------------------------------------------------------------------------
    807745
    808 BOOLEAN ringIsLocal( void )
    809 {
    810     poly    m   = pISet( 1 );
    811     poly    one = pISet( 1 );
    812 
    813     for( int i=1; i<=pVariables; i++ )
    814     {
    815         pSetExp( m,i,1 );
    816         pSetm( m );
    817 
    818         if( pComp( m,one )>0 )
    819         {
    820             return  FALSE;
    821         }
    822         pSetExp( m,i,0 );
    823     }
    824 
    825     pDelete( &m );
    826     pDelete( &one );
    827 
    828     return  TRUE;
    829 }
    830 
     746static BOOLEAN ringIsLocal( void )
     747{
     748  poly    m   = pOne();
     749  poly    one = pOne();
     750
     751  for( int i=pVariables; i>0; i-- )
     752  {
     753    pSetExp( m,i,1 );
     754    pSetm( m );
     755
     756    if( pComp( m,one )>0 )
     757    {
     758      return  FALSE;
     759    }
     760    pSetExp( m,i,0 );
     761  }
     762
     763  pDelete( &m );
     764  pDelete( &one );
     765
     766  return  TRUE;
     767}
     768
     769// ----------------------------------------------------------------------------
     770// print error message corresponding to spectrumState state:
     771// ----------------------------------------------------------------------------
     772static void spectrumPrintError(spectrumState state)
     773{
     774  switch( state )
     775  {
     776    case spectrumZero:
     777      WerrorS( "polynomial is zero" );
     778      break;
     779    case spectrumBadPoly:
     780      WerrorS( "polynomial has constant term" );
     781      break;
     782    case spectrumNoSingularity:
     783      WerrorS( "not a singularity" );
     784      break;
     785    case spectrumNotIsolated:
     786      WerrorS( "the singularity is not isolated" );
     787      break;
     788    case spectrumNoHC:
     789      WerrorS( "highest corner cannot be computed" );
     790      break;
     791    case spectrumDegenerate:
     792      WerrorS( "principal part is degenerate" );
     793      break;
     794    case spectrumOK:
     795      break;
     796
     797    default:
     798      WerrorS( "unknown error occurred" );
     799      break;
     800  }
     801}
    831802// ----------------------------------------------------------------------------
    832803//  this procedure is called from the interpreter
     
    838809BOOLEAN spectrumProc( leftv result,leftv first )
    839810{
    840     spectrumState state = spectrumOK;
    841 
    842     // -------------------
    843     //  check consistency
    844     // -------------------
    845 
    846     //  check for a local ring
    847 
    848     if( !ringIsLocal( ) )
    849     {
    850         WerrorS( "only works for local orderings" );
    851         state = spectrumWrongRing;
    852     }
    853 
    854     //  no quotient rings are allowed
    855 
    856     else if( currRing->qideal != NULL )
    857     {
    858         WerrorS( "does not work in quotient rings" );
    859         state = spectrumWrongRing;
     811  spectrumState state = spectrumOK;
     812
     813  // -------------------
     814  //  check consistency
     815  // -------------------
     816
     817  //  check for a local ring
     818
     819  if( !ringIsLocal( ) )
     820  {
     821    WerrorS( "only works for local orderings" );
     822    state = spectrumWrongRing;
     823  }
     824
     825  //  no quotient rings are allowed
     826
     827  else if( currRing->qideal != NULL )
     828  {
     829    WerrorS( "does not work in quotient rings" );
     830    state = spectrumWrongRing;
     831  }
     832  else
     833  {
     834    lists   L    = (lists)NULL;
     835    int     flag = 1; // weight corner optimization is safe
     836
     837    state = spectrumCompute( (poly)first->Data( ),&L,flag );
     838
     839    if( state==spectrumOK )
     840    {
     841      result->rtyp = LIST_CMD;
     842      result->data = (char*)L;
    860843    }
    861844    else
    862845    {
    863         lists   L    = (lists)NULL;
    864         int     flag = 1; // weight corner optimization is safe
    865 
    866         state = spectrumCompute( (poly)first->Data( ),&L,flag );
    867 
    868         if( state==spectrumOK )
    869         {
    870             result->rtyp = LIST_CMD;
    871             result->data = (char*)L;
    872         }
    873         else
    874         {
    875             switch( state )
    876             {
    877                    case spectrumZero:
    878                     WerrorS( "polynomial is zero" );
    879                     break;
    880                    case spectrumBadPoly:
    881                     WerrorS( "polynomial has constant term" );
    882                     break;
    883                    case spectrumNoSingularity:
    884                     WerrorS( "not a singularity" );
    885                     break;
    886                    case spectrumNotIsolated:
    887                     WerrorS( "the singularity is not isolated" );
    888                     break;
    889                    case spectrumNoHC:
    890                     WerrorS( "highest corner cannot be computed" );
    891                     break;
    892                 case spectrumDegenerate:
    893                     WerrorS( "principal part is degenerate" );
    894                     break;
    895 
    896                 default:
    897                     WerrorS( "unknown error occurred" );
    898                     break;
    899             }
    900         }
    901     }
    902 
    903     return  (state!=spectrumOK);
     846      spectrumPrintError(state);
     847    }
     848  }
     849
     850  return  (state!=spectrumOK);
    904851}
    905852
     
    913860BOOLEAN spectrumfProc( leftv result,leftv first )
    914861{
    915     spectrumState state = spectrumOK;
    916 
    917     // -------------------
    918     //  check consistency
    919     // -------------------
    920 
    921     //  check for a local polynomial ring
    922 
    923     if( currRing->OrdSgn != -1 )
    924     {
    925         WerrorS( "only works for local orderings" );
    926         state = spectrumWrongRing;
    927     }
    928     else if( currRing->qideal != NULL )
    929     {
    930         WerrorS( "does not work in quotient rings" );
    931         state = spectrumWrongRing;
     862  spectrumState state = spectrumOK;
     863
     864  // -------------------
     865  //  check consistency
     866  // -------------------
     867
     868  //  check for a local polynomial ring
     869
     870  if( currRing->OrdSgn != -1 )
     871  // ?? HS: the test above is also true for k[x][[y]], k[[x]][y]
     872  //if( !ringIsLocal( ) )
     873  {
     874    WerrorS( "only works for local orderings" );
     875    state = spectrumWrongRing;
     876  }
     877  else if( currRing->qideal != NULL )
     878  {
     879    WerrorS( "does not work in quotient rings" );
     880    state = spectrumWrongRing;
     881  }
     882  else
     883  {
     884    lists   L    = (lists)NULL;
     885    int     flag = 2; // symmetric optimization
     886
     887    state = spectrumCompute( (poly)first->Data( ),&L,flag );
     888
     889    if( state==spectrumOK )
     890    {
     891      result->rtyp = LIST_CMD;
     892      result->data = (char*)L;
    932893    }
    933894    else
    934895    {
    935         lists   L    = (lists)NULL;
    936         int     flag = 2; // symmetric optimization
    937 
    938         state = spectrumCompute( (poly)first->Data( ),&L,flag );
    939 
    940         if( state==spectrumOK )
    941         {
    942             result->rtyp = LIST_CMD;
    943             result->data = (char*)L;
    944         }
    945         else
    946         {
    947             switch( state )
    948             {
    949                    case spectrumZero:
    950                     WerrorS( "polynomial is zero" );
    951                     break;
    952                    case spectrumBadPoly:
    953                     WerrorS( "polynomial has constant term" );
    954                     break;
    955                    case spectrumNoSingularity:
    956                     WerrorS( "not a singularity" );
    957                     break;
    958                    case spectrumNotIsolated:
    959                     WerrorS( "the singularity is not isolated" );
    960                     break;
    961                    case spectrumNoHC:
    962                     WerrorS( "highest corner cannot be computed" );
    963                     break;
    964                 case spectrumDegenerate:
    965                     WerrorS( "principal part is degenerate" );
    966                     break;
    967 
    968                 default:
    969                     WerrorS( "unknown error occurred" );
    970                     break;
    971             }
    972         }
    973     }
    974 
    975     return  (state!=spectrumOK);
     896      spectrumPrintError(state);
     897    }
     898  }
     899
     900  return  (state!=spectrumOK);
    976901}
    977902
     
    14241349}
    14251350
     1351BOOLEAN    spectrumProc2 ( leftv res,leftv u, leftv v)
     1352{
     1353  if (((int)v->Data())==1)
     1354    return spectrumfProc(res,u);
     1355  return spectrumProc(res,u);
     1356}
     1357
     1358BOOLEAN    spectrumOp3  ( leftv res, leftv u, leftv v, leftv w )
     1359{
     1360  char v_str=(char *)v->Data();
     1361  if (strcmp(v_str, "+")==0)
     1362    return spaddProc(res,u,w);
     1363  else if (strcmp(v_str, "*")==0)
     1364    return spmulProc(res,u,w);
     1365  Werror("unknown operation '%s' for spectrum",v_str);
     1366  return TRUE;
     1367}
    14261368#endif /* HAVE_SPECTRUM */
    14271369// ----------------------------------------------------------------------------
  • Singular/spectrum.h

    r578da3 refaa52  
    99#define SPECTRUM_H
    1010
     11BOOLEAN    spectrumProc2 ( leftv,leftv, leftv );
     12BOOLEAN    spectrumOp3  ( leftv,leftv, leftv, leftv );
    1113BOOLEAN    spectrumProc ( leftv,leftv );
    1214BOOLEAN    spectrumfProc( leftv,leftv );
  • Singular/tok.h

    r578da3 refaa52  
    77* ABSTRACT: tokens, types for interpreter; general macros
    88*/
    9 /* $Id: tok.h,v 1.25 1999-08-12 10:57:37 Singular Exp $ */
     9/* $Id: tok.h,v 1.26 1999-08-13 17:12:22 Singular Exp $ */
    1010
    1111#ifndef MYYSTYPE
     
    105105#ifdef HAVE_SPECTRUM
    106106  SPECTRUM_CMD,
    107   SPECTRUMF_CMD,
    108   SPADD_CMD,
    109   SPMUL_CMD,
    110107  SEMIC_CMD,
    111108  SEMICH_CMD,
Note: See TracChangeset for help on using the changeset viewer.