Context Navigation

classify.lib @ dc062fe

spielwiese

Last change on this file since dc062fe was dc062fe, checked in by Hans Schönemann <hannes@…>, 18 years ago
*hannes: Thomas changes git-svn-id: file:///usr/local/Singular/svn/trunk@9370 2c84dea3-7e68-4137-9b89-c4e89433aadc
Property mode set to `100644`
File size: 88.9 KB

Line
1	// KK,GMG last modified: 17.12.00
2	///////////////////////////////////////////////////////////////////////////////
3	version = "$Id: classify.lib,v 1.53 2006-07-31 09:11:10 Singular Exp $";
4	category="Singularities";
5	info="
6	LIBRARY: classify.lib Arnold Classifier of Singularities
7	AUTHOR: Kai Krueger, krueger@mathematik.uni-kl.de
8
9	OVERVIEW:
10	A library for classifying isolated hypersurface singularities w.r.t. right
11	equivalence, based on the determinator of singularities by V.I. Arnold.
12
13	PROCEDURES:
14	basicinvariants(f); computes Milnor number, determinacy-bound and corank of
15	classify(f); normal form of poly f determined with Arnold's method
16	corank(f); computes the corank of f (i.e. of the Hessian of f)
17	Hcode(v); coding of intvec v acoording to the number repetitions
18	init_debug([n]); print trace and debugging information depending on int n
19	internalfunctions();display names of internal procedures of this library
20	milnorcode(f[,e]); Hilbert poly of [e-th] Milnor algebra coded with Hcode
21	morsesplit(f); residual part of f after applying the splitting lemma
22	quickclass(f) normal form of f determined by invariants (milnorcode)
23	singularity(s,[]); normal form of singularity given by its name s and index
24	A_L(s/f); shortcut for quickclass(f) or normalform(s)
25	normalform(s); normal form of singularity given by its name s
26	debug_log(lev,[]); print trace and debugging information w.r.t level>@DeBug
27	swap(a,b); swaps the arguments
28	(parameters in square brackets [] are optional)
29	";
30
31	LIB "inout.lib";
32	LIB "elim.lib";
33	LIB "sing.lib";
34	LIB "makedbm.lib";
35
36	///////////////////////////////////////////////////////////////////////////////
37	proc classify_init
38	{
39	string s;
40	link l="DBM:r NFlist";
41	s = read(l,"VERSION");
42	if (s == "" ) {
43	if (printlevel)
44	{
45	"(classify.lib): Need to create database...";
46	}
47	create_sing_dbm();
48	}
49	close(l);
50	l="DBM:r NFlist";
51	s = read(l,"VERSION");
52	//"(classify.lib): Creation done. Current version:", s;
53	}
54	///////////////////////////////////////////////////////////////////////////////
55
56	proc classify (poly f_in)
57	"USAGE: classify(f); f=poly
58	COMPUTE: normal form and singularity type of f with respect to right
59	equivalence, as given in the book \"Singularities of differentiable
60	maps, Volume I\" by V.I. Arnold, S.M. Gusein-Zade, A.N. Varchenko
61	RETURN: normal form of f, of type poly
62	REMARK: This version of classify is only beta. Please send bugs and
63	comments to: \"Kai Krueger\" <krueger@mathematik.uni-kl.de> @*
64	Be sure to have at least Singular version 1.0.1. Updates can be
65	found at: @*
66	URL=http://www.mathematik.uni-kl.de/~krueger/Singular/
67	NOTE: type init_debug(n); (0 <= n <= 10) in order to get intermediate
68	information, higher values of n give more information.
69	The proc creates several global objects with names all starting
70	with @, hence there should be no name conflicts.
71	EXAMPLE: example classify; shows an example"
72	{
73	//---------------------------- initialisation ---------------------------------
74	poly f_out;
75	int n, i, corank, show_nf;
76	string s2;
77	list v;
78	def ring_ext = basering;
79
80	if(system("with","Namespaces")) {
81	if(defined(A_Z)==0) { proc A_Z = general_lib::A_Z; export A_Z; }
82	}
83	init_debug(); // initialize trace/debug mode
84	if(checkring()) { return(f_in); }
85	n = nvars(basering);
86	show_nf = 1; // return normal form if set to '1'
87
88	// define new ring
89	ring ring_top=char(basering),(x(1..n)),(c,ds);
90
91	map conv_ext2top=ring_ext,maxideal(1);
92
93	if(defined(@ringdisplay)!=0) { kill @ringdisplay; }
94	string @ringdisplay = "ring_ext";
95	export @ringdisplay;
96
97	v = Klassifiziere(conv_ext2top(f_in));
98	poly f_out = v[1];
99	s2 = v[2]; // s2: Typ des Polynoms f z.b: E[18]
100	corank = v[3];
101
102	//---------------- collect results and create return-value --------------------
103	if( s2=="error!" \|\| s2=="NoClass") {
104	setring ring_ext;
105	return(f_in);
106	}
107
108	if(show_nf==1) {
109	poly f_nf = normalform(s2);
110	for(i=corank+1;i<=n;i=i+1) { f_nf = f_nf + x(i)^2; }
111	debug_log(2, "Normal form NF(f)=", f_nf);
112	}
113	if(corank>1) { for(i=corank+1;i<=n;i=i+1) { f_out = f_out + x(i)^2; } }
114	setring ring_ext;
115	map conv_top2ext=ring_top,maxideal(1);
116
117	if(show_nf == 1) { return(conv_top2ext(f_nf)); }
118	else { return(conv_top2ext(f_out)); }
119	}
120	example
121	{"EXAMPLE"; echo=2;
122	ring r=0,(x,y,z),ds;
123	poly f=(x2+3y-2z)^2+xyz-(x-y3+x2*z3)^3;
124	classify(f);
125	init_debug(3);
126	classify(f);
127	}
128
129	///////////////////////////////////////////////////////////////////////////////
130	static proc Klassifiziere (poly f)
131	{
132	//--------------------------- initialisation ----------------------------------
133	string s1;
134	int n, cnt, corank_f, K, Mu;
135	list v, cstn;
136	map PhiG;
137	def ring_top = basering;
138
139	n = nvars(basering); // Zahl der Variablen des aktuellen Rings.
140	PhiG = ring_top, maxideal(1);
141	cstn[4] = PhiG;
142	if( defined(@ringdisplay) == 0)
143	{
144	string @ringdisplay; // Define always 'ringdisplay' to be
145	export @ringdisplay; // able to run 'Show(f)'
146	}
147	@ringdisplay = "setring RingB;";
148	if(defined(RingB)!=0) { kill RingB; }
149	execute ("ring RingB="+charstr(basering)+",("+A_Z("x", n)+"),(c,ds);");
150	export RingB;
151	setring ring_top;
152
153	//---------------------- compute basciinvariants ------------------------------
154	if(jet(f,0) != 0 )
155	{
156	cstn[1] = corank(f); cstn[2]=-1; cstn[3]=1;
157	return(printresult(1, f, "a unit", cstn, -1));
158	}
159
160	debug_log(1, "Computing Basicinvariants of f ...");
161	K, Mu, corank_f = basicinvariants(f);
162	debug_log(0, "About the singularity :");
163	debug_log(0, " Milnor number(f) = "+string(Mu));
164	debug_log(0, " Corank(f) = "+string(corank_f));
165	debug_log(0, " Determinacy <= "+string(K));
166	cstn[1] = corank_f;
167	cstn[2] = Mu;
168	cstn[3] = K;
169
170	if( Mu == 0)
171	{
172	cstn[1]=1; cstn[3]=1;
173	return(printresult(1, f, "A[0]", cstn, 0));
174	}
175
176	if(Mu<0)
177	{
178	debug_log(0, "The Milnor number of the function is infinite.");
179	debug_log(0, "The singularity is not in Arnolds list.");
180	return(printresult(1, 1, "error!", cstn, -1));
181	}
182
183	f = jet(f, K);
184	v = HKclass(milnorcode(f));
185	if(v[2]>0) { debug_log(0, "Guessing type via Milnorcode: ", v[1]);}
186	else
187	{
188	debug_log(0, "Hilbert polynomial not recognised. Milnor code = ",
189	milnorcode(f));
190	}
191	debug_log(0, "");
192	debug_log(0, "Computing normal form ...");
193
194	//------------ step 1, classification according to corank(f) ------------------
195	if(corank_f == 0)
196	{
197	return(printresult(2, f, "A["+string(Mu)+"]", cstn, 0));
198	}
199	if(corank_f == 1)
200	{
201	return(printresult(2, f, "A["+string(Mu)+"]", cstn, 0));
202	}
203	cstn[4] = 0;
204	if(corank_f == 2) { return(Funktion1bis(f, cstn)); }
205	if(corank_f == 3) { return(Funktion1bis(f, cstn)); }
206	return(printresult(105, f, "NoClass", cstn, -1));
207	}
208
209	///////////////////////////////////////////////////////////////////////////////
210	static proc Funktion1bis (poly f, list cstn)
211	{
212	//---------------------------- initialisation ---------------------------------
213	def ring_top=basering;
214	poly g;
215	int n, corank, K;
216	map conv, PhiG;
217	string s1;
218	list v_res, v_class, v, iv;
219
220	corank = cstn[1];
221	K = cstn[3];
222	n = nvars(basering);
223
224	//-------------------- apply morsesplit if n>corank ---------------------------
225	if( n > corank) {
226	debug_log(0,
227	"I have to apply the splitting lemma. This will take some time....:-)");
228	v_res = Morse(f, K, corank, 0);
229	g = v_res[1];
230	PhiG = v_res[2];
231
232	conv = ReOrder(g);
233	g = conv(g);
234	PhiG = conv(PhiG);
235
236	if(defined(RingB) != 0 ) { kill RingB; }
237	ring ring_rest=char(basering),(x(1..corank)),(c,ds);
238
239	map MapReduce=ring_top,maxideal(1);
240	poly G = MapReduce(g);
241	map PhiG=ring_top,maxideal(1);// Konstruiere Id auf r
242	PhiG = MapReduce(PhiG);
243
244	execute("ring RingB="+charstr(basering)+",("+A_Z("x",corank)+"),(c,ds);");
245	export RingB;
246	setring ring_rest;
247	}
248	else {
249	ring ring_rest=char(basering),(x(1..corank)),(c,ds);
250	map PhiG=ring_top,maxideal(1);
251	poly G = PhiG(f);
252	}
253
254	//--------------------- step 1 of Arnold now finished -------------------------
255	map phi = ring_rest,maxideal(1);
256	cstn[4] = phi;
257	if(corank == 2) { v_class = Funktion3(G, cstn); }
258	else {
259	if(corank == 3) { v_class = Funktion50(G, cstn); }
260	else { v_class = printresult(1, f, "error!", cstn, -1); }
261	}
262	//-------------------------- classification done ------------------------------
263	poly f_result = v_class[1];
264	v[2] = v_class[2];
265	v[3] = v_class[3];
266	map Phi = v_class[4];
267	PhiG = Phi(PhiG);
268	setring ring_top;
269	if(defined(conv)!=0) { kill conv; }
270	map conv=ring_rest,maxideal(1);
271	v[1] = conv(f_result);
272	return(v);
273	}
274
275	///////////////////////////////////////////////////////////////////////////////
276	static proc Funktion3 (poly f, list cstn)
277	{
278	//---------------------------- initialisation ---------------------------------
279	poly f3 = jet(f, 3);
280	ideal Jf;
281	int Dim, Mult, Mu;
282	debug_log(1, "Step 3");
283
284	if( f3 == 0 ) { return(Funktion13(f, cstn)); }
285
286	// f3 ~ x3 , x2y+y3 , x2y
287	Jf = std(jacob(f3));
288	Dim = dim(Jf);
289	if(Dim == 0) { return(printresult(4, f, "D[4]", cstn, 0)); }
290
291	Mult = mult(Jf);
292	Mu = cstn[2];
293	if(Dim == 1) {
294	if( Mult == 1) { return(printresult(5,f,"D["+string(Mu)+"]", cstn, 0)); }
295	if( Mult == 2) { return(Funktion6(f, cstn));} // series E,J
296	debug_log(0, "dimension 1 und deg != 1, 2 => error, ",
297	"this should never occur");
298	return(printresult(3, f, "error!", cstn, -1));
299	// Should never reach this line
300	}
301	// Should never reach this line
302	return(printresult(3, f, "error!", cstn, -1));
303	}
304
305	///////////////////////////////////////////////////////////////////////////////
306	static proc Funktion6 (poly f, list cstn)
307	{ // Arnold's steps 6-12
308	//---------------------------- initialisation ---------------------------------
309	poly f3, fk;
310	ideal JetId, Jf;
311	int k, Dim, Mult, n, Mu;
312	map PhiG, Phi;
313	intvec RFlg;
314	list v;
315
316	def ring_top=basering;
317	f3 = jet(f, 3); // 3-Jet von f
318	n = nvars(basering);
319	Mu = cstn[2];
320	PhiG = cstn[4];
321	k = 1;
322	debug_log(1, " Step 6");
323
324	RFlg = GetRf(f, n);
325	v = Faktorisiere(f, f3, 3, 1, RFlg);
326	f = v[1];
327	Phi = v[2];
328	PhiG = Phi(PhiG);
329	cstn[4] = PhiG;
330
331	//---------------------------- begin of loop ----------------------------------
332	while( (6*k) <= Mu ) {
333	JetId = x(1)^3+x(2)^(3*k);
334	fk = jet(f, 3*k, weight(JetId));
335	//--------------------------- step 6(k) ---------------------------------
336	if( fk == Coeff(fk,x(1), x(1)^3)*x(1)^3 ) {
337	JetId = x(1)^3+x(2)^(3*k+1); // check jet x3,y3k+1 : E[6k]
338	fk = jet(f, 3(3k+1), weight(JetId));
339	if( Coeff(fk,x(2),x(2)^(3*k+1)) != 0 ) {
340	return(printresult(7, f, "E["+string(6*k)+"]", cstn, k-1));
341	}
342
343	JetId = x(1)^3+x(1)x(2)^(2k+1); // check jet x3,xy2k+1 : E[6k+1]
344	fk = jet(f, 3(2k+1), weight(JetId));
345	if( Coeff(fk, x(1)x(2), x(1)x(2)^(2*k+1)) != 0 ) {
346	return(printresult(8, f,"E["+string(6*k+1)+"]", cstn, k-1));
347	}
348
349	JetId = x(1)^3+x(2)^(3*k+2); // check jet x3,y3k+1 : E[6k+2]
350	fk = jet(f, 3(3k+2), weight(JetId));
351	if( Coeff(fk,x(2),x(2)^(3*k+2)) != 0 ) {
352	return(printresult(9, f,"E["+string(6*k+2)+"]", cstn, k-1));
353	}
354
355	//------------------------- step 10(k) --------------------------------
356	k++;
357	JetId = x(1)^3+x(2)^(3*k);
358	fk = jet(f, 3*k, weight(JetId));
359	Jf = std(jacob(fk));
360	Dim = dim(Jf);
361
362	if(Dim==0) { return(printresult(11,f,"J["+string(k)+",0]",cstn,k-1)); }
363	Mult = mult(Jf);
364	if( Dim ==1 && Mult==1) {
365	return(printresult(12,f,"J["+string(k)+","+string(Mu - 6*k +2)+"]",
366	cstn, k-1));
367	}
368	if( Dim == 1 && Mult == 2) {
369	if(Coeff(fk, x(2), x(2)^(3*k)) != 0) {
370	v = Faktorisiere(f, fk, 3, k, RFlg);
371	f = v[1];
372	Phi = v[2];
373	PhiG = Phi(PhiG);
374	cstn[4] = PhiG;
375	}
376	}
377	}
378	}
379	// Should never reach this line
380	return(printresult(6, f, "error!", cstn, -1));
381	}
382
383	///////////////////////////////////////////////////////////////////////////////
384	static proc Funktion13 (poly f, list cstn)
385	{
386	//---------------------------- initialisation ---------------------------------
387	poly f4;
388	ideal Jf;
389	int Dim, Mult, Mu;
390
391	debug_log(1, " Step 13");
392	Mu = cstn[2];
393	f4 = jet(f, 4);
394	if( f4 == 0 ) { return(Funktion47(f, cstn)); }
395
396	// f4 ~ x4+ax2y2+y4, x4+x2y2, x2y2, x3y, x4
397	Jf = std(jacob(f4));
398	Dim = dim(Jf);
399
400	if(Dim==0) { return(printresult(14,f,"X[9] = X[1,0] = T[2,4,4]",cstn,1)); }
401	Mult = mult(Jf);
402	if( Dim == 1) {
403	if( Mult == 1 ) {
404	return(printresult(15, f,
405	"X[1,"+string(Mu-9)+"] = T[2,4,"+string(Mu-5)+"]", cstn, 1));
406	}
407	if( Mult == 2 ) {
408	Jf = Jf, jacob(Jf);
409	Jf = std(Jf);
410	Dim = dim(Jf);
411	if(Dim==0){return(printresult(16,f,"Y[1,p,q] = T[2,4+p,4+q]",cstn,1));}
412	if( Dim == 1 ) { return(Funktion17(f, cstn)); }
413	}
414	if( Mult == 3 ) { return(Funktion25(f, cstn)); }
415	}
416	// Should never reach this line
417	return(printresult(13, f, "error!", cstn, -1));
418	}
419
420	///////////////////////////////////////////////////////////////////////////////
421	static proc Funktion17 (poly f, list cstn)
422	{ // Analog zu Fumktion 6, Kombination 17-24
423	//---------------------------- initialisation ---------------------------------
424	poly fk, ft;
425	ideal JetId, Jf;
426	int p, Dim, Mult, Mu;
427	list v;
428	map PhiG, Phi;
429
430	def ring_top=basering;
431	debug_log(1, " Step 17");
432	Mu = cstn[2];
433	PhiG = cstn[4];
434	fk = jet(f, 4);
435	p = 1;
436
437	//---------------------------- begin of loop ----------------------------------
438	while( 3*p<= Mu) {
439	debug_log(1, " Step 18("+string(p)+")");
440	v = Isomorphie_s17(f, fk, p, 1);
441	f, Phi = v[1..2];
442	PhiG = Phi(PhiG);
443	cstn[4] = PhiG;
444
445	if ( p>1) {
446	JetId = x(1)^3x(2) + x(2)^(3p);
447	fk = jet(f, 3*p, weight(JetId));
448	}
449	//--------------------------- step 18(p) --------------------------------
450	JetId = x(1)^3x(2) + x(2)^(3p+2); // check jet x3y,y3k+2 : Z[6p+5]
451	fk = jet(f, 3(3p+2), weight(JetId));
452	if( Coeff(fk, x(2), x(2)^(3*p+2)) != 0) {
453	return(printresult(19,f, "Z["+string(6*p+5)+"]", cstn, p));
454	}
455
456	JetId = x(1)^3x(2)+x(1)x(2)^(2*p+2); // check jet x3y,xy2k+2 : Z[6p+6]
457	fk = jet(f, 2(3p+2)+1, weight(JetId));
458	if( Coeff(fk, x(1)x(2), x(1)x(2)^(2*p+2)) != 0) {
459	return(printresult(20, f, "Z["+string(6*p+6)+"]", cstn, p));
460	}
461
462	JetId = x(1)^3x(2) + x(2)^(3p+3); // check jet x3y,y3k+3 : Z[6p+7]
463	fk = jet(f, 3(3p+3), weight(JetId));
464	if( Coeff(fk, x(2), x(2)^(3*p+3)) != 0) {
465	return(printresult(21, f, "Z["+string(6*p+7)+"]", cstn, p));
466	}
467
468	//---------------------------- step 22 ----------------------------------
469	p = p+1;
470	JetId = x(1)^3x(2) + x(2)^(3p+1);
471	fk = jet(f, 3*p+1, weight(JetId));
472	ft = Teile(fk, x(2));
473	Jf = std(jacob(ft));
474	Dim = dim(Jf);
475	Mult = mult(Jf);
476	if(Dim==0) { return(printresult(23,f,"Z["+string(p-1)+",0]", cstn, p)); }
477	if( Mult == 1 ) {
478	return(printresult(24, f, "Z["+string(p-1)+","+string(Mu-3-6*p)+"]",
479	cstn, p));
480	}
481	}
482	// Should never reach this line
483	return(printresult(17, f, "error!", cstn, -1));
484	}
485
486	///////////////////////////////////////////////////////////////////////////////
487	static proc Funktion25 (poly f, list cstn)
488	{ // Analog zu Fumktion 6, Kombination 25-46
489	//---------------------------- initialisation ---------------------------------
490	poly fk, ft;
491	ideal JetId, Jf;
492	int k, Dim, Mult, Mu;
493	map PhiG, Phi;
494	intvec RFlg;
495	list v;
496	def ring_top=basering;
497
498	debug_log(1, " Step 25");
499	Mu = cstn[2];
500	PhiG = cstn[4];
501	fk = jet(f, 4);
502	k = 1;
503	RFlg = GetRf(f, 2);
504
505	//---------------------------- begin of loop ----------------------------------
506	while (k<Mu) {
507	v = Faktorisiere(f, fk, 4 , k, RFlg);
508	f, Phi = v[1..2];
509	PhiG = Phi(PhiG);
510	cstn[4] = PhiG;
511
512	//--------------------------- step 26(k) --------------------------------
513	JetId = x(1)^4 + x(2)^(4*k+1); // check jet x4,y4k+1 : W[12k]
514	fk = jet(f, 4(4k+1), weight(JetId));
515	if( Coeff(fk, x(2), x(2)^(4*k+1)) != 0) {
516	return(printresult(27, f, "W["+string(12k)+"]", cstn, 3k-2));
517	}
518
519	JetId = x(1)^4 + x(1)x(2)^(3k+1); // check jet x4,xy3k+1 : W[12k+1]
520	fk = jet(f, 4(3k+1), weight(JetId));
521	if( Coeff(fk, x(1)x(2), x(1)x(2)^(3*k+1)) != 0) {
522	return(printresult(28, f, "W["+string(12k+1)+"]", cstn, 3k-2));
523	}
524
525	//--------------------------- step 29(k) --------------------------------
526	JetId = x(1)^4 + x(2)^(4*k+2);
527	fk = jet(f, 2(4k+2), weight(JetId));
528	if( Coeff(fk, x(2), x(2)^(4*k+2)) != 0) {
529	Jf = std(jacob(fk));
530	Dim = dim(Jf);
531	if(Dim==0) {return(printresult(30,f,"W["+string(k)+",0]",cstn,3*k-1));}
532	if(Dim==1) {
533	return(printresult(32, f,
534	"W#["+string(k)+","+string(Mu-12k-2)+"]", cstn, 3k-1));
535	}
536	return(printresult(29, f, "error!", cstn, -1));
537	}
538	else {
539	// x^4 oder x^2(x^2+x(2)^2k+1)
540	ft = Teile(fk, x(1)^2);
541	Jf = std(jacob(ft));
542	Dim = dim(Jf);
543	if( Dim == 0 ) {
544	return(printresult(31, f, "W["+string(k)+","+string(Mu-12*k-3)+"]",
545	cstn, 3*k-1));
546	}
547	if( Dim != 1 ) { return(printresult(29, f, "error!", cstn, -1)); }
548
549	//-------------------------- step 33(k) -------------------------------
550	JetId = x(1)^4 + x(1)x(2)^(3k+2); // check jet x4,xy3k+2 : W[12k+5]
551	fk = jet(f, 4(3k+2), weight(JetId));
552	if( Coeff(fk, x(1)x(2), x(1)x(2)^(3*k+2)) != 0) {
553	return(printresult(34, f,"W["+string(12k+5)+"]", cstn, 3k-1));
554	}
555
556	JetId = x(1)^4 + x(2)^(4*k+3); // check jet x4,y4k+3 : W[12k+6]
557	fk = jet(f, 4(4k+3), weight(JetId));
558	if( Coeff(fk, x(2), x(2)^(4*k+3)) != 0){
559	return(printresult(35, f,"W["+string(12k+6)+"]", cstn, 3k-1));
560	}
561
562	//-------------------------- step 36(k) -------------------------------
563	k = k+1;
564	JetId = x(1)^4 + x(2)^(4*k);
565	fk = jet(f, (4*k), weight(JetId));
566	Jf = std(jacob(fk));
567	Dim = dim(Jf);
568	Mult = mult(Jf);
569	if(Dim==0) {return(printresult(37,f,"X["+string(k)+",0]",cstn,3*k-1));}
570	if(Dim==1) {
571	if(Mult==1) {
572	return(printresult(38, f,"X["+string(k)+","+string(Mu-12*k+3)+"]",
573	cstn, 3*k-1));
574	}
575	if(Mult==2) {
576	ft = Teile(fk, x(1)^2);
577	Jf = std(jacob(ft));
578	Dim = dim(Jf);
579	if( Dim == 0) { return(Funktion40(f, cstn, k)); }
580	if( Dim == 1) {
581	return(printresult(39, f, "Y["+string(k)+",r,s]", cstn,3*k-1));
582	}
583	}
584	if(Mult!=3) {
585	return(printresult(36, f, "error!", cstn, -1)); }
586	}
587	else { return(printresult(36, f, "error!", cstn, -1)); }
588	}
589	} // Ende der While-Schleife
590	// Should never reach this line
591	return(printresult(25, f, "error!", cstn, -1));
592	}
593
594	///////////////////////////////////////////////////////////////////////////////
595	static proc Funktion40 (poly f, list cstn, int k)
596	{
597	//---------------------------- initialisation ---------------------------------
598	int r, kr, rr, sr, oldDebug;
599	poly fk, f2, a, b, c;
600	ideal JetId, Jfsyz;
601	string Typ, RestRing, s1;
602	list v1, v2;
603	def ring_top=basering;
604
605	debug_log(1, " Step 40" );
606
607	//------------------------------ compute f2 -----------------------------------
608	JetId = x(1)^4 + x(2)^(4*k);
609	fk = jet(f, (4*k), weight(JetId));
610	f2 = -fk / x(1)^3;
611	Jfsyz = f - fk, x(1)^3, f2;
612	matrix Mat = matrix(syz(Jfsyz));
613	a = Mat[2,1] / Mat[1,1] - Mat[2,2];
614	b = - Mat[3,1] / Mat[1,1] + Mat[3,2];
615	ring tmp_ring=char(basering), (x(1),x(2)),(c,ds);
616	map map_top2tmp=ring_top,maxideal(1);
617	oldDebug = @DeBug;
618	init_debug(-1);
619	//------------------------------ classify f2 ----------------------------------
620	v1=Klassifiziere(map_top2tmp(b));
621	init_debug(oldDebug);
622	Typ = v1[2];
623	v2 = DecodeNormalFormString(Typ);
624	Typ, kr, rr, sr = v2[1..4];
625	r = kr-k;
626	setring ring_top;
627	if( Typ == "E[6k]" ) {
628	return(printresult(42, f, "Z["+string(k)+","+string(12k+6r-1)+"]",
629	cstn, 3*k+r-2));
630	}
631	if( Typ == "E[6k+1]" ) {
632	return(printresult(43, f, "Z["+string(k)+","+string(12k+6r)+"]",
633	cstn, 3*k+r-2));
634	}
635	if( Typ == "E[6k+2]" ) {
636	return(printresult(44, f, "Z["+string(k)+","+string(12k+6r+1)+"]",
637	cstn, 3*k+r-2));
638	}
639	if( Typ == "J[k,0]" ) {
640	return(printresult(45, f, "Z["+string(k)+","+string(r)+",0]",
641	cstn, 3*k+r-2));
642	}
643	if( Typ == "J[k,r]" ) {
644	return(printresult(46,f,"Z["+string(k)+","+string(r)+","+string(rr)+"]",
645	cstn, 3*k+r-2));
646	}
647	// Should never reach this line
648	return(printresult(40, f, "error!", cstn, -1));
649	}
650
651	///////////////////////////////////////////////////////////////////////////////
652	static proc Funktion50 (poly f, list cstn)
653	{
654	//---------------------------- initialisation ---------------------------------
655	poly f3;
656	ideal Jf, Jf1, Jf2;
657	int Dim, Mult, Mu;
658	debug_log(1, "Step 50");
659
660	f3 = jet(f, 3);
661	if( f3 == 0 ) { return(printresult(104, f, "NoClass", cstn, -1)); }
662
663	// f3 ~
664	Jf1 = jacob(f3);
665	Jf = std(Jf1);
666	Dim = dim(Jf);
667	Mult = mult(Jf);
668	Mu = cstn[2];
669
670	if(Dim == 0) {
671	return(printresult(51, f, "P[8] = T[3,3,3]", cstn, 1));
672	} // x3 + y3 + z3 + axyz
673	if(Dim == 1) {
674	if (Mult == 1) {
675	return(printresult(52, f,"P["+string(Mu)+"] = T[3,3,"+string(Mu-5)+"]",
676	cstn, 1));
677	} // x3 + y3 + xyz
678	if(Mult == 2) {
679	Jf2 = wedge(jacob(Jf1),3-Dim), Jf1;
680	Jf2 = std(Jf2);
681	Dim = dim(Jf2);
682	if (Dim==0) { return(printresult(54,f,"R[p,q] = T[3,p,q]", cstn, 1)); }
683	if (Dim==1) { return(Funktion58(f, cstn)); } // x3 + yz2
684	}
685	if(Mult == 3) {
686	Jf2 = wedge(jacob(Jf1),3-Dim), Jf1;
687	Jf2 = std(Jf2);
688	Dim = dim(Jf2);
689	if(Dim == 0) { return(printresult(56, f, "T[p,q,r]", cstn, 1)); }
690	if(Dim == 1) { return(Funktion66(f, cstn)); } // x2z + yz2
691	}
692	if(Mult == 4) { return(Funktion82(f, cstn)); } // x3 + xz2
693	}
694	if(Dim == 2) {
695	if(Mult == 1) { return(Funktion97(f, cstn)); } // x2y
696	if(Mult == 2) { return(printresult(103,f,"NoClass", cstn, -1));}
697	}
698
699	// Should never reach this line
700	return(printresult(50, f, "error!", cstn, -1));
701	}
702
703	///////////////////////////////////////////////////////////////////////////////
704	static proc Funktion58 (poly fin, list cstn)
705	{
706	//---------------------------- initialisation ---------------------------------
707	poly f, f3, a, b, phi, b1, b2, b3, fa, fb, fc;
708	ideal B, Jf3, J1, J2, S, Jfsyz;
709	int kx, ky, kz;
710	string tp;
711	matrix M[2][3];
712	matrix C[2][3], D;
713	list v;
714	map PhiG, VERT;
715	def ring_top=basering;
716	debug_log(1, " Step 58");
717
718	f = fin;
719	f3 = jet(f, 3);
720	PhiG = cstn[4];
721	tp = "Nix";
722	kx = 1; // Koordinate x
723	ky = 2; // Koordinate y
724	kz = 3; // Koordinate z
725	B = maxideal(1); // ideal fuer Abbildungen
726	Jf3 = jacob(f3);
727	S = sat(Jf3, maxideal(1))[1];
728	J1 = diff(S[1], x(kx)), diff(S[1], x(ky)), diff(S[1], x(kz)),
729	diff(S[2], x(kx)), diff(S[2], x(ky)), diff(S[2], x(kz));
730	M = J1;
731	J2 = minor(M, 2), S;
732
733	//------------------ determine coordinate named 'x' -----------------------
734	S = sat(J2, maxideal(1))[1];
735	J1 = Coeff(S[1], x(1), x(1)), Coeff(S[1], x(2), x(2)),
736	Coeff(S[1], x(3), x(3)), Coeff(S[2], x(1), x(1)),
737	Coeff(S[2], x(2), x(2)), Coeff(S[2], x(3), x(3));
738	C = J1;
739	D = syz(C);
740	kill C;
741
742	b1 = D[1,1];
743	b2 = D[2,1];
744	b3 = D[3,1];
745
746	debug_log(6, "f3,s1=", Show(f3));
747	if( b1 != 0) {
748	VERT=ring_top,-1b1x(1), -1b2x(1)+x(2), -1b3x(1) + x(3);
749	kx=1; ky=2; kz=3;
750	}
751	else {
752	if( b2 != 0) {
753	VERT=ring_top, x(1) + -1b1x(2), -1b2x(2), -1b3x(2) + x(3);
754	kx=2; ky=1; kz=3;
755	}
756	else {
757	if( b3 != 0) {
758	VERT=ring_top,x(1) + -1b1x(3), x(2) + -1b2x(3), -1b3x(3);
759	kx=3; ky=1; kz=2;
760	}
761	}
762	}
763	f = VERT(f);
764	PhiG = VERT(PhiG);
765	cstn[4] = PhiG;
766	debug_log(6, VERT);
767	f3 = jet(f,3);
768	debug_log(6, "f3,s2=", Show(f3));
769
770	//---------------- compute f_2 such that j3f = xf_2+f_3 -------------------
771	debug_log(6, "1) x=", kx, " y=", ky, " z=", kz );
772	matrix C = Coeffs(f3, x(kx));
773	fb=C[2,1]; // Coeff von x^1
774	fc=C[1,1]; // Coeff von x^0
775	if(diff(fb, x(ky)) != 0) {
776	Jfsyz = fb, diff(fb, x(ky));
777	matrix Mat = matrix(syz(Jfsyz));
778	B = maxideal(1); // setze Abbildungs-ideal
779	if( nrows(Mat) == 2) {
780	poly Relation = -2 * Mat[2,1] / Mat[1,1];
781	b = Coeff(Relation, x(kz), x(kz));
782	B[rvar(x(ky))] = x(ky)-b*x(kz);
783	}
784	else {
785	Jfsyz = fb, diff(fb, x(kz));
786	Mat = matrix(syz(Jfsyz));
787	poly Relation = -2 * Mat[2,1];
788	a = Coeff(Relation, x(ky), x(ky));
789	B[rvar(x(kz))] = x(kz)-a*x(kz);
790	ky, kz = swap(ky, kz);
791	}
792	VERT = ring_top, B;
793	f = VERT(f);
794	PhiG = VERT(PhiG);
795	cstn[4] = PhiG;
796	f3 = jet(f,3);
797	kill Mat;
798	}
799	else { ky,kz = swap(ky,kz); }
800
801	//------- compute tschirnhaus for 'z' and get f3=f_1(x,y,z)y^2+z^3 --------
802	C = Coeffs(f3, x(kx));
803	fb = C[2,1]; // Coeff von x^1
804	fc = C[1,1]; // Coeff von x^0
805	v = tschirnhaus(fc, x(kz));
806	fc, VERT = v[1..2];
807	f = VERT(f);
808	PhiG = VERT(PhiG);
809	cstn[4] = PhiG;
810	f3 = jet(f,3);
811
812	//------------------- compute f_1 and get f3=xy^2+z^3 ---------------------
813	fb = (f3 - 1(Coeffs(f3, x(kz))[4,1])x(kz)^3)/(x(ky)^2);
814	fc=(x(kx)-1(Coeffs(fb,x(ky))[2,1])x(ky)-1(Coeffs(fb,x(kz))[2,1])x(kz));
815	fa = Coeffs(fb, x(kx))[2,1];
816	if ( fa != 0 ) {
817	B = maxideal(1);
818	B[rvar(x(kx))] = fc / fa;
819	VERT = ring_top, B;
820	f = VERT(f);
821	PhiG = VERT(PhiG);
822	cstn[4] = PhiG;
823	f3 = jet(f,3);
824	}
825
826	//--------------------- permutation of x,y,z -----------------------------
827	if(Coeffs(f3, x(1))[4,1]!=0) {
828	kx=1;
829	if(Coeffs(f3, x(2))[3,1]==0) { ky=2; kz=3; }
830	else { ky=3; kz=2; }
831	}
832	else {
833	if(Coeffs(f3, x(2))[4,1]!=0) {
834	kx=2;
835	if(Coeffs(f3, x(3))[3,1]==0) { ky=3; kz=1; }
836	else { ky=1; kz=3; }
837	}
838	else {
839	kx=3;
840	if(Coeffs(f3, x(1))[3,1]==0) { ky=1; kz=2; }
841	else { ky=2; kz=1; }
842	}
843	}
844	VERT = ring_top, x(kx), x(ky), x(kz);
845	f = VERT(f);
846	PhiG = VERT(PhiG);
847	cstn[4] = PhiG;
848	f3 = jet(f,3);
849	return(Funktion59(f, cstn));
850	}
851
852	///////////////////////////////////////////////////////////////////////////////
853	static proc Funktion59 (poly f, list cstn)
854	{
855	//---------------------------- initialisation ---------------------------------
856	poly phi, fr, fk, alpha, beta, f_tmp;
857	ideal JetId;
858	int p, Dim, Mult, Mu;
859	string tp;
860	list v;
861	map PhiG, Phi;
862	intvec w, RFlg;
863	def ring_top=basering;
864	debug_log(1, " Step 59");
865
866	Mu = cstn[2];
867	PhiG = cstn[4];
868	tp = "Nix";
869	p = 1;
870	phi = jet(f,3);
871	fr = f - phi;
872	RFlg = 1,2,3;
873	alpha = coeffs(fr, x(1))[1,1];
874	beta = (fr - alpha) / x(1);
875	debug_log(3, "f = ", Show(f));
876	debug_log(3, "fr = ", Show(fr));
877	debug_log(3, "alpha= ", Show(alpha));
878	debug_log(3, "beta = ", Show(beta));
879
880	//---------------------------- begin of loop ----------------------------------
881	while(6*p<Mu) {
882	JetId = x(2)^(3*p+1);
883	JetId = phi + x(2)^(3*p+1);
884	//--------------------------- step 59(k) --------------------------------
885	w = weight(JetId);
886	fk = jet(fr, 3*w[1], w);
887	if(fk!=0) { return(printresult(60,f, "Q["+string(6*p+4)+"]", cstn, p)); }
888
889	JetId = phi + x(1)x(2)^(2p+1);
890	w = weight(JetId);
891	fk = jet(fr, 3*w[1], w);
892	if(fk!=0) { return(printresult(61,f, "Q["+string(6*p+5)+"]", cstn, p)); }
893
894	JetId = phi + x(2)^(3*p+2);
895	w = weight(JetId);
896	fk = jet(fr, 3*w[1], w);
897	if(fk!=0) { return(printresult(62,f, "Q["+string(6*p+6)+"]", cstn, p)); }
898
899	//--------------------------- step 63(k) --------------------------------
900	p++;
901	JetId = phi + x(2)^(3*p);
902	w = weight(JetId);
903	fk = jet(f, 3*w[1], w);
904	JetId = std(jacob(fk));
905	Dim = dim(JetId);
906	Mult = mult(JetId);
907	if(Dim==0) { return(printresult(64, f, "Q["+string(p)+",0]", cstn, p)); }
908	if(Dim==1) {
909	if(Mult == 1) {
910	return(printresult(65, f, "Q["+string(p)+","+string(Mu-(6*p+2))+"]",
911	cstn, p));
912	}
913	if(Mult == 2) {
914	fk = jet(fr, 3*w[1], w);
915	f_tmp = Coeffs(phi, x(1))[4,1] *x(1)^3+fk;
916	v = Faktorisiere(f, f_tmp, 3 , p, RFlg);
917	f = v[1];
918	Phi = v[2];
919	PhiG = Phi(PhiG);
920	cstn[4] = PhiG;
921	fr = f - phi;
922	}
923	}
924	}
925	// Should never reach this line
926	return(printresult(59, f, "error!", cstn, -1));
927	}
928
929	///////////////////////////////////////////////////////////////////////////////
930	static proc Funktion66 (poly f, list cstn)
931	{
932	//---------------------------- initialisation ---------------------------------
933	int kx = 1; // Koordinate x
934	int ky = 2; // Koordinate y
935	int kz = 3; // Koordinate z
936	poly f3 = jet(f, 3);
937	ideal JetId;
938
939	debug_log(1, " Step 66");
940	debug_log(2, "F3=", Show(f3));
941	poly fx = diff(f3, x(kx));
942	JetId = jacob(fx);
943	JetId = std(JetId);
944	"nach x:",Show(fx), " Id=", JetId, " Dim=", dim(JetId);
945
946	poly fy = diff(f3, x(ky));
947	JetId = jacob(fx);
948	JetId = std(JetId);
949	"nach y:",Show(fy), " Id=", JetId, " Dim=", dim(JetId);
950
951	poly fz = diff(f3, x(kz));
952	JetId = jacob(fx);
953	JetId = std(JetId);
954	"nach z:",Show(fz), " Id=", JetId, " Dim=", dim(JetId);
955	return(printresult(1, 66, "error!", cstn, -1));
956	}
957
958	///////////////////////////////////////////////////////////////////////////////
959	static proc Funktion82 (poly f, list cstn)
960	{
961	//---------------------------- initialisation ---------------------------------
962	poly f3, b1, b2, b3;
963	int i, kx, ky, kz, Fall;
964	ideal Jfsyz, B;
965	intvec kv = 1,2,3;
966	matrix Mat;
967	map PhiG, VERT;
968	list v;
969	def ring_top=basering;
970	debug_log(1, " Step 82");
971
972	f3 = jet(f,3);
973	kx = 1; // Koordinate x
974	ky = 2; // Koordinate y
975	kz = 3; // Koordinate z
976	B = maxideal(1);
977	Jfsyz = jacob(f3);
978	PhiG = cstn[4];
979	Fall = 2;
980
981	//------------------ find coordinatechange that f3 ~ g(x,z) -------------------
982	if (diff(f3, x(1)) == 0) { kx, ky = swap(kx, ky); }
983	if (diff(f3, x(2)) == 0) { }
984	if (diff(f3, x(3)) == 0) { kz, ky = swap(kz, ky); }
985	if ( (diff(f3, x(1)) != 0) && (diff(f3, x(2)) != 0) &&
986	(diff(f3, x(3)) != 0) ) {
987	Mat = matrix(syz(Jfsyz));
988	b1 = Mat[1,1];
989	b2 = Mat[2,1];
990	b3 = Mat[3,1];
991
992	if( b1 != 0) {
993	VERT = ring_top,b1x(kx), b2x(kx)+x(ky), b3*x(kx) + x(kz);
994	kx, ky = swap(kx, ky);
995	}
996	else {
997	if(b2!=0) { VERT = ring_top,x(kx)+b1x(ky),b2x(ky),b3*x(ky)+x(kz); }
998	else {
999	if(b3!=0) { VERT = ring_top,x(kx)+b1x(kz),x(ky)+b2x(kz),b3*x(kz); }
1000	else { VERT = ring_top,B; }
1001	}
1002	}
1003	f = VERT(f);
1004	PhiG = VERT(PhiG);
1005	cstn[4] = PhiG;
1006	}
1007	VERT = ring_top,x(kx),x(ky),x(kz);
1008	f = VERT(f);
1009	PhiG = VERT(PhiG);
1010	cstn[4] = PhiG;
1011	f3 = jet(f,3);
1012
1013	if( (f3-subst(f3, x(kx), 0)) == 0) { kx, ky = swap(kx, ky); }
1014	if( (f3-subst(f3, x(kz), 0)) == 0) { kz, ky = swap(kz, ky); }
1015
1016	//------------ find coordinatechange for f3 ~ x3+xz2, if possible ------------
1017	matrix C = coeffs(f3, x(kx));
1018	if(size(C) == 3) { C = coeffs(f3, x(kz)); kx,kz=swap(kx, kz); }
1019	if(C[1,1] == 0 && C[3,1] == 0) { Fall = 1; }
1020	if(C[1,1] != 0 && C[3,1] != 0 ) { Fall = 3; }
1021	if(C[1,1] == 0 && C[3,1] != 0 ) { Fall = 2; }
1022	if(C[1,1] != 0 && C[3,1] == 0 ) { Fall = 2; kx,kz=swap(kx, kz); }
1023
1024	if(Fall == 1) {
1025	VERT=ring_top,x(kx),x(ky),x(kz);
1026	}
1027	if(Fall == 2) {
1028	v = tschirnhaus(f3/x(kz), x(kx));
1029	b1, VERT = v[1..2];
1030	}
1031	if(Fall == 3) {
1032	v = tschirnhaus(f3/x(kx), x(kx));
1033	b1, VERT = [1..2];
1034	debug_log(2, "B1=", Show(jet(VERT(f),3)));
1035	v = tschirnhaus(f3/x(kz), x(kz));
1036	b2, VERT = [1..2];
1037	debug_log(2, "B2=", Show(jet(VERT(f),3)));
1038	}
1039	f = VERT(f);
1040	PhiG = VERT(PhiG);
1041	cstn[4] = PhiG;
1042	f3 = jet(f,3);
1043
1044	//------------- if f3 ~ x3+xz2 then continue with classification -------------
1045	C = coeffs(f3, x(1));
1046	if( C[1,1] == 0 && C[2,1] != 0 && C[3,1] == 0 && C[4,1] != 0 ) {
1047	return(Funktion83(f, cstn));
1048	}
1049	return(printresult(82, f, "error!", cstn, -1));
1050	}
1051
1052	///////////////////////////////////////////////////////////////////////////////
1053	static proc Isomorphie_s82_z (poly f, poly fk, int p)
1054	{
1055	//---------------------------- initialisation ---------------------------------
1056	poly Relation, a, b;
1057	ideal Jfsyz, B;
1058	matrix Mat;
1059	map VERT;
1060	list v;
1061	def ring_top=basering;
1062
1063	debug_log(1, " Isomorphie 82/90 z");
1064	debug_log(2, "tt=", Show(fk));
1065	Jfsyz = fk, diff(fk, x(3));
1066	Mat = matrix(syz(Jfsyz));
1067	Relation = -2 * Mat[2,1] / Mat[1,1];
1068	a = Coeff(Relation, x(3), x(3));
1069	b = Coeff(Relation, x(2), x(2)^p);
1070	B = maxideal(1);
1071	B[rvar(x(3))] = x(3)-b*x(2)^p;
1072	VERT = ring_top,B;
1073	v = VERT(f), VERT;
1074	debug_log(2, VERT);
1075	debug_log(2, " z res=", Show(VERT(fk)));
1076	return(v);
1077	}
1078
1079	///////////////////////////////////////////////////////////////////////////////
1080	static proc Isomorphie_s82_x (poly f, poly fk, int p)
1081	{
1082	//---------------------------- initialisation ---------------------------------
1083	matrix Mat;
1084	poly Relation, a, b;
1085	ideal Jfsyz, B;
1086	map VERT;
1087	list v;
1088	def ring_top=basering;
1089
1090	debug_log(1, " Isomorphie 82/90 x");
1091	debug_log(2, "tt=", Show(fk));
1092	Jfsyz = fk, diff(fk, x(1));
1093	Mat = matrix(syz(Jfsyz));
1094	Relation = -3 * Mat[2,1] / Mat[1,1];
1095	a = Coeff(Relation, x(1), x(1));
1096	b = Coeff(Relation, x(2), x(2)^p);
1097	B = maxideal(1);
1098	B[rvar(x(1))] = x(1)-b*x(2)^p;
1099	VERT = ring_top,B;
1100	v = VERT(f), VERT;
1101	debug_log(2, VERT);
1102	debug_log(2, " x res=", Show(VERT(fk)));
1103
1104	return(v);
1105	}
1106
1107	///////////////////////////////////////////////////////////////////////////////
1108	static proc Funktion83 (poly f, list cstn)
1109	{
1110	//---------------------------- initialisation ---------------------------------
1111	poly fk, phi, a, b;
1112	ideal JetId, Jf, B;
1113	int k, Dim, Mult;
1114	intvec w;
1115	map PhiG, Phi;
1116	list v;
1117	matrix Mat;
1118	def ring_top=basering;
1119	debug_log(1, " Step 83");
1120
1121	k = 1;
1122	PhiG = cstn[4];
1123	//---------------------------- begin of loop ----------------------------------
1124	while(k<10) {
1125	phi = jet(f, 3);
1126	//--------------------------- step 83(k) --------------------------------
1127	JetId = x(1)^3 + x(3)^3 + x(2)^(3*k+1);
1128	fk = jet(f- phi, 3*w[1], weight(JetId)) ;
1129	if(fk!=0) { return(printresult(84,f,"U["+string(12k)+"]",cstn,4k-3)); }
1130
1131	JetId = x(1)^3 + x(3)^3 + x(1)x(2)^(2k+1);
1132	fk = jet(f, 3*w[1], weight(JetId)) ;
1133	//--------------------------- step 85(k) --------------------------------
1134	if ( fk != phi ) {
1135	Jf = std(jacob(fk));
1136	Dim = dim(Jf);
1137	if(Dim==0) {return(printresult(86,f,"U["+string(k)+",0]",cstn,4*k-2));}
1138	if(Dim==1) {return(printresult(87,f,"U["+string(k)+",p]",cstn,4*k-2));}
1139	}
1140
1141	//--------------------------- step 88(k) --------------------------------
1142	JetId = x(1)^3 + x(3)^3 + x(2)^(3*k+2);
1143	fk = jet(f- phi, 3*w[1], weight(JetId)) ;
1144	if(fk!=0) {return(printresult(89,f,"U["+string(12k+4)+"]",cstn,4k-2));}
1145
1146	//--------------------------- step 90(k) --------------------------------
1147	k++;
1148	JetId = x(1)^3 + x(3)^3 + x(2)^(3*k);
1149	fk = jet(f, 3*w[1], weight(JetId)) ;
1150	Jf = std(jacob(fk));
1151	Dim = dim(Jf);
1152	Mult = mult(Jf);
1153	if ( Dim == 0 ) { return(printresult(83, f, "NoClass", cstn, -1)); }
1154	if ( Dim == 1 ) {
1155	if ( Mult == 4 ) {
1156	if( fk - phi != 0) { // b!=0 und/oder b'!=0
1157	if( Coeff(fk,x(1)x(2), x(1)^2x(2)^k) == 0 ) { // b=0 und b'!=0
1158	a = (fk - Coeff(fk, x(1), x(1)^3)*x(1)^3) / x(1);
1159	v = Isomorphie_s82_z(f, a, k);
1160	}
1161	else {
1162	if( Coeff(fk,x(1)x(2)x(3), x(1)x(2)^kx(3)) == 0 ){
1163	// b!=0 und b'=0
1164	a = subst(fk, x(3), 0);
1165	v = Isomorphie_s82_x(f, a, k);
1166	}
1167	else {
1168	a = Coeff(fk,x(1)x(2)x(3), x(1)x(2)^kx(3));
1169	b = Coeff(fk,x(2)x(3), x(2)^(2k)*x(3));
1170	B = maxideal(1);
1171	B[rvar(x(1))] = x(1)-b/a*x(2)^k;
1172	Phi = ring_top,B;
1173	f = Phi(f);
1174	PhiG = Phi(PhiG);
1175	cstn[4] = PhiG;
1176	fk = jet(f, 3*w[1], w) ;
1177	a = (fk - Coeff(fk, x(1), x(1)^3)*x(1)^3) / x(1);
1178	v = Isomorphie_s82_z(f, a, k);
1179	} // ende else b!=0 und b'=0
1180	} // ende else b=0 und b'!=0
1181	f, Phi = v[1..2];
1182	PhiG = Phi(PhiG);
1183	cstn[4] = PhiG;
1184	} //ende fk-phi!=0
1185	} // ende mult=4
1186	else { return(printresult(83, f, "NoClass", cstn, -1)); }
1187	} // ende dim=1
1188	else { return(printresult(83, f, "NoClass", cstn, -1)); }
1189	} // ENDE While
1190	return(printresult(83, f, "error!", cstn, -1));
1191	}
1192
1193	///////////////////////////////////////////////////////////////////////////////
1194	static proc Funktion97 (poly f, list cstn)
1195	{
1196	//---------------------------- initialisation ---------------------------------
1197	poly f3, l1, l2, a, b, c, prod;
1198	ideal Jfsyz, Jf, B;
1199	int k, i, pt, kx, ky, kz, Dim, Mult, Mu;
1200	matrix Mat;
1201	map PhiG, VERT;
1202	def ring_top=basering;
1203	debug_log(1, " Step 97");
1204
1205	Mu = cstn[2];
1206	PhiG = cstn[4];
1207	kx = 1; // Koordinate x
1208	ky = 2; // Koordinate y
1209	kz = 3; // Koordinate z
1210	B = maxideal(1); // Abbildungs-ideal
1211	pt = 2;
1212	f3 = jet(f, 3);
1213	k = 1;
1214
1215	//--------------------------- compute f3 ~ x2y --------------------------------
1216	// vertausche 2 Koordinaten sodass d2f/dx2 <>0 ist.
1217	for(i=1;i<4;i=i+1) {
1218	if(diff(diff(f3, x(i)), x(i)) != 0) { kx = i; i=4; }
1219	}
1220	if(kx == 2) { ky = 1; kz = 3; }
1221	if(kx == 3) { ky = 2; kz = 1; }
1222	//-------------------------- compute -l1*l2 -------------------------------
1223	f3 = jet(f, 3);
1224	Jfsyz = f3, diff(f3, x(kx));
1225	Mat = matrix(syz(Jfsyz));
1226	if(deg(Mat[2,1])>1) {
1227	Jfsyz = f3, Mat[2,1];
1228	Mat = matrix(syz(Jfsyz));
1229
1230	// berechen Abb. sodass f=x2*l2
1231	l1 = Mat[2,1];
1232	a = Coeff(l1, x(kx), x(kx));
1233	l1 = l1 / number(a);
1234	b = Coeff(l1, x(ky), x(ky));
1235	c = Coeff(l1, x(kz), x(kz));
1236	B[rvar(x(kx))] = x(kx) - b * x(ky) - c * x(kz);
1237	VERT = ring_top, B;
1238	f = VERT(f);
1239	PhiG = VERT(PhiG);
1240	cstn[4] = PhiG;
1241
1242	f3 = jet(f, 3);
1243	l2 = f3 / x(kx)^2;
1244
1245	// sorge dafuer, dass b<>0 ist.
1246	b = Coeff(l2, x(ky), x(ky));
1247	if( b== 0) { ky, kz = swap(ky, kz); }
1248
1249	// Koordinaten-Transf. s.d. f=x2y
1250	b = Coeff(l2, x(ky), x(ky));
1251	l2 = l2 / number(b);
1252	a = Coeff(l2, x(kx), x(kx));
1253	c = Coeff(l2, x(kz), x(kz));
1254	B = maxideal(1);
1255	B[rvar(x(ky))] = -a * x(kx) + x(ky) - c * x(kz);
1256	VERT = ring_top, B;
1257	f = VERT(f);
1258	PhiG = VERT(PhiG);
1259	cstn[4] = PhiG;
1260	}
1261
1262	//------------------------------- step 98 ---------------------------------
1263	Jfsyz = x(kx)^2*x(ky) + x(ky)^4 + x(kz)^4;
1264	a = jet(f, 8, weight(Jfsyz));
1265	Jf = std(jacob(a));
1266	Dim = dim(Jf);
1267	Mult = mult(Jf);
1268	if( Dim == 0) { return(printresult(99, f, "V[1,0]", cstn, 3)); }
1269	if( Dim == 1) {
1270	if(Mult==1) {return(printresult(100,f,"V[1,"+string(Mu-15)+"]",cstn,3));}
1271	if(Mult==2){return(printresult(101,f,"V#[1,"+string(Mu-15)+"]",cstn,3));}
1272	}
1273	" Dim=",Dim," Mult=",Mult;
1274	return(printresult(102, f, "NoClass", cstn, -1));
1275	}
1276
1277	///////////////////////////////////////////////////////////////////////////////
1278	proc tschirnhaus (poly f, poly x)
1279	"USAGE: tschirnhaus();"
1280	{
1281	//---------------------------- initialisation ---------------------------------
1282	poly b;
1283	ideal B;
1284	int n, j, hc;
1285	matrix cf;
1286	intvec z;
1287	string s;
1288	list v;
1289	map Phi, EH;
1290	def ring_top=basering;
1291
1292	n = nvars(basering);
1293	cf = coeffs(f, x);
1294	hc = nrows(cf) - 1; // hoechster exponent von x_i
1295	b = cf[hc+1,1]; // koeffizient von x_i^hc
1296	B = maxideal(1);
1297	z[n] = 0;
1298	EH = ring_top, z;
1299	Phi = ring_top, B;
1300	v[1] = f;
1301	if ( EH(b) != 0) // pruefe ob der Koeff von x_i^hc
1302	{ B[rvar(x)] = x -1(cf[hc,1]/(hcb));
1303	v[1] = Phi(f);
1304	}
1305	v[2] = Phi;
1306	return(v);
1307	}
1308
1309	///////////////////////////////////////////////////////////////////////////////
1310	static proc Isomorphie_s17 (poly f, poly fk, int k, int ct, list #)
1311	{
1312	//---------------------------- initialisation ---------------------------------
1313	ideal Jfsyz, JetId, bb;
1314	poly a, b, c, d, Relation, an, bn;
1315	int B,C, alpha, beta, gamma, g;
1316	matrix Matx, Maty;
1317	map PhiG, VERT;
1318	list v;
1319	def ring_top=basering;
1320
1321	if(size(#)==1) { PhiG = #[1]; }
1322	else { PhiG = ring_top,maxideal(1); }
1323	bb = maxideal(1);
1324
1325	// Ziel: bestimme a,b,c,d sodass fk = (ax+by^k)^3(cx+dy) gilt.
1326	debug_log(2, "Isomorphie_s17:");
1327	debug_log(2, "Faktor: f=",Show(f)," Jet=",Show(fk)," k=",k, "cnt=", ct);
1328
1329	if( k == 1) {
1330	Jfsyz = fk, diff(fk, x(1));
1331	Matx = matrix(syz(Jfsyz));
1332	Jfsyz = fk, diff(fk, x(2));
1333	Maty = matrix(syz(Jfsyz));
1334
1335	a = Coeff(fk, x(1), x(1)^4);
1336	b = Coeff(fk, x(2), x(2)^4);
1337	c = Coeff(fk, x(1)x(2), x(1)^3x(2));
1338	d = Coeff(fk, x(1)x(2), x(1)x(2)^3);
1339
1340	if( (a != 0) && (b != 0) ) {
1341	B = -int(Coeff(Matx[1,1], x(2), x(2)));
1342	C = -int(Coeff(Maty[1,1], x(1), x(1)));
1343	alpha = int(Coeff(Matx[2,1], x(1), x(1)^2));
1344	beta = int(Coeff(Matx[2,1], x(1)x(2), x(1)x(2)));
1345	gamma = int(Coeff(Matx[2,1], x(2), x(2)^2));
1346
1347	bb[rvar(x(1))] = x(1) - (2gamma / (B - beta))x(2);
1348	bb[rvar(x(2))] = x(2) - ((C - beta) / (2gamma))x(1);
1349	VERT = ring_top,bb;
1350	Relation = VERT(f);
1351	fk = jet(Relation, 4);
1352
1353	an = Coeff(fk, x(1), x(1)^4);
1354	bn = Coeff(fk, x(2), x(2)^4);
1355	if( (an != 0) & (bn != 0) ) { VERT=ring_top,x(1),(x(2) + a*x(1))/ b; }
1356	f = VERT(f);
1357	fk = jet(f, 4);
1358	PhiG = VERT(PhiG);
1359
1360	a = Coeff(fk, x(1), x(1)^4);
1361	b = Coeff(fk, x(2), x(2)^4);
1362	c = Coeff(fk, x(1)x(2), x(1)^3x(2));
1363	d = Coeff(fk, x(1)x(2), x(1)x(2)^3);
1364	Jfsyz = fk, diff(fk, x(1));
1365	Matx = matrix(syz(Jfsyz));
1366	Jfsyz = fk, diff(fk, x(2));
1367	Maty = matrix(syz(Jfsyz));
1368	}
1369
1370	if( (a == 0) \|\| (b == 0) ) {
1371	if( a == 0) {
1372	if( c == 0) { // y3(ax+by)
1373	Relation = - Matx[2,1] / Matx[1,1];
1374	a = Coeff(Relation, x(1), x(1));
1375	b = Coeff(Relation, x(2), x(2));
1376	VERT=ring_top,ax(2)^k - bx(1), x(1);
1377	}
1378	else { // (ax+by)^3y
1379	Relation = - 3*Matx[2,1] / Matx[1,1];
1380	a = Coeff(Relation, x(1), x(1));
1381	b = Coeff(Relation, x(2), x(2));
1382	VERT=ring_top,ax(1) - bx(2), x(2);
1383	}
1384	}
1385	else {
1386	if( d == 0) { // x3(ax+by)
1387	Relation = - Maty[2,1] / Maty[1,1];
1388	a = Coeff(Relation, x(1), x(1));
1389	b = Coeff(Relation, x(2), x(2));
1390	VERT=ring_top,x(1), bx(2)^k - ax(1);
1391	}
1392	else { // x(ax+by)^3
1393	Relation = - 3*Maty[2,1] / Maty[1,1];
1394	a = Coeff(Relation, x(1), x(1));
1395	b = Coeff(Relation, x(2), x(2));
1396	VERT=ring_top,x(2), bx(1) - ax(2);
1397	}
1398	}
1399	f = VERT(f);
1400	PhiG = VERT(PhiG);
1401	}
1402	else { // "Weder b noch a sind 0";
1403	if(ct > 5) { v[1]=f; v[2]=PhiG; return(v); }
1404	fk = jet(f, 4);
1405	return(Isomorphie_s17(f, fk, k, ct+1, PhiG));
1406	}
1407	}
1408	else { // k >1
1409	a = fk/x(2);
1410	Jfsyz = a, diff(a, x(1));
1411	Matx = matrix(syz(Jfsyz));
1412	Relation = -3 * Matx[2,1] / Matx[1,1];
1413	a = Coeff(Relation, x(1), x(1));
1414	b = Coeff(Relation, x(2), x(2)^k);
1415	VERT = basering,x(1)-b*x(2)^k,x(2);
1416	f = VERT(f);
1417	PhiG = VERT(PhiG);
1418	JetId = x(1)^3x(2) + x(2)^(3k+1);
1419	fk = jet(f, 3*k+1, weight(JetId));
1420	}
1421	v = f, PhiG;
1422	debug_log(2, "Isomorphie_s17: done");
1423	debug_log(2, "Faktor: f=",Show(f)," Jet=",Show(fk)," k=",k);
1424
1425	return(v);
1426	}
1427
1428	///////////////////////////////////////////////////////////////////////////////
1429	static proc printresult (int step, poly f, string typ, list cstn, int m)
1430	{
1431	//---------------------------- initialisation ---------------------------------
1432	int corank, Mu, K;
1433	list v;
1434
1435	corank, Mu, K = cstn[1..3];
1436	debug_log(0," Arnold step number "+string(step));
1437	if( @DeBug >= 0 )
1438	{
1439	"The singularity";
1440	" "+Show(jet(f, K))+"";
1441	if( typ != "error!" && typ != "NoClass" )
1442	{
1443	"is R-equivalent to "+typ+".";
1444	}
1445	if ( typ == "NoClass" )
1446	{
1447	"is not in Arnolds list.";
1448	}
1449	// if(K>=0) { " det = "+string(K); }
1450	if(Mu>=0) { " Milnor number = "+string(Mu); }
1451	if(m>=0) { " modality = "+string(m); }
1452	}
1453	v = f, typ, corank, cstn[4];
1454	return(v);
1455	}
1456
1457	///////////////////////////////////////////////////////////////////////////////
1458	static proc Funktion47 (poly f, list cstn)
1459	{
1460	int corank = cstn[1];
1461	int Mu = cstn[2];
1462	int K = cstn[3];
1463	string s = "The Singularity ";+Show(jet(f, K), corank, K);
1464	string tp="";
1465	// return(printresult(47, f, tp, cstn, -1));
1466
1467	s = s +" has 4-jet equal to zero. (F47), mu="+string(Mu);
1468
1469	s; // +" ("+SG_Typ+")";
1470	return(Show(f), tp, corank);
1471	}
1472
1473	///////////////////////////////////////////////////////////////////////////////
1474	static proc Funktion91 (poly f, list cstn, int k)
1475	{
1476	string tp = "U*[k,0]";
1477	return(printresult(91, f, tp, cstn, -1));
1478	}
1479
1480	///////////////////////////////////////////////////////////////////////////////
1481	static proc Funktion92 (poly f, list cstn, int k)
1482	{
1483	string tp = "UP[k]";
1484	return(printresult(92, f, tp, cstn, -1));
1485	}
1486
1487	///////////////////////////////////////////////////////////////////////////////
1488	static proc Funktion93 (poly f, list cstn, int k)
1489	{
1490	string tp = "UQ[k]";
1491	return(printresult(93, f, tp, cstn, -1));
1492	}
1493
1494	///////////////////////////////////////////////////////////////////////////////
1495	static proc Funktion94 (poly f, list cstn, int k)
1496	{
1497	string tp = "UR[k]";
1498	return(printresult(94, f, tp, cstn, -1));
1499	}
1500
1501	///////////////////////////////////////////////////////////////////////////////
1502	static proc Funktion95 (poly f, list cstn, int k)
1503	{
1504	string tp = "US[k]";
1505	return(printresult(95, f, tp, cstn, -1));
1506	}
1507
1508	///////////////////////////////////////////////////////////////////////////////
1509	static proc Funktion96 (poly f, list cstn, int k)
1510	{
1511	string tp = "UT[k]";
1512	return(printresult(96, f, tp, cstn, -1));
1513	}
1514
1515	///////////////////////////////////////////////////////////////////////////////
1516	proc morsesplit(poly f)
1517	"
1518	USAGE: morsesplit(f); f=poly
1519	RETURN: Normal form of f in M^3 after application of the splitting lemma
1520	COMPUTE: apply the splitting lemma (generalized Morse lemma) to f
1521	EXAMPLE: example morsesplit; shows an example"
1522	{
1523	//---------------------------- initialisation ---------------------------------
1524	poly f_out;
1525	int n, K, Mu, corank;
1526	list v;
1527
1528	if(defined(@ringdisplay) != 0 ) { kill @ringdisplay; }
1529	string @ringdisplay = "setring "+nameof(basering);
1530	export @ringdisplay;
1531
1532	def ring_ext=basering;
1533
1534	n = nvars(basering);
1535
1536	// if trace/debug mode not set, do it!
1537	init_debug();
1538	K, Mu, corank = basicinvariants(f);
1539	ring ring_top=char(basering),(x(1..n)),(c,ds);
1540
1541	map Conv=ring_ext,maxideal(1);
1542	setring ring_top;
1543	v = Morse(jet(Conv(f),K), K, corank, 0);
1544	poly f_out = v[1];
1545	setring ring_ext;
1546	map ConvUp = ring_top, maxideal(1);
1547	return(ConvUp(f_out));
1548	}
1549	example
1550	{ "EXAMPLE"; echo=2;
1551	ring r=0,(x,y,z),ds;
1552	export r;
1553	init_debug(1);
1554	poly f=(x2+3y-2z)^2+xyz-(x-y3+x2*z3)^3;
1555	poly g=morsesplit(f);
1556	g;
1557	}
1558
1559	///////////////////////////////////////////////////////////////////////////////
1560	static proc Coeffs (list #)
1561	{
1562	matrix m=matrix(coeffs(#[1],#[2]), deg(#[1])+1, 1);
1563	return(m);
1564	}
1565
1566	///////////////////////////////////////////////////////////////////////////////
1567	static proc Morse(poly f, int K, int corank, int ShowPhi)
1568	{
1569	//---------------------------- initialisation ---------------------------------
1570	poly fc, f2, a, P, Beta, fi;
1571	ideal Jfx, B;
1572	int n, i, j, k, Rang, Done;
1573	matrix Mat;
1574	map Id, Psi, Phi, PhiG;
1575	intvec Abb, RFlg;
1576	list v;
1577
1578	fi = f;
1579	n = nvars(basering);
1580	init_debug();
1581
1582	def ring_top=basering;
1583
1584	debug_log(3, "Split the polynomial below using determinacy: ", string(K));
1585	debug_log(3, Show(fi));
1586
1587	for( j=1; j<=n ; j++) { Abb[j] = 0; }
1588
1589	RFlg = GetRf(fi, n);
1590	debug_log(2, "Permutations:", RFlg );
1591	PhiG=ring_top,maxideal(1);
1592
1593	//----------------- find quadratic term, if there is only one -----------------
1594	B = maxideal(1);
1595	if(corank == (n-1))
1596	{
1597	Done = 0;
1598	f2 = jet(fi, 2);
1599	j = 1;
1600	Jfx = f2, diff(f2, x(j));
1601	while(j<=n && (diff(f2, x(j))==0))
1602	{
1603	j++;
1604	Jfx = f2, diff(f2, x(j));
1605	}
1606	Mat = matrix(syz(Jfx));
1607	Beta = 2*Mat[2,1]/Mat[1,1];
1608	for( j=1; j<=n ; j++)
1609	{
1610	f2 = Coeff(Beta, x(RFlg[j]), x(RFlg[j]));
1611	if(f2!=0)
1612	{
1613	k = RFlg[j];
1614	break;
1615	}
1616	}
1617	for( j=1; j<=n ; j=j+1)
1618	{
1619	f2 = Coeff(Beta, x(j), x(j));
1620	if(j == k) { B[rvar(x(j))] = (2f2x(j)-Beta) / number(f2); }
1621	}
1622	Phi =ring_top,B;
1623	fi = Phi(fi);
1624	PhiG = Phi(PhiG);
1625	}
1626	if( ShowPhi > 1) { PhiG; }
1627
1628	//------------------------ compute spliting lemma -----------------------------
1629	fc = fi;
1630	i = 1; // Index fuer Variablen wird bearbeitet
1631	while( i <= n)
1632	{
1633	Phi=ring_top,maxideal(1);
1634	debug_log(6, "Pruefe Variable x(" +string(RFlg[i]) + ")");
1635	debug_log(6, "--------------------");
1636	j = i + 1; // setze j fuer evtle Verschiebung
1637	f2 = jet(fc,2);
1638	debug_log(6, "Rechne 2-Jet =" , string(f2));
1639	if( (f2 - subst(f2, x(RFlg[i]), 0)) == 0 ) { Abb[RFlg[i]] = 1; }
1640	if( (f2 - subst(f2, x(RFlg[i]), 0)) != 0 )
1641	{
1642	while( (j<=n) \|\| (i==n) )
1643	{
1644	debug_log(6, "Pruefe 2-Jet mit Wert : " + string(jet(fc,2)));
1645	a = Coeff(jet(fc,2), x(RFlg[i]), x(RFlg[i])^2);
1646	debug_log(6,"Koeffizient von x(" + string(RFlg[i]) + ")^2 ist:", a);
1647	if( (a != 0) \|\| (i==n) )
1648	{
1649	debug_log(6, "BREAK!!!!!!!!!!!!!!");
1650	break;
1651	}
1652	debug_log(6,"Verschiebe evtl Variable x(",string(RFlg[j]),") um x(",
1653	string(RFlg[i]), ")");
1654	B = maxideal(1);
1655	for( k=1; k<=n ; k=k+1)
1656	{
1657	if(k==RFlg[j]) { B[rvar(x(k))] = x(k) + x(RFlg[i]); }
1658	}
1659	Phi = ring_top,B;
1660	fc = Phi(fi);
1661	j++;
1662	} // Ende while( (j<=n) \|\| (i==n) )
1663
1664	debug_log(6, "Moegliche Verschiebung fertig!");
1665	PhiG = Phi(PhiG);
1666	if( ShowPhi > 1) { "NachVersch.:"; Phi; }
1667
1668	if( (j<=n) \|\| (i==n))
1669	{
1670	P = Coeff(fc, x(RFlg[i]), x(RFlg[i]));
1671	debug_log(6, "Koeffizient von x("+string(RFlg[i])+") ist: ",
1672	string(P));
1673	if(P != 0)
1674	{
1675	debug_log(6, "1 Koeffizient von x("+string(RFlg[i])+") ist: ",
1676	string(P));
1677	debug_log(6, "a=" + string(a));
1678	P = P / number (2 * a);
1679	debug_log(6, "2 Koeffizient von x("+string(RFlg[i])+") ist: ",
1680	string(P));
1681	B = maxideal(1);
1682	for( k=1; k<=n ; k=k+1) {if(k==RFlg[i]) {B[rvar(x(k))]=x(k)-P;}}
1683	Phi =ring_top,B;
1684	debug_log(6, "Quadratische-Ergaenzung durch:", Phi);
1685	fi = Phi(fc);
1686	PhiG = Phi(PhiG);
1687	fc = jet(fi,K);
1688	P = Coeff(fc, x(RFlg[i]), x(RFlg[i]));
1689	if( P != 0)
1690	{
1691	fi = fc;
1692	continue;
1693	}
1694	} // Ende if(P != 0)
1695	// Fertig mit Quadratischer-Ergaenzung
1696	} // Ende if( (j<=n) \|\| (i==n))
1697	} // Ende if( (f2 - subst(f2, x(RFlg[i]), 0)) != 0 )
1698
1699	fi = fc;
1700	i++;
1701	debug_log(6, "++++++++++++++++++++++++++++++++++++++++++++++++++++++++");
1702	}
1703	debug_log(6, "Ende ---------------------------------------------------");
1704
1705	//--------------------------- collect results ---------------------------------
1706	if( ShowPhi > 0 )
1707	{
1708	"Abbildung innerhalb des Morse-Lemmas:";
1709	PhiG;
1710	"Vergleich:";
1711	"PhiG(f)= " + Show(jet(PhiG(f), K));
1712	"fi = " + Show(fi);
1713	}
1714
1715	Rang = 0;
1716	B = maxideal(1);
1717	for( i=1; i<=n ; i++) { if(Abb[i] != 1) { Rang ++; B[rvar(x(i))] = 0; } }
1718	Phi = ring_top,B;
1719	PhiG = Phi(PhiG);
1720	fi = Phi(fi);
1721	v = fi, PhiG;
1722	debug_log(2, "rank determined with Morse rg=", Rang);
1723	debug_log(1, "Residual singularity f=",Show(fi));
1724	return(v);
1725	}
1726
1727	///////////////////////////////////////////////////////////////////////////////
1728	static proc Coeff(poly f, list #)
1729	{
1730	//---------------------------- initialisation ---------------------------------
1731	poly a, term;
1732	int n, i;
1733	matrix K;
1734
1735	n = nvars(basering);
1736	i = 1;
1737	term = #[2];
1738	K = coef(f, #[1]);
1739
1740	while( (i<=ncols(K)) && (K[1,i] != term) )
1741	{ i++;
1742	if(i>ncols(K)) { break; }
1743	}
1744	if(i<=ncols(K)) { a = K[2,i]; }
1745	if(i>ncols(K)) { a = 0; }
1746
1747	return(a);
1748	}
1749
1750	///////////////////////////////////////////////////////////////////////////////
1751	static proc ReOrder(poly f)
1752	{
1753	//---------------------------- initialisation ---------------------------------
1754	poly result;
1755	ideal B = maxideal(1);
1756	int i, n, Ctn, Ctv;
1757	map conv;
1758
1759	n = nvars(basering);
1760	Ctv = 1; // Zahl der Vorhandenen Variablen
1761	Ctn = n; // Zahl der Nicht-Vorhandenen Variablen
1762	def ring_top=basering;
1763
1764	for( i=1; i<=n; i=i+1)
1765	{ result = subst(f,x(i), 0) - f;
1766	if( result != 0 ) { B[rvar(x(i))] = x(Ctv); Ctv++; }
1767	else { B[rvar(x(i))] = x(Ctn); Ctn--; }
1768	}
1769
1770	conv = ring_top,B;
1771	return(conv);
1772	}
1773
1774	///////////////////////////////////////////////////////////////////////////////
1775	proc quickclass(poly f)
1776	"
1777	USAGE: quickclass(f); f=poly
1778	RETURN: Normal form of f in Arnold's list
1779	REMARK: try to determine the normal form of f by invariants, mainly by
1780	computing the Hilbert function of the Milnor algebra,
1781	no coordinate change is needed (see also proc 'milnorcode').
1782	EXAMPLE: example quickclass; shows an example"
1783	{
1784	//---------------------------- initialisation ---------------------------------
1785	string Typ;
1786	int cnt, K, Mu, corank;
1787	list v;
1788	def ring_top=basering;
1789	// check basic condition on the basering.
1790	if(checkring()) { return(f); }
1791	if( f==0 ) {
1792	"Normal form : 0";
1793	return(f);
1794	}
1795	if( jet(f,0)!=0 ) {
1796	"Normal form : 1";
1797	return(f);
1798	}
1799	K, Mu, corank = basicinvariants(f);
1800	if(Mu<0) {
1801	debug_log(0, "The Milnor number of the function is infinite.");
1802	return(f);
1803	}
1804
1805	// Do the classification of f
1806	// typ: list of types matching the milnorcode
1807	// cnt: number of matches found
1808	v = HKclass(milnorcode(f));
1809	Typ, cnt = v[1..2];
1810	"Singularity R-equivalent to :",Typ;
1811	if(cnt==0) {
1812	"Hilbert polynomial not recognised. Milnor code = ", milnorcode(f);
1813	return();
1814	}
1815	if(cnt==1) {
1816	debug_log(1,"Getting normal form from database.");
1817	"normal form :",A_L(Typ);
1818	return(A_L(Typ));
1819	}
1820	// Hier nun der Fall cnt>1
1821	"Hilbert-Code of Jf^2";
1822	"We have ", cnt, "cases to test";
1823	Cubic(f);
1824	return(v);
1825	}
1826	example
1827	{ "EXAMPLE:"; echo=2;
1828	ring r=0,(x,y,z),ds;
1829	poly f=(x2+3y-2z)^2+xyz-(x-y3+x2*z3)^3;
1830	quickclass(f);
1831	}
1832
1833	///////////////////////////////////////////////////////////////////////////////
1834	proc milnorcode (poly f, list #)
1835	"USAGE: milnorcode(f[,e]); f=poly, e=int
1836	RETURN: intvec, coding the Hilbert function of the e-th Milnor algebra
1837	of f, i.e. of basering/(jacob(f)^e) (default e=1), according
1838	to proc Hcode
1839	EXAMPLE: example milnorcode; shows an example"
1840	{
1841	int e=1;
1842	if(size(#)==1) { e=#[1]; }
1843	ideal jf=std(jacob(f)^e);
1844	intvec v=hilb(jf,2);
1845
1846	return(Hcode(v));
1847	}
1848	example
1849	{ "EXAMPLE:"; echo=2;
1850	ring r=0,(x,y,z),ds;
1851	poly f=x2y+y3+z2;
1852	milnorcode(f);
1853	milnorcode(f,2); // a big second argument may result in memory overflow
1854	}
1855
1856	///////////////////////////////////////////////////////////////////////////////
1857	proc Hcode (intvec v)
1858	"USAGE: Hcode(v); v=intvec
1859	RETURN: intvec, coding v according to the number of successive
1860	repetitions of an entry
1861	EXAMPLE: example Hcode; shows an example."
1862	{
1863	int col, len, i, cur, cnt, maxcoef, nlen;
1864	intvec hil1, hil2;
1865
1866	col = 1;
1867	len = size(v);
1868	v[len+1] = 0;
1869
1870	init_debug();
1871	debug_log(1, "Hcode:", v );
1872
1873	for(i=1; i<=len; i++) { if( v[i] > maxcoef) { maxcoef = v[i]; } }
1874
1875	nlen = 2*maxcoef - 1;
1876	hil1[nlen] = 0;
1877	hil2[nlen] = 0;
1878
1879	for(i=1; i<=nlen; i++)
1880	{ if( i > maxcoef) { hil2[i] = 2*maxcoef-i; }
1881	else { hil2[i] = i; }
1882	}
1883
1884	for(i=1; i<=nlen; i++)
1885	{ cnt=0;
1886	while( (col<=len) && (v[col] == hil2[i]) )
1887	{ cnt++; col++; }
1888	hil1[i] = cnt;
1889	}
1890	return(hil1);
1891	}
1892	example
1893	{ "EXAMPLE:"; echo=2;
1894	intvec v1 = 1,3,5,5,2;
1895	Hcode(v1);
1896	intvec v2 = 1,2,3,4,4,4,4,4,4,4,3,2,1;
1897	Hcode(v2);
1898	}
1899
1900	///////////////////////////////////////////////////////////////////////////////
1901	static proc Cubic (poly f)
1902	{
1903	//---------------------------- initialisation ---------------------------------
1904	poly f3;
1905	ideal Jf, Jf1, Jf2;
1906	int Dim, Mult;
1907
1908	f3 = jet(f, 3);
1909	if( jet(f,2) != 0) { return("2-jet non zero"); }
1910	if( f3 == 0 ) { return("null form"); }
1911
1912	Jf1 = jacob(f3);
1913	Jf = std(Jf1);
1914	Dim = dim(Jf);
1915	Mult = mult(Jf);
1916
1917	if(Dim == 0) { return("P[8]:smooth cubic"); } // x3 + y3 + z3 + axyz
1918	if(Dim == 1) {
1919	if(Mult == 2) {
1920	Jf2 = wedge(jacob(Jf1),3-Dim), Jf1;
1921	Jf2 = std(Jf2);
1922	Dim = dim(Jf2);
1923	if (Dim == 0) { return("R:conic + line"); } // x3 + xyz
1924	if (Dim == 1) { return("Q:cuspidal cubic"); } // x3 + yz2
1925	}
1926	if(Mult == 3) {
1927	Jf2 = wedge(jacob(Jf1),3-Dim), Jf1;
1928	Jf2 = std(Jf2);
1929	Dim = dim(Jf2);
1930	if(Dim == 0) { return("T:three lines"); } // xyz
1931	if(Dim == 1) { return("S:conic + tangent"); } // x2z + yz2
1932	}
1933	if(Mult == 4) { return("U:three concurrent lines"); } // x3 + xz2
1934	}
1935	if(Dim == 2) {
1936	if(Mult == 1) { return("V:doubleline + line"); } // x2y
1937	if(Mult == 2) { return("V': tripple line"); } // x3
1938	}
1939	if(Dim == 3) { return("P[9]:nodal cubic"); } // x3 + y3 + xyz
1940
1941	return("");
1942	}
1943
1944	///////////////////////////////////////////////////////////////////////////////
1945	static proc parity (int e)
1946	"USAGE: parity()"
1947	{
1948	int r = e/2;
1949	if( 2*r == e ) { return(0); }
1950	return(1);
1951	}
1952
1953	///////////////////////////////////////////////////////////////////////////////
1954	static proc HKclass (intvec sg)
1955	{
1956	//---------------------------- initialisation ---------------------------------
1957	int cnt = 0;
1958	string SG_Typ = "";
1959	list v;
1960
1961	// if trace/debug mode not set, do it!
1962	init_debug();
1963	debug_log(1, "Milnor code : ", sg );
1964	if(size(sg) == 1) { v ="A["+string(sg[1])+"]", 1; return(v); }
1965	if(size(sg) == 3) { return(HKclass3(sg, SG_Typ, cnt)); }
1966	if(size(sg) == 5) { return(HKclass5(sg, SG_Typ, cnt)); }
1967	if(size(sg) == 7) { return(HKclass7(sg, SG_Typ, cnt)); }
1968	debug_log(1, "No solution found." );
1969	v = "", 0;
1970	return(v);
1971	}
1972
1973	///////////////////////////////////////////////////////////////////////////////
1974	static proc HKclass3 (intvec sg, string SG_Typ, int cnt)
1975	{
1976	list v;
1977
1978	if(sg[1] == 1) { v = HKclass3_teil_1(sg, SG_Typ, cnt); }
1979	debug_log(6, "HKclass3: ", v[1], " cnt=", v[2]);
1980	return(v);
1981	}
1982
1983	///////////////////////////////////////////////////////////////////////////////
1984	static proc HKclass3_teil_1 (intvec sg, string SG_Typ, int cnt)
1985	{
1986	int k, r, s;
1987	list v;
1988
1989	debug_log(2, "entering HKclass3_teil_1", sg);
1990	if(sg[2]==1) { SG_Typ=SG_Typ+" D[k]=D["+string(sg[3]+3)+"]";cnt++;} // D[k]
1991	if(sg[2]>=1) {
1992	if( parity(sg[2])) { // sg[2] ist ungerade
1993	if(sg[2]<=sg[3]) {
1994	k = (sg[2]+1)/2;
1995	if(k>1) {
1996	cnt++;
1997	SG_Typ=SG_Typ+" J[k,r]=J["+string(k)+","+string(sg[3]+1-2*k)+"]";
1998	}// J[k,r]
1999	}
2000	if(sg[2]==sg[3]+2) { // E[6k+2]
2001	k = (sg[2]-1)/2;
2002	if(k>0) {cnt++; SG_Typ=SG_Typ+" E[6k+2]=E[" + string(6*k+2) + "]"; }
2003	}
2004	}
2005	else { // sg[2] ist gerade
2006	if( sg[2] == sg[3]+1) { // E[6k]
2007	k = sg[2]/2; cnt++; SG_Typ=SG_Typ + " E[6k]=E[" + string(6*k) + "]"; }
2008	if( sg[2] == sg[3]) { // E[6k+1]
2009	k=sg[2]/2; cnt++; SG_Typ=SG_Typ+" E[6k+1]=E["+string(6*k+1)+"]"; }
2010	}
2011	}
2012
2013	debug_log(2, "finishing HKclass3_teil_1");
2014	debug_log(6, "HKclass3: ", SG_Typ, " cnt=", cnt);
2015	v = SG_Typ, cnt;
2016	return(v);
2017	}
2018
2019	///////////////////////////////////////////////////////////////////////////////
2020	static proc HKclass5 (intvec sg, string SG_Typ, int cnt)
2021	{
2022	list v;
2023
2024	if(sg[1] == 1 && sg[2] == 1) { v = HKclass5_teil_1(sg, SG_Typ,cnt); }
2025	if(sg[1] == 1 && sg[2] == 0) { v = HKclass5_teil_2(sg, SG_Typ,cnt); }
2026	debug_log(6, "HKclass5: ", v[1], " cnt=", v[2]);
2027	return(v);
2028	}
2029
2030	///////////////////////////////////////////////////////////////////////////////
2031	static proc HKclass5_teil_1 (intvec sg, string SG_Typ, int cnt)
2032	{
2033	int k, r, s;
2034	list v;
2035
2036	debug_log(2, "entering HKclass5_teil_1", sg);
2037	if(parity(sg[3])) { // Dritte Stelle soll ungerade sein
2038	k = (sg[3]+1)/2;
2039	if(sg[3] > sg[4]) {
2040	k--;
2041	if( (sg[4]==sg[5]) && (sg[3] == sg[4]+1) && k>0 ) { // W[12k+6]
2042	SG_Typ = SG_Typ + " W[12k+6]=W["+string(12*k+6)+"]"; cnt++;
2043	}
2044	if( (sg[3]==sg[5]) && (sg[3] == sg[4]+2) && k>0 ) { // W[12k+5]
2045	SG_Typ = SG_Typ + " W[12k+5]=W["+string(12*k+5)+"]"; cnt++;
2046	}
2047	}
2048	else { // sg[3] <= sg[4]
2049	if( (sg[3]==sg[4]) && (sg[5] >= sg[3]) ) {
2050	r = sg[5] - sg[4];
2051	SG_Typ=SG_Typ +" X[k,r]=X["+string(k)+","+string(r)+"]"; cnt++;
2052	}
2053	if( (sg[3]==1) && (sg[4]==3) && (sg[5]>=sg[4])){ // Z[1,r]
2054	r = sg[5] - sg[4];
2055	SG_Typ = SG_Typ + " Z[1,r]=Z[1,"+string(r)+"]"; cnt++;
2056	}
2057
2058	if( sg[4] == sg[5]) {
2059	if(parity(sg[4])) { // Z[k,r,0]
2060	r = (sg[4] - sg[3])/2;
2061	if( r>0 ) { cnt++;
2062	SG_Typ = SG_Typ + " Z[k,r,0]=Z["+string(k)+","+string(r)+",0]";
2063	}
2064	}
2065	else { // Z[k,12k+6r]
2066	r = (sg[4] - 2*k)/2; cnt++;
2067	SG_Typ = SG_Typ+" Z[k,12k+6r]=Z["+string(k)+","+string(12k+6r)+"]";
2068	}
2069	}
2070
2071	if( parity(sg[4]) ) { // 4. Stelle ist ungerade
2072	if(sg[4] == sg[5]+2) { // Z[k,12k+6r+1]
2073	r = (sg[4]-2*k-1)/2; cnt++;
2074	SG_Typ=SG_Typ+" Z[k,12k+6r+1]=Z["+string(k)+",";
2075	SG_Typ=SG_Typ+string(12k+6r+1)+"]";
2076	}
2077	if( (sg[5]>sg[4]) && (sg[4]>sg[3]) ) { // Z[k,r,s]
2078	r = (sg[4] - sg[3])/2; cnt++;
2079	s = sg[5] - sg[4];
2080	SG_Typ = SG_Typ + " Z[k,r,s]=";
2081	SG_Typ = SG_Typ + "Z["+string(k)+","+string(r)+","+string(s)+"]";
2082	}
2083	}
2084	else { // 4. Stelle ist gerade
2085	if( sg[4] == sg[5]+1) { // Z[k,12k+6r-1]
2086	r = (sg[4] - 2*k)/2; cnt++;
2087	SG_Typ=SG_Typ+" Z[k,12k+6r-1]=Z["+string(k)+",";
2088	SG_Typ=SG_Typ+string(12k+6r-1)+"]";
2089	}
2090	}
2091
2092	if(sg[4]>sg[3]) { // Y[k,r,s]
2093	r = sg[4] - sg[3];
2094	s = sg[5] - sg[3] + r;
2095	if( s<0 ) { s = -s; }
2096	SG_Typ = SG_Typ + " Y[k,r,s]="; cnt++;
2097	SG_Typ = SG_Typ + "Y["+string(k)+","+string(r)+","+string(s)+"]";
2098	}
2099	}
2100	}
2101	else { // Dritte Stelle soll gerade sein
2102	k = sg[3]/2;
2103	// sortiere verschiedene W's
2104	if(k>0) {
2105	if( (sg[4]==2*k-1) && (sg[4]==sg[5]) ) { // W[12k]
2106	SG_Typ = SG_Typ + " W[12k]=W["+string(12*k)+"]"; cnt++;
2107	}
2108	if( (sg[4]==2*k-1) && (sg[3]==sg[5]) ) { // W[12k+1]
2109	SG_Typ = SG_Typ + " W[12k+1]=W["+string(12*k+1)+"]"; cnt++;
2110	}
2111	if( (sg[4]==2*k) && (sg[5]>=sg[4]) ) { // W[k,r]
2112	r = sg[5] - sg[4];
2113	SG_Typ=SG_Typ+" W[k,r]=W["+string(k)+","+string(r)+"]"; cnt++;
2114	}
2115	if( (sg[5]==2*k-1) && (sg[4]>sg[3]) ) { // W#[k,2r-1]
2116	r = sg[4] - sg[3]; cnt++;
2117	SG_Typ = SG_Typ + " W#[k,2r-1]=W["+string(k)+","+string(2*r-1)+"]";
2118	}
2119	if( (sg[5]==2*k) && (sg[4]>sg[3]) ) { // W#[k,2r]
2120	r = sg[4] - sg[3]; cnt++;
2121	SG_Typ = SG_Typ + " W#[k,2r]=W["+string(k)+","+string(2*r)+"]";
2122	}
2123	} // ENDIF k>0
2124	}
2125	debug_log(2, "finishing HKclass5_teil_1");
2126	debug_log(6, "HKclass5_teil_1: ", SG_Typ, " cnt=", cnt);
2127	v = SG_Typ, cnt;
2128	return(v);
2129	}
2130
2131	///////////////////////////////////////////////////////////////////////////////
2132	static proc HKclass5_teil_2 (intvec sg, string SG_Typ, int cnt)
2133	{
2134	int k, r, s;
2135	list v;
2136
2137	debug_log(2, "entering HKclass5_teil_2", sg);
2138	// finde T[p,q,r]
2139	k = sg[3] + 1;
2140	r = sg[4] + k;
2141	s = sg[5] + r - 1;
2142	if(k>2 && r>2 && s>2) { // T[k,r,s]
2143	cnt++;
2144	SG_Typ = SG_Typ + " T[k,r,s]=T["+string(k)+","+string(r)+","+string(s)+"]";
2145	}
2146
2147	// finde weitere Moeglicjkeiten.
2148	if(sg[3]==2) { // Q[...]
2149	if(parity(sg[4])) { // 4. Stelle ist ungerade.
2150	if(sg[4]==sg[5]) { // Q[6k+4]
2151	k=(sg[4]+1)/2; cnt++; SG_Typ=SG_Typ+" Q[6k+4]=Q["+string(6*k+4)+"]";
2152	}
2153	if(sg[4]+1==sg[5]) { // Q[6k+5]
2154	k=sg[5]/2; cnt++; SG_Typ=SG_Typ+" Q[6k+5]=Q["+string(6*k+5)+"]";
2155	}
2156	}
2157	else { // 4. Stelle ist gerade.
2158	if(sg[4]==sg[5]+1) { // Q[6k+6]
2159	k=sg[4]/2; cnt++; SG_Typ=SG_Typ+" Q[6k+6]=Q["+string(6*k+6)+"]";
2160	}
2161	if(sg[4]<sg[5]) { // Q[k,r]
2162	k = (sg[4]+2)/2;
2163	if(k>=2) {
2164	r=sg[5]+1-2*k; cnt++;
2165	SG_Typ=SG_Typ+" Q[k,r]=Q["+string(k)+","+string(r)+"]";
2166	}
2167	}
2168	}
2169	}
2170	else { // S[...]
2171	if(parity(sg[3])) { // 3. Stelle ist ungerade.
2172	k = (sg[3]-1)/2;
2173	if(sg[3]==sg[4]+3 && sg[3]==sg[5]+2) { // S[12k-1]
2174	cnt++; SG_Typ = SG_Typ + " S[12k-1]=S["+string(12*k-1)+"]";
2175	}
2176	if(sg[3]==sg[4]+3 && sg[3]==sg[5]+1) { // s[12k]
2177	cnt++; SG_Typ = SG_Typ + " S[12k]=S["+string(12*k)+"]";
2178	}
2179	if(sg[3]==sg[4]+2 && sg[5]>=sg[4]+1) { // S[k,r]
2180	r = sg[5] - 2*k; cnt++;
2181	SG_Typ = SG_Typ + " S[k,r]=S["+string(k)+","+string(r)+"]";
2182	}
2183	if(sg[3]==sg[5]+2 && sg[4]>=sg[5]) { // S#[k,2r-1]
2184	r = sg[4] - 2*k + 1; cnt++;
2185	SG_Typ = SG_Typ + " S#[k,2r-1]=S#["+string(k)+","+string(2*r-1)+"]";
2186	}
2187	if(sg[3]==sg[5]+1 && sg[4]>=sg[5]) { // S#[k,2r]
2188	r = sg[4] - 2*k + 1; cnt++;
2189	SG_Typ = SG_Typ + " S#[k,2r]=S#["+string(k)+","+string(2*r)+"]";
2190	}
2191	}
2192	else { // 3. Stelle ist gerade.
2193	if(sg[3]==sg[5]+1 && sg[5]==sg[4]+3) { // S[12k+4]
2194	k = (sg[3]-2)/2; cnt++;
2195	SG_Typ = SG_Typ + " S[12k+4]=S["+string(12*k+4)+"]";
2196	}
2197	if(sg[3]==sg[5]+2 && sg[5]==sg[4]+1) { // S[12k+5]
2198	k = (sg[3]-2)/2; cnt++;
2199	SG_Typ = SG_Typ + " S[12k+5]=S["+string(12*k+5)+"]";
2200	}
2201	}
2202	}
2203	debug_log(2, "finishing HKclass5_teil_2");
2204	debug_log(6, "HKclass5_teil_2: ", SG_Typ, " cnt=", cnt);
2205	v = SG_Typ, cnt;
2206	return(v);
2207	}
2208
2209	///////////////////////////////////////////////////////////////////////////////
2210	static proc HKclass7 (intvec sg, string SG_Typ, int cnt)
2211	{
2212	list v;
2213
2214	if(sg[1]==1 && sg[2]==0 && sg[3]==1) { v=HKclass7_teil_1(sg, SG_Typ, cnt); }
2215	debug_log(6, "HKclass7: ", v[1], " cnt=", v[2]);
2216	return(v);
2217	}
2218
2219	///////////////////////////////////////////////////////////////////////////////
2220	static proc HKclass7_teil_1 (intvec sg, string SG_Typ, int cnt)
2221	{
2222	int k, r, s;
2223	list v;
2224
2225	debug_log(2, "entering HKclass7_teil_1", sg);
2226	if(sg[4] == 2) { // V[...]
2227	if(sg[5] == 0 && sg[6] == 1 && sg[7]>0) { // V[1,r]
2228	r = sg[7] - 1; cnt++; SG_Typ = SG_Typ + " V[1,r]=V[1,"+string(r)+"]";
2229	}
2230	if(sg[5] == 1 && sg[7] == 1) { // V#[1,2r-1]
2231	r=sg[6]+1; cnt++; SG_Typ=SG_Typ+" V#[1,2r-1]=V#[1,"+string(2*r-1)+"]";
2232	}
2233	if(sg[5] == 1 && sg[7] == 2) { // V#[1,2r]
2234	r=sg[6]+1; cnt++; SG_Typ=SG_Typ+" V#[1,2r]=V#[1,"+string(2*r)+"]";
2235	}
2236	}
2237	// Moegliche U[...]'s
2238	k = sg[4];
2239	if(sg[5]==2*k-1 && sg[6]==0 && sg[7]==sg[5]) { // U[12k]
2240	cnt++;SG_Typ = SG_Typ + " U[12k]=U["+string(12*k)+"]";
2241	}
2242	if(sg[5]==2*k && sg[6]==0 && sg[7]==sg[5]) { // U[12k+4]
2243	cnt++;SG_Typ = SG_Typ + " U[12k+4]=U["+string(12*k+4)+"]";
2244	}
2245	if(sg[5]==2*k-1 && sg[6]>0 && sg[7]==sg[5]) { // U[k,2r-1]
2246	r=sg[6]-1; cnt++;
2247	SG_Typ=SG_Typ+" U[k,2r-1]=U["+string(k)+","+string(2*r-1)+"]";
2248	}
2249	if(sg[5]==2k-1 && sg[6]>0 && sg[7]==2k) { // U[k,2r]
2250	r = sg[6]; cnt++;
2251	SG_Typ = SG_Typ + " U[k,2r]=U["+string(k)+","+string(2*r)+"]";
2252	}
2253	debug_log(2, "finishing HKclass7_teil_1");
2254	debug_log(6, "HKclass7_teil_1: ", SG_Typ, " cnt=", cnt);
2255	v = SG_Typ, cnt;
2256	return(v);
2257	}
2258
2259	///////////////////////////////////////////////////////////////////////////////
2260	proc singularity(string typ, list #)
2261	"USAGE: singularity(t, l); t=string (name of singularity),
2262	l=list of integers/polynomials (indices/parmeters of singularity)
2263	COMPUTE: get the singularity named by type t from the database.
2264	list l is as follows: @*
2265	l= k [,r [,s [,a [,b [,c [,d]..]: k,r,s=int a,b,c,d=poly. @*
2266	The name of the dbm-databasefile is: NFlist.[dir,pag].
2267	The file is found in the current directory. If it does not
2268	exist, please run the script MakeDBM first.
2269	RETURN: Normal form and corank of the singularity named by type t and its
2270	index (indices) l.
2271	EXAMPLE: example singularity; shows an example"
2272	{
2273	poly a1, a2, a3, a4, f;
2274	int k, r, s;
2275	int len = size(#);
2276	list v, ret;
2277
2278	classify_init();
2279	ret = 0, 0;
2280	k = #[1];
2281	if(len>=2) { r = #[2]; }
2282	else { r = 0; }
2283	if(len>=3) { s = #[3]; }
2284	else { s = 0; }
2285	if( k<0 \|\| r<0 \|\| s<0) {
2286	"Initial condition failed: k>=0; r>=0; s>=0";
2287	"k="+string(k)+" r="+string(r)+" s="+string(s);
2288	return(ret);
2289	}
2290	int crk;
2291
2292	init_debug();
2293	def ring_top=basering;
2294
2295	if(len>=4) { a1 = #[4]; }
2296	else { a1=1; }
2297	if(len>=5) { a2 = #[5]; }
2298	else { a2=1; }
2299	if(len>=6) { a3 = #[6]; }
2300	else { a3=1; }
2301	if(len>=7) { a4 = #[7]; }
2302	else { a4=1; }
2303
2304	debug_log(4, "Values: len=", string(len), " k=", string(k), " r=",
2305	string(r));
2306	if(defined(RingNF) != 0 ) { kill RingNF; }
2307	ring RingNF=char(basering),(x,y,z),(c,ds);
2308	poly f;
2309	map Conv=ring_top,maxideal(1);
2310	v = Singularitaet(typ, k, r, s, Conv(a1), Conv(a2), Conv(a4), Conv(a4));
2311	f, crk = v[1..2];
2312	debug_log(2, "Info=", f );
2313	setring ring_top;
2314	if(defined(Phi) != 0 ) { kill Phi; }
2315	map Phi=RingNF,maxideal(1);
2316
2317	ret = Phi(f), crk;
2318	return(ret);
2319	}
2320	example
2321	{ "EXAMPLE"; echo=2;
2322	ring r=0,(x,y,z),(c,ds);
2323	init_debug(0);
2324	singularity("E[6k]",6);
2325	singularity("T[k,r,s]", 3, 7, 5);
2326	poly f=y;
2327	singularity("J[k,r]", 4, 0, 0, f);
2328	}
2329
2330	///////////////////////////////////////////////////////////////////////////////
2331	static proc Singularitaet (string typ,int k,int r,int s,poly a,poly b,poly c,poly d)
2332	{
2333	list v;
2334	string DBMPATH=system("getenv","DBMPATH");
2335	string DatabasePath, Database, S, Text, Tp;
2336	poly f, f1;
2337	int crk, Mu, ret;
2338	intvec MlnCd;
2339
2340	if( DBMPATH != "" ) { DatabasePath = DBMPATH+"/NFlist"; }
2341	else { DatabasePath = "NFlist"; }
2342	Database="DBM: ",DatabasePath;
2343
2344	link dbmLink=Database;
2345	debug_log(2, "Opening Singalarity-database: ", newline, Database);
2346	Tp = read(dbmLink, typ);
2347	debug_log(2,"DBMread(", typ, ")=", Tp, ".");
2348	if( Tp != "(null)" && Tp !="" ) {
2349	string Key = "I_", typ;
2350	S = "f = ", Tp, ";";
2351	debug_log(2,"S=", S, " Tp=", Tp, "Key=", Key);
2352	execute(S);
2353	execute(read(dbmLink, Key)+";");
2354	debug_log(1, "Polynom f=", f, " crk=", crk, " Mu=", Mu,
2355	" MlnCd=", MlnCd);
2356	v = f, crk, Mu, MlnCd;
2357	}
2358	else {
2359	v = 0, 0, 0, 0;
2360	}
2361	close(dbmLink);
2362	return(v);
2363	}
2364
2365	///////////////////////////////////////////////////////////////////////////////
2366	proc RandomPolyK (int M, string Typ)
2367	"USAGE: RandomPolyK(M, Typ)"
2368	{
2369	//---------------------------- initialisation ---------------------------------
2370	int n, b, i, k, r, s, crk;
2371	ideal B;
2372	map Phi;
2373	string txt, Tp;
2374	list v;
2375
2376	def ring_ext=basering;
2377	n=4;
2378	if(M<5) { M = 5; }
2379
2380	k = random(1, M);
2381	r = random(-5, 2*M);
2382	s = random(-5, 2*M);
2383	if(r<0) { r = 0; }
2384	if(s<0) { s = 0; }
2385
2386	ring RgAnf=char(basering),(x,y,z,t),(c,ds);
2387	poly f;
2388
2389	v = singularity(Typ, k, r, s);
2390	f, crk = v[1..2];
2391	// f = f +t2;
2392	if(crk==1) { f = f + y2 + z2; }
2393	if(crk==2) { f = f + z2; }
2394	txt="RandomPoly-Series: gewaehlt fall "+Typ+" mit";
2395	txt=txt+" f="+string(f);
2396	txt;
2397	setring ring_ext;
2398	B = maxideal(1);
2399
2400	r=1;
2401	for(i=n; i>0; i--,r++) {
2402	// for(i=1; i<=n; i=i+1)
2403	B[rvar(x(r))] = x(i);
2404	if(i>2 && random(1,10)<3) { B[rvar(x(r))] = B[rvar(x(r))] + x(i-1); }
2405	// if(i==1 && random(1,10)<4) { B[rvar(x(r))] = B[rvar(x(r))]- x(n); }
2406	if(i>0) {
2407	for(b=3; b<5; b=b+1) {
2408	// B[rvar(x(r))] = B[rvar(x(r))] + random(0,9) * x(i)^(b+2);
2409	if(random(1,20)<3) {
2410	B[rvar(x(r))] = B[rvar(x(r))] - random(-2,2)*x(b)^2;
2411	}
2412	}
2413	}
2414	}
2415	Phi=RgAnf, B;
2416	Phi;
2417	poly fr=Phi(f);
2418	fr = fr+(x(1)+x(2))^2;
2419	// return(Phi(f));
2420	return(fr);
2421	}
2422
2423	///////////////////////////////////////////////////////////////////////////////
2424	proc debug_log (int level, list #)
2425	"USAGE: debug_log(level,li); level=int, li=comma separated \"message\" list
2426	COMPUTE: print \"messages\" if level>=@DeBug.
2427	useful for user-defined trace messages.
2428	EXAMPLE: example debug_log; shows an example
2429	SEE ALSO: init_debug
2430	"
2431	{
2432	int len = size(#);
2433	// int printresult = printlevel - level +1;
2434	// if(level>1) {
2435	// dbprint(printresult, "Debug:("+ string(level)+ "): ", #[2..len]);
2436	// }
2437	// else { dbprint(printresult, #[1..len]); }
2438	if( defined(@DeBug) == 0 ) { init_debug(); }
2439	if(@DeBug>=level) {
2440	if(level>1) { "Debug:("+ string(level)+ "): ", #[1..len]; }
2441	else { #[1..len]; }
2442	}
2443	}
2444	example
2445	{ "EXAMPLE:"; echo=2;
2446	example init_debug;
2447	}
2448
2449	///////////////////////////////////////////////////////////////////////////////
2450	proc init_debug(list #)
2451	"USAGE: init_debug([level]); level=int
2452	COMPUTE: Set the global variable @DeBug to level. The variable @DeBug is
2453	used by the function debug_log(level, list of strings) to know
2454	when to print the list of strings. init_debug() reports only
2455	changes of @DeBug.
2456	NOTE: The procedure init_debug(n); is usefull as trace-mode. n may
2457	range from 0 to 10, higher values of n give more information.
2458	EXAMPLE: example init_debug; shows an example"
2459	{
2460	int newDebug=0;
2461	if( defined(@DeBug) != 0 ) { newDebug = @DeBug; }
2462
2463	if( size(#) > 0 ) {
2464	newDebug=#[1];
2465	}
2466	else {
2467	string s=system("getenv", "SG_DEBUG");
2468	if( s != "" && defined(@DeBug)==0) {
2469	s="newDebug="+s;
2470	execute(s);
2471	}
2472	}
2473	if( defined(@DeBug) == 0) {
2474	int @DeBug = newDebug;
2475	export @DeBug;
2476	if(@DeBug>0) { "Debugging level is set to ", @DeBug; }
2477	}
2478	else {
2479	if( (size(#) == 0) && (newDebug < @DeBug) ) { return(); }
2480	if( @DeBug != newDebug) {
2481	int oldDebug = @DeBug;
2482	@DeBug = newDebug;
2483	if(@DeBug>0) { "Debugging level change from ", oldDebug, " to ",@DeBug; }
2484	else {
2485	if( @DeBug==0 && oldDebug>0 ) { "Debugging switched off."; }
2486	}
2487	}
2488	}
2489	printlevel = @DeBug;
2490	}
2491	example
2492	{ "EXAMPLE:"; echo=2;
2493	init_debug();
2494	debug_log(1,"no trace information printed");
2495	init_debug(1);
2496	debug_log(1,"some trace information");
2497	init_debug(2);
2498	debug_log(2,"nice for debugging scripts");
2499	init_debug(0);
2500	}
2501
2502	///////////////////////////////////////////////////////////////////////////////
2503	proc basicinvariants(poly f)
2504	"USAGE: basicinvariants(f); f = poly
2505	COMPUTE: Compute basic invariants of f: an upper bound d for the
2506	determinacy, the milnor number mu and the corank c of f
2507	RETURN: intvec: d, mu, c
2508	EXAMPLE: example basicinvariants; shows an example"
2509	{
2510	intvec v;
2511	ideal Jfs = std(jacob(f));
2512	v[1] = system("HC")+1;
2513	v[2] = vdim(Jfs);
2514	v[3] = corank(f);
2515	if( v[2]<v[1] ) { v[1] = v[2]+1; }
2516	return(v);
2517	}
2518	example
2519	{ "EXAMPLE:"; echo=2;
2520	ring r=0,(x,y,z),ds;
2521	basicinvariants((x2+3y-2z)^2+xyz-(x-y3+x2*z3)^3);
2522	}
2523
2524	///////////////////////////////////////////////////////////////////////////////
2525	proc corank(poly f)
2526	"USAGE: corank(f); f=poly
2527	RETURN: the corank of the Hessian matrix of f, of type int
2528	REMARK: corank(f) is the number of variables occuring in the residual
2529	singularity after applying 'morsesplit' to f
2530	EXAMPLE: example corank; shows an example"
2531	{
2532	matrix M = jacob(jacob(jet(f,2)));
2533	list lba = bareiss(M);
2534	int cr = nvars(basering) - size(module(lba[1]));
2535	return(cr);
2536	}
2537	example
2538	{ "EXAMPLE:"; echo=2;
2539	ring r=0,(x,y,z),ds;
2540	poly f=(x2+3y-2z)^2+xyz-(x-y3+x2*z3)^3;
2541	corank(f);
2542	}
2543	///////////////////////////////////////////////////////////////////////////////
2544	static proc Faktorisiere(poly f, poly fk, int pt, int k, intvec RFlg)
2545	{
2546	//---------------------------- initialisation ---------------------------------
2547	poly a, b, Relation;
2548	ideal B, Jfsyz;
2549	map PhiG, VERT;
2550	matrix Mat;
2551	list v;
2552	def ring_top=basering;
2553
2554	// Ziel: bestimme a,b sodass fk = (ax+by^k)^pt gilt.
2555	B = maxideal(1);
2556	PhiG = ring_top,B;
2557	debug_log(2, "Faktor: f=",Show(f)," Jet=",Show(fk)," k=",k," exp=",pt);
2558
2559	//----------------------- compute role of x and y -----------------------------
2560	Jfsyz = fk, diff(fk, x(1));
2561	Mat = matrix(syz(Jfsyz));
2562	if( (fk-subst(fk,x(1),0)) != 0 && (fk-subst(fk,x(2),0)) != 0 ) {
2563	// Wenn k>0 ist die Wahl fuer x & y bereits getroffen
2564	// sonst bestimmen x und y
2565	Jfsyz = fk, diff(fk, x(1));
2566	Mat = matrix(syz(Jfsyz));
2567	Relation = -pt * Mat[2,1] / Mat[1,1];
2568	a = Coeff(Relation, x(1), x(1));
2569	b = Coeff(Relation, x(2), x(2)^k);
2570	B = maxideal(1);
2571	if( (RFlg[1]==1 && k==1) \|\| k>1) { B[rvar(x(1))] = x(1)-b*x(2)^k; }
2572	else { B[rvar(x(2))] = x(2)-b*x(1)^k; }
2573	VERT = basering,B;
2574	f = VERT(f);
2575	PhiG = VERT(PhiG);
2576	}
2577
2578	//------------------- permutation of x and y, if needed -----------------------
2579	if( k==1 ) {
2580	debug_log(2, "Fak-7:",Show(f)," jet=",Show(fk));
2581	if(Coeff(jet(f, pt), x(1), x(1)^pt) == 0) {
2582	VERT = basering,x(2),x(1);
2583	f = VERT(f);
2584	PhiG = VERT(PhiG);
2585	}
2586	}
2587	debug_log(2, "Fak-8:",Show(f)," jet=",Show(fk));
2588	debug_log(6, "Faktorisiere liefert: f=", Show(f));
2589	v[1] = f;
2590	v[2] = PhiG;
2591	return(v);
2592	}
2593
2594	///////////////////////////////////////////////////////////////////////////////
2595	static proc Teile(poly f, poly fk)
2596	{
2597	ideal Jfsyz = f, fk;
2598	poly Relation;
2599	matrix Mat = matrix(syz(Jfsyz));
2600	Relation = -1 * Mat[2,1]/Mat[1,1];
2601	return(Relation);
2602	}
2603
2604	///////////////////////////////////////////////////////////////////////////////
2605	static proc GetRf (poly fi, int n)
2606	"USAGE: GetRf();"
2607	{
2608	//---------------------------- initialisation ---------------------------------
2609	int j, k, l1, l1w;
2610	matrix Koef;
2611	intvec RFlg;
2612
2613	RFlg[n] = 0;
2614	intvec Haeufigkeit = RFlg;
2615
2616	for( j=1; j<=n ; j=j+1) {
2617	Koef = coef(fi, x(j));
2618	Haeufigkeit[j] = ncols(Koef);
2619	if(Coeff(fi, x(j),0) == 0) { Haeufigkeit[j] = Haeufigkeit[j] + 1;}
2620	}
2621	for( j=n; j>0 ; j=j-1) {
2622	l1 = 0;
2623	l1w = 0;
2624	for(k=1;k<=n;k=k+1) { if(Haeufigkeit[k]>l1w) { l1=k; l1w=Haeufigkeit[k];}}
2625	RFlg[j] = l1;
2626	Haeufigkeit[l1] = 0;
2627	}
2628	debug_log(2, "Permutations:", RFlg);
2629	return(RFlg);
2630	}
2631
2632	///////////////////////////////////////////////////////////////////////////////
2633	static proc Show(poly g)
2634	{
2635	string s;
2636	def ring_save=basering;
2637
2638	execute(@ringdisplay);
2639	map showpoly=ring_save,maxideal(1);
2640	s = string(showpoly(g));
2641	setring ring_save;
2642	return (s);
2643	}
2644
2645	///////////////////////////////////////////////////////////////////////////////
2646	static proc checkring
2647	{
2648	int CH = char(basering);
2649	if(CH >= 2 && CH<=13) {
2650	"Ring has characteristic ",CH;
2651	"Characteristic other than 0 or 0<char<13 is not yet implemented";
2652	return(1);
2653	}
2654	return(0); // characteristic of ring is OK, return (0)
2655	}
2656
2657	///////////////////////////////////////////////////////////////////////////////
2658	static proc DecodeNormalFormString (string S_in)
2659	"USAGE: DecodeNormalFormString"
2660	{
2661	//---------------------------- initialisation ---------------------------------
2662	int C_eq, a, b, i, t, k, r, s;
2663	string s1, s2, s3, s4, s_in, Typ;
2664	list v = "Error", 0, 0, 0;
2665
2666	C_eq = find(S_in, "=")+1;
2667	s_in = S_in[C_eq,30];
2668	debug_log(2, "Decode:");
2669
2670	debug_log(2, "S_in=", S_in," s_in=",s_in );
2671	a = find(s_in, "[")+1;
2672	b = find(s_in, "]")-1;
2673	t = 1;
2674	k = 0;
2675	r = 0;
2676	s = 0;
2677
2678	if(a<0 \|\| b<0) { return(v); }
2679	Typ = s_in[1..a-1];
2680	s1 = s_in[a..b];
2681	debug_log(6, "Suche Type:", Typ);
2682	//---------------------- decode between brackets ----------------------------
2683	if( find(s1, ",") == 0) {
2684	debug_log(8, " Number of columns: 0");
2685	s2 = "k = "+s1+";";
2686	execute(s2);
2687	if( (Typ=="A[") \|\| (Typ=="D[") ) { s3 = "k"; }
2688	if( Typ == "E[") { t = 6; }
2689	if( Typ == "W[") { t = 12; }
2690	if( Typ == "Q[") { t = 6; }
2691	if( Typ == "Z[") { t = 6; }
2692	if( Typ == "U[") { t = 12; }
2693	if( t > 1 ) {
2694	i = k;
2695	k = k/t;
2696	b = i - t*k;
2697	if( (s1 == "Q[") && (b==0) ) { k=k-1; b=6; }
2698	if(Typ == "Z[") {
2699	if(b==0) { k=k-1; b=6; }
2700	if(b==1) { k=k-1; b=7; }
2701	}
2702	if( b == 0 ) { s3 = string(t)+"k"; }
2703	else { s3 = string(t)+"k+"+string(b); }
2704	}
2705	if( Typ == "S[") {
2706	i = k+1;
2707	k = i/12;
2708	b = i - 12*k;
2709	if( b == 1 ) { s3 = "k"; }
2710	else {
2711	if(b==0) { s3 = "12k"+string(b-1); }
2712	else { s3 = "12k+"+string(b-1); }
2713	}
2714	}
2715	s2 = Typ + s3 +"]";
2716	} // es kommt mindestens ein komma vor...
2717	//----------------------- more than 1 paramater -----------------------------
2718	else {
2719	b = find(s1, ",");
2720	s2 = "k = ",s1[1..b-1],";";
2721	execute(s2);
2722	s1 = s1[b+1..size(s1)];
2723	if(find(s1, ",") == 0) {
2724	debug_log(8, " Number of columns 1");
2725	s2 = "r = "+s1+";";
2726	execute(s2);
2727	s4 = "r";
2728	s3 = "k";
2729	if(r==0) { s4 = string(0); }
2730	if(k==0 && Typ=="Z[") { s3 = string(1); }
2731	if(Typ[2] == "#") {
2732	i = r+1;
2733	r = i/2;
2734	b = i - 2*r;
2735	if( b == 1 ) { s4 = "2r"; }
2736	else { s4 = "2r-1"; }
2737	}
2738	s2 = Typ + s3 + "," + s4 +"]";
2739	} // es kommt mindestens zwei komma vor...
2740	//----------------------- get third parameter -----------------------------
2741	else {
2742	debug_log(8, " Number of columns >=2");
2743	debug_log(2, "Y[k,r,s] / Z[k,r,s] / T[k,r,s]");
2744	b = find(s1, ",");
2745	s2 = "r = ",s1[1..b-1],";";
2746	execute(s2);
2747	s2 = "s = ",s1[b+1..size(s1)],";";
2748	execute(s2);
2749	if(Typ=="Y[") { s2 = "Y[k,r,s]"; }
2750	if(Typ=="Z[") { s2 = "Z[k,r,s]"; }
2751	if(Typ=="T[") { s2 = "T[k,r,s]"; }
2752	}
2753	}
2754	debug_log(2, "Looking for Normalform of ",s2,"with (k,r,s) = (",
2755	k,",",r,",", s, ")" );
2756	v = s2, k, r, s;
2757	return(v);
2758	}
2759
2760	///////////////////////////////////////////////////////////////////////////////
2761	proc A_L
2762	"USAGE: A_L(f); f poly
2763	A_L(s); s string, the name of the singularity
2764	COMPUTE: the normal form of f in Arnold's list of singularities in case 1,
2765	in case 2 nothing has to be computed.
2766	RETURN: A_L(f): compute via 'milnorcode' the class of f and return the normal
2767	form of f found in the database.
2768	A_L(\"name\"): get the normal form from the database for the
2769	singularity given by its name.
2770	EXAMPLE: example A_L; shows an example"
2771	{
2772	// if trace/debug mode not set, do it!
2773	init_debug();
2774
2775	if( typeof(#[1]) == "string" ) {
2776	if(checkring()) { return(#[1]); }
2777	return(normalform(#[1]));
2778	}
2779	if( typeof(#[1]) == "poly" ) {
2780	if(checkring()) { return(#[1]); }
2781	return(quickclass(#[1]));
2782	}
2783
2784	}
2785	example
2786	{ "EXAMPLE:"; echo=2;
2787	ring r=0,(a,b,c),ds;
2788	poly f=A_L("E[13]");
2789	f;
2790	A_L(f);
2791	}
2792
2793	///////////////////////////////////////////////////////////////////////////////
2794	proc normalform(string s_in)
2795	"USAGE: normalform(s); s=string
2796	RETURN: Arnold's normal form of singularity with name s
2797	EXAMPLE: example normalform; shows an example."
2798	{
2799	string Typ;
2800	int k, r, s, crk;
2801	int n, i;
2802	poly f;
2803	list v;
2804	def ring_ext = basering;
2805	n = nvars(basering);
2806	ring ring_top=char(basering),(x(1..n)),(c,ds);
2807
2808	if(checkring()) { return(s_in); }
2809	if(nvars(basering)<=1) {
2810	"We need at least 2 variables in basering, you have",nvars(basering),".";
2811	return();
2812	}
2813	// if trace/debug mode not set, do it!
2814	init_debug();
2815
2816	v = DecodeNormalFormString(s_in);
2817	Typ, k, r, s = v[1..4];
2818	if(Typ=="Error") { return(0); }
2819	v = singularity(Typ, k, r, s);
2820	poly f_out;
2821	f_out, crk = v[1..2];
2822	if(crk>1) { for(i=crk+1;i<=n;i=i+1) { f_out = f_out + x(i)^2; } }
2823	setring ring_ext;
2824	map conv_top2ext=ring_top,maxideal(1);
2825	f = conv_top2ext(f_out);
2826	// f, crk = v[1..2];
2827	return(f);
2828	}
2829	example
2830	{ "EXAMPLE:"; echo=2;
2831	ring r=0,(a,b,c),ds;
2832	normalform("E[13]");
2833	}
2834
2835	///////////////////////////////////////////////////////////////////////////////
2836	proc swap
2837	"USAGE: swap(a,b);
2838	RETURN: b,a if a,b is the input (any type)"
2839	{
2840	return(#[2],#[1]);
2841	}
2842	example
2843	{ "EXAMPLE:"; echo=2;
2844	swap("variable1","variable2");
2845	}
2846
2847	///////////////////////////////////////////////////////////////////////////////
2848	proc Setring(int c, string name)
2849	"USAGE: "
2850	{
2851	string s="ring "+name+"=0,(x(1.."+ string(c) +")),(c,ds);";
2852	return(s);
2853	}
2854
2855	///////////////////////////////////////////////////////////////////////////////
2856	proc internalfunctions()
2857	"USAGE: internalfunctions();
2858	RETURN: nothing, display names of internal procedures of classify.lib
2859	EXAMPLE: no example"
2860	{ " Internal functions for the classification using Arnold's method,";
2861	" the function numbers correspond to numbers in Arnold's classifier:";
2862	"Klassifiziere(poly f); //determine the type of the singularity f
2863	Funktion1bis (poly f, list cstn)
2864	Funktion3 (poly f, list cstn)
2865	Funktion6 (poly f, list cstn)
2866	Funktion13 (poly f, list cstn)
2867	Funktion17 (poly f, list cstn)
2868	Funktion25 (poly f, list cstn)
2869	Funktion40 (poly f, list cstn, int k)
2870	Funktion47 (poly f, list cstn)
2871	Funktion50 (poly f, list cstn)
2872	Funktion58 (poly fin, list cstn)
2873	Funktion59 (poly f, list cstn)
2874	Funktion66 (poly f, list cstn)
2875	Funktion82 (poly f, list cstn)
2876	Funktion83 (poly f, list cstn)
2877	Funktion91 (poly f, list cstn, int k)
2878	Funktion92 (poly f, list cstn, int k)
2879	Funktion93 (poly f, list cstn, int k)
2880	Funktion94 (poly f, list cstn, int k)
2881	Funktion95 (poly f, list cstn, int k)
2882	Funktion96 (poly f, list cstn, int k)
2883	Funktion97 (poly f, list cstn)
2884	Isomorphie_s82_x (poly f, poly fk, int k)
2885	Isomorphie_s82_z (poly f, poly fk, int k)
2886	Isomorphie_s17 (poly f, poly fk, int k, int ct)
2887	printresult (string f,string typ,int Mu,int m,int corank,int K)
2888	";
2889	" Internal functions for the classifcation by invariants:
2890	Cubic (poly f)
2891	parity (int e) //return the parity of e
2892	HKclass (intvec i)
2893	HKclass3( intvec i, string SG_Typ, int cnt)
2894	HKclass3_teil_1 (intvec i, string SG_Typ, int cnt)
2895	HKclass5 (intvec i, string SG_Typ, int cnt)
2896	HKclass5_teil_1 (intvec i, string SG_Typ, int cnt)
2897	HKclass5_teil_2 (intvec i, string SG_Typ, int cnt)
2898	HKclass7 (intvec i, string SG_Typ, int cnt)
2899	HKclass7_teil_1 (intvec i, string SG_Typ, int cnt)
2900	";
2901	" Internal functions for the Morse-splitting lemma:
2902	Morse(poly fi, int K, int corank) //splitting lemma itself
2903	Coeffs (list #)
2904	Coeff
2905	";
2906	" Internal functions providing tools:
2907	ReOrder(poly f)
2908	Singularitaet(string typ,int k,int r,int s,poly a,poly b,poly c,poly d)
2909	RandomPolyK
2910	Faktorisiere(poly f, poly g, int p, int k) compute g = (ax+by^k)^p
2911	Teile(poly f, poly g); //divides f by g
2912	GetRf(poly f, int n);
2913	Show(poly f);
2914	checkring();
2915	DecodeNormalFormString(string s);
2916	Setring(int n, string ringname);
2917	";
2918	}
2919	example
2920	{
2921	"EXAMPLE"; echo=2;
2922	internalfunctions();
2923	}
2924
2925	///////////////////////////////////////////////////////////////////////////////
2926	// E n d O f F i l e
2927

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

Download in other formats:

Original Format