mathguide.py

# -*- coding: UTF-8 -*-

# mathGUIde Version 2.1"
# Copyright © 2004-2011 Hartmut Ring"
# http://www.math.uni-siegen.de/ring/mathGUIde/
# This file is for internal use only by the mathGUIde GUI.
#
# mathGUIde is free software; you can redistribute it and/or modify it under the
# terms of the GNU General Public License as published by the Free Software
# Foundation; either version 3 of the License, or (at your option) any later
# version.
# This program is distributed in the hope that it will be useful, but WITHOUT
# ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
# FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
# You should have received a copy of the GNU General Public License along with
# this program; if not, see http://www.gnu.org/licenses/.""")
#............................................................

import operator, random, math, re, inspect, functools

class _MathGUIdeError(Exception):
 def __init__(self, value):
 self.value = value
 def __str__(self):
 return repr(self.value)

_pow = pow # import everything from math except pow:
from math import *
pow = _pow

# make original Python functions floor and ceil invisible:
_floor, _ceil = floor, ceil

def _cmp(a,b):
 return int(a > b) - int(a < b)

#............................................................
#begin

# Change floor and ceil so that they map from R-->Z
# instead of R-->R as given by Python.

def floor(x):
 """ floor(x) 
 highest integer n <= x 
 ----de----
 floor(x) 
 größte ganze Zahl n <= x (Gauß-Klammer) 
 """
 return int(_floor(x))

def ceil(x):
 """ floor(x) 
 smallest integer n >= x 
 ----de----
 ceil(x) 
 kleinste ganze Zahl n >= x 
 """
 return int(_ceil(x))

#............................................................

def copy(x):
 """ copy(x) 
 independent copy ("clone") of x
 ----de----
 copy(x) 
 unabhängige Kopie
 """
 try:
 c = x.copy()
 except:
 c = x
 return c

#===============================================================

def fromTo(a, b, step=1):
 """ fromTo(a, b, step=1) 
 Iterator over equidistant integers from <var>a</var> to <var>b</var> 
 If no <var>step</var> is given, <var>step</var> is set to 1. 
 Examples:<pre>
 list(fromTo(3, 6)) # ---> [3,4,5,6]
 list(fromTo(3, 6, 2)) # ---> [3,5]
 list(fromTo(6, 3)) # ---> []
 list(fromTo(6, 3, -2)) # ---> [6,4]
 [2*i for i in fromTo(6, 3, -2)] # ---> [12,8]
 sum(fromTo(1, 10)) # ---> 55</pre>
 ----de----
 fromTo(a, b, step=1) 
 Iterator über die ganzen Zahlen von a bis einschl. b im Abstand step 
 Wenn <var>step</var> nicht übergeben wird, wird <var>step</var> auf 1 gesetzt. 
 Beispiele:<pre>
 list(fromTo(3, 6)) # ---> [3,4,5,6]
 list(fromTo(3, 6, 2)) # ---> [3,5]
 list(fromTo(6, 3)) # ---> []
 list(fromTo(6, 3, -2)) # ---> [6,4]
 [2*i for i in fromTo(6, 3, -2)] # ---> [12,8]
 sum(fromTo(1, 10)) # ---> 55</pre>
 """
 if step > 0:
 return range(a, b+1, step)
 else:
 return range(a, b-1, step)

#===============================================================
# Elementary Functions
#---------------------------------------------------------------

def isInteger(n):
 """ isInteger(n) 
 returns <code>True</code> if and only if n is of type <code>int</code>.
 ----de----
 isInteger(n) 
 gibt <code>True</code> zurück genau dann, wenn n vom Typ <code>int</code> ist.
 """
 return type(n) == int or type(n) == Rational and n.q == 1

def isNatural(n):
 return isInteger(n) and n > 0

#............................................................

def sum(obj, start=None, end=None, step=1):
 """ sum(obj, start=None, end=None, step=1) 
 The sum of the elements of iterable obj or of all obj(i) for i in fromTo(start,end,step) 
 The element type must define the operator + (e.g. numbers or strings). 
 Examples: 
 <code>sum([i^2 for i n fromTo(1,10)])</code> 
 <code>sum(fibonacci, 0, 10)</code> 
 <code>sum(lambda i:i^2, 2, 10, 2)</code>
 ----de----
 sum(obj, start=None, end=None, step=1) 
 Summe der Elemente des iterierbaren Objekts a bzw. aller obj(i) für i in fromTo(start,end,step) 
 Die Elemente müssen mit dem Operator + verknüpfbar sein
 (z. B. Zahlen oder Strings). 
 Beispiele: 
 <code>sum([i^2 for i n fromTo(1,10)])</code> 
 <code>sum(fibonacci, 0, 10)</code> 
 <code>sum(lambda i:i^2, 2, 10, 2)</code>
 """
 s = 0 # for an empty iteration the element type ist assumed to be numeric.
 if str(type(obj)) != "<class 'function'>":
 # obj is sequence or iterator or generator
 for i in obj:
 if s == 0: # adding to 0 would not work for types with neutral element != 0 (e.g. for a list of strings).
 s = i
 else:
 s += i
 else: # obj is a function
 for i in fromTo(start, end, step):
 if s == None:
 s = obj(i)
 else:
 s += obj(i)
 return s

#............................................................

def product(obj, start=None, end=None, step=1):
 """ product(obj, start=None, end=None, step=1) 
 The product of the elements of iterable obj or of all obj(i) for i in fromTo(start,end,step) 
 The element type must define the operator *. 
 Examples: 
 <code>product([i^2 for i n fromTo(1,10)])</code> 
 <code>product(fibonacci, 0, 10)</code> 
 <code>product(lambda i:i^2, 2, 10, 2)</code>
 ----de----
 product(obj, start=None, end=None, step=1) 
 Product der Elemente des iterierbaren Objekts a bzw. aller obj(i) für i in fromTo(start,end,step) 
 Die Elemente müssen mit dem Operator * verknüpfbar sein 
 Beispiele: 
 <code>product([i^2 for i n fromTo(1,10)])</code> 
 <code>product(fibonacci, 0, 10)</code> 
 <code>product(lambda i:i^2, 2, 10, 2)</code>
 """
 s = 1 # for an empty iteration the element type ist assumed to be numeric.
 if str(type(obj)) != "<class 'function'>":
 # obj is sequence or iterator or generator
 for i in obj:
 if s == 1: # multiplying to 1 would not work for types with neutral element != 1
 s = i
 else:
 s *= i
 else: # obj is a function
 for i in fromTo(start, end, step):
 if s == None:
 s = obj(i)
 else:
 s *= obj(i)
 return s

#............................................................

def closure(elements, operation):
 """ closure(elements, operation) 
 The closure of the elements under the binary operation 
 operation must be given as a function with two parameters.
 ----de----
 closure(elements, operation) 
 Der Abschluss der Menge elements unter der binären Operation operation 
 operation muss eine function mit zwei Parametern sein.
 """
 closureList = list(elements)
 stringList = [str(x) for x in closureList]
 more = True
 while more:
 more = False
 for x in closureList:
 for y in closureList:
 e = operation(x,y)
 s = str(e)
 if s not in stringList:
 closureList.append(e)
 stringList.append(s)
 more = True
 closureList.sort(key=str)
 return closureList

#............................................................

def join(l, filler=""):
 """ join(l, filler="") 
 Joins the elements of <code>l</code> to a string 
 Between each two adjacent elements <code>filler</code> is inserted. 
 This is a replacement for the string method <code>join()</code>. 
 In addition to the string method the list elements are converted to strings. 
 Example: 
 <code>join(fromTo(1,5), "-")</code>
 returns <code>"1-2-3-4-5"</code>.
 ----de----
 join(l, filler="") 
 verbindet die Elemente von <code>l</code> zu einem String 
 Zwischen je zwei Elemente wird <code>filler</code> eingefügt. 
 Ersatz für die schwer lesbare String-Methode join()
 bzw. für die Funktion string.join(). 
 Im Gegensatz zu diesen Funktionen werden zusätzlich
 die Listenelemente zu Strings konvertiert. 
 Beispiel: 
 <code>join(fromTo(1,5), "-")</code>
 liefert <code>"1-2-3-4-5"</code>.
 """
 return filler.join([str(x) for x in l])

#===============================================================
# Algorithms of Elementary Number Theory
#---------------------------------------------------------------

#---------------------------------------------------------------
# Large random numbers

def rand(n):
 """ rand(n) 
 Integer random number from the range 0 ... <var>n</var>-1 
 <var>n</var> must be a natural number (int).
 ----de----
 rand(n) 
 Zufallszahl im Bereich 0 ... <var>n</var>-1 
 <var>n</var> kann eine beliebige natürliche Zahl sein (int).
 """
 return random.randrange(n)

#---------------------------------------------------------------
# Conversion between positional notation systems

class NumeralSystem: # conversions between numeral systems # Umrechnung zwischen Stellenwertsystemen
 """ Class <code>NumeralSystem</code>: Numeral Systems 
 This is a static class (only static methods available).
 ----de----
 Klasse <code>NumeralSystem</code>: Stellenwertsysteme 
 Statische Klasse: nur Klassenmethoden.
 """
 #............................................................
 _staticMethods = ("toBase", "binary", "fromBase")

 @staticmethod
 def toBase(n, b, ascending=False):
 """ toBase(n, b) 
 <var>n</var> in the base <var>b</var> numeral system 
 The result is returned as a list of the digits. 
 See also <a href="#binary"><code>binary</code></a>,
 <a href="#fromBase"><code>fromBase</code></a>.
 ----de----
 toBase(n, b) 
 <var>n</var> im Stellenwertsystem zur Basis <var>b</var> 
 Das Ergebnis wird als Liste der Stellen ausgegeben. 
 Siehe auch <a href="#binary"><code>binary</code></a>,
 <a href="#fromBase"><code>fromBase</code></a>.
 """
 v = []
 while n > 0:
 if ascending:
 v.append(n % b) # insert n % b at the end
 else:
 v.insert(0, n % b) # insert n % b at the beginning
 n //= b
 return v

 @staticmethod
 def binary(n, ascending=False):
 """ binary(n) 
 Binary number representation of <var>n</var> as list 
 Abbreviation for <code>toBase(n,2)</code> 
 Example: 
 <code>NumeralSystem.binary(11)</code> returns <code>[1, 0, 1, 1]</code>. 
 See also <a href="#toBase"><code>toBase</code></a>,
 <a href="#fromBase"><code>fromBase</code></a>.
 ----de----
 binary(n) 
 Binärdarstellung von <var>n</var> als Liste 
 Kurzform für <code>toBase(n,2)</code> 
 Beispiel: 
 <code>NumeralSystem.binary(11)</code> gibt <code>[1, 0, 1, 1]</code> zurück. 
 Siehe auch <a href="#toBase"><code>toBase</code></a>,
 <a href="#fromBase"><code>fromBase</code></a>.
 """
 return NumeralSystem.toBase(n, 2, ascending)

 @staticmethod
 def fromBase(v, b, ascending=False):
 """ fromBase(v, b) 
 Converts <var>v</var> from the base <var>b</var> numeral system to a number 
 <var>v</var> is expected as list of digits. 
 See also <a href="#toBase"><code>toBase</code></a>,
 <a href="#binary"><code>binary</code></a>.
 ----de----
 fromBase(v, b) 
 <var>v</var> aus dem Stellenwertsystem zur Basis <var>b</var> als Zahl 
 <var>v</var> wird als Liste der Stellen erwartet. 
 Siehe auch <a href="#toBase"><code>toBase</code></a>,
 <a href="#binary"><code>binary</code></a>.
 """
 n = 0
 if ascending:
 for i in fromTo(len(v)-1, 0, -1):
 n = b * n + i
 else:
 for i in v:
 n = b * n + i
 return n

#---------------------------------------------------------------
# Integer functions

def binomial(n, k):
 """ binomial(n, k) 
 Binomial coefficient 
 n must be a nonnegative integer and k in the range 0 .. n.
 ----de----
 binomial(n, k) 
 Binomialkoeffizient (n über k) 
 n muss nichtnegative ganze Zahl sein, k im Bereich 0 .. n.
 """
 if k > n//2:
 k = n - k
 b = 1
 for i in fromTo(1, k):
 b = b * (n+1-i) // i
 return b

#............................................................

def factorial(n):
 """ factorial(n) 
 n! = 1*2*3* ... *n (factorial of n)
 ----de----
 factorial(n) 
 n! = 1*2*3* ... *n (n Fakultät)
 """
 f = 1
 for i in fromTo(2, n):
 f *= i
 return f

#............................................................

def fibonacci(n):
 """ fibonacci(n) 
 <var>n</var>-th element of the Fibonacci series
 fibonacci(0) = 0 
 fibonacci(1) = 1 
 fibonacci(n) = fibonacci(n-1) + fibonacci(n-2) für n >= 2
 ----de----
 fibonacci(n) 
 <var>n</var>-tes Element der Fibonacci-Folge
 fibonacci(0) = 0 
 fibonacci(1) = 1 
 fibonacci(n) = fibonacci(n-1) + fibonacci(n-2) für n >= 2
 """
 a, b = 0, 1
 while n > 0:
 a, b = b, a+b
 n -= 1
 return a

#---------------------------------------------------------------
# Elementary number theoretic functions
# Funktionen der elementaren Zahlentheorie

def gcd(a,b):
 """ gcd(a,b) 
 greatest common divisor of a and b 
 Calculated using the Euclidean algorithm 
 See also <a href="#lcm"><code>lcm</code></a>
 ----de----
 gcd(a,b) 
 Größter gemeinsamer Teiler von a und b 
 Berechnung mit dem euklidischen Algorithmus 
 Siehe auch <a href="#lcm"><code>lcm</code></a>
 """
 while b != 0:
 a, b = b, a % b
 return a

#............................................................

def _doc_gcdBinary(a,b):
 """ Binary Euclidean algorithm
 For optimization replace:
 x //= 2 by x >>= 1
 x *= 2 by x <<= 1
 x % 2 == 0 by x & 1 == 0
 ----de----
 Binärer Euklidischer Algorithmus
 Zur Optimierung sollte ersetzt werden:
 x //= 2 durch x >>= 1
 x *= 2 durch x <<= 1
 x % 2 == 0 durch x & 1 == 0
 """
 assert a >= b
 g = 1
 while a % 2 == 0 and b % 2 == 0:
 # größte gemeinsame Zweierpotenz:
 a //= 2
 b //= 2
 g *= 2
 while a != 0:
 while a % 2 == 0:
 a //= 2
 while b % 2 == 0:
 b //= 2
 t = abs(a-b) // 2
 if a >= b:
 a = t
 else:
 b = t
 return g * b

#............................................................

def _doc_gcdRecursive(a,b):
 """ recursive version of gcd
 ----de----
 rekursive Variante von gcd
 """
 if b == 0:
 return a
 else:
 return _doc_gcdRecursive(b, a % b)

#............................................................

def gcdExt(a,b):
 """ gcdExt(a,b) 
 Extended Euclidean algorithm. 
 Returns (d, x, y) with d == gcd(a,b) and x*a + y*b == d
 ----de----
 gcdExt(a,b) 
 Erweiterter Euklidischer Algorithmus. 
 Gibt (d, x, y) zurück mit d == gcd(a,b) und x*a + y*b == d
 """
 if b == 0: # Basis for Math. Induction over b
 return a, 1, 0 # gcd(a,b) == 1*a + 0*b
 else:
 q, r = divmod(a, b) # r < b
 d, z, x = gcdExt(b, r) # d == z*b + x*r (Induction assumption)
 # == z*b + x*(a - b*q)
 # == x*a + (z - x*q)*b
 y = z - x * q # d == x*a + y*b
 return d, x, y

#............................................................

def lcm(a, b):
 """ lcm(a,b) 
 Least common multiple of <var>a</var> and <var>b</var> 
 See also <a href="#gcd"><code>gcd</code></a>
 ----de----
 lcm(a,b) 
 Kleinstes gemeinsames Vielfaches von <var>a</var> und <var>b</var> 
 Siehe auch <a href="#gcd"><code>gcd</code></a>
 """
 assert a > 0 and b > 0
 return (a * b) // gcd(a, b)

#............................................................

def _chinese2(mod1, mod2):
 m, n = mod1.m, mod2.m
 a, b = mod1.n, mod2.n
 d, x, y = gcdExt(m,n)
 assert d == 1 # m*x + n*y == 1
 x1 = x % n # m*x1 % n == 1
 y1 = y % m # n*y1 % m == 1
 c = n * y1 # c % m == 1
 # c % n == 0
 d = m * x1 # d % n == 1
 # d % m == 0
 # (a*c + b*d) % m == a
 # (a*c + b*d) % n == b
 return Mod(a*c + b*d, m*n)

def chinese(mod1, *mod2):
 """ chinese(mod1, mod2,...) 
 Mod object from mod1, mod2, ... according to the Chinese remainder theorem 
 Requires that mod1, mod2, ... are Mod objects with pairwise coprime moduli. 
 <code>chinese(Mod(a1, m1), Mod(a2, m2), ...)</code> 
 returns Mod(x, m1*m2*...) with 0 <= <code>x</code> < m1*m2*... and 
 <code>x % m1 == a1</code> 
 <code>x % m2 == a2</code> , ... 
 Example:<pre>
 chinese(Mod(1,5), Mod(2,7)) # ---> Mod(16,35)</pre>
 (16 mod 5 = 1, 16 mod 7 = 2). 
 See also <a href="Mod.html"><code>Mod</code></a>
 ----de----
 chinese(mod1, mod2,...) 
 Mod-Objekt zu mod1 mod2, ... nach Chinesischem Restsatz 
 Voraussetzung: mod1, mod2, ... sind Mod-Objekte mit paarweise
 teilerfremden Moduln. 
 <code>chinese(Mod(a1, m1), Mod(a2, m2, ...)</code> 
 berechnet Mod(x, m1*m2*...) mit 0 <= <code>x</code> < m1*m2*... und 
 <code>x % m1 == a1</code> 
 <code>x % m2 == a2</code> , ... 
 Beispiel:<pre>
 chinese(Mod(1,5), Mod(2,7)) # ---> Mod(16,35)</pre>
 (16 mod 5 = 1, 16 mod 7 = 2). 
 Siehe auch <a href="Mod.html"><code>Mod</code></a>
 """
 r = mod1.copy()
 for m in mod2:
 r = _chinese2(r, m)
 return r

#............................................................

def eulerPhi(n):
 """ eulerPhi(n) 
 Euler's totient function 
 Number of positive integers less than or equal to <var>n</var> that are coprime to <var>n</var>
 ----de----
 eulerPhi(n) 
 Eulersche Phi-Funktion 
 Anzahl der zu <var>n</var> teilerfremden Zahlen im Bereich von 1 bis <var>n</var>
 """
 #--- directly following the definition:
 # return sum([int(gcd(i,n) == 1) for i in Range(n)])

 #--- more efficient:
 for p in factors(n):
 n = n*(p-1)//p
 return n

#............................................................

def expMod(a, x, m):
 """ expMod(a, x, m) 
 Calculates (a**x) % m efficiently. 
 Equivalent to the Python standard function <code>pow</code>.
 ----de----
 expMod(a, x, m) 
 Berechnet (a**x) % m effizient. 
 Synonym für die Python-Standardfunktion <code>pow</code>.
 """
 return pow(a, x, m) # Python standard function pow

#............................................................
# The following two implementations are for documentation only:
# Die folgenden beiden Implementierungen dienen nur der Dokumentation:

def _doc_expMod_recursive(a, x, n):
 """ recursive implementation of expMod
 (Stack overflow possible!)
 ----de----
 rekursive Implementierung von expMod
 (Stack-Überlauf möglich!)
 """
 if x == 1:
 return a % n
 elif x % 2 == 1:
 return (a * _doc_expMod_recursive(a,x-1,n)) % n
 else:
 t = _doc_expMod_recursive(a, x//2, n)
 return (t * t) % n

def _doc_expMod_iter(a,x,n):
 """ iterative implementation of expMod
 ----de----
 iterative Implementierung von expMod
 """
 b, c = 1, a # TODO: verification
 while x > 1:
 if x % 2 == 1:
 b = (b*c) % n
 c = (c*c) % n
 x = x // 2
 return (b*c) % n

#............................................................

def _divisor(n):
 """ Pollard-Rho Algorithm (see Cormen, p. 845)
 """
 i, k = 1, 2
 x = rand(n)
 y = x
 while i < n:
 i += 1
 x = (x*x-1) % n
 d = gcd(y-x, n)
 if d != 1 and d != n:
 return d
 if i == k:
 y = x
 k = 2 * k

#............................................................

def isPrime(n,s=20):
 """ isPrime(n) 
 Returns bool (is <var>n</var> a prime number?) 
 Implementation using Miller Rabin test 
 (see. Cormen, Leiserson, Rivest: Algorithms, p. 841) 
 See also <a href="#prime"><code>prime</code></a>,
 <a href="#nextPrime"><code>nextPrime</code></a>.
 ----de----
 isPrime(n) 
 Rückgabewert bool (ist <var>n</var> Primzahl?) 
 Implementierung mit Miller-Rabin-Test 
 (vgl. Cormen, Leiserson, Rivest: Algorithms, S. 841) 
 Siehe auch <a href="#prime"><code>prime</code></a>,
 <a href="#nextPrime"><code>nextPrime</code></a>.
 """
 def _witness(a,n):
 """ witness for n being a product (see. Cormen p. 840)
 Local helper function
 ----de----
 Zeuge dafür, dass n zusammengesetzt ist (vgl. Cormen S. 840)
 Lokale Hilfsfunktion
 """
 b = NumeralSystem.binary(n-1)
 d = 1
 for i in b:
 x = d
 d = (d*d) % n
 if d == 1 and x != 1 and x != n-1:
 return True
 if i == 1:
 d = (d*a) % n
 return d != 1

 if n < 4:
 return n > 1
 for i in fromTo(1, s):
 a = 1 + rand(n-1)
 if _witness(a,n):
 return False
 return True

#............................................................
# nextPrime

def nextPrime(n):
 """ nextPrime(n) 
 next prime larger than or equal to <var>n</var> 
 See also <a href="#isPrime"><code>isPrime</code></a>,
 <a href="#prime"><code>prime</code></a>.
 ----de----
 nextPrime(n) 
 nächste Primzahl ab (einschließlich) <var>n</var> 
 Siehe auch <a href="#isPrime"><code>isPrime</code></a>,
 <a href="#prime"><code>prime</code></a>.
 """
 if n < 3:
 return 2
 if n % 2 == 0:
 n += 1
 while not isPrime(n):
 n += 2
 return n

#............................................................

# pre-calculated primes:
# prime(0), prime(1000), prime(2000), ..., prime(99000),
_P1000 = [ 0, 7919, 17389, 27449, 37813, 48611, 59359, 70657, 81799, 93179,
 104729, 116447, 128189, 139901, 151703, 163841, 176081, 187963, 200183, 212369,
 224737, 237203, 249439, 262139, 274529, 287117, 300023, 312583, 324949, 337541,
 350377, 363269, 376127, 389171, 401987, 414977, 427991, 440723, 453889, 467213,
 479909, 493127, 506131, 519227, 532333, 545747, 559081, 572311, 585493, 598687,
 611953, 625187, 638977, 652429, 665659, 679277, 692543, 706019, 719639, 732923,
 746773, 760267, 773603, 787207, 800573, 814279, 827719, 841459, 855359, 868771,
 882377, 896009, 909683, 923591, 937379, 951161, 965113, 978947, 993107, 1006721,
 1020379, 1034221, 1048129, 1062511, 1076143, 1090373, 1103923, 1117579, 1131617, 1145689,
 1159523, 1173301, 1187003, 1200949, 1215133, 1229269, 1243709, 1257517, 1271293, 1285517]

def prime(n):
 """ prime(n) 
 <var>n</var>-th prime number 
 See also <a href="#isPrime"><code>isPrime</code></a>,
 <a href="#nextPrime"><code>nextPrime</code></a>.
 ----de----
 prime(n) 
 <var>n</var>-te Primzahl 
 Siehe auch <a href="#isPrime"><code>isPrime</code></a>,
 <a href="#nextPrime"><code>nextPrime</code></a>.
 """
 if n < 100000:
 t, e = divmod(n, 1000)
 p = _P1000[t]
 for i in range(e):
 p = nextPrime(p+1)
 return p
 else:
 raise _MathGUIdeError("Argument is to large (must be < 100000)!")

#............................................................

def factor(n):
 """ factor(n) 
 Prime factorization 
 Returns a list of pairs [factor, multiplicity] 
 Example: 
 <code>factor(24)</code> 
 returns <code>[[2,3],[3,1]]</code>. 
 See also <a href="#factors"><code>factors</code></a>,
 <a href="#allFactors"><code>allFactors</code></a>.
 ----de----
 factor(n) 
 Primfaktorenzerlegung 
 Rückgabewert: Liste mit Paaren [Primfaktor, Exponent] 
 Beispiel: 
 <code>factor(24)</code> 
 gibt <code>[[2,3],[3,1]]</code> zurück. 
 Siehe auch <a href="#factors"><code>factors</code></a>,
 <a href="#allFactors"><code>allFactors</code></a>.
 """
 p = 2
 i = 0
 r = []
 while p*p <= n:
 while n % p == 0:
 i += 1
 n //= p
 if i > 0:
 r.append([p,i])
 p = nextPrime(p+1)
 i = 0
 if n > 1:
 r.append([n,1])
 return r

#............................................................

def factors(n):
 """ factors(n) 
 List of prime factors of <var>n</var> (without multiplicity) 
 Example: 
 <code>factors(24)</code>
 returns <code>[2, 3]</code>. 
 See also <a href="#factor"><code>factor</code></a>,
 <a href="#allFactors"><code>allFactors</code></a>.
 ----de----
 factors(n) 
 Liste der Primfaktoren von <var>n</var> (jeder nur einmal) 
 Beispiel: 
 <code>factors(24)</code>
 gibt <code>[2, 3]</code> zurück. 
 Siehe auch <a href="#factor"><code>factor</code></a>,
 <a href="#allFactors"><code>allFactors</code></a>.
 """
 return [f[0] for f in factor(n)]

#............................................................

def allFactors(n):
 """ allFactors(n) 
 List of all prime factors of <var>n</var> 
 Example: <code>allFactors(24)</code>
 returns <code>[2, 2, 2, 3]</code>. 
 See also <a href="#factor"><code>factor</code></a>,
 <a href="#factors"><code>factors</code></a>.
 ----de----
 allFactors(n) 
 Liste aller Primfaktoren von <var>n</var> 
 Beispiel: <code>allFactors(24)</code>
 gibt <code>[2, 2, 2, 3]</code> zurück. 
 Siehe auch <a href="#factor"><code>factor</code></a>,
 <a href="#factors"><code>factors</code></a>.
 """
 p = 2
 r = []
 while p*p <= n:
 while n % p == 0:
 r.append(p)
 n //= p
 p = nextPrime(p+1)
 if n > 1:
 r.append(n)
 return r

#............................................................

def divisors(n):
 """ divisors(n) 
 List of all divisors of <var>n</var>
 ----de----
 divisors(n) 
 Liste der Teiler von <var>n</var>
 """
 # TODO: more efficient implementation!
 d = []
 for i in fromTo(1, n//2):
 if n % i == 0:
 d.append(i)
 d.append(n)
 return d

#===============================================================
# Rational numbers and continued fractions
# Rationale Zahlen, Kettenbrüche
#---------------------------------------------------------------

class Rational: # Rational numbers # Rationale Zahlen
 """ Class <code>Rational</code>: Rational numbers 
 Rational(p, q) Example: <code>Rational(3,17)</code> 
 Rational(n) Example: <code>Rational(3)</code> 
 Rational(s) Example: <code>Rational("3/17")</code> 
 For <code>Rational</code> objects the arithmetical operators
 +, -, *, / are defined.
 ----de----
 Klasse <code>Rational</code>: Rationale Zahl 
 Rational(p, q) Beispiel: <code>Rational(3,17)</code> 
 Rational(n) Beispiel: <code>Rational(3)</code> 
 Rational(s) Beispiel: <code>Rational("3/17")</code> 
 <code>Rational</code>-Objekte können mit den
 arithmetischen Operatoren +, -, *, / verknüpft werden.
 """
 def operators(self):
 """ Rational.operators() 
 For documentation only 
 The following operators are defined in the class Rational: 
 <table>
 <tr><th>Op.</th><th>Function</th><th>Examples</th></tr>
 <tr><th>+</th><td>Addition</td><td><code>a + b; a += b</code></td></tr>
 <tr><th>-</th><td>Subtraction</td><td><code>a - b; a -= b</code></td></tr>
 <tr><th>*</th><td>Multiplication</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>/</th><td>Division</td><td><code>a / b; a /= b</code></td></tr>
 <tr><th>-</th><td>Unary minus</td><td><code>-a</code></td></tr>
 <tr><th>< > <= >=</th><td>Comparision operators</td><td></td></tr>
 </table>
 This method is for documentation only. It has no effect.
 ----de----
 Rational.operators() 
 Nur zur Dokumentation 
 Die Klasse Rational erlaubt folgende Operatoren: 
 <table>
 <tr><th>Op.</th><th>Funktion</th><th>Beispiele</th></tr>
 <tr><th>+</th><td>Addition</td><td><code>a + b; a += b</code></td></tr>
 <tr><th>-</th><td>Subtraktion</td><td><code>a - b; a -= b</code></td></tr>
 <tr><th>*</th><td>Multiplikation</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>/</th><td>Division</td><td><code>a / b; a /= b</code></td></tr>
 <tr><th>-</th><td>unäres Minus</td><td><code>-a</code></td></tr>
 <tr><th>< > <= >=</th><td>Vergleichsoperatoren</td><td></td></tr>
 </table>
 Diese Methode dient nur zur Dokumentation. Sie hat keine Wirkung.
 """
 pass

 def __init__(self, p, q=1):
 """ Constructor
 ----de----
 Konstruktor
 """
 if type(p) == float:
 r = ContFrac(p,12).toRational()
 self.p, self.q = r.p, r.q
 return
 if type(p) == Rational:
 self.p, self.q = p.p, p.q
 return
 if type(p) == str: # e.g. "3/17"
 l = p.split("/")
 p = int(l[0])
 if len(l) > 1:
 q = int(l[1])
 if q < 0:
 p,q = -p, -q # denominator always positive
 d = gcd(abs(p), q)
 if d > 1: # reduce fraction if possible
 p //= d
 q //= d
 self.p = p
 self.q = q

 def copy(self):
 """ r.copy() 
 independent copy ("clone") of r
 ----de----
 r.copy() 
 unabhängige Kopie von r
 """
 return Rational(self.p, self.q)

 def __abs__(self):
 """ absolute value
 ----de----
 Absolutbetrag
 """
 return Rational(abs(self.p), self.q)

 def abs(self):
 """ r.abs() 
 absolute value 
 Also as global function: abs(r) 
 ----de----
 r.abs() 
 Absolutbetrag 
 Auch globale Funktion: abs(r) 
 """
 return Rational(abs(self.p), self.q)

 def __float__(self):
 """ Conversion to floating point number
 ----de----
 Umwandlung in Gleitkommazahl
 """
 return float(self.p) / float(self.q)

 def toFloat(self):
 """ r.toFloat() 
 Conversion to floating point number 
 Also as global function: float(r) 
 ----de----
 r.toFloat() 
 Konvertierung in Gleitkommazahl 
 Auch globale Funktion: float(r) 
 """
 return float(self.p) / float(self.q)

 def __int__(self):
 """ Round to the nearest integer
 ----de----
 Auf-/Abrundung zur nächsten Ganzzahl
 """
 assert self.q > 0
 return (self.p + self.q/2) // self.q

 def isInteger(self):
 return self.q == 1

 def toInt(self):
 """ r.toInt() 
 Round to the nearest integer 
 Also as global function: int(r) 
 ----de----
 r.toInt() 
 Rundung zur nächsten ganzen Zahl 
 Auch globale Funktion: int(r) 
 """

 def __long__(self):
 """ Round to the nearest integer
 ----de----
 Auf-/Abrundung zur nächsten Ganzzahl
 """
 assert self.q > 0
 return (self.p + self.q/2) // self.q

 def _unify(self, other):
 return (self, Rational(other))

 def __repr__(self):
 """ Object representation as string
 """
 if self.q == 1:
 return str(self.p)
 else:
 return "{0}/{1}".format(self.p, self.q)

 def toStr(self):
 """ r.toStr() 
 Representation as string 
 Also as global function: str(r) 
 ----de----
 r.toStr() 
 Darstellung als String 
 Auch globale Funktion: str(r) 
 """

 def _cmp(a, b):
 if type(b) == float:
 return _cmp(float(a), b)
 if not isinstance(b, Rational):
 a,b = a._unify(b)
 return _cmp(a.p*b.q, b.p*a.q)

 def __lt__(a,b): return a._cmp(b) < 0
 def __le__(a,b): return a._cmp(b) <= 0
 def __eq__(a,b): return a._cmp(b) == 0
 def __ge__(a,b): return a._cmp(b) >= 0
 def __gt__(a,b): return a._cmp(b) > 0
 def __ne__(a,b): return a._cmp(b) != 0

 def __neg__(self):
 """ Unary operator -
 -a --> a.__neg__()
 """
 return Rational(-self.p, self.q)

 def __add__(a, b):
 """ a.__add__(b) 
 Addition 
 Also operator notation: a + b 
 Instead of <code>a = a + b</code> you can write <code>a += b</code>.
 """
 if type(b) == float:
 return float(a) + b
 if not isinstance(b, Rational):
 a, b = a._unify(b)
 return Rational(a.p*b.q + b.p*a.q, a.q*b.q)

 def __radd__(a, b):
 """ Addition, see __add__
 """
 return a + b

 def __sub__(a, b):
 """ Subtraction:
 a - b --> a.__sub__(b), if a is of type Rational
 b.__rsub__(a), otherwise
 """
 if type(b) == float:
 return float(a) - b
 if not isinstance(b, Rational):
 a,b = a._unify(b)
 return Rational(a.p*b.q - b.p*a.q, a.q*b.q)

 def __rsub__(b, a):
 """ Subtraction, see __sub__
 """
 if type(a) == float:
 return a - float(b)
 if not isinstance(a, Rational):
 b,a = b._unify(a)
 return a - b

 def __mul__(a, b):
 """ Multiplication:
 a * b --> a.__mul__(b), if a is of type Rational
 b.__rmul__(a), otherwise
 """
 if type(b) == float:
 return float(a) * b
 if type(b) in (Matrix, Poly):
 return b * a
 if not isinstance(b, Rational):
 a,b = a._unify(b)
 return Rational(a.p*b.p, a.q*b.q)

 def __rmul__(b, a):
 """ Multiplication, see __mul__
 """
 return b * a

 def __div__(a, b):
 """ exact Division
 a / b --> a.__div__(b), if a is of type Rational
 b.__rdiv__(a), otherwise
 """
 if type(b) == float:
 return float(a) / b
 if not isinstance(b, Rational):
 a,b = a._unify(b)
 return Rational(a.p*b.q, a.q*b.p)

 def __rdiv__(b, a):
 """ exact Division, see __div__
 """
 if type(a) == float:
 return a / float(b)
 if not isinstance(a, Rational):
 b,a = b._unify(a)
 return Rational(a.p*b.q, a.q*b.p)

 def __truediv__(a, b): return a.__div__(b)
 def __rtruediv__(a, b): return a.__rdiv__(b)
 def __floordiv__(a, b): return a.__div__(b)
 def __rfloordiv__(a, b): return a.__rdiv__(b)

def toNumber(s):
 """ toNumber(s) 
 converts the string s into a rational or floating point number 
 Examples: 
 <code>toNumber("3.1")</code> 
 <code>toNumber("7/3")</code>
 ----de----
 toNumber(s) 
 verwandelt den String s in eine rationale oder Gleitkommazahl 
 Beispiele: 
 <code>toNumber("3.1")</code> 
 <code>toNumber("7/3")</code>
 """
 if "." in s or "e" in s or "E" in s:
 return float(s)
 else:
 return Rational(s)

#---------------------------------------------------------------

class ContFrac: # Continued fractions # Kettenbrüche
 """ ContFrac(number, steps=5) 
 Continued fraction from number with maximal length steps 
 Example: 
 <code>ContFrac(pi, 3)</code> returns
 <code>ContFrac(3, 7, 15)</code>. 
 This represent the number 3 + 1/(7 + 1/15).
 ----de----
 ContFrac(number, steps=5) 
 Kettenbruch aus number mit Maximallänge steps 
 Beispiel: 
 <code>ContFrac(pi, 3)</code> ergibt
 <code>ContFrac(3, 7, 15)</code>. 
 Das steht für die Zahl 3 + 1/(7 + 1/15).
 """
 def __init__(self, number, steps=5):
 """ Constructor
 ----de----
 Konstruktor
 """
 self._l = []
 for i in range(steps):
 i = floor(number)
 self._l.append(i)
 if abs(number-i) < 0.00000001:
 break
 number = 1 / (number-i)

 def __repr__(self):
 """ Object representation as string
 """
 return "ContFrac({0})".format(join(self._l, ","))

 def toRational(self, n=0):
 """ toRational(n=0) 
 Rational number defined by the first n list elements 
 Converts the continued fraction int a rational number. 
 Example: 
 <code>ContFrac(pi).toRational(2)</code> returns 22/7.
 ----de----
 toRational(n=0) 
 Rationale Zahl aus den ersten n Gliedern 
 Wandelt den Kettenbruch in eine rationale Zahl um. 
 Beispiel: 
 <code>ContFrac(pi).toRational(2)</code> ergibt 22/7.
 """
 def contFracListToRat(cf):
 if len(cf) == 0: return Rational(0)
 if len(cf) == 1: return Rational(cf[0])
 return cf[0] + Rational(1) / contFracListToRat(cf[1:])

 if n == 0 or n > len(self._l):
 n = len(self._l)
 return contFracListToRat(self._l[0:n])

#===============================================================
# Lineare Algebra
# Klassen Vector und Matrix
#---------------------------------------------------------------

class Vector (list): # Vectors # Vektoren
 """ Class <code>Vector</code> (with arithmetic operators) 
 Vector(v) (v Liste or Tuple) 
 The Indices of the elements of Vector are counted from 0 (as in Python lists etc.).
 Example: 
 <code>Vector([1,2,3])</code>
 ----de----
 Klasse <code>Vector</code> (mit arithmetischen Operatoren) 
 Vector(v) (v Liste oder Tupel) 
 Die Indizes der Elemente von Vector werden (wie in Python-Listen etc.) ab 0 gezählt.
 Beispiel: 
 <code>Vector([1,2,3])</code>
 """
 #............................................................
 _staticMethods = ("fromFunction", "fromString")

 @staticmethod
 def fromFunction(n, fn, offset=0):
 """ Vector.fromFunction(n, fn, offset=0) 
 Vector, the elements of which are to be calculated using the function fn 
 Returns: <code>Vector([fn(i+offset) for i in range(n)])</code> 
 Example: 
 <code>Vector.fromFunction(5, sqrt, 1)</code> returns Vector([1.0, 1.414, 1.732, 2.0, 2.236])
 ----de----
 Vector.fromFunction(n, fn, offset=0) 
 Vector, dessen Elemente mit der Funktion fn berechnet werden 
 Rückgabewert: <code>Vector([fn(i+offset) for i in range(n)])</code> 
 Beispiel: 
 <code>Vector.fromFunction(5, sqrt, 1)</code> ergibt Vector([1.0, 1.414, 1.732, 2.0, 2.236])
 """
 return Vector([fn(i+offset) for i in range(n)])

 @staticmethod
 def fromString(s):
 """ Vector.fromString(s) 
 Vector from String 
 Elements (may also bbe Rationals) must be given as comma separated string. 
 Example: 
 <code>Vector.fromString("1, 2/3, 2")</code>
 ----de----
 Vector.fromString(s) 
 Vector aus String 
 Elemente (auch rationale Zahlen) müssen im String durch Kommas getrennt werden. 
 Beispiel: 
 <code>Vector.fromString("1, 2/3, 2")</code>
 """
 return Vector([toNumber(x) for x in s.split(",")])
 #............................................................

 def operators(self):
 """ Vector.operators() 
 For documentation only 
 The following operators are defined in the class Vector: 
 <table>
 <tr><th>Op.</th><th>Function</th><th>Examples</th></tr>
 <tr><th>+</th><td>elementwise Addition</td><td><code>a + b; a += b</code></td></tr>
 <tr><th>-</th><td>elementwise Subtraction</td><td><code>a - b; a -= b</code></td></tr>
 <tr><th>*</th><td>if both operators are Vectors: inner product,
 otherwise: scalar multiplication (elementwise)</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>-</th><td>Unary minus</td><td><code>-a</code></td></tr>
 <tr><th>[ ]</th><td>Index operator</td><td><code>A[i]</code> (Element)</td></tr>
 </table>
 This method is for documentation only. It has no effect.
 ----de----
 Vector.operators() 
 Nur zur Dokumentation 
 Die Klasse Vector erlaubt folgende Operatoren: 
 <table>
 <tr><th>Op.</th><th>Funktion</th><th>Beispiele</th></tr>
 <tr><th>+</th><td>Elementweise Addition</td><td><code>a + b; a += b</code></td></tr>
 <tr><th>-</th><td>Elementweise Subtraktion</td><td><code>a - b; a -= b</code></td></tr>
 <tr><th>*</th><td>wenn beide Operanden Vektoren sind: Skalarprodukt,
 sonst Skalarmultiplikation (elementweise)</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>-</th><td>unäres Minus</td><td><code>-a</code></td></tr>
 <tr><th>[ ]</th><td>Indexoperator</td><td><code>A[i]</code> (Element)</td></tr>
 </table>
 Diese Methode dient nur zur Dokumentation. Sie hat keine Wirkung.
 """
 pass

 def __init__(v, n):
 """ Constructor
 ----de----
 Konstruktor
 """
 if isinstance(n, list):
 list.__init__(v, n)
 elif isinstance(n, tuple):
 list.__init__(v, list(n))

 def __repr__(A):
 """ Object representation as string
 """
 return "Vector({0})".format(list.__repr__(A))

 def __invert__(v):
 """ Operator ~ (transposition)
 """
 return v.transp()

 def transp(v):
 """ v.transp() 
 berechnet die transponierte Matrix (Spaltenvektor) 
 kürzer: ~v
 ----de----
 v.transp() 
 calculates the transposed Matrix (column vektor) 
 shorter: ~v
 """
 return Matrix([[a] for a in v])

 def join(x, y):
 """ v.join(w) 
 Concatenation with Vector w
 ----de----
 v.join(w) 
 Verkettung mit dem Vektor w
 """
 return Vector(list(x)+list(y))

 def __neg__(v):
 """ unary oparator - (negative Vektor)
 """
 return Vector([-x for x in v])

 def __add__(x, y):
 """ operators + and += (Elementwise addition)
 """
 # a + b --> a.__add__ (b), if a is of type Vector
 # b.__radd__(a), otherwise
 n = len(x)
 assert n == len(y)
 return Vector([x[i]+y[i] for i in range(n)])

 def __sub__(x, y):
 """ operators - and -= (elementwise subtraction)
 """
 return x + (-y)

 def __mul__(x, y):
 """ Operators * and *=
 Multiplikation
 if both operands are vektors sind: inner product,
 otherwise scalar multiplication (elementwise)
 """
 if isinstance(y, Vector): # Skalarprodukt
 n = len(x)
 assert n == len(y)
 return sum([x[i]*y[i] for i in range(n)])
 else:
 return Vector([y*x[i] for i in range(len(x))])

 def __rmul__(x, c):
 return x * c

 def __imul__(x, c):
 for i in range(len(x)):
 x[i] *= c
 return x

 def norm(x):
 """ v.norm() 
 Norm des Vektors v
 ----de----
 v.norm() 
 Norm of the Vector v
 """
 return math.sqrt(x*x)

class Matrix (Vector): # Matrices # Matrizen
 """ class <code>Matrix</code> 
 see also: Menu: Insert -- Matrix 
 Matrix(v): v list of lists with uniform length. 
 Example: 
 <code>Matrix([[11,12,13],[21,22,23]])</code> 
 The class <code>Matrix</code> ist derived from <code>Vector</code>: 
 A <code>Matrix</code> object ist a <code>Vector</code>
 of <code>Vector</code> objects with uniform length. 
 The elements can be addressed using double indices (e.g.: <code>A[i,k]</code>).
 Indices are counted from 0 (as in Python lists etc.). 
 see also. class methods with <code>Matrix.</code>
 ----de----
 Klasse <code>Matrix</code> 
 vgl. Menü: Einfügen -- Matrix 
 Matrix(v): v Liste von gleichlangen Listen, 
 Beispiel: 
 <code>Matrix([[11,12,13],[21,22,23]])</code> 
 Die Klasse <code>Matrix</code> ist abgeleitet von <code>Vector</code>: 
 Ein <code>Matrix</code>-Objekt ist ein <code>Vector</code>
 von gleichlangen <code>Vector</code>-Objekten. 
 Die Elemente können mit Doppelindizes (z.B.: <code>A[i,k]</code>) angesprochen werden.
 Indizes werden (wie in Python-Listen etc.) ab 0 gezählt. 
 vgl. Klassenmethoden mit <code>Matrix.</code>
 """
 #............................................................
 _staticMethods = ("null", "identity", "fromFunction", "fromString", "random")

 @staticmethod
 def null(m, n=0):
 """ Matrix.null(m, n=0) 
 m*n-Null matrix 
 Square Matrix, if n is omitted. 
 ----de----
 Matrix.null(m, n=0) 
 m*n-Nullmatrix 
 Quadratische Matrix, falls n weggelassen wird. 
 """
 if n == 0:
 n = m
 return Matrix([ [0 for k in range(n)]
 for i in range(m)])

 @staticmethod
 def identity(n):
 """ Matrix.identity(n) 
 n*n Identity Matrix 
 ----de----
 Matrix.identity(n) 
 n*n-Einheitsmatrix 
 """
 return Matrix([ [Rational(int(i==k)) for k in range(n)]
 for i in range(n)])

 @staticmethod
 def fromFunction(m, n, fn, offset=0):
 """ Matrix.fromFunction(m,n, fn, offset=0) 
 <code>m*n</code>-Matrix defined by function <code>fn</code> 
 <code>fn</code> must be a binary function. 
 The Matrix element <code>A[i,k]</code> is defined as <code>fn(i+offset,k+offset)</code>. 
 Example: 
 <code>Matrix.fromFunction(2,4, pow, 1)</code> 
 returns the Matrix<pre>
 / 1 1 1 1 \
 | |
 \ 2 4 8 16 /</pre>
 ----de----
 Matrix.fromFunction(m,n, fn, offset=0) 
 Durch Funktion <code>fn</code> berechnete <code>m*n</code>-Matrix 
 <code>fn</code> muss eine zweistellige Funktion sein. 
 Das allgemeine Element <code>A[i,k]</code> wird mit <code>fn(i+offset,k+offset)</code> berechnet. 
 Beispiel: 
 <code>Matrix.fromFunction(2,4, pow, 1)</code> 
 liefert die Matrix<pre>
 / 1 1 1 1 \
 | |
 \ 2 4 8 16 /</pre>
 """
 return Matrix([ [fn(i+offset,k+offset) for k in range(n)]
 for i in range(m)])

 @staticmethod
 def fromString(s):
 """ Matrix.fromString(s) 
 Matrix given by a string 
 The string s must contains the elementa row by row. 
 The rows are divided by semicolon, the elemente within a row by comma. 
 Elemente may also be rational (notated with slash). 
 Example: 
 <code>Matrix.fromString("1, 2, 3.14; 4/5, 5, 6")</code>
 returns the matrix<pre>
 / 1 2 3.14 \
 | |
 \ 4/5 5 6 /</pre>
 ----de----
 Matrix.fromString(s) 
 Matrix aus String 
 Im String s werden die Elemente der Matrix zeilenweise angegeben. 
 Die Zeilen werden durch Semikolon getrennt, die Elemente innerhalb der Zeilen durch Komma. 
 Elemente können auch rational sein (mit Schrägstrich). 
 Beispiel: 
 <code>Matrix.fromString("1, 2, 3.14; 4/5, 5, 6")</code>
 liefert die Matrix<pre>
 / 1 2 3.14 \
 | |
 \ 4/5 5 6 /</pre>
 """
 return Matrix([ [toNumber(a) for a in row.split(",")]
 for row in s.split(";")])

 @staticmethod
 def random(m,n=0, r=100):
 """ Matrix.random(m,n=0, r=100) 
 m*n matrix with random elements in 0..r-1
 ----de----
 Matrix.random(m,n=0, r=100) 
 m*n-Matrix mit Zufallswerten 0..r-1
 """
 if n == 0:
 n = m
 return Matrix([ [Rational(rand(r)) for k in range(n)]
 for i in range(m)])
 #............................................................

 def operators(self):
 """ Matrix.operators() 
 For documentation only 
 The following operators are defined in the class Matrix: 
 <table>
 <tr><th>Op.</th><th>Function</th><th width="30%">Examples</th></tr>
 <tr><th>+</th><td>Addition</td><td><code>a + b; a += b</code></td></tr>
 <tr><th>-</th><td>Subtraction</td><td><code>a - b; a -= b</code></td></tr>
 <tr><th>*</th><td>Matrizenmultiplikation if both operands are Matrices, otherwise: Scalar multiplication</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>-</th><td>Unary minus</td><td><code>-a</code></td></tr>
 <tr><th>~</th><td>Transposed Matrix</td><td><code>~A</code></td></tr>
 <tr><th>|</th><td>Concatenation</td><td><code>A | B; A |= B</code></td></tr>
 <tr><th>[ ]</th><td>Index operator</td><td><code>A[i,k]</code> <code>A[i]</code> (Row Vector)</td></tr>
 </table>
 This method is for documentation only. It has no effect.
 ----de----
 Matrix.operators() 
 Nur zur Dokumentation 
 Die Klasse Matrix erlaubt folgende Operatoren: 
 <table>
 <tr><th>Op.</th><th>Funktion</th><th width="30%">Beispiele</th></tr>
 <tr><th>+</th><td>Addition</td><td><code>A + B; A += B</code></td></tr>
 <tr><th>-</th><td>Subtraktion</td><td><code>A - B; A -= B</code></td></tr>
 <tr><th>*</th><td>Matrizenmultiplikation (wenn beide Operanden Matrizen sind) bzw. Skalarmultiplikation (sonst)</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>-</th><td>unäres Minus</td><td><code>-A</code></td></tr>
 <tr><th>~</th><td>Transponierte Matrix</td><td><code>~A</code></td></tr>
 <tr><th>|</th><td>Verkettung</td><td><code>A | B; A |= B</code></td></tr>
 <tr><th>[ ]</th><td>Indexoperator</td><td><code>A[i,k]</code> <code>A[i]</code> (Zeilenvektor)</td></tr>
 </table>
 Diese Methode dient nur zur Dokumentation. Sie hat keine Wirkung.
 """
 pass

 def __init__(A, v):
 """ Constructor
 ----de----
 Konstruktor
 """
 Vector.__init__(A, [Vector(v[i]) for i in range(len(v))])
 A._makeRational()

 def _makeRational(A):
 for i in A.rowRange():
 for k in A.colRange():
 if isInteger(A[i,k]):
 A[i,k] = Rational(A[i,k])

 def float(A):
 """ A.float() 
 returns the Matrix A with all elements converted to <code>float</code>
 ----de----
 A.float() 
 gibt die Matrix A zurück, in der alle Elemente nach <code>float</code> konvertiert sind
 """
 return Matrix([ [float(A[i,k]) for k in A.colRange()]
 for i in A.rowRange()])

 def copy(A):
 """ A.copy() 
 independent copy of the matrix ("clone") A
 ----de----
 A.copy() 
 unabhängige Kopie der Matrix A
 """
 return Matrix([ [copy(A[i,k]) for k in A.colRange()]
 for i in A.rowRange()])

 def isSquare(A):
 """ A.isSquare() 
 Truth value (A is square matrix)
 ----de----
 A.isSquare() 
 Wahrheitswert (A quadratisch)
 """
 return A.height() == A.width()

 def str(A, i,k):
 """ A.str(i,k) 
 string representation of the element [i,k]
 ----de----
 A.str(i,k) 
 String-Darstellung des Elements [i,k]
 """
 if type(A[i,k]) == float:
 return "%g" % A[i,k]
 else:
 return str(A[i,k])

 def _pp(A, name=""):
 n = A.height()
 s = """{{html}}<table><tr valign="middle">
 <td width="30">{0}</td>
 <td><table cellspacing="0" cellpadding="-1">""".format(name)
 for i in range(n):
 left, right = "21", "22"
 if n > 1:
 if i == 0: left, right = "11", "12"
 elif i == n-1: left, right = "31", "32"
 else:
 left, right = "01", "02"
 s += """<tr><td align="center" valign="middle"><img src=":/img/matrix/{0}.png"/></td>""".format(left)
 spc = {False: " ", True: ""}
 for k in range(A.width()):
 s += """<td align="center" valign="middle">{1}{0}{2}</td>""".format(A[i,k], spc[k==0], spc[k==A.width()-1])
 s += """<td align="center" valign="middle"><img src=":/img/matrix/{0}.png"/></td></tr>""".format(right)
 s += """</table></td></table>{/html}"""
 return s

 def __repr__(A):
 """ Object representation as string
 """
 return A._pp()

 def __getitem__(A, i):
 """ Definition of the Index operator []
 Counts from 0
 A[i] --> i-th row vector
 A[i,k] --> element in i-th row, k-th column
 """
 if isinstance(i, tuple):
 return A[i[0]][i[1]]
 else:
 return Vector.__getitem__(A, i)

 def __setitem__(A, i, x):
 """ Definition of the Index operator [] for assignments
 Counts from 0
 A[i] --> i-th row vector
 A[i,k] --> element in i-th row, k-th column
 """
 if isinstance(i, tuple):
 A[i[0]][i[1]] = x
 else:
 Vector.__setitem__(A, i, x)

 def height(A):
 """ A.height() 
 Number of rows of the Matrix A
 ----de----
 A.height() 
 Anzahl der Zeilen der Matrix A
 """
 return len(A)

 def width(A):
 """ A.height() 
 Number of columns of the Matrix A
 ----de----
 A.width() 
 Anzahl der Spalten der Matrix A
 """
 return len(A[0])

 def rowRange(A):
 """ A.rowRange() 
 List of all row indices
 ----de----
 A.rowRange() 
 Liste aller Zeilenindizes
 """
 return range(A.height())

 def colRange(A):
 """ A.rowRange() 
 List of all column indices
 ----de----
 A.rowRange() 
 Liste aller Spaltenindizes
 """
 return range(A.width())

 def __neg__(A):
 """ unary minus operator
 """
 return Matrix([ [-A[i,k] for k in A.colRange()]
 for i in A.rowRange()])

 def joinRight(A, B):
 """ A.joinRight(B) 
 Horizontal concatenation of A with B 
 Requires: Both Matrices must have the same row count.
 ----de----
 A.concat(B) 
 Horizontale Verkettung von A mit B 
 Voraussetzung: Die beiden Matrizen müssen gleich viele Zeilen haben.
 """
 assert A.height() == B.height()
 return Matrix([ list(A[i])+list(B[i]) for i in A.rowRange()])

 def joinBottom(A, B):
 """ A.joinBottom(B) 
 Vertical concatenation of A with B 
 Requires: Both Matrices must have the same column count.
 ----de----
 A.concat(B) 
 Vertikale Verkettung von A mit B 
 Voraussetzung: Die beiden Matrizen müssen gleich viele Spalten haben.
 """
 assert A.width() == B.width()
 def el(i,k):
 if i < A.height(): return A[i,k]
 else: return B[i-A.height(), k]
 return Matrix.fromFunction(A.height()+B.height(), A.width(), el)

 def __or__(A, B):
 """ A.__or__(B) 
 Horizontal concatenation with B 
 A | B is the same as A.joinRight(B) 
 Instead of A = A | B you may write A |= B.
 Requires: Both matrices must have the same row count.
 """
 assert A.height() == B.height()
 return Matrix([ list(A[i])+list(B[i]) for i in A.rowRange()])

 def __add__(A, B):
 """ A.__add__(B) 
 Matrix addition operator 
 Instead of A = A + B you may write A += B.
 """
 # a + b --> a.__add__ (b), wenn a vom Typ Matrix ist
 # b.__radd__(a), sonst
 assert A.height() == B.height() and A.width() == B.width()
 return Matrix([ [A[i,k]+B[i,k] for k in A.colRange()]
 for i in A.rowRange()
 ])

 def __sub__(A, B):
 """ A.__sub__(B) 
 Matrix subtraction operator 
 Instead of A = A - B you may write A -= B.
 """
 return A + (-B)

 def __mul__(A, B):
 """ A.__sub__(B) 
 Scalar or matrix multiplication 
 Instead of A = A * B you may write A *= B.
 """
 # If both operands are matrices, matrix multiplication is applied,
 # otherwise scalar multiplikation
 if isinstance(B, Matrix): # matrix multiplication
 assert A.width() == B.height()
 return Matrix(
 [ [ sum([A[i,j] * B[j,k] for j in A.colRange()])
 for k in B.colRange()]
 for i in A.rowRange()])
 else: # scalar multiplikation
 return Matrix([ [B*A[i,k] for k in A.colRange()]
 for i in A.rowRange()])

 def __rmul__(A, B):
 assert not isinstance(B, Matrix)
 return A * B

 def transp(A):
 """ A.transp() 
 returns the transposed Matrix 
 short form: ~A
 ----de----
 A.transp() 
 berechnet die Transponierte der Matrix A 
 kürzer: ~A
 """
 return Matrix([ [A[i,k] for i in A.rowRange()]
 for k in A.colRange()])

 def __invert__(A):
 """ Operator ~ (transposition)
 """
 return A.transp()

 def submatrix(A, i, k, m, n):
 """ A.submatrix(i, k, m, n) 
 m*n Submatrix (rows i..i+m-1, columns k..k+n-1) 
 ----de----
 A.submatrix(i, k, m, n) 
 m*n-Untermatrix (Zeilen i..i+m-1, Spalten k..k+n-1) 
 """
 assert i+m-1 < A.height() and k+n-1 < A.width()
 return Matrix([ [A[i1,k1] for k1 in range(k,k+n)]
 for i1 in range(i,i+m)])

 def complement(A, i, k):
 """ A.complement(i, k) 
 Matrix without row i and column k 
 (Algebraic complement)
 ----de----
 A.complement(i, k) 
 Matrix ohne die Zeile i und die Spalte k 
 (Algebraisches Komplement)
 """
 return Matrix([ [A[i1,k1] for k1 in A.colRange() if k1 != k]
 for i1 in A.rowRange() if i1 != i])

 def minor(A, i, k):
 """ A.minor(i, k) 
 Minor for row i, column k 
 (Determinant of the algebraic complement)
 ----de----
 A.minor(i, k) 
 Minor zur Zeile i und der Spalte k 
 (Determinante des Algebraischen Komplements)
 """
 return A.complement(i,k).det()

 def cofactor(A, i, k):
 """ A.cofactor(i, k) 
 Cofactor for row i, columns k 
 (signed minor)
 ----de----
 A.cofactor(i, k) 
 Kofaktor zur Zeile i und der Spalte k 
 (Minor mit Vorzeichen)
 """
 return (-1)**(i+k) * A.minor(i,k)

 def adjoint(A):
 """ A.adjoint() 
 Adjoint matrix 
 ----de----
 A.adjoint() 
 Adjungierte Matrix 
 """
 assert A.isSquare()
 return Matrix([ [A.cofactor(k,i) for k in A.rowRange()]
 for i in A.colRange()])
 def det(A):
 """ A.det() 
 Determinant of the matrix A 
 ----de----
 A.det() 
 Determinante der Matrix A 
 """
 if not A.isSquare():
 raise _MathGUIdeError("A is not a square matrix")
 n = A.height()
 if n == 1:
 return A[0,0]
 else:
 return sum([A[0,k] * A.cofactor(0,k)
 for k in A.colRange()])

 def _gaussElim0(A, jordan=True):
 """ unsichere Vorstufe zu gaussElim
 """
 m, n = A.height(), A.width()
 k = 0
 for i0 in A.rowRange():
 k = i0
 # oberste Zeile für führende Eins passend durchmultiplizieren
 A[i0] *= 1/A[i0,k]

 # Passende Vielfache zu anderen Zeilen addieren, so dass
 # unterhalb der führenden Eins Nullen entstehen
 for i in range(i0+1, m):
 A[i] -= A[i0] * A[i,k]
 if jordan:
 # Gauß-Jordan: dto. auch oberhalb der führenden Eins
 for i in range(i0):
 A[i] -= A[i0] * A[i,k]

 def gaussElim(self, jordan=True):
 """ A.gaussElim(jordan=True) 
 Row echelon form of the matrix A 
 If <var>jordan</var> is <code>True</code>, the reduced row echelon form
 is returned (Gaussian elimination), otherwise the normal row echelon form.
 ----de----
 A.gaussElim(jordan=True) 
 Zeilenstufenform der Matrix A. 
 Wenn <var>jordan</var> <code>True</code> ist, wird die reduzierte Zeilenstufenform
 (Gauß-Jordan-Elimination) zurückgegeben, sonst die einfache Zeilenstufenform (Gauß-Elimination).
 """
 A = self.copy()
 m, n = A.height(), A.width()
 k = -1
 for i0 in A.rowRange():
 # Bestimme die am weitesten links stehende Spalte k, die
 # (ab Zeile i0) von Null verschiedene Elemente enthält:
 k += 1
 aMax, iMax = 0, i0
 while aMax == 0 and k < n:
 for i in range(i0,m):
 if abs(A[i,k]) > aMax:
 aMax, iMax = abs(A[i,k]), i
 if aMax == 0:
 k += 1
 if k < n:
 # Oberste Zeile mit der vertauschen, die das (dem Betrag nach)
 # größte Element in Spalte i0 enthält
 A[iMax], A[i0] = A[i0], A[iMax]

 # oberste Zeile für führende Eins passend durchmultiplizieren
 A[i0] *= 1/A[i0,k]

 # Passende Vielfache zu anderen Zeilen addieren, so dass
 # unterhalb der führenden Eins Nullen entstehen
 for i in range(i0+1, m):
 A[i] -= A[i0] * A[i,k]
 if jordan:
 # Gauß-Jordan: dto. auch oberhalb der führenden Eins
 for i in range(i0):
 A[i] -= A[i0] * A[i,k]
 return A

 def rank(A):
 """ A.rank() 
 Rank of the matrix A 
 ----de----
 A.rank() 
 Rang der Matrix A 
 """
 B = A.gaussElim()
 r = 0
 while r < B.height() and B[r].norm() != 0:
 r += 1
 return r

 def inverse(A):
 """ A.inverse() 
 Inverse og the matrix A 
 ----de----
 A.inverse() 
 Inverse Matrix A 
 """
 if not A.isSquare():
 raise _MathGUIdeError("A is not a square matrix")
 n = A.height()
 I = Matrix.identity(n)
 B = A.joinRight(I).gaussElim()
 if B.submatrix(0,0,n,n) != I:
 raise _MathGUIdeError("Matrix is not invertible")
 return B.submatrix(0,n,n,n)

 def gramSchmidt(self, normalize=False):
 """ A.gramSchmidt(normalize=False) 
 Orthogonalizing using the Gram–Schmidt process 
 <code>A</code> must be the matrix of the (linearly independent) vectors
 to be orthogonalized (given as rows). 
 If <var>normalize</var> is True, the vectors are normalized. 
 return value: Matrix consisting of the orthogonalen vectors. 
 Example: <pre>
 A = Matrix([[1,1],[2,0]])
 A.gramSchmidt()</pre>
 returns the matrix<pre>
 / 1 1 \
 | |
 \ 1 -1 /</pre>
 ----de----
 A.gramSchmidt(normalize=False) 
 Gram-Schmidt'sches Orthogonalisierungsverfahren. 
 <code>A</code> muss Matrix aus den (linear unabhängigen)
 zu orthogonalisierenden Vektoren (als Zeilen) sein. 
 Wenn <var>normalize</var> True ist, werden die Vektoren normalisiert. 
 Rückgabewert: Matrix mit den orthogonalen Vektoren. 
 Beispiel<pre>
 A = Matrix([[1,1],[2,0]])
 A.gramSchmidt()</pre>
 liefert die Matrix<pre>
 / 1 1 \
 | |
 \ 1 -1 /</pre>
 """
 A = self.copy()
 n = A.height()
 for i in range(1,n):
 for k in range(i):
 A[i] -= A[k] * ((A[i]*A[k])/(A[k]*A[k]))
 if normalize:
 for v in A:
 v *= 1/v.norm()
 return A

 def _luRecursive(A):
 """ Vorstufe: funktioniert nicht für jede reguläre Matrix!
 """
 assert A.isSquare()
 n = A.height()
 if n == 1: # Rekursionsabbruch: 1·1-Matrix
 return Matrix.identity(1), A.copy()

 # /a11 : w \ / 1 : 0 \ /a11: w \
 # |..........| |.............| |................ |
 # A = | : | = | 1 : | * | : 1 |
 # | ~v : A1 | | --- ~v : I | | 0 : A1- --- ~v*w|
 # \ : / \ a11 : / \ : a11 /
 v = A.submatrix(1,0, n-1,1)
 w = A.submatrix(0,1, 1,n-1)
 A1 = A.submatrix(1,1, n-1,n-1)
 B = A1 - (1/A[0,0]) * v * w
 L1, U1 = B._luRecursive() # Rekursiver Aufruf für (n-1)·(n-1)-Matrix
 L_top = Matrix.identity(1).joinRight(Matrix.null(1,n-1))
 L_btm = ((1/A[0,0]) * v).joinRight(L1)
 L = L_top.joinBottom(L_btm)
 U_top = A.submatrix(0,0, 1,1).joinRight(w)
 U_btm = Matrix.null(n-1,1).joinRight(U1)
 U = U_top.joinBottom(U_btm)
 return L,U

 def _luIter(self):
 A = self.copy()
 assert A.isSquare()
 n = A.height()
 for k in range(n-1):
 assert A[k,k] != 0 # sonst scheitert der Algorithmus!
 for i in range(k+1,n):
 A[i,k] /= A[k,k]
 for j in range(k+1,n):
 A[i,j] -= A[i,k] * A[k,j]
 # A in linke und rechte Dreiecksmatrix zerlegen
 L = Matrix.identity(n)
 U = Matrix.identity(n)
 for i in range(n):
 for k in range(n):
 if k < i: L[i,k] = A[i,k]
 else: U[i,k] = A[i,k]
 return L, U

 def lup(self):
 """ A.lup() 
 LUP decomposition of the matrix 
 Returns a triplet (L, U, P) with P*A = L*R, where: 
 P is a Permutation matrix, 
 L is a Lower triangular matrix where all diagonal elements are one. 
 R is an upper triangular matrix. 
 <bRequires: Matrix ist invertible.
 ----de----
 .lup() 
 LRP-Zerlegung (LUP decomposition) der Matrix 
 Rückgabewert: Tripel (L, U, P) mit P*A = L*R, wobei: 
 P Permutationsmatrix, 
 L linke (lower) Dreiecksmatrix mit Einsen auf Diagonale, 
 U rechte (upper) Dreiecksmatrix ist. 
 Voraussetzung: Matrix ist invertierbar.
 """
 A = self.copy()
 assert A.isSquare()
 n = A.height()
 P = Matrix.identity(n)
 for k in range(n-1):
 # Suche das dem Betrag nach größte Element
 # der k-ten Spalte ab der k-ten Zeile und seinen Zeilenindex
 absMax, kMax = max([(abs(A[i,k]), -i) for i in range(k,n)])
 kMax = -kMax # -i im Paar --> bei mehreren gleichen erstes
 assert absMax > 0

 # Vertausche die Zeile des größten Elements
 # mit der k-ten Zeile (Permutation und Matrix)
 P[k], P[kMax] = P[kMax], P[k]
 A[k], A[kMax] = A[kMax], A[k]

 for i in range(k+1,n):
 A[i,k] /= A[k,k]
 for j in range(k+1,n):
 A[i,j] -= A[i,k] * A[k,j]

 # A in linke und rechte Dreiecksmatrix zerlegen
 L = Matrix.identity(n)
 U = Matrix.identity(n)
 for i in range(n):
 for k in range(n):
 if k < i: L[i,k] = A[i,k]
 else: U[i,k] = A[i,k]
 return L, U, P

 def solve(A, b):
 """ A.solve(b) 
 Solution of the linear equation system Ax = b 
 b may be Vektor or column Matrix.
 ----de----
 A.solve(b) 
 Lösung x des linearen Gleichungssystems Ax = b 
 b kann Vektor oder Spaltenmatrix sein.
 """
 if not isinstance(b, Matrix):
 assert isinstance(b, Vector)
 b = ~b
 return A.inverse() * b

 def leastSquares(A, b):
 """ A.leastSquares(b) 
 Solution of the normal system 
 Returns: <code>(~A*A).solve(~A*b)</code>
 ----de----
 A.leastSquares(b) 
 Lösung der Normalgleichung 
 Rückgabewert: <code>(~A*A).solve(~A*b)</code>
 """
 return (~A*A).solve(~A*b)

def fit(data, values):
 """ fit(data, values) 
 Interpolation using the Gaussian method of least squares 
 <code>fit(data, values)</code>
 calculates (using the method of least squares)
 a linear combination of the functions mapping x
 to the terms in <code>values</code>. 
 (corresponds to the Mathematica function <code>Fit</code>). 
 Example:<pre>
fit([(-1,2), (1,1), (2,1), (3,0), (5,3)],
 ["1", "x", "x**2"])</pre>
 calculates <code>[c1, c2, c3]</code> such that the function 
 <code>f(x) = c1*1 + c2*x + c3 * x**2</code> 
 approximates the conditions f(-1)=2, f(1)=1, f(2)=1, f(3)=0, f(5)=3 . 
 Result: "(6/5 * 1) + (-53/70 * x) + (3/14 * x**2)"
 ----de----
 fit(data, values) 
 Interpolation mit der Gaußschen Methode der kleinsten Quadrate 
 <code>fit(data, values)</code>
 berechnet nach der Methode der kleinsten Quadrate
 eine Linearkombination der Funktionen, die x auf die
 Terme in <code>values</code> abbilden. 
 (entspricht der Mathematica-Funktion <code>Fit</code>). 
 Beispiel:<pre>
fit([(-1,2), (1,1), (2,1), (3,0), (5,3)],
 ["1", "x", "x**2"])</pre>
 berechnet <code>[c1, c2, c3]</code> so dass die Funktion 
 <code>f(x) = c1*1 + c2*x + c3 * x**2</code> 
 die Bedingungen f(-1)=2, f(1)=1, f(2)=1, f(3)=0, f(5)=3 approximiert. 
 Ergebnis: "(6/5 * 1) + (-53/70 * x) + (3/14 * x**2)"
 """
 # from the example:
 m = len(data) # m = 5
 n = len(values) # n = 3
 functions = [eval("lambda x:{}".format(val)) for val in values]
 A = Matrix.fromFunction(m,n, lambda i,k: functions[k](data[i][0]))
 b = ~Vector.fromFunction(m, lambda i: data[i][1])
 # / 1 -1 1 \ / 2 \
 # | 1 1 1 | | 1 |
 # A = | 1 2 4 | b = | 1 |
 # | 1 3 9 | | 0 |
 # \ 1 5 25 / \ 3 /
 y = A.leastSquares(b) # Solution of the normal system (column matrix):
 s = "" # y = [[6/5], [-53/70], [3/14]]
 for i in y.rowRange(): # Create linear combination:
 s += "({0} * {1})".format(str(y[i,0]), values[i])
 if i < y.height()-1:
 s += " + "
 return s # "(6/5 * 1) + (-53/70 * x) + (3/14 * x**2)"

#---------------------------------------------------------------
# RSA-Chiffrierung

class Rsa: # RSA Encryption # RSA-Verschlüsselung
 """ Class <code>Rsa</code>: RSA Encryption 
 This is a static class (only static methods available).
 ----de----
 Klasse <code>Rsa</code>: RSA-Verschlüsselung 
 Statische Klasse: nur Klassenmethoden.
 """

 _staticMethods = ("createKeys", "encrypt")

 @staticmethod
 def createKeys(bits=768):
 """ createKeys(bits=768) 
 returns an RSA key pair 
 The pair is returned as tuple <code>((e,n),(d,n))</code>. 
 Example: 
 <code>public, private = Rsa.createKeys(24) 
 m = 12345 # message 
 c = Rsa.encrypt(m, public) # encrypted message 
 assert Rsa.encrypt(c, private) == m</code> 
 See also <a href="#Rsa.encrypt"><code>Rsa.encrypt</code></a>.
 ----de----
 createKeys(bits=768) 
 liefert ein RSA-Schlüsselpaar 
 Das Paar wird als Tupel <code>((e,n),(d,n))</code> zurückgegeben. 
 Beispiel: 
 <code>public, private = Rsa.createKeys(24) 
 m = 12345 # Nachricht 
 c = Rsa.encrypt(m, public) # verschlüsselte Nachricht 
 assert Rsa.encrypt(c, private) == m</code> 
 Siehe auch <a href="#Rsa.encrypt"><code>Rsa.encrypt</code></a>.
 """
 m = 2 ** (bits//2 - 1)
 # wähle zwei Primzahlen p, q
 p = nextPrime(m + rand(m))
 q = nextPrime(m + rand(m))
 n = p * q
 phi = (p-1)*(q-1) # phi(n) (Eulersche Phi-Funktion)
 d = phi # wähle d aus 1..phi(n)
 while gcd(d,phi) != 1: # mit d, phi(n) teilerfremd
 d = rand(phi)
 t,x,y = gcdExt(d,phi)
 e = x % phi
 publicKey = (e,n)
 privateKey = (d,n)
 return (publicKey, privateKey)

 @staticmethod
 def encrypt(m, key):
 """ encrypt(m, key) 
 RSA encryption 
 encrypts m using the mit dem key <code>key = (e,n)</code> 
 See also <a href="#Rsa.createKeys"><code>Rsa.createKeys</code></a>.
 ----de----
 encrypt(m, key) 
 RSA-Verschlüsselung 
 Verschlüsselt m mit dem Schlüssel key = (e,n) 
 Siehe auch <a href="#Rsa.createKeys"><code>Rsa.createKeys</code></a>.
 """
 return pow(m, key[0], key[1])

#===============================================================
# Permutationen, Permutationsgruppen
#---------------------------------------------------------------

def permute(v, fn, prefix=[]):
 """ permute(v, fn) 
 applies procedure <code>fn</code> to all permutations of liste v. 
 If v is sorted, the permutations are also generated in sorted order.
 ----de----
 permute(v, fn) 
 wendet Prozedur <code>fn</code> auf alle Permutationen von Liste v an. 
 Wenn v sortiert ist, werden auch die Permutationen in
 sortierter Reihenfolge erzeugt.
 """
 if len(v) == 0: # nichts zu permutieren!
 fn(prefix)
 else:
 for i in v:
 p = prefix + [i] # i an Praefix anhängen
 w = list(v) # Duplikat von w, ...
 w.remove(i) # ... ohne i
 permute(w, fn, p)

def permutations(v, prefix=[]):
 """ permutations(v) 
 List of all permutations of v 
 If v is sorted, the permutations are also generated in sorted order.
 ----de----
 permutations(v) 
 Liste aller Permutationen von v 
 Wenn v sortiert ist, werden auch die Permutationen in
 sortierter Reihenfolge erzeugt.
 """
 if len(v) == 0: # nichts zu permutieren!
 yield prefix
 else:
 for i in v:
 p = prefix + [i] # i an Praefix anhängen
 w = list(v) # Duplikat von w, ...
 w.remove(i) # ... ohne i
 for i in permutations(w, p):
 yield i

class Perm: # Permutations (from dictionary or cycle representation) # Permutationen (aus Dictionary oder Zyklendarstellung)
 """ Class <code>Perm</code> 
 Permutations (from dictionary or cycle representation)
 Examples (three times the same permutation: 
 <code>Perm({1:4, 2:3, 3:2, 4:1})</code> 
 <code>Perm("(1,4)(2,3)")</code> 
 <code>Perm([[1,4],[2,3]])</code>
 ----de----
 Klasse <code>Perm</code> 
 Permutationen (aus Dictionary oder Zyklendarstellung)
 Beispiele (dreimal die gleiche Permutation: 
 <code>Perm({1:4, 2:3, 3:2, 4:1})</code> 
 <code>Perm("(1,4)(2,3)")</code> 
 <code>Perm([[1,4],[2,3]])</code>
 """
 _staticMethods = []

 def __init__(self, p):
 """ Constructor
 ----de----
 Konstruktor
 """
 if type(p) == dict:
 # Zuordnungstafel
 # Beispiel: Perm({1:4, 2:3, 3:2, 4:1})
 self._d = p
 elif type(p) == str:
 # Zyklendarstellung als String
 # Der String darf nur Zahlen enthalten
 # Beispiel: Perm("(1,4)(2,3)")
 v = []
 p = p.replace(" ", "") # Leerzeichen entfernen
 assert p[0] == "(" and p[-1] == ")"
 p = p[1:-1] # Klammern vorn und hinten entfernen
 for cycle in p.split(")("):
 v.append([eval(x) for x in cycle.split(",")])
 self.fromVector(v)
 else:
 # Zyklendarstellung
 # Beispiel: Perm([[1,4],[2,3]])
 self.fromVector(p)
 self.reduce()

 def fromVector(self, v):
 """ p.fromVector(v) 
 Generation from vector v with argument/value pairs. 
 Example: <code>p.fromVector([[1,4],[2,3]])</code>
 ----de----
 p.fromVector(v) 
 Erzeugung aus Vektor v mit Argument-Wert-Paaren. 
 Beispiel: <code>p.fromVector([[1,4],[2,3]])</code>
 """
 assert type(v) == list
 self._d = {}
 for cycle in v:
 assert type(cycle) == list
 for i in range(len(cycle)):
 if i < len(cycle) - 1:
 self._d[cycle[i]] = cycle[i+1]
 else:
 self._d[cycle[i]] = cycle[0]

 def __repr__(self):
 """ Object representation as string (as cycles)
 """
 if len(self._d) == 0:
 return "identity"
 else:
 return self.cycles()

 def __getitem__(self, i):
 return self._d.get(i, i)

 def __mul__(p1, p2):
 """ product of two Permutations (left first, then right!)
 """
 p = {}
 m = p1.members()
 for i in m:
 p[i] = p2[p1[i]]
 for i in p2.members():
 if i not in m:
 p[i] = p2[i]
 return Perm(p)

 def __eq__(p1,p2): return p1.cycles() == p2.cycles()
 def __ne__(p1,p2): return p1.cycles() != p2.cycles()

 def members(self):
 """ p.members() 
 List of the permuted elements
 ----de----
 p.members() 
 Liste der permutierten Elemente
 """
 m = list(self._d.keys())
 m.sort()
 return m

 def inverse(self):
 """ p.inverse() 
 inverse permutation
 ----de----
 p.inverse() 
 Umkehrpermutation
 """
 d = {}
 for i in self._d:
 d[self._d[i]] = i
 return Perm(d)

 def cycles(self):
 """ p.cycles() 
 Cycle representation
 ----de----
 p.cycles() 
 Zyklendarstellung
 """
 s = ""
 m = self.members()
 while len(m) > 0:
 first = m[0]
 m.remove(first)
 i = first
 s += "({0}".format(str(first))
 while self._d[i] != first:
 i = self._d[i]
 m.remove(i)
 s += ",{0}".format(str(i))
 s += ")"
 return s

 def check(self):
 """ p.check() 
 check for consistency
 ----de----
 p.check() 
 Konsistenzprüfung
 """
 m = self.members()
 v = self._d.values()
 v.sort()
 assert(k == v)

 def reduce(self):
 """ p.reduce() 
 removes fixed points
 ----de----
 p.reduce() 
 entfernt Fixpunkte
 """
 id = []
 for i in self._d:
 if self._d[i] == i:
 id.append(i)
 for i in id:
 del self._d[i]

 def inversions(self):
 """ p.inversions() 
 Inversion number (count of "swapped pairs")
 ----de----
 p.inversions() 
 Inversionszahl (Anzahl der "vertauschten Paare")
 """
 m = self.members()
 s = 0
 n = len(m)
 for i in range(n):
 for j in range(i+1, n):
 if self[m[i]] > self[m[j]]:
 s += 1
 return s

 def sign(self):
 """ p.sign() 
 signum 
 1 for even, -1 for odd permutations
 ----de----
 p.sign() 
 Signum 
 1 bzw. -1 für gerade bzw. ungerade Permutation
 """
 return 1 - 2*(self.inversions() % 2)

 def odd(self):
 """ p.odd() 
 odd permutation? (bool)
 ----de----
 p.odd() 
 ungerade Permutation? (bool)
 """
 return self.inversions() % 2 == 1

#===============================================================
# Modulo-Arithmetik und Rechnen in Polynomringen
#---------------------------------------------------------------

class Mod: # Elements of the residue class rings Z[m] # Elemente der Restklassenringe Z[m]
 """ Class <code>Mod</code> 
 Elements of the residue class rings Z[m] 
 <code>Mod(n, m)</code> returns element n des of the residue class ring Z[m] 
 Arithmetic operators may be applied to the objects. 
 Some operations are defined only,
 if Z[m] is a field (<==> m is prime number) 
 See also <a href="mathguide.html#mod"><code>Mod</code></a>
 ----de----
 Klasse <code>Mod</code> 
 Elemente der Restklassenringe Z[m] 
 <code>Mod(n, m)</code> liefert Element n des Restklassenrings Z[m] 
 Die Objekte können mit arithmetischen Operatoren verknüpft werden. 
 Einige Operationen sind nur definiert,
 wenn Z[m] ein Körper ist (<==> m ist Primzahl) 
 Siehe auch <a href="mathguide.html#mod"><code>Mod</code></a>
 """
 _staticMethods = []

 def operators(self):
 """ Mod.operators() 
 For documentation only 
 The following operators are defined in the class Mod: 
 <table>
 <tr><th>Op.</th><th>Function</th><th>Examples</th></tr>
 <tr><th>+</th><td>Addition</td><td><code>a + b; a += b</code></td></tr>
 <tr><th>-</th><td>Subtraction</td><td><code>a - b; a -= b</code></td></tr>
 <tr><th>*</th><td>Multiplication</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>/</th><td>Division</td><td><code>a / b; a /= b</code></td></tr>
 <tr><th>%</th><td>Modulo (in fields only , i.e. modulus is prime)</td><td><code>a & b; a &= b</code></td></tr>
 <tr><th>-</th><td>Unary minus</td><td><code>-a</code></td></tr>
 </table>
 This method is for documentation only. It has no effect.
 ----de----
 Mod.operators() 
 Nur zur Dokumentation 
 Die Klasse Mod erlaubt folgende Operatoren: 
 <table>
 <tr><th>Op.</th><th>Funktion</th><th>Beispiele</th></tr>
 <tr><th>+</th><td>Addition</td><td><code>a + b; a += b</code></td></tr>
 <tr><th>-</th><td>Subtraktion</td><td><code>a - b; a -= b</code></td></tr>
 <tr><th>*</th><td>Multiplikation</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>/</th><td>Division (nur in Körpern, d.h. Primzahlmodul)</td><td><code>a / b; a /= b</code></td></tr>
 <tr><th>%</th><td>Modulo (nur in Körpern, d.h. Primzahlmodul)</td><td><code>a & b; a &= b</code></td></tr>
 <tr><th>-</th><td>unäres Minus</td><td>-a</td></tr>
 </table>
 Diese Methode dient nur zur Dokumentation. Sie hat keine Wirkung.
 """
 pass

 def __init__(self, n, m):
 """ Constructor
 ----de----
 Konstruktor
 """
 #assert isInteger(n) and isNatural(m)
 self.n = n % m
 self.m = m

 def copy(self):
 """ A.copy() 
 independent copy ("clone") of A
 ----de----
 A.copy() 
 unabhängige Kopie von A
 """
 return Mod(self.n, self.m)

 def __int__(self):
 return self.n

 def __repr__(self):
 """ Object representation as string
 """
 return "Mod({0}, {1})".format(str(self.n), str(self.m))

 def __cmp__(A, B):
 if not isinstance(B, Mod):
 B = Mod(B, A.m)
 return _cmp(A.n, B.n)

 def inverse(self):
 """ p.inverse() 
 Multiplicative inverse 
 defined in fields only. Therefore modulus must be prime number!
 ----de----
 p.inverse() 
 Multiplikatives Inverses 
 Nur in Körpern definiert. Deshalb muss der Modul Primzahl sein!
 """
 def gcdExt(a,b):
 if b == 0:
 return a, 1, 0
 else:
 q, r = divmod(a, b)
 d, z, x = gcdExt(b, r)
 return d, x, z - x * q

 d, inv, y = gcdExt(self.n, self.m)
 assert d == 1
 # inv * n + y * m == 1 ==> (inv * n) % m == 1
 assert (inv * self.n) % self.m == 1
 return Mod(inv, self.m)

 def __eq__(a,b):
 """ Test for equality (operator ==)
 with other Mod object or 0
 """
 if type(b) == int and b == 0:
 return a.n == 0
 return (a.n == b.n and a.m == b.m)

 def __abs__(self):
 """ Absolute value
 """
 return self

 def __ne__(a,b):
 """ Test for unequality (Operator !=)
 with other Mod object or 0
 """
 return not a.__eq__(b)

 def __neg__(self):
 """ Unary Operator -
 -a --> a.__neg__()
 """
 return Mod((self.m - self.n) % self.m, self.m)

 def __add__(a, b):
 """ a.__add__(b) 
 Addition 
 Also operator notation: a + b 
 Instead of <code>a = a + b</code> you can write <code>a += b</code>.
 """
 # a + b --> a.__add__(b), wenn a vom Typ Mod ist
 # b.__radd__(a), sonst
 if isinstance(b, Mod):
 d = gcd(a.m, b.m)
 return Mod(a.n + b.n, d)
 else:
 return Mod(a.n + b, a.m)

 def __radd__(b, a):
 """ Addition, see __add__ """
 return b.__add__(a)

 def __sub__(a, b):
 """ Subtraction:
 a - b --> a.__sub__(b), if a is of type Mod
 b.__rsub__(a), otherwise
 """
 if isinstance(b, Mod):
 d = gcd(a.m, b.m)
 return Mod(a.n - b.n, d)
 else:
 return Mod(a.n - b, a.m)

 def __rsub__(b, a):
 """ Subtraction, see __sub__
 """
 return Mod(a - b.m, b.m)

 def __mul__(a, b):
 """ Multiplication:
 a * b --> a.__mul__(b), if a is of type Mod
 b.__rmul__(a), otherwise
 """
 if isinstance(b, Mod):
 d = gcd(a.m, b.m)
 return Mod(a.n * b.n, d)
 else:
 return Mod(a.n * b, a.m)

 def __rmul__(b, a):
 """ Multiplikation, see __mul__
 """
 return b.__mul__(a)

 def __truediv__(a, b):
 """ exact division:
 a / b --> a.__div__(b), if a is of type Mod
 b.__rdiv__(a), otherwise
 """
 if not isinstance(b, Mod):
 b = Mod(b, a.m)
 assert a.m == b.m
 return a * b.inverse()

 def __rtruediv__(a, b):
 """ exact division, see __div__
 """
 assert not isinstance(b, Mod)
 b = Mod(b, a.m)
 return b.__div__(a)

 def __div__(a, b): return a.__truediv__(b)
 def __rdiv__(a, b): return a.__rtruediv__(b)

#===============================================================
# Poly
#---------------------------------------------------------------

class Poly: # Polynomials over arbitrary rings # Polynome über beliebigen Ringen
 """ Class <code>Poly</code> 
 Polynomials 
 <code>Poly(list)</code> returns a polynomial with the coefficients from <code>list</code> (ascending). 
 The coefficients must comply with the ring axioms, 
 for division and modulo operator (//, %, divmod) also the field axioms, 
 and implement the equality with 0. 
 Short form for polynomial with Mod objektes: 
 Example 
 <code>Poly([Mod(1,2),Mod(0,2),Mod(1,2)])</code> 
 Short form: 
 <code>Poly([1,0,1,1],2)</code>
 ----de----
 Klasse <code>Poly</code> 
 Polynome 
 <code>Poly(list)</code> liefert Polynom mit den Koeffizienten aus <code>list</code> (aufsteigend). 
 Die Koeffizienten müssen die Ringaxiome erfüllen, 
 für Division und Modulo-Operator (//, %, divmod) auch Körperaxiome, 
 und die Abfrage auf Gleichheit mit 0 implementieren. 
 Kurzform für Polynom mit Mod-Objekten: 
 Beispiel 
 <code>Poly([Mod(1,2),Mod(0,2),Mod(1,2)])</code> 
 Dafür auch Kurzschreibweise: 
 <code>Poly([1,0,1,1],2)</code>
 """
 #............................................................
 _staticMethods = ("gcd", "gcdExt", "irreducibles")

 @staticmethod
 def gcd(a,b):
 """ gcd(a,b) 
 greatest common divisor of a and b 
 Calculated using the Euclidean algorithm
 ----de----
 gcd(a,b) 
 Größter gemeinsamer Teiler von a und b 
 Berechnung mit dem Euklidischen Algorithmus
 """
 p = a.copy()
 q = b.copy()
 while not q.isZero():
 p, q = q, p % q
 return p.normalized()

 @staticmethod
 def gcdExt(a,b):
 """ gcdExt(a,b) 
 Extended Euclidean algorithm. 
 Returns (d, x, y) with d == gcd(a,b) and x*a + y*b == d
 ----de----
 gcdExt(a,b) 
 Erweiterter Euklidischer Algorithmus. 
 Gibt (d, x, y) zurück mit d == gcd(a,b) und x*a + y*b == d
 """
 assert type(a[0]) == Mod
 m = a[0].m
 if b.isZero():
 return a.copy(), Poly([1],m), Poly([0],m)
 else:
 q, r = divmod(a, b)
 d, z, x = Poly.gcdExt(b, r)
 y = z - x * q
 return d, x, y

 @staticmethod
 def irreducibles(maxDegree, mod=2):
 """ Poly.irreducibles(maxDegree, mod=2) 
 Irreducible polynomials 
 Returns: List of irreducible polynomials
 over Z[mod] up to maxDegree
 ----de----
 Poly.irreducibles(maxDegree, mod=2) 
 Irreduzible Polynome 
 Rückgabewert: Liste der irreduziblen Polynome
 über Z[mod] bis zum Grad maxDegree
 """
 primeList = [Poly([0,1],mod)]
 for i in fromTo(mod+1, mod**(maxDegree+1)-1):
 f = Poly(NumeralSystem.toBase(i, mod, True), mod)
 divisible = False
 for p in primeList:
 if f % p == 0:
 divisible = True
 break
 if not divisible:
 primeList.append(f)
 return primeList
 #............................................................

 def operators(self):
 """ Poly.operators() 
 For documentation only 
 The following operators are defined in the class Poly: 
 <table>
 <tr><th>Op.</th><th>Function</th><th>Examples</th></tr>
 <tr><th>+</th><td>Addition</td><td><code>a + b; a += b</code></td></tr>
 <tr><th>-</th><td>Subtraction</td><td><code>a - b; a -= b</code></td></tr>
 <tr><th>*</th><td>Multiplication</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>/</th><td>Polynomial division (with remainder)</td><td><code>a / b; a /= b</code></td></tr>
 <tr><th>divmod</th><td>Division and remainder simultaneously</td><td><code>d,m = divmod(a,b)</code></td></tr>
 <tr><th>-</th><td>Unary minus</td><td><code>-a</code></td></tr>
 </table>
 This method is for documentation only. It has no effect.
 ----de----
 Poly.operators() 
 Nur zur Dokumentation 
 Die Klasse Poly erlaubt folgende Operatoren: 
 <table>
 <tr><th>Op.</th><th>Funktion</th><th>Beispiele</th></tr>
 <tr><th>+</th><td>Addition</td><td><code>a + b; a += b</code></td></tr>
 <tr><th>-</th><td>Subtraktion</td><td><code>a - b; a -= b</code></td></tr>
 <tr><th>*</th><td>Multiplikation</td><td><code>a * b; a *= b</code></td></tr>
 <tr><th>//</th><td>Polynomdivision (mit Rest)</td><td><code>a // b; a //= b</code></td></tr>
 <tr><th>divmod</th><td>Division und Modulo simultan</td><td><code>d,m = divmod(a,b)</code></td></tr>
 <tr><th>-</th><td>unäres Minus</td><td>-a</td></tr>
 </table>
 Diese Methode dient nur zur Dokumentation. Sie hat keine Wirkung.
 """
 pass
 def __init__(self, coeff=[], mod=0):
 """ Constructor
 ----de----
 Konstruktor
 """
 if mod != 0:
 self.__p = [Mod(i,mod) for i in coeff]
 if len(coeff) == 0: # zero polynomial
 self.__p = [Mod(0, mod)]
 else:
 # make all coefficients rational numbers
 self.__p = [Rational(c) for c in coeff]
 if len(coeff) == 0: # zero polynomial
 self.__p = [Rational(0)]

 # drop leading zero coefficients
 n = len(self.__p)
 while n > 0 and self[n-1] == 0:
 n -= 1
 self.__p = self.__p[0:n]

 def copy(self):
 """ p.copy() 
 independent copy ("clone") of p
 ----de----
 p.copy() 
 unabhängige Kopie von p
 """
 c = Poly()
 c.__p = [copy(x) for x in self.__p]
 return c

 def isZero(self):
 return (len(self.__p) == 0
 or len(self.__p) == 1 and self.__p[0] == 0)

 def isOne(self):
 if len(self.__p) != 1:
 return False
 if type(self.__p[0]) == Mod:
 return self.__p[0].n == 1
 return self.__p[0] == 1

 def degree(self):
 """ p.degree() 
 Grad des Polynoms
 ----de----
 p.degree() 
 degree of the polynomial
 """
 if self.isZero():
 return -1
 else:
 return len(self.__p) - 1

 def value(self, x):
 """ p.value(x) 
 Value of the polynomial for argument x
 ----de----
 p.value(x) 
 Wert des Polynoms für Argument x
 """
 val = 0
 n = self.degree()
 for i in fromTo(0, self.degree()):
 val = (x * val) + self.__p[n-i]
 return val

 def __normalize(self):
 """ reduce to true degree
 """
 n = len(self.__p)
 while n > 0 and self[n-1] == 0:
 n -= 1
 self.__p = self.__p[0:n]

 def toDigits(f):
 """ p.toDigits() 
 String of coefficients 
 Short form, nur für singled digit Mod coefficients: 
 Concatenation of the coefficient digits (without mod). 
 Highest coeffizient left! 
 Example 
 <code>p=Poly([1,0,1,1],2)</code> 
 <code>p.toDigits()</code> 
 Result: <code>"1101"</code>
 ----de----
 p.toDigits() 
 String aus Koeffizienten 
 Kurzform, nur für einstellige Mod-Koeffizienten: 
 Verkettung der Koeffizienten-Ziffern (ohne mod-Angabe). 
 Höchster Koeffizient links! 
 Beispiel 
 <code>p=Poly([1,0,1,1],2)</code> 
 <code>p.toDigits()</code> 
 Ergebnis: <code>"1101"</code>
 """
 assert isinstance(f[0], Mod) and f[0].n < 10
 l = [str(x.n) for x in f.__p]
 l.reverse()
 return join(l)

 def __repr__(f):
 """ Object representation as string
 """
 # prüfen, ob mind. ein Koeffizien vom Typ Mod ist:
 m = 0
 for x in f.__p:
 if isinstance(x, Mod):
 m = x.m
 break
 if m > 0:
 # Voraussetzung: Mod-Objekte kompatibel!
 return "Poly({0}, {1})".format(str([int(x) for x in f.__p]), m)
 else:
 return "Poly([{0}])".format(join(f.__p, ", "))

 def expand(f, x="x"):
 """ p.expand(x="x") 
 Polynomial written for the variable x 
 ----de----
 p.expand(x="x") 
 ausgeschriebenes Polynom für die Variable x 
 """
 def pot(i):
 if i == 0:
 return 1
 elif i == 1:
 return x
 else:
 return "{0}^{1}".format(x, i)
 l = ["({0} * {1})".format(str(f[i]), pot(i)) for i in fromTo(0, f.degree())]
 return join(l, " + ")

 def pretty(f, width=0):
 """ p.pretty(width=0) 
 Human readable written representation 
 ----de----
 p.pretty(width=0) 
 Lesbare Darstellung des Polynoms 
 """
 if f.isZero():
 pp = "0"
 else:
 x = "x"
 def pot(i):
 if i == 0: return ""
 elif i == 1: return x
 elif i == 2: return "{0}²".format(x)
 elif i == 3: return "{0}³".format(x)
 else: return "{0}^{1}".format(x, i)
 def sign(i):
 if f[i] < 0: return "- "
 elif i > 0: return "+ "
 else: return ""
 pp = ""
 for i in fromTo(f.degree(), 0, -1):
 # sign
 sign = ""
 if type(f[i]) != Mod and f[i] < 0: sign = " - "
 elif i < f.degree(): sign = " + "

 # abs value
 a = str(abs(f[i]))
 if type(f[i]) == Mod:
 a = f[i].n
 if i > 0 and type(f[i]) == Mod and f[i].m == 2:
 a = ""

 if f[i] != 0 and pot != "":
 pp += "{0}{1}{2}".format(sign, a, pot(i))
 if len(pp) < width:
 width -= len(pp)
 pp = "{0}{1}{2}".format((width//2)*" ", pp, (width-width//2)*" ")
 return pp

 def __getitem__(f, i):
 """ Definition of the index operator []
 f[i] returns the i-th coefficient of the polynomial.
 """
 if i <= f.degree():
 return f.__p[i]
 else:
 return 0 * f.__p[i]

 #--- unäre Operatoren: ---

 def __neg__(f):
 """ unärer Operator -
 """
 r = Poly()
 r.__p = [-copy(x) for x in f.coeff]
 return r

 #--- binäre Operatoren: ---

 def __eq__(f,g):
 """ Test for equality (operator ==)
 with other polynomial or ring element
 """
 if not isinstance(g, Poly):
 g = Poly([g])
 return f.__p == g.__p

 def __ne__(f,g):
 """ Test for unequality (operator !=)
 with other polynomial or ring element
 """
 return not f.__eq__(g)

 def __add__(f, g):
 """ a.__add__(b) 
 Addition 
 Also operator notation: a + b 
 Instead of <code>a = a + b</code> you can write <code>a += b</code>.
 """
 # a + b --> a.__add__ (b), wenn a vom Typ Poly ist
 # b.__radd__(a), sonst
 if not isinstance(g, Poly):
 g = Poly([g])
 fDeg, gDeg = f.degree(), g.degree()
 minDeg = min(fDeg, gDeg)
 maxDeg = max(fDeg, gDeg)
 s = Poly()
 s.__p = (maxDeg + 1) * [0] # list with (maxDeg + 1) zeroes
 for i in fromTo(0, maxDeg):
 if i > fDeg:
 s.__p[i] = copy(g.__p[i])
 elif i > gDeg:
 s.__p[i] = copy(f.__p[i])
 else:
 s.__p[i] = f.__p[i] + g.__p[i]
 s.__normalize()
 return s

 def __radd__(f, g):
 """ Addition, see __add__
 """
 return f + g

 def __sub__(f, g):
 """ Subtraction (operators - and -=):
 a - b --> a.__sub__(b), if a is of type Poly
 b.__rsub__(a), otherwise
 """
 if not isinstance(g, Poly):
 g = Poly([g])
 fDeg, gDeg = f.degree(), g.degree()
 minDeg = min(fDeg, gDeg)
 maxDeg = max(fDeg, gDeg)
 s = Poly()
 s.__p = (maxDeg + 1) * [0] # list with (maxDeg + 1) zeroes
 for i in fromTo(0, maxDeg):
 if i > fDeg:
 s.__p[i] = -g.__p[i]
 elif i > gDeg:
 s.__p[i] = copy(f.__p[i])
 else:
 s.__p[i] = f.__p[i] - g.__p[i]
 s.__normalize()
 return s

 def __rsub__(f, g):
 """ Subtraktion, see __sub__
 """
 return Poly([g]) - f

 def __mul__(f, g):
 """ Multiplication (operators * and *=):
 a * b --> a.__mul__(b), if a is of type Poly
 b.__rmul__(a), otherwise
 """
 if not isinstance(g, Poly):
 g = Poly([g])
 if f.isZero():
 return f.copy()
 if g.isZero():
 return g.copy()
 n, m = f.degree(), g.degree()

 p = Poly()
 p.__p = (n + m + 1) * [0] # list with (n + m + 1) zeroes
 for i in fromTo(0, n):
 for j in fromTo(0, m):
 p.__p[i+j] += f.__p[i] * g.__p[j]
 p.__normalize()
 return p

 def __rmul__(f, g):
 """ Multiplikation, see __mul__
 """
 return f * g

 def __rshift__(f, n):
 """ f >> n -- Multiplication of Polynomial by x**n
 """
 r = Poly()
 r.__p = n * [0] + [copy(x) for x in f.__p]
 return r

 def __divmod__(f, g):
 """ Polynomial division with remainder. (global function divmod)
 divmod(f, g) --> f.__divmod__(g), if f is of type Poly
 g.__rdivmod__(f), otherwise

 Returns the pair (q,r)
 with: f == g * q + r
 and: r.degree() < g.degree() or r == 0
 """
 if not isinstance(g, Poly):
 g = Poly([g])
 r = f.copy() # Kopie: Zählerpolynom wird zu Restpolynom
 q = Poly() # Divisionsergebnis: beginnen mit Nullpolynom
 while r.degree() >= g.degree() and not r.isZero():
 assert f == g * q + r
 qCoeff = r[r.degree()] / g[g.degree()] # Quotient der hoechsten Koeff.
 qPoly = Poly()
 qPoly.__p = [qCoeff]
 qPoly >>= (r.degree() - g.degree())
 q += qPoly
 assert f == g * (q-qPoly) + r
 deg0 = r.degree()
 r -= qPoly * g # grad(r) muss jetzt kleiner werden
 assert r.degree() < deg0
 assert f == g * (q-qPoly) + (r + qPoly*g)
 assert f == g * q + r
 return q, r

 def __rdivmod__(f, g):
 """ Polynomial division with remainder. See __divmod__
 """
 assert not isinstance(g, Poly)
 g = Poly([g])
 return g.__divmod__(f)

 def __div__(f, g):
 """ Polynomial division with remainder (operators // und //=):
 f // g --> f.__floordiv__ (g), if f is of type Poly
 g.__rfloordiv__(f), otherwise
 """
 return divmod(f, g)[0]

 def __rdiv__(f, g):
 """ Polynomial division with remainder. See __floordiv__
 """
 assert not isinstance(g, Poly)
 g = Poly([g])
 return divmod(g, f)[0]

 def __truediv__(a, b): return a.__div__(b)
 def __rtruediv__(a, b): return a.__rdiv__(b)
 def __floordiv__(a, b): return a.__div__(b)
 def __rfloordiv__(a, b): return a.__rdiv__(b)

 #-- Modulo-Operator ---
 # a & b --> __mod__ (a, b), wenn a ein Poly-Objekt ist
 # __rmod__(b, a), sonst

 def __mod__(f, d):
 """ Modulo operator for polynomials (operators % und %=):
 a % b --> a.__mod__ (b), if f is of type Poly
 b.__rmod__(a), otherwise
 """
 return divmod(f, d)[1]

 def __rmod__(f, g):
 """ Modulo operator for polynomials. See __mod__
 """
 assert not isinstance(g, Poly)
 g = Poly([g])
 return divmod(g, f)[1]

 def __pow__(a, x):
 assert isNatural(x)
 if x == 1:
 return a
 elif x % 2 == 1:
 return a * a**(x-1)
 else:
 t = a ** (x//2)
 return t * t

 def normalized(self):
 """ Divides the polynomial by the highest coefficient
 ----de----
 Dividiert das Polynom durch den höchsten Koeffizienten
 """
 if self.degree() >= 0:
 r = self.copy()
 r.__p = [c/r.__p[-1] for c in r.__p]
 return r
 else:
 return self.copy()

 def inverse(self, mod):
 """ p.inverse(m) 
 Multiplicative inverse modulo polynomial m 
 The polynomial mod must be relative prime to p!
 ----de----
 p.inverse(m) 
 Multiplikatives Inverses modulo Polynom m 
 Das Polynom mod muss teilerfremd zu p sein!
 """
 assert type(self.__p[0]) == Mod
 m = self.__p[0].m
 d, inv, y = Poly.gcdExt(self, mod)
 if not d.isOne():
 raise _MathGUIdeError("Polynomial is not invertible")
 # inv * self + y * mod == d == 1
 # ==> (inv * self) % mod == 1
 return inv % mod

#===============================================================
# Set
#---------------------------------------------------------------

class Set (list): # Polymorphic sets # Polymorphe Mengen
 """ Class <code>Set</code> 
 Polymorphic sets 
 Example: <code>A = Set([1,2,3])</code> 
 The class Set allows the following operations: 
 <ul>
 <li>Intersection (operator &)</li>
 <li>Union (operator +)</li>
 <li>Set difference (operator -)</li>
 <li>Symmetric difference (operator ^)</li>
 <li>Power set (method powerset)</li>
 <li>Cartesian product (Operator *)</li>
 <li>Element test (Operator <code>in</code>)</li>
 <li>Cardinality (Python function <code>len</code>)</li>
 </ul>
 ----de----
 Klasse <code>Set</code> 
 Polymorphe Mengen 
 Beispiel: <code>A = Set([1,2,3])</code> 
 Die Klasse Set erlaubt folgende Operationen: 
 <ul>
 <li>Durchschnitt (Operator &)</li>
 <li>Vereinigung (Operator +)</li>
 <li>Mengendifferenz (Operator -)</li>
 <li>Symmetrische Differenz (Operator ^)</li>
 <li>Potenzmenge (powerset)</li>
 <li>Kartesisches Produkt (Operator *)</li>
 <li>Elementtest (Operator <code>in</code>)</li>
 <li>Kardinalzahl (Python-Funktion <code>len</code>)</li>
 </ul>
 """
 _staticMethods = []

 def operators(self):
 """ Set.operators() 
 For documentation only 
 The following operators are defined in the class Set: 
 <table>
 <tr><th>Op.</th><th>Function</th><th>Examples</th></tr>
 <tr><th>&</th><td>Intersection</td><td><code>A & B; A &= B</code></td></tr>
 <tr><th>+</th><td>Union</td><td><code>A + B; A += B</code></td></tr>
 <tr><th>-</th><td>Set difference</td><td><code>A - B; A -= B</code></td></tr>
 <tr><th>^</th><td>Symmetric difference</td><td><code>A ^ B; A ^= B</code></td></tr>
 <tr><th>*</th><td>Cartesian product</td><td><code>A * B; A *= B</code></td></tr>
 <tr><th>in</th><td>is element of</td><td><code>x in A</code></td></tr>
 </table>
 This method is for documentation only. It has no effect.
 ----de----
 Set.operators() 
 Nur zur Dokumentation 
 Die Klasse Set erlaubt folgende Operatoren: 
 <table>
 <tr><th>Op.</th><th>Funktion</th><th>Beispiele</th></tr>
 <tr><th>&</th><td>Durchschnitt</td><td><code>A & B; A &= B</code></td></tr>
 <tr><th>+</th><td>Vereinigung</td><td><code>A + B; A += B</code></td></tr>
 <tr><th>-</th><td>Mengendifferenz</td><td><code>A - B; A -= B</code></td></tr>
 <tr><th>^</th><td>Symmetrische Differenz</td><td><code>A ^ B; A ^= B</code></td></tr>
 <tr><th>*</th><td>Kartesisches Produkt</td><td><code>A * B; A *= B</code></td></tr>
 <tr><th>in</th><td>ist Element von</td><td><code>x in A</code></td></tr>
 </table>
 Diese Methode dient nur zur Dokumentation. Sie hat keine Wirkung.
 """
 pass

 def __init__(self, el):
 """ Constructor
 ----de----
 Konstruktor
 """
 list.__init__(self, list(el))
 self.__normalize()

 def __normalize(self):
 self.sort()
 for i in range(len(self)-1, 0, -1):
 if self[i] == self[i-1]:
 del self[i]

 def __repr__(self):
 """ Object representation as string (in curly braces)
 """
 r = "{"
 for e in self:
 r += str(e) + ","
 return (r+"}").replace(",}","}")

 # In der Basisklasse list sind die "rich comparison"
 # Methoden definiert und haben Vorrang vor __cmp___
 # Deshalb werden sie hier wieder auf __cmp__ geleitet:
 def __lt__(A, B):
 return A.__cmp__(B) < 0

 def __le__(A, B):
 return A.__cmp__(B) <= 0

 def __eq__(A, B):
 return A.__cmp__(B) == 0

 def __ne__(A, B):
 return A.__cmp__(B) != 0

 def __gt__(A, B):
 return A.__cmp__(B) > 0

 def __ge__(A, B):
 return A.__cmp__(B) >= 0

 def __cmp__(A, B):
 """ For ordering sets
 Elements, wihch are sets themself come first and are sorted by length (cardinality)

 simultaneous overloading of < <=, >, >=, == and !=

 If A is a Set and b not, A.__cmp__(b) is called for comparing A with b (or vice versa)!

 """
 if not isinstance(B, Set):
 return -1 # erst Mengen, dann Nichtmengen
 elif len(A) != len(B): # zwei Mengen: nach Größe ordnen
 return _cmp(len(A), len(B))
 elif len(A) == 0:
 return 0 # zwei leere Mengen
 else: # zwei gleich große nichtleere Mengen
 if A[0] != B[0]:
 return _cmp(A[0], B[0])
 else:
 return _cmp(Set(A[1:]), Set(B[1:]))

 def insert(self, e):
 """ A.insert(e) 
 Inserts an additional element ti the set 
 A.insert(e) is equivalent to A += Set(e)
 ----de----
 A.insert(e) 
 Fügt ein zusätzliches Element in die Menge ein 
 A.insert(e) ist äquivalent zu A += Set(e)
 """
 self.append(e)
 self.__normalize()

 def issubset(A, B):
 """ A.issubset(B) 
 Is A subset of (oder equal to) B? 
 Return value: bool
 ----de----
 A.issubset(B) 
 Ist A Teilmenge (oder gleich) B? 
 Rückgabewert: bool
 """
 for e in A:
 if e not in B:
 return False
 return True

 def issuperset(A, B):
 """ A.issuperset(B) 
 Is A superset of (oder equal to) B? 
 Return value: bool
 ----de----
 A.issuperset(B) 
 Ist A Obermenge (oder gleich) B 
 Rückgabewert: bool
 """
 return B.issubset(A)

 def __add__(A, B):
 """ Union (operators + and +=):
 a + b --> __add__ (a, b), if a is of type Set
 __radd__(b, a), otherwise
 """
 return Set(list(A) + list(B))

 def union(A,B):
 """ A.union(B) 
 Union set 
 Also operator notation: A + B 
 Instead of <code>A = A + B</code> you can write <code>A += B</code>.
 ----de----
 A.union(B) 
 Vereinigungsmenge 
 Auch Operatorschreibweise: A + B 
 Statt A = A + B kann auch A += B geschrieben werden.
 """
 return A + B

 def __sub__(A, B):
 """ Set difference (operators - and -=)
 """
 result = list(A)
 for e in B:
 if e in result:
 result.remove(e)
 return Set(result)

 def difference(A,B):
 """ A.difference(B) 
 Set difference 
 Also operator notation: A - B 
 Instead of <code>A = A - B</code> you can write <code>A -= B</code>.
 ----de----
 A.difference(B) 
 Mengendifferenz 
 Auch Operatorschreibweise: A - B 
 Statt A = A - B kann auch A -= B geschrieben werden.
 """
 return A - B

 def __and__(A, B):
 """ Intersection (operators & and &=)
 """
 result = list(A)
 for e in A:
 if e not in B:
 result.remove(e)
 return Set(result)

 def intersection(A,B):
 """ A.intersection(B) 
 Intersection 
 Also operator notation: A - B 
 Instead of <code>A = A & B</code> you can write <code>A &= B</code>.
 ----de----
 A.intersection(B) 
 Schnittmenge 
 Auch Operatorschreibweise: A & B 
 Statt A = A & B kann auch A &= B geschrieben werden.
 """
 return A & B

 def __xor__(A, B):
 """ Symmetric difference (operators ^ and ^=)
 """
 return (A + B) - (A & B)

 def __pow__(A, B):
 """ Symmetric difference (operators ^ and ^=)
 """
 return (A + B) - (A & B)

 def symmDiff(A,B):
 """ A.symmDiff(B) 
 Symmetric difference 
 Also operator notation: A ^ B 
 Instead of <code>A = A ^ B</code> you can write <code>A ^= B</code>.
 ----de----
 A.symmDiff(B) 
 Symmetrische Differenz 
 Auch Operatorschreibweise: A ^ B 
 Statt A = A ^ B kann auch A ^= B geschrieben werden.
 """
 return A^B

 def __mul__(A, B):
 """ Cartesian product (operator *)
 """
 return Set([(x,y) for x in A for y in B])

 def product(A,B):
 """ A.product(B) 
 Cartesian product 
 Also operator notation: A * B 
 Instead of <code>A = A * B</code> you can write <code>A *= B</code>.
 ----de----
 A.product(B) 
 Kartesisches Produkt 
 Auch Operatorschreibweise: A * B 
 Statt A = A * B kann auch A *= B geschrieben werden.
 """
 return A*B

 def powerset(self):
 """ A.powerset() 
 Powerset (set of alle subsets) 
 ----de----
 A.powerset() 
 Potenzmenge (Menge aller Teilmengen) 
 """
 if len(self) == 0:
 return Set([Set([])])
 else:
 first = self[0]
 rest = self[1:]
 powRest = Set(list(rest)).powerset()
 pl = [Set([first])+x for x in powRest]
 return powRest + Set(pl)

 def card(self):
 """ A.card() 
 Cardinality 
 Same as len(A)
 ----de----
 A.card() 
 Kardinalzahl 
 Synonym zu len(A)
 """
 return len(self)

###"Repräsentative Mengen"
##
## # 1 2 12
##A1 = Set(1,12) # + - +
##A2 = Set(2,12) # - + +
##
## # 1 2 3 12 13 23 123
##B1 = Set(1,12,13,123) # + - - + + - +
##B2 = Set(2,12,23,123) # - + - + - + +
##B3 = Set(3,13,23,123) # - - + - + + +
##
##C1 = Set(1,12,13,14, 123, 124, 134, 1234)
##C2 = Set(2,12,23,24, 123, 124, 234, 1234)
##C3 = Set(3,13,23,34, 123, 134, 234, 1234)
##C4 = Set(4,14,24,34, 124, 134, 234, 1234)
##
##D1 = Set(1, 12,13,14,15, 123,124,125,134,135,145, 1234,1235,1245,1345, 12345)
##D2 = Set(2, 12,23,24,25, 123,124,125,234,235,245, 1234,1235,1245,2345, 12345)
##D3 = Set(3, 13,23,34,35, 123,134,135,234,234,345, 1234,1235,1345,2345, 12345)
##D4 = Set(4, 14,24,34,45, 124,134,145,234,245,345, 1234,1245,1345,2345, 12345)
##D5 = Set(5, 15,25,35,45, 125,135,145,235,245,345, 1235,1245,1345,2345, 12345)

#===============================================================
# Propositional Logic
#---------------------------------------------------------------

class Logic: # propositional logic # Aussagenlogik
 """ Class <code>Logic</code>: Propositional Logic 
 This is a static class (only static methods available) 
 Formulae of propositional logic, passed to the static methods,
 must be composed as follows: 
 Variables must be single letters. 
 Junctors: "<code>not</code>", "<code>and</code>", "<code>or</code>", "<code>xor</code>", "<code>-></code>", "<code><-></code>".
 ----de----
 Klasse <code>Logic</code>: Aussagenlogik 
 Statische Klasse: nur Klassenmethoden 
 Aussagenlogische Formeln, die den Klassenmethoden übergeben werden,
 müssen wie folgt aufgebaut sein: 
 Variablen müssen einzelne Buchstaben sein. 
 Junktoren: "<code>not</code>", "<code>and</code>", "<code>or</code>", "<code>xor</code>", "<code>-></code>", "<code><-></code>".
 """
 _staticMethods = ("truthTables", "satisfies", "valid", "satisfiable", "implies", "equivalent")

 @staticmethod
 def truthTables(formula):
 """ truthTables(formula) 
 All possible truth tables for variables in formula 
 The result is returned as a generator object (for a single iteration). 
 Example 
 <code>for t in Logic.truthTables("A and B"): print(t)</code> 
 returns: 
 <code>{'A': False, 'B': False} 
 {'A': False, 'B': True} 
 {'A': True, 'B': False} 
 {'A': True, 'B': True}</code>
 ----de----
 truthTables(formula) 
 Alle möglichen Wahrheitswerttabellen für im Ausdruck
 formula vorkommenden Variablen 
 Das Ergebnis wird als Generatorobjekt (zur einmaligen Iteration) ausgegeben. 
 Beispiel 
 <code>for t in Logic.truthTables("A and B"): print(t)</code> 
 gibt folgendes aus: 
 <code>{'A': False, 'B': False} 
 {'A': False, 'B': True} 
 {'A': True, 'B': False} 
 {'A': True, 'B': True}</code>
 """
 table = {}
 for op in ("not", "and", "or", "xor"):
 formula = formula.replace(op, " ")
 for c in formula:
 if c.isalpha(): # Variable!
 table[c] = False
 n = len(table)
 vars = list(table.keys())
 while True:
 yield table
 k = n-1
 table[vars[k]] = not table[vars[k]]
 while not table[vars[k]]: # "Übertrag"
 k -= 1
 if k < 0:
 return
 table[vars[k]] = not table[vars[k]]

 @staticmethod
 def satisfies(f, truthTable):
 """ satisfies(f, truthTable) 
 Is the formula f true under the truthTable? 
 Example 
 <code>Logic.satisfies("A or B", {"A":True, "B":False})</code>
 returns <code>True<code>.
 ----de----
 satisfies(f, truthTable) 
 Wird die Formel f unter der Belegung truthTable wahr? 
 Beispiel 
 <code>Logic.satisfies("A or B", {"A":True, "B":False})</code>
 gibt <code>True</code> zurück.
 """
 f = f.replace("not", "@0").replace("and", "@1").replace("or", "@2").replace("xor", "@3")
 for var in truthTable:
 f = f.replace(var, str(truthTable[var]))
 f = f.replace("@0", "not").replace("@1", "and").replace("@2", "or").replace("@3", "xor")
 return eval(f)

 @staticmethod
 def valid(f):
 """ satisfies(f, truthTable) 
 Is the formula f valid? 
 Example 
 <code>Logic.valid("A or not A")</code>
 returns <code>True<code>.
 ----de----
 satisfies(f, truthTable) 
 Ist die Formel f allgemeingültig? 
 Beispiel 
 <code>Logic.valid("A or not A")</code>
 gibt <code>True</code> zurück.
 """
 f = f.replace("->", "<=").replace("<->", "==").replace("xor", "!=").replace("<<=", "==")
 for t in Logic.truthTables(f):
 if not Logic.satisfies(f, t):
 return False
 return True

 @staticmethod
 def satisfiable(f):
 """ satisfiable(f) 
 Is the formula f satisfiable? 
 Example 
 <code>Logic.satisfiable("A and not A")</code>
 returns <code>False<code>.
 ----de----
 satisfiable(f) 
 Ist die Formel f erfüllbar? 
 Beispiel 
 <code>Logic.satisfiable("A and not A")</code>
 gibt <code>False</code> zurück.
 """
 return not Logic.valid("not({0})".format(f))

 @staticmethod
 def implies(f1, f2):
 """ implies(f1, f2) 
 implies the formula f1 the formula f2?? 
 Example 
 <code>Logic.implies("A and B", "A or B")</code>
 returns <code>True<code>.
 ----de----
 implies(f1, f2) 
 Folgt aus der Formel f1 die Formel f2? 
 Beispiel 
 <code>Logic.implies("A and B", "A or B")</code>
 gibt <code>True</code> zurück.
 """
 return Logic.valid("not({0}) or ({1})".format(f1, f2))

 @staticmethod
 def equivalent(f1, f2):
 """ equivalent(f1, f2) 
 Are the formulae f1 and f2 equivalent?? 
 Example 
 <code>Logic.equivalent("not (A and B)", "not A or not B"))</code>
 returns <code>True<code>.
 ----de----
 equivalent(f1, f2) 
 Sind die Formeln f1 und f2 äquivalent? 
 Beispiel 
 <code>Logic.equivalent("not (A and B)", "not A or not B"))</code>
 gibt <code>True</code> zurück.
 """
 return Logic.valid("({0}) == ({1})".format(f1, f2))

#===============================================================
# Graphs and Trees
#---------------------------------------------------------------

class Edge: # Edges of a Graph # Kanten eines Graphen
 """ Class <code>Edge</code>: Edge of a <a href="Graph.html">Graph</a> 
 ----de----
 Klasse <code>Edge</code>: Kante eines Graphen (Klasse <a href="Graph.html">Graph</a>) 
 """
 _staticMethods = []

 def __init__(self, source, target, weight=0, highlight=False):
 """ Constructor
 ----de----
 Konstruktor
 """
 self.source = source
 self.target = target
 self.weight = weight
 self.highlight = highlight

 def __repr__(self):
 """ Object representation as string
 """
 if self.weight == 0:
 return "{{{0},{1}}}".format(self.source, self.target)
 else:
 return "({{{0},{1}}}, w={2})".format(self.source, self.target, self.weight)

 @staticmethod
 def _fromShortRepr(s):
 p = s.split(",")
 n0, n1 = p[0:2]
 w = 0
 m = len(p) >= 4
 if len(p) >= 3:
 w = p[2]
 return Edge(int(n0), int(n1), int(w), m)

 def isIncidentNode(self, node):
 return node in (self.source, self.target)

 def adjacentNode(self, node):
 if node == self.source:
 return self.target
 elif node == self.target:
 return self.source
 else:
 return None

 def _shortRepr(self):
 """ short representation as string
 """
 if self.highlight:
 return "{0},{1},{2},1".format(self.source, self.target, self.weight)
 elif self.weight != 0:
 return "{0},{1},{2}".format(self.source, self.target, self.weight)
 else:
 return "{0},{1}".format(self.source, self.target)

class Graph: # Graphs # Graphen
 """ Class <code>Graph</code>: Graph (in the sense of graph theory) 
 A Graph consists of nodes and edges.
 An edge (class <a href="Edge.html">Edge</a>) connects two nodes. 
 Edges optionally may have a weight. 
 Although a pure Graph does not carry geometrical data, mathGUIde allows to attach
 positions to nodes corresponding to a Graph representation in the Graph window.
 ----de----
 Klasse <code>Graph</code>: Graph (im Sinne der Graphentheorie) 
 Ein Graph besteht aus Knoten und Kanten.
 Eine Kante (Klasse <a href="Edge.html">Edge</a>) verbindet zwei Knoten. 
 Kanten können auch mit einem Gewicht versehen werden. 
 Obwohl zur Graphendefinition keine geometrischen Informationen gehören,
 kann mathGUIde die Positionen von im Graphenfenster angeordneten Knoten speichern.
 """
 #............................................................
 _staticMethods = ("fromShortRepr", "example1", "complete", "completeBipartite", "petersen")

 @staticmethod
 def fromShortRepr(s):
 """ Graph.fromShortRepr(s) 
 Makes the Graph from the short representation s
 ----de----
 Graph.fromShortRepr(s) 
 Erzeugt den Graphen aus der Kurzdarstellung s
 """
 g = Graph()
 parts = s.split("|")
 if len(parts[0]) > 0:
 # nodes
 for n in parts[0].split(";"):
 p = n.split(",")
 g.nodes.append(int(p[0]))
 if len(p) >= 3:
 g.locations[int(p[0])] = float(p[1]), float(p[2])
 if len(parts[1]) > 0:
 # edges
 for e in parts[1].split(";"):
 g.edges.append(Edge._fromShortRepr(e))
 return g

 @staticmethod
 def example1():
 """ Graph.example1()
 ----de----
 Graph.example1()
 """
 return Graph.fromShortRepr("1,80,20;2,240,20;3,0,100;4,160,100;5,320,100;6,0,180;7,160,180;8,320,180;9,80,260;10,240,260|1,2,36;1,4,15;2,4,69;2,5,1;4,5,53;3,4,92;1,3,99;3,6,37;6,7,21;7,8,88;4,6,27;4,7,32;5,7,10;5,8,66;6,9,3;7,9,18;7,10,49;8,10,43;9,10,74")

 @staticmethod
 def complete(n, r=130, margin=20):
 """ Graph.complete(n, r=130, margin=20) 
 Complete Graph of order n, displayed on a circle with radius r.
 ----de----
 Graph.complete(n, r=130, margin=20) 
 Vollständiger Graph der Ordnung n angeordnet auf einem Kreis mit Radius r.
 """
 c = margin + r # center x and y
 nodes = ";".join(["{0},{1},{2}".format(i+1, c+r*sin(i*2*pi/n), c-r*cos(i*2*pi/n)) for i in range(n)])
 edges = ";".join(["{0},{1},{2}".format(i+1,j+1,0) for i in range(n) for j in range(n) if i < j])
 return Graph.fromShortRepr("{0}|{1}".format(nodes, edges))

 @staticmethod
 def completeBipartite(n1, n2):
 """ Graph.completeBipartite(n1,n2) 
 Complete bipartite Graph with partition sizes n1 and n2.
 ----de----
 Graph.completeBipartite(n1,n2) 
 Vollständiger bipartiter Graph mit den Partitionsgrößen n1 und n2.
 """
 dx, dy, margin = 60, 160, 20
 margin1 = margin2 = margin
 if n1 < n2:
 margin1 += (n2-n1) * dx / 2
 elif n2 < n1:
 margin2 += (n1-n2) * dx / 2
 nodes1 = ";".join(["{0},{1},{2}".format(i+1, margin1+i*dx, margin) for i in range(n1)])
 nodes2 = ";".join(["{0},{1},{2}".format(n1+i+1, margin2+i*dx, margin+dy) for i in range(n2)])
 edges = ";".join(["{0},{1},{2}".format(i+1,n1+j+1,0) for i in range(n1) for j in range(n2)])
 return Graph.fromShortRepr("{0};{1}|{2}".format(nodes1, nodes2, edges))

 @staticmethod
 def petersen(r=150, margin=20):
 """ Graph.petersen() 
 The Peterson Graph 
 It serves as an example and counterexample for many problems in graph theory.
 ----de----
 Graph.completeBipartite() 
 Der Peterson-Graph 
 Er dient als Beispiel oder Gegenbeispiel für viele Probleme der Graphentheorie.
 """
 c = margin + r # center x and y
 r0 = r/2
 nodes1 = ";".join(["{0},{1},{2}".format(i+1, c+r*sin(i*2*pi/5), c-r*cos(i*2*pi/5)) for i in range(5)])
 nodes2 = ";".join(["{0},{1},{2}".format(i+6, c+r0*sin(i*2*pi/5), c-r0*cos(i*2*pi/5)) for i in range(5)])
 edges = "1,2;1,5;1,6;2,3;2,7;3,4;3,8;4,5;4,9;5,10;6,8;6,9;7,9;7,10;8,10"
 return Graph.fromShortRepr("{0};{1}|{2}".format(nodes1, nodes2, edges))

 #............................................................

 def __init__(self):
 """ Constructor
 ----de----
 Konstruktor
 """
 self.nodes = []
 self.edges = []
 self.locations = {}

 def __repr__(self):
 """ Object representation as string
 """
 return "Graph(Nodes = {{{0}}}, Edges = {{{1}}})".format(
 join([str(n) for n in self.nodes], ","),
 join([str(e) for e in self.edges], ","))

 def shortRepr(self):
 """ G.shortRepr() 
 short representation of Graph G as string 
 ----de----
 .shortRepr() 
 Kurzdarstellung des Graphen G (als String)
 """
 # example: "1;2;3;4|1,2,36;1,3,99;1,4,15;2,4,69"
 return "{0}|{1}".format(join([self._shortNodeRepr(n) for n in self.nodes], ";"),
 join([e._shortRepr() for e in self.edges], ";"))

 def _shortNodeRepr(self, node):
 """ short representation of a Node as string
 """
 if node in self.locations:
 return "{0},{1},{2}".format(node, self.locations[node][0], self.locations[node][1])
 else:
 return "{0}".format(node)

 def order(self):
 """ G.order() 
 Order (number of vertices) of the Graph G.
 ----de----
 G.order() 
 Ordnung (Anzahl der Knoten) des Graphen G.
 """
 return len(self.nodes)

 def complement(self):
 """ G.complement() 
 Returns the graph with the same vertices as G, such that two vertices are adjacent if and only if they are not adjacent in G.
 ----de----
 G.complement() 
 Gibt den Graphen mit den gleichen Knoten wie G zurück, in dem sich zwischen zwei Knoten genau dann eine Kante befindet, wenn G dort keine Kante hat.
 """
 g = Graph()
 g.nodes = list(self.nodes) # independent copy
 g.edges = []
 g.locations = dict(self.locations)
 for n1 in g.nodes:
 for n2 in g.nodes:
 if n1 < n2:
 if not self.adjacent(n1, n2):
 g.addEdges(Edge(n1, n2))
 return g

 def connectedNodes(self, node):
 """ G.componentsCount() 
 Returns a list of all nodes of G which are connected with node (including node).
 ----de----
 G.componentsCount() 
 Gibt eine Liste aller mit node durch einen Kantenzug verbundenen Knoten in G zurück (einschließlich node).
 """
 checked = 0
 nodes = [node]
 while checked < len(nodes):
 n = nodes[checked]
 for m in self.adjacentNodes(n):
 if (m not in nodes):
 nodes.append(m)
 checked += 1
 nodes.sort()
 return nodes

 def componentsCount(self):
 """ G.componentsCount() 
 The number connectivity components of the Graph G
 ----de----
 G.componentsCount() 
 Anzahl der Zusammenhangskomponenten des Graphen G
 """
 count = 0
 uncheckedNodes = list(self.nodes) # independent copy
 while len(uncheckedNodes) > 0:
 count += 1
 n = uncheckedNodes[0]
 for m in self.connectedNodes(n):
 uncheckedNodes.remove(m)
 return count

 def isConnected(self):
 """ G.isConnected() 
 Returns a truth value indicating wether the Graph G is connected
 ----de----
 G.isConnected() 
 Gibt als Wahrheitswert zurück, ob der Graph G zusammenhängend ist
 """
 if len(self.nodes) == 0: # an empty graph is connected!
 return True
 else:
 return (self.componentsCount() == 1)

 def isTree(self):
 """ G.isTree() 
 Returns a truth value indicating wether the Graph G is a tree
 ----de----
 G.isTree() 
 Gibt als Wahrheitswert zurück, ob der Graph G ein Baum ist
 """
 return self.isConnected() and len(self.edges) == len(self.nodes) - 1

 def adjacentNodes(self, node):
 """ G.adjacentNodes(node) 
 Returns a sorted list of all nodes of G adjacent to node.
 ----de----
 G.adjacentNodes(node) 
 Gibt eine sortierte Liste aller zu node adjazenten Knoten von G zurück.
 """
 adjNodes = []
 for e in self.edges:
 if node == e.source:
 adjNodes.append(e.target)
 elif node == e.target:
 adjNodes.append(e.source)
 adjNodes.sort()
 return adjNodes

 def adjacent(self, node1, node2):
 """ G.adjacentNodes(node1, node2) 
 Determines wether node1 and node2 are adjacent (bool return value) in G.
 ----de----
 G.adjacentNodes(node1, node2) 
 Gibt an, ob node1 und node2 in G adjazent sind (bool).
 """
 return node2 in self.adjacentNodes(node1)

 def incidentEdges(self, node):
 incEdges = []
 for e in self.edges:
 if e.isIncidentNode(node):
 incEdges.append(e)
 return incEdges

 def degree(self, node):
 """ G.degree(node) 
 Returns the degree of node (count of adjacent nodes)
 ----de----
 G.degree(node) 
 Gibt den Grad des Knotens node (Anzahl der adjazenten Knoten) in G zurück
 """
 return len(self.incidentEdges(node))

 def edge(self, node1, node2):
 """ G.edge(node1, node2) 
 The connecting edge between the two nodes 
 If node1 and node2 are adjacent, the connecting edge is returned, otherwise <code>None</code>.
 ----de----
 G.edge(node1, node2) 
 Die Verbindungskante zwischen den zwei Knoten 
 Wenn node1 und node2 adjazent sind, wird die the Verbindungskante zurückgegeben, sonst <code>None</code>.
 """
 for e in self.edges:
 if (e.source == node1 and e.target == node2) or (e.source == node2 and e.target == node1):
 return e
 return None

 def display(self, n=0):
 """ G.display(n=0) 
 Displays the Graph G in the n-th Graph tab. 
 if n == 0, the Graph is displayed in the active Graph window.
 ----de----
 G.display(n=0) 
 Zeigt den Graphen G im n-ten Graphen-Fenster an.
 Ist n = 0, wird der Graph im aktiven Graphen-Fenster angezeigt.
 """
 return "@G{0}:{1}".format(n, self.shortRepr())

 def addNodes(self, *nodes):
 self.nodes += nodes

 def addEdges(self, *edges):
 for e in edges:
 if e.source not in self.nodes:
 self.nodes.append(e.source)
 if e.target not in self.nodes:
 self.nodes.append(e.target)
 self.edges += edges

 def copy(self):
 """ G.copy() 
 independent copy ("clone") of G
 ----de----
 G.copy() 
 unabhängige Kopie von G
 """
 return Graph.fromShortRepr(self.shortRepr())

 def reduceToHighlightedEdges(self):
 """ G.reduceToHighlightedEdges() 
 Removes all not highlighted edges from the Graph G. 
 The highlighting of the remaining edges is turned off.
 ----de----
 G.reduceToHighlightedEdges() 
 Entfernt alle nicht hervorgehobenen Kanten aus dem Graphen G. 
 Die Hervorhebung der verbleibenden Kanten wird zurückgenommen.
 """
 highlightedEdges = [e for e in self.edges if e.highlight]
 otherEdges = [e for e in self.edges if not e.highlight]
 for e in highlightedEdges:
 e.highlight = False
 for e in otherEdges:
 self.edges.remove(e)

 def cancelHighlight(self):
 """ G.cancelHighlight() 
 Resets all highlighted edges to normal state.
 ----de----
 G.cancelHighlight() 
 Setzt alle hervorgehobenen Kanten in den Normalzustand zurück.
 """
 for e in self.edges:
 e.highlight = False

 def highlightShortestPath(self, startNode, goalNode):
 """ G.highlightShortestPath(startNode, goalNode) 
 Highlightes a shortest path (if existent) from startNode to goalNode using the Dijkstra algorithm. 
 Returns the length of the path (or inf if there is no connection). 
 Use the method <a href="Graph.display"><code>display</code></a> in order
 to make the highlighting visible.
 ----de----
 G.highlightShortestPath(startNode, goalNode) 
 Hebt einen kürzesten Pfad von startNode nach goalNode hervor (mit dem Dijkstra-Algorithmus). 
 Die Länge des Pfades wird zurückgegeben (oder inf, falls keine Verbindung vorhanden). 
 Mit der Methode <a href="Graph.display"><code>display</code></a> können Sie
 die Hervorhebung sichtbar machen.
 """
 self.cancelHighlight()
 predecessorNode = {}
 unfinishedNodes = list(self.nodes) # independent copy

 upperBound = {} # dictionary with upper bounds for path length from startNode to node
 infinite = float("inf")
 for n in self.nodes:
 upperBound[n] = infinite
 upperBound[startNode] = 0

 activeNode = startNode
 while True:
 unfinishedNodes.remove(activeNode)
 fail = True
 for e in self.incidentEdges(activeNode):
 n = e.adjacentNode(activeNode)
 if n in unfinishedNodes:
 fail = False
 # add 1 to weights in order to count path lengths
 # for unweighted Graphs (default weight = 0):
 dist = upperBound[activeNode] + e.weight + 1
 if dist < upperBound[n]: # a shorter path to n discovered
 upperBound[n] = dist
 predecessorNode[n] = activeNode
 if fail:
 break
 best, newNode = min([(upperBound[n], n) for n in unfinishedNodes])
 activeNode = newNode

 if goalNode in unfinishedNodes:
 return infinite
 # go back from goalNode to startNode:
 n = goalNode
 length = 0
 while n != startNode:
 p = predecessorNode[n]
 e = self.edge(p, n)
 e.highlight = True
 length += e.weight
 n = p
 return length

 def highlightMinimalSpanTree(g, trace=False):
 """ G.minimalSpanTree(trace=False) 
 Highlightes a minimal spanning tree using the Kruskal algorithm. 
 If <code>trace</code> is <code>True</code>, all steps are printed.
 ----de----
 G.minimalSpanTree(trace=False) 
 Hebt einen minimalen aufspannenden Baum hervor (mit dem Kruskal-Algorithmus). 
 Wenn <code>trace</code> den Wert <code>True</code> hat, werden alle Schritte angezeigt.
 """

 nodeColor = {}
 marked = {}
 for n in g.nodes:
 nodeColor[n] = 0
 def weight(edge):
 return edge.weight
 g.edges.sort(key=weight)

 curColor = 1
 spanTreeSize = 0

 for edge in g.edges:
 marked[edge] = False

 for edge in g.edges:
 choose = True
 n0, n1 = edge.source, edge.target
 if nodeColor[n0] == nodeColor[n1]:
 # both nodes have tha same color
 if nodeColor[n0] == 0:
 nodeColor[n0] = curColor # both nodes uncolored (no edge):
 nodeColor[n1] = curColor # use current color,
 curColor += 1 # new acurrent color
 else:
 choose = False # coonecting would result in a circle ==> skip edge
 else:
 # both nodes have different colors
 if min(nodeColor[n0], nodeColor[n1]) == 0:
 # one node uncolored: use color of other node
 if nodeColor[n0] == 0:
 nodeColor[n0] = nodeColor[n1]
 else:
 nodeColor[n1] = nodeColor[n0]
 else:
 # both nodes colored: join two components (change all nodes with lower color to higher color)
 lowColor = min(nodeColor[n0], nodeColor[n1])
 highColor = max(nodeColor[n0], nodeColor[n1])
 for n in nodeColor:
 if nodeColor[n] == lowColor:
 nodeColor[n] = highColor
 marked[edge] = choose
 if trace:
 if choose:
 spanTreeSize += 1
 print("{0} Choose".format(spanTreeSize), end=" ")
 else:
 print("Skip", end=" ")
 print("Edge {0}-{1} (Weight {2})".format(edge.source, edge.target, edge.weight))
 if spanTreeSize == len(nodeColor) - 1:
 break
 for edge in marked:
 edge.highlight = marked[edge]

 def isomorphic(G,H):
 """ G.isomorphic(H) 
 Determines wether G and H are isomorphic Graphs. 
 ----de----
 G.isomorphic(H) 
 Entscheidet, ob G und H isomorphe Graphen sind. 
 """
 # Very inefficient implementation!
 # As a first optimization we could permute only nodes of equal degree
 g, h = G.nodes, H.nodes
 n = len(g)
 if n != len(h):
 return False
 if n > 8:
 raise _MathGUIdeError("max. 8 nodes!")
 for p in permutations(g):
 match = True
 for i in fromTo(0, n-2):
 for k in fromTo(i+1, n-1):
 if G.adjacent(p[i], p[k]) != H.adjacent(h[i], h[k]):
 match = False
 break
 if not match:
 break
 if match:
 return True
 return False

#===============================================================

def printTable(height, width, cellFunction, **options):
 """ printTable(height, width, cellFunction, rowHeadFn=None, colHeadFn=None, corner="") 
 displays a table of <var>height</var> rows and <var>width</var> columns, whose cells are calculated by <var>cellFunction</var>. 
 <var>cellFunction</var> must be a function with two parameters (row, column). 
 The following optional named parameters are supported: 
 <var>rowHeadFn</var>: if given, row header cells are added and calculated using this function 
 <var>colHeadFn</var>: if given, column header cells are added and calculated using this function 
 <var>headAlign</var>: "left" (default), "right" or "center" 
 <var>cellAlign</var>: "left", "right" (default) or "center" 
 <var>title</var>: if given, a title line is added, using this text 
 <var>corner</var>: if given (and rowHeadFn and colHeadFn are given), this text is displayed in the top left table corner. 
 This function is used by the menu command Insert -- Table. 
 Example:<pre>
printTable(7, 7, lambda i,k: i*k % 7,
 rowHeadFn = lambda i: i,
 colHeadFn = lambda k: k,
 title = "Multiplication mod 7",
 corner = "*")</pre>
 prints a multiplication table modulo 7.
 ----de----
 printTable(height, width, cellFunction, rowHeadFn=None, colHeadFn=None, corner="") 
 gibt eine Tabelle mit <var>height</var> Zeilen und <var>width</var> Spalten aus, deren Zellen mit <var>cellFunction</var> berechnet werden. 
 <var>cellFunction</var> muss eine Funktion mit zwei Parametern (Zeile, Spalte) sein. 
 Folgende optionalen benannten Parameter werden unterstützt: 
 <var>rowHeadFn</var>: Falls vorhanden, werden Zeilenköpfe angezeigt und mit dieser Funktion berechnet. 
 <var>colHeadFn</var>: Falls vorhanden, werden Spaltenköpfe angezeigt und mit dieser Funktion berechnet. 
 <var>headAlign</var>: "left" (Standard), "right" oder "center" 
 <var>cellAlign</var>: "left", "right" (Standard) oder "center" 
 <var>title</var>: Falls vorhanden, wird eine Titelzeile mit diesem Text über der Tabelle angezeigt 
 <var>corner</var>: Falls vorhanden (und auch rowHeadFn und colHeadFn), wird dieser Text
 in der oberen linken Tabellenecke angezeigt. 
 Diese Funktion wird vom Menübefehl Einfügen -- Tabelle verwendet. 
 Beispiel:<pre>
printTable(7, 7, lambda i,k: i*k % 7,
 rowHeadFn = lambda i: i,
 colHeadFn = lambda k: k,
 title = "Multiplikation mod 7",
 corner = "*")</pre>
 Gibt eine Multiplikationstabelle modulo 7 aus.
 """
 rowHeadFn = options.get("rowHeadFn", None)
 colHeadFn = options.get("colHeadFn", None)
 headAlign = options.get("headAlign", "left")
 cellAlign = options.get("cellAlign", "right")
 title = options.get("title", "")
 titleAlign = options.get("titleAlign", "center")
 corner = options.get("corner", "")

 tdHead = '<td bgcolor="#C6D9F1" align="{0}">'.format(headAlign)
 tdValue = '<td bgcolor="#F2F7FC" align="{0}">'.format(cellAlign)
 print('{html}<table cellspacing="6">') # don't convert to HTML again !
 if len(title) > 0:
 colspan = width
 if rowHeadFn:
 colspan = width + 1
 print('<tr><td colspan="{0}" align="{1}" bgcolor="#376092">{2}</td></tr>'.format(colspan, titleAlign, title))
 if colHeadFn:
 print("<tr>", end="") # heading
 if rowHeadFn:
 print("<td>{0}</td>".format(corner), end="") # top left corner
 for k in range(width):
 print(tdHead, colHeadFn(k), "</td>", end="")
 print("</tr>")
 for i in range(height):
 print("<tr>", end="") # heading
 if rowHeadFn:
 print(tdHead, rowHeadFn(i), "</td>", end="")
 for k in range(width):
 print(tdValue, cellFunction(i,k), "</td>", end="")
 print("</tr>")
 print("</table>{/html}")

def printValueTable(values, variable="n", start=0, end=10, increment=1, vertical=True):
 """ printValueTable(values, variable="n", start=0, end=10, increment=1, vertical=True) 
 Displays a value table. 
 This function is used by the menu command Insert -- Value Table. 
 Example:<pre>
printValueTable([(lambda n:n^2, "n²"), (lambda n:2^n, "2n")],
 "n", 0, 10, True)</pre>
 displays the values of n² and 2n for n from 0 to 10.
 ----de----
 printValueTable(values, variable="n", start=0, end=10, increment=1, vertical=True) 
 Gibt eine Wertetafel aus. 
 Diese Funktion wird vom Menübefehl Einfügen -- Wertetafel verwendet. 
 Beispiel:<pre>
printValueTable([(lambda n:n^2, "n²"), (lambda n:2^n, "2n")],
 "n", 0, 10, True)</pre>
 Listet die Werte n² und 2n für n von 0 bis 10 auf.
 """

 titles = [variable] + [v[1] for v in values]
 functions = [lambda n:n] + [v[0] for v in values]

 tdHead = '<td bgcolor="#C6D9F1" align="left">'
 tdValue = '<td bgcolor="#F2F7FC" align="right">'

 def args():
 n = start
 if increment > 0:
 while n <= end:
 yield n
 n += increment
 else:
 while n >= end:
 yield n
 n += increment

 print('{html}<table cellspacing="6">') # don't convert to HTML again !
 if vertical:
 print("<tr>", end="") # heading
 for t in titles:
 print(tdHead, t, "</td>", end="")
 print("</tr>")

 for n in args():
 print("<tr>", end="") # heading
 for f in functions:
 print(tdValue, f(n), "</td>", end="")
 print("</tr>")
 else:
 for i in range(len(functions)):
 print("<tr>", end="")
 print(tdHead, titles[i], "</td>", end="")
 for n in args():
 print(tdValue, functions[i](n), "</td>", end="")
 print("</table>{/html}")

def plot(terms, x0, x1, **options):
 termList = [terms]
 if type(terms) in (tuple, list):
 termList = terms
 nSteps = options.get("nSteps", 0)
 dots = options.get("dots", False)
 if nSteps == 0 and isInteger(x0) and isInteger(x1) and x0 < x1:
 bars = options.get("bars", True)
 lines = options.get("lines", False)
 nSteps = x1 - x0
 else:
 nSteps = options.get("nSteps", 200)
 bars = options.get("bars", False)
 lines = options.get("lines", True)
 optStr = ""
 if dots: optStr += "d"
 if bars: optStr += "b"
 if not lines: optStr += "n"
 if nSteps == x1 - x0:
 d = 1
 else:
 d = (x1-x0)/float(nSteps)
 val = []
 for i in range(nSteps+1):
 x = x0 + i*d
 try:
 s = str(float(eval(termList[0])))
 except:
 s = "0"
 for j in range(1, len(termList)):
 try:
 s += "&" + str(float(eval(termList[j])))
 except:
 s += "&0"
 val.append(s)
 mark = ["{0},{1}".format(p[0], p[1]) for p in options.get("mark", [])]
 return "@P{0}|{1}|{2}|{3}|{4}".format(optStr, str(x0), str(x1), ";".join(val), ";".join(mark))

def plotFunctions(functions, x0, x1, **options):
 if type(functions) not in (tuple, list):
 functions = [functions]
 nSteps = options.get("nSteps", 0)
 dots = options.get("dots", False)
 if nSteps == 0 and isInteger(x0) and isInteger(x1) and x0 < x1:
 bars = options.get("bars", True)
 lines = options.get("lines", False)
 nSteps = x1 - x0
 else:
 nSteps = options.get("nSteps", 200)
 bars = options.get("bars", False)
 lines = options.get("lines", True)
 optStr = ""
 if dots: optStr += "d"
 if bars: optStr += "b"
 if not lines: optStr += "n"
 if nSteps == x1 - x0:
 d = 1
 else:
 d = (x1-x0)/float(nSteps)
 val = []
 for i in range(nSteps+1):
 x = x0 + i*d
 try:
 s = str(float(functions[0](x)))
 except:
 s = "0"
 for j in range(1, len(functions)):
 try:
 s += "&" + str(float(functions[j](x)))
 except:
 s += "&0"
 val.append(s)
 mark = ["{0},{1}".format(p[0], p[1]) for p in options.get("mark", [])]
 return "@P{0}|{1}|{2}|{3}|{4}".format(optStr, str(x0), str(x1), ";".join(val), ";".join(mark))

#===============================================================
# Introspection
# (internal functions for the mathGUIde GUI)
#---------------------------------------------------------------

def _x_getFunc(glob, loc):
 all = [x.strip() for x in glob.keys() if not "_" in x]
 classes = [x+"C" for x in all if str(glob[x]).startswith("<class") ]
 functions = [x+"F" for x in all if str(glob[x]).startswith("<function") ]
 all = classes + functions
 all.sort(key=str.lower)
 return "L" + ",".join(all)

def _x_getMeth(c, glob, loc):
 s_static = [x.strip() for x in dir(eval(c, glob, loc)) if x[0] == "_"]
 all = [x.strip() for x in dir(eval(c, glob, loc)) if x[0] != "_"]
 typeDescr = str(type(eval(c, glob, loc)))
 elements = []
 if typeDescr == "<class 'type'>":
 # class: display static methods
 try:
 elements = [m+"F" for m in eval(c+"._staticMethods")]
 except:
 pass
 elif typeDescr[:6] == "<class":
 # object: display members
 for x in all:
 typeStr = str(type(eval(c+"."+x, glob, loc)))
 if "function" not in typeStr:
 if "'method'" in typeStr:
 elements.append(x+"F")
 else:
 elements.append(x+"E")
 else:
 # global functions
 elements = [x+"()" for x in all if "function" in str(type(eval(c+"."+x, glob, loc)))]
 elements.sort(key=str.lower)
 return "L" + ",".join(elements)

def _x_getFnPar(fn, glob, loc):
 fn = eval(fn, glob, loc)
 descr = str(fn)
 s = inspect.getdoc(fn)
 return "S" + s

def _x_getMethPar(obj, fn, glob, loc):
 fn = eval(obj + "." + fn, glob, loc)
 s = inspect.getdoc(fn)
 return "S" + s

def _x_parse(s):
 op = s[-1]
 s = s[:-1]
 if op == "(": # \w --> word characters [a-zA-Z0-9_]
 m = re.search(r"\w+\s*\.\s*\w+\s*$", s) # \s --> whitespace
 if m:
 return [x.strip() for x in s[m.start():].split(".")] + ["("]
 m = re.search(r"\w+\s*$", s)
 if m:
 return ["", s[m.start():].strip(), op]

def _contextInfo(s, glob, loc):
 """ Erwartet den Inhalt der Cursorzeile links vom Cursor.
 Rückgabewert:
 - erstes Zeichen "S": danach Info-String für Statuszeile
 - erstes Zeichen "L": danach durch Komma getrennte Methoden
 """
 if s == "@@@": # globale Klassen und Funktionen:
 try:
 r = _x_getFunc(glob, loc)
 except:
 r = "***"
 return r

 p = _x_parse(s)
 try:
 if p[-1] == ".":
 r = _x_getMeth(p[-2], glob, loc)
 elif p[-1] == "(":
 if p[-3] != "":
 r = _x_getMethPar(p[-3], p[-2], glob, loc)
 else:
 r = _x_getFnPar(p[-2], glob, loc)
 except:
 r = "***"
 return r

def _x_classMethInfo(s, glob, loc):
 c, m, p = _x_parse(s)
 assert p == "("
 if c == "": # global
 return "global.{0}".format(m)
 else:
 t = str(type(eval(c, glob, loc)))
 if not t == "<class 'type'>":
 c = t[18:-2] # e.g. "<class 'mathguide.Matrix'>" --> "Matrix"
 return "{0}.{1}".format(c, m)

#end