Context Navigation

source: git/Singular/ipshell.cc @ 45c9bd3

spielwiese

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

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

Download in other formats: