Das Java Programm
In diesem Kapitel soll das Gaußsche Kreisteilungsverfahren als Java-Programm implementiert werden. Der vollständige Quellcode und die jar-Files können hier heruntergeladen werden. In diesem Kapitel werden nur Auszüge des Programms erläutert und zwar hauptsächlich die, die mit dem Gaußschen Verfahren zu tun haben.
Eigentlich handelt es sich um zwei Programme Hermes und Kompressor. Das Programm Kompressor wird hier aber nicht weiter behandelt, denn es komprimiert, wie der Name schon sagt, den Output des ersten Programms. Insbesondere werden alle Wurzelausdrücke gelöscht, die in den endgültigen Ausdruck nicht eingesetzt werden. Beim 257-Eck kann man auf diese Weise 89 von 129 Ausdrücken weglassen, beim 65537-Eck sogar 8009 von 10114 (die gar nicht erst berechneten Ausdrücke nicht mitgezählt). Ferner teilt das Programm die beim 65537-Eck doch sehr länglichen Wurzeln auf mehrere Zeilen auf.
Da man mit Hermes gleich LaTeX-Notation erzeugen kann, um sie unverändert in das Kap.8 bzw. in das voluminöse PDF zu übernehmen, wurde Kompressor nur für das LaTeX-Format angepasst.
Beide Programme sind Java-Konsol-Anwendungen, es gibt also keine aufwändige grafische Oberfläche. Beide Programme wurden bewusst spartanisch gehalten, damit das Wesentliche, nämlich der Algorithmus, nicht in einem Wust von Trivialitäten untergeht. So verzichten beide Programme sogar auf die Ausgabe in eine Datei und geben die Listen auf System.out aus (Kompressor kommt natürlich um eine Datei-Eingabe nicht herum). Es sollte daher tunlichst der Ausgabestrom auf eine Datei umgelenkt werden, wenn man die Berechnung für das 65537-Eck durchführen möchte.
Installation
Zur Installation kopiert man die Dateien Hermes.jar und Kompressor.jar einfach in ein beliebiges Verzeichnis. Zum Ablauf ist das Java® SE Runtime Environment (JRE) Version 7 (oder höher) erforderlich.
Benutzeroberfläche
Der Aufruf erfolgt unter der Eingabeaufforderung bzw. von einer Linux-Konsole aus mit
java -jar Hermes.jar [num] [opt]bei Auswertung für das 65537-Eck empfiehlt sich der Aufruf mit
java -jar -Xms512m -Xmx768m -XX:MaxPermSize=256m Hermes.jar [num] [opt]
sonst könnte der Platz in der Heap knapp werden. Dabei können für [num] die Werte 5, 17, 257, 65537 angegeben werden. Auf eine Berechnung für das Dreieck wurde verzichtet, weil dann an einigen Stellen unschöne Ausnahmebehandlungen hätten eingebaut werden müssen.
Das Argument [opt] dient zur Steuerung des Ausgabe-Layouts:
-tex : Ausgabe der Wurzelausdrücke im LaTeX-Format -ascii : Ausgabe der Wurzelausdrücke im ascii-Format -num : Ausgabe der Wurzelausdrücke mit eingesetzten numerischen Werten -val : Ausgabe des numerischen Werts der Periode -ref : Ausgabe des numerischen Werts der Periode plus Referenzwert und Abweichung -roots : Ausgabe der Exponenten der Einheitswurzeln in der Periode -sorted : Sortierte Ausgabe der Exponenten der Einheitswurzeln in der Periode
Die einzelnen Optionen können kombiniert werden, dann werden die unterschiedlichen Ausgabeformate gemischt. Die Ausgabe des Programms erfolgt, wie oben erwähnt, auf System.out.
Der Aufruf des Nachbereitungsprogramms Kompressor erfordert zwei Argumente:
java -jar Kompressor.jar [num] [file]
Für [num]
können nur die Werte 17,
257,
65537
angegeben werden – beim Fünfeck gibt es nichts zu komprimieren. Die Angabe wird gebraucht, um den Startpunkt für das Einsetzen
der Wurzelausdrücke zu finden. Beim 17-Eck ist dies die Periode $p_{3,0},$ beim
257-Eck die Periode $p_{7,0}$ und beim 65537-Eck die Periode $p_{15,0}.$ In einem ersten Durchlauf wird nämlich z.B. beim
65537-Eck zunächst ermittelt, welche Perioden kleineren Levels zur Darstellung der Periode $p_{15,0}$ benötigt werden,
anschließend werden rekursiv rückwärts die für die Darstellung der dabei gefundenen Perioden benötigten Perioden kleineren
Levels ermittelt, bis man schließlich bei der Periode $0$-ten Levels angekommen ist. Alle gefundenen Perioden werden als
notwendig
gekennzeichnet. In einem zweiten Durchlauf werden nur die so gekennzeichneten Perioden – diesmal
in der richtigen Reihenfolge – ausgegeben und alle nicht benötigten Perioden weggelassen.
Das Argument [file] übergibt an den Kompressor den Namen der Eingabedatei, die ihrerseits eine Ausgabe von Hermes enthalten muss, die mit der Option -tex erzeugt wurde. Wie bei Hermes erfolgt die Ausgabe auf System.out.
Hermes
im Schnelldurchlauf
In Prosa lässt sich die Arbeitsweise des Programms wie folgt zusammenfassen: Zunächst wird
die nullte Gaußsche Periode durch fortlaufendes Potenzieren des
primitiven Elements (im Programm grundsätzlich $3$) erzeugt (siehe fill).
Die Perioden der Stufe $1$ werden durch Aufspalten (siehe split) der nullten
Periode erzeugt und allgemein die Perioden der Stufe $n+1$ durch Aufspalten der
Perioden der $n$-ten Stufe.
Die beiden Tochterperioden einer aufgespalteten Periode werden miteinander
multipliziert – um die Addition muss man sich nicht kümmern, denn dabei entsteht, wie wir wissen,
nur die schon bekannte Elternperiode. Bei der Multiplikation (siehe mult) wird einfach, analog wie
wir es für das Siebzehneck gemacht haben, eine Multiplikationstabelle aufgebaut.
Dabei sind lediglich die Exponenten der primitiven Elemente zu addieren und gleich
$\pmod p$ zu reduzieren. Es wird hier gleich etwas Buchführung betrieben, indem
man für jeden der dabei entstehenden Exponenten einen Zähler versorgt, der die
Vielfachheit des jeweiligen Exponenten angibt. Das ist ganz ähnlich wie beim
Siebzehneck, wo wir ja ebenfalls in der Multiplikationstabelle gezählt haben, wie oft
jeder Exponent vorkommt.
Mit diesen Werten kommt man zum Kern des ganzen Algorithmus. Es wird im Elternlevel nach den Perioden gesucht , aus denen sich das Multiplikationsergebnis als Linearkombination zusammensetzen lässt. Findet man einen Exponenten (aus der Multiplikationstabelle) in einer Periode des Elternlevels, wird die gesamte gefundene Periode zu der Linearkombination hinzugenommen und alle Exponenten aus dieser Periode werden in der Multiplikationstafel in ihrer Vielfachheit um eins vermindert. Dass dies funktionieren muss , ohne dass negative Vielfachheiten auftreten oder Exponenten in der Multiplikationstafel übrig bleiben, haben wir im theoretischen Teil mit viel Mühe bewiesen. Natürlich wird das im Programm zusätzlich noch überprüft.
Nach Ende dieser Prozedur hat man die Koeffizienten für die Linearkombination gefunden, und es wird klar, warum wir soviel Wert darauf gelegt haben, dass diese Koeffizienten natürliche Zahlen sind. Andernfalls wäre nämlich der Ansatz, die Koeffizienten durch Zählen zu bestimmen, von vornherein zum Scheitern verurteilt.
Summe und Produkt der aktuell betrachteten Perioden sind nun als Linearkombination aus
Perioden des Vorgängerlevels dargestellt und können in dieser Form in die Formel für die
Lösung der quadratischen Gleichung eingesetzt werden, siehe compute. Anschließend kann man zum nächsten
Tochterlevel übergehen, indem man die aktuellen Perioden wieder aufspaltet. Dieses Verfahren
(Aufspalten - Multiplizieren - Koeffizienten suchen - Wurzelausdruck berechnen) wird
fortgesetzt, bis die Periodenlänge zwei erreicht wird.
Das Programm erzeugt eine fortlaufende Ausgabe von Wurzelausdrücken, wobei in jedem Ausruck nur Summen und Produkte von vorher ausgegebenen Wurzelausdrücken vorkommen. Da der Wert der nullten Gaußschen Periode immer $-1$ ist, hat man somit für jede Periode eine Darstellung aus Summen und Produkten rationaler Zahlen und ineinander geschachtelter Quadratwurzeln. Das ist aber nichts anderes als die Konstruktionsbeschreibung für das gerade betrachtete regelmäßige Vieleck.
Soweit der Überblick, nun zu den Einzelheiten:
Grundstruktur des Programms Hermes
Das Programm Hermes ist auf oberstem Level wie folgt gegliedert:
package hermes;
class Period {
...
}
class Hermes {
public static void main(String[] args) {
...
}
}
Die main-Methode der Klasse
Hermes wertet die
Argumente aus und baut eine verkettete Liste von Instanzen der
Klasse Period auf, die – natürlich – die Gaußschen
Perioden repräsentieren.
Die Verkettung erfolgt dabei einmal horizontal innerhalb aller
Perioden der gleichen Stufe und einmal vertikal, so dass man von
jeder Periode zu der Periode zurückgehen kann, aus der sie durch
Aufspaltung entstanden ist.
Die Klasse Period ist das Kernstück des ganzen
Algorithmus, daher wollen wir damit beginnen, sie zu analysieren.
Fangen wir mit den Datenstrukturen an:
class Period {
int len; // length of this period
int modulus; // it's either 5, 17, 257 or 65537
int level; // level of the period
int twoPowLevel; // 2^level
int seqnr; // sequent. number within level
int nr; // number of period within level
// all exponents are 3^k with
// k = nr (mod 2^level)
int primroot = 3; // primitive root 3
int[] exponents = null; // exponents of the roots of unity
int[] exponents_sorted = null; // the same sorted
Period leftSplit = null; // used for splitting
Period rightSplit = null; // the same
Period father = null; // vertical link
Period neighbour = null; // horizontal link
Period mostleft = null; // to 1st period in level
Period[] msummands = null; // representation as lin. comb. of periods
int[] mfactors = null; // with these factors
char sign; // sign of the squareroot
double theP; // used for solving the
double theQ; // quadratic equation
double numericValue; // value of period by sqareroots
double referenceValue; // value of period by cosinus
String name; // name of period p_{m,n}
Hier treffen wir viele Größen wieder, die bei der Behandlung der
Gaußschen Perioden eine Rolle gespielt haben. Die Länge len der
Periode, die bei der Aufteilung jeweils halbiert wird und in der
Stufe (dem level) $0$ den Wert $p-1$ bzw. $2^N$ hat.
Es folgen einige Größen, die bei der Nummerierung der Perioden, wie wir
sie im vorigen Kapitel erläutert haben, eine Rolle spielen.
Die wichtigste Größe ist aber das Array exponents, denn dies
ist ja praktisch die Gaußsche Periode. Wir hinterlegen dort
natürlich nur die Exponenten der Einheitswurzeln und zwar in
der Reihenfolge, wie sie durch fortlaufende
Potenzierung der primroot bzw. durch den Aufspaltungsprozess
entstanden sind. Das Array exponents_sorted wird zum einen für die
sortierte Ausgabe der Exponenten benötigt, zum anderen bei der Suche
nach einem bestimmten Exponenten in dieser Periode.
Die folgenden Pointer realisieren die vertikale und horizontale
Verkettung unter den Perioden. Man kann also ein level
von links (mostleft)
via neighbour nach rechts
durchsuchen. Das ist z. B. erforderlich, wenn für das Produkt zweier
Perioden die Darstellung als Linearkombination ermittelt werden soll.
Diese Linearkombination aus Perioden der Vorgängerstufe wird in einem Array namens
msummands abgespeichert, wobei die Vielfachheit der beteiligten
Perioden, m.a.W. der Koeffizient, mit dem jede Periode in der Summe
zu multiplizieren ist, parallel in
mfactors mitgeführt wird.
Die folgenden double-Werte werden für die Berechnung des
numerischen Werts der Periode benötigt. In numericValue wird
dieser Wert durch Lösen der quadratischen Gleichung ermittelt, in
referenceValue geschieht dies durch direktes Aufsummieren der
entsprechenden Cosinus-Werte.
Schließlich bekommt jede Periode noch einen name, der nach den
Konventionen gebildet wird, wie sie im vorigen Kapitel erläutert wurden.
Der Name wird gleich in LaTeX-Notation hinterlegt.
Konstruktor und Initialisierungen
Der Konstruktor für die Klasse Period ist klein gehalten, weil
fast alle Instanzen von Period durch Aufspaltung zustande kommen
und ihre Daten vom father
erben. Deswegen ist initialize
vom Konstruktor abgespalten.
public Period(String name, int len, int modulus, int primroot) {
this.len = len;
this.modulus = modulus;
this.name = name;
this.primroot = primroot;
exponents = new int[len];
}
public void initialize(Period leftest, double numV, int level, int nr) {
this.level = level;
this.seqnr = nr;
this.twoPowLevel = 1;
for (int i=0;i < level;i++)
this.twoPowLevel*=2;
this.mostleft = leftest;
this.numericValue = -1.0d;
}
Exponenten berechnen
Interessanter ist die folgende Methode, die die Gaußsche Periode erzeugt:
public void fill(int primitive) {
exponents_sorted = null;
exponents[0] = 1;
for (int i = 1; i < len; i++) {
exponents[i] = (exponents[i-1]*primitive) % modulus;
}
setReferenceValue();
}
Hier wird das Array exponents durch fortlaufendes Potenzieren
des primitiven Elements primitive gefüllt. Dabei wird bei jedem Schritt
sofort die Restklasse nach dem modulus gebildet.
Zum Abschluss wird der Referenzwert durch Aufruf von setReferenceValue
berechnet. Diese Methode ist wie folgt implementiert:
private void setReferenceValue() {
double re = 0.0d;
for (int i = 0; i < this.len; i++) {
int exp = this.exponents[i];
re += Math.cos(2 * exp * Math.PI / modulus);
}
this.referenceValue = re;
}
Es wird der Realteil der jeweiligen Einheitswurzel berechnet und die
Summe gebildet. Es wird m.a.W. der Ausdruck
\[ \sum_{p_{n,m}} \Re(\zeta^{3^k}) = \sum_{p_{n,m}} \cos(\frac{2\pi 3^k}{p})
\]
berechnet, wobei exp $=3^k$ durch alle Exponenten der jeweiligen
Gaußschen Periode läuft.
Der Referenzwert dient vor allem dazu, für den
Quadratwurzelausdruck das richtige Vorzeichen zu finden.
Aufspalten der Periode
Die folgende Methode gehört zu den Kernmethoden des Algorithmus. Sie realisiert das Aufspalten einer Gaußschen Periode in zwei neue Perioden von jeweils halber Länge.
public Period split(int which, int level, int seqnr) {
if (leftSplit == null) {
int l = this.len / 2;
int index = this.nr;
leftSplit = new Period("p_{" + level + "," + index + "}",l, modulus,primroot);
leftSplit.father = this;
leftSplit.level = level;
leftSplit.seqnr = seqnr;
leftSplit.nr = index;
leftSplit.twoPowLevel = this.twoPowLevel * 2;
seqnr++;
index += leftSplit.twoPowLevel/2;
rightSplit = new Period("p_{" + level + "," + index + "}",l, modulus,primroot);
rightSplit.father = this;
rightSplit.level = level;
rightSplit.seqnr = seqnr;
rightSplit.nr = index;
rightSplit.twoPowLevel = leftSplit.twoPowLevel;
int j = 0;
for (int i = 0; i < l; i++) {
leftSplit.exponents[i] = this.exponents[j];
j++;
rightSplit.exponents[i] = this.exponents[j];
j++;
}
leftSplit.setReferenceValue();
rightSplit.setReferenceValue();
}
if (which == 0) {
return leftSplit;
}
return rightSplit;
}
Die Methode liefert je nach dem Wert des Parameters which die
linke , d.h. die die den ersten Exponenten der alten Periode
enthält, oder die rechte neue Periode zurück. Es werden aber schon beide
neuen Perioden gebildet, die rechte sozusagen auf Vorrat .
Etwas kompliziert ist die Namensbildung der rechten Periode, denn nach
der Konvention des vorigen Kapitels wird ja die Periode $p_{n-1,m}$ in
$p_{n,m}$ und $p_{n,m+2^{n-1}}$ aufgespalten. Die entsprechende Berechnung
findet in der Zeile index += leftSplit.twoPowLevel/2 statt.
Die Methode trägt ihre eigene
Instanz this
als father für die beiden neuen Perioden
ein. Außerdem wird noch eine laufende Nummer seqnr mitgeführt. Diese
dient hauptsächlich dazu, bei den großen Moduli, 257 oder 65537
die Berechnung einer Periode abzubrechen, wenn klar ist, dass die Werte
später nicht mehr gebraucht werden. Das eigentliche Aufteilen des
father auf
seine Kinder in der for-Schleife erfolgt
ganz kanonisch.
Multiplizieren
Es gibt noch eine Methode add, die zwei Perioden addiert und die
nur im Rahmen eines Plausibilitätschecks verwendet wird. Diese Methode ist
nicht besonders interessant und soll hier nicht weiter erläutert werden.
Komplizierter und wichtiger ist die Multiplikation zweier Perioden.
public boolean mult(Period other) {
int[] temp = new int[modulus];
for (int i = 0; i < this.len; i++) {
for (int j = 0; j < other.len; j++) {
int result = (this.exponents[i] + other.exponents[j]) % modulus;
temp[result]++;
}
}
int le = modulus / (this.len + this.len);
msummands = new Period[le];
mfactors = new int[le];
Period start = this.father.mostleft;
int k = 0;
while (start != null) {
msummands[k] = start;
k++;
start = start.neighbour;
}
for (int i = 0; i < modulus; i++) {
while (temp[i] > 0) {
int j;
// search exponent i in msummands
for (j = 0; j < msummands.length; j++) {
if (msummands[j] != null && msummands[j].contains(i)) {
break;
}
}
if (j < msummands.length) {
for (int l = 0; l < msummands[j].len; l++) {
temp[msummands[j].exponents[l]]--;
}
mfactors[j]++;
}
else
{
System.out.println("% Root "+i+" = "+temp[i]+"not resolved");
return false;
}
}
}
// plausibilty check
for (int i = 0; i < modulus; i++) {
if (temp[i] != 0) {
return false;
}
}
simplify();
other.msummands = this.msummands;
other.mfactors = this.mfactors;
compute(other);
return true;
}
Zunächst wird ein Hilfs-Array temp angelegt, das das Ergebnis
der Multiplikation vorübergehend aufnimmt. Die Multiplikation selbst
erfolgt in der Schleife für jedes Paar von Summanden der beiden
Perioden this und
other, wobei die Exponenten addiert werden.
Die Restklassenbildung mit dem modulus erfolgt sofort nach der
Addition. Im Array temp wird mitgezählt, wie oft jeder Exponent
als Ergebnis vorkommt.
Nun kommt der wesentliche Schritt, denn jetzt muss das Ergebnis der
Multiplikation als Linearkombination von Perioden des Vorgänger-Levels
dargestellt werden. Dazu dienen die Arrays msummands und
mfactors,
die in der Länge modulus/2*len
angelegt werden. In der while-Schleife werden die potenziellen
Summanden in das Array msummands eingetragen. Dazu wird
start auf die erste,
linkeste Periode des Vorgänger-Levels
gesetzt, nämlich this.father.mostleft und über den Link zum
neighbour das ganze Level abgearbeitet.
Die nächsten beiden geschachtelten Schleifen, die äußere for-Schleife
und die innere while-Schleife realisieren die Suche nach den
Exponenten, die durch die Multiplikation entstanden sind, in den Perioden
des Vorgänger-Levels. Die for-Schleife arbeitet sukzessive jeden
Exponenten aus temp ab,
die while-Schleife berücksichtigt
die Vielfachheit, mit der der Exponent auftritt.
Die mit dem Kommentar search exponent i in msummands gekennzeichnete
Schleife tut genau das, sie sucht den Exponenten, der ja Ergebnis der Multiplikation
ist, in den Perioden des Vorgänger-Levels, die wir alle in msummands
hinterlegt haben. Dazu wird die Methode contains benutzt,
die wir später kurz vorstellen werden. Wenn die Suche erfolgreich ist, und sie muss aufgrund
unseres Beweises erfolgreich sein, wird die innere Schleife mit break
verlassen. Die Sicherheitsabfrage msummands[j] != null ist
übrigens erforderlich, weil wir manche Level nicht vollständig berechnen, dazu später mehr.
Wenn die Schleife verlassen wird, muss j < msummands.length sein.
Der else-Zweig der entsprechenden Abfrage wird nur dann durchlaufen,
wenn wir beim Abbruch der Rechnung innerhalb eines Levels einen Fehler
gemacht haben, d.h. wenn Perioden nicht berechnet wurden, obwohl sie
später noch gebraucht werden. Nach dem Debuggen des Programms sollte das
nicht mehr vorkommen. Im if-Zweig vertrauen wir nun wieder auf unsere
Theorie: Wir gehen davon aus, dass das Multiplikationsergebnis, wenn es den
Exponenten i
aus der Periode msummands[j]
enthält, alle
Exponenten aus msummands[j] enthält. Wir zählen für alle
diese Exponenten in dem Array temp die Vielfachheit um eins herunter
und notieren in mfactors, dass wir die gesamte Periode mit der
Vielfachheit eins als Bestandteil der Linearkombination erkannt haben.
Die abschließende, als plausibilty check gekennzeichnete Schleife
überprüft, ob tatsächlich alle Exponenten aus temp einer Periode
zugeordnet werden konnten. Mit Aufruf der Methode simplify wird
versucht, in der gerade gefundenen Linearkombination Summanden zusammenzufassen.
Dann wird die vereinfachte Linearkombination auch der zweiten an der
Multiplikation beteiligten Periode other zugewiesen. Die Darstellung
als Linearkombination gehört ja zum Produkt von this
und other,
ebenso wie als Summe von this und
other in beiden Perioden der
gleiche father eingetragen ist.
Schließlich wird noch compute aufgerufen, um die quadratische Gleichung
wirklich numerisch zu lösen – dazu gleich mehr.
Die Methode mult stellt somit den im vorigen Kapitel erwähnten
experimentellen Beweis für die Tatsache dar, dass das Ergebnis der
Multiplikation der beiden Perioden $p_{n,m} \cdot p_{n,m+2^{n-1}}$ eine Linearkombination aus Perioden des
Vorgänger-Levels ist.
Es soll nicht verschwiegen werden, dass man mit etwas mehr Aufwand, sowohl bei der zugrunde liegenden Theorie als auch bei der Implementierung, für die Linearkombinationen einfachere Ausdrücke erhalten kann, bei denen z.T. negative Koeffizienten vorkommen. Näheres dazu findet sich in [Bracke 2011]. Unser Anliegen war jedoch, sowohl die benutzte Theorie wie die Implementierung so einfach und durchsichtig wie möglich zu halten.
Vereinfachen und Vorzeichen setzen
In dem Programm wird mit der Methode simplify versucht,
den Ausdruck, der die Linearkombination darstellt, zu vereinfachen,
allerdings nur auf sehr primitivem Niveau:
private void simplify() {
for (int i = 0; i < msummands.length - 1; i++) {
if (mfactors[i] != 0 && msummands[i].seqnr % 2 == 0 && msummands[i + 1].seqnr == msummands[i].seqnr + 1) {
if (mfactors[i] == mfactors[i + 1]) {
if (msummands[i].father != null) {
mfactors[i + 1] = 0;
msummands[i + 1] = null;
msummands[i] = msummands[i].father;
}
}
}
}
}
Es wird lediglich überprüft, ob zwei unmittelbar aufeinanderfolgende
Perioden vom gleichen father abstammen und mit der gleichen
Vielfachheit vorkommen, dann wird die Summe dieser Perioden durch
ihren father, der ja die Summe darstellt, ersetzt. Die
Abfrage erfolgt hier durch Vergleich der seqnr, man könnte zwar
gleich die Links zum father vergleichen, müsste dann aber
die Sicherheitsabfrage father != null zusätzlich für den Fall
vornehmen, dass wegen ungleicher Vielfachheit eine Vereinfachung gar nicht
infrage kommt.
Es kann kein Zweifel bestehen, dass die beschriebene Methode bestenfalls einen
ersten Ansatz zur Vereinfachung der Ausdrücke darstellt. Wenn jemand Lust hat, hier
einen ausgefeilteren Algorithmus zu implementieren, nur zu! Im Erfolgsfall, wäre ich
für eine kurze Mitteilung sehr dankbar. Immerhin schafft es aber auch das simple
simplify beim 65537-Eck ca. 3400 Paare von Perioden zusammenzufassen.
Die Methode compute rechnet den Wert der beiden Perioden, die
wir gerade multipliziert haben,
als Wurzelausdruck aus. Wir erinnern uns, dass sich der Wert der Perioden
als Lösung einer quadratischen Gleichung ergibt:
\[
x_{0;1} = \frac{(x_0+x_1) \pm \sqrt{(x_0+x_1)^2-4x_0 x_1}}{2}\eqndot
\]
wobei $x_0+x_1$ der Wert der Summe der Perioden, also der Wert des father
ist und wir $x_0 x_1$ gerade mit einiger Mühe als Linearkombination aus Perioden dargestellt
haben, deren Werte ebenfalls bekannt sind. In der Methode compute
wird der Wert $x_0+x_1$ mit $p$ bezeichnet und der Wert des Produktes $x_0 x_1$ mit
$q,$ wie in der Lösungsformel für quadratische Gleichungen üblich.
public void compute(Period other){
double p = father.numericValue;
double q = 0.0d;
for (int i = 0; i < msummands.length; i++) {
if (mfactors[i] != 0) {
q = q + mfactors[i] * msummands[i].numericValue;
}
}
double D = p * p - 4 * q;
double root1 = (p + Math.sqrt(D)) / 2;
double root2 = (p - Math.sqrt(D)) / 2;
this.theP = p;
this.theQ = q;
if (Math.abs(root1 - this.referenceValue) < Math.abs(root2 - this.referenceValue)) {
this.numericValue = root1;
other.numericValue = root2;
this.sign = '+';
other.sign = '-';
} else {
this.numericValue = root2;
other.numericValue = root1;
this.sign = '-';
other.sign = '+';
}
}
Der Wert von p ergibt sich unmittelbar aus dem
father.numericValue.
Der Wert von q wird durch Summieren
der Werte der in der Linearkombination enthaltenen Perioden, jeweils
in der richtigen Vielfachheit, gebildet. Für die beiden Wurzeln root1
und root2 muss nur noch die Lösungsformel für die quadratische
Gleichung angewendet werden.
Anschließend wird der schon erwähnte Trick zur Bestimmung des richtigen
Vorzeichens angewendet, indem der soeben aus dem Wurzelausdruck errechnete
Wert mit dem über Cosinus vorher berechneten Referenzwert verglichen
wird. Je nach Ergebnis des Vergleichs bekommt entweder die
Periode this
oder die Periode other das positive bzw.
negative Wurzelvorzeichen zugeordnet. Der Trick wurde übrigens bereits von
Gauß selbst empfohlen (siehe [Gauß 1801]):
Hat man diese aus den Sinustafeln auf nur wenige Stellen berechnet, namlich soweit dass man entscheiden kann, welche größer, welche kleiner sind, so kann kein Zweifel mehr übrig sein, durch welche Bezeichnungen die einzelnen Wurzeln der Gleichung (A) zu unterscheiden sind.
Die Differenzen zwischen
root1
bzw. root2 und
dem referenceValue liegen
am Anfang der Berechnung in der Größenordnung von $10^{-12},$
beim 65537-Eck schaukeln sich mit der Zeit aber Rundungsfehler auf,
so dass bei den höheren Leveln Differenzen von ca. $10^{-4}$ auftreten.
Dazu später mehr.
Hilfs- und Druckroutinen
Die kleine Hilfsroutine contains stellt fest, ob ein gesuchter Exponent
unter den Exponenten dieser Periode vorhanden ist. Die Exponenten werden
dazu zunächst sortiert. Das mag aufwendig erscheinen, aber da contains
evtl. sehr oft aufgerufen wird und die Routine sort nur
sortiert, wenn exponents_sorted
noch null ist,
lohnt sich dieser Aufwand doch.
}
public boolean contains(int what) {
this.sort();
for (int i : this.exponents_sorted) {
if (i == what) {
return true;
}
if (i > what) {
break;
}
}
return false;
}
Eine weitere kleine Hilfsroutine, die ebenfals die sortierten Exponenten voraussetzt,
ist equals. Sie stellt fest, ob zwei Perioden gleich sind. Die
Routine wird nur für Plausibilitätschecks benutzt.
public boolean equals(Period other) {
if (this.len != other.len) {
return false;
}
if (this.modulus != other.modulus) {
return false;
}
this.sort();
other.sort();
for (int i = 0; i < this.len; i++) {
if (this.exponents_sorted[i] != other.exponents_sorted[i]) {
return false;
}
}
return true;
}
Unter den Druckroutinen soll hier nur der Zweig für die Ausgabe von LaTeX-Ausdrücken genauer betrachtet werden:
public void print(String arg) {
switch (arg) {
case "-tex":
String texstring0 = "\\[%s = %.12f\\]%n";
String texstring = "\\[%s = \\tfrac{%s %s \\sqrt{%s^2 - 4(%s)}}{2}\\]%n";
if (this.father == null){
System.out.printf(texstring0,this.name, this.numericValue);
break;
}
System.out.printf(texstring, this.name, father.name, this.sign, father.name, this.toString());
break;
...
}
}
Der Formatstring texstring formatiert den Quadratwurzelausdruck
dieser Periode im LaTeX-Format. Man sieht, dass für den Summanden
vor der sqrt sowie für den zu quadrierenden Summanden
unter der sqrt
jeweils der Name des father eingesetzt wird,
denn der father
ist die Summe der beiden Perioden, von denen eine
hier betrachtet wird. Das Produkt der beiden Perioden wird durch
die Linearkombination dargestellt, die in der Routine mult aufgebaut
wurde. Die Ausgabe dieser Linearkombination erfolgt durch Aufruf der
toString-Methode für diese Periode:
public @Override
String toString() {
StringBuilder buf = new StringBuilder();
if (msummands != null) {
for (int i = 0; i < msummands.length; i++) {
if (mfactors[i] != 0) {
if (mfactors[i] > 1) {
buf.append(Integer.toString(mfactors[i]));
}
buf.append(msummands[i].name);
buf.append("+");
}
}
}
buf.setLength(buf.length() - 1);
return buf.toString();
}
Der Aufbau des String erfolgt ebenfalls völlig kanonisch: Aus
msummands werden die Namen der an der Linearkombination beteiligten
Perioden ermittelt und mit dem Faktor aus mfactors versehen in
den Summenausdruck aufgenommen. In der vorletzten Zeile wird das überflüssige
letzte Pluszeichen beseitigt.
Beim 65537-Eck kann der Ausdruck sehr lang werden. Glücklicherweise
sind wir aber nie an die Grenzen der Kapazität eines StringBuilder
gestoßen. Für die Druckaufbereitung wird dieser lange String durch das
schon erwähnte separate Programm Kompressor aufgeteilt.
Nicht beschriebene Methoden von Period
Einige Methoden der Klasse Period werden hier nicht näher beschrieben.
Dazu zählen die Standardalgorithmen, die in sort
und quicksort
Anwendung finden, sowie die Routine cleanUp, die ursprünglich dazu
gedacht war beim 65537-Eck Heap-Speicher freizugeben, falls es
Engpässe geben sollte. Es hat sich aber gezeigt, dass bei Aufruf mit
-Xms512m -Xmx768m -XX:MaxPermSize=256m, wie oben beschrieben,
keine Probleme auftreten. Die Laufzeit hält sich ebenfalls in Grenzen: Das Programm wurde ursprünglich auf einem recht
schwachbrüstigen Rechner
(Baujahr 2008; 3 GHz Doppelkern-CPU mit 2GB RAM) entwickelt. Die in
Kap.7 bzw. im Original-Output des Programms zu findende
CPU-Zeit von ca. 33 Sek. beim 65537-Eck bezieht sich auf diesen Rechner. Auf einem etwas
aktuelleren Rechner (Baujahr 2016; 3.3 GHz Core i5-6500 mit 8 GB RAM) benötigt die Berechnung für das 65537-Eck rund 10 Sek.
Ebenfalls nicht beschrieben wird die Routine refactor. Sie ist eine Kopie von
compute,
benutzt aber zur Berechnung statt der Werte numericValue
des father
bzw. der Perioden aus msummands,
sondern referenceValue. Wie erwähnt, liegen die Rundungsfehler
bei den letzten Perioden des 65537-Eck in der Größenordnung von
$10^{-4},$ was normalerweise nicht tragisch wäre, doch ausgerechnet bei der
allerletzten benötigten Periode $p_{15,0}$ führt dies dazu, dass der Wert
unter der Wurzel negativ wird und als Ergebnis daher NaN herauskommt.
Daher muss für die allerletzte Periode mittels refactor neu aufgesetzt werden.
Näheres dazu wird ausführlich im Kap.8 zum 65537-Eck geschildert.
Das Hauptprogramm
Das Hauptprogramm analysiert zunächst die mitgelieferten Argumente, das muss
hier nicht weiter erläutert werden. Dann wird die
Gaußsche Periode der Stufe $0$ (hier p_base genannt)
initialisiert und durch Aufruf von p_base.fill(3) mit
den Exponenten gefüllt. Hier ist das primitive Element $3$
standardmäßig vorgegeben. Nun folgt die
Haupt-Schleife, in der die vertikal und horizontal verlinkten Perioden
aufgebaut werden:
public static void main(String[] args) {
long startTimeNano = System.nanoTime();
System.err.println("% Job started.");
System.out.print("% Running with arguments: ");
for (String arg : args) {
System.out.print(arg + " ");
}
...
Hermes cyc = new Hermes();
int base = fermatPrime - 1;
Period p_base = new Period("p_{0,0}", base, fermatPrime, 3);
p_base.initialize(p_base, -1.0d, 0, 0);
p_base.fill(3);
int periodLen = base;
int level = 0;
Period actV = p_base;
p_base.print(args);
while (periodLen > 2) {
periodLen /= 2;
Period actH = actV;
Period baseH = null;
Period leftest = null;
level++;
int nr = 0;
while (actH != null) {
Period left = actH.split(0, level, nr);
if (baseH != null) {
baseH.neighbour = left;
} else {
leftest = left;
}
baseH = left;
left.mostleft = leftest;
Period right = actH.split(1, level, nr);
if (!actH.equals(left.add(right))) {
System.err.println("% Error when splitting!");
}
baseH.neighbour = right;
baseH = right;
right.mostleft = leftest;
// the kernel of the whole program
if (left.mult(right)) {
left.print(args);
right.print(args);
} else {
System.err.println("Error in multiplication!");
}
if (periodLen == 2) {
finalPrint(left, right, fermatPrime, args);
break;
}
actH = actH.neighbour;
nr += 2;
if (level == 14 && nr > 16) {
break;
}
if (level == 13 && nr > 1900) {
break;
}
}
actV = leftest;
}
long taskTimeNano = System.nanoTime() - startTimeNano;
System.out.println("% Time used: "
+ taskTimeNano / 1000000000.0d + "sec");
System.err.println("% Job finished.");
}
Die beiden Perioden actV für
vertikal und
actH für
horizontal werden benutzt, um durch die im Aufbau befindlichen Listen
zu navigieren. Die äußere while-Schleife bewegt sich vertikal, also von Stufe zu Stufe,
die innere while-Schleife erzeugt innerhalb der Stufe durch
Halbieren der Perioden aus der Vorgänger-Stufe die neuen Perioden und
verkettet sie über neighbour horizontal. Dabei ist der direkte
Zugriff auf neighbour
und mostleft, die sich ja
in einer fremden Klasse befinden, sicher nicht lupenreines Java, aber hier
setter-Methoden bereitzustellen würde den Code nur aufblähen und
vom Wesentlichen ablenken.
Das Wesentliche wird etwas pathetisch durch den Kommentar
the kernel of the whole program gekennzeichnet. Hier werden die
soeben durch Aufspalten erzeugten Perioden left
und right
miteinander multipliziert, wobei wir gesehen haben, was alles
hinter dieser Multiplikation steckt.
Anschließend wird gleich die Ausgabe für die beiden Perioden angestoßen. Das liegt einerseits nahe, weil alle benötigten Informationen zu diesem Zeitpunkt vollständig vorhanden sind – es wird halt sukzessiv gearbeitet, andererseits bestand zunächst die Sorge, dass beim 65537-Eck der Heap-Speicher nicht ausreichen könnte. Die Information sollte rasch aufs Papier , damit man sie so schnell wie möglich dem Garbage-Collector vorwerfen kann.
Damit ist die wesentliche Arbeit geleistet, es kommen nur noch einige
Sonderbehandlungen. Für die letzte Stufe mit der Periodenlänge $2$
müssen wir nur die erste Periode berechnen, denn deren Wert ist, wie wir
gesehen haben, $2\cos(2\pi/p).$ Also wird bei periodLen == 2
nur noch die finale Ausgabe angestoßen und die innere Schleife
verlassen. Weil die Periodenlänge $2$ erreicht ist, führt das zum
Verlassen der äußeren Schleife und damit zum Ende des Programms.
Die beiden if-Abfragen auf
level == 14 bzw.
level == 13
sollen die Datenflut beim 65537-Eck etwas
eindämmen. Diese Begrenzungen wurden erst nach dem ersten Probelauf des
Programms eingebaut, nachdem man feststellen konnte, dass die Perioden, die
hier weggelassen werden, für die Berechnung von $p_{15,0}$ nicht gebraucht
werden. Man spart damit immerhin in Level 13 die Ausgabe von
mehr als 6000 Zeilen, in Level 14 sogar mehr als 16000 Zeilen.
Wir haben schon erwähnt, dass im Programm Überprüfungen stattfinden, die
fälschlicherweise entfernte Perioden entdecken.
Man beachte, dass als Kriterium für den Abbruch die laufenden Nummer nr
herangezogen wird. Die letzten berechneten Wurzeln im jeweiligen Level sind nach
dem verwendeten Nummerierungs-Schema $p_{13,3164}$ und $p_{14,3072}.$
Bleibt noch die Routine finalPrint nachzutragen:
public static void finalPrint(Period left, Period right,
int fermatPrime, String[] args){
System.out.println("% 1/2 * " + left.name + " = " + left.numericValue / 2);
System.out.println("% cos(2*pi/" + fermatPrime + "): " + Math.cos(2 * Math.PI / fermatPrime));
if (Double.isNaN(left.numericValue)) {
System.out.println("% Taking reference values!");
left.refactor(right);
left.print(args);
System.out.println("% 1/2 * " + left.name + " = " + left.numericValue / 2);
}
}
Es wird der über Wurzelausdrücke berechnete Wert der Periode ausgegeben und
zum Vergleich der über Cosinus berechnete Wert. Anschließend
erfolgt die Sonderbehandlung, falls die Rundungsfehler dazu geführt haben,
dass der Wert der letzten Periode NaN ist. In diesem Fall wird
die Methode refactor aufgerufen und der Wert noch einmal
ausgegeben. Es sei nochmals bemerkt, dass diese Sonderbehandlung nur bei
der allerletzten Periode des 65537-Ecks stattfindet.