source: git/Singular/fglm.cc @ 18dd47

spielwiese
Last change on this file since 18dd47 was d720e3, checked in by Hans Schönemann <hannes@…>, 27 years ago
* hannes: corrected check for NOSTREAMIO in fglm.cc, fglmzero.cc git-svn-id: file:///usr/local/Singular/svn/trunk@594 2c84dea3-7e68-4137-9b89-c4e89433aadc
  • Property mode set to 100644
File size: 11.0 KB
Line 
1// emacs edit mode for this file is -*- C++ -*-
2// $Id: fglm.cc,v 1.9 1997-08-01 10:53:05 Singular Exp $
3
4/****************************************
5*  Computer Algebra System SINGULAR     *
6****************************************/
7/*
8* ABSTRACT - The FGLM-Algorithm
9*   Calculate a reduced groebner basis for one ordering, given a
10*   reduced groebner basis for another ordering.
11*   In this file the input is checked. Furthermore we decide, if
12*   the input is 0-dimensional ( then fglmzero.cc is used ) or
13*   if the input is homogeneous ( then fglmhom.cc is used. Yet
14*   not implemented ).
15*/
16
17#include "mod2.h"
18
19#ifndef NOSTREAMIO
20#include <iostream.h>
21#endif
22
23#ifdef HAVE_FGLM
24#include "tok.h"
25#include "structs.h"
26#include "polys.h"
27#include "ideals.h"
28#include "ring.h"
29#include "ipid.h"
30#include "ipshell.h"
31#include "febase.h"
32#include "maps.h"
33#include "mmemory.h"
34#include "kstd1.h"   
35#include "fglm.h"
36
37//     enumeration to handle the various errors to occour.
38enum FglmState{ 
39    FglmOk, 
40    FglmHasOne, 
41    FglmNoIdeal,
42    FglmNotReduced,
43    FglmNotZeroDim, 
44    FglmIncompatibleRings
45};
46
47// Has to be called, if currQuotient != NULL. ( i.e. qring-case )
48// Then a new ideal is build, consisting of the generators of sourceIdeal
49// and the generators of currQuotient, which are completely reduced by
50// the sourceIdeal. This means: If sourceIdeal is reduced, then the new
51// ideal will be reduced as well.
52// Assumes that currRing == sourceRing
53ideal fglmUpdatesource( const ideal sourceIdeal ) 
54{
55    int k, l, offset;
56    BOOLEAN found;
57    ideal newSource= idInit( IDELEMS( sourceIdeal ) + IDELEMS( currQuotient ), 1 );
58    for ( k= IDELEMS( sourceIdeal )-1; k >=0; k-- )
59        (newSource->m)[k]= pCopy( (sourceIdeal->m)[k] );
60    offset= IDELEMS( sourceIdeal );
61    for ( l= IDELEMS( currQuotient )-1; l >= 0; l-- ) {
62        if ( (currQuotient->m)[l] != NULL ) {
63            found= FALSE;
64            for ( k= IDELEMS( sourceIdeal )-1; (k >= 0) && (found == FALSE); k-- )
65                if ( pDivisibleBy( (sourceIdeal->m)[k], (currQuotient->m)[l] ) )
66                    found= TRUE;
67            if ( ! found ) {
68                (newSource->m)[offset]= pCopy( (currQuotient->m)[l] );
69                offset++;
70            }
71        }
72    }
73    idSkipZeroes( newSource );
74    return newSource;
75}
76
77// Has to be called, if currQuotient != NULL, i.e. in qring-case.
78// Gets rid of the elements of result which are contained in
79// currQuotient and skips Zeroes.
80// Assumes that currRing == destRing
81void
82fglmUpdateresult( ideal & result ) 
83{
84    int k, l;
85    BOOLEAN found;
86    for ( k= IDELEMS( result )-1; k >=0; k-- ) {
87        if ( (result->m)[k] != NULL ) {
88            found= FALSE;
89            for ( l= IDELEMS( currQuotient )-1; (l >= 0) && ( found == FALSE ); l-- ) 
90                if ( pDivisibleBy( (currQuotient->m)[l], (result->m)[k] ) )
91                    found= TRUE;
92            if ( found ) pDelete( & ((result->m)[k]) );
93        }
94    }
95    idSkipZeroes( result );
96}
97
98// Checks if the two rings sringHdl and dringHdl are compatible enough to
99// be used for the fglm. This means:
100//  1) Same Characteristic, 2) globalOrderings in both rings,
101//  3) Same number of variables, 4) same number of parameters
102//  5) variables in one ring are permutated variables of the other one
103//  6) parameters in one ring are permutated parameters of the other one
104//  7) either both rings are rings or both rings are qrings
105//  8) if they are qrings, the quotientIdeals of both must coincide.
106// vperm must be a vector of length pVariables+1, initialized by 0.
107// If both rings are compatible, it stores the permutation of the
108// variables if mapped from sringHdl to dringHdl.
109// if the rings are compatible, it returns FglmOk.
110// Should be called with currRing= IDRING( sringHdl );
111FglmState
112fglmConsistency( idhdl sringHdl, idhdl dringHdl, int * vperm ) 
113{
114    int k;
115    FglmState state= FglmOk;
116    ring dring = IDRING( dringHdl );
117    ring sring = IDRING( sringHdl );
118   
119    if ( sring->ch != dring->ch ) {
120        WerrorS( "rings must have same characteristic" );
121        state= FglmIncompatibleRings;
122    }
123    if ( (sring->OrdSgn != 1) || (dring->OrdSgn != 1) ) {
124        WerrorS( "only works for global orderings" );
125        state= FglmIncompatibleRings;
126    }
127    if ( sring->N != dring->N ) {
128        WerrorS( "rings must have same number of variables" );
129        state= FglmIncompatibleRings;
130    }
131    if ( sring->P != dring->P ) {
132        WerrorS( "rings must have same number of parameters" );
133        state= FglmIncompatibleRings;
134    }
135    if ( state != FglmOk ) return state;
136    // now the rings have the same number of variables resp. parameters.
137    // check if the names of the variables resp. parameters do agree:
138    int nvar = sring->N;
139    int npar = sring->P;
140    int * pperm;
141    if ( npar > 0 ) 
142        pperm= (int *)Alloc0( (npar+1)*sizeof( int ) );
143    else
144        pperm= NULL;
145    maFindPerm( sring->names, nvar, sring->parameter, npar, dring->names, nvar, dring->parameter, npar, vperm, pperm );
146    for ( k= nvar; (k > 0) && (state == FglmOk); k-- )
147        if ( vperm[k] <= 0 ) {
148            WerrorS( "variable names do not agree" );
149            state= FglmIncompatibleRings;
150        }
151    for ( k= npar-1; (k >= 0) && (state == FglmOk); k-- )
152        if ( pperm[k] >= 0 ) {
153            WerrorS( "paramater names do not agree" );
154            state= FglmIncompatibleRings;
155        }
156    Free( (ADDRESS)pperm, (npar+1)*sizeof( int ) );
157    if ( state != FglmOk ) return state;
158    // check if both rings are qrings or not
159    if ( sring->qideal != NULL ) {
160        if ( dring->qideal == NULL ) {
161            Werror( "%s is a qring, current ring not", sringHdl->id );
162            return FglmIncompatibleRings;
163        }
164        // both rings are qrings, now check if both quotients define the same ideal.
165        // check if sring->qideal is contained in dring->qideal:
166        rSetHdl( dringHdl, TRUE );
167        nSetMap( sring->ch, sring->parameter, npar, sring->minpoly );
168        ideal sqind = idInit( IDELEMS( sring->qideal ), 1 );
169        for ( k= IDELEMS( sring->qideal )-1; k >= 0; k-- )
170            (sqind->m)[k]= pPermPoly( (sring->qideal->m)[k], vperm, nvar );
171        ideal sqindred = kNF( dring->qideal, NULL, sqind );
172        if ( ! idIs0( sqindred ) ) {
173            WerrorS( "the quotients do not agree" );
174            state= FglmIncompatibleRings;
175        }
176        idDelete( & sqind );
177        idDelete( & sqindred );
178        rSetHdl( sringHdl, TRUE );
179        if ( state != FglmOk ) return state;
180        // check if dring->qideal is contained in sring->qideal:
181        int * dsvperm = (int *)Alloc0( (nvar+1)*sizeof( int ) );
182        maFindPerm( dring->names, nvar, NULL, 0, sring->names, nvar, NULL, 0, dsvperm, NULL );
183        nSetMap( dring->ch, dring->parameter, npar, dring->minpoly );
184        ideal dqins = idInit( IDELEMS( dring->qideal ), 1 );
185        for ( k= IDELEMS( dring->qideal )-1; k >= 0; k-- ) 
186            (dqins->m)[k]= pPermPoly( (dring->qideal->m)[k], dsvperm, nvar );
187        ideal dqinsred = kNF( sring->qideal, NULL, dqins );
188        if ( ! idIs0( dqinsred ) ) {
189            WerrorS( "the quotients do not agree" );
190            state= FglmIncompatibleRings;
191        }
192        idDelete( & dqins );
193        idDelete( & dqinsred );
194        Free( (ADDRESS)dsvperm, (nvar+1)*sizeof( int ) );
195        if ( state != FglmOk ) return state;
196    } 
197    else {
198        if ( dring->qideal != NULL ) {
199            Werror( "current ring is a qring, %s not", sringHdl->id );
200            return FglmIncompatibleRings;
201        }
202    }
203    return FglmOk;
204}
205
206//     Checks if the ideal "theIdeal" is zero-dimensional and minimal. It does
207//      not check, if it is reduced.
208//     returns FglmOk if we can use theIdeal for CalculateFunctionals (this
209//                     function reports an error if theIdeal is not reduced,
210//                     so this need not to be tested here)
211//             FglmNotReduced if theIdeal is not minimal
212//             FglmNotZeroDim if it is not zero-dimensional
213//             FglmHasOne if 1 belongs to theIdeal
214FglmState
215fglmIdealcheck( const ideal theIdeal )
216{
217    FglmState state = FglmOk;
218    int power;
219    int k; 
220    BOOLEAN * purePowers = (BOOLEAN *)Alloc( pVariables*sizeof( BOOLEAN ) );
221    for ( k= pVariables-1; k >= 0; k-- ) 
222        purePowers[k]= FALSE;
223
224    for ( k= IDELEMS( theIdeal ) - 1; (state == FglmOk) && (k >= 0); k-- ) {
225        poly p = (theIdeal->m)[k];
226        if ( pIsConstant( p ) ) state= FglmHasOne;
227        else if ( (power= pIsPurePower( p )) > 0 ) {
228            fglmASSERT( 0 < power && power <= pVariables, "illegal power" );
229            if ( purePowers[power-1] == TRUE  ) state= FglmNotReduced;
230            else purePowers[power-1]= TRUE; 
231        }
232        for ( int l = IDELEMS( theIdeal ) - 1; state == FglmOk && l >= 0; l-- ) 
233            if ( (k != l) && pDivisibleBy( p, (theIdeal->m)[l] ) )
234                state= FglmNotReduced;
235    }
236    if ( state == FglmOk ) {
237        for ( k= pVariables-1 ; (state == FglmOk) && (k >= 0); k-- ) 
238            if ( purePowers[k] == FALSE ) state= FglmNotZeroDim;
239    }
240    Free( (ADDRESS)purePowers, pVariables*sizeof( BOOLEAN ) );
241    return state;
242}
243
244//     the main function for the fglm-Algorithm.
245//     Checks the input-data, calls CalculateFunctionals, handles change
246//     of ring-vars and finaly calls GroebnerViaFunctionals.
247//     returns the new groebnerbasis or 0 if an error occoured.
248BOOLEAN
249fglmProc( leftv result, leftv first, leftv second ) 
250{
251    FglmState state = FglmOk;
252
253    idhdl destRingHdl = currRingHdl;
254    ring destRing = currRing;
255    ideal destIdeal = NULL;
256    idhdl sourceRingHdl = (idhdl)first->data;
257    rSetHdl( sourceRingHdl, TRUE );
258    ring sourceRing = currRing;
259
260    int * vperm = (int *)Alloc0( (pVariables+1)*sizeof( int ) );
261    state= fglmConsistency( sourceRingHdl, destRingHdl, vperm );
262    Free( (ADDRESS)vperm, (pVariables+1)*sizeof(int) );
263
264    if ( state == FglmOk ) {
265        idhdl ih = currRing->idroot->get( second->Name(), myynest );
266        if ( (ih != NULL) && (IDTYP(ih)==IDEAL_CMD) ) {
267            ideal sourceIdeal;
268            if ( currQuotient != NULL ) 
269                sourceIdeal= fglmUpdatesource( IDIDEAL( ih ) );
270            else
271                sourceIdeal = IDIDEAL( ih );
272            state= fglmIdealcheck( sourceIdeal );
273            if ( state == FglmOk ) { 
274                // Now the settings are compatible with FGLM
275                assumeStdFlag( (leftv)ih );
276                if ( fglmzero( sourceRingHdl, sourceIdeal, destRingHdl, destIdeal, FALSE, (currQuotient != NULL) ) == FALSE )
277                    state= FglmNotReduced;
278            }
279        } else state= FglmNoIdeal;
280    }
281    if ( currRingHdl != destRingHdl )
282        rSetHdl( destRingHdl, TRUE );
283    switch (state) {
284        case FglmOk:
285            if ( currQuotient != NULL ) fglmUpdateresult( destIdeal );
286            break;
287        case FglmHasOne:
288            destIdeal= idInit(1,1);
289            (destIdeal->m)[0]= pOne();
290            state= FglmOk;
291            break;
292        case FglmIncompatibleRings:
293            Werror( "ring %s and current ring are incompatible", first->Name() );
294            destIdeal= idInit(0,0);
295            break;
296        case FglmNoIdeal:
297            Werror( "Can't find ideal %s in ring %s", second->Name(), first->Name() );
298            destIdeal= idInit(0,0);
299            break;
300        case FglmNotZeroDim:
301            Werror( "The ideal %s has to be 0-dimensional", second->Name() );
302            destIdeal= idInit(0,0);
303            break;
304        case FglmNotReduced:
305            Werror( "The ideal %s has to be reduced", second->Name() );
306            destIdeal= idInit(0,0);
307            break;
308        default:
309            destIdeal= idInit(1,1);
310    }
311
312    result->rtyp = IDEAL_CMD;
313    result->data= (void *)destIdeal;
314    setFlag( result, FLAG_STD );
315    if ( state == FglmOk )
316        return FALSE;
317    else
318        return TRUE;
319}
320
321#endif
322// ----------------------------------------------------------------------------
323// Local Variables: ***
324// compile-command: "make Singular" ***
325// page-delimiter: "^\\(\\|//!\\)" ***
326// fold-internal-margins: nil ***
327// End: ***
Note: See TracBrowser for help on using the repository browser.