17 #include "factoryconf.h" 23 #include "kernel/subexpr.h" 27 #include "kernel/ipid.h" 28 #include "kernel/ipshell.h" 41 #define PROT(msg) if (BTEST1(OPT_PROT)) Print(msg) 42 #define STICKYPROT(msg) if (BTEST1(OPT_PROT)) Print(msg) 43 #define PROT2(msg,arg) if (BTEST1(OPT_PROT)) Print(msg,arg) 44 #define STICKYPROT2(msg,arg) if (BTEST1(OPT_PROT)) Print(msg,arg) 45 #define fglmASSERT(ignore1,ignore2) 63 homogElem() : v(), dv(), basis(0), destbasis(0), inDest(
FALSE) {}
70 #ifndef HAVE_EXPLICIT_CONSTR 91 doublepoly * sourceHeads;
101 int numberofdestbasismonoms;
107 hfglmNextdegree(
intvec * source, ideal current,
int & deg )
114 if ( deg < newhilb->
length() )
116 if ( deg < source->
length() )
117 numelems= (*newhilb)[deg]-(*source)[deg];
119 numelems= (*newhilb)[deg];
123 if (deg < source->
length())
124 numelems= -(*source)[deg];
139 generateMonoms( poly
m,
int var,
int deg, homogData * dat )
144 poly mon =
pCopy( m );
149 for ( i= dat->numSourceHeads - 1; (i >= 0) && (inSource==
FALSE); i-- ) {
154 for ( i= dat->numDestPolys - 1; (i >= 0) && (inDest==
FALSE); i-- ) {
159 if ( (!inSource) || (!inDest) ) {
162 basis= ++(dat->basisSize);
164 ++dat->numberofdestbasismonoms;
165 if ( dat->numMonoms == dat->monlistmax ) {
167 #ifdef HAVE_EXPLICIT_CONSTR 171 (dat->monlistmax)*
sizeof( homogElem ),
172 (dat->monlistmax+dat->monlistblock) *
sizeof( homogElem ) );
173 for ( k= dat->monlistmax; k < (dat->monlistmax+dat->monlistblock); k++ )
174 dat->monlist[k].homogElem();
177 int newsize = dat->monlistmax + dat->monlistblock;
178 homogElem * tempelem =
new homogElem[ newsize ];
180 for ( k= dat->monlistmax - 1; k >= 0; k-- )
181 tempelem[
k].initialize( dat->monlist[k] );
183 homogElem = tempelem;
185 dat->monlistmax+= dat->monlistblock;
187 #ifdef HAVE_EXPLICIT_CONSTR 188 dat->monlist[dat->numMonoms]= homogElem( mon, basis, inDest );
190 dat->monlist[dat->numMonoms].initialize( mon, basis, inDest );
193 if ( inSource && ! inDest )
PROT(
"\\" );
194 if ( ! inSource && inDest )
PROT(
"/" );
195 if ( ! inSource && ! inDest )
PROT(
"." );
203 poly newm =
pCopy( m );
205 generateMonoms( newm, var+1, deg, dat );
216 mapMonoms( ring oldRing, homogData & dat )
223 for ( s= dat.numMonoms - 1; s >= 0; s-- ) {
227 dat.monlist[
s].mon.sm= pPermPoly( dat.monlist[s].mon.dm, vperm, oldRing,
NULL, 0 );
232 getVectorRep( homogData & dat )
236 for ( s= 0; s < dat.numMonoms; s++ ) {
237 if ( dat.monlist[s].inDest ==
FALSE ) {
239 if ( dat.monlist[s].basis == 0 ) {
243 poly
nf =
kNF( dat.sourceIdeal,
NULL, dat.monlist[s].mon.sm );
246 while (temp !=
NULL ) {
248 for ( t= dat.numMonoms - 1; t >= 0; t-- ) {
249 if ( dat.monlist[t].basis > 0 ) {
250 if (
pLmEqual( dat.monlist[t].mon.sm, temp ) ) {
252 v.
setelem( dat.monlist[t].basis, coeff );
262 v=
fglmVector( dat.basisSize, dat.monlist[s].basis );
270 remapVectors( ring oldring, homogData & dat )
275 for ( s= dat.numMonoms - 1; s >= 0; s-- ) {
276 if ( dat.monlist[s].inDest ==
FALSE ) {
279 for ( k= dat.basisSize; k > 0; k-- ){
280 number newnum= nMap( dat.monlist[s].v.getelem( k ) );
281 newv.setelem( k, newnum );
283 dat.monlist[
s].dv= newv;
289 gaussreduce( homogData & dat,
int maxnum,
int BS )
294 int destbasisSize = 0;
297 for ( s= 0; (s < dat.numMonoms) && (found < maxnum); s++ ) {
298 if ( dat.monlist[s].inDest ==
FALSE ) {
299 if ( gauss.reduce( dat.monlist[s].dv ) ==
FALSE ) {
301 dat.monlist[
s].destbasis= destbasisSize;
311 for ( k= 1; k < p.
size(); k++ ) {
313 while ( dat.monlist[l].destbasis != k )
315 poly temp =
pCopy( dat.monlist[l].mon.dm );
317 result=
pAdd( result, temp );
324 (dat.destIdeal->m)[dat.numDestPolys]= result;
326 if (
IDELEMS(dat.destIdeal) == dat.numDestPolys ) {
336 PROT2(
"/%i)", dat.numberofdestbasismonoms );
341 fglmhomog( ring sourceRing, ideal sourceIdeal, ring destRing, ideal & destIdeal )
343 #define groebnerBS 16 354 dat.sourceIdeal= sourceIdeal;
355 dat.sourceHeads= (doublepoly *)
omAlloc(
IDELEMS( sourceIdeal ) *
sizeof( doublepoly ) );
356 for ( s=
IDELEMS( sourceIdeal ) - 1; s >= 0; s-- )
358 dat.sourceHeads[
s].sm=
pHead( (sourceIdeal->m)[s] );
360 dat.numSourceHeads=
IDELEMS( sourceIdeal );
364 int * vperm = (
int *)
omAlloc( (sourceRing->N + 1)*
sizeof(int) );
369 for ( s=
IDELEMS( sourceIdeal ) - 1; s >= 0; s-- )
371 dat.sourceHeads[
s].dm= pPermPoly( dat.sourceHeads[s].sm, vperm, sourceRing,
NULL, 0 );
374 dat.destIdeal=
idInit( groebnerBS, 1 );
377 while ( (numGBelems= hfglmNextdegree( hilb, dat.destIdeal, deg )) != 0 )
380 PROT2(
"deg= %i ", deg );
381 PROT2(
"num= %i\ngen>", numGBelems );
382 dat.monlistblock= 512;
383 dat.monlistmax= dat.monlistblock;
384 #ifdef HAVE_EXPLICIT_CONSTR 385 dat.monlist= (homogElem *)
omAlloc( dat.monlistmax*
sizeof( homogElem ) );
387 for ( j= dat.monlistmax - 1; j >= 0; j-- ) dat.monlist[
j].homogElem();
389 dat.monlist =
new homogElem[ dat.monlistmax ];
394 dat.numberofdestbasismonoms= 0;
397 generateMonoms( start, 1, deg, &dat );
400 PROT2(
"(%i/", dat.basisSize );
401 PROT2(
"%i)\nvec>", dat.overall );
404 mapMonoms( destRing, dat );
409 remapVectors( sourceRing, dat );
413 gaussreduce( dat, numGBelems, groebnerBS );
415 #ifdef HAVE_EXPLICIT_CONSTR 418 delete [] dat.monlist;
423 destIdeal= dat.destIdeal;
const CanonicalForm int s
poly kNF(ideal F, ideal Q, poly p, int syzComp, int lazyReduce)
Compatiblity layer for legacy polynomial operations (over currRing)
number getconstelem(int i) const
#define omFreeSize(addr, size)
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy ...
intvec * hHstdSeries(ideal S, intvec *modulweight, intvec *wdegree, ideal Q, ring tailRing)
VAR ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
#define omReallocSize(addr, o_size, size)
#define pHead(p)
returns newly allocated copy of Lm(p), coef is copied, next=NULL, p might be NULL ...
void idSkipZeroes(ideal ide)
gives an ideal/module the minimal possible size
void rChangeCurrRing(ring r)
ideal idInit(int idsize, int rank)
initialise an ideal / module
const Variable & v
< [in] a sqrfree bivariate poly
void maFindPerm(char const *const *const preim_names, int preim_n, char const *const *const preim_par, int preim_p, char const *const *const names, int n, char const *const *const par, int nop, int *perm, int *par_perm, n_coeffType ch)
#define pDivisibleBy(a, b)
returns TRUE, if leading monom of a divides leading monom of b i.e., if there exists a expvector c > ...
void pEnlargeSet(poly **p, int l, int increment)
void setelem(int i, number &n)
static BOOLEAN length(leftv result, leftv arg)
#define pSetCoeff(p, n)
deletes old coeff before setting the new one
#define pCopy(p)
return a copy of the poly