17 #include "factoryconf.h"
23 #include "kernel/subexpr.h"
27 #include "kernel/ipid.h"
28 #include "kernel/ipshell.h"
42 #define PROT(msg) if (BTEST1(OPT_PROT)) Print(msg)
43 #define STICKYPROT(msg) if (BTEST1(OPT_PROT)) Print(msg)
44 #define PROT2(msg,arg) if (BTEST1(OPT_PROT)) Print(msg,arg)
45 #define STICKYPROT2(msg,arg) if (BTEST1(OPT_PROT)) Print(msg,arg)
46 #define fglmASSERT(ignore1,ignore2)
64 homogElem() :
v(), dv(), basis(0), destbasis(0), inDest(
FALSE) {}
71 #ifndef HAVE_EXPLICIT_CONSTR
92 doublepoly * sourceHeads;
102 int numberofdestbasismonoms;
108 hfglmNextdegree(
intvec * source, ideal current,
int & deg )
115 if ( deg < newhilb->
length() )
117 if ( deg < source->
length() )
118 numelems= (*newhilb)[deg]-(*source)[deg];
120 numelems= (*newhilb)[deg];
124 if (deg < source->
length())
125 numelems= -(*source)[deg];
140 generateMonoms( poly
m,
int var,
int deg, homogData * dat )
150 for (
i= dat->numSourceHeads - 1; (
i >= 0) && (inSource==
FALSE);
i-- ) {
155 for (
i= dat->numDestPolys - 1; (
i >= 0) && (inDest==
FALSE);
i-- ) {
160 if ( (!inSource) || (!inDest) ) {
163 basis= ++(dat->basisSize);
165 ++dat->numberofdestbasismonoms;
166 if ( dat->numMonoms == dat->monlistmax ) {
168 #ifdef HAVE_EXPLICIT_CONSTR
172 (dat->monlistmax)*
sizeof( homogElem ),
173 (dat->monlistmax+dat->monlistblock) *
sizeof( homogElem ) );
174 for (
k= dat->monlistmax;
k < (dat->monlistmax+dat->monlistblock);
k++ )
175 dat->monlist[
k].homogElem();
178 int newsize = dat->monlistmax + dat->monlistblock;
179 homogElem * tempelem =
new homogElem[ newsize ];
181 for (
k= dat->monlistmax - 1;
k >= 0;
k-- )
182 tempelem[
k].initialize( dat->monlist[
k] );
184 homogElem = tempelem;
186 dat->monlistmax+= dat->monlistblock;
188 #ifdef HAVE_EXPLICIT_CONSTR
189 dat->monlist[dat->numMonoms]= homogElem( mon, basis, inDest );
191 dat->monlist[dat->numMonoms].initialize( mon, basis, inDest );
194 if ( inSource && ! inDest )
PROT(
"\\" );
195 if ( ! inSource && inDest )
PROT(
"/" );
196 if ( ! inSource && ! inDest )
PROT(
"." );
206 generateMonoms( newm, var+1, deg, dat );
217 mapMonoms( ring oldRing, homogData & dat )
224 for (
s= dat.numMonoms - 1;
s >= 0;
s-- ) {
228 dat.monlist[
s].mon.sm= pPermPoly( dat.monlist[
s].mon.dm, vperm, oldRing,
NULL, 0 );
233 getVectorRep( homogData & dat )
237 for (
s= 0;
s < dat.numMonoms;
s++ ) {
238 if ( dat.monlist[
s].inDest ==
FALSE ) {
240 if ( dat.monlist[
s].basis == 0 ) {
244 poly
nf =
kNF( dat.sourceIdeal,
NULL, dat.monlist[
s].mon.sm );
247 while (temp !=
NULL ) {
249 for ( t= dat.numMonoms - 1; t >= 0; t-- ) {
250 if ( dat.monlist[t].basis > 0 ) {
251 if (
pLmEqual( dat.monlist[t].mon.sm, temp ) ) {
253 v.setelem( dat.monlist[t].basis, coeff );
271 remapVectors( ring oldring, homogData & dat )
276 for (
s= dat.numMonoms - 1;
s >= 0;
s-- ) {
277 if ( dat.monlist[
s].inDest ==
FALSE ) {
280 for (
k= dat.basisSize;
k > 0;
k-- ){
281 number newnum= nMap( dat.monlist[
s].v.getelem(
k ) );
282 newv.setelem(
k, newnum );
284 dat.monlist[
s].dv= newv;
290 gaussreduce( homogData & dat,
int maxnum,
int BS )
295 int destbasisSize = 0;
298 for (
s= 0; (
s < dat.numMonoms) && (
found < maxnum);
s++ ) {
299 if ( dat.monlist[
s].inDest ==
FALSE ) {
300 if ( gauss.reduce( dat.monlist[
s].dv ) ==
FALSE ) {
302 dat.monlist[
s].destbasis= destbasisSize;
312 for (
k= 1;
k <
p.size();
k++ ) {
313 if ( !
p.elemIsZero(
k ) ) {
314 while ( dat.monlist[
l].destbasis !=
k )
316 poly temp =
pCopy( dat.monlist[
l].mon.dm );
325 (dat.destIdeal->m)[dat.numDestPolys]=
result;
327 if (
IDELEMS(dat.destIdeal) == dat.numDestPolys ) {
337 PROT2(
"/%i)", dat.numberofdestbasismonoms );
342 fglmhomog( ring sourceRing, ideal sourceIdeal, ring destRing, ideal & destIdeal )
344 #define groebnerBS 16
355 dat.sourceIdeal= sourceIdeal;
356 dat.sourceHeads= (doublepoly *)
omAlloc(
IDELEMS( sourceIdeal ) *
sizeof( doublepoly ) );
357 for (
s=
IDELEMS( sourceIdeal ) - 1;
s >= 0;
s-- )
359 dat.sourceHeads[
s].sm=
pHead( (sourceIdeal->m)[
s] );
361 dat.numSourceHeads=
IDELEMS( sourceIdeal );
365 int * vperm = (
int *)
omAlloc( (sourceRing->N + 1)*
sizeof(int) );
370 for (
s=
IDELEMS( sourceIdeal ) - 1;
s >= 0;
s-- )
372 dat.sourceHeads[
s].dm= pPermPoly( dat.sourceHeads[
s].sm, vperm, sourceRing,
NULL, 0 );
375 dat.destIdeal=
idInit( groebnerBS, 1 );
378 while ( (numGBelems= hfglmNextdegree( hilb, dat.destIdeal, deg )) != 0 )
381 PROT2(
"deg= %i ", deg );
382 PROT2(
"num= %i\ngen>", numGBelems );
383 dat.monlistblock= 512;
384 dat.monlistmax= dat.monlistblock;
385 #ifdef HAVE_EXPLICIT_CONSTR
386 dat.monlist= (homogElem *)
omAlloc( dat.monlistmax*
sizeof( homogElem ) );
388 for (
j= dat.monlistmax - 1;
j >= 0;
j-- ) dat.monlist[
j].homogElem();
390 dat.monlist =
new homogElem[ dat.monlistmax ];
395 dat.numberofdestbasismonoms= 0;
398 generateMonoms( start, 1, deg, &dat );
401 PROT2(
"(%i/", dat.basisSize );
402 PROT2(
"%i)\nvec>", dat.overall );
405 mapMonoms( destRing, dat );
410 remapVectors( sourceRing, dat );
414 gaussreduce( dat, numGBelems, groebnerBS );
416 #ifdef HAVE_EXPLICIT_CONSTR
419 delete [] dat.monlist;
424 destIdeal= dat.destIdeal;
const CanonicalForm int s
const Variable & v
< [in] a sqrfree bivariate poly
intvec * hHstdSeries(ideal S, intvec *modulweight, intvec *wdegree, ideal Q, ring tailRing)
static BOOLEAN length(leftv result, leftv arg)
poly kNF(ideal F, ideal Q, poly p, int syzComp, int lazyReduce)
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)
static number & pGetCoeff(poly p)
return an alias to the leading coefficient of p assumes that p != NULL NOTE: not copy
#define omFreeSize(addr, size)
#define omReallocSize(addr, o_size, size)
void pEnlargeSet(poly **p, int l, int increment)
ring currRing
Widely used global variable which specifies the current polynomial ring for Singular interpreter and ...
void rChangeCurrRing(ring r)
Compatiblity layer for legacy polynomial operations (over currRing)
#define pHead(p)
returns newly allocated copy of Lm(p), coef is copied, next=NULL, p might be NULL
#define pSetCoeff(p, n)
deletes old coeff before setting the new one
#define pDivisibleBy(a, b)
returns TRUE, if leading monom of a divides leading monom of b i.e., if there exists a expvector c > ...
#define pCopy(p)
return a copy of the poly
ideal idInit(int idsize, int rank)
initialise an ideal / module
void idSkipZeroes(ideal ide)
gives an ideal/module the minimal possible size