Context Navigation

source: git/Singular/ipshell.cc @ 999726

spielwiese

Last change on this file since 999726 was 999726, checked in by Hans Schoenemann <hannes@…>, 14 years ago
debug stuff git-svn-id: file:///usr/local/Singular/svn/trunk@12836 2c84dea3-7e68-4137-9b89-c4e89433aadc
Property mode set to `100644`
File size: 125.2 KB

Line
1	/****************************************
2	* Computer Algebra System SINGULAR *
3	****************************************/
4	/* $Id$ */
5	/*
6	* ABSTRACT:
7	*/
8
9	//#include <stdlib.h>
10	#include <stdio.h>
11	#include <string.h>
12	#include <ctype.h>
13	#include <math.h>
14
15	#include "mod2.h"
16	#include "tok.h"
17	#include "options.h"
18	#include "ipid.h"
19	#include "intvec.h"
20	#include "omalloc.h"
21	#include "febase.h"
22	#include "polys.h"
23	#include "prCopy.h"
24	#include "ideals.h"
25	#include "matpol.h"
26	#include "kstd1.h"
27	#include "ring.h"
28	#include "subexpr.h"
29	#include "maps.h"
30	#include "syz.h"
31	#include "numbers.h"
32	#include "modulop.h"
33	#include "longalg.h"
34	#include "lists.h"
35	#include "attrib.h"
36	#include "ipconv.h"
37	#include "silink.h"
38	#include "stairc.h"
39	#include "weight.h"
40	#include "semic.h"
41	#include "splist.h"
42	#include "spectrum.h"
43	#include "gnumpfl.h"
44	#include "mpr_base.h"
45	#include "ffields.h"
46	#include "clapsing.h"
47	#include "hutil.h"
48	#include "gring.h"
49	#include "ipshell.h"
50	#ifdef HAVE_FACTORY
51	#define SI_DONT_HAVE_GLOBAL_VARS
52	#include <factory.h>
53	#endif
54
55	// define this if you want to use the fast_map routine for mapping ideals
56	#define FAST_MAP
57
58	#ifdef FAST_MAP
59	#include "fast_maps.h"
60	#endif
61
62	leftv iiCurrArgs=NULL;
63	idhdl iiCurrProc=NULL;
64	int traceit = 0;
65	const char *lastreserved=NULL;
66
67	int myynest = -1;
68
69	static BOOLEAN iiNoKeepRing=TRUE;
70
71	/0 implementation/
72
73	const char * iiTwoOps(int t)
74	{
75	if (t<127)
76	{
77	static char ch[2];
78	switch (t)
79	{
80	case '&':
81	return "and";
82	case '\|':
83	return "or";
84	default:
85	ch[0]=t;
86	ch[1]='\0';
87	return ch;
88	}
89	}
90	switch (t)
91	{
92	case COLONCOLON: return "::";
93	case DOTDOT: return "..";
94	//case PLUSEQUAL: return "+=";
95	//case MINUSEQUAL: return "-=";
96	case MINUSMINUS: return "--";
97	case PLUSPLUS: return "++";
98	case EQUAL_EQUAL: return "==";
99	case LE: return "<=";
100	case GE: return ">=";
101	case NOTEQUAL: return "<>";
102	default: return Tok2Cmdname(t);
103	}
104	}
105
106	static void list1(const char* s, idhdl h,BOOLEAN c, BOOLEAN fullname)
107	{
108	char buffer[22];
109	int l;
110	char buf2[128];
111
112	if(fullname) sprintf(buf2, "%s::%s", "", IDID(h));
113	else sprintf(buf2, "%s", IDID(h));
114
115	Print("%s%-20.20s [%d] ",s,buf2,IDLEV(h));
116	if (h == currRingHdl) PrintS("*");
117	PrintS(Tok2Cmdname((int)IDTYP(h)));
118
119	ipListFlag(h);
120	switch(IDTYP(h))
121	{
122	case INT_CMD: Print(" %d",IDINT(h)); break;
123	case INTVEC_CMD:Print(" (%d)",IDINTVEC(h)->length()); break;
124	case INTMAT_CMD:Print(" %d x %d",IDINTVEC(h)->rows(),IDINTVEC(h)->cols());
125	break;
126	case POLY_CMD:
127	case VECTOR_CMD:if (c)
128	{
129	PrintS(" ");wrp(IDPOLY(h));
130	if(IDPOLY(h) != NULL)
131	{
132	Print(", %d monomial(s)",pLength(IDPOLY(h)));
133	}
134	}
135	break;
136	case MODUL_CMD: Print(", rk %d", (int)(IDIDEAL(h)->rank));
137	case IDEAL_CMD: Print(", %u generator(s)",
138	IDELEMS(IDIDEAL(h))); break;
139	case MAP_CMD:
140	Print(" from %s",IDMAP(h)->preimage); break;
141	case MATRIX_CMD:Print(" %u x %u"
142	,MATROWS(IDMATRIX(h))
143	,MATCOLS(IDMATRIX(h))
144	);
145	break;
146	case PACKAGE_CMD:
147	PrintS(" (");
148	switch (IDPACKAGE(h)->language)
149	{
150	case LANG_SINGULAR: PrintS("S"); break;
151	case LANG_C: PrintS("C"); break;
152	case LANG_TOP: PrintS("T"); break;
153	case LANG_NONE: PrintS("N"); break;
154	default: PrintS("U");
155	}
156	if(IDPACKAGE(h)->libname!=NULL)
157	Print(",%s", IDPACKAGE(h)->libname);
158	PrintS(")");
159	break;
160	case PROC_CMD: if(strlen(IDPROC(h)->libname)>0)
161	Print(" from %s",IDPROC(h)->libname);
162	if(IDPROC(h)->is_static)
163	PrintS(" (static)");
164	break;
165	case STRING_CMD:
166	{
167	char *s;
168	l=strlen(IDSTRING(h));
169	memset(buffer,0,22);
170	strncpy(buffer,IDSTRING(h),si_min(l,20));
171	if ((s=strchr(buffer,'\n'))!=NULL)
172	{
173	*s='\0';
174	}
175	PrintS(" ");
176	PrintS(buffer);
177	if((s!=NULL) \|\|(l>20))
178	{
179	Print("..., %d char(s)",l);
180	}
181	break;
182	}
183	case LIST_CMD: Print(", size: %d",IDLIST(h)->nr+1);
184	break;
185	case QRING_CMD:
186	case RING_CMD:
187	if ((IDRING(h)==currRing) && (currRingHdl!=h))
188	PrintS("()"); / this is an alias to currRing */
189	#ifdef RDEBUG
190	if (traceit &TRACE_SHOW_RINGS)
191	Print(" <%lx>",(long)(IDRING(h)));
192	#endif
193	break;
194	/default: break;/
195	}
196	PrintLn();
197	}
198
199	void type_cmd(idhdl h)
200	{
201	BOOLEAN oldShortOut = FALSE;
202
203	if (currRing != NULL)
204	{
205	oldShortOut = currRing->ShortOut;
206	currRing->ShortOut = 1;
207	}
208	list1("// ",h,FALSE,FALSE);
209	if (IDTYP(h)!=INT_CMD)
210	{
211	sleftv expr;
212	memset(&expr,0,sizeof(expr));
213	expr.rtyp=IDHDL;
214	expr.name=IDID(h);
215	expr.data=(void *)h;
216	expr.Print();
217	}
218	if (currRing != NULL)
219	currRing->ShortOut = oldShortOut;
220	}
221
222	static void killlocals0(int v, idhdl * localhdl, const ring r)
223	{
224	idhdl h = *localhdl;
225	while (h!=NULL)
226	{
227	int vv;
228	//Print("consider %s, lev: %d:",IDID(h),IDLEV(h));
229	if ((vv=IDLEV(h))>0)
230	{
231	if (vv < v)
232	{
233	if (iiNoKeepRing)
234	{
235	//PrintS(" break\n");
236	return;
237	}
238	h = IDNEXT(h);
239	//PrintLn();
240	}
241	else //if (vv >= v)
242	{
243	idhdl nexth = IDNEXT(h);
244	killhdl2(h,localhdl,r);
245	h = nexth;
246	//PrintS("kill\n");
247	}
248	}
249	else
250	{
251	h = IDNEXT(h);
252	//PrintLn();
253	}
254	}
255	}
256	void killlocals_list(lists l,int v)
257	{
258	int i;
259	for(i=l->nr; i>=0; i--)
260	{
261	if (l->m[i].rtyp == LIST_CMD)
262	killlocals_list((lists)l->m[i].data,v);
263	else if ((l->m[i].rtyp == RING_CMD) \|\| (l->m[i].rtyp == QRING_CMD))
264	killlocals0(v,&(((ring)(l->m[i].data))->idroot),currRing);
265	}
266	}
267	void killlocals_rec(idhdl *root,int v, ring r)
268	{
269	idhdl h=*root;
270	while (h!=NULL)
271	{
272	if (IDLEV(h)>=v)
273	{
274	// Print("kill %s, lev %d for lev %d\n",IDID(h),IDLEV(h),v);
275	idhdl n=IDNEXT(h);
276	killhdl2(h,root,r);
277	h=n;
278	}
279	else if (IDTYP(h)==PACKAGE_CMD)
280	{
281	// Print("into pack %s, lev %d for lev %d\n",IDID(h),IDLEV(h),v);
282	if (IDPACKAGE(h)!=basePack)
283	killlocals_rec(&(IDRING(h)->idroot),v,r);
284	h=IDNEXT(h);
285	}
286	else if ((IDTYP(h)==RING_CMD)
287	\|\|(IDTYP(h)==QRING_CMD))
288	{
289	if ((IDRING(h)!=NULL) && (IDRING(h)->idroot!=NULL))
290	// we have to test IDRING(h)!=NULL: qring Q=groebner(...): killlocals
291	{
292	// Print("into ring %s, lev %d for lev %d\n",IDID(h),IDLEV(h),v);
293	killlocals_rec(&(IDRING(h)->idroot),v,IDRING(h));
294	}
295	h=IDNEXT(h);
296	}
297	else
298	{
299	// Print("skip %s lev %d for lev %d\n",IDID(h),IDLEV(h),v);
300	h=IDNEXT(h);
301	}
302	}
303	}
304	BOOLEAN killlocals_list(int v, lists L)
305	{
306	if (L==NULL) return FALSE;
307	BOOLEAN changed=FALSE;
308	int n=L->nr;
309	for(;n>=0;n--)
310	{
311	leftv h=&(L->m[n]);
312	void *d=h->data;
313	if (((h->rtyp==RING_CMD) \|\| (h->rtyp==QRING_CMD))
314	&& (((ring)d)->idroot!=NULL))
315	{
316	if (d!=currRing) {changed=TRUE;rChangeCurrRing((ring)d);}
317	killlocals0(v,&(((ring)h->data)->idroot),(ring)h->data);
318	}
319	else if (h->rtyp==LIST_CMD)
320	changed\|=killlocals_list(v,(lists)d);
321	}
322	return changed;
323	}
324	void killlocals(int v)
325	{
326	BOOLEAN changed=FALSE;
327	idhdl sh=currRingHdl;
328	ring cr=currRing;
329	if (sh!=NULL) changed=((IDLEV(sh)<v) \|\| (IDRING(sh)->ref>0));
330	//if (changed) Print("currRing=%s(%x), lev=%d,ref=%d\n",IDID(sh),IDRING(sh),IDLEV(sh),IDRING(sh)->ref);
331
332	killlocals_rec(&(basePack->idroot),v,currRing);
333
334	if (iiRETURNEXPR_len > myynest)
335	{
336	int t=iiRETURNEXPR[myynest].Typ();
337	if ((/iiRETURNEXPR[myynest].Typ()/ t==RING_CMD)
338	\|\| (/iiRETURNEXPR[myynest].Typ()/ t==QRING_CMD))
339	{
340	leftv h=&iiRETURNEXPR[myynest];
341	if (((ring)h->data)->idroot!=NULL)
342	killlocals0(v,&(((ring)h->data)->idroot),(ring)h->data);
343	}
344	else if (/iiRETURNEXPR[myynest].Typ()/ t==LIST_CMD)
345	{
346	leftv h=&iiRETURNEXPR[myynest];
347	changed \|=killlocals_list(v,(lists)h->data);
348	}
349	}
350	if (changed)
351	{
352	currRingHdl=rFindHdl(cr,NULL,NULL);
353	if (currRingHdl==NULL)
354	currRing=NULL;
355	else
356	rChangeCurrRing(cr);
357	}
358
359	if (myynest<=1) iiNoKeepRing=TRUE;
360	//Print("end killlocals >= %d\n",v);
361	//listall();
362	}
363
364	void list_cmd(int typ, const char* what, const char *prefix,BOOLEAN iterate, BOOLEAN fullname)
365	{
366	idhdl h,start;
367	BOOLEAN all = typ<0;
368	BOOLEAN really_all=FALSE;
369	BOOLEAN do_packages=FALSE;
370
371	if ( typ == -1 ) do_packages=TRUE;
372	if ( typ==0 )
373	{
374	if (strcmp(what,"all")==0)
375	{
376	really_all=TRUE;
377	h=basePack->idroot;
378	}
379	else
380	{
381	h = ggetid(what);
382	if (h!=NULL)
383	{
384	if (iterate) list1(prefix,h,TRUE,fullname);
385	if ((IDTYP(h)==RING_CMD)
386	\|\| (IDTYP(h)==QRING_CMD)
387	//\|\| (IDTYP(h)==PACKE_CMD)
388	)
389	{
390	h=IDRING(h)->idroot;
391	}
392	else if((IDTYP(h)==PACKAGE_CMD) \|\| (IDTYP(h)==POINTER_CMD))
393	{
394	//Print("list_cmd:package or pointer\n");
395	all=TRUE;typ=PROC_CMD;fullname=TRUE;really_all=TRUE;
396	h=IDPACKAGE(h)->idroot;
397	}
398	else
399	return;
400	}
401	else
402	{
403	Werror("%s is undefined",what);
404	return;
405	}
406	}
407	all=TRUE;
408	}
409	else if (RingDependend(typ))
410	{
411	h = currRing->idroot;
412	}
413	else
414	h = IDROOT;
415	start=h;
416	while (h!=NULL)
417	{
418	if ((all && (IDTYP(h)!=PROC_CMD) &&(IDTYP(h)!=PACKAGE_CMD))
419	\|\| (typ == IDTYP(h))
420	\|\| ((IDTYP(h)==QRING_CMD) && (typ==RING_CMD)))
421	{
422	list1(prefix,h,start==currRingHdl, fullname);
423	if (((IDTYP(h)==RING_CMD)\|\|(IDTYP(h)==QRING_CMD))
424	&& (really_all \|\| (all && (h==currRingHdl)))
425	&& ((IDLEV(h)==0)\|\|(IDLEV(h)==myynest)))
426	{
427	list_cmd(0,IDID(h),"// ",FALSE);
428	}
429	if (IDTYP(h)==PACKAGE_CMD && really_all)
430	{
431	package save_p=currPack;
432	currPack=IDPACKAGE(h);
433	list_cmd(0,IDID(h),"// ",FALSE);
434	currPack=save_p;
435	}
436	}
437	h = IDNEXT(h);
438	}
439	}
440
441	void test_cmd(int i)
442	{
443	int ii;
444
445	if (i<0)
446	{
447	ii= -i;
448	if (ii < 32)
449	{
450	test &= ~Sy_bit(ii);
451	}
452	else if (ii < 64)
453	{
454	verbose &= ~Sy_bit(ii-32);
455	}
456	else
457	WerrorS("out of bounds\n");
458	}
459	else if (i<32)
460	{
461	ii=i;
462	if (Sy_bit(ii) & kOptions)
463	{
464	Warn("Gerhard, use the option command");
465	test \|= Sy_bit(ii);
466	}
467	else if (Sy_bit(ii) & validOpts)
468	test \|= Sy_bit(ii);
469	}
470	else if (i<64)
471	{
472	ii=i-32;
473	verbose \|= Sy_bit(ii);
474	}
475	else
476	WerrorS("out of bounds\n");
477	}
478
479	int exprlist_length(leftv v)
480	{
481	int rc = 0;
482	while (v!=NULL)
483	{
484	switch (v->Typ())
485	{
486	case INT_CMD:
487	case POLY_CMD:
488	case VECTOR_CMD:
489	case NUMBER_CMD:
490	rc++;
491	break;
492	case INTVEC_CMD:
493	case INTMAT_CMD:
494	rc += ((intvec *)(v->Data()))->length();
495	break;
496	case MATRIX_CMD:
497	case IDEAL_CMD:
498	case MODUL_CMD:
499	{
500	matrix mm = (matrix)(v->Data());
501	rc += mm->rows() * mm->cols();
502	}
503	break;
504	case LIST_CMD:
505	rc+=((lists)v->Data())->nr+1;
506	break;
507	default:
508	rc++;
509	}
510	v = v->next;
511	}
512	return rc;
513	}
514
515	int iiIsPrime0(unsigned p) /* brute force !!!! */
516	{
517	unsigned i,j=0 /only to avoid compiler warnings/;
518	#ifdef HAVE_FACTORY
519	if (p<=32749) // max. small prime in factory
520	{
521	int a=0;
522	int e=cf_getNumSmallPrimes()-1;
523	i=e/2;
524	do
525	{
526	j=cf_getSmallPrime(i);
527	if (p==j) return p;
528	if (p<j) e=i-1;
529	else a=i+1;
530	i=a+(e-a)/2;
531	} while ( a<= e);
532	if (p>j) return j;
533	else return cf_getSmallPrime(i-1);
534	}
535	#endif
536	#ifdef HAVE_FACTORY
537	unsigned end_i=cf_getNumSmallPrimes()-1;
538	#else
539	unsigned end_i=p/2;
540	#endif
541	unsigned end_p=(unsigned)sqrt((double)p);
542	restart:
543	for (i=0; i<end_i; i++)
544	{
545	#ifdef HAVE_FACTORY
546	j=cf_getSmallPrime(i);
547	#else
548	if (i==0) j=2;
549	else j=2*i-1;
550	#endif
551	if ((p%j) == 0)
552	{
553	#ifdef HAVE_FACTORY
554	if (p<=32751) return iiIsPrime0(p-2);
555	#endif
556	p-=2;
557	goto restart;
558	}
559	if (j > end_p) return p;
560	}
561	#ifdef HAVE_FACTORY
562	if (i>=end_i)
563	{
564	while(j<=end_p)
565	{
566	j+=2;
567	if ((p%j) == 0)
568	{
569	if (p<=32751) return iiIsPrime0(p-2);
570	p-=2;
571	goto restart;
572	}
573	}
574	}
575	#endif
576	return p;
577	}
578	int IsPrime(int p) /* brute force !!!! */
579	{
580	int i,j;
581	if (p == 0) return 0;
582	else if (p == 1) return 1/1/;
583	else if ((p == 2)\|\|(p==3)) return p;
584	else if (p < 0) return 2; //(iiIsPrime0((unsigned)(-p)));
585	else if ((p & 1)==0) return iiIsPrime0((unsigned)(p-1));
586	return iiIsPrime0((unsigned)(p));
587	}
588
589	BOOLEAN iiWRITE(leftv res,leftv v)
590	{
591	sleftv vf;
592	if (iiConvert(v->Typ(),LINK_CMD,iiTestConvert(v->Typ(),LINK_CMD),v,&vf))
593	{
594	WerrorS("link expected");
595	return TRUE;
596	}
597	si_link l=(si_link)vf.Data();
598	if (vf.next == NULL)
599	{
600	WerrorS("write: need at least two arguments");
601	return TRUE;
602	}
603
604	BOOLEAN b=slWrite(l,vf.next); /* iiConvert preserves next */
605	if (b)
606	{
607	const char *s;
608	if ((l!=NULL)&&(l->name!=NULL)) s=l->name;
609	else s=sNoName;
610	Werror("cannot write to %s",s);
611	}
612	vf.CleanUp();
613	return b;
614	}
615
616	leftv iiMap(map theMap, const char * what)
617	{
618	idhdl w,r;
619	leftv v;
620	int i;
621	nMapFunc nMap;
622
623	r=IDROOT->get(theMap->preimage,myynest);
624	if ((currPack!=basePack)
625	&&((r==NULL) \|\| ((r->typ != RING_CMD) && (r->typ != QRING_CMD))))
626	r=basePack->idroot->get(theMap->preimage,myynest);
627	if ((r==NULL) && (currRingHdl!=NULL)
628	&& (strcmp(theMap->preimage,IDID(currRingHdl))==0))
629	{
630	r=currRingHdl;
631	}
632	if ((r!=NULL) && ((r->typ == RING_CMD) \|\| (r->typ== QRING_CMD)))
633	{
634	//if ((nMap=nSetMap(rInternalChar(IDRING(r)),
635	// IDRING(r)->parameter,
636	// rPar(IDRING(r)),
637	// IDRING(r)->minpoly)))
638	if ((nMap=nSetMap(IDRING(r)))==NULL)
639	{
640	if (rEqual(IDRING(r),currRing))
641	{
642	nMap=nCopy;
643	}
644	else
645	{
646	Werror("can not map from ground field of %s to current ground field",
647	theMap->preimage);
648	return NULL;
649	}
650	}
651	if (IDELEMS(theMap)<IDRING(r)->N)
652	{
653	theMap->m=(polyset)omReallocSize((ADDRESS)theMap->m,
654	IDELEMS(theMap)*sizeof(poly),
655	(IDRING(r)->N)*sizeof(poly));
656	for(i=IDELEMS(theMap);i<IDRING(r)->N;i++)
657	theMap->m[i]=NULL;
658	IDELEMS(theMap)=IDRING(r)->N;
659	}
660	if (what==NULL)
661	{
662	WerrorS("argument of a map must have a name");
663	}
664	else if ((w=IDRING(r)->idroot->get(what,myynest))!=NULL)
665	{
666	char *save_r=NULL;
667	v=(leftv)omAlloc0Bin(sleftv_bin);
668	sleftv tmpW;
669	memset(&tmpW,0,sizeof(sleftv));
670	tmpW.rtyp=IDTYP(w);
671	if (tmpW.rtyp==MAP_CMD)
672	{
673	tmpW.rtyp=IDEAL_CMD;
674	save_r=IDMAP(w)->preimage;
675	IDMAP(w)->preimage=0;
676	}
677	tmpW.data=IDDATA(w);
678	#if 0
679	if (((tmpW.rtyp==IDEAL_CMD)\|\|(tmpW.rtyp==MODUL_CMD)) && idIs0(IDIDEAL(w)))
680	{
681	v->rtyp=tmpW.rtyp;
682	v->data=idInit(IDELEMS(IDIDEAL(w)),IDIDEAL(w)->rank);
683	}
684	else
685	#endif
686	{
687	#ifdef FAST_MAP
688	if ((tmpW.rtyp==IDEAL_CMD) && (nMap==nCopy)
689	#ifdef HAVE_PLURAL
690	&& (!rIsPluralRing(currRing))
691	#endif
692	)
693	{
694	v->rtyp=IDEAL_CMD;
695	v->data=fast_map(IDIDEAL(w), IDRING(r), (ideal)theMap, currRing);
696	}
697	else
698	#endif
699	if (maApplyFetch(MAP_CMD,theMap,v,&tmpW,IDRING(r),NULL,NULL,0,nMap))
700	{
701	Werror("cannot map %s(%d)",Tok2Cmdname(w->typ),w->typ);
702	omFreeBin((ADDRESS)v, sleftv_bin);
703	if (save_r!=NULL) IDMAP(w)->preimage=save_r;
704	return NULL;
705	}
706	}
707	if (save_r!=NULL)
708	{
709	IDMAP(w)->preimage=save_r;
710	IDMAP((idhdl)v)->preimage=omStrDup(save_r);
711	v->rtyp=MAP_CMD;
712	}
713	return v;
714	}
715	else
716	{
717	Werror("%s undefined in %s",what,theMap->preimage);
718	}
719	}
720	else
721	{
722	Werror("cannot find preimage %s",theMap->preimage);
723	}
724	return NULL;
725	}
726
727	#ifdef OLD_RES
728	void iiMakeResolv(resolvente r, int length, int rlen, char * name, int typ0,
729	intvec ** weights)
730	{
731	lists L=liMakeResolv(r,length,rlen,typ0,weights);
732	int i=0;
733	idhdl h;
734	char * s=(char *)omAlloc(strlen(name)+5);
735
736	while (i<=L->nr)
737	{
738	sprintf(s,"%s(%d)",name,i+1);
739	if (i==0)
740	h=enterid(s,myynest,typ0,&(currRing->idroot), FALSE);
741	else
742	h=enterid(s,myynest,MODUL_CMD,&(currRing->idroot), FALSE);
743	if (h!=NULL)
744	{
745	h->data.uideal=(ideal)L->m[i].data;
746	h->attribute=L->m[i].attribute;
747	if (BVERBOSE(V_DEF_RES))
748	Print("//defining: %s as %d-th syzygy module\n",s,i+1);
749	}
750	else
751	{
752	idDelete((ideal *)&(L->m[i].data));
753	Warn("cannot define %s",s);
754	}
755	//L->m[i].data=NULL;
756	//L->m[i].rtyp=0;
757	//L->m[i].attribute=NULL;
758	i++;
759	}
760	omFreeSize((ADDRESS)L->m,(L->nr+1)*sizeof(sleftv));
761	omFreeBin((ADDRESS)L, slists_bin);
762	omFreeSize((ADDRESS)s,strlen(name)+5);
763	}
764	#endif
765
766	//resolvente iiFindRes(char * name, int * len, int *typ0)
767	//{
768	// char s=(char )omAlloc(strlen(name)+5);
769	// int i=-1;
770	// resolvente r;
771	// idhdl h;
772	//
773	// do
774	// {
775	// i++;
776	// sprintf(s,"%s(%d)",name,i+1);
777	// h=currRing->idroot->get(s,myynest);
778	// } while (h!=NULL);
779	// *len=i-1;
780	// if (*len<=0)
781	// {
782	// Werror("no objects %s(1),.. found",name);
783	// omFreeSize((ADDRESS)s,strlen(name)+5);
784	// return NULL;
785	// }
786	// r=(ideal )omAlloc(/(len+1)/ isizeof(ideal));
787	// memset(r,0,(len)sizeof(ideal));
788	// i=-1;
789	// *typ0=MODUL_CMD;
790	// while (i<(*len))
791	// {
792	// i++;
793	// sprintf(s,"%s(%d)",name,i+1);
794	// h=currRing->idroot->get(s,myynest);
795	// if (h->typ != MODUL_CMD)
796	// {
797	// if ((i!=0) \|\| (h->typ!=IDEAL_CMD))
798	// {
799	// Werror("%s is not of type module",s);
800	// omFreeSize((ADDRESS)r,(len)sizeof(ideal));
801	// omFreeSize((ADDRESS)s,strlen(name)+5);
802	// return NULL;
803	// }
804	// *typ0=IDEAL_CMD;
805	// }
806	// if ((i>0) && (idIs0(r[i-1])))
807	// {
808	// *len=i-1;
809	// break;
810	// }
811	// r[i]=IDIDEAL(h);
812	// }
813	// omFreeSize((ADDRESS)s,strlen(name)+5);
814	// return r;
815	//}
816
817	static resolvente iiCopyRes(resolvente r, int l)
818	{
819	int i;
820	resolvente res=(ideal )omAlloc0((l+1)sizeof(ideal));
821
822	for (i=0; i<l; i++)
823	res[i]=idCopy(r[i]);
824	return res;
825	}
826
827	BOOLEAN jjMINRES(leftv res, leftv v)
828	{
829	int len=0;
830	int typ0;
831	lists L=(lists)v->Data();
832	intvec weights=(intvec)atGet(v,"isHomog",INTVEC_CMD);
833	int add_row_shift = 0;
834	if (weights==NULL)
835	weights=(intvec*)atGet(&(L->m[0]),"isHomog",INTVEC_CMD);
836	if (weights!=NULL) add_row_shift=weights->min_in();
837	resolvente rr=liFindRes(L,&len,&typ0);
838	if (rr==NULL) return TRUE;
839	resolvente r=iiCopyRes(rr,len);
840
841	syMinimizeResolvente(r,len,0);
842	omFreeSize((ADDRESS)rr,len*sizeof(ideal));
843	len++;
844	res->data=(char *)liMakeResolv(r,len,-1,typ0,NULL,add_row_shift);
845	return FALSE;
846	}
847
848	BOOLEAN jjBETTI(leftv res, leftv u)
849	{
850	sleftv tmp;
851	memset(&tmp,0,sizeof(tmp));
852	tmp.rtyp=INT_CMD;
853	tmp.data=(void *)1;
854	if ((u->Typ()==IDEAL_CMD)
855	\|\| (u->Typ()==MODUL_CMD))
856	return jjBETTI2_ID(res,u,&tmp);
857	else
858	return jjBETTI2(res,u,&tmp);
859	}
860
861	BOOLEAN jjBETTI2_ID(leftv res, leftv u, leftv v)
862	{
863	lists l=(lists) omAllocBin(slists_bin);
864	l->Init(1);
865	l->m[0].rtyp=u->Typ();
866	l->m[0].data=u->Data();
867	l->m[0].attribute=u->attribute;
868	sleftv tmp2;
869	memset(&tmp2,0,sizeof(tmp2));
870	tmp2.rtyp=LIST_CMD;
871	tmp2.data=(void *)l;
872	BOOLEAN r=jjBETTI2(res,&tmp2,v);
873	l->m[0].data=NULL;
874	l->m[0].attribute=NULL;
875	l->m[0].rtyp=DEF_CMD;
876	l->Clean();
877	return r;
878	}
879
880	BOOLEAN jjBETTI2(leftv res, leftv u, leftv v)
881	{
882	resolvente r;
883	int len;
884	int reg,typ0;
885	lists l=(lists)u->Data();
886
887	intvec *weights=NULL;
888	int add_row_shift=0;
889	intvec ww=(intvec )atGet(&(l->m[0]),"isHomog",INTVEC_CMD);
890	if (ww!=NULL)
891	{
892	weights=ivCopy(ww);
893	add_row_shift = ww->min_in();
894	(*weights) -= add_row_shift;
895	}
896	//Print("attr:%x\n",weights);
897
898	r=liFindRes(l,&len,&typ0);
899	if (r==NULL) return TRUE;
900	res->data=(char *)syBetti(r,len,&reg,weights,(int)(long)v->Data());
901	omFreeSize((ADDRESS)r,(len)*sizeof(ideal));
902	atSet(res,omStrDup("rowShift"),(void*)add_row_shift,INT_CMD);
903	if (weights!=NULL) delete weights;
904	return FALSE;
905	}
906
907	int iiRegularity(lists L)
908	{
909	int len,reg,typ0;
910
911	resolvente r=liFindRes(L,&len,&typ0);
912
913	if (r==NULL)
914	return -2;
915	intvec *weights=NULL;
916	int add_row_shift=0;
917	intvec ww=(intvec )atGet(&(L->m[0]),"isHomog",INTVEC_CMD);
918	if (ww!=NULL)
919	{
920	weights=ivCopy(ww);
921	add_row_shift = ww->min_in();
922	(*weights) -= add_row_shift;
923	}
924	//Print("attr:%x\n",weights);
925
926	intvec *dummy=syBetti(r,len,&reg,weights);
927	if (weights!=NULL) delete weights;
928	delete dummy;
929	omFreeSize((ADDRESS)r,len*sizeof(ideal));
930	return reg+1+add_row_shift;
931	}
932
933	BOOLEAN iiDebugMarker=TRUE;
934	#define BREAK_LINE_LENGTH 80
935	void iiDebug()
936	{
937	Print("\n-- break point in %s --\n",VoiceName());
938	if (iiDebugMarker) VoiceBackTrack();
939	char * s;
940	iiDebugMarker=FALSE;
941	s = (char *)omAlloc(BREAK_LINE_LENGTH+4);
942	loop
943	{
944	memset(s,0,80);
945	fe_fgets_stdin("",s,BREAK_LINE_LENGTH);
946	if (s[BREAK_LINE_LENGTH-1]!='\0')
947	{
948	Print("line too long, max is %d chars\n",BREAK_LINE_LENGTH);
949	}
950	else
951	break;
952	}
953	if (*s=='\n')
954	{
955	iiDebugMarker=TRUE;
956	}
957	#if MDEBUG
958	else if(strncmp(s,"cont;",5)==0)
959	{
960	iiDebugMarker=TRUE;
961	}
962	#endif /* MDEBUG */
963	else
964	{
965	strcat( s, "\n;~\n");
966	newBuffer(s,BT_execute);
967	}
968	}
969
970	lists scIndIndset(ideal S, BOOLEAN all, ideal Q)
971	{
972	int i;
973	indset save;
974	lists res=(lists)omAlloc0Bin(slists_bin);
975
976	hexist = hInit(S, Q, &hNexist);
977	if ((hNexist == 0) \|\| (hisModule!=0))
978	{
979	res->Init(0);
980	return res;
981	}
982	save = ISet = (indset)omAlloc0Bin(indlist_bin);
983	hMu = 0;
984	hwork = (scfmon)omAlloc(hNexist * sizeof(scmon));
985	hvar = (varset)omAlloc((pVariables + 1) * sizeof(int));
986	hpure = (scmon)omAlloc((1 + (pVariables * pVariables)) * sizeof(long));
987	hrad = hexist;
988	hNrad = hNexist;
989	radmem = hCreate(pVariables - 1);
990	hCo = pVariables + 1;
991	hNvar = pVariables;
992	hRadical(hrad, &hNrad, hNvar);
993	hSupp(hrad, hNrad, hvar, &hNvar);
994	if (hNvar)
995	{
996	hCo = hNvar;
997	memset(hpure, 0, (pVariables + 1) * sizeof(long));
998	hPure(hrad, 0, &hNrad, hvar, hNvar, hpure, &hNpure);
999	hLexR(hrad, hNrad, hvar, hNvar);
1000	hDimSolve(hpure, hNpure, hrad, hNrad, hvar, hNvar);
1001	}
1002	if (hCo && (hCo < pVariables))
1003	{
1004	hIndMult(hpure, hNpure, hrad, hNrad, hvar, hNvar);
1005	}
1006	if (hMu!=0)
1007	{
1008	ISet = save;
1009	hMu2 = 0;
1010	if (all && (hCo+1 < pVariables))
1011	{
1012	JSet = (indset)omAlloc0Bin(indlist_bin);
1013	hIndAllMult(hpure, hNpure, hrad, hNrad, hvar, hNvar);
1014	i=hMu+hMu2;
1015	res->Init(i);
1016	if (hMu2 == 0)
1017	{
1018	omFreeBin((ADDRESS)JSet, indlist_bin);
1019	}
1020	}
1021	else
1022	{
1023	res->Init(hMu);
1024	}
1025	for (i=0;i<hMu;i++)
1026	{
1027	res->m[i].data = (void *)save->set;
1028	res->m[i].rtyp = INTVEC_CMD;
1029	ISet = save;
1030	save = save->nx;
1031	omFreeBin((ADDRESS)ISet, indlist_bin);
1032	}
1033	omFreeBin((ADDRESS)save, indlist_bin);
1034	if (hMu2 != 0)
1035	{
1036	save = JSet;
1037	for (i=hMu;i<hMu+hMu2;i++)
1038	{
1039	res->m[i].data = (void *)save->set;
1040	res->m[i].rtyp = INTVEC_CMD;
1041	JSet = save;
1042	save = save->nx;
1043	omFreeBin((ADDRESS)JSet, indlist_bin);
1044	}
1045	omFreeBin((ADDRESS)save, indlist_bin);
1046	}
1047	}
1048	else
1049	{
1050	res->Init(0);
1051	omFreeBin((ADDRESS)ISet, indlist_bin);
1052	}
1053	hKill(radmem, pVariables - 1);
1054	omFreeSize((ADDRESS)hpure, (1 + (pVariables * pVariables)) * sizeof(long));
1055	omFreeSize((ADDRESS)hvar, (pVariables + 1) * sizeof(int));
1056	omFreeSize((ADDRESS)hwork, hNexist * sizeof(scmon));
1057	hDelete(hexist, hNexist);
1058	return res;
1059	}
1060
1061	int iiDeclCommand(leftv sy, leftv name, int lev,int t, idhdl* root,BOOLEAN isring, BOOLEAN init_b)
1062	{
1063	BOOLEAN res=FALSE;
1064	const char *id = name->name;
1065
1066	memset(sy,0,sizeof(sleftv));
1067	if ((name->name==NULL)\|\|(isdigit(name->name[0])))
1068	{
1069	WerrorS("object to declare is not a name");
1070	res=TRUE;
1071	}
1072	else
1073	{
1074	//if (name->rtyp!=0)
1075	//{
1076	// Warn("`%s` is already in use",name->name);
1077	//}
1078	{
1079	sy->data = (char *)enterid(id,lev,t,root,init_b);
1080	}
1081	if (sy->data!=NULL)
1082	{
1083	sy->rtyp=IDHDL;
1084	currid=sy->name=IDID((idhdl)sy->data);
1085	// name->name=NULL; /* used in enterid */
1086	//sy->e = NULL;
1087	if (name->next!=NULL)
1088	{
1089	sy->next=(leftv)omAllocBin(sleftv_bin);
1090	res=iiDeclCommand(sy->next,name->next,lev,t,root, isring);
1091	}
1092	}
1093	else res=TRUE;
1094	}
1095	name->CleanUp();
1096	return res;
1097	}
1098
1099	BOOLEAN iiDefaultParameter(leftv p)
1100	{
1101	attr at=NULL;
1102	if (iiCurrProc!=NULL)
1103	at=iiCurrProc->attribute->get("default_arg");
1104	if (at==NULL)
1105	return FALSE;
1106	sleftv tmp;
1107	memset(&tmp,0,sizeof(sleftv));
1108	tmp.rtyp=at->atyp;
1109	tmp.data=at->CopyA();
1110	return iiAssign(p,&tmp);
1111	}
1112	BOOLEAN iiParameter(leftv p)
1113	{
1114	if (iiCurrArgs==NULL)
1115	{
1116	if (strcmp(p->name,"#")==0)
1117	return iiDefaultParameter(p);
1118	Werror("not enough arguments for proc %s",VoiceName());
1119	p->CleanUp();
1120	return TRUE;
1121	}
1122	leftv h=iiCurrArgs;
1123	if (strcmp(p->name,"#")==0)
1124	{
1125	iiCurrArgs=NULL;
1126	}
1127	else
1128	{
1129	iiCurrArgs=h->next;
1130	h->next=NULL;
1131	}
1132	BOOLEAN res=iiAssign(p,h);
1133	h->CleanUp();
1134	omFreeBin((ADDRESS)h, sleftv_bin);
1135	return res;
1136	}
1137
1138	static BOOLEAN iiInternalExport (leftv v, int toLev)
1139	{
1140	idhdl h=(idhdl)v->data;
1141	//Print("iiInternalExport('%s',%d)%s\n", v->name, toLev,"");
1142	if (IDLEV(h)==0) Warn("`%s` is already global",IDID(h));
1143	else
1144	{
1145	h=IDROOT->get(v->name,toLev);
1146	idhdl *root=&IDROOT;
1147	if ((h==NULL)&&(currRing!=NULL))
1148	{
1149	h=currRing->idroot->get(v->name,toLev);
1150	root=&currRing->idroot;
1151	}
1152	BOOLEAN keepring=FALSE;
1153	if ((h!=NULL)&&(IDLEV(h)==toLev))
1154	{
1155	if (IDTYP(h)==v->Typ())
1156	{
1157	if (((IDTYP(h)==RING_CMD)\|\|(IDTYP(h)==QRING_CMD))
1158	&& (v->Data()==IDDATA(h)))
1159	{
1160	IDRING(h)->ref++;
1161	keepring=TRUE;
1162	IDLEV(h)=toLev;
1163	//WarnS("keepring");
1164	return FALSE;
1165	}
1166	if (BVERBOSE(V_REDEFINE))
1167	{
1168	Warn("redefining %s",IDID(h));
1169	}
1170	#ifdef USE_IILOCALRING
1171	if (iiLocalRing[0]==IDRING(h) && (!keepring)) iiLocalRing[0]=NULL;
1172	#else
1173	proclevel *p=procstack;
1174	while (p->next!=NULL) p=p->next;
1175	if ((p->cRing==IDRING(h)) && (!keepring))
1176	{
1177	p->cRing=NULL;
1178	p->cRingHdl=NULL;
1179	}
1180	#endif
1181	killhdl2(h,root,currRing);
1182	}
1183	else
1184	{
1185	return TRUE;
1186	}
1187	}
1188	h=(idhdl)v->data;
1189	IDLEV(h)=toLev;
1190	if (keepring) IDRING(h)->ref--;
1191	iiNoKeepRing=FALSE;
1192	//Print("export %s\n",IDID(h));
1193	}
1194	return FALSE;
1195	}
1196
1197	BOOLEAN iiInternalExport (leftv v, int toLev, idhdl roothdl)
1198	{
1199	idhdl h=(idhdl)v->data;
1200	if(h==NULL)
1201	{
1202	Warn("'%s': no such identifier\n", v->name);
1203	return FALSE;
1204	}
1205	package frompack=v->req_packhdl;
1206	if (frompack==NULL) frompack=currPack;
1207	package rootpack = IDPACKAGE(roothdl);
1208	// Print("iiInternalExport('%s',%d,%s->%s) typ:%d\n", v->name, toLev, IDID(currPackHdl),IDID(roothdl),v->Typ());
1209	if ((RingDependend(IDTYP(h)))
1210	\|\| ((IDTYP(h)==LIST_CMD)
1211	&& (lRingDependend(IDLIST(h)))
1212	)
1213	)
1214	{
1215	//Print("// ==> Ringdependent set nesting to 0\n");
1216	return (iiInternalExport(v, toLev));
1217	}
1218	else
1219	{
1220	IDLEV(h)=toLev;
1221	v->req_packhdl=rootpack;
1222	if (h==frompack->idroot)
1223	{
1224	frompack->idroot=h->next;
1225	}
1226	else
1227	{
1228	idhdl hh=frompack->idroot;
1229	while ((hh!=NULL) && (hh->next!=h))
1230	hh=hh->next;
1231	if ((hh!=NULL) && (hh->next==h))
1232	hh->next=h->next;
1233	else
1234	{
1235	Werror("`%s` not found",v->Name());
1236	return TRUE;
1237	}
1238	}
1239	h->next=rootpack->idroot;
1240	rootpack->idroot=h;
1241	}
1242	return FALSE;
1243	}
1244
1245	BOOLEAN iiExport (leftv v, int toLev)
1246	{
1247	#ifndef NDEBUG
1248	checkall();
1249	#endif
1250	BOOLEAN nok=FALSE;
1251	leftv r=v;
1252	while (v!=NULL)
1253	{
1254	if ((v->name==NULL)\|\|(v->rtyp==0)\|\|(v->e!=NULL))
1255	{
1256	WerrorS("cannot export");
1257	nok=TRUE;
1258	}
1259	else
1260	{
1261	if(iiInternalExport(v, toLev))
1262	{
1263	r->CleanUp();
1264	return TRUE;
1265	}
1266	}
1267	v=v->next;
1268	}
1269	r->CleanUp();
1270	#ifndef NDEBUG
1271	checkall();
1272	#endif
1273	return nok;
1274	}
1275
1276	/assume root!=idroot/
1277	BOOLEAN iiExport (leftv v, int toLev, idhdl root)
1278	{
1279	#ifndef NDEBUG
1280	checkall();
1281	#endif
1282	// Print("iiExport1: pack=%s\n",IDID(root));
1283	package pack=IDPACKAGE(root);
1284	BOOLEAN nok=FALSE;
1285	leftv rv=v;
1286	while (v!=NULL)
1287	{
1288	if ((v->name==NULL)\|\|(v->rtyp==0)\|\|(v->e!=NULL)
1289	)
1290	{
1291	WerrorS("cannot export");
1292	nok=TRUE;
1293	}
1294	else
1295	{
1296	idhdl old=pack->idroot->get( v->name,toLev);
1297	if (old!=NULL)
1298	{
1299	if ((pack==currPack) && (old==(idhdl)v->data))
1300	{
1301	Warn("`%s` is already global",IDID(old));
1302	break;
1303	}
1304	else if (IDTYP(old)==v->Typ())
1305	{
1306	if (BVERBOSE(V_REDEFINE))
1307	{
1308	Warn("redefining %s",IDID(old));
1309	}
1310	v->name=omStrDup(v->name);
1311	killhdl2(old,&(pack->idroot),currRing);
1312	}
1313	else
1314	{
1315	rv->CleanUp();
1316	return TRUE;
1317	}
1318	}
1319	//Print("iiExport: pack=%s\n",IDID(root));
1320	if(iiInternalExport(v, toLev, root))
1321	{
1322	rv->CleanUp();
1323	return TRUE;
1324	}
1325	}
1326	v=v->next;
1327	}
1328	rv->CleanUp();
1329	#ifndef NDEBUG
1330	checkall();
1331	#endif
1332	return nok;
1333	}
1334
1335	BOOLEAN iiCheckRing(int i)
1336	{
1337	if (currRingHdl==NULL)
1338	{
1339	#ifdef SIQ
1340	if (siq<=0)
1341	{
1342	#endif
1343	if (RingDependend(i))
1344	{
1345	WerrorS("no ring active");
1346	return TRUE;
1347	}
1348	#ifdef SIQ
1349	}
1350	#endif
1351	}
1352	return FALSE;
1353	}
1354
1355	poly iiHighCorner(ideal I, int ak)
1356	{
1357	int i;
1358	if(!idIsZeroDim(I)) return NULL; // not zero-dim.
1359	poly po=NULL;
1360	if (rHasLocalOrMixedOrdering_currRing())
1361	{
1362	scComputeHC(I,currQuotient,ak,po);
1363	if (po!=NULL)
1364	{
1365	pGetCoeff(po)=nInit(1);
1366	for (i=pVariables; i>0; i--)
1367	{
1368	if (pGetExp(po, i) > 0) pDecrExp(po,i);
1369	}
1370	pSetComp(po,ak);
1371	pSetm(po);
1372	}
1373	}
1374	else
1375	po=pOne();
1376	return po;
1377	}
1378
1379	void iiCheckPack(package &p)
1380	{
1381	if (p==basePack) return;
1382
1383	idhdl t=basePack->idroot;
1384
1385	while ((t!=NULL) && (IDTYP(t)!=PACKAGE_CMD) && (IDPACKAGE(t)!=p)) t=t->next;
1386
1387	if (t==NULL)
1388	{
1389	WarnS("package not found\n");
1390	p=basePack;
1391	}
1392	return;
1393	}
1394
1395	idhdl rDefault(const char *s)
1396	{
1397	idhdl tmp=NULL;
1398
1399	if (s!=NULL) tmp = enterid(s, myynest, RING_CMD, &IDROOT);
1400	if (tmp==NULL) return NULL;
1401
1402	if (ppNoether!=NULL) pDelete(&ppNoether);
1403	if (sLastPrinted.RingDependend())
1404	{
1405	sLastPrinted.CleanUp();
1406	memset(&sLastPrinted,0,sizeof(sleftv));
1407	}
1408
1409	ring r = IDRING(tmp);
1410
1411	r->ch = 32003;
1412	r->N = 3;
1413	/r->P = 0; Alloc0 in idhdl::set, ipid.cc/
1414	/names/
1415	r->names = (char *) omAlloc0(3 sizeof(char_ptr));
1416	r->names[0] = omStrDup("x");
1417	r->names[1] = omStrDup("y");
1418	r->names[2] = omStrDup("z");
1419	/weights: entries for 3 blocks: NULL/
1420	r->wvhdl = (int *)omAlloc0(3 sizeof(int_ptr));
1421	/order: dp,C,0/
1422	r->order = (int ) omAlloc(3 sizeof(int *));
1423	r->block0 = (int )omAlloc0(3 sizeof(int *));
1424	r->block1 = (int )omAlloc0(3 sizeof(int *));
1425	/* ringorder dp for the first block: var 1..3 */
1426	r->order[0] = ringorder_dp;
1427	r->block0[0] = 1;
1428	r->block1[0] = 3;
1429	/* ringorder C for the second block: no vars */
1430	r->order[1] = ringorder_C;
1431	/* the last block: everything is 0 */
1432	r->order[2] = 0;
1433	/polynomial ring/
1434	r->OrdSgn = 1;
1435
1436	/* complete ring intializations */
1437	rComplete(r);
1438	rSetHdl(tmp);
1439	return currRingHdl;
1440	}
1441
1442	idhdl rFindHdl(ring r, idhdl n, idhdl w)
1443	{
1444	idhdl h=rSimpleFindHdl(r,IDROOT,n);
1445	if (h!=NULL) return h;
1446	if (IDROOT!=basePack->idroot) h=rSimpleFindHdl(r,basePack->idroot,n);
1447	if (h!=NULL) return h;
1448	proclevel *p=procstack;
1449	while(p!=NULL)
1450	{
1451	if ((p->cPack!=basePack)
1452	&& (p->cPack!=currPack))
1453	h=rSimpleFindHdl(r,p->cPack->idroot,n);
1454	if (h!=NULL) return h;
1455	p=p->next;
1456	}
1457	idhdl tmp=basePack->idroot;
1458	while (tmp!=NULL)
1459	{
1460	if (IDTYP(tmp)==PACKAGE_CMD)
1461	h=rSimpleFindHdl(r,IDPACKAGE(tmp)->idroot,n);
1462	if (h!=NULL) return h;
1463	tmp=IDNEXT(tmp);
1464	}
1465	return NULL;
1466	}
1467
1468	void rDecomposeCF(leftv h,const ring r,const ring R)
1469	{
1470	lists L=(lists)omAlloc0Bin(slists_bin);
1471	L->Init(4);
1472	h->rtyp=LIST_CMD;
1473	h->data=(void *)L;
1474	// 0: char/ cf - ring
1475	// 1: list (var)
1476	// 2: list (ord)
1477	// 3: qideal
1478	// ----------------------------------------
1479	// 0: char/ cf - ring
1480	L->m[0].rtyp=INT_CMD;
1481	L->m[0].data=(void *)r->ch;
1482	// ----------------------------------------
1483	// 1: list (var)
1484	lists LL=(lists)omAlloc0Bin(slists_bin);
1485	LL->Init(r->N);
1486	int i;
1487	for(i=0; i<r->N; i++)
1488	{
1489	LL->m[i].rtyp=STRING_CMD;
1490	LL->m[i].data=(void *)omStrDup(r->names[i]);
1491	}
1492	L->m[1].rtyp=LIST_CMD;
1493	L->m[1].data=(void *)LL;
1494	// ----------------------------------------
1495	// 2: list (ord)
1496	LL=(lists)omAlloc0Bin(slists_bin);
1497	i=rBlocks(r)-1;
1498	LL->Init(i);
1499	i--;
1500	lists LLL;
1501	for(; i>=0; i--)
1502	{
1503	intvec *iv;
1504	int j;
1505	LL->m[i].rtyp=LIST_CMD;
1506	LLL=(lists)omAlloc0Bin(slists_bin);
1507	LLL->Init(2);
1508	LLL->m[0].rtyp=STRING_CMD;
1509	LLL->m[0].data=(void *)omStrDup(rSimpleOrdStr(r->order[i]));
1510	if (r->block1[i]-r->block0[i] >=0 )
1511	{
1512	j=r->block1[i]-r->block0[i];
1513	if(r->order[i]==ringorder_M) j=(j+1)*(j+1)-1;
1514	iv=new intvec(j+1);
1515	if ((r->wvhdl!=NULL) && (r->wvhdl[i]!=NULL))
1516	{
1517	for(;j>=0; j--) (*iv)[j]=r->wvhdl[i][j];
1518	}
1519	else switch (r->order[i])
1520	{
1521	case ringorder_dp:
1522	case ringorder_Dp:
1523	case ringorder_ds:
1524	case ringorder_Ds:
1525	case ringorder_lp:
1526	for(;j>=0; j--) (*iv)[j]=1;
1527	break;
1528	default: /* do nothing */;
1529	}
1530	}
1531	else
1532	{
1533	iv=new intvec(1);
1534	}
1535	LLL->m[1].rtyp=INTVEC_CMD;
1536	LLL->m[1].data=(void *)iv;
1537	LL->m[i].data=(void *)LLL;
1538	}
1539	L->m[2].rtyp=LIST_CMD;
1540	L->m[2].data=(void *)LL;
1541	// ----------------------------------------
1542	// 3: qideal
1543	L->m[3].rtyp=IDEAL_CMD;
1544	if (R->minpoly==NULL)
1545	L->m[3].data=(void *)idInit(1,1);
1546	else
1547	{
1548	ideal I=idInit(1,1);
1549	L->m[3].data=(void *)I;
1550	I->m[0]=pNSet(R->minpoly);
1551	}
1552	// ----------------------------------------
1553	}
1554	void rDecomposeC(leftv h,const ring R)
1555	/* field is R or C */
1556	{
1557	lists L=(lists)omAlloc0Bin(slists_bin);
1558	if (rField_is_long_C(R)) L->Init(3);
1559	else L->Init(2);
1560	h->rtyp=LIST_CMD;
1561	h->data=(void *)L;
1562	// 0: char/ cf - ring
1563	// 1: list (var)
1564	// 2: list (ord)
1565	// ----------------------------------------
1566	// 0: char/ cf - ring
1567	L->m[0].rtyp=INT_CMD;
1568	L->m[0].data=(void *)0;
1569	// ----------------------------------------
1570	// 1:
1571	lists LL=(lists)omAlloc0Bin(slists_bin);
1572	LL->Init(2);
1573	LL->m[0].rtyp=INT_CMD;
1574	LL->m[0].data=(void *)si_max(R->float_len,SHORT_REAL_LENGTH/2);
1575	LL->m[1].rtyp=INT_CMD;
1576	LL->m[1].data=(void *)si_max(R->float_len2,SHORT_REAL_LENGTH);
1577	L->m[1].rtyp=LIST_CMD;
1578	L->m[1].data=(void *)LL;
1579	// ----------------------------------------
1580	// 2: list (par)
1581	if (rField_is_long_C(R))
1582	{
1583	L->m[2].rtyp=STRING_CMD;
1584	L->m[2].data=(void *)omStrDup(R->parameter[0]);
1585	}
1586	// ----------------------------------------
1587	}
1588
1589	#ifdef HAVE_RINGS
1590	void rDecomposeRing(leftv h,const ring R)
1591	/* field is R or C */
1592	{
1593	lists L=(lists)omAlloc0Bin(slists_bin);
1594	if (rField_is_Ring_Z(R)) L->Init(1);
1595	else L->Init(2);
1596	h->rtyp=LIST_CMD;
1597	h->data=(void *)L;
1598	// 0: char/ cf - ring
1599	// 1: list (module)
1600	// ----------------------------------------
1601	// 0: char/ cf - ring
1602	L->m[0].rtyp=STRING_CMD;
1603	L->m[0].data=(void *)omStrDup("integer");
1604	// ----------------------------------------
1605	// 1: module
1606	if (rField_is_Ring_Z(R)) return;
1607	lists LL=(lists)omAlloc0Bin(slists_bin);
1608	LL->Init(2);
1609	LL->m[0].rtyp=BIGINT_CMD;
1610	LL->m[0].data=nlMapGMP((number) R->ringflaga);
1611	LL->m[1].rtyp=INT_CMD;
1612	LL->m[1].data=(void *) R->ringflagb;
1613	L->m[1].rtyp=LIST_CMD;
1614	L->m[1].data=(void *)LL;
1615	}
1616	#endif
1617
1618
1619	lists rDecompose(const ring r)
1620	{
1621	// sanity check: require currRing==r for rings with polynomial data
1622	if ((r!=currRing)
1623	&& ((r->minpoly!=NULL) \|\| (r->qideal!=NULL) \|\| (r->minideal!=NULL)
1624	#ifdef HAVE_PLURAL
1625	\|\| (rIsPluralRing(r))
1626	#endif
1627	))
1628	{
1629	WerrorS("ring with polynomial data must be the base ring or compatible");
1630	return NULL;
1631	}
1632	// 0: char/ cf - ring
1633	// 1: list (var)
1634	// 2: list (ord)
1635	// 3: qideal
1636	// possibly:
1637	// 4: C
1638	// 5: D
1639	lists L=(lists)omAlloc0Bin(slists_bin);
1640	if (rIsPluralRing(r))
1641	L->Init(6);
1642	else
1643	L->Init(4);
1644	// ----------------------------------------
1645	// 0: char/ cf - ring
1646	#if 1 /* TODO */
1647	if (rField_is_numeric(r))
1648	{
1649	rDecomposeC(&(L->m[0]),r);
1650	}
1651	#ifdef HAVE_RINGS
1652	else if (rField_is_Ring(r))
1653	{
1654	rDecomposeRing(&(L->m[0]),r);
1655	}
1656	#endif
1657	else if (rIsExtension(r))
1658	{
1659	if (r->algring!=NULL)
1660	rDecomposeCF(&(L->m[0]),r->algring,r);
1661	else
1662	{
1663	lists Lc=(lists)omAlloc0Bin(slists_bin);
1664	Lc->Init(4);
1665	// char:
1666	Lc->m[0].rtyp=INT_CMD;
1667	Lc->m[0].data=(void*)r->ch;
1668	// var:
1669	lists Lv=(lists)omAlloc0Bin(slists_bin);
1670	Lv->Init(1);
1671	Lv->m[0].rtyp=STRING_CMD;
1672	Lv->m[0].data=(void *)omStrDup(r->parameter[0]);
1673	Lc->m[1].rtyp=LIST_CMD;
1674	Lc->m[1].data=(void*)Lv;
1675	// ord:
1676	lists Lo=(lists)omAlloc0Bin(slists_bin);
1677	Lo->Init(1);
1678	lists Loo=(lists)omAlloc0Bin(slists_bin);
1679	Loo->Init(2);
1680	Loo->m[0].rtyp=STRING_CMD;
1681	Loo->m[0].data=(void *)omStrDup(rSimpleOrdStr(ringorder_lp));
1682	intvec iv=new intvec(1); (iv)[0]=1;
1683	Loo->m[1].rtyp=INTVEC_CMD;
1684	Loo->m[1].data=(void *)iv;
1685	Lo->m[0].rtyp=LIST_CMD;
1686	Lo->m[0].data=(void*)Loo;
1687
1688	Lc->m[2].rtyp=LIST_CMD;
1689	Lc->m[2].data=(void*)Lo;
1690	// q-ideal:
1691	Lc->m[3].rtyp=IDEAL_CMD;
1692	Lc->m[3].data=(void *)idInit(1,1);
1693	// ----------------------
1694	L->m[0].rtyp=LIST_CMD;
1695	L->m[0].data=(void*)Lc;
1696	}
1697	if (L->m[0].rtyp==0)
1698	{
1699	//omFreeBin(slists_bin,(void *)L);
1700	return NULL;
1701	}
1702	}
1703	else
1704	#endif
1705	{
1706	L->m[0].rtyp=INT_CMD;
1707	L->m[0].data=(void *)r->ch;
1708	}
1709	// ----------------------------------------
1710	// 1: list (var)
1711	lists LL=(lists)omAlloc0Bin(slists_bin);
1712	LL->Init(r->N);
1713	int i;
1714	for(i=0; i<r->N; i++)
1715	{
1716	LL->m[i].rtyp=STRING_CMD;
1717	LL->m[i].data=(void *)omStrDup(r->names[i]);
1718	}
1719	L->m[1].rtyp=LIST_CMD;
1720	L->m[1].data=(void *)LL;
1721	// ----------------------------------------
1722	// 2: list (ord)
1723	LL=(lists)omAlloc0Bin(slists_bin);
1724	i=rBlocks(r)-1;
1725	LL->Init(i);
1726	i--;
1727	lists LLL;
1728	for(; i>=0; i--)
1729	{
1730	intvec *iv;
1731	int j;
1732	LL->m[i].rtyp=LIST_CMD;
1733	LLL=(lists)omAlloc0Bin(slists_bin);
1734	LLL->Init(2);
1735	LLL->m[0].rtyp=STRING_CMD;
1736	LLL->m[0].data=(void *)omStrDup(rSimpleOrdStr(r->order[i]));
1737	if (r->block1[i]-r->block0[i] >=0 )
1738	{
1739	j=r->block1[i]-r->block0[i];
1740	if (r->order[i]==ringorder_M) j=(j+1)*(j+1)-1;
1741	iv=new intvec(j+1);
1742	if ((r->wvhdl!=NULL) && (r->wvhdl[i]!=NULL))
1743	{
1744	for(;j>=0; j--) (*iv)[j]=r->wvhdl[i][j];
1745	}
1746	else switch (r->order[i])
1747	{
1748	case ringorder_dp:
1749	case ringorder_Dp:
1750	case ringorder_ds:
1751	case ringorder_Ds:
1752	case ringorder_lp:
1753	for(;j>=0; j--) (*iv)[j]=1;
1754	break;
1755	default: /* do nothing */;
1756	}
1757	}
1758	else
1759	{
1760	iv=new intvec(1);
1761	}
1762	LLL->m[1].rtyp=INTVEC_CMD;
1763	LLL->m[1].data=(void *)iv;
1764	LL->m[i].data=(void *)LLL;
1765	}
1766	L->m[2].rtyp=LIST_CMD;
1767	L->m[2].data=(void *)LL;
1768	// ----------------------------------------
1769	// 3: qideal
1770	L->m[3].rtyp=IDEAL_CMD;
1771	if (r->qideal==NULL)
1772	L->m[3].data=(void *)idInit(1,1);
1773	else
1774	L->m[3].data=(void *)idCopy(r->qideal);
1775	// ----------------------------------------
1776	#ifdef HAVE_PLURAL // NC! in rDecompose
1777	if (rIsPluralRing(r))
1778	{
1779	L->m[4].rtyp=MATRIX_CMD;
1780	L->m[4].data=(void *)mpCopy(r->GetNC()->C);
1781	L->m[5].rtyp=MATRIX_CMD;
1782	L->m[5].data=(void *)mpCopy(r->GetNC()->D);
1783	}
1784	#endif
1785	return L;
1786	}
1787
1788	void rComposeC(lists L, ring R)
1789	/* field is R or C */
1790	{
1791	// ----------------------------------------
1792	// 0: char/ cf - ring
1793	if ((L->m[0].rtyp!=INT_CMD) \|\| (L->m[0].data!=(char *)0))
1794	{
1795	Werror("invald coeff. field description, expecting 0");
1796	return;
1797	}
1798	R->ch=-1;
1799	// ----------------------------------------
1800	// 1:
1801	if (L->m[1].rtyp!=LIST_CMD)
1802	Werror("invald coeff. field description, expecting precision list");
1803	lists LL=(lists)L->m[1].data;
1804	int r1=(int)(long)LL->m[0].data;
1805	int r2=(int)(long)LL->m[1].data;
1806	if ((r1<=SHORT_REAL_LENGTH)
1807	&& (r2=SHORT_REAL_LENGTH))
1808	{
1809	R->float_len=SHORT_REAL_LENGTH/2;
1810	R->float_len2=SHORT_REAL_LENGTH;
1811	}
1812	else
1813	{
1814	R->float_len=si_min(r1,32767);
1815	R->float_len2=si_min(r2,32767);
1816	}
1817	// ----------------------------------------
1818	// 2: list (par)
1819	if (L->nr==2)
1820	{
1821	R->P=1;
1822	if (L->m[2].rtyp!=STRING_CMD)
1823	{
1824	Werror("invald coeff. field description, expecting parameter name");
1825	return;
1826	}
1827	R->parameter=(char*)omAlloc0(R->Psizeof(char_ptr));
1828	R->parameter[0]=omStrDup((char *)L->m[2].data);
1829	}
1830	// ----------------------------------------
1831	}
1832
1833	#ifdef HAVE_RINGS
1834	void rComposeRing(lists L, ring R)
1835	/* field is R or C */
1836	{
1837	// ----------------------------------------
1838	// 0: string: integer
1839	// no further entries --> Z
1840	R->ringflaga = (int_number) omAlloc(sizeof(mpz_t));
1841	if (L->nr == 0)
1842	{
1843	mpz_init_set_ui(R->ringflaga,0);
1844	R->ringflagb = 1;
1845	}
1846	// ----------------------------------------
1847	// 1:
1848	else
1849	{
1850	if (L->m[1].rtyp!=LIST_CMD) Werror("invald data, expecting list of numbers");
1851	lists LL=(lists)L->m[1].data;
1852	mpz_init(R->ringflaga);
1853	if ((LL->nr >= 0) && LL->m[0].rtyp == BIGINT_CMD)
1854	{
1855	number ringflaga = (number) LL->m[0].data;
1856	nlGMP(ringflaga, (number) R->ringflaga);
1857	LL->m[0].data = (void *)ringflaga;
1858	}
1859	else if ((LL->nr >= 0) && LL->m[0].rtyp == INT_CMD)
1860	{
1861	mpz_set_ui(R->ringflaga,(unsigned long) LL->m[0].data);
1862	}
1863	else
1864	{
1865	mpz_set_ui(R->ringflaga,0);
1866	}
1867	if (LL->nr >= 1)
1868	{
1869	R->ringflagb = (unsigned long) LL->m[1].data;
1870	}
1871	else
1872	{
1873	R->ringflagb = 1;
1874	}
1875	}
1876	// ----------------------------------------
1877	if ((mpz_cmp_ui(R->ringflaga, 1) == 0) && (mpz_cmp_ui(R->ringflaga, 0) < 0))
1878	{
1879	Werror("Wrong ground ring specification (module is 1)");
1880	return;
1881	}
1882	if (R->ringflagb < 1)
1883	{
1884	Werror("Wrong ground ring specification (exponent smaller than 1");
1885	return;
1886	}
1887	// module is 0 ---> integers
1888	if (mpz_cmp_ui(R->ringflaga, 0) == 0)
1889	{
1890	R->ch = 0;
1891	R->ringtype = 4;
1892	}
1893	// we have an exponent
1894	else if (R->ringflagb > 1)
1895	{
1896	R->ch = R->ringflagb;
1897	if ((mpz_cmp_ui(R->ringflaga, 2) == 0) && (R->ringflagb + 2 <= 8*sizeof(NATNUMBER)))
1898	{
1899	R->ringtype = 1; // Use Z/2^ch
1900	}
1901	else
1902	{
1903	R->ringtype = 3;
1904	}
1905	}
1906	// just a module m > 1
1907	else
1908	{
1909	R->ringtype = 2;
1910	R->ch = mpz_get_ui(R->ringflaga);
1911	}
1912	}
1913	#endif
1914
1915	ring rCompose(const lists L)
1916	{
1917	if ((L->nr!=3)
1918	#ifdef HAVE_PLURAL
1919	&&(L->nr!=5)
1920	#endif
1921	)
1922	return NULL;
1923	int is_gf_char=0;
1924	// 0: char/ cf - ring
1925	// 1: list (var)
1926	// 2: list (ord)
1927	// 3: qideal
1928	// possibly:
1929	// 4: C
1930	// 5: D
1931	ring R=(ring) omAlloc0Bin(sip_sring_bin);
1932	// ------------------------- VARS ---------------------------
1933	if (L->m[1].Typ()==LIST_CMD)
1934	{
1935	lists v=(lists)L->m[1].Data();
1936	R->N = v->nr+1;
1937	R->names = (char *)omAlloc0(R->N sizeof(char_ptr));
1938	int i;
1939	for(i=0;i<R->N;i++)
1940	{
1941	if (v->m[i].Typ()==STRING_CMD)
1942	R->names[i]=omStrDup((char *)v->m[i].Data());
1943	else if (v->m[i].Typ()==POLY_CMD)
1944	{
1945	poly p=(poly)v->m[i].Data();
1946	int nr=pIsPurePower(p);
1947	if (nr>0)
1948	R->names[i]=omStrDup(currRing->names[nr-1]);
1949	else
1950	{
1951	Werror("var name %d must be a string or a ring variable",i+1);
1952	goto rCompose_err;
1953	}
1954	}
1955	else
1956	{
1957	Werror("var name %d must be `string`",i+1);
1958	goto rCompose_err;
1959	}
1960	}
1961	}
1962	else
1963	{
1964	WerrorS("variable must be given as `list`");
1965	goto rCompose_err;
1966	}
1967	// ------------------------ ORDER ------------------------------
1968	if (L->m[2].Typ()==LIST_CMD)
1969	{
1970	lists v=(lists)L->m[2].Data();
1971	int n= v->nr+2;
1972	int j;
1973	// initialize fields of R
1974	R->order=(int )omAlloc0(nsizeof(int));
1975	R->block0=(int )omAlloc0(nsizeof(int));
1976	R->block1=(int )omAlloc0(nsizeof(int));
1977	R->wvhdl=(int*)omAlloc0(nsizeof(int_ptr));
1978	// init order, so that rBlocks works correctly
1979	for (j=0; j < n-1; j++)
1980	R->order[j] = (int) ringorder_unspec;
1981	// orderings
1982	R->OrdSgn=1;
1983	for(j=0;j<n-1;j++)
1984	{
1985	// todo: a(..), M
1986	if (v->m[j].Typ()!=LIST_CMD)
1987	{
1988	WerrorS("ordering must be list of lists");
1989	goto rCompose_err;
1990	}
1991	lists vv=(lists)v->m[j].Data();
1992	if ((vv->nr!=1)
1993	\|\| (vv->m[0].Typ()!=STRING_CMD)
1994	\|\| ((vv->m[1].Typ()!=INTVEC_CMD) && (vv->m[1].Typ()!=INT_CMD)))
1995	{
1996	WerrorS("ordering name must be a (string,intvec)");
1997	goto rCompose_err;
1998	}
1999	R->order[j]=rOrderName(omStrDup((char*)vv->m[0].Data())); // assume STRING
2000	if (j==0) R->block0[0]=1;
2001	else
2002	{
2003	int jj=j-1;
2004	while((jj>=0)
2005	&& ((R->order[jj]== ringorder_a)
2006	\|\| (R->order[jj]== ringorder_aa)
2007	\|\| (R->order[jj]== ringorder_c)
2008	\|\| (R->order[jj]== ringorder_C)
2009	))
2010	{
2011	//Print("jj=%, skip %s\n",rSimpleOrdStr(R->order[jj]));
2012	jj--;
2013	}
2014	if (jj<0) R->block0[j]=1;
2015	else R->block0[j]=R->block1[jj]+1;
2016	}
2017	intvec *iv;
2018	if (vv->m[1].Typ()==INT_CMD)
2019	iv=new intvec((int)(long)vv->m[1].Data(),(int)(long)vv->m[1].Data());
2020	else
2021	iv=ivCopy((intvec*)vv->m[1].Data()); //assume INTVEC
2022	R->block1[j]=si_max(R->block0[j],R->block0[j]+iv->length()-1);
2023	//Print("block %d from %d to %d\n",j,R->block0[j], R->block1[j]);
2024	int i;
2025	switch (R->order[j])
2026	{
2027	case ringorder_ws:
2028	case ringorder_Ws:
2029	R->OrdSgn=-1;
2030	case ringorder_aa:
2031	case ringorder_a:
2032	case ringorder_wp:
2033	case ringorder_Wp:
2034	R->wvhdl[j] =( int )omAlloc((iv->length())sizeof(int));
2035	for (i=0; i<iv->length();i++)
2036	{
2037	R->wvhdl[j][i]=(*iv)[i];
2038	}
2039	break;
2040	case ringorder_M:
2041	R->wvhdl[j] =( int )omAlloc((iv->length())sizeof(int));
2042	for (i=0; i<iv->length();i++) R->wvhdl[j][i]=(*iv)[i];
2043	R->block1[j]=si_max(R->block0[j],R->block0[j]+(int)sqrt((double)(iv->length()-1)));
2044	break;
2045	case ringorder_ls:
2046	case ringorder_ds:
2047	case ringorder_Ds:
2048	case ringorder_rs:
2049	R->OrdSgn=-1;
2050	case ringorder_lp:
2051	case ringorder_dp:
2052	case ringorder_Dp:
2053	case ringorder_rp:
2054	break;
2055	case ringorder_S:
2056	break;
2057	case ringorder_c:
2058	case ringorder_C:
2059	R->block1[j]=R->block0[j]=0;
2060	break;
2061	case 0:
2062	case ringorder_unspec:
2063	break;
2064	}
2065	delete iv;
2066	}
2067	// sanity check
2068	j=n-2;
2069	if ((R->order[j]==ringorder_c)
2070	\|\| (R->order[j]==ringorder_C)
2071	\|\| (R->order[j]==ringorder_unspec)) j--;
2072	if (R->block1[j] != R->N)
2073	{
2074	if (((R->order[j]==ringorder_dp) \|\|
2075	(R->order[j]==ringorder_ds) \|\|
2076	(R->order[j]==ringorder_Dp) \|\|
2077	(R->order[j]==ringorder_Ds) \|\|
2078	(R->order[j]==ringorder_rp) \|\|
2079	(R->order[j]==ringorder_rs) \|\|
2080	(R->order[j]==ringorder_lp) \|\|
2081	(R->order[j]==ringorder_ls))
2082	&&
2083	R->block0[j] <= R->N)
2084	{
2085	R->block1[j] = R->N;
2086	}
2087	else
2088	{
2089	Werror("ordering incomplete: size (%d) should be %d",R->block1[j],R->N);
2090	goto rCompose_err;
2091	}
2092	}
2093	}
2094	else
2095	{
2096	WerrorS("ordering must be given as `list`");
2097	goto rCompose_err;
2098	}
2099	// ------------------------------------------------------------------
2100	// 0: char:
2101	if (L->m[0].Typ()==INT_CMD)
2102	{
2103	R->ch=(int)(long)L->m[0].Data();
2104	if (R->ch!=-1)
2105	{
2106	int l=0;
2107	if (((R->ch!=0) && (R->ch<2) && (is_gf_char=-1))
2108	#ifndef NV_OPS
2109	\|\| (R->ch > 32003)
2110	#endif
2111	\|\| ((l=IsPrime(R->ch))!=R->ch)
2112	)
2113	{
2114	Warn("%d is invalid characteristic of ground field. %d is used.", R->ch,l);
2115	R->ch=l;
2116	}
2117	}
2118	}
2119	else if (L->m[0].Typ()==LIST_CMD)
2120	{
2121	lists LL=(lists)L->m[0].Data();
2122	#ifdef HAVE_RINGS
2123	if (LL->m[0].Typ() == STRING_CMD)
2124	{
2125	rComposeRing(LL,R); /* Ring */
2126	}
2127	else
2128	#endif
2129	if (LL->nr<3)
2130	rComposeC(LL,R); /* R, long_R, long_C */
2131	else
2132	{
2133	if (LL->m[0].Typ()==INT_CMD)
2134	{
2135	int ch=(int)(long)LL->m[0].Data();
2136	while ((ch!=fftable[is_gf_char]) && (fftable[is_gf_char])) is_gf_char++;
2137	if (fftable[is_gf_char]==0) is_gf_char=-1;
2138	}
2139	if (is_gf_char==-1)
2140	{
2141	R->algring=rCompose((lists)L->m[0].Data());
2142	if (R->algring==NULL)
2143	{
2144	WerrorS("could not create rational function coefficient field");
2145	goto rCompose_err;
2146	}
2147	if (R->algring->ch>0)
2148	R->ch= -R->algring->ch;
2149	else
2150	R->ch=1;
2151	R->P=R->algring->N;
2152	R->parameter=(char*)omAlloc0(R->Psizeof(char_ptr));
2153	int i;
2154	for(i=R->P-1;i>=0;i--)
2155	R->parameter[i]=omStrDup(R->algring->names[i]);
2156	if (R->algring->qideal!=NULL)
2157	{
2158	if (IDELEMS(R->algring->qideal)==1)
2159	{
2160	R->minpoly=naInit(1,R);
2161	lnumber n=(lnumber)R->minpoly;
2162	n->z=R->algring->qideal->m[0];
2163	R->algring->qideal->m[0]=NULL;
2164	idDelete(&(R->algring->qideal));
2165	}
2166	else
2167	{
2168	WerrorS("not implemented yet.");
2169	}
2170	}
2171	}
2172	else
2173	{ // gf-char
2174	R->ch=fftable[is_gf_char];
2175	R->P=1;
2176	R->parameter=(char*)omAlloc0(1sizeof(char_ptr));
2177	R->parameter[0]=omStrDup((char*)((lists)(LL->m[1].Data()))->m[0].Data());
2178	}
2179	}
2180	}
2181	else
2182	{
2183	WerrorS("coefficient field must be described by `int` or `list`");
2184	goto rCompose_err;
2185	}
2186	rNameCheck(R);
2187	// ------------------------ Q-IDEAL ------------------------
2188	rComplete(R);
2189
2190	if (L->m[3].Typ()==IDEAL_CMD)
2191	{
2192	ideal q=(ideal)L->m[3].Data();
2193	if (q->m[0]!=NULL)
2194	{
2195	if (R->ch!=currRing->ch)
2196	{
2197	#if 0
2198	WerrorS("coefficient fields must be equal if q-ideal !=0");
2199	goto rCompose_err;
2200	#else
2201	ring orig_ring=currRing;
2202	rChangeCurrRing(R);
2203	int *perm=NULL;
2204	int *par_perm=NULL;
2205	int par_perm_size=0;
2206	nMapFunc nMap;
2207	BOOLEAN bo;
2208
2209	if ((nMap=nSetMap(orig_ring))==NULL)
2210	{
2211	if (rEqual(orig_ring,currRing))
2212	{
2213	nMap=nCopy;
2214	}
2215	else
2216	// Allow imap/fetch to be make an exception only for:
2217	if ( (rField_is_Q_a(orig_ring) && // Q(a..) -> Q(a..) \|\| Q \|\| Zp \|\| Zp(a)
2218	(rField_is_Q() \|\| rField_is_Q_a() \|\|
2219	(rField_is_Zp() \|\| rField_is_Zp_a())))
2220	\|\|
2221	(rField_is_Zp_a(orig_ring) && // Zp(a..) -> Zp(a..) \|\| Zp
2222	(rField_is_Zp(currRing, rInternalChar(orig_ring)) \|\|
2223	rField_is_Zp_a(currRing, rInternalChar(orig_ring)))) )
2224	{
2225	par_perm_size=rPar(orig_ring);
2226	BITSET save_test=test;
2227	naSetChar(rInternalChar(orig_ring),orig_ring);
2228	nSetChar(currRing);
2229	test=save_test;
2230	}
2231	else
2232	{
2233	WerrorS("coefficient fields must be equal if q-ideal !=0");
2234	goto rCompose_err;
2235	}
2236	}
2237	perm=(int )omAlloc0((orig_ring->N+1)sizeof(int));
2238	if (par_perm_size!=0)
2239	par_perm=(int )omAlloc0(par_perm_sizesizeof(int));
2240	int i;
2241	#if 0
2242	// use imap:
2243	maFindPerm(orig_ring->names,orig_ring->N,orig_ring->parameter,orig_ring->P,
2244	currRing->names,currRing->N,currRing->parameter, currRing->P,
2245	perm,par_perm, currRing->ch);
2246	#else
2247	// use fetch
2248	if ((rPar(orig_ring)>0) && (rPar(currRing)==0))
2249	{
2250	for(i=si_min(rPar(orig_ring),rVar(currRing))-1;i>=0;i--) par_perm[i]=i+1;
2251	}
2252	else if (par_perm_size!=0)
2253	for(i=si_min(rPar(orig_ring),rPar(currRing))-1;i>=0;i--) par_perm[i]=-(i+1);
2254	for(i=si_min(orig_ring->N,pVariables);i>0;i--) perm[i]=i;
2255	#endif
2256	ideal dest_id=idInit(IDELEMS(q),1);
2257	for(i=IDELEMS(q)-1; i>=0; i--)
2258	{
2259	dest_id->m[i]=pPermPoly(q->m[i],perm,orig_ring,nMap,
2260	par_perm,par_perm_size);
2261	// PrintS("map:");pWrite(dest_id->m[i]);PrintLn();
2262	pTest(dest_id->m[i]);
2263	}
2264	R->qideal=dest_id;
2265	if (perm!=NULL)
2266	omFreeSize((ADDRESS)perm,(orig_ring->N+1)*sizeof(int));
2267	if (par_perm!=NULL)
2268	omFreeSize((ADDRESS)par_perm,par_perm_size*sizeof(int));
2269	rChangeCurrRing(orig_ring);
2270	#endif
2271	}
2272	else
2273	R->qideal=idrCopyR(q,currRing,R);
2274	}
2275	}
2276	else
2277	{
2278	WerrorS("q-ideal must be given as `ideal`");
2279	goto rCompose_err;
2280	}
2281
2282
2283	// ---------------------------------------------------------------
2284	#ifdef HAVE_PLURAL
2285	if (L->nr==5)
2286	{
2287	if (nc_CallPlural((matrix)L->m[4].Data(),(matrix)L->m[5].Data(),NULL,NULL,R, true)) goto rCompose_err;
2288	// takes care about non-comm. quotient! i.e. calls "nc_SetupQuotient" due to last true
2289	}
2290	#endif
2291	return R;
2292
2293	rCompose_err:
2294	if (R->N>0)
2295	{
2296	int i;
2297	if (R->names!=NULL)
2298	{
2299	i=R->N-1;
2300	while (i>=0) { if (R->names[i]!=NULL) omFree(R->names[i]); i--; }
2301	omFree(R->names);
2302	}
2303	}
2304	if (R->order!=NULL) omFree(R->order);
2305	if (R->block0!=NULL) omFree(R->block0);
2306	if (R->block1!=NULL) omFree(R->block1);
2307	if (R->wvhdl!=NULL) omFree(R->wvhdl);
2308	omFree(R);
2309	return NULL;
2310	}
2311
2312	// from matpol.cc
2313
2314	/*2
2315	* compute the jacobi matrix of an ideal
2316	*/
2317	BOOLEAN mpJacobi(leftv res,leftv a)
2318	{
2319	int i,j;
2320	matrix result;
2321	ideal id=(ideal)a->Data();
2322
2323	result =mpNew(IDELEMS(id),pVariables);
2324	for (i=1; i<=IDELEMS(id); i++)
2325	{
2326	for (j=1; j<=pVariables; j++)
2327	{
2328	MATELEM(result,i,j) = pDiff(id->m[i-1],j);
2329	}
2330	}
2331	res->data=(char *)result;
2332	return FALSE;
2333	}
2334
2335	/*2
2336	* returns the Koszul-matrix of degree d of a vectorspace with dimension n
2337	* uses the first n entrees of id, if id <> NULL
2338	*/
2339	BOOLEAN mpKoszul(leftv res,leftv c/ip/, leftv b/in/, leftv id)
2340	{
2341	int n=(int)(long)b->Data();
2342	int d=(int)(long)c->Data();
2343	int k,l,sign,row,col;
2344	matrix result;
2345	ideal temp;
2346	BOOLEAN bo;
2347	poly p;
2348
2349	if ((d>n) \|\| (d<1) \|\| (n<1))
2350	{
2351	res->data=(char *)mpNew(1,1);
2352	return FALSE;
2353	}
2354	int choise = (int)omAlloc(d*sizeof(int));
2355	if (id==NULL)
2356	temp=idMaxIdeal(1);
2357	else
2358	temp=(ideal)id->Data();
2359
2360	k = binom(n,d);
2361	l = k*d;
2362	l /= n-d+1;
2363	result =mpNew(l,k);
2364	col = 1;
2365	idInitChoise(d,1,n,&bo,choise);
2366	while (!bo)
2367	{
2368	sign = 1;
2369	for (l=1;l<=d;l++)
2370	{
2371	if (choise[l-1]<=IDELEMS(temp))
2372	{
2373	p = pCopy(temp->m[choise[l-1]-1]);
2374	if (sign == -1) p = pNeg(p);
2375	sign *= -1;
2376	row = idGetNumberOfChoise(l-1,d,1,n,choise);
2377	MATELEM(result,row,col) = p;
2378	}
2379	}
2380	col++;
2381	idGetNextChoise(d,n,&bo,choise);
2382	}
2383	if (id==NULL) idDelete(&temp);
2384
2385	res->data=(char *)result;
2386	return FALSE;
2387	}
2388
2389	// from syz1.cc
2390	/*2
2391	* read out the Betti numbers from resolution
2392	* (interpreter interface)
2393	*/
2394	BOOLEAN syBetti2(leftv res, leftv u, leftv w)
2395	{
2396	syStrategy syzstr=(syStrategy)u->Data();
2397	BOOLEAN minim=(int)(long)w->Data();
2398	int row_shift=0;
2399	int add_row_shift=0;
2400	intvec *weights=NULL;
2401	intvec ww=(intvec )atGet(u,"isHomog",INTVEC_CMD);
2402	if (ww!=NULL)
2403	{
2404	weights=ivCopy(ww);
2405	add_row_shift = ww->min_in();
2406	(*weights) -= add_row_shift;
2407	}
2408	res->data=(void *)syBettiOfComputation(syzstr,minim,&row_shift,weights);
2409	//row_shift += add_row_shift;
2410	//Print("row_shift=%d, add_row_shift=%d\n",row_shift,add_row_shift);
2411	atSet(res,omStrDup("rowShift"),(void*)add_row_shift,INT_CMD);
2412	return FALSE;
2413	}
2414	BOOLEAN syBetti1(leftv res, leftv u)
2415	{
2416	sleftv tmp;
2417	memset(&tmp,0,sizeof(tmp));
2418	tmp.rtyp=INT_CMD;
2419	tmp.data=(void *)1;
2420	return syBetti2(res,u,&tmp);
2421	}
2422
2423	/*3
2424	* converts a resolution into a list of modules
2425	*/
2426	lists syConvRes(syStrategy syzstr,BOOLEAN toDel,int add_row_shift)
2427	{
2428	if ((syzstr->fullres==NULL) && (syzstr->minres==NULL))
2429	{
2430	if (syzstr->hilb_coeffs==NULL)
2431	{
2432	syzstr->fullres = syReorder(syzstr->res,syzstr->length,syzstr);
2433	}
2434	else
2435	{
2436	syzstr->minres = syReorder(syzstr->orderedRes,syzstr->length,syzstr);
2437	syKillEmptyEntres(syzstr->minres,syzstr->length);
2438	}
2439	}
2440	resolvente tr;
2441	int typ0=IDEAL_CMD;
2442	if (syzstr->minres!=NULL)
2443	tr = syzstr->minres;
2444	else
2445	tr = syzstr->fullres;
2446	resolvente trueres=NULL;
2447	intvec ** w=NULL;
2448	if (syzstr->length>0)
2449	{
2450	trueres=(resolvente)omAlloc0((syzstr->length)*sizeof(ideal));
2451	for (int i=(syzstr->length)-1;i>=0;i--)
2452	{
2453	if (tr[i]!=NULL)
2454	{
2455	trueres[i] = idCopy(tr[i]);
2456	}
2457	}
2458	if (idRankFreeModule(trueres[0]) > 0)
2459	typ0 = MODUL_CMD;
2460	if (syzstr->weights!=NULL)
2461	{
2462	w = (intvec*)omAlloc0((syzstr->length)sizeof(intvec*));
2463	for (int i=(syzstr->length)-1;i>=0;i--)
2464	{
2465	if (syzstr->weights[i]!=NULL) w[i] = ivCopy(syzstr->weights[i]);
2466	}
2467	}
2468	}
2469	lists li = liMakeResolv(trueres,syzstr->length,syzstr->list_length,typ0,
2470	w,add_row_shift);
2471	if (w != NULL) omFreeSize(w, (syzstr->length)sizeof(intvec));
2472	if (toDel) syKillComputation(syzstr);
2473	return li;
2474	}
2475
2476	/*3
2477	* converts a list of modules into a resolution
2478	*/
2479	syStrategy syConvList(lists li,BOOLEAN toDel)
2480	{
2481	int typ0;
2482	syStrategy result=(syStrategy)omAlloc0(sizeof(ssyStrategy));
2483
2484	resolvente fr = liFindRes(li,&(result->length),&typ0,&(result->weights));
2485	if (fr != NULL)
2486	{
2487
2488	result->fullres = (resolvente)omAlloc0((result->length+1)*sizeof(ideal));
2489	for (int i=result->length-1;i>=0;i--)
2490	{
2491	if (fr[i]!=NULL)
2492	result->fullres[i] = idCopy(fr[i]);
2493	}
2494	result->list_length=result->length;
2495	omFreeSize((ADDRESS)fr,(result->length)*sizeof(ideal));
2496	}
2497	else
2498	{
2499	omFreeSize(result, sizeof(ssyStrategy));
2500	result = NULL;
2501	}
2502	if (toDel) li->Clean();
2503	return result;
2504	}
2505
2506	/*3
2507	* converts a list of modules into a minimal resolution
2508	*/
2509	syStrategy syForceMin(lists li)
2510	{
2511	int typ0;
2512	syStrategy result=(syStrategy)omAlloc0(sizeof(ssyStrategy));
2513
2514	resolvente fr = liFindRes(li,&(result->length),&typ0);
2515	result->minres = (resolvente)omAlloc0((result->length+1)*sizeof(ideal));
2516	for (int i=result->length-1;i>=0;i--)
2517	{
2518	if (fr[i]!=NULL)
2519	result->minres[i] = idCopy(fr[i]);
2520	}
2521	omFreeSize((ADDRESS)fr,(result->length)*sizeof(ideal));
2522	return result;
2523	}
2524	// from weight.cc
2525	BOOLEAN kWeight(leftv res,leftv id)
2526	{
2527	ideal F=(ideal)id->Data();
2528	intvec * iv = new intvec(pVariables);
2529	polyset s;
2530	int sl, n, i;
2531	int *x;
2532
2533	res->data=(char *)iv;
2534	s = F->m;
2535	sl = IDELEMS(F) - 1;
2536	n = pVariables;
2537	double wNsqr = (double)2.0 / (double)n;
2538	wFunctional = wFunctionalBuch;
2539	x = (int * )omAlloc(2 * (n + 1) * sizeof(int));
2540	wCall(s, sl, x, wNsqr);
2541	for (i = n; i!=0; i--)
2542	(*iv)[i-1] = x[i + n + 1];
2543	omFreeSize((ADDRESS)x, 2 * (n + 1) * sizeof(int));
2544	return FALSE;
2545	}
2546
2547	BOOLEAN kQHWeight(leftv res,leftv v)
2548	{
2549	res->data=(char *)idQHomWeight((ideal)v->Data());
2550	if (res->data==NULL)
2551	res->data=(char *)new intvec(pVariables);
2552	return FALSE;
2553	}
2554	/==============================================================/
2555	// from clapsing.cc
2556	#if 0
2557	BOOLEAN jjIS_SQR_FREE(leftv res, leftv u)
2558	{
2559	BOOLEAN b=singclap_factorize((poly)(u->CopyD()), &v, 0);
2560	res->data=(void *)b;
2561	}
2562	#endif
2563
2564	BOOLEAN jjRESULTANT(leftv res, leftv u, leftv v, leftv w)
2565	{
2566	res->data=singclap_resultant((poly)u->Data(),(poly)v->Data(), (poly)w->Data());
2567	return errorreported;
2568	}
2569	BOOLEAN jjCHARSERIES(leftv res, leftv u)
2570	{
2571	res->data=singclap_irrCharSeries((ideal)u->Data());
2572	return (res->data==NULL);
2573	}
2574
2575	// from semic.cc
2576	#ifdef HAVE_SPECTRUM
2577
2578	// ----------------------------------------------------------------------------
2579	// Initialize a spectrum deep from another spectrum
2580	// ----------------------------------------------------------------------------
2581
2582	void spectrum::copy_deep( const spectrum &spec )
2583	{
2584	mu = spec.mu;
2585	pg = spec.pg;
2586	n = spec.n;
2587
2588	copy_new( n );
2589
2590	for( int i=0; i<n; i++ )
2591	{
2592	s[i] = spec.s[i];
2593	w[i] = spec.w[i];
2594	}
2595	}
2596
2597	// ----------------------------------------------------------------------------
2598	// Initialize a spectrum deep from a singular lists
2599	// ----------------------------------------------------------------------------
2600
2601	void spectrum::copy_deep( lists l )
2602	{
2603	mu = (int)(long)(l->m[0].Data( ));
2604	pg = (int)(long)(l->m[1].Data( ));
2605	n = (int)(long)(l->m[2].Data( ));
2606
2607	copy_new( n );
2608
2609	intvec num = (intvec)l->m[3].Data( );
2610	intvec den = (intvec)l->m[4].Data( );
2611	intvec mul = (intvec)l->m[5].Data( );
2612
2613	for( int i=0; i<n; i++ )
2614	{
2615	s[i] = (Rational)((num)[i])/(Rational)((den)[i]);
2616	w[i] = (*mul)[i];
2617	}
2618	}
2619
2620	// ----------------------------------------------------------------------------
2621	// singular lists constructor for spectrum
2622	// ----------------------------------------------------------------------------
2623
2624	spectrum::spectrum( lists l )
2625	{
2626	copy_deep( l );
2627	}
2628
2629	// ----------------------------------------------------------------------------
2630	// generate a Singular lists from a spectrum
2631	// ----------------------------------------------------------------------------
2632
2633	lists spectrum::thelist( void )
2634	{
2635	lists L = (lists)omAllocBin( slists_bin);
2636
2637	L->Init( 6 );
2638
2639	intvec *num = new intvec( n );
2640	intvec *den = new intvec( n );
2641	intvec *mult = new intvec( n );
2642
2643	for( int i=0; i<n; i++ )
2644	{
2645	(*num) [i] = s[i].get_num_si( );
2646	(*den) [i] = s[i].get_den_si( );
2647	(*mult)[i] = w[i];
2648	}
2649
2650	L->m[0].rtyp = INT_CMD; // milnor number
2651	L->m[1].rtyp = INT_CMD; // geometrical genus
2652	L->m[2].rtyp = INT_CMD; // # of spectrum numbers
2653	L->m[3].rtyp = INTVEC_CMD; // numerators
2654	L->m[4].rtyp = INTVEC_CMD; // denomiantors
2655	L->m[5].rtyp = INTVEC_CMD; // multiplicities
2656
2657	L->m[0].data = (void*)mu;
2658	L->m[1].data = (void*)pg;
2659	L->m[2].data = (void*)n;
2660	L->m[3].data = (void*)num;
2661	L->m[4].data = (void*)den;
2662	L->m[5].data = (void*)mult;
2663
2664	return L;
2665	}
2666	// from spectrum.cc
2667	// ----------------------------------------------------------------------------
2668	// print out an error message for a spectrum list
2669	// ----------------------------------------------------------------------------
2670
2671	void list_error( semicState state )
2672	{
2673	switch( state )
2674	{
2675	case semicListTooShort:
2676	WerrorS( "the list is too short" );
2677	break;
2678	case semicListTooLong:
2679	WerrorS( "the list is too long" );
2680	break;
2681
2682	case semicListFirstElementWrongType:
2683	WerrorS( "first element of the list should be int" );
2684	break;
2685	case semicListSecondElementWrongType:
2686	WerrorS( "second element of the list should be int" );
2687	break;
2688	case semicListThirdElementWrongType:
2689	WerrorS( "third element of the list should be int" );
2690	break;
2691	case semicListFourthElementWrongType:
2692	WerrorS( "fourth element of the list should be intvec" );
2693	break;
2694	case semicListFifthElementWrongType:
2695	WerrorS( "fifth element of the list should be intvec" );
2696	break;
2697	case semicListSixthElementWrongType:
2698	WerrorS( "sixth element of the list should be intvec" );
2699	break;
2700
2701	case semicListNNegative:
2702	WerrorS( "first element of the list should be positive" );
2703	break;
2704	case semicListWrongNumberOfNumerators:
2705	WerrorS( "wrong number of numerators" );
2706	break;
2707	case semicListWrongNumberOfDenominators:
2708	WerrorS( "wrong number of denominators" );
2709	break;
2710	case semicListWrongNumberOfMultiplicities:
2711	WerrorS( "wrong number of multiplicities" );
2712	break;
2713
2714	case semicListMuNegative:
2715	WerrorS( "the Milnor number should be positive" );
2716	break;
2717	case semicListPgNegative:
2718	WerrorS( "the geometrical genus should be nonnegative" );
2719	break;
2720	case semicListNumNegative:
2721	WerrorS( "all numerators should be positive" );
2722	break;
2723	case semicListDenNegative:
2724	WerrorS( "all denominators should be positive" );
2725	break;
2726	case semicListMulNegative:
2727	WerrorS( "all multiplicities should be positive" );
2728	break;
2729
2730	case semicListNotSymmetric:
2731	WerrorS( "it is not symmetric" );
2732	break;
2733	case semicListNotMonotonous:
2734	WerrorS( "it is not monotonous" );
2735	break;
2736
2737	case semicListMilnorWrong:
2738	WerrorS( "the Milnor number is wrong" );
2739	break;
2740	case semicListPGWrong:
2741	WerrorS( "the geometrical genus is wrong" );
2742	break;
2743
2744	default:
2745	WerrorS( "unspecific error" );
2746	break;
2747	}
2748	}
2749	// ----------------------------------------------------------------------------
2750	// this is the main spectrum computation function
2751	// ----------------------------------------------------------------------------
2752
2753	spectrumState spectrumCompute( poly h,lists *L,int fast )
2754	{
2755	int i,j;
2756
2757	#ifdef SPECTRUM_DEBUG
2758	#ifdef SPECTRUM_PRINT
2759	#ifdef SPECTRUM_IOSTREAM
2760	cout << "spectrumCompute\n";
2761	if( fast==0 ) cout << " no optimization" << endl;
2762	if( fast==1 ) cout << " weight optimization" << endl;
2763	if( fast==2 ) cout << " symmetry optimization" << endl;
2764	#else
2765	fprintf( stdout,"spectrumCompute\n" );
2766	if( fast==0 ) fprintf( stdout," no optimization\n" );
2767	if( fast==1 ) fprintf( stdout," weight optimization\n" );
2768	if( fast==2 ) fprintf( stdout," symmetry optimization\n" );
2769	#endif
2770	#endif
2771	#endif
2772
2773	// ----------------------
2774	// check if h is zero
2775	// ----------------------
2776
2777	if( h==(poly)NULL )
2778	{
2779	return spectrumZero;
2780	}
2781
2782	// ----------------------------------
2783	// check if h has a constant term
2784	// ----------------------------------
2785
2786	if( hasConstTerm( h ) )
2787	{
2788	return spectrumBadPoly;
2789	}
2790
2791	// --------------------------------
2792	// check if h has a linear term
2793	// --------------------------------
2794
2795	if( hasLinearTerm( h ) )
2796	{
2797	*L = (lists)omAllocBin( slists_bin);
2798	(*L)->Init( 1 );
2799	(*L)->m[0].rtyp = INT_CMD; // milnor number
2800	/* (L)->m[0].data = (void)0;a -- done by Init */
2801
2802	return spectrumNoSingularity;
2803	}
2804
2805	// ----------------------------------
2806	// compute the jacobi ideal of (h)
2807	// ----------------------------------
2808
2809	ideal J = NULL;
2810	J = idInit( pVariables,1 );
2811
2812	#ifdef SPECTRUM_DEBUG
2813	#ifdef SPECTRUM_PRINT
2814	#ifdef SPECTRUM_IOSTREAM
2815	cout << "\n computing the Jacobi ideal...\n";
2816	#else
2817	fprintf( stdout,"\n computing the Jacobi ideal...\n" );
2818	#endif
2819	#endif
2820	#endif
2821
2822	for( i=0; i<pVariables; i++ )
2823	{
2824	J->m[i] = pDiff( h,i+1); //j );
2825
2826	#ifdef SPECTRUM_DEBUG
2827	#ifdef SPECTRUM_PRINT
2828	#ifdef SPECTRUM_IOSTREAM
2829	cout << " ";
2830	#else
2831	fprintf( stdout," " );
2832	#endif
2833	pWrite( J->m[i] );
2834	#endif
2835	#endif
2836	}
2837
2838	// --------------------------------------------
2839	// compute a standard basis stdJ of jac(h)
2840	// --------------------------------------------
2841
2842	#ifdef SPECTRUM_DEBUG
2843	#ifdef SPECTRUM_PRINT
2844	#ifdef SPECTRUM_IOSTREAM
2845	cout << endl;
2846	cout << " computing a standard basis..." << endl;
2847	#else
2848	fprintf( stdout,"\n" );
2849	fprintf( stdout," computing a standard basis...\n" );
2850	#endif
2851	#endif
2852	#endif
2853
2854	ideal stdJ = kStd(J,currQuotient,isNotHomog,NULL);
2855	idSkipZeroes( stdJ );
2856
2857	#ifdef SPECTRUM_DEBUG
2858	#ifdef SPECTRUM_PRINT
2859	for( i=0; i<IDELEMS(stdJ); i++ )
2860	{
2861	#ifdef SPECTRUM_IOSTREAM
2862	cout << " ";
2863	#else
2864	fprintf( stdout," " );
2865	#endif
2866
2867	pWrite( stdJ->m[i] );
2868	}
2869	#endif
2870	#endif
2871
2872	idDelete( &J );
2873
2874	// ------------------------------------------
2875	// check if the h has a singularity
2876	// ------------------------------------------
2877
2878	if( hasOne( stdJ ) )
2879	{
2880	// -------------------------------
2881	// h is smooth in the origin
2882	// return only the Milnor number
2883	// -------------------------------
2884
2885	*L = (lists)omAllocBin( slists_bin);
2886	(*L)->Init( 1 );
2887	(*L)->m[0].rtyp = INT_CMD; // milnor number
2888	/* (L)->m[0].data = (void)0;a -- done by Init */
2889
2890	return spectrumNoSingularity;
2891	}
2892
2893	// ------------------------------------------
2894	// check if the singularity h is isolated
2895	// ------------------------------------------
2896
2897	for( i=pVariables; i>0; i-- )
2898	{
2899	if( hasAxis( stdJ,i )==FALSE )
2900	{
2901	return spectrumNotIsolated;
2902	}
2903	}
2904
2905	// ------------------------------------------
2906	// compute the highest corner hc of stdJ
2907	// ------------------------------------------
2908
2909	#ifdef SPECTRUM_DEBUG
2910	#ifdef SPECTRUM_PRINT
2911	#ifdef SPECTRUM_IOSTREAM
2912	cout << "\n computing the highest corner...\n";
2913	#else
2914	fprintf( stdout,"\n computing the highest corner...\n" );
2915	#endif
2916	#endif
2917	#endif
2918
2919	poly hc = (poly)NULL;
2920
2921	scComputeHC( stdJ,currQuotient, 0,hc );
2922
2923	if( hc!=(poly)NULL )
2924	{
2925	pGetCoeff(hc) = nInit(1);
2926
2927	for( i=pVariables; i>0; i-- )
2928	{
2929	if( pGetExp( hc,i )>0 ) pDecrExp( hc,i );
2930	}
2931	pSetm( hc );
2932	}
2933	else
2934	{
2935	return spectrumNoHC;
2936	}
2937
2938	#ifdef SPECTRUM_DEBUG
2939	#ifdef SPECTRUM_PRINT
2940	#ifdef SPECTRUM_IOSTREAM
2941	cout << " ";
2942	#else
2943	fprintf( stdout," " );
2944	#endif
2945	pWrite( hc );
2946	#endif
2947	#endif
2948
2949	// ----------------------------------------
2950	// compute the Newton polygon nph of h
2951	// ----------------------------------------
2952
2953	#ifdef SPECTRUM_DEBUG
2954	#ifdef SPECTRUM_PRINT
2955	#ifdef SPECTRUM_IOSTREAM
2956	cout << "\n computing the newton polygon...\n";
2957	#else
2958	fprintf( stdout,"\n computing the newton polygon...\n" );
2959	#endif
2960	#endif
2961	#endif
2962
2963	newtonPolygon nph( h );
2964
2965	#ifdef SPECTRUM_DEBUG
2966	#ifdef SPECTRUM_PRINT
2967	cout << nph;
2968	#endif
2969	#endif
2970
2971	// -----------------------------------------------
2972	// compute the weight corner wc of (stdj,nph)
2973	// -----------------------------------------------
2974
2975	#ifdef SPECTRUM_DEBUG
2976	#ifdef SPECTRUM_PRINT
2977	#ifdef SPECTRUM_IOSTREAM
2978	cout << "\n computing the weight corner...\n";
2979	#else
2980	fprintf( stdout,"\n computing the weight corner...\n" );
2981	#endif
2982	#endif
2983	#endif
2984
2985	poly wc = ( fast==0 ? pCopy( hc ) :
2986	( fast==1 ? computeWC( nph,(Rational)pVariables ) :
2987	/* fast==2 */computeWC( nph,((Rational)pVariables)/(Rational)2 ) ) );
2988
2989	#ifdef SPECTRUM_DEBUG
2990	#ifdef SPECTRUM_PRINT
2991	#ifdef SPECTRUM_IOSTREAM
2992	cout << " ";
2993	#else
2994	fprintf( stdout," " );
2995	#endif
2996	pWrite( wc );
2997	#endif
2998	#endif
2999
3000	// -------------
3001	// compute NF
3002	// -------------
3003
3004	#ifdef SPECTRUM_DEBUG
3005	#ifdef SPECTRUM_PRINT
3006	#ifdef SPECTRUM_IOSTREAM
3007	cout << "\n computing NF...\n" << endl;
3008	#else
3009	fprintf( stdout,"\n computing NF...\n" );
3010	#endif
3011	#endif
3012	#endif
3013
3014	spectrumPolyList NF( &nph );
3015
3016	computeNF( stdJ,hc,wc,&NF );
3017
3018	#ifdef SPECTRUM_DEBUG
3019	#ifdef SPECTRUM_PRINT
3020	cout << NF;
3021	#ifdef SPECTRUM_IOSTREAM
3022	cout << endl;
3023	#else
3024	fprintf( stdout,"\n" );
3025	#endif
3026	#endif
3027	#endif
3028
3029	// ----------------------------
3030	// compute the spectrum of h
3031	// ----------------------------
3032
3033	return NF.spectrum( L,fast );
3034	}
3035
3036	// ----------------------------------------------------------------------------
3037	// this procedure is called from the interpreter
3038	// ----------------------------------------------------------------------------
3039	// first = polynomial
3040	// result = list of spectrum numbers
3041	// ----------------------------------------------------------------------------
3042
3043	BOOLEAN spectrumProc( leftv result,leftv first )
3044	{
3045	spectrumState state = spectrumOK;
3046
3047	// -------------------
3048	// check consistency
3049	// -------------------
3050
3051	// check for a local ring
3052
3053	if( !ringIsLocal( ) )
3054	{
3055	WerrorS( "only works for local orderings" );
3056	state = spectrumWrongRing;
3057	}
3058
3059	// no quotient rings are allowed
3060
3061	else if( currRing->qideal != NULL )
3062	{
3063	WerrorS( "does not work in quotient rings" );
3064	state = spectrumWrongRing;
3065	}
3066	else
3067	{
3068	lists L = (lists)NULL;
3069	int flag = 1; // weight corner optimization is safe
3070
3071	state = spectrumCompute( (poly)first->Data( ),&L,flag );
3072
3073	if( state==spectrumOK )
3074	{
3075	result->rtyp = LIST_CMD;
3076	result->data = (char*)L;
3077	}
3078	else
3079	{
3080	spectrumPrintError(state);
3081	}
3082	}
3083
3084	return (state!=spectrumOK);
3085	}
3086
3087	// ----------------------------------------------------------------------------
3088	// this procedure is called from the interpreter
3089	// ----------------------------------------------------------------------------
3090	// first = polynomial
3091	// result = list of spectrum numbers
3092	// ----------------------------------------------------------------------------
3093
3094	BOOLEAN spectrumfProc( leftv result,leftv first )
3095	{
3096	spectrumState state = spectrumOK;
3097
3098	// -------------------
3099	// check consistency
3100	// -------------------
3101
3102	// check for a local polynomial ring
3103
3104	if( currRing->OrdSgn != -1 )
3105	// ?? HS: the test above is also true for k[x][[y]], k[[x]][y]
3106	// or should we use:
3107	//if( !ringIsLocal( ) )
3108	{
3109	WerrorS( "only works for local orderings" );
3110	state = spectrumWrongRing;
3111	}
3112	else if( currRing->qideal != NULL )
3113	{
3114	WerrorS( "does not work in quotient rings" );
3115	state = spectrumWrongRing;
3116	}
3117	else
3118	{
3119	lists L = (lists)NULL;
3120	int flag = 2; // symmetric optimization
3121
3122	state = spectrumCompute( (poly)first->Data( ),&L,flag );
3123
3124	if( state==spectrumOK )
3125	{
3126	result->rtyp = LIST_CMD;
3127	result->data = (char*)L;
3128	}
3129	else
3130	{
3131	spectrumPrintError(state);
3132	}
3133	}
3134
3135	return (state!=spectrumOK);
3136	}
3137
3138	// ----------------------------------------------------------------------------
3139	// check if a list is a spectrum
3140	// check for:
3141	// list has 6 elements
3142	// 1st element is int (mu=Milnor number)
3143	// 2nd element is int (pg=geometrical genus)
3144	// 3rd element is int (n =number of different spectrum numbers)
3145	// 4th element is intvec (num=numerators)
3146	// 5th element is intvec (den=denomiantors)
3147	// 6th element is intvec (mul=multiplicities)
3148	// exactly n numerators
3149	// exactly n denominators
3150	// exactly n multiplicities
3151	// mu>0
3152	// pg>=0
3153	// n>0
3154	// num>0
3155	// den>0
3156	// mul>0
3157	// symmetriy with respect to numberofvariables/2
3158	// monotony
3159	// mu = sum of all multiplicities
3160	// pg = sum of all multiplicities where num/den<=1
3161	// ----------------------------------------------------------------------------
3162
3163	semicState list_is_spectrum( lists l )
3164	{
3165	// -------------------
3166	// check list length
3167	// -------------------
3168
3169	if( l->nr < 5 )
3170	{
3171	return semicListTooShort;
3172	}
3173	else if( l->nr > 5 )
3174	{
3175	return semicListTooLong;
3176	}
3177
3178	// -------------
3179	// check types
3180	// -------------
3181
3182	if( l->m[0].rtyp != INT_CMD )
3183	{
3184	return semicListFirstElementWrongType;
3185	}
3186	else if( l->m[1].rtyp != INT_CMD )
3187	{
3188	return semicListSecondElementWrongType;
3189	}
3190	else if( l->m[2].rtyp != INT_CMD )
3191	{
3192	return semicListThirdElementWrongType;
3193	}
3194	else if( l->m[3].rtyp != INTVEC_CMD )
3195	{
3196	return semicListFourthElementWrongType;
3197	}
3198	else if( l->m[4].rtyp != INTVEC_CMD )
3199	{
3200	return semicListFifthElementWrongType;
3201	}
3202	else if( l->m[5].rtyp != INTVEC_CMD )
3203	{
3204	return semicListSixthElementWrongType;
3205	}
3206
3207	// -------------------------
3208	// check number of entries
3209	// -------------------------
3210
3211	int mu = (int)(long)(l->m[0].Data( ));
3212	int pg = (int)(long)(l->m[1].Data( ));
3213	int n = (int)(long)(l->m[2].Data( ));
3214
3215	if( n <= 0 )
3216	{
3217	return semicListNNegative;
3218	}
3219
3220	intvec num = (intvec)l->m[3].Data( );
3221	intvec den = (intvec)l->m[4].Data( );
3222	intvec mul = (intvec)l->m[5].Data( );
3223
3224	if( n != num->length( ) )
3225	{
3226	return semicListWrongNumberOfNumerators;
3227	}
3228	else if( n != den->length( ) )
3229	{
3230	return semicListWrongNumberOfDenominators;
3231	}
3232	else if( n != mul->length( ) )
3233	{
3234	return semicListWrongNumberOfMultiplicities;
3235	}
3236
3237	// --------
3238	// values
3239	// --------
3240
3241	if( mu <= 0 )
3242	{
3243	return semicListMuNegative;
3244	}
3245	if( pg < 0 )
3246	{
3247	return semicListPgNegative;
3248	}
3249
3250	int i;
3251
3252	for( i=0; i<n; i++ )
3253	{
3254	if( (*num)[i] <= 0 )
3255	{
3256	return semicListNumNegative;
3257	}
3258	if( (*den)[i] <= 0 )
3259	{
3260	return semicListDenNegative;
3261	}
3262	if( (*mul)[i] <= 0 )
3263	{
3264	return semicListMulNegative;
3265	}
3266	}
3267
3268	// ----------------
3269	// check symmetry
3270	// ----------------
3271
3272	int j;
3273
3274	for( i=0, j=n-1; i<=j; i++,j-- )
3275	{
3276	if( (num)[i] != pVariables((den)[i]) - (num)[j] \|\|
3277	(den)[i] != (den)[j] \|\|
3278	(mul)[i] != (mul)[j] )
3279	{
3280	return semicListNotSymmetric;
3281	}
3282	}
3283
3284	// ----------------
3285	// check monotony
3286	// ----------------
3287
3288	for( i=0, j=1; i<n/2; i++,j++ )
3289	{
3290	if( (num)[i](den)[j] >= (num)[j](den)[i] )
3291	{
3292	return semicListNotMonotonous;
3293	}
3294	}
3295
3296	// ---------------------
3297	// check Milnor number
3298	// ---------------------
3299
3300	for( mu=0, i=0; i<n; i++ )
3301	{
3302	mu += (*mul)[i];
3303	}
3304
3305	if( mu != (int)(long)(l->m[0].Data( )) )
3306	{
3307	return semicListMilnorWrong;
3308	}
3309
3310	// -------------------------
3311	// check geometrical genus
3312	// -------------------------
3313
3314	for( pg=0, i=0; i<n; i++ )
3315	{
3316	if( (num)[i]<=(den)[i] )
3317	{
3318	pg += (*mul)[i];
3319	}
3320	}
3321
3322	if( pg != (int)(long)(l->m[1].Data( )) )
3323	{
3324	return semicListPGWrong;
3325	}
3326
3327	return semicOK;
3328	}
3329
3330	// ----------------------------------------------------------------------------
3331	// this procedure is called from the interpreter
3332	// ----------------------------------------------------------------------------
3333	// first = list of spectrum numbers
3334	// second = list of spectrum numbers
3335	// result = sum of the two lists
3336	// ----------------------------------------------------------------------------
3337
3338	BOOLEAN spaddProc( leftv result,leftv first,leftv second )
3339	{
3340	semicState state;
3341
3342	// -----------------
3343	// check arguments
3344	// -----------------
3345
3346	lists l1 = (lists)first->Data( );
3347	lists l2 = (lists)second->Data( );
3348
3349	if( (state=list_is_spectrum( l1 )) != semicOK )
3350	{
3351	WerrorS( "first argument is not a spectrum:" );
3352	list_error( state );
3353	}
3354	else if( (state=list_is_spectrum( l2 )) != semicOK )
3355	{
3356	WerrorS( "second argument is not a spectrum:" );
3357	list_error( state );
3358	}
3359	else
3360	{
3361	spectrum s1( l1 );
3362	spectrum s2( l2 );
3363	spectrum sum( s1+s2 );
3364
3365	result->rtyp = LIST_CMD;
3366	result->data = (char*)(sum.thelist( ));
3367	}
3368
3369	return (state!=semicOK);
3370	}
3371
3372	// ----------------------------------------------------------------------------
3373	// this procedure is called from the interpreter
3374	// ----------------------------------------------------------------------------
3375	// first = list of spectrum numbers
3376	// second = integer
3377	// result = the multiple of the first list by the second factor
3378	// ----------------------------------------------------------------------------
3379
3380	BOOLEAN spmulProc( leftv result,leftv first,leftv second )
3381	{
3382	semicState state;
3383
3384	// -----------------
3385	// check arguments
3386	// -----------------
3387
3388	lists l = (lists)first->Data( );
3389	int k = (int)(long)second->Data( );
3390
3391	if( (state=list_is_spectrum( l ))!=semicOK )
3392	{
3393	WerrorS( "first argument is not a spectrum" );
3394	list_error( state );
3395	}
3396	else if( k < 0 )
3397	{
3398	WerrorS( "second argument should be positive" );
3399	state = semicMulNegative;
3400	}
3401	else
3402	{
3403	spectrum s( l );
3404	spectrum product( k*s );
3405
3406	result->rtyp = LIST_CMD;
3407	result->data = (char*)product.thelist( );
3408	}
3409
3410	return (state!=semicOK);
3411	}
3412
3413	// ----------------------------------------------------------------------------
3414	// this procedure is called from the interpreter
3415	// ----------------------------------------------------------------------------
3416	// first = list of spectrum numbers
3417	// second = list of spectrum numbers
3418	// result = semicontinuity index
3419	// ----------------------------------------------------------------------------
3420
3421	BOOLEAN semicProc3 ( leftv res,leftv u,leftv v,leftv w )
3422	{
3423	semicState state;
3424	BOOLEAN qh=(((int)(long)w->Data())==1);
3425
3426	// -----------------
3427	// check arguments
3428	// -----------------
3429
3430	lists l1 = (lists)u->Data( );
3431	lists l2 = (lists)v->Data( );
3432
3433	if( (state=list_is_spectrum( l1 ))!=semicOK )
3434	{
3435	WerrorS( "first argument is not a spectrum" );
3436	list_error( state );
3437	}
3438	else if( (state=list_is_spectrum( l2 ))!=semicOK )
3439	{
3440	WerrorS( "second argument is not a spectrum" );
3441	list_error( state );
3442	}
3443	else
3444	{
3445	spectrum s1( l1 );
3446	spectrum s2( l2 );
3447
3448	res->rtyp = INT_CMD;
3449	if (qh)
3450	res->data = (void*)(s1.mult_spectrumh( s2 ));
3451	else
3452	res->data = (void*)(s1.mult_spectrum( s2 ));
3453	}
3454
3455	// -----------------
3456	// check status
3457	// -----------------
3458
3459	return (state!=semicOK);
3460	}
3461	BOOLEAN semicProc ( leftv res,leftv u,leftv v )
3462	{
3463	sleftv tmp;
3464	memset(&tmp,0,sizeof(tmp));
3465	tmp.rtyp=INT_CMD;
3466	/* tmp.data = (void )0; -- done by memset /
3467
3468	return semicProc3(res,u,v,&tmp);
3469	}
3470	// from splist.cc
3471	// ----------------------------------------------------------------------------
3472	// Compute the spectrum of a spectrumPolyList
3473	// ----------------------------------------------------------------------------
3474
3475	spectrumState spectrumPolyList::spectrum( lists *L,int fast )
3476	{
3477	spectrumPolyNode **node = &root;
3478	spectrumPolyNode *search;
3479
3480	poly f,tmp;
3481	int found,cmp;
3482
3483	Rational smax( ( fast==0 ? 0 : pVariables ),
3484	( fast==2 ? 2 : 1 ) );
3485
3486	Rational weight_prev( 0,1 );
3487
3488	int mu = 0; // the milnor number
3489	int pg = 0; // the geometrical genus
3490	int n = 0; // number of different spectral numbers
3491	int z = 0; // number of spectral number equal to smax
3492
3493	int k = 0;
3494
3495	while( (node)!=(spectrumPolyNode)NULL &&
3496	( fast==0 \|\| (*node)->weight<=smax ) )
3497	{
3498	// ---------------------------------------
3499	// determine the first normal form which
3500	// contains the monomial node->mon
3501	// ---------------------------------------
3502
3503	found = FALSE;
3504	search = *node;
3505
3506	while( search!=(spectrumPolyNode*)NULL && found==FALSE )
3507	{
3508	if( search->nf!=(poly)NULL )
3509	{
3510	f = search->nf;
3511
3512	do
3513	{
3514	// --------------------------------
3515	// look for (*node)->mon in f
3516	// --------------------------------
3517
3518	cmp = pCmp( (*node)->mon,f );
3519
3520	if( cmp<0 )
3521	{
3522	f = pNext( f );
3523	}
3524	else if( cmp==0 )
3525	{
3526	// -----------------------------
3527	// we have found a normal form
3528	// -----------------------------
3529
3530	found = TRUE;
3531
3532	// normalize coefficient
3533
3534	number inv = nInvers( pGetCoeff( f ) );
3535	pMult_nn( search->nf,inv );
3536	nDelete( &inv );
3537
3538	// exchange normal forms
3539
3540	tmp = (*node)->nf;
3541	(*node)->nf = search->nf;
3542	search->nf = tmp;
3543	}
3544	}
3545	while( cmp<0 && f!=(poly)NULL );
3546	}
3547	search = search->next;
3548	}
3549
3550	if( found==FALSE )
3551	{
3552	// ------------------------------------------------
3553	// the weight of node->mon is a spectrum number
3554	// ------------------------------------------------
3555
3556	mu++;
3557
3558	if( (*node)->weight<=(Rational)1 ) pg++;
3559	if( (*node)->weight==smax ) z++;
3560	if( (*node)->weight>weight_prev ) n++;
3561
3562	weight_prev = (*node)->weight;
3563	node = &((*node)->next);
3564	}
3565	else
3566	{
3567	// -----------------------------------------------
3568	// determine all other normal form which contain
3569	// the monomial node->mon
3570	// replace for node->mon its normal form
3571	// -----------------------------------------------
3572
3573	while( search!=(spectrumPolyNode*)NULL )
3574	{
3575	if( search->nf!=(poly)NULL )
3576	{
3577	f = search->nf;
3578
3579	do
3580	{
3581	// --------------------------------
3582	// look for (*node)->mon in f
3583	// --------------------------------
3584
3585	cmp = pCmp( (*node)->mon,f );
3586
3587	if( cmp<0 )
3588	{
3589	f = pNext( f );
3590	}
3591	else if( cmp==0 )
3592	{
3593	search->nf = pSub( search->nf,
3594	ppMult_nn( (*node)->nf,pGetCoeff( f ) ) );
3595	pNorm( search->nf );
3596	}
3597	}
3598	while( cmp<0 && f!=(poly)NULL );
3599	}
3600	search = search->next;
3601	}
3602	delete_node( node );
3603	}
3604
3605	}
3606
3607	// --------------------------------------------------------
3608	// fast computation exploits the symmetry of the spectrum
3609	// --------------------------------------------------------
3610
3611	if( fast==2 )
3612	{
3613	mu = 2*mu - z;
3614	n = ( z > 0 ? 2n - 1 : 2n );
3615	}
3616
3617	// --------------------------------------------------------
3618	// compute the spectrum numbers with their multiplicities
3619	// --------------------------------------------------------
3620
3621	intvec *nom = new intvec( n );
3622	intvec *den = new intvec( n );
3623	intvec *mult = new intvec( n );
3624
3625	int count = 0;
3626	int multiplicity = 1;
3627
3628	for( search=root; search!=(spectrumPolyNode*)NULL &&
3629	( fast==0 \|\| search->weight<=smax );
3630	search=search->next )
3631	{
3632	if( search->next==(spectrumPolyNode*)NULL \|\|
3633	search->weight<search->next->weight )
3634	{
3635	(*nom) [count] = search->weight.get_num_si( );
3636	(*den) [count] = search->weight.get_den_si( );
3637	(*mult)[count] = multiplicity;
3638
3639	multiplicity=1;
3640	count++;
3641	}
3642	else
3643	{
3644	multiplicity++;
3645	}
3646	}
3647
3648	// --------------------------------------------------------
3649	// fast computation exploits the symmetry of the spectrum
3650	// --------------------------------------------------------
3651
3652	if( fast==2 )
3653	{
3654	int n1,n2;
3655	for( n1=0, n2=n-1; n1<n2; n1++, n2-- )
3656	{
3657	(nom) [n2] = pVariables(den)[n1]-(nom)[n1];
3658	(den) [n2] = (den)[n1];
3659	(mult)[n2] = (mult)[n1];
3660	}
3661	}
3662
3663	// -----------------------------------
3664	// test if the spectrum is symmetric
3665	// -----------------------------------
3666
3667	if( fast==0 \|\| fast==1 )
3668	{
3669	int symmetric=TRUE;
3670
3671	for( int n1=0, n2=n-1 ; n1<n2 && symmetric==TRUE; n1++, n2-- )
3672	{
3673	if( (mult)[n1]!=(mult)[n2] \|\|
3674	(den) [n1]!= (den)[n2] \|\|
3675	(nom)[n1]+(nom)[n2]!=pVariables(den) [n1] )
3676	{
3677	symmetric = FALSE;
3678	}
3679	}
3680
3681	if( symmetric==FALSE )
3682	{
3683	// ---------------------------------------------
3684	// the spectrum is not symmetric => degenerate
3685	// principal part
3686	// ---------------------------------------------
3687
3688	*L = (lists)omAllocBin( slists_bin);
3689	(*L)->Init( 1 );
3690	(*L)->m[0].rtyp = INT_CMD; // milnor number
3691	(L)->m[0].data = (void)mu;
3692
3693	return spectrumDegenerate;
3694	}
3695	}
3696
3697	*L = (lists)omAllocBin( slists_bin);
3698
3699	(*L)->Init( 6 );
3700
3701	(*L)->m[0].rtyp = INT_CMD; // milnor number
3702	(*L)->m[1].rtyp = INT_CMD; // geometrical genus
3703	(*L)->m[2].rtyp = INT_CMD; // number of spectrum values
3704	(*L)->m[3].rtyp = INTVEC_CMD; // nominators
3705	(*L)->m[4].rtyp = INTVEC_CMD; // denomiantors
3706	(*L)->m[5].rtyp = INTVEC_CMD; // multiplicities
3707
3708	(L)->m[0].data = (void)mu;
3709	(L)->m[1].data = (void)pg;
3710	(L)->m[2].data = (void)n;
3711	(L)->m[3].data = (void)nom;
3712	(L)->m[4].data = (void)den;
3713	(L)->m[5].data = (void)mult;
3714
3715	return spectrumOK;
3716	}
3717
3718	#endif
3719
3720	//from mpr_inout.cc
3721	extern void nPrint(number n);
3722
3723	BOOLEAN loNewtonP( leftv res, leftv arg1 )
3724	{
3725	res->data= (void*)loNewtonPolytope( (ideal)arg1->Data() );
3726	return FALSE;
3727	}
3728
3729	BOOLEAN loSimplex( leftv res, leftv args )
3730	{
3731	if ( !(rField_is_long_R()) )
3732	{
3733	WerrorS("Ground field not implemented!");
3734	return TRUE;
3735	}
3736
3737	simplex * LP;
3738	matrix m;
3739
3740	leftv v= args;
3741	if ( v->Typ() != MATRIX_CMD ) // 1: matrix
3742	return TRUE;
3743	else
3744	m= (matrix)(v->CopyD());
3745
3746	LP = new simplex(MATROWS(m),MATCOLS(m));
3747	LP->mapFromMatrix(m);
3748
3749	v= v->next;
3750	if ( v->Typ() != INT_CMD ) // 2: m = number of constraints
3751	return TRUE;
3752	else
3753	LP->m= (int)(long)(v->Data());
3754
3755	v= v->next;
3756	if ( v->Typ() != INT_CMD ) // 3: n = number of variables
3757	return TRUE;
3758	else
3759	LP->n= (int)(long)(v->Data());
3760
3761	v= v->next;
3762	if ( v->Typ() != INT_CMD ) // 4: m1 = number of <= constraints
3763	return TRUE;
3764	else
3765	LP->m1= (int)(long)(v->Data());
3766
3767	v= v->next;
3768	if ( v->Typ() != INT_CMD ) // 5: m2 = number of >= constraints
3769	return TRUE;
3770	else
3771	LP->m2= (int)(long)(v->Data());
3772
3773	v= v->next;
3774	if ( v->Typ() != INT_CMD ) // 6: m3 = number of == constraints
3775	return TRUE;
3776	else
3777	LP->m3= (int)(long)(v->Data());
3778
3779	#ifdef mprDEBUG_PROT
3780	Print("m (constraints) %d\n",LP->m);
3781	Print("n (columns) %d\n",LP->n);
3782	Print("m1 (<=) %d\n",LP->m1);
3783	Print("m2 (>=) %d\n",LP->m2);
3784	Print("m3 (==) %d\n",LP->m3);
3785	#endif
3786
3787	LP->compute();
3788
3789	lists lres= (lists)omAlloc( sizeof(slists) );
3790	lres->Init( 6 );
3791
3792	lres->m[0].rtyp= MATRIX_CMD; // output matrix
3793	lres->m[0].data=(void*)LP->mapToMatrix(m);
3794
3795	lres->m[1].rtyp= INT_CMD; // found a solution?
3796	lres->m[1].data=(void*)LP->icase;
3797
3798	lres->m[2].rtyp= INTVEC_CMD;
3799	lres->m[2].data=(void*)LP->posvToIV();
3800
3801	lres->m[3].rtyp= INTVEC_CMD;
3802	lres->m[3].data=(void*)LP->zrovToIV();
3803
3804	lres->m[4].rtyp= INT_CMD;
3805	lres->m[4].data=(void*)LP->m;
3806
3807	lres->m[5].rtyp= INT_CMD;
3808	lres->m[5].data=(void*)LP->n;
3809
3810	res->data= (void*)lres;
3811
3812	return FALSE;
3813	}
3814
3815	BOOLEAN nuMPResMat( leftv res, leftv arg1, leftv arg2 )
3816	{
3817	ideal gls = (ideal)(arg1->Data());
3818	int imtype= (int)(long)arg2->Data();
3819
3820	uResultant::resMatType mtype= determineMType( imtype );
3821
3822	// check input ideal ( = polynomial system )
3823	if ( mprIdealCheck( gls, arg1->Name(), mtype, true ) != mprOk )
3824	{
3825	return TRUE;
3826	}
3827
3828	uResultant *resMat= new uResultant( gls, mtype, false );
3829	if (resMat!=NULL)
3830	{
3831	res->rtyp = MODUL_CMD;
3832	res->data= (void*)resMat->accessResMat()->getMatrix();
3833	if (!errorreported) delete resMat;
3834	}
3835	return errorreported;
3836	}
3837
3838	BOOLEAN nuLagSolve( leftv res, leftv arg1, leftv arg2, leftv arg3 )
3839	{
3840
3841	poly gls;
3842	gls= (poly)(arg1->Data());
3843	int howclean= (int)(long)arg3->Data();
3844
3845	if ( !(rField_is_R() \|\|
3846	rField_is_Q() \|\|
3847	rField_is_long_R() \|\|
3848	rField_is_long_C()) )
3849	{
3850	WerrorS("Ground field not implemented!");
3851	return TRUE;
3852	}
3853
3854	if ( !(rField_is_R()\|\|rField_is_long_R()\|\|rField_is_long_C()) )
3855	{
3856	unsigned long int ii = (unsigned long int)arg2->Data();
3857	setGMPFloatDigits( ii, ii );
3858	}
3859
3860	if ( gls == NULL \|\| pIsConstant( gls ) )
3861	{
3862	WerrorS("Input polynomial is constant!");
3863	return TRUE;
3864	}
3865
3866	int ldummy;
3867	int deg= pLDeg( gls, &ldummy, currRing );
3868	// int deg= pDeg( gls );
3869	int len= pLength( gls );
3870	int i,vpos=0;
3871	poly piter;
3872	lists elist;
3873	lists rlist;
3874
3875	elist= (lists)omAlloc( sizeof(slists) );
3876	elist->Init( 0 );
3877
3878	if ( pVariables > 1 )
3879	{
3880	piter= gls;
3881	for ( i= 1; i <= pVariables; i++ )
3882	if ( pGetExp( piter, i ) )
3883	{
3884	vpos= i;
3885	break;
3886	}
3887	while ( piter )
3888	{
3889	for ( i= 1; i <= pVariables; i++ )
3890	if ( (vpos != i) && (pGetExp( piter, i ) != 0) )
3891	{
3892	WerrorS("The input polynomial must be univariate!");
3893	return TRUE;
3894	}
3895	pIter( piter );
3896	}
3897	}
3898
3899	rootContainer * roots= new rootContainer();
3900	number * pcoeffs= (number )omAlloc( (deg+1) sizeof( number ) );
3901	piter= gls;
3902	for ( i= deg; i >= 0; i-- )
3903	{
3904	//if ( piter ) Print("deg %d, pDeg(piter) %d\n",i,pTotaldegree(piter));
3905	if ( piter && pTotaldegree(piter) == i )
3906	{
3907	pcoeffs[i]= nCopy( pGetCoeff( piter ) );
3908	//nPrint( pcoeffs[i] );PrintS(" ");
3909	pIter( piter );
3910	}
3911	else
3912	{
3913	pcoeffs[i]= nInit(0);
3914	}
3915	}
3916
3917	#ifdef mprDEBUG_PROT
3918	for (i=deg; i >= 0; i--)
3919	{
3920	nPrint( pcoeffs[i] );PrintS(" ");
3921	}
3922	PrintLn();
3923	#endif
3924
3925	roots->fillContainer( pcoeffs, NULL, 1, deg, rootContainer::onepoly, 1 );
3926	roots->solver( howclean );
3927
3928	int elem= roots->getAnzRoots();
3929	char *out;
3930	char *dummy;
3931	int j;
3932
3933	rlist= (lists)omAlloc( sizeof(slists) );
3934	rlist->Init( elem );
3935
3936	if (rField_is_long_C())
3937	{
3938	for ( j= 0; j < elem; j++ )
3939	{
3940	rlist->m[j].rtyp=NUMBER_CMD;
3941	rlist->m[j].data=(void *)nCopy((number)(roots->getRoot(j)));
3942	//rlist->m[j].data=(void *)(number)(roots->getRoot(j));
3943	}
3944	}
3945	else
3946	{
3947	for ( j= 0; j < elem; j++ )
3948	{
3949	dummy = complexToStr( (*roots)[j], gmp_output_digits );
3950	rlist->m[j].rtyp=STRING_CMD;
3951	rlist->m[j].data=(void *)dummy;
3952	}
3953	}
3954
3955	elist->Clean();
3956	//omFreeSize( (ADDRESS) elist, sizeof(slists) );
3957
3958	// this is (via fillContainer) the same data as in root
3959	//for ( i= deg; i >= 0; i-- ) nDelete( &pcoeffs[i] );
3960	//omFreeSize( (ADDRESS) pcoeffs, (deg+1) * sizeof( number ) );
3961
3962	delete roots;
3963
3964	res->rtyp= LIST_CMD;
3965	res->data= (void*)rlist;
3966
3967	return FALSE;
3968	}
3969
3970	BOOLEAN nuVanderSys( leftv res, leftv arg1, leftv arg2, leftv arg3)
3971	{
3972	int i;
3973	ideal p,w;
3974	p= (ideal)arg1->Data();
3975	w= (ideal)arg2->Data();
3976
3977	// w[0] = f(p^0)
3978	// w[1] = f(p^1)
3979	// ...
3980	// p can be a vector of numbers (multivariate polynom)
3981	// or one number (univariate polynom)
3982	// tdg = deg(f)
3983
3984	int n= IDELEMS( p );
3985	int m= IDELEMS( w );
3986	int tdg= (int)(long)arg3->Data();
3987
3988	res->data= (void*)NULL;
3989
3990	// check the input
3991	if ( tdg < 1 )
3992	{
3993	WerrorS("Last input parameter must be > 0!");
3994	return TRUE;
3995	}
3996	if ( n != pVariables )
3997	{
3998	Werror("Size of first input ideal must be equal to %d!",pVariables);
3999	return TRUE;
4000	}
4001	if ( m != (int)pow((double)tdg+1,(double)n) )
4002	{
4003	Werror("Size of second input ideal must be equal to %d!",
4004	(int)pow((double)tdg+1,(double)n));
4005	return TRUE;
4006	}
4007	if ( !(rField_is_Q() /* \|\|
4008	rField_is_R() \|\| rField_is_long_R() \|\|
4009	rField_is_long_C()*/ ) )
4010	{
4011	WerrorS("Ground field not implemented!");
4012	return TRUE;
4013	}
4014
4015	number tmp;
4016	number pevpoint= (number )omAlloc( n * sizeof( number ) );
4017	for ( i= 0; i < n; i++ )
4018	{
4019	pevpoint[i]=nInit(0);
4020	if ( (p->m)[i] )
4021	{
4022	tmp = pGetCoeff( (p->m)[i] );
4023	if ( nIsZero(tmp) \|\| nIsOne(tmp) \|\| nIsMOne(tmp) )
4024	{
4025	omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4026	WerrorS("Elements of first input ideal must not be equal to -1, 0, 1!");
4027	return TRUE;
4028	}
4029	} else tmp= NULL;
4030	if ( !nIsZero(tmp) )
4031	{
4032	if ( !pIsConstant((p->m)[i]))
4033	{
4034	omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4035	WerrorS("Elements of first input ideal must be numbers!");
4036	return TRUE;
4037	}
4038	pevpoint[i]= nCopy( tmp );
4039	}
4040	}
4041
4042	number wresults= (number )omAlloc( m * sizeof( number ) );
4043	for ( i= 0; i < m; i++ )
4044	{
4045	wresults[i]= nInit(0);
4046	if ( (w->m)[i] && !nIsZero(pGetCoeff((w->m)[i])) )
4047	{
4048	if ( !pIsConstant((w->m)[i]))
4049	{
4050	omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4051	omFreeSize( (ADDRESS)wresults, m * sizeof( number ) );
4052	WerrorS("Elements of second input ideal must be numbers!");
4053	return TRUE;
4054	}
4055	wresults[i]= nCopy(pGetCoeff((w->m)[i]));
4056	}
4057	}
4058
4059	vandermonde vm( m, n, tdg, pevpoint, FALSE );
4060	number *ncpoly= vm.interpolateDense( wresults );
4061	// do not free ncpoly[]!!
4062	poly rpoly= vm.numvec2poly( ncpoly );
4063
4064	omFreeSize( (ADDRESS)pevpoint, n * sizeof( number ) );
4065	omFreeSize( (ADDRESS)wresults, m * sizeof( number ) );
4066
4067	res->data= (void*)rpoly;
4068	return FALSE;
4069	}
4070
4071	BOOLEAN nuUResSolve( leftv res, leftv args )
4072	{
4073	leftv v= args;
4074
4075	ideal gls;
4076	int imtype;
4077	int howclean;
4078
4079	// get ideal
4080	if ( v->Typ() != IDEAL_CMD )
4081	return TRUE;
4082	else gls= (ideal)(v->Data());
4083	v= v->next;
4084
4085	// get resultant matrix type to use (0,1)
4086	if ( v->Typ() != INT_CMD )
4087	return TRUE;
4088	else imtype= (int)(long)v->Data();
4089	v= v->next;
4090
4091	if (imtype==0)
4092	{
4093	ideal test_id=idInit(1,1);
4094	int j;
4095	for(j=IDELEMS(gls)-1;j>=0;j--)
4096	{
4097	if (gls->m[j]!=NULL)
4098	{
4099	test_id->m[0]=gls->m[j];
4100	intvec *dummy_w=idQHomWeight(test_id);
4101	if (dummy_w!=NULL)
4102	{
4103	WerrorS("Newton polytope not of expected dimension");
4104	delete dummy_w;
4105	return TRUE;
4106	}
4107	}
4108	}
4109	}
4110
4111	// get and set precision in digits ( > 0 )
4112	if ( v->Typ() != INT_CMD )
4113	return TRUE;
4114	else if ( !(rField_is_R()\|\|rField_is_long_R()\|\|rField_is_long_C()) )
4115	{
4116	unsigned long int ii=(unsigned long int)v->Data();
4117	setGMPFloatDigits( ii, ii );
4118	}
4119	v= v->next;
4120
4121	// get interpolation steps (0,1,2)
4122	if ( v->Typ() != INT_CMD )
4123	return TRUE;
4124	else howclean= (int)(long)v->Data();
4125
4126	uResultant::resMatType mtype= determineMType( imtype );
4127	int i,c,count;
4128	lists listofroots= NULL;
4129	lists emptylist;
4130	number smv= NULL;
4131	BOOLEAN interpolate_det= (mtype==uResultant::denseResMat)?TRUE:FALSE;
4132
4133	//emptylist= (lists)omAlloc( sizeof(slists) );
4134	//emptylist->Init( 0 );
4135
4136	//res->rtyp = LIST_CMD;
4137	//res->data= (void *)emptylist;
4138
4139	// check input ideal ( = polynomial system )
4140	if ( mprIdealCheck( gls, args->Name(), mtype ) != mprOk )
4141	{
4142	return TRUE;
4143	}
4144
4145	uResultant * ures;
4146	rootContainer ** iproots;
4147	rootContainer ** muiproots;
4148	rootArranger * arranger;
4149
4150	// main task 1: setup of resultant matrix
4151	ures= new uResultant( gls, mtype );
4152	if ( ures->accessResMat()->initState() != resMatrixBase::ready )
4153	{
4154	WerrorS("Error occurred during matrix setup!");
4155	return TRUE;
4156	}
4157
4158	// if dense resultant, check if minor nonsingular
4159	if ( mtype == uResultant::denseResMat )
4160	{
4161	smv= ures->accessResMat()->getSubDet();
4162	#ifdef mprDEBUG_PROT
4163	PrintS("// Determinant of submatrix: ");nPrint(smv);PrintLn();
4164	#endif
4165	if ( nIsZero(smv) )
4166	{
4167	WerrorS("Unsuitable input ideal: Minor of resultant matrix is singular!");
4168	return TRUE;
4169	}
4170	}
4171
4172	// main task 2: Interpolate specialized resultant polynomials
4173	if ( interpolate_det )
4174	iproots= ures->interpolateDenseSP( false, smv );
4175	else
4176	iproots= ures->specializeInU( false, smv );
4177
4178	// main task 3: Interpolate specialized resultant polynomials
4179	if ( interpolate_det )
4180	muiproots= ures->interpolateDenseSP( true, smv );
4181	else
4182	muiproots= ures->specializeInU( true, smv );
4183
4184	#ifdef mprDEBUG_PROT
4185	c= iproots[0]->getAnzElems();
4186	for (i=0; i < c; i++) pWrite(iproots[i]->getPoly());
4187	c= muiproots[0]->getAnzElems();
4188	for (i=0; i < c; i++) pWrite(muiproots[i]->getPoly());
4189	#endif
4190
4191	// main task 4: Compute roots of specialized polys and match them up
4192	arranger= new rootArranger( iproots, muiproots, howclean );
4193	arranger->solve_all();
4194
4195	// get list of roots
4196	if ( arranger->success() )
4197	{
4198	arranger->arrange();
4199	listofroots= arranger->listOfRoots( gmp_output_digits );
4200	}
4201	else
4202	{
4203	WerrorS("Solver was unable to find any roots!");
4204	return TRUE;
4205	}
4206
4207	// free everything
4208	count= iproots[0]->getAnzElems();
4209	for (i=0; i < count; i++) delete iproots[i];
4210	omFreeSize( (ADDRESS) iproots, count * sizeof(rootContainer*) );
4211	count= muiproots[0]->getAnzElems();
4212	for (i=0; i < count; i++) delete muiproots[i];
4213	omFreeSize( (ADDRESS) muiproots, count * sizeof(rootContainer*) );
4214
4215	delete ures;
4216	delete arranger;
4217	nDelete( &smv );
4218
4219	res->data= (void *)listofroots;
4220
4221	//emptylist->Clean();
4222	// omFreeSize( (ADDRESS) emptylist, sizeof(slists) );
4223
4224	return FALSE;
4225	}
4226
4227	// from mpr_numeric.cc
4228	lists rootArranger::listOfRoots( const unsigned int oprec )
4229	{
4230	int i,j,tr;
4231	int count= roots[0]->getAnzRoots(); // number of roots
4232	int elem= roots[0]->getAnzElems(); // number of koordinates per root
4233
4234	lists listofroots= (lists)omAlloc( sizeof(slists) ); // must be done this way!
4235
4236	if ( found_roots )
4237	{
4238	listofroots->Init( count );
4239
4240	for (i=0; i < count; i++)
4241	{
4242	lists onepoint= (lists)omAlloc(sizeof(slists)); // must be done this way!
4243	onepoint->Init(elem);
4244	for ( j= 0; j < elem; j++ )
4245	{
4246	if ( !rField_is_long_C() )
4247	{
4248	onepoint->m[j].rtyp=STRING_CMD;
4249	onepoint->m[j].data=(void )complexToStr((roots[j])[i],oprec);
4250	}
4251	else
4252	{
4253	onepoint->m[j].rtyp=NUMBER_CMD;
4254	onepoint->m[j].data=(void *)nCopy((number)(roots[j]->getRoot(i)));
4255	}
4256	onepoint->m[j].next= NULL;
4257	onepoint->m[j].name= NULL;
4258	}
4259	listofroots->m[i].rtyp=LIST_CMD;
4260	listofroots->m[i].data=(void *)onepoint;
4261	listofroots->m[j].next= NULL;
4262	listofroots->m[j].name= NULL;
4263	}
4264
4265	}
4266	else
4267	{
4268	listofroots->Init( 0 );
4269	}
4270
4271	return listofroots;
4272	}
4273
4274	// from ring.cc
4275	void rSetHdl(idhdl h)
4276	{
4277	int i;
4278	ring rg = NULL;
4279	if (h!=NULL)
4280	{
4281	// Print(" new ring:%s (l:%d)\n",IDID(h),IDLEV(h));
4282	rg = IDRING(h);
4283	if (rg==NULL) return; //id <>NULL, ring==NULL
4284	omCheckAddrSize((ADDRESS)h,sizeof(idrec));
4285	if (IDID(h)) // OB: ????
4286	omCheckAddr((ADDRESS)IDID(h));
4287	rTest(rg);
4288	}
4289
4290	// clean up history
4291	if (sLastPrinted.RingDependend())
4292	{
4293	sLastPrinted.CleanUp();
4294	memset(&sLastPrinted,0,sizeof(sleftv));
4295	}
4296
4297	// test for valid "currRing":
4298	if ((rg!=NULL) && (rg->idroot==NULL))
4299	{
4300	ring old=rg;
4301	rg=rAssure_HasComp(rg);
4302	if (old!=rg)
4303	{
4304	rKill(old);
4305	IDRING(h)=rg;
4306	}
4307	}
4308	/------------ change the global ring -----------------------/
4309	rChangeCurrRing(rg);
4310	currRingHdl = h;
4311	}
4312
4313	BOOLEAN rSleftvOrdering2Ordering(sleftv *ord, ring R)
4314	{
4315	int last = 0, o=0, n = 1, i=0, typ = 1, j;
4316	sleftv *sl = ord;
4317
4318	// determine nBlocks
4319	while (sl!=NULL)
4320	{
4321	intvec iv = (intvec )(sl->data);
4322	if (((iv)[1]==ringorder_c)\|\|((iv)[1]==ringorder_C))
4323	i++;
4324	else if ((*iv)[1]==ringorder_L)
4325	{
4326	R->bitmask=(*iv)[2];
4327	n--;
4328	}
4329	else if (((*iv)[1]!=ringorder_a)
4330	&& ((*iv)[1]!=ringorder_a64))
4331	o++;
4332	n++;
4333	sl=sl->next;
4334	}
4335	// check whether at least one real ordering
4336	if (o==0)
4337	{
4338	WerrorS("invalid combination of orderings");
4339	return TRUE;
4340	}
4341	// if no c/C ordering is given, increment n
4342	if (i==0) n++;
4343	else if (i != 1)
4344	{
4345	// throw error if more than one is given
4346	WerrorS("more than one ordering c/C specified");
4347	return TRUE;
4348	}
4349
4350	// initialize fields of R
4351	R->order=(int )omAlloc0(nsizeof(int));
4352	R->block0=(int )omAlloc0(nsizeof(int));
4353	R->block1=(int )omAlloc0(nsizeof(int));
4354	R->wvhdl=(int*)omAlloc0(nsizeof(int_ptr));
4355
4356	int weights=(int)omAlloc0((R->N+1)*sizeof(int));
4357
4358	// init order, so that rBlocks works correctly
4359	for (j=0; j < n-1; j++)
4360	R->order[j] = (int) ringorder_unspec;
4361	// set last _C order, if no c/C order was given
4362	if (i == 0) R->order[n-2] = ringorder_C;
4363
4364	/* init orders */
4365	sl=ord;
4366	n=-1;
4367	while (sl!=NULL)
4368	{
4369	intvec *iv;
4370	iv = (intvec *)(sl->data);
4371	if ((*iv)[1]!=ringorder_L)
4372	{
4373	n++;
4374
4375	/* the format of an ordering:
4376	* iv[0]: factor
4377	* iv[1]: ordering
4378	* iv[2..end]: weights
4379	*/
4380	R->order[n] = (*iv)[1];
4381	typ=1;
4382	switch ((*iv)[1])
4383	{
4384	case ringorder_ws:
4385	case ringorder_Ws:
4386	typ=-1;
4387	case ringorder_wp:
4388	case ringorder_Wp:
4389	R->wvhdl[n]=(int)omAlloc((iv->length()-1)sizeof(int));
4390	R->block0[n] = last+1;
4391	for (i=2; i<iv->length(); i++)
4392	{
4393	R->wvhdl[n][i-2] = (*iv)[i];
4394	last++;
4395	if (weights[last]==0) weights[last]=(iv)[i]typ;
4396	}
4397	R->block1[n] = last;
4398	break;
4399	case ringorder_ls:
4400	case ringorder_ds:
4401	case ringorder_Ds:
4402	case ringorder_rs:
4403	typ=-1;
4404	case ringorder_lp:
4405	case ringorder_dp:
4406	case ringorder_Dp:
4407	case ringorder_rp:
4408	R->block0[n] = last+1;
4409	if (iv->length() == 3) last+=(*iv)[2];
4410	else last += (*iv)[0];
4411	R->block1[n] = last;
4412	//if ((R->block0[n]>R->block1[n])
4413	//\|\| (R->block1[n]>rVar(R)))
4414	//{
4415	// R->block1[n]=rVar(R);
4416	// //WerrorS("ordering larger than number of variables");
4417	// break;
4418	//}
4419	if (rCheckIV(iv)) return TRUE;
4420	for(i=si_min(rVar(R),R->block1[n]);i>=R->block0[n];i--)
4421	{
4422	if (weights[i]==0) weights[i]=typ;
4423	}
4424	break;
4425	case ringorder_S:
4426	case ringorder_c:
4427	case ringorder_C:
4428	if (rCheckIV(iv)) return TRUE;
4429	break;
4430	case ringorder_aa:
4431	case ringorder_a:
4432	R->block0[n] = last+1;
4433	R->block1[n] = si_min(last+iv->length()-2 , rVar(R));
4434	R->wvhdl[n] = (int)omAlloc((iv->length()-1)sizeof(int));
4435	for (i=2; i<iv->length(); i++)
4436	{
4437	R->wvhdl[n][i-2]=(*iv)[i];
4438	last++;
4439	if (weights[last]==0) weights[last]=(iv)[i]typ;
4440	}
4441	last=R->block0[n]-1;
4442	break;
4443	case ringorder_a64:
4444	{
4445	R->block0[n] = last+1;
4446	R->block1[n] = si_min(last+iv->length()-2 , rVar(R));
4447	R->wvhdl[n] = (int)omAlloc((iv->length()-1)sizeof(int64));
4448	int64 w=(int64 )R->wvhdl[n];
4449	for (i=2; i<iv->length(); i++)
4450	{
4451	w[i-2]=(*iv)[i];
4452	last++;
4453	if (weights[last]==0) weights[last]=(iv)[i]typ;
4454	}
4455	last=R->block0[n]-1;
4456	break;
4457	}
4458	case ringorder_M:
4459	{
4460	int Mtyp=rTypeOfMatrixOrder(iv);
4461	if (Mtyp==0) return TRUE;
4462	if (Mtyp==-1) typ = -1;
4463
4464	R->wvhdl[n] =( int )omAlloc((iv->length()-1)sizeof(int));
4465	for (i=2; i<iv->length();i++)
4466	R->wvhdl[n][i-2]=(*iv)[i];
4467
4468	R->block0[n] = last+1;
4469	last += (int)sqrt((double)(iv->length()-2));
4470	R->block1[n] = last;
4471	for(i=si_min(rVar(R),R->block1[n]);i>=R->block0[n];i--)
4472	{
4473	if (weights[i]==0) weights[i]=typ;
4474	}
4475	break;
4476	}
4477
4478	case ringorder_no:
4479	R->order[n] = ringorder_unspec;
4480	return TRUE;
4481
4482	default:
4483	Werror("Internal Error: Unknown ordering %d", (*iv)[1]);
4484	R->order[n] = ringorder_unspec;
4485	return TRUE;
4486	}
4487	}
4488	sl=sl->next;
4489	}
4490
4491	// check for complete coverage
4492	if ((R->order[n]==ringorder_c) \|\| (R->order[n]==ringorder_C)) n--;
4493	if (R->block1[n] != R->N)
4494	{
4495	if (((R->order[n]==ringorder_dp) \|\|
4496	(R->order[n]==ringorder_ds) \|\|
4497	(R->order[n]==ringorder_Dp) \|\|
4498	(R->order[n]==ringorder_Ds) \|\|
4499	(R->order[n]==ringorder_rp) \|\|
4500	(R->order[n]==ringorder_rs) \|\|
4501	(R->order[n]==ringorder_lp) \|\|
4502	(R->order[n]==ringorder_ls))
4503	&&
4504	R->block0[n] <= R->N)
4505	{
4506	R->block1[n] = R->N;
4507	}
4508	else
4509	{
4510	Werror("mismatch of number of vars (%d) and ordering (%d vars)",
4511	R->N,R->block1[n]);
4512	return TRUE;
4513	}
4514	}
4515	// find OrdSgn:
4516	R->OrdSgn = 1;
4517	for(i=1;i<=R->N;i++)
4518	{ if (weights[i]<0) { R->OrdSgn=-1;break; }}
4519	omFree(weights);
4520	return FALSE;
4521	}
4522
4523	BOOLEAN rSleftvList2StringArray(sleftv* sl, char** p)
4524	{
4525
4526	while(sl!=NULL)
4527	{
4528	if (sl->Name() == sNoName)
4529	{
4530	if (sl->Typ()==POLY_CMD)
4531	{
4532	sleftv s_sl;
4533	iiConvert(POLY_CMD,ANY_TYPE,-1,sl,&s_sl);
4534	if (s_sl.Name() != sNoName)
4535	*p = omStrDup(s_sl.Name());
4536	else
4537	*p = NULL;
4538	sl->next = s_sl.next;
4539	s_sl.next = NULL;
4540	s_sl.CleanUp();
4541	if (*p == NULL) return TRUE;
4542	}
4543	else
4544	return TRUE;
4545	}
4546	else
4547	*p = omStrDup(sl->Name());
4548	p++;
4549	sl=sl->next;
4550	}
4551	return FALSE;
4552	}
4553
4554	////////////////////
4555	//
4556	// rInit itself:
4557	//
4558	// INPUT: s: name, pn: ch & parameter (names), rv: variable (names)
4559	// ord: ordering
4560	// RETURN: currRingHdl on success
4561	// NULL on error
4562	// NOTE: * makes new ring to current ring, on success
4563	// * considers input sleftv's as read-only
4564	//idhdl rInit(char s, sleftv pn, sleftv* rv, sleftv* ord)
4565	ring rInit(sleftv* pn, sleftv* rv, sleftv* ord)
4566	{
4567	int ch;
4568	#ifdef HAVE_RINGS
4569	unsigned int ringtype = 0;
4570	int_number ringflaga = NULL;
4571	unsigned int ringflagb = 1;
4572	#endif
4573	int float_len=0;
4574	int float_len2=0;
4575	ring R = NULL;
4576	idhdl tmp = NULL;
4577	BOOLEAN ffChar=FALSE;
4578	int typ = 1;
4579
4580	/* ch -------------------------------------------------------*/
4581	// get ch of ground field
4582	int numberOfAllocatedBlocks;
4583
4584	if (pn->Typ()==INT_CMD)
4585	{
4586	ch=(int)(long)pn->Data();
4587	}
4588	else if ((pn->name != NULL)
4589	&& ((strcmp(pn->name,"real")==0) \|\| (strcmp(pn->name,"complex")==0)))
4590	{
4591	BOOLEAN complex_flag=(strcmp(pn->name,"complex")==0);
4592	ch=-1;
4593	if ((pn->next!=NULL) && (pn->next->Typ()==INT_CMD))
4594	{
4595	float_len=(int)(long)pn->next->Data();
4596	float_len2=float_len;
4597	pn=pn->next;
4598	if ((pn->next!=NULL) && (pn->next->Typ()==INT_CMD))
4599	{
4600	float_len2=(int)(long)pn->next->Data();
4601	pn=pn->next;
4602	}
4603	}
4604	if ((pn->next==NULL) && complex_flag)
4605	{
4606	pn->next=(leftv)omAlloc0Bin(sleftv_bin);
4607	pn->next->name=omStrDup("i");
4608	}
4609	}
4610	#ifdef HAVE_RINGS
4611	else if ((pn->name != NULL) && (strcmp(pn->name, "integer") == 0))
4612	{
4613	ringflaga = (int_number) omAlloc(sizeof(mpz_t));
4614	mpz_init_set_si(ringflaga, 0);
4615	if ((pn->next!=NULL) && (pn->next->Typ()==INT_CMD))
4616	{
4617	mpz_set_ui(ringflaga, (int)(long) pn->next->Data());
4618	pn=pn->next;
4619	if ((pn->next!=NULL) && (pn->next->Typ()==INT_CMD))
4620	{
4621	ringflagb = (long) pn->next->Data();
4622	pn=pn->next;
4623	}
4624	while ((pn->next!=NULL) && (pn->next->Typ()==INT_CMD))
4625	{
4626	mpz_mul_ui(ringflaga, ringflaga, (int)(long) pn->next->Data());
4627	pn=pn->next;
4628	}
4629	}
4630	if ((mpz_cmp_ui(ringflaga, 1) == 0) && (mpz_cmp_ui(ringflaga, 0) < 0))
4631	{
4632	Werror("Wrong ground ring specification (module is 1)");
4633	goto rInitError;
4634	}
4635	if (ringflagb < 1)
4636	{
4637	Werror("Wrong ground ring specification (exponent smaller than 1");
4638	goto rInitError;
4639	}
4640	// module is 0 ---> integers
4641	if (mpz_cmp_ui(ringflaga, 0) == 0)
4642	{
4643	ch = 0;
4644	ringtype = 4;
4645	}
4646	// we have an exponent
4647	else if (ringflagb > 1)
4648	{
4649	ch = ringflagb;
4650	if ((mpz_cmp_ui(ringflaga, 2) == 0) && (ringflagb + 2 <= 8*sizeof(NATNUMBER)))
4651	{
4652	ringtype = 1; // Use Z/2^ch
4653	}
4654	else
4655	{
4656	ringtype = 3;
4657	}
4658	}
4659	// just a module m > 1
4660	else
4661	{
4662	ringtype = 2;
4663	ch = mpz_get_ui(ringflaga);
4664	}
4665	}
4666	#endif
4667	else
4668	{
4669	Werror("Wrong ground field specification");
4670	goto rInitError;
4671	}
4672	pn=pn->next;
4673
4674	#ifdef HAVE_RINGS
4675	if (ringtype > 0)
4676	{
4677	WarnS("You are using coefficient rings which are not fields.");
4678	WarnS("Please note that only limited functionality is available");
4679	WarnS("for these coefficients.");
4680	WarnS("");
4681	WarnS("The following commands are meant to work:");
4682	WarnS("- basic polynomial arithmetic");
4683	WarnS("- std");
4684	WarnS("- syz");
4685	WarnS("- lift");
4686	WarnS("- reduce");
4687	}
4688	#endif
4689
4690	int l, last;
4691	sleftv * sl;
4692	/every entry in the new ring is initialized to 0/
4693
4694	/* characteristic -----------------------------------------------*/
4695	/* input: 0 ch=0 : Q parameter=NULL ffChar=FALSE float_len
4696	* 0 1 : Q(a,...) *names FALSE
4697	* 0 -1 : R NULL FALSE 0
4698	* 0 -1 : R NULL FALSE prec. >6
4699	* 0 -1 : C *names FALSE prec. 0..?
4700	* p p : Fp NULL FALSE
4701	* p -p : Fp(a) *names FALSE
4702	* q q : GF(q=p^n) *names TRUE
4703	*/
4704	if ((ch!=-1)
4705	#ifdef HAVE_RINGS
4706	&& (ringtype == 0)
4707	#endif
4708	)
4709	{
4710	int l = 0;
4711
4712	if (((ch!=0) && (ch<2))
4713	#ifndef NV_OPS
4714	\|\| (ch > 32003)
4715	#endif
4716	)
4717	{
4718	Warn("%d is invalid characteristic of ground field. 32003 is used.", ch);
4719	ch=32003;
4720	}
4721	// load fftable, if necessary
4722	if (pn!=NULL)
4723	{
4724	while ((ch!=fftable[l]) && (fftable[l])) l++;
4725	if (fftable[l]==0) ch = IsPrime(ch);
4726	else
4727	{
4728	char m[1]={(char )sNoName};
4729	nfSetChar(ch,m);
4730	if (errorreported) goto rInitError;
4731	else ffChar=TRUE;
4732	}
4733	}
4734	else
4735	{
4736	ch = IsPrime(ch);
4737	}
4738	}
4739	// allocated ring and set ch
4740	R = (ring) omAlloc0Bin(sip_sring_bin);
4741	R->ch = ch;
4742	#ifdef HAVE_RINGS
4743	R->ringtype = ringtype;
4744	R->ringflaga = ringflaga;
4745	R->ringflagb = ringflagb;
4746	#endif
4747	if (ch == -1)
4748	{
4749	R->float_len= si_min(float_len,32767);
4750	R->float_len2= si_min(float_len2,32767);
4751	}
4752
4753	/* parameter -------------------------------------------------------*/
4754	if (pn!=NULL)
4755	{
4756	R->P=pn->listLength();
4757	//if ((ffChar\|\| (ch == 1)) && (R->P > 1))
4758	if ((R->P > 1) && (ffChar \|\| (ch == -1)))
4759	{
4760	WerrorS("too many parameters");
4761	goto rInitError;
4762	}
4763	R->parameter=(char*)omAlloc0(R->Psizeof(char_ptr));
4764	if (rSleftvList2StringArray(pn, R->parameter))
4765	{
4766	WerrorS("parameter expected");
4767	goto rInitError;
4768	}
4769	if (ch>1 && !ffChar) R->ch=-ch;
4770	else if (ch==0) R->ch=1;
4771	}
4772	else if (ffChar)
4773	{
4774	WerrorS("need one parameter");
4775	goto rInitError;
4776	}
4777	/* post-processing of field description */
4778	// we have short reals, but no short complex
4779	if ((R->ch == - 1)
4780	&& (R->parameter !=NULL)
4781	&& (R->float_len < SHORT_REAL_LENGTH))
4782	{
4783	R->float_len = SHORT_REAL_LENGTH;
4784	R->float_len2 = SHORT_REAL_LENGTH;
4785	}
4786
4787	/* names and number of variables-------------------------------------*/
4788	{
4789	int l=rv->listLength();
4790	#if SIZEOF_SHORT == 2
4791	#define MAX_SHORT 0x7fff
4792	#endif
4793	if (l>MAX_SHORT)
4794	{
4795	Werror("too many ring variables(%d), max is %d",l,MAX_SHORT);
4796	goto rInitError;
4797	}
4798	R->N = l; /rv->listLength();/
4799	}
4800	R->names = (char *)omAlloc0(R->N sizeof(char_ptr));
4801	if (rSleftvList2StringArray(rv, R->names))
4802	{
4803	WerrorS("name of ring variable expected");
4804	goto rInitError;
4805	}
4806
4807	/* check names and parameters for conflicts ------------------------- */
4808	rNameCheck(R); // conflicting variables will be renamed
4809	/* ordering -------------------------------------------------------------*/
4810	if (rSleftvOrdering2Ordering(ord, R))
4811	goto rInitError;
4812
4813	// Complete the initialization
4814	if (rComplete(R,1))
4815	goto rInitError;
4816
4817	rTest(R);
4818
4819	// try to enter the ring into the name list
4820	// need to clean up sleftv here, before this ring can be set to
4821	// new currRing or currRing can be killed beacuse new ring has
4822	// same name
4823	if (pn != NULL) pn->CleanUp();
4824	if (rv != NULL) rv->CleanUp();
4825	if (ord != NULL) ord->CleanUp();
4826	//if ((tmp = enterid(s, myynest, RING_CMD, &IDROOT))==NULL)
4827	// goto rInitError;
4828
4829	//memcpy(IDRING(tmp),R,sizeof(*R));
4830	// set current ring
4831	//omFreeBin(R, ip_sring_bin);
4832	//return tmp;
4833	return R;
4834
4835	// error case:
4836	rInitError:
4837	if (R != NULL) rDelete(R);
4838	if (pn != NULL) pn->CleanUp();
4839	if (rv != NULL) rv->CleanUp();
4840	if (ord != NULL) ord->CleanUp();
4841	return NULL;
4842	}
4843
4844	ring rSubring(ring org_ring, sleftv* rv)
4845	{
4846	ring R = rCopy0(org_ring);
4847	int perm=(int )omAlloc0((org_ring->N+1)*sizeof(int));
4848	int last = 0, o=0, n = rBlocks(org_ring), i=0, typ = 1, j;
4849
4850	/* names and number of variables-------------------------------------*/
4851	{
4852	int l=rv->listLength();
4853	if (l>MAX_SHORT)
4854	{
4855	Werror("too many ring variables(%d), max is %d",l,MAX_SHORT);
4856	goto rInitError;
4857	}
4858	R->N = l; /rv->listLength();/
4859	}
4860	omFree(R->names);
4861	R->names = (char *)omAlloc0(R->N sizeof(char_ptr));
4862	if (rSleftvList2StringArray(rv, R->names))
4863	{
4864	WerrorS("name of ring variable expected");
4865	goto rInitError;
4866	}
4867
4868	/* check names for subring in org_ring ------------------------- */
4869	{
4870	i=0;
4871
4872	for(j=0;j<R->N;j++)
4873	{
4874	for(;i<org_ring->N;i++)
4875	{
4876	if (strcmp(org_ring->names[i],R->names[j])==0)
4877	{
4878	perm[i+1]=j+1;
4879	break;
4880	}
4881	}
4882	if (i>org_ring->N)
4883	{
4884	Werror("variable %d (%s) not in basering",j+1,R->names[j]);
4885	break;
4886	}
4887	}
4888	}
4889	//Print("perm=");
4890	//for(i=1;i<org_ring->N;i++) Print("v%d -> v%d\n",i,perm[i]);
4891	/* ordering -------------------------------------------------------------*/
4892
4893	for(i=0;i<n;i++)
4894	{
4895	int min_var=-1;
4896	int max_var=-1;
4897	for(j=R->block0[i];j<=R->block1[i];j++)
4898	{
4899	if (perm[j]>0)
4900	{
4901	if (min_var==-1) min_var=perm[j];
4902	max_var=perm[j];
4903	}
4904	}
4905	if (min_var!=-1)
4906	{
4907	//Print("block %d: old %d..%d, now:%d..%d\n",
4908	// i,R->block0[i],R->block1[i],min_var,max_var);
4909	R->block0[i]=min_var;
4910	R->block1[i]=max_var;
4911	if (R->wvhdl[i]!=NULL)
4912	{
4913	omFree(R->wvhdl[i]);
4914	R->wvhdl[i]=(int)omAlloc0((max_var-min_var+1)sizeof(int));
4915	for(j=org_ring->block0[i];j<=org_ring->block1[i];j++)
4916	{
4917	if (perm[j]>0)
4918	{
4919	R->wvhdl[i][perm[j]-R->block0[i]]=
4920	org_ring->wvhdl[i][j-org_ring->block0[i]];
4921	//Print("w%d=%d (orig_w%d)\n",perm[j],R->wvhdl[i][perm[j]-R->block0[i]],j);
4922	}
4923	}
4924	}
4925	}
4926	else
4927	{
4928	if(R->block0[i]>0)
4929	{
4930	//Print("skip block %d\n",i);
4931	R->order[i]=ringorder_unspec;
4932	if (R->wvhdl[i] !=NULL) omFree(R->wvhdl[i]);
4933	R->wvhdl[i]=NULL;
4934	}
4935	//else Print("keep block %d\n",i);
4936	}
4937	}
4938	i=n-1;
4939	while(i>0)
4940	{
4941	// removed unneded blocks
4942	if(R->order[i-1]==ringorder_unspec)
4943	{
4944	for(j=i;j<=n;j++)
4945	{
4946	R->order[j-1]=R->order[j];
4947	R->block0[j-1]=R->block0[j];
4948	R->block1[j-1]=R->block1[j];
4949	if (R->wvhdl[j-1] !=NULL) omFree(R->wvhdl[j-1]);
4950	R->wvhdl[j-1]=R->wvhdl[j];
4951	}
4952	R->order[n]=ringorder_unspec;
4953	n--;
4954	}
4955	i--;
4956	}
4957	n=rBlocks(org_ring)-1;
4958	while (R->order[n]==0) n--;
4959	while (R->order[n]==ringorder_unspec) n--;
4960	if ((R->order[n]==ringorder_c) \|\| (R->order[n]==ringorder_C)) n--;
4961	if (R->block1[n] != R->N)
4962	{
4963	if (((R->order[n]==ringorder_dp) \|\|
4964	(R->order[n]==ringorder_ds) \|\|
4965	(R->order[n]==ringorder_Dp) \|\|
4966	(R->order[n]==ringorder_Ds) \|\|
4967	(R->order[n]==ringorder_rp) \|\|
4968	(R->order[n]==ringorder_rs) \|\|
4969	(R->order[n]==ringorder_lp) \|\|
4970	(R->order[n]==ringorder_ls))
4971	&&
4972	R->block0[n] <= R->N)
4973	{
4974	R->block1[n] = R->N;
4975	}
4976	else
4977	{
4978	Werror("mismatch of number of vars (%d) and ordering (%d vars) in block %d",
4979	R->N,R->block1[n],n);
4980	return NULL;
4981	}
4982	}
4983	omFree(perm);
4984	// find OrdSgn:
4985	R->OrdSgn = org_ring->OrdSgn; // IMPROVE!
4986	//for(i=1;i<=R->N;i++)
4987	//{ if (weights[i]<0) { R->OrdSgn=-1;break; }}
4988	//omFree(weights);
4989	// Complete the initialization
4990	if (rComplete(R,1))
4991	goto rInitError;
4992
4993	rTest(R);
4994
4995	if (rv != NULL) rv->CleanUp();
4996
4997	return R;
4998
4999	// error case:
5000	rInitError:
5001	if (R != NULL) rDelete(R);
5002	if (rv != NULL) rv->CleanUp();
5003	return NULL;
5004	}
5005
5006	void rKill(ring r)
5007	{
5008	if ((r->ref<=0)&&(r->order!=NULL))
5009	{
5010	#ifdef RDEBUG
5011	if (traceit &TRACE_SHOW_RINGS) Print("kill ring %lx\n",(long)r);
5012	#endif
5013	if (r->qideal!=NULL)
5014	{
5015	id_Delete(&r->qideal, r);
5016	r->qideal = NULL;
5017	}
5018	int i=1;
5019	int j;
5020	int *pi=r->order;
5021	#ifdef USE_IILOCALRING
5022	for (j=0;j<iiRETURNEXPR_len;j++)
5023	{
5024	if (iiLocalRing[j]==r)
5025	{
5026	if (j<myynest) Warn("killing the basering for level %d",j);
5027	iiLocalRing[j]=NULL;
5028	}
5029	}
5030	#else /* USE_IILOCALRING */
5031	//#endif /* USE_IILOCALRING */
5032	{
5033	proclevel * nshdl = procstack;
5034	int lev=myynest-1;
5035
5036	for(; nshdl != NULL; nshdl = nshdl->next)
5037	{
5038	if (nshdl->cRing==r)
5039	{
5040	Warn("killing the basering for level %d",lev);
5041	nshdl->cRing=NULL;
5042	nshdl->cRingHdl=NULL;
5043	}
5044	}
5045	}
5046	#endif /* USE_IILOCALRING */
5047	// any variables depending on r ?
5048	while (r->idroot!=NULL)
5049	{
5050	killhdl2(r->idroot,&(r->idroot),r);
5051	}
5052	if (r==currRing)
5053	{
5054	// all dependend stuff is done, clean global vars:
5055	if (r->qideal!=NULL)
5056	{
5057	currQuotient=NULL;
5058	}
5059	if (ppNoether!=NULL) pDelete(&ppNoether);
5060	if (sLastPrinted.RingDependend())
5061	{
5062	sLastPrinted.CleanUp();
5063	}
5064	if ((myynest>0) && (iiRETURNEXPR[myynest].RingDependend()))
5065	{
5066	WerrorS("return value depends on local ring variable (export missing ?)");
5067	iiRETURNEXPR[myynest].CleanUp();
5068	}
5069	currRing=NULL;
5070	currRingHdl=NULL;
5071	}
5072
5073	/* nKillChar(r); will be called from inside of rDelete */
5074	rDelete(r);
5075	return;
5076	}
5077	r->ref--;
5078	}
5079
5080	void rKill(idhdl h)
5081	{
5082	ring r = IDRING(h);
5083	int ref=0;
5084	if (r!=NULL)
5085	{
5086	ref=r->ref;
5087	rKill(r);
5088	}
5089	if (h==currRingHdl)
5090	{
5091	if (ref<=0) { currRing=NULL; currRingHdl=NULL;}
5092	else
5093	{
5094	currRingHdl=rFindHdl(r,currRingHdl,NULL);
5095	}
5096	}
5097	}
5098
5099	idhdl rSimpleFindHdl(ring r, idhdl root, idhdl n=NULL)
5100	{
5101	//idhdl next_best=NULL;
5102	idhdl h=root;
5103	while (h!=NULL)
5104	{
5105	if (((IDTYP(h)==RING_CMD)\|\|(IDTYP(h)==QRING_CMD))
5106	&& (h!=n)
5107	&& (IDRING(h)==r)
5108	)
5109	{
5110	// if (IDLEV(h)==myynest)
5111	// return h;
5112	// if ((IDLEV(h)==0) \|\| (next_best==NULL))
5113	// next_best=h;
5114	// else if (IDLEV(next_best)<IDLEV(h))
5115	// next_best=h;
5116	return h;
5117	}
5118	h=IDNEXT(h);
5119	}
5120	//return next_best;
5121	return NULL;
5122	}
5123
5124	extern BOOLEAN jjPROC(leftv res, leftv u, leftv v);
5125	ideal kGroebner(ideal F, ideal Q)
5126	{
5127	//test\|=Sy_bit(OPT_PROT);
5128	idhdl save_ringhdl=currRingHdl;
5129	ideal resid;
5130	idhdl new_ring=NULL;
5131	if ((currRingHdl==NULL) \|\| (IDRING(currRingHdl)!=currRing))
5132	{
5133	currRingHdl=enterid(omStrDup(" GROEBNERring"),0,RING_CMD,&IDROOT,FALSE);
5134	new_ring=currRingHdl;
5135	IDRING(currRingHdl)=currRing;
5136	}
5137	sleftv v; memset(&v,0,sizeof(v)); v.rtyp=IDEAL_CMD; v.data=(char *) F;
5138	idhdl h=ggetid("groebner");
5139	sleftv u; memset(&u,0,sizeof(u)); u.rtyp=IDHDL; u.data=(char *) h;
5140	u.name=IDID(h);
5141
5142	sleftv res; memset(&res,0,sizeof(res));
5143	if(jjPROC(&res,&u,&v))
5144	{
5145	resid=kStd(F,Q,testHomog,NULL);
5146	}
5147	else
5148	{
5149	//printf("typ:%d\n",res.rtyp);
5150	resid=(ideal)(res.data);
5151	}
5152	// cleanup GROEBNERring, save_ringhdl, u,v,(res )
5153	if (new_ring!=NULL)
5154	{
5155	idhdl h=IDROOT;
5156	if (h==new_ring) IDROOT=h->next;
5157	else
5158	{
5159	while ((h!=NULL) &&(h->next!=new_ring)) h=h->next;
5160	if (h!=NULL) h->next=h->next->next;
5161	}
5162	if (h!=NULL) omFreeSize(h,sizeof(*h));
5163	}
5164	currRingHdl=save_ringhdl;
5165	u.CleanUp();
5166	v.CleanUp();
5167	return resid;
5168	}
5169
5170	void jjINT_S_TO_ID(int n,int *e, leftv res)
5171	{
5172	if (n==0) n=1;
5173	ideal l=idInit(n,1);
5174	int i;
5175	poly p;
5176	for(i=pVariables;i>0;i--)
5177	{
5178	if (e[i]>0)
5179	{
5180	n--;
5181	p=pOne();
5182	pSetExp(p,i,1);
5183	pSetm(p);
5184	l->m[n]=p;
5185	if (n==0) break;
5186	}
5187	}
5188	res->data=(char*)l;
5189	omFreeSize((ADDRESS)e,(pVariables+1)*sizeof(int));
5190	}
5191	BOOLEAN jjVARIABLES_P(leftv res, leftv u)
5192	{
5193	int e=(int )omAlloc0((pVariables+1)*sizeof(int));
5194	int n=pGetVariables((poly)u->Data(),e);
5195	jjINT_S_TO_ID(n,e,res);
5196	return FALSE;
5197	}
5198
5199	BOOLEAN jjVARIABLES_ID(leftv res, leftv u)
5200	{
5201	int e=(int )omAlloc0((pVariables+1)*sizeof(int));
5202	ideal I=(ideal)u->Data();
5203	int i;
5204	int n=0;
5205	for(i=I->nrows*I->ncols-1;i>=0;i--)
5206	{
5207	n=pGetVariables(I->m[i],e);
5208	}
5209	jjINT_S_TO_ID(n,e,res);
5210	return FALSE;
5211	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats: