Klasse zur Primfactor-zerlegung und zum primzahlen errechnen

[nuss]

so hab noch'n paar Sachen verbessert:
http://paste.pocoo.org/show/89434/

@EyDu:
Hab jetzt in der Funktion get_new_primes proof und cache auf None gesetzt.
Die anderen Funktionen können ruhig die errechneten Primzahlen weiterverwenden.


Wiie kriegst du das raus, obs Tabs oder Spaces sind ?
Wenn ich python -tt ausführe, das Modul exportiere,
und alle Funktionen einmal aufrufe, gibt er mir jedenfalls keine Fehler aus.

man könnte die funktionen get_factors und isprime so vereinen:

Code: Alles auswählen

def _prime_ops(number, split, cache=[2]):
    
    proof = _get_first(cache)
    factors = []
    counter = 0

    while True:
        try:
            prime = cache[counter]
        except IndexError:

            cache, proof = get_new_primes(proof+1, proof, cache)
            continue

        if number % prime is 0:
            if split is False:
                return (False, cache)
            else:
        	    factors.append(prime)
        	    number /= prime
        elif number / prime < prime:
            if split is False:
                return (True, cache)
            else:
                if number is not 1:
            	    factors.append(number)
                return (factors, cache)
        else:
        	counter += 1

def get_factors(number, cache=[2]):
    return _prime_ops(number, True, cache)

def isprime(number, cache=[2]):
    return _prime_ops(number, False, cache)

[name]

Schmerz! Need a Medic!
Du ueberpruefst Werte mit is, das ist ganz ganz boese. "is" prueft Identitaet, "==" den Wert.
PS: Nicht das ich euren Spass verderben will, aber wie waer es damit:

Code: Alles auswählen

def get_primes(max_):
    sieve = [True] * (max_ + 1)
    sieve[0] = sieve[1] = False
    for cur_number in xrange(2, max_ + 1):
        if sieve[cur_number]:
            j = 2 * cur_number
            while j <= max_:
                sieve[j] = False
                j += cur_number
            yield cur_number

[Qubit]

Schmerz! Need a Medic!
Du ueberpruefst Werte mit is, das ist ganz ganz boese. "is" prueft Identitaet, "==" den Wert.
PS: Nicht das ich euren Spass verderben will, aber wie waer es damit:

Ein Thema: Primzahlen finden.
Anderes Thema: Primfaktoren einer Zahl finden

Das erste ist simpel und dafür hast du eine mögliche Lösung geboten.
Das andere ist alles andere als simple. Dafür hast du keine Lösung geboten und alle anderen hier sind ziemlich "unperformant".
Um das zu sehen, nehme man mal so "einfache" Zahlen wie:
348234823482389479
oder
3482348234823894793453453455


Natürlich kann man darüber streiten, ob Python überhaupt für solche Kategorien numerischer Probleme geeignet ist. Aber mit einem halbwegs vernünftigen Algorithmus lässt sich auch hier etwas rausholen (mal mit den obigen Zahlen testen)

Code: Alles auswählen

import math
def primtest(z):
        """
        Primfaktorenzerlegung
        (c)2008 Qubit
        hannover.post@gmail.com
        """
        
        primfak=[]
        m = 0
        found = 0
        limit = int((math.sqrt(z)-1)/2)
        while z % 2 == 0:
                z = z/2
                primfak.append((2,z))
                print "found:",primfak[-1]
        i = (z-1)/2
        while (i-1) % 3 == 0:
                i = (i - 1)/3
                primfak.append((3,2*i+1))
                print "found:",primfak[-1]
        if (i-2) % 3 == 0:
                u = (i-2)/3
                s = 2
        if (i-3) % 3 == 0:
                u = (i-3)/3
                s = 3
        while m < u and m < limit:
                if s == 2:
                        if not found and (u-m) % (6*m+5) == 0:
                                n = (u-m) / (6*m+5)
                                m0 = 3*m + 2
                                n0 = 3*n
                                primfak.append((2*m0+1,2*n0+1))
                                print "found:",primfak[-1]
                                found = 1
                        if not found and (u-5*m-5) % (6*m+7) == 0:
                                n = (u-5*m-5) / (6*m+7)
                                m0 = 3*m + 3
                                n0 = 3*n + 2
                                primfak.append((2*m0+1,2*n0+1))
                                print "found:",primfak[-1]
                                found = 1
                if s == 3:
                        if not found and (u-5*m-3) % (6*m+5) == 0:
                                n = (u-5*m-3) / (6*m+5)
                                m0 = 3*m + 2
                                n0 = 3*n + 2
                                primfak.append((2*m0+1,2*n0+1))
                                print "found:",primfak[-1]
                                found = 1
                        if not found and (u-m) % (6*m+7) == 0:
                                n = (u-m) / (6*m+7)
                                m0 = 3*m + 3
                                n0 = 3*n
                                primfak.append((2*m0+1,2*n0+1))
                                print "found:",primfak[-1]
                                found = 1
                if found == 1:
                        i = n0
                        while (i-1) % 3 == 0:
                                i = (i - 1)/3
                                primfak.append((3,2*i+1))
                                print "found:",primfak[-1]
                        if (i-2) % 3 == 0:
                                u = (i-2)/3
                                s = 2
                        if (i-3) % 3 == 0:
                                u = (i-3)/3
                                s = 3
                        limit = int((math.sqrt(2*i+1)-1)/2)
                        found = 0
                        m -= 1
                m += 1

        if len(primfak) > 0:
                x = [x[0] for x in primfak]+[primfak[-1][1]]
                if x[-1] == 1: x.pop()
        else: x=[z]
        print "Primfaktoren: ",sorted(x)

[name]

Code: Alles auswählen

....
....

Scheint zu funktionieren, schaut aber eher, naja, nicht sehr uebersichtlich aus. Koenntest du den Algorithmus bitte einfach so, bzw mit Pseudocode erklaeren?

[Qubit]

Code: Alles auswählen

....
....

Scheint zu funktionieren, schaut aber eher, naja, nicht sehr uebersichtlich aus. Koenntest du den Algorithmus bitte einfach so, bzw mit Pseudocode erklaeren?

Nun, die mathematischen Hintergründe halte ich in diesem Forum für "beyond the scope", nur soviel, es ist etwas "Selbstgestricktes" und soll leicht portierbar sein

Habe es noch mal als Klasse implementiert..

Code: Alles auswählen

import math,time

class primetest(object):
        """
        Primfaktorenzerlegung, Primzahlentest, Primzahlensuche
        (c)2008 Qubit
        hannover.post@gmail.com
        """
        __time = {
                'fak': 0,
                'is_prime': 0,
                'primes': 0
                }

        def __getattr__(self,attr):
                attrs = {
                        't_fak': self.__time['fak'],
                        't_is_prime': self.__time['is_prime'],
                        't_primes': self.__time['primes']
                        }
                get = attrs.get(attr,None)
                if not get == None:
                        return tuple(str('%.3f,sec' %(get)).split(','))
                else:
                        object.__getattribute__(self,attr)
                
        def __fak(self,z,mode):
                primfak=[]
                m = 0
                found = 0
                limit = int((math.sqrt(z)-1)/2)
                while z % 2 == 0:
                        if mode == 1 and not z == 2: return [z,False]
                        z = z/2
                        primfak.append((2,z))
                i = (z-1)/2
                while (i-1) % 3 == 0:
                        if mode == 1 and not i == 1: return [z,False]
                        i = (i - 1)/3
                        primfak.append((3,2*i+1))
                if (i-2) % 3 == 0:
                        u = (i-2)/3
                        s = 2
                if (i-3) % 3 == 0:
                        u = (i-3)/3
                        s = 3
                while m < u and m < limit:
                  if s == 2:
                          if (u-m) % (6*m+5) == 0:
                                n = (u-m) / (6*m+5)
                                m0 = 3*m + 2
                                n0 = 3*n
                                primfak.append((2*m0+1,2*n0+1))
                                found = 1
                          if not found and (u-5*m-5) % (6*m+7) == 0:
                                n = (u-5*m-5) / (6*m+7)
                                m0 = 3*m + 3
                                n0 = 3*n + 2
                                primfak.append((2*m0+1,2*n0+1))
                                found = 1
                  if s == 3:
                          if (u-5*m-3) % (6*m+5) == 0:
                                n = (u-5*m-3) / (6*m+5)
                                m0 = 3*m + 2
                                n0 = 3*n + 2
                                primfak.append((2*m0+1,2*n0+1))
                                found = 1
                          if not found and (u-m) % (6*m+7) == 0:
                                n = (u-m) / (6*m+7)
                                m0 = 3*m + 3
                                n0 = 3*n
                                primfak.append((2*m0+1,2*n0+1))
                                found = 1
                  if found == 1:
                          if mode == 1 and not 2*m0+1 == z:
                                  return [z,False]
                          i = n0
                          while (i-1) % 3 == 0:
                                i = (i - 1)/3
                                primfak.append((3,2*i+1))
                          if (i-2) % 3 == 0:
                                u = (i-2)/3
                                s = 2
                          if (i-3) % 3 == 0:
                                u = (i-3)/3
                                s = 3
                          limit = int((math.sqrt(2*i+1)-1)/2)
                          found = 0
                          m -= 1
                  m += 1

                if len(primfak) > 0:
                        x = [x[0] for x in primfak]+[primfak[-1][1]]
                        if x[-1] == 1: x.pop()
                else: x=[z]
                return x

        def __is_prime(self,z): return len(self.__fak(z,mode=1)) == 1
        def __primes(self,bis,von): return filter(self.__is_prime, range(von,bis+1))

        def fak(self,z=2):
                t_start = time.clock()
                res = self.__fak(z,mode=0)
                t_elapsed = time.clock() - t_start
                self.__time['fak'] = t_elapsed
                return res
        def is_prime(self,z=2):
                t_start = time.clock()
                res = self.__is_prime(z)
                t_elapsed = time.clock() - t_start
                self.__time['is_prime'] = t_elapsed
                return res
        def primes(self,bis=2,von=2):
                t_start = time.clock()
                res = self.__primes(bis,von)
                t_elapsed = time.clock() - t_start
                self.__time['primes'] = t_elapsed
                return res


>>> prim=primetest()
>>> prim.is_prime(348234823482389479)
False
>>> prim.t_is_prime
('0.350', 'sec')
>>> prim.fak(348234823482389479)
[413243, 842687773253L]
>>> prim.t_fak
('2.057', 'sec')
>>> prim.is_prime(123456789)
False
>>> prim.t_is_prime
('0.000', 'sec')
>>> prim.fak(123456789)
[3, 3, 3607, 3803]
>>> prim.t_fak
('0.001', 'sec')
>>> res = prim.primes(10000)
>>> prim.t_primes
('0.142', 'sec')
>>> res[-10:]
[9887, 9901, 9907, 9923, 9929, 9931, 9941, 9949, 9967, 9973]
>>> len(res)
1229
>>> 

[audax]

Ich erzittere vor deinen Mathkenntnissen. Trotzdem ist der Code hässlich wie dar Tag und die Klasse mehr als überflüssig.

Vor allem ist in dem Code irre viel Wiederholung und keine einzige Erklärung.

[name]

Ich erzittere vor deinen Mathkenntnissen. Trotzdem ist der Code hässlich wie dar Tag und die Klasse mehr als überflüssig.

Vor allem ist in dem Code irre viel Wiederholung und keine einzige Erklärung.

Amen. Wuerde mich echt voll intressieren wie das so schnell geht.

[Qubit]

Ich erzittere vor deinen Mathkenntnissen. Trotzdem ist der Code hässlich wie dar Tag und die Klasse mehr als überflüssig.

Warte erst mal ab, bis die richtigen Zahlentheoretiker kommen
Ja, der Code ist nicht hübsch, aber er ist portable und erfüllt seinen Zweck.
Die Klasse ist auch eher "standalone", mir fallen nicht allzu viele (praktische) Anwendungen ein, wo man Primzahlen bräuchte.
Aber zB bei Brüchen zum Kürzen :

Code: Alles auswählen

class ratio(object):
    def __init__(self,zaehler,nenner):
        if nenner == 0:
            raise ZeroDivisionError
        else:
            abs_n = abs(int(nenner))
            sign_n = int(nenner) /abs_n            
        if zaehler == 0:
            sign_z = 1
            r =(0,abs_n)
        else:
            abs_z = abs(int(zaehler))
            sign_z = int(zaehler)/abs_z
            r = self.__reduce(abs_z,abs_n)
        self.z = sign_z*divint(r[0])
        self.n = sign_n*divint(r[1])
    def __add__(x,y):
        return divint(x.z*y.n+y.z*x.n)/divint(x.n*y.n)
    def __sub__(x,y):
        return divint(x.z*y.n-y.z*x.n)/divint(x.n*y.n)
    def __mul__(x,y):
        return divint(x.z*y.z)/divint(x.n*y.n)
    def __rmul__(x,y):
        return divint(x.z*y)/divint(x.n)
    def __pow__(x,p):
        return divint(x.z**p)/divint(x.n**p)
    def __div__(x,y):
        return divint(x.z*y.n)/divint(x.n*y.z)
    def __repr__(self):
        return str(self.z) + '/' + str(self.n)
    def __mini(self,x,y):
        if len(x) <= len(y):
            return x[:],x,y
        else: return y[:],y,x
    def __reduce(self,a,b):
        a = primetest().fak(a)
        b = primetest().fak(b)
        x=self.__mini(a,b)
        for i in x[0]:
            if i in x[2]:
                x[1].remove(i)
                x[2].remove(i)
        x[1].append(1)
        x[2].append(1)
        return reduce(int.__mul__,a),reduce(int.__mul__,b)

class divint(int):
    def __init__(self,x):
        self.z = x
        self.n = 1
    def __div__(x,y):
            return ratio(x,y)

>>> r1 = ratio(2,4)
>>> r2 = ratio(8,12)
>>> r3 = divint(12)/divint(15)
>>> r1,r2,r3
(1/2, 2/3, 4/5)
>>> r1+r2+r3
59/30
>>> r1-r2+r3
19/30
>>> r2+r1*r3
16/15
>>> r1*r2**2
2/9
>>> 

[Leonidas]

Python 2.6 bringt einen Bruch-Datentyp in der Stdlib mit.

[BlackJack]

@Qubit: Beim Anwenden sieht man ganz gut, das die Klasse überflüssig ist. Wenn man ein Exemplar nur erzeugt, um darauf eine Methode auf zu rufen und es dann einfach wieder verwirft, ist das ein deutlicher Hinweis, dass die Klasse überflüssig ist.

[HWK]

@Qubit: Die Zeiten sind ja erst einmal beeindruckend. Sie hängen aber doch ganz erheblich von den Ausgangszahlen ab. Bei z.B. '3165771122567177' sieht das Ergebnis deutlich schlechter aus.
Ich habe als Vergleich einmal die Pollard-Rho-Methode implementiert. Die Ergebnisse sind natürlich auch hier von den Ausgangszahlen abhängig, insgesamt aber sicher nicht schlechter. Die Implementierung ist aber wohl verständlicher.

Code: Alles auswählen

from math import sqrt
from random import randint

def ggt(a, b):
    a, b = abs(a), abs(b)
    while b:
        a, b = b, a % b
    return a

def rho(n, x0=None):
    if x0 is None:
        x0 = randint(2, int(sqrt(n)))
    while True:
        d = randint(1, 1000)
        x = y = x0
        nn = n
        result = []
        while nn > 1:
            while True:
                x = (x * x + d) % nn
                t = (y * y + d) % nn
                y = (t * t + d) % nn
                g = ggt(abs(y - x), nn)
                if g > 1:
                    result.append(g)
                    nn /= g
                    break
        if result[0] != n:
            result.sort()
            return result

if __name__ == '__main__':
    from time import clock
    for n in (34823482348238947, 348234823482389479, 3165771122567177):
        t = clock()
        print rho(n)
        print clock() - t

Ergebnis:

Code: Alles auswählen

[11L, 29576089L, 107038193L]
0.0694772659654
[413243L, 842687773253L]
23.8080616391
[29576089L, 107038193L]
0.258115486745

Ergebnis bei Deiner Variante mit denselben Zahlen:

Code: Alles auswählen

[11, 29576089, 107038193L]
('19.252', 'sec')
[413243, 842687773253L]
('1.840', 'sec')
[29576089, 107038193L]
('19.322', 'sec')

MfG
HWK

[HWK]

Das war wohl etwas voreilig. Die Implementierung funktioniert so leider nicht. Ich habe es jetzt etwas angepasst. Dazu war allerdings eine Überprüfung auf Primzahleigenschaften nötig. Dazu habe ich den Miller-Rabin-Test verwendet, auch wenn er eine gewisse Fehlerquote lässt. Mit 9.0949470177292824e-013 ist diese aber sicher zu verschmerzen. Die Geschwindigkeit wurde insgesamt durch die Änderungen erheblich gesteigert.
Hier der Code:

Code: Alles auswählen

from math import sqrt
from random import randint

TIMES = 20

def ggt(a, b):
    a, b = abs(a), abs(b)
    while b:
        a, b = b, a % b
    return a

def is_prime(n):
    if n in [1, 2, 3]:
        return True
    s = 0
    t = n - 1
    while not (t % 2):
        s += 1
        t //= 2
    for _ in xrange(TIMES):
        b = randint(2, n-2)
        z = pow(b, t, n)
        if z != 1 and z != n-1:
            for _ in xrange(TIMES):
                z = pow(z, 2, n)
                if z == n-1:
                    break
            else:
                return False
    return True

def rho(n, x0=None):
    if n == 1:
        return []
    if n == 4:
        return [2, 2]
    if is_prime(n):
        return [n]
    if x0 is None:
        x0 = randint(2, int(sqrt(n)))
    while True:
        x = y = x0
        d = randint(1, 1000)
        while True:
            x = (x * x + d) % n
            t = (y * y + d) % n
            y = (t * t + d) % n
            g = ggt(abs(y - x), n)
            if g > 1:
                if g != n:
                    return rho(n / g) + rho(g)
                break

if __name__ == '__main__':
    from time import clock
    for n in (34823482348238947, 348234823482389479, 3165771122567177):
##     for _ in xrange(5):
##         n = randint(10000000, 1000000000000000000)
        print n
        t = clock()
        print rho(n)
        print clock() - t

Und das Ergebnis:

Code: Alles auswählen

34823482348238947
[107038193L, 29576089L, 11L]
0.0919432497706
348234823482389479
[842687773253L, 413243L]
0.0252710889233
3165771122567177
[107038193L, 29576089L]
0.191376557635

MfG
HWK

[Qubit]

@Qubit: Die Zeiten sind ja erst einmal beeindruckend. Sie hängen aber doch ganz erheblich von den Ausgangszahlen ab. Bei z.B. '3165771122567177' sieht das Ergebnis deutlich schlechter aus.

Hallo HWK,
bin mir schon bewusst, dass der von mir gezeigte Algorithmus i.a. gegen probabilistische Tests, AKS, etc. schlechter aussieht. Zumindestens ist mein Algorithmus aber 100% sicher
Beschäftige mich mit Primzahlen nur nebensächlich, bin da kein Profi noch ist das ein Hobby von mir. Habe dabei aber ein paar Zusammenhänge erkannt, die mE so noch nicht publik sind (aber wie gesagt, bin auf dem Gebiet kein Profi) und im Algorithmus umgesetzt wurden. Dieser implementiert allerdings noch nicht alle diese Zusammenhänge, so dass ich ihn in absehbarer Zeit in dieser Richtung noch verbessern werde. Dann schauen wir weiter