Algorithmensammlung: Zahlentheorie: Euklidischer Algorithmus

Euklidischer Algorithmus

Mit dem Euklidischen Algorithmus kann man den größten gemeinsamen Teiler zweier Zahlen herausfinden. Für eine genaue Beschreibung siehe Euklidischer Algorithmus.

Laufzeit

Man kann zeigen, dass der moderne, iterative Euklidische Algorithmus (Division und Modulo-Funktion statt Subtraktion) aus den ganzen Zahlen a und b nach spätestens n Schritten^[1]

\textstyle n\leq 2\cdot \log _{2}(\max {(a,b)})

den größten gemeinsamen Teiler $\textstyle ggT(a,b)$ errechnet hat.

Die längste Laufzeit für ein Zahlenpaar ergibt sich, wenn es sich um zwei aufeinander folgende Zahlen aus der Fibonacci-Folge handelt, damit lässt sich die Laufzeitschätzung weiter nach unten abschätzen.

Implementierungen

Basic X 11

 rem Groesster gemeinsamer Teiler zweier ganzer Zahlen a, b
 print "Dieses Programm berechnet den groessten gemeinsamen Teiler zweier Zahlen."
 input a
 input b
 rem absolute Beträge bilden, weil ggT auch für negative Zahlen definiert ist
 rem ggT(a, b) = ggT(|a|, |b|)
 a = Abs(a)
 b = Abs(b)
 While (a > 0) And (b > 0)
 rem Weiter, solange a·b > 0
 rem Berechnung für den allgemeinen Fall a·b <> 0
        If a > b
          rem Der größere Wert wird überschrieben
            a = a Mod b
        Else
            b = b Mod a
        Endif
 Wend
 ggT = a + b
 print ggT

C

a und b sind die beiden Zahlen von denen man den größten Teiler haben möchte.

int Euklid(int a, int b)
{
    if (a == 0)                          //Wenn a=0 ist b der größte gemeinsame Teiler laut Definition
    {
    	return b;
    }
    else {
    while(b != 0)                        //So lange wiederholen, wie b nicht 0 ist.
    {
    	if (a > b)
    	{
    		a = a - b;               //Wenn a größer als b, subtrahiere b von a.
    	}
        else
    	{
    		b = b - a;               //In jedem anderen Fall subtrahiere a von b.
    	}
    }
    return a;                            //In a steht jetzt der größte gemeinsame Teiler von a und b.
}
}

C#

Klassische Variante

Rekursiv

public int euklid(int a, int b)
{
	if (b == 0) {
		return a;
	} else if (a == 0) {
		return b;
	} else if (a > b) {
		return euklid(a - b, b);
	} else {
		return euklid(a, b - a);
	}
}

Iterativ

public int euklid(int a, int b)
{
	if (a == 0) {
		return b;
	} else {
		while (b != 0) {
			if (a > b) {
				a -= b;
			} else {
				b -= a;
			}
		}
	}
	return a;
}

Moderne Variante

Rekursiv

public int euklid(int a, int b)
{
	if (b == 0)
		return a;
	else
		return euklid(b, a % b);
	
}

Iterativ

public int euklid(int a, int b)
{
	int tmp = 0;
	while (b != 0) {
		tmp = a % b;
		a = b;
		b = tmp;
	}
	return a;
}

FreePascal

klassischer Algorithmus

Folgender Quelltext kombiniert Pascal und X86-64-assembler:

type
	naturalNumber = 1..high(qword);

function gcdEuclid(n, m: naturalNumber): naturalNumber;
{$ifdef CPUx86_64} // - - - - - - - - - optimized implementation
assembler; register;
{$push}
{$goto on}
label
	gcdEuclid_iterate, gcdEuclid_subtract, gcdEuclid_condition;
{$asmMode intel}
asm
	// use cmp and jump instructions in loop condition block
	jmp gcdEuclid_condition // goto condition
	
gcdEuclid_iterate:
	// definition: n always contains the larger number,
	// so the sub-instruction operates on the same regs
	jae gcdEuclid_subtract // if n >= m then goto subtract
	xchg n, m // swap values in n and m
	
gcdEuclid_subtract:
	// sub modifies registers => put it here inbetween
	sub n, m // n := n - m
	
gcdEuclid_condition:
	// iterate, while n <> m
	cmp n, m // n = m
	jne gcdEuclid_iterate // if n <> m then goto iterate
	
	mov rax, n // result := n
end;
{$pop}
{$else} // - - - - - - - - - - - - - - -  default implementation
var
	x: qword;
begin
	while n <> m do
	begin
		if n < m then
		begin
			x := n;
			n := m;
			m := x;
		end;
		n := n - m;
	end;
	gcdEuclid := n;
end;
{$endif}

Der FreePascal Compiler generiert mit -O3 (höchste Optimierungsstufe) bereits sehr schnellen und kompakten Code. Lediglich n <> m und n < m, welches durch zwei cmp-Instruktion umgesetzt wird, läßt sich – wie oben im $ifdef-Block geschehen – in eine einzige cmp-Instruktion auflösen, da die relevanten Flaggen sich zwischenzeitlich nicht ändern. Zudem kann man den Dreieckstausch mit xchg schreiben. Drei konsekutive mov werden aber genauso gut von der Datenfluß-analyse-Einheit erkannt und entsprechend behandelt.

Es sei angemerkt, daß obige Funktion zur Laufzeit auch 0 akzeptiert, aber alleine nicht abfängt. Dies wurde weggelassen, da dies nicht zum eigentlichen Euklidischen Algorithmus zugehört.

Die Schnelligkeit läßt sich (ungefähr) mit dem Werkzeug time(1) ermitteln:

$ echo '4294967291\n9223372036854775783' | time ./greatestCommonDivisor

Hier wurde mit Primzahlen getestet, um zum einen (auszugsweise) die Richtigkeit des Algorithmuses bestätigt zu sehen (Ergebnis muß 1 sein), und zum anderen möglichst viele Iteration zu erreichen (meßbare Zeitdauern erreichen), da schließlich bis 1 runtersubtrahiert werden muß.

Java

/**
	gcd: greatest common divisor
	according to modern euclidean algorithm
*/
public static long gcdEuclid (long x, long y) {
	while ( y != 0 ) {
		long z = x % y;
		x = y;
		y = z;
	}
	return x;
}

Julia

Moderne Variante

Rekursiv

function euklid(a, b)
    if b == 0
        return a
    else
        return euklid(b, mod(a, b))
    end
end

Iterativ

function euklid(a, b)
    while b != 0
        a, b = b, mod(a,b)
    end
    return a
end

Lua

Klassische Variante

Rekursiv

function euklid(a, b)
    if b == 0 then
        return a
    elseif a == 0 then
        return b
    elseif a > b then
        return euklid((a - b),b)
    else
        return euklid(a,(b - a))
    end
end

Iterativ

function euklid(a, b)
    if a == 0 then
        return b
    else
        while b ~= 0 do
            if a > b then
                a = (a - b)
            else
                b = (b - a)
            end
        end
    end
    
    return a
end

PHP

Klassische Variante

Rekursiv

function euklid($a, $b) {
    if($b == 0) {
        return $a;
    } elseif($a == 0) {
        return $b;
    } elseif($a > $b) {
        return euklid(($a - $b),$b);
    } else {
        return euklid($a,($b - $a));
    }
}

Iterativ

function euklid($a, $b) {
    if($a == 0) {
        return $b;
    } else {
        while($b != 0) {
            if($a > $b) {
                $a -= $b;
            } else {
                $b -= $a;
            }
        }
    }

    return $a;
}

Moderne Variante

Rekursiv

function euklid($a, $b) {
    if($b == 0) {
        return $a;
    } else {
        return euklid($b,($a % $b));
    }
}

Iterativ

function euklid($a, $b) {
    while ($b != 0) {
        $r = $a % $b;
        $a = $b;
        $b = $r;
    }
    return $a;
}

Steinscher Algorithmus

Die Ausnutzung einiger Rechenregeln führt zum steinschen Algorithmus:

function stein($a, $b) {
    $k = 0;

    while(!($a & 1) && !($b & 1)) {
        $a >>= 1;
        $b >>= 1;

        ++$k;
    }

    if($a & 1) {
        $t = -$b;
    } else {
        $t = $a;
    }

    while($t != 0) {
        while(!($t & 1)) {
            $t >>= 1;
        }

        if($t > 0) {
            $a = $t;
        } else {
            $b = -$t;
        }

        $t = ($a - $b);
    }

    return ($a << $k);
}

Python

Der Algorithmus in moderner Form:

# Getestet unter Python 3.4, sollte aber unter Python Version 2 und 3 laufen

def euklid(a, b):
    while b != 0:
        a, b = b, a % b
    return a

Ruby

Klassische Variante

Rekursiv

def euklid(a, b)
  if b == 0
    a
  elsif a == 0
    b
  elsif a > b
    euklid (a - b), b
  else
    euklid a, (b - a)
  end
end

Visual Basic .NET

Klassische Variante

Rekursiv

    Function euklid(ByVal a As Integer, ByVal b As Integer) As Integer
        If b = 0 Then
            Return a
        ElseIf a = 0 Then
            Return b
        ElseIf a > b Then
            Return euklid(a - b, b)
        Else
            Return euklid(a, b - a)
        End If
    End Function

Iterativ

 Function euklid(ByVal a As Integer, ByVal b As Integer) As Integer
     If a = 0 Then
            Return b
        Else
            While b <> 0
                If a > b Then
                    a -= b
                Else
                    b -= a
                End If
            End While
        End If
        Return a
    End Function

Moderne Variante

Rekursiv

Function euklid(ByVal a As Integer, ByVal b As Integer) As Integer
   Return If(b = 0, a, euklid(b, a Mod b))
End Function

Iterativ

 Function euklid(ByVal a As Integer, ByVal b As Integer) As Integer
        Dim tmp As Integer = 0
        While b <> 0
            tmp = a Mod b
            a = b
            b = tmp
        End While
        Return a
 End Function

Visual Basic for Applications

Allgemeine Hinweise zu VBA-Code:

Wenn ein Argument einer Funktion als Datentyp Long deklariert, aber Dezimalzahlen übergeben wurden, dann Rundet sie VBA bei der Übergabe an die Funktion auf ganzzahlige Werte (impliziter Typecast durch VBA, dabei wird nach IEEE-754 ( mathematische Rundung) gerundet. Dadurch erhält die Funktion beispielsweise bei Übergabe von x = 0,5 den abgerundeten Wert x = 0.
Die Operation \ ist in VBA die ganzzahlige Division (div).
Die Bereichsgrenzen von Datentypen beruhen auf der 32-Bit-VBA Umgebung.

Größter gemeinsamer Teiler (ggT)

Weil in Excel ("ziehen" von Formeln in Tabellen) und Access (VBA-Funktionen in Abfragen) je nach Aufgabe sehr viele Berechnungen erforderlich sein können, sollten nur performante und numerisch stabile Algorithmen verwendet werden. Daher wird hier nur die iterativ-moderne Variante vorgestellt:

' *** Größter gemeinsamer Teiler zweier ganzer Zahlen a, b
Public Function ggT(ByVal a As Long, ByVal b As Long) As Long
    ' Beträge bilden, weil ggT auch für negative Zahlen
    ' ggT(a, b) = ggT(|a|, |b|) definiert ist:
    a = Abs(a)
    b = Abs(b)
    
    While (a > 0) And (b > 0) ' Weiter, solange a·b > 0
        '*** Berechnung für den allgemeinen Fall a·b <> 0
        If a > b Then ' Der größere Wert wird überschrieben
            a = a Mod b
        Else
            b = b Mod a
        End If
    Wend
    
    ggT = a + b ' Ein Summand ist sicher 0
End Function

Hinweise zum Algorithmus: Der Algorithmus kann keine Zwischen- oder Endergebnisse erzeugen, die zum Überlauf oder Fehler führen. Die logischen Testbedingungen in der While-Anweisung sollten nicht durch a·b ersetzt werden, da sonst ein Überlauf bei zu großem a·b auftreten könnte. Die Argumente a, b müssen ganzzahlig sein (sonst werden sie von VBA gerundet) und im Datentyp Long darstellbar sein (2.147.483.647 ≥ max(a, b) ≥ min(a, b) ≥ -2.147.483.648), die Reihenfolge der Eingabe ist egal, es muss nicht nach Größe sortiert sein.

Ergebnis: Das Ergebnis liegt stets zwischen $\textstyle \min(|a|,|b|)\geq ggT(a,b)\geq 0$ .

Kleinstes gemeinsames Vielfaches (kgV)

Mit dem dem Euklidischen Algorithmus eng verwandt ist das kleinste gemeinsame Vielfache (kgV), der Zusammenhang ergibt sich aus der Gleichung $a\cdot b=ggT(a,b)\cdot kgV(a,b)$ .

Public Function kgV(a As Long, b As Long) As Variant
    
    ' kgV ist für a·b = 0 nicht definiert:
    If (a = 0) Or (b = 0) Then Exit Function
    
    ' Typumwandlung in Double verhindern
    kgV = a \ ggT(a, b) ' oder: kgV = CLng(a / ggT(a, b))
    
    kgV = Abs(kgV * b) ' Eventuell Typumwandlung in Double

End Function

Hinweise zum Algorithmus: Da es möglich ist, dass das $\textstyle kgV(a,b)\geq \max(a,b)$ ist, könnte das Programm einen Überlauf produzieren, wenn der Rückgabewert vom Datentyp Long wäre, daher wurde für die Ausgabe der Datentyp Variant gewählt. Die Division wird als erste Operation ausgeführt, um mit möglichst kleinen Werten zu rechnen, andererseits bewirkt sie unabhängig vom Ergebnis eine automatische Typumwandlung zu Double, obwohl das Ergebnis ganzzahlig sein muss. Dies wird entweder durch Verwendung von \ oder einen Typecast mit CLng() verhindert. Durch den Aufbau der Funktion wird im zweiten Schritt das Ergebnis von kgV durch VBA nur dann automatisch in den Typ Double umgewandelt, wenn es nicht mehr als Typ Long darstellbar ist.

Der Rückgabewert ist kgV(a, b) = 0 falls a·b = 0, sonst ist kgV(a, b) ≥ max(a, b)

Einzelnachweise

↑ Ziegenbalg, "Algorithmen", 2. Auflage 2007, Harri Deutsch Verlag Frankfurt/M, ISBN 978-3-8171-1814-4, siehe Kap. 5.5 mit Verweis auf Donald E. Knuths "Art of Computer Programming" Volume 2, Abschnitt 4.5.3

[1] Ziegenbalg, "Algorithmen", 2. Auflage 2007, Harri Deutsch Verlag Frankfurt/M, ISBN 978-3-8171-1814-4, siehe Kap. 5.5 mit Verweis auf Donald E. Knuths "Art of Computer Programming" Volume 2, Abschnitt 4.5.3

[1]