4 Algorithmen

4.1 GCD-Algorithmen

4.1.1 EUKLID’s Algorithmus

Der Algorithmus von EUKLID berechnet den größten gemeinsamen Teiler 𝑑 = gcd(𝑎,𝑏) zweier natürlicher Zahlen 𝑎,𝑏 ∈ N. Die grundlegende, iterativ angewendete Rechenoperation dabei ist modulare Division.³⁵ Algorithmus 1 beschreibt das klassische Verfahren [16, 31, 7, 22, 2]. Es beginnt unter der Voraussetzung 𝑎 > 𝑏 mit einer Modulo-Division 𝑐₂ = 𝑎 mod 𝑏. Im nächsten Schritt wird 𝑐₂ als Modul verwendet und 𝑐₃ = 𝑏 mod 𝑐₂ berechnet, wonach 𝑐₄ = 𝑐₂ mod 𝑐₃ folgt usw. . Diese (wegen 𝑐_𝑘+1 < 𝑐_𝑘) absteigende Sequenz endet wenn 𝑐_𝑘+1 = 0 wird³⁶ – das Ergebnis 𝑑 befindet sich dann in 𝑐_𝑘.

Beweis Jeder Iterationsschritt 𝑐_𝑖+1 = 𝑐_𝑖−1 mod 𝑐_𝑖 kann in der Umkehrung (vgl. Restklassenbeziehung 22 in Abschnitt 3) als

geschrieben werden. So gesehen wird durch den Algorithmus der folgende Abstieg vorgenommen:

bis 𝑐_𝑘+1 verschwindet.

Warum 𝑑 = 𝑐_𝑘 ein gemeinsamer Faktor von 𝑎 und 𝑏 ist, wird klar wenn man die Folge rückwärts betrachtet. In der letzten Zeile steht 𝑐_𝑘−1 = 𝑞_𝑘+1𝑑, also teilt 𝑑 den Rest 𝑐_𝑘−1 (oder kürzer 𝑑∣𝑐_𝑘−1). Wenn 𝑑 aber als Faktor in 𝑐_𝑘−1 enthalten ist, dann kann man 𝑑 auf der rechten Seite der vorletzten Gleichung ausklammern, weshalb es auch als Faktor in 𝑐_𝑘−2 vorkommen muß. Dies setzt sich bis in die erste Gleichung fort (sämtliche Divisionsreste 𝑐₀,𝑐₁,…,𝑐_𝑘 enthalten folglich 𝑑 als Teiler), in der die Startbedingung 𝑐₀ = 𝑎 und 𝑐₁ = 𝑏 verankert ist. Deshalb ist der gemeinsame Teiler 𝑑 sowohl in 𝑎 als auch 𝑏 enthalten.

Viel kürzer kann man damit argumentieren, daß die Modulo-Division 𝑐_𝑖+1 = 𝑐_𝑖−1 mod 𝑐_𝑖 eine GCD erhaltende Operation ist. Denn mit der Zerlegung 𝑐_𝑖 = 𝑐 _𝑖𝑑 gilt ausgehend von Formel 22:

d. h. der gemeinsame Teiler 𝑑 in 𝑐_𝑖−1 und 𝑐_𝑖 ist auch in 𝑐_𝑖+1 wieder enthalten.³⁷

(36)

Damit 𝑑 wirklich den größten gemeinsamen Teiler stellt, muß es überhaupt alle gemeinsamen Teiler enthalten. Mit dem Ziel dies nachzuweisen betrachten wir nochmals die Folge beginnend mit der vorletzten Gleichung, welche nach 𝑐_𝑘 = 𝑑 umgestellt wird. Ersetzt man darin 𝑐_𝑘−1 mit Hilfe der vorvorletzten Gleichung und fährt aufsteigend fort, so erhält man eine lineare Darstellung für 𝑑 = gcd(𝑎,𝑏),

welche auf ganzen Zahlen 𝛼,𝛽 ∈ Z sowie 𝑐₀ = 𝑎 und 𝑐₁ = 𝑏 beruht.

Mit Hilfe von Formel 37, welche auch Satz von BéZOUT genannt wird,³⁸ kann jetzt relativ einfach bewiesen werden, daß 𝑑 wirklich der größte gemeinsame Teiler ist. Nehmen wir dazu an, es gäbe einen weiteren gemeinsamen Teiler 𝑑′. Dann würde dieser (auf der rechten Seite von Gleichung 37) als Faktor von 𝑎 und 𝑏 auszuklammern sein und deshalb (wenn man die linke Seite betrachtet) als Teiler von 𝑑 auftreten. Mit anderen Worten stecken alle weiteren Teiler (schon) in 𝑑, weshalb nur dieser der größte gemeinsame Teiler sein kann.

Hinweis Es existieren unendlich viele lineare Darstellungen für 𝑑 als Linearkombination von 𝑎 und 𝑏 (vgl. auch die kurze Betrachtung zu diophantischen Gleichungen auf Seite §). Jede Substitution der Art 𝛼 := 𝛼 +𝑛𝑏 und 𝛽 := 𝛽 −𝑛𝑎, mit

erfüllt BéZOUT’s Identität 37 ebenso.

Die kleinsten Werte für |𝛼| und |𝛽| zeichnen sich folglich durch die Relationen |𝛼| < 𝑏 und |𝛽| < 𝑎 aus. Um sie zu ermitteln kann man entweder 𝑎 und 𝑏 sukzessive von 𝛼 und 𝛽 subtrahieren/addieren oder man reduziert die BéZOUT-Kofaktoren mit Hilfe von 𝑞 = bzw. 𝑞 = und folgender Formeln:³⁹

4.1.2 Erweiterter euklidischer Algorithmus

Der erweiterte euklidische Algorithmus erlaubt eine effiziente Berechnung der BéZOUT-Kofaktoren 𝛼 und 𝛽 zusammen (und gleichzeitig) mit dem größten gemeinsamen Teiler (siehe auch [16, 31, 7, 2]). Dazu berücksichtigt er einige Erkenntnisse aus dem vorigen Abschnitt, insbesondere daß:

• der größte gemeinsame Teiler 𝑑 = gcd(𝑎,𝑏) für 𝛼,𝛽 ∈ Z als Linearkombination 𝛼𝑎 +𝛽𝑏 darstellbar ist;
• 𝑑 für 0 ≤ 𝑖 ≤ 𝑘 jeden Rest 𝑐_𝑖 teilt und damit die Darstellung 𝑐_𝑖 = 𝑑𝑐 _𝑖 = 𝛼_𝑖𝑎 + 𝛽_𝑖𝑏 rechtfertigt;
• die Berechnung von 𝑐_𝑖+1 mit Hilfe von 𝑞_𝑖+1 = und 𝑐_𝑖+1 = 𝑐_𝑖−1−𝑞_𝑖+1𝑐_𝑖 erfolgen kann.

Durch Induktion ergibt sich aus

(39)

wenn man 𝑐_𝑖−1 = 𝛼_𝑖−1𝑎+𝛽_𝑖−1𝑏 und 𝑐_𝑖 = 𝛼_𝑖𝑎+𝛽_𝑖𝑏 berücksichtigt:

Vergleich mit Ausgangsformel 39 erlaubt die Bestimmung von 𝛼_𝑖+1 und 𝛽_𝑖+1.

Mit den Startwerten 𝛼₀ = 1 und 𝛽₀ = 0 (gewährleistet 𝑐₀ = 𝑎) bzw. 𝛼₁ = 0 und 𝛽₁ = 1 (ebenso für 𝑐₁ = 𝑏) kann man den erweiterten euklidischen Algorithmus 2 formulieren.

Require:   𝑎 ≥ 𝑏
  (𝑐₀,𝛼₀,𝛽₀) ⇐ (𝑎,1,0)
  (𝑐₁,𝛼₁,𝛽₁) ⇐ (𝑏,0,1)
  𝑘 ⇐ 0
  repeat
   𝑘 ⇐ 𝑘 +1
   𝑞 ⇐ {Integer-Division}
   (𝑐_𝑘+1,𝛼_𝑘+1,𝛽_𝑘+1) ⇐ (𝑐_𝑘−1,𝛼_𝑘−1,𝛽_𝑘−1) −𝑞(𝑐_𝑘,𝛼_𝑘,𝛽_𝑘) {𝑐_𝑘+1 = 𝑐_𝑘−1 mod 𝑐_𝑘}
  until  𝑐_𝑘+1 = 0
  (𝑑,𝛼,𝛽) ⇐ (𝑐_𝑘,𝛼_𝑘,𝛽_𝑘)

Er endet ganz genauso wie Algorithmus 1 bei 𝑐_𝑘+1 = 0, ermittelt aber zusätzlich die BéZOUT-Kofaktoren 𝛼 = 𝛼_𝑘 und 𝛽 = 𝛽_𝑘 sowie die teilerfremden Anteile 𝑎 = 𝑎/𝑑 und 𝑏 = 𝑏/𝑑. Letztere findet man wegen 𝑐_𝑘+1 = 𝛼_𝑘+1𝑎+𝛽_𝑘+1𝑏 = 0 am Ende in 𝛽_𝑘+1 und 𝛼_𝑘+1. Dies wird relativ schnell aus

ersichtlich, wenn man auf beiden Seiten der letzten Gleichung das kleinste gemeinsame Vielfache anvisiert.

Hinweis Als Ergänzung zum Hinweis von Seite § kann man sogar feststellen, daß die „Wahl“ von 𝑞 nicht unbedingt mit einer Modulo-Division verbunden sein muß. Etwas anders könnte man 𝑞 einfach so wählen, daß

• die Folge der 𝑐_𝑘 weiterhin eine absteigende Sequenz ist;
• und in jedem Iterationsschritt 𝑞𝑐_𝑘 ≤ 𝑐_𝑘−1 gilt (wodurch 𝑐_𝑘+1 ≥ 0 gewährleistet ist).

Den steilsten Abstieg 𝛥 = 𝑐_𝑘+1 − 𝑐_𝑘 erreicht man allerdings, wenn 𝑞 im jeweiligen Iterationsschritt maximal ist. Dies ist aber unter der Voraussetzung 𝑐_𝑘+1 ≥ 0 genau bei einer Modulo-Division 𝑐_𝑘+1 = 𝑐_𝑘−1 mod 𝑐_𝑘 der Fall. Ansonsten ist jeder Wert 𝑞 = 0… geeignet — einzig die Konvergenzgeschwindigkeit nimmt mit kleineren Werten von 𝑞 ab.⁴⁰

4.1.3 EUKLID’s Algorithmus für Polynome über Z₂

Für Polynome mit Koeffizienten aus GF(2) kann man ganz ähnlich verfahren [11, S. 2.3.6],[31, S. 17.3], [7, S. 2.2.1], [2, S. 2.1]. Mit Rückblick auf den Hinweis von Seite § muß man nur folgendes beachten:

1. Addition und Subtraktion sind gleichwertig (und die Addition wird zu einem Exklusiv-Oder ⊕).
2. An die Stelle eines wertbehafteten Vergleichs 𝑐_𝑘−1 ≥ 𝑐_𝑘 tritt ein Vergleich des Polynomgrades deg 𝑐_𝑘−1 ≥ deg 𝑐_𝑘.
3. Die Vorbedingung deg 𝑎 ≥ deg 𝑏 kann entfallen, da der Algorithmus selbst die Vertauschung der Argumente vornimmt.

4. Die Ermittlung von 𝑞(𝑥) läßt sich in folgende Einzelfälle zerlegen:

   • deg 𝑐_𝑘−1(𝑥) < deg 𝑐_𝑘(𝑥) führt zu 𝑞(𝑥) = 0, was nach Algorithmus 3 einer Vertauschung (𝑐_𝑘,𝛼_𝑘,𝛽_𝑘)⟺(𝑐_𝑘−1,𝛼_𝑘−1,𝛽_𝑘−1) entspricht.
   • deg 𝑐_𝑘−1(𝑥) = deg 𝑐_𝑘(𝑥) ergibt 𝑞(𝑥) = 1, was mit (𝑐_𝑘+1,𝛼_𝑘+1,𝛽_𝑘+1) ⇐ (𝑐_𝑘−1,𝛼_𝑘−1,𝛽_𝑘−1) +(𝑐_𝑘,𝛼_𝑘,𝛽_𝑘) im nachfolgenden Fall aufgeht.
   • Für deg 𝑐_𝑘−1(𝑥) > deg 𝑐_𝑘(𝑥) müßte man 𝑞(𝑥) = berechnen. Mit der Option auch ein kleineres 𝑞(𝑥) zu verwenden (siehe Hinweis auf Seite §), kann man sich eine vollständige Polynomdivision sparen. Stattdessen ist eine sinnvolle Vorgehensweise sich bei der Divison auf das höchstwertige Bit von 𝑐_𝑘−1(𝑥) und 𝑐_𝑘(𝑥) beschränken. In einer monischen Polynomrepräsentation (der höchstwertige Koeffizient ist 1) führt die Division so zu: 𝑥^{deg 𝑐_𝑘−1−deg 𝑐_𝑘}.

4.1.4 Binärer GCD-Algorithmus

Der binäre GCD-Algorithmus kommt im Gegensatz zum klassischen euklidischen Algorithmus ohne Divisionen aus [7, 16, 32, 34]. Statt dessen wird (beginnend mit 𝑎₀ = 𝑎 und 𝑏₀ = 𝑏) durch Subtraktion und Halbierung eine stetige Reduktion der Argumente 𝑎 und 𝑏 vorgenommen, welche letztlich zum größten gemeinsamen Teiler 𝑑 = gcd(𝑎,𝑏) führt.⁴¹

Der Algorithmus endet spätestens dann, wenn im Schritt 𝑖 = 𝑘 eines der beiden Argumente 𝑎_𝑖 oder 𝑏_𝑖 verschwindet (im obigen Fall beispielhaft für 𝑏_𝑘 = 0). Das Reduktionsschema zeigt Tabelle 2.

𝑎𝑖	𝑏𝑖	gcd(𝑎𝑖+1,𝑏𝑖+1) =	Begründung
gerade		2 gcd	gemeinsamer Teiler ist 2
gerade	ungerade	gcd	2 ist kein gemeinsamer Teiler
ungerade	gerade	gcd	2 ist kein gemeinsamer Teiler
ungerade, 𝑎𝑖 ≥ 𝑏𝑖		gcd	𝑎𝑖 −𝑏𝑖, wie auch 𝑎𝑖 +𝑏𝑖, enthalten den Teiler gcd(𝑎𝑖,𝑏𝑖) = 𝑑, denn 𝑎𝑖 = 𝑎𝑖𝑑, 𝑏𝑖 = 𝑏𝑖𝑑 führt zu gcd(𝑎𝑖 −𝑏𝑖,𝑏𝑖) = gcd[𝑑(𝑎𝑖 −𝑏𝑖),𝑏𝑖𝑑]. Außerdem sind sowohl Summe als auch Differenz gerade Zahlen (𝑎𝑖 = 2𝜈 +1, 𝑏𝑖 = 2𝜇 +1 ergibt 𝑎𝑖 −𝑏𝑖 = 2(𝜈 −𝜇)), weshalb 2 als gemeinsamer Teiler wieder ausgeschlossen werden kann.
ungerade, 𝑏𝑖 ≥ 𝑎𝑖		gcd(𝑎𝑖,)	siehe Fall 𝑎𝑖 > 𝑏𝑖
𝑎𝑖 > 0	0	𝑎𝑖	𝑎𝑖 ist größter gemeinsamer Teiler mit 0
0	𝑏𝑖 > 0	𝑏𝑖	𝑏𝑖 ist größter gemeinsamer Teiler mit 0

Anmerkungen

1. Bezüglich 𝑎_𝑘−1 gibt es genau eine Situation, die das Verschwinden von 𝑎_𝑘 hervorruft (genauso bezüglich 𝑏_𝑘−1 und 𝑏_𝑘) . Betrachtet man dazu Tabelle 2, so wird 𝑎_𝑘−1 in folgenden Fällen reduziert:

   1. 𝑎_𝑘−1 und 𝑏_𝑘−1 sind gerade: Halbierung von 𝑎_𝑘−1 = 1 ergibt 0, bedeutet aber ungerades 𝑎_𝑘−1 (Widerspruch).
   2. 𝑎_𝑘−1 gerade, 𝑏_𝑘−1 ungerade: Halbierung von 𝑎_𝑘−1 = 1 ergibt 0, was ebenfalls ungerades 𝑎_𝑘−1 bedeutet (Widerspruch).
   3. 𝑎_𝑘−1 und 𝑏_𝑘−1 sind ungerade: Halbierung von 𝑏_𝑘−1 − 𝑎_𝑘−1 = 1 ergibt 0, aber unter der Voraussetzung 𝑏_𝑘−1 = 𝑎_𝑘−1 + 1 können niemals beide ungerade sein (Widerspruch).
   4. 𝑎_𝑘−1 und 𝑏_𝑘−1 sind ungerade: Halbierung von 𝑏_𝑘−1 = 𝑎_𝑘−1 führt zu 𝑎_𝑘 = 0 und stellt damit den einzig möglichen Fall im Schritt 𝑘 −1 dar.
2. In jedem Iterationszyklus wird entweder 𝑎_𝑖 oder 𝑏_𝑖 um mindestens ein Bit reduziert.⁴² Aus diesem Umstand läßt sich (im Fall 𝑎₀ > 𝑏₀) für die Anzahl der Iterationsschritte ≤ 𝑘 ≤ 2 schlußfolgern (vgl. auch [15, Theorem 2]).
3. Die Halbierung des Arguments für den Fall, daß 𝑎_𝑖 und 𝑏_𝑖 ungerade sind, muß man nicht unbedingt im Iterationsschritt 𝑖 ausführen. Es ist durchaus legitim, falls beispielsweise 𝑎_𝑖 ≥ 𝑏_𝑖 gilt, einfach nur 𝑎_𝑖+1 = 𝑎_𝑖 −𝑏_𝑖 zu berechnen. Die Differenz ergibt ja bekanntlich wieder eine gerade Zahl (dann ist 𝑎_𝑖+1 gerade, 𝑏_𝑖+1 ungerade) und es kommt im nächsten Iterationsschritt zu der gewünschten Halbierung.

4. Den Fall, daß sowohl 𝑎_𝑖 als auch 𝑏_𝑖 gerade sind, kann man grundsätzlich aus der Haupt-Iterationsschleife herausziehen. Denn wird einmal eines der beiden Argumente ungerade, dann kann dieser Fall niemals wieder eintreten (genau das ungerade Argument ist in allen anderen Reduktionsfällen unveränderlich). Solange beide Argumente gerade sind, kann man sie also kontinuierlich reduzieren, bis nach 𝑛 Schritten entweder 𝑎_𝑖 oder 𝑏_𝑖 ungerade geworden ist. Danach kann auf 𝑎_𝑛 = 𝑎/2^𝑛 und 𝑏_𝑛 = 𝑏/2^𝑛 irgendein binärer (erweiterter) GCD-Algorithmus angewendet werden, der gcd(𝑎_𝑛,𝑏_𝑛) = 𝛼_𝑛𝑎_𝑛 +𝛽_𝑛𝑏_𝑛 ermittelt.

   

   Der so ermittelte größte gemeinsame Teiler muß am Schluß nur noch mit 2^𝑛 multipliziert werden (also um 𝑛 Bit verschoben), was auch Algorithmus 4 entsprechend wiedergibt.⁴³ Wir können deshalb in den Betrachtungen zum erweiterten binären GCD-Algorithmus (siehe nächste Abschnitte) immer voraussetzen, daß 𝑎 oder 𝑏 ungerade ist.

     𝑛 ⇐ 0
     while 𝑎 gerade ∧ 𝑏 gerade do
      𝑎 ⇐ 𝑎/2
      𝑏 ⇐ 𝑏/2
      𝑛 ⇐ 𝑛+1
     end while
     𝑑 ⇐ 2^𝑛 gcd(𝑎,𝑏) {GCD-Algorithmus mit Voraussetzung: 𝑎 oder 𝑏 ungerade}

5. Für den Fall, daß entweder 𝑎 oder 𝑏 ungerade ist (oder beide), kann man ausgehend von 𝑎 = 𝑑𝑎 und 𝑏 = 𝑑𝑏 feststellen:⁴⁴

   • Der größte gemeinsame Teiler 𝑑 ist ungerade.
   • Sollte 𝑎 ungerade sein, dann ist es dessen teilerfremde Anteil 𝑎 ebenfalls.
   • Sollte 𝑏 ungerade sein, dann ist es dessen teilerfremde Anteil 𝑏 auch.

4.1.5 Erweiterter binärer GCD-Algorithmus

Algorithmus nach KALISKI⁴⁵ [15] Betrachtet man das Reduktionsschema des binären GCD-Algorithmus, so kann man verschiedenste lineare Transformationen der Art

definieren, wobei sich die einzelnen Betrachtungsweisen durch unterschiedliche Werte in den Koeffizienten unterscheiden. Bei jedem Schritt sind so 𝑎_𝑖 und 𝑏_𝑖 als lineare Funktion von 𝑎_𝑖+1 und 𝑏_𝑖+1 darstellbar, die Ausgangsgrößen 𝑎₀ und 𝑏₀ deshalb als lineare Funktionen von 𝑎_𝑖 und 𝑏_𝑖.

Durch Einsetzen von 𝑎_𝑖 und 𝑏_𝑖 in folgende Gleichung,

gefolgt von einem Koeffizientenvergleich, kann man für die konkreten Reduktionsfälle die jeweilige Transformationen der Linearfaktoren ableiten (siehe Tabelle 3). Am Beispiel 𝑎_𝑖+1 = (𝑎_𝑖 −𝑏_𝑖)/2, 𝑏_𝑖+1 = 𝑏_𝑖 soll das Vorgehen exemplarisch verdeutlicht werden.

Aus den Transformationen nach Tabelle 3 kann man außerdem die Determinante

bestimmen.

𝑎𝑖	𝑏𝑖	𝑎𝑖+1 =	𝑏𝑖+1 =	𝑢𝑖+1 =	𝑣𝑖+1 =	𝑠𝑖+1 =	𝑡𝑖+1 =	Di+1 =
gerade				2𝑢𝑖	2𝑣𝑖	2𝑠𝑖	2𝑡𝑖	4Di
gerade	ungerade		𝑏𝑖	2𝑢𝑖	𝑣𝑖	2𝑠𝑖	𝑡𝑖	2Di
ungerade	gerade	𝑎𝑖		𝑢𝑖	2𝑣𝑖	𝑠𝑖	2𝑡𝑖	2Di
ungerade, 𝑎𝑖 ≥ 𝑏𝑖			𝑏𝑖	2𝑢𝑖	𝑢𝑖 +𝑣𝑖	2𝑠𝑖	𝑠𝑖 +𝑡𝑖	2Di
ungerade, 𝑏𝑖 ≥ 𝑎𝑖		𝑎𝑖		𝑢𝑖 +𝑣𝑖	2𝑣𝑖	𝑠𝑖 +𝑡𝑖	2𝑡𝑖	2Di

Um den Zusammenhang mit den BéZOUT-Koeffizienten 𝛼,𝛽 ∈ Z in gcd(𝑎,𝑏) = 𝛼𝑎₀ +𝛽𝑏₀ herzustellen, stellen wir die Formeln 40 und 41 nach 𝑎_𝑖 und 𝑏_𝑖 um. Dazu werden beide Gleichungen mit den Linearfaktoren der jeweils anderen multipliziert und dann wechselweise voneinander subtrahiert.

Schließt man den Fall aus, daß 𝑎₀ und 𝑏₀ gerade sind (D_i = 2^𝑖, vgl. Anmerkungen auf Seite §), dann kann für den letzten Iterationsschritt 𝑖 = 𝑘, in Abhängigkeit davon ob 𝑎_𝑘 oder 𝑏_𝑘 zuerst verschwindet, folgendermaßen konkretisiert werden:

An diesem Punkt stellen wir jedoch fest, daß es sich bei den Größen 𝑢_𝑘, 𝑣_𝑘, 𝑠_𝑘 und 𝑡_𝑘 nicht um die Kofaktoren von gcd(𝑎,𝑏), sondern um eine Linearfaktordarstellung von 2^𝑘 gcd(𝑎,𝑏) handelt.⁴⁶

Um aus 𝑡_𝑘, 𝑣_𝑘 und 𝑠_𝑘, 𝑢_𝑘 die BéZOUT-Koeffizienten 𝛼 und 𝛽 zu bestimmen, sind zusätzliche Korrekturschritte nötig, welche in [15] auch Korrekturphase (oder Phase II) genannt werden. Ziel ist es dabei, die rechte Seite der letzten Gleichung durch 2^𝑘 zu dividieren (oder in 𝑘 Schritten wiederholt durch 2). Da die linke Seite immer eine gerade Zahl ist, können folgende Schlußfolgerungen im Falle 𝑏_𝑘 = 0 gezogen werden (falls 𝑎_𝑘 = 0 war, äquivalent für 𝑡_𝑘 → 𝑠_𝑘 und 𝑣_𝑘 → 𝑢_𝑘):⁴⁷

• Sind 𝑣 und 𝑡 gerade, dann ist eine ganzzahlige Halbierung 𝑣 := 𝑣^(g)/2, 𝑡 := 𝑡^(g)/2 möglich (d. h. führt wieder zu einer ganzen Zahl).

• In allen anderen Fällen (𝑡 bzw. 𝑣 ungerade) kann man wegen

  

  die Reduktion 𝑡 := (𝑡+𝑏)/2 bzw. 𝑣 := (𝑣+𝑎)/2 vornehmen, ohne daß sich die Gleichung ändert. Man wandelt jedoch die ungerade Zahl 𝑡 bzw. 𝑣 in eine gerade Zahl, die dann (wie gewünscht) durch 2 teilbar ist. Die Ursache liegt einerseits darin begründet, daß die linke (und demzufolge auch rechte) Seite der Gleichung

  

  immer eine gerade Zahl repräsentiert und andererseits, entweder 𝑎 oder/und 𝑏 als ungerade vorausgesetzt wurden. Berücksichtigt man die folgenden Gesetzmäßigkeiten:

  1. die Summe zweier ungerader oder zweier gerader Zahlen ist eine gerade Zahl;
  2. die Summe einer ungeraden und einer geraden Zahl ist eine ungerade Zahl;
  3. das Produkt zweier ungerader Zahlen ist eine ungerade Zahl;
  4. alle anderen Kombinationen ergeben bei der Multiplikation eine gerade Zahl;

  dann sind genau drei Fälle konstruierbar, in denen 𝑡 oder 𝑣 ungerade ist.

  

  In allen diesen Varianten führt aber die Addition 𝑡 +𝑏 bzw. 𝑣+𝑎 zu einer geraden Zahl, womit sich das Korrekturverfahren bestätigt.

Require:   𝑎 ungerade ∨ 𝑏 ungerade {Phase I}
  (𝑎₀,𝑢₀,𝑠₀) ⇐ (𝑎,1,0)
  (𝑏₀,𝑣₀,𝑡₀) ⇐ (𝑏,0,1)
  𝑓 ⇐ 0 {Flag, daß anzeigt, ob schlußendlich 𝛼 oder 𝛽 negativ ist}
  𝑘 ⇐ 0
  while 𝑎_𝑘 > 0 do
   if 𝑏_𝑘 > 𝑎_𝑘 then
   (𝑏_𝑘,𝑣_𝑘,𝑡_𝑘)⟺(𝑎_𝑘,𝑢_𝑘,𝑠_𝑘) {gewährleistet 𝑎_𝑘 ≥ 𝑏_𝑘}
   𝑓 ⇐ 𝑓 {invertiere Flag}
   end if
   if 𝑎_𝑘 ungerade ∧ 𝑏_𝑘 ungerade then
   (𝑎_𝑘,𝑣_𝑘,𝑡_𝑘)(𝑎_𝑘 −𝑏_𝑘,𝑣_𝑘 +𝑢_𝑘,𝑡_𝑘 +𝑠_𝑘) {𝑎_𝑘 ist jetzt gerade, 𝑏_𝑘 weiterhin ungerade}
   end if
   if 𝑎_𝑘 gerade then
   (𝑎_𝑘+1,𝑢_𝑘+1,𝑠_𝑘+1)(𝑎_𝑘/2,2𝑢_𝑘,2𝑠_𝑘) {𝑎_𝑘 gerade, 𝑏_𝑘 ungerade}
   else
   (𝑏_𝑘+1,𝑣_𝑘+1,𝑡_𝑘+1)(𝑏_𝑘/2,2𝑣_𝑘,2𝑡_𝑘) {𝑎_𝑘 ungerade, 𝑏_𝑘 gerade}
   end if
   𝑘 ⇐ 𝑘 +1
  end while
  gcd(𝑎,𝑏) ⇐ 𝑏_𝑘 {Teilergebnis 𝑑 = gcd(𝑎,𝑏)} {Phase II}
  𝑖 ⇐ 𝑘
  while 𝑖 > 0 do
   if 𝑠_𝑘 ungerade ∨ 𝑢_𝑘 ungerade then
   𝑠_𝑘𝑠_𝑘 +𝑡_𝑘 {Addition von 𝑡_𝑘 = 𝑏 (vgl. Bemerkung 3)}
   𝑢_𝑘𝑢_𝑘 +𝑣_𝑘 {Addition von 𝑣_𝑘 = 𝑎 (vgl. Bemerkung 3)}
   end if
   𝑠_𝑘𝑠_𝑘/2 {Korrektur 𝛼}
   𝑢_𝑘𝑢_𝑘/2 {Korrektur 𝛽}
   𝑖 ⇐ 𝑖 −1
  end while
  (𝛼,𝛽)(𝑠_𝑘,𝑢_𝑘) {Wenn 𝑓 = 0, dann 𝛼 < 0, sonst 𝛽}

Als Ergebnis kann man Algorithmus 5 formulieren, wobei außerdem folgende Anmerkungen berücksichtigt wurden:

1. Aus den Gleichungen 40 und 41 läßt sich im letzten Iterationsschritt (von Phase I)

   

   schlußfolgern, d. h. bei 𝑢_𝑘 und 𝑠_𝑘 handelt es sich im Fall 𝑏_𝑘 = 0 (gleichermaßen für 𝑣_𝑘, 𝑡_𝑘 im Fall 𝑎_𝑘 = 0) um die teilerfremden Faktoren

   
2. Wegen der stetigen Reduktion von 𝑎_𝑖 oder 𝑏_𝑖 müssen die Linearfaktoren 𝑢_𝑖, 𝑠_𝑖, 𝑣_𝑖 und 𝑡_𝑖 in Phase I schrittweise anwachsen, denn nur so können 𝑎₀ und 𝑏₀ nach Gleichung 40 und 41 konstant bleiben. In [15, Theorem 1] wird bewiesen, daß keine dieser Größen während der Ausführung des Algorithmus den Maximalwert 2𝑎₀ −1 überschreitet (für 𝑎₀ > 𝑏₀).⁴⁸

3. Die in Phase II vorzunehmende Addition von 𝑎 bzw. 𝑏 kann (entsprechend Anmerkung 5 auf Seite §) ersetzt werden durch eine Addition von 𝑎₀ bzw. 𝑏₀.⁴⁹

   

Algorithmus nach PENK [16, Exercise 4.5.2.39], [22, S. 14.4.3] Dieser Algorithmus stellt 𝑎_𝑖+1 und 𝑏_𝑖+1 als lineare Funktion von 𝑎_𝑖 und 𝑏_𝑖 dar, letztlich wird also von (𝑎₀,𝑏₀) auf (𝑎_𝑖,𝑏_𝑖) geschlossen.

Vorteilhaft wirkt sich aus, daß im letzten Reduktionsschritt (wenn 𝑎_𝑘 oder 𝑏_𝑘 verschwindet) 𝑢_𝑘 und 𝑣_𝑘 bzw. 𝑠_𝑘 und 𝑡_𝑘 direkt die gesuchten BéZOUT-Kofaktoren 𝛼, 𝛽 darstellen.

Durch Einsetzen von 𝑎_𝑖+1 und 𝑏_𝑖+1 (entsprechend Tabelle 2) in die Gleichungen

kann man für den jeweiligen Reduktionsfall die zugehörige Transformationen der Linearfaktoren ableiten. Wieder am Beispiel 𝑎_𝑖+1 = (𝑎_𝑖 −𝑏_𝑖)/2, 𝑏_𝑖+1 = 𝑏_𝑖 soll das Vorgehen veranschaulicht werden (𝑎_𝑖 und 𝑏_𝑖 sind ungerade, 𝑎_𝑖 ≥ 𝑏_𝑖).

Vergleich von linker und rechter Seite läßt den Schluß zu:

Die anderen Kombinationen können genau nach demselben Schema abgeleitet werden. Tabelle 4 faßt die Ergebnisse in übersichtlicher Form zusammen.⁵⁰

𝑎𝑖	𝑏𝑖	𝑎𝑖+1 =	𝑏𝑖+1 =	𝑢𝑖+1 =	𝑣𝑖+1 =	𝑠𝑖+1 =	𝑡𝑖+1 =
gerade	ungerade		𝑏𝑖			𝑠𝑖	𝑡𝑖
ungerade	gerade	𝑎𝑖		𝑢𝑖	𝑣𝑖
ungerade, 𝑎𝑖 ≥ 𝑏𝑖			𝑏𝑖			𝑠𝑖	𝑡𝑖
ungerade, 𝑏𝑖 ≥ 𝑎𝑖		𝑎𝑖		𝑢𝑖	𝑣𝑖

Dabei tritt allerdings wieder das bekannte Problem auf: Wie kann man die Ganzzahligkeit des jeweiligen Linearfaktoren bei der Halbierung wahren? Die Antwort haben wir schon beim vorangegangenen Algorithmus geliefert – indem die Ausgangsgleichungen 42 und 43 folgendermaßen erweitert:

und dadurch ungerade Zahlen 𝑢_𝑖, 𝑣_𝑖 bzw. 𝑠_𝑖, 𝑡_𝑖 in gerade umwandelt.

Beschränken wir uns auf die Kombinationen nach Tabelle 4, in denen 𝑢_𝑖 und 𝑣_𝑖 verändert werden (äquivalent für 𝑢_𝑖 → 𝑠_𝑖 und 𝑣_𝑖 → 𝑡_𝑖). Für den Fall, daß 𝑎_𝑖 gerade und 𝑏_𝑖 ungerade ist (zweite Zeile in 4), können wir auf die Argumentation von Seite § zurückgreifen (Phase II der Methode nach KALISKI). Sie erlaubt uns die Subtraktion/Addition von 𝑎₀ und 𝑏₀ für den Fall, daß 𝑢_𝑖 oder 𝑣_𝑖 ungerade ist. Die Situation, daß 𝑎_𝑖 und 𝑏_𝑖 ungerade sind (vorletzte Zeile in 4), kann man durch gedankliche Verzögerung der Halbierung in den nächsten Iterationsschritt erklären. Wählt man als modifizierten Einzelschritt 𝑎_𝑖+1 = 𝑎_𝑖 −𝑏_𝑖 (und entsprechend 𝑢_𝑖+1 = 𝑢_𝑖 −𝑠_𝑖, 𝑣_𝑖+1 = 𝑣_𝑖 −𝑡_𝑖), so reduziert sich die Fragestellung wieder auf eine gerade Zahl 𝑎_𝑖+1 = 𝑢_𝑖+1𝑎₀ +𝑣_𝑖+1𝑏₀, also auf den vorangegangenen Fall.

Als Ergebnis der Ausführungen kann man Algorithmus 6 formulieren. Der Vorteil des Algorithmus (gegenüber KALISKI’s) liegt vor allem darin, daß keine Korrekturphase nötig ist. Nachteilig für eine praktische Umsetzung ist die notwendige Vorzeichen-Arithmetik.

Require:   𝑏 ungerade
  (𝑎₀,𝑢₀,𝑣₀) ⇐ (𝑎,1,0)
  (𝑏₀,𝑠₀,𝑡₀) ⇐ (𝑏,0,1)
  𝑖 ⇐ 0
  while 𝑎_𝑖 > 0 do
   if 𝑏_𝑖 > 𝑎_𝑖 then
   (𝑏_𝑖,𝑠_𝑖,𝑡_𝑖)⟺(𝑎_𝑖,𝑢_𝑖,𝑣_𝑖) {gewährleistet 𝑎_𝑖 ≥ 𝑏_𝑖}
   end if
   if 𝑎_𝑖 ungerade ∧ 𝑏_𝑖 ungerade then
   (𝑎_𝑖,𝑢_𝑖,𝑣_𝑖)(𝑎_𝑖 −𝑏_𝑖,𝑢_𝑖 −𝑠_𝑖,𝑣_𝑖 −𝑡_𝑖) {𝑎_𝑖 ist jetzt gerade, 𝑏_𝑖 weiterhin ungerade}
   end if
   if 𝑎_𝑖 gerade then
   𝑎_𝑖+1𝑎_𝑖/2 {𝑎_𝑖 gerade, 𝑏_𝑖 ungerade}
   if 𝑢_𝑖 gerade ∧ 𝑣_𝑖 gerade then
   𝑢_𝑖+1𝑢_𝑖/2
   𝑣_𝑖+1𝑣_𝑖/2
   else
   𝑢_𝑖+1(𝑢_𝑖 +𝑏)/2
   𝑣_𝑖+1(𝑣_𝑖 −𝑎)/2
   end if
   else
   𝑏_𝑖+1𝑏_𝑖/2 {𝑎_𝑖 ungerade, 𝑏_𝑖 gerade}
   if 𝑠_𝑖 gerade ∧ 𝑡_𝑖 gerade then
   𝑠_𝑖+1𝑠_𝑖/2
   𝑡_𝑖+1𝑡_𝑖/2
   else
   𝑠_𝑖+1(𝑠_𝑖 +𝑏)/2
   𝑡_𝑖+1(𝑡_𝑖 −𝑎)/2
   end if
   end if
   𝑖 ⇐ 𝑖 +1
  end while
  (𝑑,𝛼,𝛽) ⇐ (𝑏_𝑖,𝑠_𝑖,𝑡_𝑖) {𝑑 = gcd(𝑎,𝑏) = 𝛼𝑎+𝛽𝑏}

4.2 Lineare diophantische Gleichungen

Gleichungen mit ausschließlich ganzzahligen Lösungen 𝑥,𝑦,… ∈ Z nennt man diophantisch, wobei sie im linearen Fall (für zwei Variablen) die folgende Form haben:

(44)

Solche Gleichungen haben genau dann eine Lösung, wenn der größte gemeinsame Teiler 𝑑 = gcd(𝑎,𝑏) den Wert 𝑐 teilt. Im Zusammenhang mit EUKLID’s Algorithmus wurde dies praktisch schon nachgewiesen – auch daß es unendlich viele solcher Lösungen gibt, kam in Abschnitt 4.1.1 zur Sprache.

Zuerst bemerken wir, daß die Differenz zweier Lösungen (𝑥₁,𝑦₁) und (𝑥₂,𝑦₂) die homogene Gleichung 𝑎𝑥 +𝑏𝑦 = 0 erfüllt.

Mit 𝑎 = 𝑑𝑎 und 𝑏 = 𝑑𝑏 läßt sich sogar schreiben

und es liegen die Lösungen (𝑥,𝑦) = (𝑛𝑏,−𝑛𝑎), 𝑛 ∈ Z auf der Hand (Einsetzen ergibt 𝑎𝑥+𝑏𝑦 = 𝑎𝑛𝑏−𝑏𝑛𝑎 = 0).⁵¹ Findet man jetzt noch eine partikuläre Lösung (𝑥₂,𝑦₂), dann ergeben sich alle weiteren zu:

Partikuläre Lösungen (𝑥₂,𝑦₂) für 𝑎𝑥₂ +𝑏𝑦₂ = 𝑐 haben wir aber schon mittels der erweiterten GCD-Algorithmen zur Verfügung, denn mit der BéZOUT’s Identität (den Index 2 jetzt weggelassen)

gilt:

Eine partikuläre Lösung (𝑥₂𝑥, 𝑦₂𝑦) kann deshalb mit Hilfe der BéZOUT-Kofaktoren 𝛼,𝛽 gegeben werden.

Die Vielfalt aller Lösungen stellt sich dadurch wiefolgt dar:

Für den Spezialfall 𝑑 = gcd(𝑎,𝑏) = 1 entartet die Formel zu:

4.3 Chinesischer Restsatz

4.3.1 Hilfssatz für zwei Kongruenzen

Um sich dem Chinesischen Restsatz⁵² zu nähern, betrachten wir zunächst einen etwas einfacheren Fall, der die grundsätzliche Fragestellung jedoch beinhaltet: Welche natürliche Zahl 𝑧 erfüllt die folgenden beiden Kongruenzen:

wenn vorausgesetzt wird, daß 𝑛,𝑚 > 0 relativ prim zueinander sind?

Um sie zu beantworten formulieren wir den Ausgangspunkt zuerst einmal entsprechend Restklassenbeziehung 22:

Deren Darstellung als

(48)

zeigt mit Verweis auf die Form von 44, daß es sich um eine lineare diophantische Gleichung handelt. Wegen gcd(𝑛,𝑚) = 1 könnte die zugehörige Lösungsformel 46 zwar sofort zur Anwendung kommen – naheliegend (da kurz) ist jedoch auch die Anwendung von BéZOUT’s Identität. Denn multipliziert man gcd(𝑛,𝑚) = 𝛼𝑛+𝛽𝑚 = 1 mit 𝑐, dann kann aus

durch Vergleich mit 48 sofort abgelesen werden, daß 𝑥 = 𝛼𝑐 und 𝑦 = −𝛽𝑐 gelten muß.⁵³ Mit der Erkenntnis aus Formel 46, daß es sich bei den Lösungen einer linearen diophantischen Gleichungen immer um eine ganze Lösungsmenge handelt, resultiert:⁵⁴

Durch Einsetzen in Formel 47 erhält man

und so eine geschlossene Darstellung für die Lösungsmenge.

(50)

Mit der von den Restklassen bekannten Kongruenz 23 für teilerfremde Zahlen:⁵⁵

können wir die Probe machen.

Die Lösungsmenge 𝑧_𝑘 bildet demzufolge eine Restklasse 𝑚𝑛 .

(51)

4.3.2 Ein System von Kongruenzen

Nehmen wir jetzt ein ganzes System von Kongruenzen an:

wobei gcd(𝑚_𝑖,𝑚_𝑘) = 1 für 𝑖 ≠ 𝑘 gelten soll. Wie schon im Fall von zwei Kongruenzen stellt man die Frage, welche Lösung 𝑥 kongruent zu allen 𝑎_𝑖 modulo 𝑚_𝑖 ist (𝑖 = 1,…,𝑛). Ohne einen exakten mathematischen Beweis anzutreten, scheint als Schlußfolgerung aus Abschnitt 4.3.1 einleuchtend, daß die Lösung(en) 𝑥 eine Restklasse modulo 𝑚 =  _𝑛𝑚_𝑖 darstellen [16, 7].⁵⁶

Bezeichnen wir mit 𝑚_𝑖 das Produkt aller Moduli ausgeschlossen 𝑚_𝑖, also

dann gilt wegen der Teilerfremdheit der einzelnen Moduli gcd(𝑚_𝑖,𝑚_𝑖) = 1 bzw. mit dem Satz von BéZOUT:⁵⁷

Im Gegensatz dazu verschwindet 𝛼_𝑘𝑚_𝑘 mod 𝑚_𝑖 für 𝑖 ≠ 𝑘, denn 𝑚_𝑖 ist als Faktor in 𝑚_𝑘 enthalten. Die Orthogonalität beider Fälle kann mit Hilfe des KRONECKER-Symbols 𝛿_𝑖𝑘 ausgedrückt werden:

Die endgültige Lösungsidee besteht nun darin, für jede Kongruenz 𝑖 den folgenden Ausdruck zu bilden:

Wie zu sehen, ist die Summe auf der linken Seite kongruent zu 𝑎_𝑖 modulo 𝑚_𝑖, was

(52)

als finales Ergebnis rechtfertigt.

Für den einfachen Fall 𝑛 = 2 von Abschnitt 4.3.1 kann man mit 𝑚₁ = 𝑚₂, 𝑚₂ = 𝑚₁ sowie

relativ einfach Kongruenz 51 verifizieren.⁵⁸

4.4 Quadratwurzeln in F_𝑝

4.4.1 Vorbetrachtungen

Das Bestimmen der Quadratwurzel kann man im Körper F_𝑝 recht effizient mit dem TONELLI-SHANKS Algorithmus erledigen. Bevor man aber dazu übergehen kann den Algorithmus zu erläutern, muß man sich unbedingt über einige Fakten in Bezug auf das Quadrieren von Körperelementen klar werden. Ohne wesentlich an Allgemeinheit zu verlieren, betrachten wird dazu die multiplikative Gruppe Z_𝑝^∗ = (Z_𝑝 \ {0}) und formulieren zuerst die folgenden Grundaussagen in Bezug auf:⁵⁹

1. Aus der Beziehung 𝑦 = 𝑥² (mod 𝑝) ersieht man sofort, daß immer zwei Elemente 𝑥,𝑥′ genau zum selben 𝑦 führen – ein „Positives“ und ein „Negatives“: 𝑥′ ≡ −𝑥 (mod 𝑝).
2. Da Z_𝑝^∗ genau 𝑝 − 1 Elemente enthält (besitzt die Ordnung |Z_𝑝^∗| = |Z_𝑝| − 1 = 𝑝 − 1), können (𝑝 − 1)/2 Elemente 𝑦 ∈ Z_𝑝^∗ keine Quadratwurzel besitzen. Denn umgekehrt betrachtet sind ja die verfügbaren 𝑝 − 1 Elemente 𝑥 wegen der Paarbildung beim Quadrieren (welche zu (𝑝−1)/2 Elementen 𝑦 führt) schon „aufgebraucht“.

3. Für die (𝑝−1)/2 Elemente 𝑦, welche der Beziehung 𝑦 = 𝑥² (mod 𝑝) genügen, gilt wegen des kleinen Satzes von FERMAT (vgl. Abschnitt 3.3.2):

   (53)

4. Für die restlichen (𝑝−1)/2 Elemente 𝑦, welche keine zugeordnete Quadratwurzel 𝑥 besitzen, gilt hingegen:⁶⁰

   (54)

   Für einen kurzen Beweis nehmen wir uns ein bestimmtes (konstantes) Element 𝑦 heraus und betrachten die Zerlegung 𝑦 = 𝑎𝑥 (mod 𝑝), welche es in Z_𝑝^∗ immer geben muß: 𝑎 = 𝑥⁻¹𝑦 (mod 𝑝). Wenn man nun für 𝑥 alle möglichen Werte 1…𝑝 − 1 setzt, dann muß auch 𝑎 entsprechend variieren. Wegen:

   • der Abgeschlossenheit der Multiplikation;
   • der Eindeutigkeit des Inversen;
   • und der Kommutativität von Elementen

   in Z_𝑝^∗ wird 𝑎 einerseits alle möglichen Werte 1…𝑝−1 annehmen und andererseits wird jedes Produkt 𝑎𝑥 insofern doppelt vorkommen, daß 𝑎 und 𝑥 nur die Rolle getauscht haben. In Abbildung 3 sind mögliche Produkte für 𝑦 = 8 (hat keine Wurzel) und 𝑦 = 9 (Wurzeln sind ±3) beispielhaft in Z₁₁^∗ dargestellt (wobei gestrichelte Linien Produkte kennzeichnen, bei denen 𝑎 und 𝑥 nur die Rolle getauscht haben).

(a) 𝑦 = 8 (b) 𝑦 = (±3)² = 9

Multiplizieren wir alle Produkte 𝑎𝑥 miteinander und schließen dabei die Doppelungen aus (d. h. berücksichtigen nur die Produkte in Abbildung 3a, welche mit durchgezogenen Linien dargestellt sind), dann hat man alle Elemente 𝑥 ∈ Z_𝑝^∗ miteinander multipliziert.

Bei den (𝑝−1)/2 Produkten 𝑎_𝑖𝑥_𝑖 handelt es sich aber genau um die Potenz

und so kann man mit dem Satz von WILSON (siehe Formel 27) den Beweis abschließen:

Noch bevor man also versucht die Quadratwurzel 𝑥 = zu berechnen, kann man mit Hilfe des gerade eingeführten (EULER-) Kriteriums prüfen, ob diese überhaupt existiert. Die Nomenklatur für das EULER-Kriterium stammt allerdings von A.-M. LEGENDRE und wird als LEGENDRE-Symbol bezeichnet:

(55)

4.4.2 Der Spezialfall 𝑝 mod 4 = 3

Im speziellen Fall 𝑝 ≡ 3 (mod 4) kann man eine sehr einfache Lösung für die Quadratwurzel 𝑥 = erhalten, sofern sie denn existiert [28, Theorem 9.3], [1, S. D.1.4]:

(56)

mit 𝑠 = (𝑝−3)/4.

Ausgangspunkt für einen kurzen Beweis könnte die folgende Darstellung für das Modul 𝑝 sein:⁶¹

Wird sie im Zusammenhang mit dem EULER-Kriterium von Gleichung 53 angewendet, dann bestätigt sich Formel 56 wiefolgt:

4.4.3 Der TONELLI-SHANKS Algorithmus

Idee Der TONELLI-SHANKS Algorithmus nach [33, 30] reduziert, ausgehend von einem Anfangswert 𝑥₀ ∈ Z_𝑝>2 , den Exponenten 2^𝑡_𝑖 im Ausdruck

(57)

immer weiter, bis schließlich 𝑡_𝑛 = 0 wird und deshalb 𝑥_𝑛² ≡ 𝑦 (mod 𝑝) die gesuchte Lösung darstellt (siehe z. B. auch [24, S. 2.9], [7, S. 2.3.2]). Mit 2^𝑡_𝑖 soll es sich aber nicht um irgendeine Zweierpotenz handeln, sondern um die Ordnung des Elements 𝑧_𝑖 ≡ 𝑥_𝑖²𝑦⁻¹ (mod 𝑝).

Es gilt also einen Algorithmus zu beschreiben, bei dem sowohl für den Startwert 𝑧₀ als auch für alle Werte 𝑧_𝑖 in

(58)

die Ordnung |𝑧_𝑖| = 2^𝑡_𝑖 einerseits (und überhaupt) eine Zweierpotenz ist und andererseits (auch noch) abnimmt.⁶²

Startwerte Als Anfangswerte setzt man:

Dabei entspricht der Wert 2^𝑡₀′ nicht zwingendermaßen der Ordnung |𝑧₀| = 2^𝑡₀, er befriedigt aber Kongruenz 58 grundsätzlich:

1. Da es sich bei 𝑝 um eine Primzahl mit 𝑝 > 2 handeln soll (𝑝 ist zwingendermaßen ungerade), kann man

   (59)

   setzen (mit 𝑟 ≥ 1, maximal) und erhält so:

   

2. Soll die Quadratwurzel 𝑥 = wirklich existieren (wovon wir ausgehen) dann gilt mit dem EULER-Kriterium 55:

   

3. Faßt man beide Punkte zusammen, dann bestätigt sich:

   (60)

Interpretieren wir Kongruenz 60 jetzt mit dem Wissen aus Abschnitt 2.1.2, dann muß es sich bei 2^𝑡₀′ entweder um die Ordnung |𝑧₀| selbst oder aber um ein Vielfaches davon handeln: 2^𝑡₀′ = 𝐾|𝑧₀|. Deshalb können 𝐾 und |𝑧₀| nur Zweierpotenzen sein (was 𝐾 = 1 nicht ausschließt), d. h. für die Ordnung kann man ohne weitere Bedenken |𝑧₀| = 2^𝑡₀ ansetzen (mit 𝑡₀ ≤ 𝑡₀^′). Um nun aus 𝑡₀^′ den exakten Exponenten 𝑡₀ bzw. die Ordnung |𝑧₀| zu ermitteln, kann man 𝑧₀ einfach so oft quadrieren bis (nach FERMAT’s kleinem Satz) die Bedingung 𝑧₀^|𝑧₀| ≡ 1 (mod 𝑝) erfüllt ist.

Iteration Bevor wir zum iterativen Teil des Algorithmus übergehen können, benötigen wir noch die Mithilfe irgendeines (zufällig gewählten) Elements 𝛼 ∈ Z_𝑝, welches keine Wurzel in Z_𝑝^∗ haben darf. Für 𝛼 muß das EULER-Kriterium 55 dann den Wert

annehmen. Aus diesem Element 𝛼 erzeugen wir ein neues Element 𝛽 = 𝛼^𝑠 (mod 𝑝) mit den Eigenschaften:

• potenziert zu 2^𝑟 ergibt sich mit FERMAT’s kleinem Satz:

  

• aus dem EULER-Kriterium für 𝛼 läßt sich schlußfolgern:

  (61)

Nun wird die Näherung schrittweise nach folgender Vorschrift verbessert:⁶³

(62)

wobei sich der Exponent 𝑡_𝑖 jeweils um 1 verringert (𝑡_𝑖+1^′ = 𝑡_𝑖 −1, mit 2^{𝑡_𝑖+1′} einem mglw. Vielfachen der Ordnung |𝑧_𝑖+1|).

Um die Verringerung des Exponenten nachzuweisen bilden wir einfach 𝑧_𝑖+1^{2𝑡_𝑖−1} mit Hilfe von Iterationsformel 62:

Beide Faktoren auf der rechten Seite sind aber vom Wert −1 (mod 𝑝), was letztlich

bedeutet. Zur Begründung folgende Argumentation:

• 𝛽^2𝑟−1 ≡ −1 (mod 𝑝) gilt wegen des EULER-Kriteriums für 𝛼 (siehe Kongruenz 61).

• Der Beweis für 𝑧_𝑖^{2𝑡_𝑖−1} ≡ −1 (mod 𝑝) stützt sich auf die Darstellung:

  

  welche ja bei jeder Iteration als +1 (mod 𝑝) vorausgesetzt wird. Aus diesem Grund kann der Wert innerhalb des Klammerausdrucks nur +1 oder −1 sein [24, Hilfssatz 2.10]. Wenn 2^𝑡_𝑖 aber die Ordnung von 𝑧_𝑖 derstellt, dann muß 𝑡_𝑖 der kleinste Wert sein, welcher zu 𝑧_𝑖^2𝑡_𝑖 ≡ +1 (mod 𝑝) führt. Deshalb kann 2^{𝑡_𝑖−1} nicht (ebenfalls) die Ordnung von 𝑧_𝑖 sein und es bleibt nur die Schlußfolgerung: 𝑧_𝑖^{2𝑡_𝑖−1} ≡ −1 (mod 𝑝).

Aus dem letzten Punkt ergibt sich die Notwendigkeit, daß man am Anfang des Iterationszyklus jedesmal die exakte Ordnung |𝑧_𝑖| bzw. den zugehörigen Exponenten 𝑡_𝑖 ≤ 𝑡_𝑖^′ bestimmt (vgl. auch die Erläuterungen zu den Startwerten).

Algorithmus Bevor wir in Algorithmus 7 alle Erkenntnisse zusammenfassen, noch zwei Hinweise zur Umsetzung:

1. Die Berechnung des inversen Elements 𝑦⁻¹ kann man einmalig zu Beginn des Algorithmus durchführen.

2. Setzt man für die Reduktion der Ordnung in jedem Schritt 𝛥_𝑖 = 𝑡_𝑖−1 − 𝑡_𝑖 und führt in Formel 62 die Abkürzungen

   

   ein, dann kann man auch 𝛾_𝑖 rekursiv berechnen:

   

   Um den korrekten Startwert 𝛾₀ = 𝛽^{2𝑟−𝑡₀−1} zu erzielen, vergleichen wir jetzt noch die Werte für 𝑖 = 0:

   

   und erhalten:

   

Require:   𝑝 > 2
Ensure:   𝑥 =
  𝑟,𝑠 ⇐ aus 𝑝−1 = 2^𝑟𝑠, mit 𝑟 maximal
  repeat
   𝛼 ⇐ zufälliges Element aus Z_𝑝
  until 𝛼^(𝑝−1)/2 mod 𝑝 = 𝑝−1 {EULER-Kriterium: 𝛼^(𝑝−1)/2 ≡ −1 (mod 𝑝)}
  𝛾 ⇐ 𝛼^𝑠 mod 𝑝 {𝛾₋₁ = 𝛽}
  𝑥 ⇐ 𝑦^(𝑠+1)/2 mod 𝑝 {𝑥₀}
  𝑦_inv = 𝑦⁻¹ mod 𝑝
  𝑡^′ ⇐ 𝑟 −1
  𝑡 ⇐ aus 2^𝑡 = ord(𝑥²𝑦_inv mod 𝑝) {Ordnung von 𝑧₀}
  while 𝑡 > 0 do
   𝛥 ⇐ 𝑡^′ −𝑡
   𝛾 ⇐ 𝛾^2𝛥 mod 𝑝 {𝛾_𝑖 = 𝛽^{2𝑟−𝑡_𝑖−1} }
   𝑥 ⇐ 𝑥𝛾 mod 𝑝
   𝑡^′ ⇐ 𝑡
   𝑡 aus 2^𝑡 = ord(𝑥²𝑦_inv mod 𝑝) {Ordnung von 𝑧_𝑖 ≡ 𝑥_𝑖²𝑦⁻¹ (mod 𝑝)}
  end while {Probe (siehe auch [21, Appendix C]),
  anstatt zu Beginn das EULER-Kriterium (𝑦|𝑝) = 1 zu verifizieren}
  𝑦′ ⇐ 𝑥² mod 𝑝
  if 𝑦′ = 𝑦 then
   return𝑦
  else
   returnkeine Lösung
  end if

4.5 Quadratische Gleichungen in F_2^𝑛

4.5.1 Problemstellung

Quadratische Gleichungen (mit 𝑎 ≠ 0) der Form

kann man bezüglich 𝑦 in F_2^𝑛 nicht so einfach lösen wie in R. Wie noch zu sehen sein wird, hat diese Gleichung in F_2^𝑛:

• möglicherweise keine Lösung;
• genau eine Lösung für den Fall 𝑏 = 0;
• oder sogar zwei Lösungen.

Zuersteinmal kann man sie unter der Voraussetzung 𝑎,𝑏 ≠ 0 mit Hilfe der Substitution 𝑧 = 𝑎𝑦/𝑏 und der naheliegenden Abkürzung 𝛽 = 𝑎𝑐/𝑏² in eine allgemeine Form überführen:

Davon ausgehend können wir uns den mathematischen Grundlagen einer Lösung zuwenden [2, 31].

4.5.2 Trace

Zuallererst soll eine Hilfsformel eingeführt werden, welche in der einen oder anderen Form immer wieder benötigt wird.⁶⁴ Dazu sei die folgende Summe 𝜆_𝑖 für ein Element 𝛼 ∈ F_𝑝^𝑛 definiert und mit 𝑖 ∈ N indiziert (vorerst soll 𝑖 ≠ 0 gelten):

Dann berechnen wir unter Zuhilfenahme des „Anfänger-Traums“ nach Formel 33 bzw.  34 die 𝑝-te Potenz:

und stellen fest:

(63)

Definition Der Trace ist die Summe der Konjugierten eines Elements im Körper F_𝑝^𝑛.

Mit Formel 63 ist der Trace auch darstellbar als:

(64)

Eigenschaften Mit Hilfe der Polynomdarstellung für ein Element 𝛼 ∈ F_𝑝^𝑛

kann man einige interessante Eigenschaften und Beziehungen für den Trace ableiten:

1. Aus den Formeln 63 und 64 ergeben sich mit dem kleinen Satz von FERMAT (𝛼^𝑝𝑛 −𝛼 = 0) die Beziehungen:

   

   d. h. der Trace ist eine lineare Abbildung F_𝑝^𝑛 → F_𝑝. Der Grund liegt ganz einfach darin, daß Gleichung 65 nur im Körper F_𝑝 Gültigkeit hat (nicht in F_𝑝^𝑛) – und deshalb: tr(𝛼) ∈ F_𝑝.

2. Der Trace eines Elements wird nicht verändert, wenn man das Element vorher (mehrfach) zur Potenz 𝑝 potenziert:

   

   Warum dies so ist, erkennt man durch Anwendung der im Punkt 1 hergeleiteten Formel 65 sowie (wieder) des „Anfänger-Traums“ nach Formel 34:

   

3. Der Trace ist linear:

   1. tr(𝛼 +𝛽) = tr(𝛼) +tr(𝛽), mit 𝛼,𝛽 ∈ F_𝑝^𝑛

      In ausführlicher Schreibweise des Trace und (wieder) mit Formel 34 kann man auch dies schnell beweisen.

      

   2. tr(𝑐𝛼) = 𝑐tr(𝛼), mit 𝑐 ∈ F_𝑝, 𝛼 ∈ F_𝑝^𝑛

      Ausklammern von 𝑐 = 𝑐^𝑝 = 𝑐^𝑝𝑖 aus der Summendarstellung des Trace führt zu:

      

Im Körper F_2^𝑛 Speziell im Binärkörper F_2^𝑛 entsteht aus Formel 63:

(67)

und mit der Rekursion 𝜆_𝑖 = 𝜆_𝑖−1 +𝛼^2𝑖−1 die (insbesondere für einen Berechnungsalgorithmus) interessante Beziehung:⁶⁵

(68)

Insbesondere aus den Punkten 1 und 2 auf Seite § ergeben sich weitere spezialisierte Eigenschaften:

1. Quadrieren des Arguments 𝛼 in tr(𝛼) ändert nach Formel 66 nichts am Wert des Trace.

   (69)

2. Das Potenzieren eines Trace ändert ebenfalls nichts an seinem Wert.⁶⁶

   (70)

   Grundlage des Beweises bilden die folgenden beiden Reduktionsfälle:

   

   Sie erlauben es, den Exponenten iterativ bis auf 1 zu reduzieren, was letztlich zu tr ^𝑖(𝛼) = tr(𝛼) führt.

   In einer kurzen Fallbetrachtung für 𝑖 gerade:

   

   bzw. 𝑖 ungerade (und so 𝑖 −1 gerade):

   

   kann man die Begründung finden.

3. Aus der Linearität des Trace und aus Formel 66 ergibt sich direkt:

   

4. Eine bemerkenswerte Formel, insbesondere im Zusammenhang mit der Lösung quadratischer Gleichungen in F_2^𝑛, ist die folgende:

   (71)

   Ihr Beweis startet, indem man auf die äußere Summe den „Anfänger-Traum“ anwendet und danach mit 𝑖 := 𝑖 − 1 umindiziert. Mittels Formel 68 wird im Anschluß 𝜆_𝑖−1² substituiert.

   

   Läßt man die Summe bei 𝑖 = 1 starten, dann ist zwar eine Korrektur um 𝛽²(𝜆₁ +𝛼) vonnöten, letztlich ändert sich wegen 𝜆₁ = 𝛼 aber überhaupt nichts (𝜆₁ +𝛼 = 0). Jetzt noch den Summanden für 𝑖 = 𝑛 herausziehen, ausmultiplizieren und 𝛽^2𝑛 = 𝛽 sowie 𝜆_𝑛 = tr(𝛼) entsprechend Formel 64 berücksichtigen:

   

   Im letzten Schritt wird noch der Summand für 𝑖 = 0 hinzugenommen: ∑ _𝑖=1^𝑛−1𝛽^2𝑖 = 𝛽+∑_𝑖=0^𝑛−1𝛽^2𝑖 = 𝛽+tr(𝛽) und man erkennt, daß der ganz rechte Term genau 𝛾 darstellt.

   

Algorithmus Die Idee für einen Algorithmus steckt in der Darstellung 64 für den Trace bzw. hinter Rekursionsformel 68.

4.5.3 Halb-Trace

Der Halb-Trace ist eine wichtige Funktion für den Fall, daß 𝑛 für den Binärkörper F_2^𝑛 ungerade ist.

Definition

Eigenschaften

1. Der Halb-Trace ist genauso linear wie der Trace.

   

2. htr(𝛼²) = htr ²(𝛼), es gilt jedoch hier: htr(𝛼²) ≠ htr(𝛼), wehalb htr(𝛼) kein Element aus F₂ (im Gegensatz zum Trace), sondern aus F_2^𝑛 ist.

   Einsetzen der Polynomdarstellung von 𝛼 führt in ähnlicher Art und Weise wie beim Trace zu:

   

3. htr(𝛼²) +htr(𝛼) = 𝛼 +tr(𝛼)

   Etwas ausführlicher

   

   Umindizierung der beiden Summen mit 𝑙 := 2𝑖 +1 bzw. 𝑙 := 2𝑖 führt unter zusätzlicher Berücksichtigung von 𝛼^2𝑛 = 𝛼 zum Ergebnis:

   

4. Aus den beiden vorangegangenen Punkten ergibt sich:

   (72)

Algorithmus Der Algorithmus arbeitet ähnlich wie beim Trace (vgl. Algorithmus 8), außer daß der zusätzlichen Potenzierung Rechnung getragen wird.

4.5.4 Lösung

Ist 𝛽 = 0 dann sind die Lösungen 0 und 1.

Im allgemeinen gilt:

d. h. eine Voraussetzung für die Existenz der zwei Lösungen ist tr(𝛽) = 0.

Für den Spezialfall, daß 𝑛 ungerade ist, kann man mit der Grundvoraussetzung tr(𝛽) = 0 aus Formel 72 direkt eine Lösung der quadratischen Gleichung ablesen:⁶⁷

(73)

Mit htr ²(𝛽) +htr(𝛽) = [htr(𝛽) +1]htr(𝛽) sind die zwei Lösungen demzufolge 𝑧₁ = htr(𝛽) und 𝑧₂ = htr(𝛽) +1.

Ist 𝑛 jedoch gerade (bzw. im allgemeinen Fall), dann besteht eine Berechnungsmöglichkeit durch Anwendung von Formel 71. Setzt man nämlich als Grundvoraussetzung tr(𝛽) = 0, dann gilt:

d. h. findet man ein Element 𝛼 mit tr(𝛼) = 1, wodurch

dann kann man die Lösung 𝑧 aus diesem Element 𝛼 berechnen. Glücklicherweise besitzt im Mittel die Hälfte aller Elemente 𝛼 ∈ F_2^𝑛 einen Trace von 1, so daß man 𝛼 durchaus zufällig wählen kann (mit 50% Wahrscheinlichkeit für tr(𝛼) = 1).

Algorithmus Für den allgemeinen Fall (bzw. wenn 𝑛 gerade ist) kann man auf Algorithmus 10 zurückgreifen, um die quadratische Gleichung in F_2^𝑛 zu lösen [1, 13].

Beide Summen können vereint werden, wenn man 𝜆_𝑖 iterativ in der äußeren Summe mitberechnet (nach derselben Methodik wie in Algorithmus 8), was dann im letzten Durchlauf 𝜆_𝑛 = tr(𝛼) ergibt.

Und es gilt ja nach Formel 71, unter der Voraussetzung tr(𝛽) = 0 und tr(𝛼 = 1), immer:

Optimierung Schreibt man dazu unter Zuhilfenahme der Polynomdarstellung für das Element 𝛼 und noch unter Berücksichtigung von Formel 33:

Umsortieren der Summen und mit 𝛼_𝑗^𝑙 = 𝛼_𝑗 in F₂ und dann wieder mit Formel 33

Den Wert für tr(𝑥^𝑗) kann man vorausberechnen.

4.6 MONTGOMERY-Potenzierung

Bei der Methode nach MONTGOMERY handelt es sich eigentlich um eine Multiplikationsmethode, bei der die Modulo-Reduktion des (Zwischen-) Ergebnisses ohne echte Langzahldivision erfolgen kann [23], [7, S. 9.2.1], [11, S. 2.2.4]. Da allerdings zuerst eine Transformation beider Faktoren in den MONTGOMERY-“Raum“ vorgenommen werden muß (welche zwei Modulo-Division erfordert) und außerdem einmal der erweiterte GCD-Algorithmus bemüht werden muß, kommen die Geschwindigkeitsvorteile nur bei der Potenzierung wirklich zum tragen.

Um einen leicht verständlichen Zugang zu finden, konzentrieren wir uns auf eine einzelne Multiplikation 𝑐 = 𝑎𝑏 (mod 𝑚). Dazu soll eine Zahl 𝑟 vorausgesetzt werden, die keinen gemeinsamen Teiler mit dem Modul 𝑚 hat und für die 𝑟 > 𝑚 gewährleistet ist. Die MONTGOMERY-Multiplikation kann man nun, unter der Voraussetzung gcd(𝑟,𝑚) = 1 (bzw. mit dem Satz von BéZOUT 𝛼𝑟−𝛽𝑚 = 1), in folgenden Einzelschritten darstellen:

1. Berechnung der BéZOUT-Koeffizienten 𝛼,𝛽 mit Hilfe des erweiterten euklidischen Algorithmus;
2. Eingangstransformation der Faktoren 𝑎 und 𝑏 zu â = 𝑎𝑟 mod 𝑚 und = 𝑏𝑟 mod 𝑚 (eine normale Modulo-Division);

3. (Wiederholte) Multiplikation im MONTGOMERY-Bereich:

   1. Berechnung des Produkts 𝑥 = â · (man beachte, daß hierbei keine Modulo-Reduktion vorgenommen wird);

   2. MONTGOMERY-Reduktion von 𝑥 = (𝑎𝑟 mod 𝑚)·(𝑏𝑟 mod 𝑚) zu ĉ = 𝑎𝑏𝑟 mod 𝑚 = 𝑐𝑟 mod 𝑚;

      Die nötige Korrektur⁶⁸ ĉ = M(𝑥,𝑟,𝑚) = 𝑥𝑟⁻¹ mod 𝑚 führt zu einer MONTGOMERY-Darstellung für 𝑐 (als Voraussetzung für den nächsten Teilschritt).⁶⁹

   3. Eine weitere (optionale) Multiplikationen im MONTGOMERY-“Raum“, die ĉ als einen der Faktoren verwendet;
4. Rücktransformation des Ergebnisses ĉ zu 𝑐 = ĉ𝑟⁻¹ mod 𝑚, ebenfalls eine MONTGOMERY-Reduktion: 𝑐 = M(ĉ,𝑟,𝑚).

MONTGOMERY-Reduktion Zur Herleitung von MONTGOMERY’s effizienter Berechnungsmethode für M(𝑥,𝑟,𝑚) = 𝑥𝑟⁻¹ mod 𝑚 wählen wir als Ausgangspunkt:

Subtraktion beider Gleichungen ergibt

was mit gcd(𝑟,𝑚) = 𝛼𝑟 −𝛽𝑚 = 1 zu

führt.⁷⁰

Mit dem Wissen, daß es sich bei 𝛼 um das multiplikativ inverse Element von 𝑟 im Restklassensystem modulo 𝑚 handelt (𝛼 = 𝑟⁻¹ (mod 𝑚), 𝛽 = −𝑚⁻¹ (mod 𝑟) ), stellen wir noch um:

und gewinnen letztlich MONTGOMERY’s Reduktionsformel.

(74)

Praktisch benötigt man 𝑢 und 𝑣 nicht, sondern geht meist nach folgendem Algorithmus vor:

Ensure:   𝑦 = 𝑥𝑟⁻¹ mod 𝑚
  𝑡 ⇐ 𝛽𝑥 mod 𝑟
  𝑦 ⇐
  if 𝑦 ≥ 𝑚 then
   𝑦 ⇐ 𝑦−𝑚
  end if

Bei geschickter Wahl von 𝑟 zum Beispiel als 𝑟 = 2^𝑠 sind für alle (Modulo-) Divisionen in Formel 74 nur logische oder Schiebeoperationen nötig.⁷¹ Um gcd(2^𝑠,𝑚) = 1 zu garantieren ist die einfachste Bedingung die, 𝑚 als ungerade vorauszusetzen.⁷²