Primfaktorzerlegung optimieren

SMaag · Beitrag von **SMaag** » 13.03.2023 14:05

brauche ne' Primfaktorzerlegung!
Bei rosettacode gibt's das, nur nicht besonders optimal!

https://rosettacode.org/wiki/Prime_deco ... #PureBasic
das ist der original Code von Rosetta, minimal angepast, so dass auf x32 und x66 direkt läuft

Code: Alles auswählen

Procedure.q Factor(Number.q, List Factors.q())
  Protected.q I, Max, cnt
  
  While Number % 2 = 0
    AddElement(Factors())
    Factors() = 2
    Number / 2
  Wend
  
  Max = Number
  I = 3
  While I <= Max And Number > 1
    While Number % I = 0
      AddElement(Factors())
      Factors() = I
      Number/I
    Wend
    I + 2
    cnt +1 
  Wend
  ProcedureReturn cnt
EndProcedure

EnableExplicit

Define Number.q = 9007199254740991

Define.q cnt, res
Define.i time.i, str.s

NewList Factors.q()

time = ElapsedMilliseconds()
  cnt = Factor(Number, Factors())
time = ElapsedMilliseconds()-time

res = 1
ForEach Factors()
  res * Factors()    
Next


str = "Number    = " + Str(Number) + #CRLF$
str + "Result P!  = " + Str(res) + #CRLF$ 
str + "Iterations = " + cnt + #CRLF$ 
str + "Time       = " + Str(time) + "ms" + #CRLF$
str + #CRLF$

ForEach Factors()
  str + #CRLF$ + Str(Factors())
Next

ClearClipboard()     ; Clear the Clipboard
SetClipboardText(str) ; Paste text to the clipboard..

MessageRequester("", str)

Meine bisherigen Optimierungen!

1. es werden nicht alle ungeraden Zahlen als potentielle Teiler untersucht,
sondern nur Primzahlkandidatenn im Modulo30 Segment
das sind die Reste 1,7,11,13, 17,19,23,29

2. Abfrage auf Number>1 hinter die While Schleife verlegt
Das ist nur für letzten Faktor relevant

3. Dynamische Anpassung der Max Grenze für die Teilerrpüfung

Das ist schon teils Faktor 1000 schneller

Code: Alles auswählen


EnableExplicit
; Kandidaten für Primzahlen sind nur Zahlen, die in Mod30 die Reste
; 1,7,11,13, 17,19,23,29 aufweisen. Alle anderen sind auf jedenfall
; teilbare Zahlen. Grund: Primfakultät  p!  das Produkt der ersten
; Primzahlen 2*3*5 = 30. D.h. in allen vielfachen Segmenten von 30
; sind dadurch alle Vielfachen von 2,3 und 5 bereits ausgestrichen.
; die nächste Primfakultät wäre dann 2*3*5*7 = 210. D.h. in allen
; Vielfachen Segementen von 210 sind dann die Vielfachen von 
; 2,3,5 und 7 bereits ausgestrichen!
; Nur die übrigen Zahlen sind noch Kandidten für Primzahlen

Global Dim Delta.i(7) 
; Delta() enthält die Abstände bis zum jeweils nächsten Primzahlkandidaten in MOD30
; damit kann man die Multiplikationen auf die Primzahlkandidaten beschränken statt auf
; alle ungeraden Zahlen
Delta(0) = 6  ; =7-1       MOD30_Rest(1) - MOD30_Rest(0)
Delta(1) = 4  ; =11-7      MOD30_Rest(2) - MOD30_Rest(1)
Delta(2) = 2  ; =13-11
Delta(3) = 4  ; =17-13
Delta(4) = 2  ; =19-17
Delta(5) = 4  ; =23-19
Delta(6) = 6  ; =29-23
Delta(7) = 2  ; =31-29     MOD30_Rest(7+1) - MOD30_Rest(7)

Procedure.q PrimeFactors(Value.q, List lstPrimeFactors.q()) 
  Protected.q P, I, Max, cnt
  
  ; ----------------------------------------------------------------------
  ; Zuerst die MOD30 = (2*3*5) Primfaktoren testen
  ; ----------------------------------------------------------------------
  While (Value & 1) = 0  ; if lo Bit=0 => devideable by 2
    Value >> 1           ; Value /2
    AddElement(lstPrimeFactors())
    lstPrimeFactors() = 2
  Wend
  
  While (Value % 3) = 0  
    Value / 3           ; Value /3
    AddElement(lstPrimeFactors())
    lstPrimeFactors() = 3
    Debug "Faktor 3"
  Wend
  
  While (Value % 5) = 0  
    Value / 5           ; Value /3
    AddElement(lstPrimeFactors())
    lstPrimeFactors() = 5
  Wend
  
  P = 7   ; nun können wir mit 7, dem nächsten Teilerkandidaten starten
  I = 0   ; Index für das OffsetArray Delta() mit den Abständen der Primkandidaten
  Max = Value/P
  
  While P <= Max
    
    While Not (Value % P)        ; (Value % P) = 0, Rest =0 => Teiler gefunden
      ; Faktor gefunden 
      AddElement(lstPrimeFactors())
      lstPrimeFactors() = P
      Value = Value / P
    Wend
       
    I = (I + 1) & 7     ; nächster OffsetIndex für Delat()-Array
    P = P + Delta(I)    ; nächster Prim Teilerkandidat aus Mod30 Segment
    
    Max = Value / P     ; die Schranke für den neuen Teiler P
    cnt +1 ; Runden zählen ! nur Statistik
  Wend
  
  If P > Max And Value > 1
    ; wenn P über Max rausgelaufen ist, so ist Value nur durch sich selbst Teilbar
    ; also Prim. Deshalb Value als letzen Primfaktor hinzufügen 
    AddElement(lstPrimeFactors())
    lstPrimeFactors() = Value      
  EndIf 

  ProcedureReturn cnt
EndProcedure

Define.q Number, time, cnt, res
Define.s str

Number.q = 27; 90071992547409; 91
NewList Factors.q()

time = ElapsedMilliseconds()
cnt = PrimeFactors(Number, Factors())
time = ElapsedMilliseconds()-time

ResetList(Factors())

res = 1
ForEach Factors()
  res * Factors()    
Next

str = "Number   = " + Str(Number) + #CRLF$
str + "Result P!  = " + Str(res) + #CRLF$ 
str + "Iterations = " + Str(cnt) + #CRLF$
str + "Time       = " + Str(time) + "ms" + #CRLF$
str + #CRLF$
str + "Prime Factors " 
ForEach Factors()
  str + #CRLF$ + Str(Factors())
Next

ClearClipboard()  ; Clear the Clipboard
SetClipboardText(str) ; Paste text to the clipboard..

MessageRequester("Prime Factors", str, #PB_MessageRequester_Ok)

Kennt da jemand noch weiter Optimierungen?

Beitrag von **NicTheQuick** » 13.03.2023 14:29

Hast du einen speziellen Anwendungsfall?
Es gibt mehrere Algorithmen dafür:

Und für ganz spezielle Fälle wie zum Beispiel das Zerlegen von Private Keys mit Längen über 1024 Bit und mehr gibt es auch Möglichkeiten, aber wie du dir sicherlich vorstellen kannst, schaffen die es nur dann in sinnvoller Zeit, wenn dieser Schlüssel aus zwei Faktoren besteht, die heutige Sicherheitsvorkehrungen nicht mehr erfüllen. Also in der Praxis bringt dir das wenig.

SMaag · Beitrag von **SMaag** » 13.03.2023 14:55

Danke! Auf die Suche nach Prim-Sieben bin ich gar nicht gekommen!
Da geht sicher noch einiges an Geschwindigkeit!

Ich hab jetzt mal mein Code weiter getestet!
mit der 8ten und 8ten Mersenn Primzahl
2^31 -1
2^61 -1

für 2^31 -1
-----------------------------------
Number = 2147483647
Result P! = 2147483647
Iterations = 12.357
Time = 0ms
Prime Factors
2147483647

für 2^61 -1
-----------------------------------
Number = 2305843009213693951
Result P! = 2305843009213693951
Iterations = 404.933.399
Time = 1198ms
Prime Factors
2305843009213693951

falls meine Version korrekt ist, reicht das für meine Zwecke im Moment

der originale Code rechnet immer noch! Hab die Hoffung, dass der noch je fertig wird fast aufgegeben!

Beitrag von **NicTheQuick** » 13.03.2023 17:04

Eventuell könntest du das ganze auch geschickt auf Threads aufteilen und alle deine Cores nutzen, dann wird es noch schneller.

SMaag · Beitrag von **SMaag** » 13.03.2023 18:59

Ich brauch nur 32 Bit ULong, das reicht mir im Moment.

Hab noch was optimiert! Schranke = Sqr(Number)
damit für 2^61 -1 nur noch ~595ms, für alle 32 Bit ist man weit unter 1ms
und soweit ich feststellen konnte stimmen die Berechnungen

-----------------------------------
Number = 2305843009213693951
Result P! = 2305843009213693951
Iterations = 404.933.400
Time = 682ms, 595ms mit HandMade Modulo
Prime Factors
2305843009213693951

Hier der finale Code

Code: Alles auswählen


EnableExplicit
; Kandidaten für Primzahlen sind nur Zahlen, die in Mod30 die Reste
; 1,7,11,13, 17,19,23,29 aufweisen. Alle anderen sind auf jedenfall
; teilbare Zahlen. Grund: Primfakultät  p!  das Produkt der ersten
; Primzahlen 2*3*5 = 30. D.h. in allen vielfachen Segmenten von 30
; sind dadurch alle Vielfachen von 2,3 und 5 bereits ausgestrichen.
; die nächste Primfakultät wäre dann 2*3*5*7 = 210. D.h. in allen
; Vielfachen Segementen von 210 sind dann die Vielfachen von 
; 2,3,5 und 7 bereits ausgestrichen!
; Nur die übrigen Zahlen sind noch Kandidten für Primzahlen

Dim PFDelta.i(7) 
; PFDelta() enthält die Abstände bis zum jeweils nächsten Primzahlkandidaten in MOD30
; damit kann man die Multiplikationen auf die Primzahlkandidaten beschränken statt auf
; alle ungeraden Zahlen

Procedure.q PrimeFactors(Number.q, List lstPrimeFactors.q())
  
  Shared PFDelta()
  Protected.q P, val, I, Max
  
  If PFDelta(0) <> 6
    ; Init Table
    PFDelta(0) = 6  ; =7-1       MOD30_Rest(1) - MOD30_Rest(0)
    PFDelta(1) = 4  ; =11-7      MOD30_Rest(2) - MOD30_Rest(1)
    PFDelta(2) = 2  ; =13-11
    PFDelta(3) = 4  ; =17-13
    PFDelta(4) = 2  ; =19-17
    PFDelta(5) = 4  ; =23-19
    PFDelta(6) = 6  ; =29-23
    PFDelta(7) = 2  ; =31-29     MOD30_Rest(7+1) - MOD30_Rest(7)   
  EndIf
  
  ClearList(lstPrimeFactors())

  ; die max Anzahl der nötigen Iterationen = Sqr(val)/30*8
  ; also pro 30iger Segment sind es 8 Iterationen
  ; für die 2^31-1 => 12.352; das geht in unter 1ms
  
  val = Number
  ; ----------------------------------------------------------------------
  ; Zuerst die MOD30 = (2*3*5) Primfaktoren testen
  ; ----------------------------------------------------------------------
  While (val & 1) = 0  ; if lo Bit=0 => devideable by 2
    val >> 1           ; val /2
    AddElement(lstPrimeFactors())
    lstPrimeFactors() = 2
  Wend
  
  While Not (val % 3)  
    val / 3           ; val /3
    AddElement(lstPrimeFactors())
    lstPrimeFactors() = 3
    Debug "Faktor 3"
  Wend
  
  While Not (val % 5)
    val / 5           ; val /3
    AddElement(lstPrimeFactors())
    lstPrimeFactors() = 5
  Wend
  
  P = 7   ; nun können wir mit 7, dem nächsten Teilerkandidaten starten
  I = 0   ; Index für das OffsetArray PFDelta() mit den Abständen der Primkandidaten
  
  ; ACHTUNG Sqr arbeitet mit double Float. D.h. 52Bit Mantisse. Das fürt bei bei
  ; Das geht aber, da die FPU intern mit 80 oder mehr Bits arbeitet
  Max = Sqr(val)+1 ; val/P
  
  While P <= Max  
    ; While Not (val % P)        ; (val % P) = 0, Rest =0 => Teiler gefunden
     While ((val/P)*P-val) = 0   ; HandMade Modulo ist schneller

      ; Faktor gefunden 
      AddElement(lstPrimeFactors())
      lstPrimeFactors() = P
      val = val / P
      Max = Sqr(val)+1  ; die Schranke für den neuen Teiler P
    Wend      
    I = I + 1               ; Iterations; daraus lässt sich der Index mit (I & 7) bestimmen
    P = P + PFDelta(I & 7)    ; nächster Prim Teilerkandidat aus Mod30 Segment   
  Wend
  
  If P > Max And val > 1
    ; wenn P über Max rausgelaufen ist, so ist val nur durch sich selbst Teilbar
    ; also Prim. Deshalb val als letzen Primfaktor hinzufügen 
    AddElement(lstPrimeFactors())
    lstPrimeFactors() = val      
  EndIf 

  ProcedureReturn I;  Bool(Number=val) ; um als Rückgabe IsPrim(Number) zu erhalten
EndProcedure

Define.q Number, time, cnt, res
Define.s str

Number.q = 9007199254740991
;Number = (1<<31)-1  ; 8te Mersenn Primzahl 2^31-1 =            2.147.483.647
Number = (1<<61)-1   ; 9te Mersenn Primzahl 2^61-1 = 2.305.843.009.213.693.951

NewList Factors.q()

time = ElapsedMilliseconds()
cnt = PrimeFactors(Number, Factors())
time = ElapsedMilliseconds()-time

ResetList(Factors())

res = 1
ForEach Factors()
  res * Factors()    
Next

str = "Number   = " + Str(Number) + #CRLF$
str + "Result P!  = " + Str(res) + #CRLF$ 
str + "Iterations = " + Str(cnt) + #CRLF$
str + "Time       = " + Str(time) + "ms" + #CRLF$
str + #CRLF$
str + "Prime Factors " 
ForEach Factors()
  str + #CRLF$ + Str(Factors())
Next

ClearClipboard()  ; Clear the Clipboard
SetClipboardText(str) ; Paste text to the clipboard..

MessageRequester("Prime Factors", str, #PB_MessageRequester_Ok)

PureBoard

Primfaktorzerlegung optimieren

Primfaktorzerlegung optimieren

Re: Primfaktorzerlegung optimieren

Re: Primfaktorzerlegung optimieren

Re: Primfaktorzerlegung optimieren

Re: Primfaktorzerlegung optimieren