Algorithmensammlung: Zahlentheorie: Euklidischer Algorithmus

Algorithmensammlung: Zahlentheorie

Euklidischer Algorithmus

Bearbeiten

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

Laufzeit

Bearbeiten

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]

 

den größten gemeinsamen Teiler   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

Bearbeiten

Basic X 11

Bearbeiten
 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
  • 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.
}
}

Klassische Variante

Bearbeiten
Rekursiv
Bearbeiten
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
Bearbeiten
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

Bearbeiten
Rekursiv
Bearbeiten
public int euklid(int a, int b)
{
	if (b == 0)
		return a;
	else
		return euklid(b, a % b);
	
}
Iterativ
Bearbeiten
public int euklid(int a, int b)
{
	int tmp = 0;
	while (b != 0) {
		tmp = a % b;
		a = b;
		b = tmp;
	}
	return a;
}

FreePascal

Bearbeiten

klassischer Algorithmus

Bearbeiten

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ß.

/**
	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;
}

Moderne Variante

Rekursiv

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

Iterativ

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

Klassische Variante

Bearbeiten
Rekursiv
Bearbeiten
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
Bearbeiten
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

Klassische Variante

Bearbeiten
Rekursiv
Bearbeiten
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
Bearbeiten
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

Bearbeiten
Rekursiv
Bearbeiten
function euklid($a, $b) {
    if($b == 0) {
        return $a;
    } else {
        return euklid($b,($a % $b));
    }
}
Iterativ
Bearbeiten
function euklid($a, $b) {
    while ($b != 0) {
        $r = $a % $b;
        $a = $b;
        $b = $r;
    }
    return $a;
}

Steinscher Algorithmus

Bearbeiten

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);
}

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

Klassische Variante

Bearbeiten
Rekursiv
Bearbeiten
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

Bearbeiten

Klassische Variante

Bearbeiten
Rekursiv
Bearbeiten
    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
Bearbeiten
 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

Bearbeiten
Rekursiv
Bearbeiten
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
Bearbeiten
 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

Bearbeiten

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)

Bearbeiten

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  .

Kleinstes gemeinsames Vielfaches (kgV)

Bearbeiten

Mit dem dem Euklidischen Algorithmus eng verwandt ist das   kleinste gemeinsame Vielfache (kgV), der Zusammenhang ergibt sich aus der Gleichung  .

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   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

Bearbeiten
  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