"Fossies" - the Fresh Open Source Software Archive 
Member "KASH3-lib-archindep-2008-07-31/lib/unit_group_res.g" (3 Sep 2008, 23753 Bytes) of package /linux/misc/old/KASH3-lib-archindep-2008-07-31.tar.gz:
As a special service "Fossies" has tried to format the requested text file into HTML format (style:
standard) with prefixed line numbers.
Alternatively you can here
view or
download the uninterpreted source code file.
1 ############################################
2 ############################################
3 # Gesucht ist das Primideal p, das entsteht, wenn ein Primideal P aus der
4 # Maximalordnung mit o geschnitten wird:
5
6 intersection := function(ORDER, P)
7 local Bo, BP, n, M, N, H, i, Bp;
8
9 n := ORDER.n;
10 Bo := TransformationMatrix(ORDER.o,ORDER.O);
11 BP := BasisMatrix(P);
12
13 # 1.) Setze die Basismatrizen von o und P zu M uebereinander:
14 M := VerticalJoin(Matrix(ORDER.Ri,Bo),Matrix(ORDER.Ri,BP));
15 # 2.) Bestimme eine Matrix mit maximalen Rang N, so dass N*M = 0:
16 N := Matrix(ORDER.Ri,NullspaceMatrix(M));
17
18 # 3.) Streiche die letzten n Spalten von N weg:
19 for i in [1..n] do
20 N := RemoveColumn(N, n+1);
21 od;
22
23 # 4.) Bp ist die Basismatrix von o geschnitten P:
24 return Ideal(ORDER.o, N);
25 end;
26
27 ##########################################################################
28
29 # Fuer Testbeispiele in kash3, fuer die allgemein der Schnitt zweier
30 # Ideale ueber der selben Ordnung benoetigt wird:
31
32 InstallMethod(
33 rec(name := "Intersection",
34 kind := "FUNCTION",
35 sin := [[elt-ids^int/ord^num,"A"],[elt-ids^int/ord^num,"B"]],
36 sou := [[elt-ids^int/ord^num,"C"]],
37 short:= "The intersection 'C' of the ideals 'A' and 'B'. Both ideals must be in the same order.",
38 author := ["Carolin Just"],
39 ex := ["o := EquationOrder(X^2+5);\nA :=256*o;\nB := 6*o;\nIntersection(A,B);"]
40 ),
41 function(a,b)
42 local A, B, n, o, M, N, i, AB;
43
44 A := BasisMatrix(a);
45 B := BasisMatrix(b);
46 n := NumberOfRows(A);
47 o := Order(a);
48
49 if o <> Order(b) then
50 Error("Intersection: Both ideals must be in the same order");
51 fi;
52 # 1.) Setze die Basismatrizen von a und b zu M uebereinander:
53 M := VerticalJoin(Matrix(Z,A),Matrix(Z,B));
54 # 2.) Bestimme eine Matrix mit maximalen Rang N, so dass N*M = 0:
55 N := Matrix(Z,NullspaceMatrix(M));
56
57 # 3.) Streiche die letzten n Spalten von N weg:
58 for i in [1..n] do
59 N := RemoveColumn(N, n+1);
60 od;
61
62 # 4.) AB ist die Basismatrix von a geschnitten b:
63 AB := N*A;
64 return Ideal(o, AB);
65 end
66 );
67
68 ##########################################################################
69
70 # Gesucht ist das Ideal a, das entsteht, wenn ein Ideal A aus der
71 # mit einer Unterordnung o geschnitten wird:
72
73
74 InstallMethod(
75 rec(name := "Intersection",
76 kind := "FUNCTION",
77 sin := [[ord^num,"o"],[elt-ids^int/ord^num,"A"]],
78 sou := [[elt-ids^int/ord^num,"a"]],
79 short:= "The intersection 'a' of the ideal 'A' and the order 'o'. The order 'o' must be a suborder of the order of 'A'",
80 author := ["Carolin Just"],
81 ex := ["o := EquationOrder(X^2-5);\nO := MaximalOrder(o);\nA :=256*O;\na := Intersection(o,A);"]
82 ),
83 function(o, A)
84 local Ri, Bo, BA, n, M, N, H, i;
85
86 if not FieldOfFractions(o)=FieldOfFractions(Order(A)) then
87 Error("Intersection: the order of A and the order o are not compatible");
88 fi;
89
90 n := Degree(o);
91 Bo := TransformationMatrix(o,Order(A));
92
93 if not IsUnit(Bo.ext1) then
94 Error("Intersection: the order o must be a suborder of the order of A");
95 fi;
96
97 BA := BasisMatrix(A);
98 Ri := CoefficientRing(o);
99
100
101 # 1.) Setze die Basismatrizen von o und A zu M uebereinander:
102 M := VerticalJoin(Matrix(Ri,Bo),Matrix(Ri,BA));
103 # 2.) Bestimme eine Matrix mit maximalen Rang N, so dass N*M = 0:
104 N := Matrix(Ri,NullspaceMatrix(M));
105
106 # 3.) Streiche die letzten n Spalten von N weg:
107 for i in [1..n] do
108 N := RemoveColumn(N, n+1);
109 od;
110
111 # 4.) Ba ist die Basismatrix von o geschnitten A:
112 return Ideal(o, N);
113 end
114 );
115
116
117 ##########################################################################
118 ##########################################################################
119 # Gesucht sind die Primideale p, die a enthalten.
120
121 prime_ideals_of_a := function(ORDER, a)
122 local amax, prid_amax, prid_a, P, p;
123
124 # Berechne aO in der MaxOrd O:
125 amax := Ideal(ORDER.O, Generators(a));
126
127 # Zerlege aO in Primideale:
128 prid_amax := Factorization(amax);
129 prid_a := [];
130
131 # Speichere Primideale nicht doppelt ab:
132 for P in prid_amax do
133 p := intersection(ORDER, P[1]);
134 if not (p in prid_a) then
135 Add_(prid_a,p);
136 fi;
137 od;
138
139 # Gib die Liste der Primideale p zurueck:
140 return prid_a;
141 end;
142
143 ###############################################################################
144
145 # Bestimmung des Exponenten m mit p^mOp < aOp
146 # Test: Ist ein Ideal einer Ordnung in einer Primidealpotenz der selben
147 # Ordnung enthalten?
148
149 mSearch := function(o,a,p)
150 local ch_p, my, apmy, m, s, ende, pm, gm, pm_in_apmy, alpha;
151
152 ch_p := IsPrimePower(Index(o,p)).ext1;
153 my := Valuation(Index(o,a),ch_p).base;
154 apmy := a+p^my;
155 s := 0;
156 ende := Ilog(2,my)+1;
157
158 while s<=ende do
159 m := Minimum(2^s,my);
160 pm := p^m;
161 gm := Generators(pm);
162 pm_in_apmy := true;
163 for alpha in gm do
164 if not(alpha in apmy) then
165 pm_in_apmy := false;
166 fi;
167 od;
168 if pm_in_apmy then
169 s := ende+1;
170 else
171 s := s+1;
172 fi;
173 od;
174 return m;
175 end;
176
177 ################################################################################
178 ################################################################################
179
180
181 # Diese Funktion erzeugt ein record, das alle wichtigen Infos, die nicht doppelt
182 # berechnet werden sollen, enthaelt
183
184 infos_zu_p := function(ORDER,a,p,m)
185
186 local o, n, Ri, Frac, omega_vec, p_k, P_k, P_k_inv, P_invQ, tau_multi, Quot_k,
187 Relationen_k_l, ha, ende, s, k, pk, Pk, S, Q_lk, Q_lk_inv, Bas_vec, Bas_m,
188 i, Bas_apm, iota;
189
190 o := ORDER.o;
191 n := ORDER.n;
192 Ri := ORDER.Ri;
193 Frac := ORDER.Frac;
194
195 # Initialisierung der record-Elemente
196 omega_vec := Matrix(o,n,1,Basis(o));
197 p_k := []; # Liste der Ideale p^k
198 P_k := []; # Liste der TrafoMatrizen Pk von p^k
199 P_k_inv := []; # Liste der Inversen Pk^-1
200 P_invQ := []; # Liste der Matrizen Pk^-1*Q_lk. Diese werden fuer den
201 # diskreten Logarithmus in (1+p^k)/(1+p^l) gebraucht.
202 tau_multi := [];# Liste der Erzeuger (mult. Gruppe)
203 Quot_k := []; # Liste der Restklassenringe o/p^k
204 Relationen_k_l := []; # Liste der Relationenmatrizen von (1+p^k)/(1+p^l)
205 ha := []; # Liste der Erzeuger von (1+a+p^m)/(1+p^m)
206 ende := Ilog(2,m);
207
208 # Fuellen der record-Elemente
209 for s in [0..ende+1] do
210 k := Minimum(2^s,m);
211 # Das Ideal b^k, seine Trafomatrix und deren Inverses:
212 pk := p^k; # p^k = p^(2^s)
213 # Achtung: TrafoMatrizen werden bei Kash3 immer schon in HNF angegeben!
214 Pk := TransformationMatrix(pk);
215 p_k[s+1] := pk;
216 P_k[s+1] := Pk.base;
217 P_k_inv[s+1] := Matrix(Frac,Pk)^-1;
218
219 # Relationen von (1+b^k)/(1+b^l). Ausnahmsweise eingetragen an der Stelle fuer l (!!)
220 # Umgeaendert zu einer SNF
221 if (s<>0) then
222 S := SmithForm(Matrix(ORDER.Ri, Pk * P_k_inv[s]));
223 Relationen_k_l[s+1] := S.base;
224 Q_lk := S.ext2; # Zur Anpassung der Erzeugervektoren tau
225 Q_lk_inv := Matrix(Ri,Matrix(Frac,Q_lk)^(-1));
226 P_invQ[s] := P_k_inv[s]*Q_lk; # fuer diskreten Log.
227 # in (1+p^k)/(1+p^l)
228
229 # Neue Konstruktion der Basen:
230 # Die additiven Basisvektoren - an die SNF von Pl*Pk^(-1) angepasst
231 Bas_vec := Matrix(o,Q_lk_inv*P_k[s])*omega_vec;
232
233 # Die multiplikativen Erzeuger von (1+p^k)/(1+p^l). Als Vektor abgespeichert:
234 # Bas_m := Bas_vec;
235 for i in [1..n] do
236 # SetEntry(Bas_m,i,1, Bas_m[i][1]+1);
237 SetEntry(Bas_vec,i,1, Bas_vec[i][1]+1);
238 od;
239 # tau_multi[s] := Bas_m;
240 tau_multi[s] := Bas_vec;
241 fi;
242
243 # Die Restklassenringe o/p^k
244 Quot_k[s+1] := Quotient(o, pk);
245
246 # Erzeuger von (1+(a+p^m)^k)/(1+(a+p^m)^l)
247 # Brauche ich fuer die Berechnung der Matrix A in unitGroup
248 if k<>m then
249 Bas_apm := Basis((a+p^m)^k);
250 # Liste der Basiselemente von (a+p^m)^k
251 for iota in [1..n] do
252 ha[s*n+iota] := Bas_apm[iota]+1;
253 # = Basiselemente +1
254 od;
255 fi;
256 od;
257
258 return rec(p_k:=p_k,
259 P_k:=P_k,
260 P_k_inv:=P_k_inv,
261 P_invQ := P_invQ,
262 Relationen_k_l:=Relationen_k_l,
263 tau_multi:=tau_multi,
264 Quot_k:= Quot_k,
265 ha:=ha);
266 end;
267
268 #############################################################################################
269 #############################################################################################
270
271
272 ##################################################################################
273 # Der diskrete Logarithmus fuer (1+b)/(1+b^m)
274
275 # Es gibt 1 Hilfsfunktion:
276 # Die Fkt. fuer den diskreten Log. in (1+b^k)/(1+b^l) gibt einen Vektor c
277 # mit Koeffizienten fuer die Darstellung eines Elementes in (1+b^k)/(1+b^l)
278 # aus.
279 #
280 # Die Hauptfunktion gibt eine Sequenz a aus, die alle Koeffizienten enhaelt,
281 # so dass [eta] aus (1+b)/(1+b^m) darstellbar ist als [a*g], g sind die Erzeuger
282 # von (1+b)/(1+b^m)
283
284 ###################################################################################
285 # Hilfsfunktion:
286
287 # s zeigt auf k!!
288 # Berechnet wird der diskrete Logarithmus in (1+b^k)/(1+b^l):
289
290 discrete_logarithm_bk_bl := function(ORDER,IDEAL,s,eta)
291 local beta;
292
293 # 1.) Bestimme die Koeffizienten von eta-1 bzgl. omega. Sie werden fuer (2.) gebraucht:
294 beta := ElementToSequence(Coerce(ORDER.o,eta-1));
295
296 # 2.) Bestimme die Koffizienten von eta-1 bzgl. tau_k und gib sie aus:
297 return Matrix(Z,Matrix(ORDER.Frac, 1, ORDER.n, beta) * IDEAL.P_invQ[s]);
298 # Achtung, hier steht Z - zu beachten fuer R = F_q[t]!!!
299 end;
300
301 #######################################################################################
302
303 diskreter_logarithmus:= function(ORDER,IDEAL,eta,m)
304
305 local k,coeff,s,Quot_m,l,c,g;
306
307 k := 1;
308 coeff := [];
309 s := 1; # s zeigt auf die Stelle fuer b^k,
310 # s+1 zeigt auf die Stelle fuer b^l
311 Quot_m := IDEAL.Quot_k[Ilog(2,m)+2]; # = o/p^m
312 eta := Coerce(Quot_m,eta);
313
314 if k=m then # Fuer den Sonderfall m=1
315 coeff := ElementToSequence(discrete_logarithm_bk_bl(ORDER,IDEAL,s,eta));
316 fi;
317
318 while k<>m do
319 l := Minimum(2*k,m);
320 c := discrete_logarithm_bk_bl(ORDER,IDEAL,s,eta);
321 g := Transpose(Matrix(Quot_m, IDEAL.tau_multi[s]));
322 # die Eintraege der Erzeuger nach o/p^m schieben, damit PowerProduct
323 # schneller laeuft
324 eta := eta*(PowerProduct(g,c)^(-1));
325 coeff := Concatenation(coeff,ElementToSequence(c));
326 k := l;
327 s := s+1;
328 od;
329
330 return coeff;
331 end;
332
333
334 ################################################################################
335 ################################################################################
336
337
338 ##############################################################################
339 # Mg gehoeren zu (1+b)/(1+b^k).
340 # Nh gehoeren zu (1+b^k)/(1+b^l).
341 # s zeigt auf k!!
342
343 exact_sequence := function(ORDER, IDEAL, s, Mg, Nh)
344
345 local M, N, g, g_hoch_M_Seq, g_hoch_M, k, r, P_seq, i, P, c, Rel;
346
347 M := Mg[1];
348 r := NumberOfRows(M);
349 N := Nh[1];
350
351 # 1.) Den Vektor M*g berechnen: - Bis jetzt funktioniert PowerProduct nicht, wenn M eine Matrix ist,
352 # deswegen erstmal noch so umstaendlich...
353 g := Matrix(IDEAL.Quot_k[s+1],Transpose(Mg[2]));
354 g_hoch_M_Seq := [];
355 for k in [1..r] do
356 g_hoch_M_Seq[k] := PowerProduct(g,M[k]);
357 od;
358 g_hoch_M := Matrix(ORDER.o,r,1,g_hoch_M_Seq);
359
360 # 2.) Fuer jeden Eintrag von g_hoch_M den diskreten Logarithmus in (1+b^k)/(1+b^l)
361 # loesen und in einer Liste abspeichern:
362 P_seq := [];
363 for i in [1..r] do
364 P_seq[i] :=
365 discrete_logarithm_bk_bl(ORDER,IDEAL,s,g_hoch_M[i][1]);
366 od;
367 # 3.) Aus der Liste die Matrix P bauen:
368 P := VerticalJoin(P_seq);
369
370 # 4.) M und N diagonal anordnen, danach -P einfuegen:
371 c := NumberOfColumns(M);
372 Rel := InsertBlock(DiagonalJoin(M,N),Matrix(Z,-P),1,c+1);
373
374 # 5.) Ausgegeben werden die neue Relationenmatrix M und der neue Erzeugervektor g:
375 return [HermiteForm(Rel).base, VerticalJoin(Mg[2],Nh[2])];
376 end;
377
378 ##############################################################################
379 # s zeigt auf l!!!
380
381 computing_1_plus_p_nach_1_plus_p_m := function(ORDER, IDEAL,m)
382 local ende, Mg, k, s, l, Nh;
383
384 ende := Ilog(2,m);
385
386 # 1.) Berechne Erzeuger und Relationen zu (1+p)/(1+p^2):
387 Mg := [IDEAL.Relationen_k_l[2], IDEAL.tau_multi[1]];
388
389 # 2.) Von (1+p)/(1+p^k) und (1+p^k)/(1+p^l) nach (1+p)/(1+p^l).
390 k := 2;
391 s := 2;
392 while k<m do
393 l := Minimum(2*k,m);
394 Nh := [IDEAL.Relationen_k_l[s+1], IDEAL.tau_multi[s]];
395 Mg := exact_sequence(ORDER, IDEAL, s, Mg, Nh);
396 k := l;
397 s := s+1;
398 od;
399
400 return Mg;
401 end;
402
403 #############################################################################
404 #############################################################################
405 # Ein Testbeispiel:
406
407 # Hier ein Beispiel fuer eine Nicht-Maximalordnung und Nicht-Gleichungsordnung:
408 # (laeuft auch durch!!)
409
410 #oe := EquationOrder(X^2-5);
411 #o := Order(oe,Matrix(Z,2,2,[2,0,0,6]),2);
412 #n := 2;
413 #Frac := Q;
414 #Ri := Z;
415
416 #o := MaximalOrder(X^20 -5*X^19 +3*X^18 -7*X^17 +12*X^16 -3*X^15 +7*X^14 +17*X^13-6*X^12
417 #-9*X^11 +5*X^10 +11*X^9 +3*X^8 -10*X^7 -18*X^6
418 #+10*X^5 +17*X^4 +16*X^3 -14*X^2 +3*X -15);
419 #Frac := Q;
420 #Ri := Z;
421 #n := 20;
422 #b := Ideal(o,7,1+o.2);
423 #Gens := Generators(b);
424 #eta := 789*Gens[1] + 456233*Gens[2] +1;
425 #m := 15;
426
427 #o := MaximalOrder(X^14-15*X^13+14*X^12-5*X^11-17*X^10 +17*X^9 +4*X^8 +8*X^7 +16*X^6
428 #-17*X^5 +8*X^4 +12*X^3 -4*X^2 -16*X -10);
429 #Frac := Q;
430 #Ri := Z;
431 #n := 14;
432 #b := Ideal(o,5,o.2); # laut kash3 ebenfalls ganz o??
433
434 #b := Ideal(o,o.2);
435 #eta := 45629074*o.2 +1;
436 #m := 17;
437
438 #o := EquationOrder(X^6+11*X^5+20*X^4+3*X^3-18*X^2+7*X+11);
439 #O := MaximalOrder(X^6+11*X^5+20*X^4+3*X^3-18*X^2+7*X+11);
440 #Frac := Q;
441 #Ri := Z;
442 #n := 6;
443 #b := Ideal(o,47,[40,20,1,0,0,0]);
444 #Gens := Generators(b);
445 #eta := 789*Gens[1] + 456233*Gens[2] +1;
446 #m := 40;
447
448 # b := Ideal(o,3,[1,0,0,1,0,0]); Das ist bereits die ganze Ordnung!!
449
450 #o := MaximalOrder(X^3+23*X^2+88*X+51);
451 #Frac := Q;
452 #Ri := Z;
453 #n := 3;
454 #p := 7*o;
455 #b := Ideal(o,29,[18,1,0]);
456 #Gens := Generators(b);
457 #eta := 789*Gens[1] + 456233*Gens[2] +1;
458 #m := 31;
459
460 #o := MaximalOrder(X^4+14*X^3+8*X^2+4*X+11);
461 #Frac := Q;
462 #Ri := Z;
463 #n := 4;
464 #b := Ideal(o,2,[1,1,0,0]);
465 #Gens := Generators(b);
466 #eta := Element(o,[4,1,0,0]);
467 #eta := 789*Gens[1] + 456233*Gens[2] +1;
468
469 #o := EquationOrder(X^2-5);
470 #m := 5;
471 #b := 5*o;
472 #b3 := b^3;
473 #b5 := b^5;
474
475 #############################################################################
476 ################################################################################
477 # Das ist die Hauptfunktion, die alle bisher bereit gestellten Techniken
478 # zusammenfuehrt. Ihr Ablauf entspricht dem Ablauf des Pseudocodes auf S.9
479
480 #############################################################################
481 # Hilfsfunktion: Chinesischer Restsatz fuer unseren speziellen Fall
482 # Input: Ordnung o
483 # Liste der Ideale a+p^m,
484 # den i-ten Erzeugervektor g,
485 # den Zeiger i
486 # die Anzahl s der Primideale, die a enthalten
487 # Output: angepasster Erzeugervektor h
488
489 chin_rem := function(o, L_apm, g, i, s)
490 local b, j, apm, len_g, h, k;
491
492 # 1.) Multipliziere die (a+p_j^m_j), j <> i
493 b := 1;
494 for j in [1..i-1] do
495 b := b*L_apm[j];
496 od;
497 for j in [i+1..s] do
498 b := b*L_apm[j];
499 od;
500
501 # 2.) Loese fuer jeden Eintrag von g den Chinesischen Restsatz
502 apm := L_apm[i];
503 len_g := Length(g);
504 h := [];
505 for k in [1..len_g] do
506 Add_(h, ChineseRemainderTheorem(apm,b,Coerce(o,g[k][1]), Coerce(o,1)));
507 od;
508
509 # 3.) Gib h als Vektor zurueck
510 return Matrix(len_g,1,h);
511 end;
512
513 ###############################################################################
514 # Hauptfunktion:
515 # Input: Restklassenring o/a
516 # Output: Einheitengruppe (o/a)* einschliesslich der Abbildungen zwischen
517 # (o/a)* und (o/a)
518
519 unitGroup := function(oa)
520
521 local a, o, ORDER, omega, Ri, PrIds, s, m_List, i, P_IDEAL, Rel, Erz, L_apm,
522 o_by_p, p, m, K, q, UK, Kphi, UKphi, Zeta, vp, t, Xi_mod_pm, Xi,
523 RES_FIELD, Mg, M, g, len_ha, A_entry, j, A, MA, H, Erz_fuer_psi,
524 Rel_fuer_psi, S, Rel_fuer_phi, Q_smith, Erz_fuer_phi, zeiger, nrow_rel,
525 G, App, U, phi, psi, disc_log;
526
527 # 0.) Grundsaetzlich brauchen wir:
528 a := Modulus(oa);
529 o := Order(a);
530
531 # 1.) Falls a=o, gib Fehlermeldung aus:
532 if IsUnit(Determinant(BasisMatrix(a))) then
533 Error("ideal = order. This quotient ring has no unit group.");
534 fi;
535
536 # 2.) Konstruiere ein record mit den wichtigen Angaben zur Ordnung o:
537 ORDER := rec();
538 omega := Basis(o);
539 Ri := CoefficientRing(o);
540 ORDER.o := o;
541 ORDER.omega := omega;
542 ORDER.n := Length(omega);
543 ORDER.Ri := Ri;
544 ORDER.Frac := FieldOfFractions(Ri);
545 ORDER.O := MaximalOrder(o);
546
547 # 3.a) Bestimme alle Primideale p, die a enthalten
548 PrIds := prime_ideals_of_a(ORDER,a);
549 # 3.b) Bestimme die Anzahl s der Primideale
550 s := Length(PrIds);
551 # 3.c) Bestimme alle Exponenten m
552 m_List := [];
553 for i in [1..s] do
554 m_List[i] := mSearch(o,a,PrIds[i]);
555 od;
556 # 3.d) Sammle zu jedem Primideal alle noetigen Informationen
557 P_IDEAL:= [];
558 for i in [1..s] do
559 P_IDEAL[i] := infos_zu_p(ORDER,a,PrIds[i],m_List[i]);
560 od;
561
562 # 4) Berechne fuer jedes Primideal p die Einheitengruppe von (O_p/aO_p):
563 Rel := []; # Liste aller Relationen der O_p_i/aO_p_i
564 Erz := []; # Liste aller Erzeuger der O_p_i/aO_p_i
565 L_apm := []; # Liste aller a+p_i^(m_i)
566 o_by_p := []; # Liste der Abbildungen fuer die (o/p_i)*
567
568 for i in [1..s] do
569 p := PrIds[i]; # Primideal p
570 m := m_List[i]; # Exponent m zu p
571 L_apm[i] := a+p^m;
572
573 # 4.1.) Berechne die multiplikative Gruppe von (O/p):
574 # Koerper K := (o/p) und dessen Groesse q.
575 # Dabei gibt (q-1) Relation des Erzeugers von (o/p)* an
576 K := ResidueClassField(p);
577 q := Size(Base(K));
578
579 # Einheitengruppe UK := (o/p)*
580 UK := UnitGroup(K);
581
582 # Fuer den diskreten Logarithmus in (o/p)*
583 Kphi := K.ext1;
584 UKphi := UK.ext1;
585
586 # Repraesentant von (o/p)* in o
587 Zeta := Preimage(UKphi(UK.1),Kphi);
588
589 # 4.1.a) Falls m=1, sind wir mit der Berechnung von (o/p)* fertig
590 if m=1 then
591 Erz[i] := Matrix([[Zeta]]);
592 # = Erzeugervektor
593 Rel[i] := Matrix([[q-1]]);
594 # = Relationenmatrix
595 RES_FIELD := rec();
596 RES_FIELD.Kphi := Kphi;
597 RES_FIELD.UKphi := UKphi;
598 o_by_p[i] := RES_FIELD;
599 # = Record, das die Abbildungen zwischen (o/p)* und o enthaelt.
600
601 # 4.1.b) Fuer alle anderen m berechnen wir (1+p)/(1+a+p^m)
602 else
603 # Dafuer Hensellifting von Zeta:
604 vp := Valuation(Coerce(o,q),p);
605 t := 0;
606 while (Minimum(q^t,1+t*vp)<m) do
607 t := t+1;
608 od;
609 Xi_mod_pm := Coerce(Quotient(o,p^m),Zeta)^(q^t);
610
611 # Der Repraesentant von Xi_mod_pm als Matrix
612 Xi := Matrix([[Preimage(Xi_mod_pm,Quotient(o,p^m).ext1)]]);
613
614 # Die Abbildungen fuer (o/p)* in einem record sammeln:
615 RES_FIELD := rec();
616 RES_FIELD.Kphi := Kphi;
617 RES_FIELD.UKphi := UKphi;
618 o_by_p[i] := RES_FIELD;
619
620 # 4.2.) Berechne (1+p)/(1+p^m):
621 Mg := computing_1_plus_p_nach_1_plus_p_m(ORDER,P_IDEAL[i],m);
622 M := Mg[1];
623 g := Mg[2];
624
625 # 4.3.) Berechne (1+a+p^m)/(1+p^m):
626 # a) Erzeuger fuer (1+(a+p^m))/(1+(a+p^m)^2), ..., (1+(a+p^m)^(2^ende))/(1+(a+p^m)^m)
627 # sind unter IDEAL.ha bereits abgespeichert.
628 # b) Berechne fuer jeden einzelnen Erzeuger den diskreten Logarithmus bzgl. der
629 # Erzeuger g von (1+p)/(1+p^m)
630 # c) Setze die erhaltenen Sequenzen zur Matrix A zusammen:
631
632 len_ha := Length(P_IDEAL[i].ha);
633 if (len_ha>0) then
634 A_entry := [];
635 for j in [1..len_ha] do
636 A_entry[j] := Matrix(Ri, 1, len_ha,
637 diskreter_logarithmus(ORDER,P_IDEAL[i],P_IDEAL[i].ha[j],m));
638 od;
639 A := VerticalJoin(A_entry);
640
641 # 4.4.) Setze alle Relationenmatrizen zusammen:
642 # Setze M und A uebereinander =: MA, gib MA in HNF zurueck:
643 MA := HermiteForm(VerticalJoin(M,A)).base;
644
645 # Die untere Haelfte von MA ist Nullmatrix, sie wird abgeschnitten:
646 MA := ExtractBlock(MA, 1, 1, len_ha, len_ha);
647
648 else
649 MA := M;
650 fi;
651
652 # Setze (q-1) und MA zu Diagonalmatrix zusammen:
653 Rel[i] := DiagonalJoin(Matrix([[q-1]]), MA);
654
655 # 4.5.) Setze die Erzeuger Xi von (o/p)* und g zusammen:
656 Erz[i] := VerticalJoin(Xi,g);
657 fi;
658
659 od;
660
661 ## Hier muessen vier Listen herauskommen:
662 # a) Erz: Liste aller Erzeuger, also zu jedem p_i
663 # b) Rel: Liste aller Relationen
664 # c) L_apm : Liste der a+p^m
665 # d) o_by_p: Liste mit Abbildungen zu (o/p_i)*
666
667
668 # 5) Jetzt noch den Chinesischen Restsatz anwenden:
669 # (Dafuer habe ich oben die Liste L_apm der Ideale a+p^m angelegt)
670
671 if s=1 then
672 H := Erz;
673 else
674 H := [];
675 for i in [1..s] do
676 H[i] := chin_rem(o,L_apm,Erz[i],i,s);
677 od;
678 fi;
679
680 # 6.) Zusammensetzen von Erzeugern und Relationen
681 # Fuer psi brauche ich genau den zusammengesetzten Vektor aus H
682 # Fuer phi brauche ich die Basis, die entsteht, wenn die SNF angewandt
683 # wird. - Nicht die Erzeuger
684
685 # 6.1) Daten zu psi
686 Erz_fuer_psi := VerticalJoin(H);
687 Rel_fuer_psi := DiagonalJoin(Rel);
688
689 # 6.2) Daten zu phi
690 S := SmithForm(Rel_fuer_psi);
691 Rel_fuer_phi := S.base;
692 Q_smith := Matrix(Z,Matrix(Q,S.ext2)^(-1));
693 Erz_fuer_phi := [];
694 for i in [1..NumberOfRows(Erz_fuer_psi)] do
695 Erz_fuer_phi[i] :=
696 PowerProduct(Transpose(Matrix(oa,Erz_fuer_psi)), Q_smith[i]);
697 od;
698 zeiger := 1;
699 while (Rel_fuer_phi[zeiger][zeiger] = 1) do
700 Remove_(Erz_fuer_phi,1);
701 zeiger := zeiger+1;
702 od;
703
704
705 # 7.) Konstruktion der Einheitengruppe U := (o/a)*
706
707 nrow_rel := NumberOfRows(Rel_fuer_phi);
708 G := FreeAbelianGroup(nrow_rel);
709 App := Apply([1..nrow_rel], x-> Coerce(G,List(Rel_fuer_phi[x])));
710 U := Quotient(G,App);
711
712
713 # 8.) Diskreter Logarithmus:
714
715 # Input: Ein Element u aus U := (o/a)*
716 # Output: Eine "Map" in kash. Das ist:
717 # eine Abbildung phi mit:
718 # Urbildbereich (o/a)*,
719 # Bildbereich (o/a),
720 # Umkehrabbildung psi, die in kash durch "Preimage" aufgerufen wird
721
722 # 8.1.) phi: Von (o/a)* nach o/a
723 # Input: Element u aus (o/a)*
724 # Output: Repraesentant in o/a (?)
725
726 # Ich brauche hierfuer den per SNF veraenderten Erzeugervektor Erz_fuer_phi,
727 # am besten bereits nach (o/a) geschoben und als Sequenz.
728
729 phi := function(u)
730 local c;
731
732 c := ElementToSequence(u);
733 return PowerProduct(Erz_fuer_phi,c);
734 end;
735
736 # 8.2.) psi: Von o/a nach (o/a)*
737 # Input: Repraesentant alpha aus o
738 # Output: Koeffizientenvektor von alpha bzgl. der Erzeuger von (o/a)*
739
740 # Ich brauche hierfuer den unveraenderten Erzeugervektor Erz_fuer_psi
741
742 psi := function(alpha)
743 local dl, i, a0, a1, dl_smith;
744
745 # I.) Wenn [alpha] kein Element der Einheitengruppe ist, wirf Fehlermeldung
746 if (not IsUnit(alpha)) then
747 Error("Element is no unit.");
748 fi;
749
750 # II.) Setze den diskreten Logarithmus aus den diskreten Logarithmen der opi/aopi
751 # zusammen:
752 dl := [];
753 for i in [1..s] do
754 # diskreter Logarithmus von alpha in (o/p)*
755 a0 := ElementToSequence(Preimage(o_by_p[i].Kphi(Coerce(o,alpha)),o_by_p[i].UKphi));
756 alpha := alpha*(Coerce(oa,H[i][1][1])^(-a0[1])); # so ungefaehr...
757 dl := Concatenation(dl,a0);
758 # diskreter Logarithmus von alpha in (1+p)/(1+a+p^m)
759 g := Remove(ElementToSequence(Matrix(oa,H[i])),1);
760 # (1.Eintrag von H[i] wegstreichen liefert Erzeugervektor von (1+p_i)/(1+p_i^m_i))
761 # (H[i] nach oa schieben verkuerzt Rechenaufwand beim Potenzieren)
762 if Length(g)>0 then
763 a1 := diskreter_logarithmus(ORDER,P_IDEAL[i],alpha,m_List[i]);
764 alpha := alpha*PowerProduct(g,a1);
765 dl := Concatenation(dl,a1);
766 fi;
767 od;
768
769 # III.) Passe dl an die Basis (Erz_fuer_psi) an:
770 dl_smith := ElementToSequence(Matrix(Z,1,Length(dl),dl)*S.ext2);
771 zeiger := 1;
772 while (Rel_fuer_phi[zeiger][zeiger] = 1) do
773 Remove_(dl_smith,1);
774 zeiger := zeiger+1;
775 od;
776
777 # IV.) Wandle dl_smith in kanonische Repraesentanten um und gib zurueck:
778 # return ElementToSequence(Coerce(U,dl_smith));
779 return Coerce(U,dl_smith);
780 end;
781
782 # 8.3.) disc_log gibt phi als Umkehrabbildung von psi aus.
783 disc_log := Map(U,oa,phi,psi);
784
785 # 9.) Rueckgabe eines records, das als Basis die Einheitengruppe U = (o/a)* und als
786 # Erweiterung die Map disc_log enthaelt
787 return rec(base := U.base, ext1 := disc_log);
788 end;
789
790
791 #InstallMethod(
792 #rec(name:="UnitGroup",
793 # kind:="FUNCTION",
794 # sin :=[[res^num,"oa"]],
795 # sou :=[[grp^abl,"U"],[map(grp^abl,res^num),"phi"]],
796 # short := "carolin weiss wie es geht",
797 # ex := [],
798 # see := [],
799 # author :=["Carolin Just"]
800 # ),
801 #unitGroup);
802
803 _UnitGroup := UnitGroup;
804 UnitGroup := function(arg)
805 if Is(Type(arg[1]),res^num) then
806 if not IsMaximal(Order(Modulus(arg[1]))) then
807 return unitGroup(arg[1]);
808 fi;
809 fi;
810 return _CallFunction(_UnitGroup,arg);
811 end;
812
813