FFT-Calculator, AVX2-und AVX-Version

Hier könnt Ihr gute, von Euch geschriebene Codes posten. Sie müssen auf jeden Fall funktionieren und sollten möglichst effizient, elegant und beispielhaft oder einfach nur cool sein.
Benutzeravatar
Helle
Beiträge: 566
Registriert: 11.11.2004 16:13
Wohnort: Magdeburg

FFT-Calculator, AVX2-und AVX-Version

Beitrag von Helle »

Zur Feier der "Findung" der 50(?).Mersenne-Zahl (s.Heise heute) habe ich meinen FFT-Calculator aktualisiert und bin mutig (besoffen?) genug, den Code hier zu zeigen. Der war eigentlich nur für heimische Tests gedacht (deshalb auch z.B. kaum Prüfung auf Eingabe-Fehler), also bitte Nachsicht :lol: . Benötigt wird allerdings eine CPU mit AVX2-Unterstützung.
Meine 5GHz-CPU spuckt die 23.249.425 Ziffern von Mersenne M50 übrigens nach 44s aus.
Hier der (noch verbesserungsbedürftige) Code:

Code: Alles auswählen

;- Schnelle Unsigned Integer-Multiplication mit FFT
;- see https://en.wikipedia.org/wiki/Butterfly_diagram (radix-2 decimation-in-time FFT algorithm)
;- Windows 64-Bit, getestet mit PureBasic 5.46 LTS (x64) / PureBasic 5.61 (x64)
;- Ascii und Unicode
;- Benötigt CPU mit AVX2-Unterstützung
;- "Helle" Klaus Helbing, 05.01.2018: Mersenne M50 hinzugefügt
;- Um in den Genuss der Info-Anzeige zu kommen, sollten Programme wie der Taskmanager nicht gestartet sein!
;- Version vom 11.01.2018

Global.q Asc_Uni, Bits_M, CalcPot, FakA, GadgetZeile, Potenz, ZeitA, ZeitE
Global.s Basis$, Result$

Declare.s FFT_Mul(Para1, Para2) 
Declare Calc()
Declare RangeCharSSE42(Para1, Para2)
Declare FakInfo()

;Test auf AVX2
!xor ecx,ecx
!mov eax,7         ;Test auf 7 geschenkt 
!cpuid
!test ebx,20h      ;Bit5 für AVX2
!jnz @f            ;Yeah!
MessageRequester("Schade!", "Diese CPU unterstützt kein AVX2!" + #LFCR$ + "Ende")
End 
!@@:

FensterTitel$ = "Helles FFT-Calculator AVX2-Version"
OpenWindow(0, 0, 0, 800, 600, FensterTitel$, #PB_Window_SystemMenu | #PB_Window_ScreenCentered)

CompilerIf #PB_Compiler_Unicode
  Asc_Uni = 2 : #Asc_Uni = 2 ;mit Konstante für Step usw.
 CompilerElse
  Asc_Uni = 1 : #Asc_Uni = 1
CompilerEndIf

FontHigh = Int(9.0 / (GetDeviceCaps_(GetDC_(WindowID(0)), #LOGPIXELSY) / 96.0))     ;GetDeviceCaps kann bei der Bildschirmanpassung nützlich sein
;LoadFont(0, "Arial", FontHigh)
LoadFont(0, "Trebuchet MS Fett", FontHigh)       ;in Windows 7 verfügbar; evtl.testen

Repeat
  TextGadget(1, 20, 20, 300, 15, "Was soll berechnet werden:")
  OptionGadget(2, 20, 45, 780, 15, "Multiplikation zweier positiver Integer-Faktoren")
  OptionGadget(3, 20, 65, 780, 15, "Berechnung einer positiven ganzzahligen Potenz")
  OptionGadget(4, 20, 85, 780, 15, "Berechnung einer Mersenne-Prim-Zahl (M25...M50)")
  OptionGadget(5, 20, 105, 780, 15, "Fakultät")
  OptionGadget(6, 20, 125, 780, 15, "Neuner-Test")
  ButtonGadget(7, 300, 560, 200, 30, "W e i t e r")
  SetGadgetState(2, 1) : FFT = 2       ;Voreinstellung Multiplikation, falls nichts ausgewählt wird

  For i = 1 To 7
    SetGadgetFont(i, FontID(0))
  Next

  Repeat
    Event = WaitWindowEvent()
    If Event = #PB_Event_CloseWindow
      CloseWindow(0)
      End
    EndIf
    If Event = #PB_Event_Gadget
      Select EventGadget()
        Case 7
          Break          
        Default
          FFT = EventGadget()
      EndSelect
    EndIf
  ForEver

  For i = 1 To 7
    FreeGadget(i)
  Next

  Select FFT
;=========================================================================
    Case 2                   ;Multiplikation 
      Info0$ = "Multiplikation "
      SetWindowTitle(0, Info0$)

      CharRange$ = Chr(1) + Chr(47) + Chr(58) + Chr(255)   ;für Suche nach Chars, die keine Dezimal-Ziffer sind. Chr(0) geht nicht!

      TextGadget(1, 50, 0, 300, 20, "Faktor1 eintragen (z.B. einen String mit Strg+V)")
      TextGadget(2, 450, 0, 300, 20, "Faktor2 eintragen (z.B. einen String mit Strg+V)")
      ButtonGadget(3, 325, 560, 150, 30, "Berechnung starten")
      ButtonGadget(4, 125, 560, 150, 30, "Faktor1 überprüfen")
      ButtonGadget(5, 525, 560, 150, 30, "Faktor2 überprüfen")

      For i = 1 To 5
        SetGadgetFont(i, FontID(0))
      Next

      EditorGadget(6, 10, 20, 385, 530)
      EditorGadget(7, 405, 20, 385, 530)
      SetActiveGadget(6)

      Repeat
        Event = WaitWindowEvent()
        If Event = #PB_Event_CloseWindow
          CloseWindow(0)
          End 
        EndIf
        If Event = #PB_Event_Gadget
          Gadget = EventGadget()
          Select Gadget
            Case 3
              Break
            Case 4, 5
              If Gadget = 4
                Faktor$ = GetGadgetText(6)
               Else
                Faktor$ = GetGadgetText(7)
              EndIf
              RangeChar = RangeCharSSE42(@Faktor$, @CharRange$)
              If RangeChar
                Ascii = Asc(Mid(Faktor$, RangeChar, 1)) : Char$ = Mid(Faktor$, RangeChar, 1) 
                MessageRequester("Faktor" + Str(Gadget - 3) + ": Fehler!", "Digit an Stelle " + Str(RangeChar) + " ist keine Dezimal-Ziffer! (ASCII=" + Ascii + "), Char=" + Char$, #MB_ICONERROR)
               Else
                MessageRequester("Faktor" + Str(Gadget - 3) + ": O.K.!", "Alle Digits sind Dezimal-Ziffern!")
              EndIf
          EndSelect
        EndIf
      ForEver

      Faktor1$ = GetGadgetText(6)
      Faktor2$ = GetGadgetText(7)

      For i = 1 To 7
        FreeGadget(i)
      Next

      ZeitA = ElapsedMilliseconds()
      Result$ = FFT_Mul(@Faktor1$, @Faktor2$)
      ZeitE = ElapsedMilliseconds() - ZeitA
;=========================================================================
    Case 3                   ;Potenzen
      Info0$ = "Berechnung von X^Y "
      SetWindowTitle(0, Info0$)
      TextGadget(1, 50, 10, 300, 20, "Basis X eintragen für X^Y")
      TextGadget(2, 400, 10, 300, 20, "Exponent Y eintragen für X^Y")
      StringGadget(3, 50, 50, 300, 20, "2", #PB_String_Numeric) 
      StringGadget(4, 400, 50, 300, 20, "1", #PB_String_Numeric)
      ButtonGadget(5, 300, 560, 200, 30, "Berechnung starten")
      SetActiveGadget(2)

      For i = 1 To 5         ;das StringGadget einfach mal drin lassen
        SetGadgetFont(i, FontID(0))
      Next

      Repeat 
        Repeat
          Event = WaitWindowEvent()
          If Event = #PB_Event_CloseWindow
            CloseWindow(0)
            End 
          EndIf
        Until EventGadget() = 5

        Potenz = Val(GetGadgetText(4)) ;Exponent
        If (Len(GetGadgetText(4)) > 17) Or (Potenz > 1E16) ;willkürlich gewählte Werte, sichere Seite
          MessageRequester("Abbruch!", "Exponent zu groß! Bitte neue Eingabe!")
         Else
          Break
        EndIf     
      ForEver

      Basis$ = GetGadgetText(3)

      Basis = 0
      If Len(Basis$) < 2
        Basis = Len(Basis$)
      EndIf

      For i = 1 To 5
        FreeGadget(i)
      Next

      Info0$ = Basis$ + "^" + Str(Potenz) + ": "
      SetWindowTitle(0, Info0$)

      If Len(Basis$) > 0 Or Basis > 0
        If Potenz > 0
          Result$ = Basis$
          Calc()
         Else
          Result$ = "1"
        EndIf
       Else
        Result$ = "0"        ;0^Y
      EndIf
;=========================================================================
    Case 4                   ;Mersenne
      Info0$ = "Berechnung von Mersenne "
      SetWindowTitle(0, Info0$ + "Mxx")
      ;Werte für Mersenne M25-M49, die Strings sind nur zur Info
      M25.q = 21701    : M25$ = "(2^21701)-1"    : Z25$ = "44867916611904333479...57410828353511882751"
      M26.q = 23209    : M26$ = "(2^23209)-1"    : Z26$ = "40287411577898877818...36743355523779264511"
      M27.q = 44497    : M27$ = "(2^44497)-1"    : Z27$ = "85450982430363380319...44867686961011228671"
      M28.q = 86243    : M28$ = "(2^86243)-1"    : Z28$ = "53692799550275632152...99857021709433438207"
      M29.q = 110503   : M29$ = "(2^110503)-1"   : Z29$ = "52192831334175505976...69951621083465515007"
      M30.q = 132049   : M30$ = "(2^132049)-1"   : Z30$ = "51274027626932072381...52138578455730061311"
      M31.q = 216091   : M31$ = "(2^216091)-1"   : Z31$ = "74609310306466134368...91336204103815528447"
      M32.q = 756839   : M32$ = "(2^756839)-1"   : Z32$ = "17413590682008709732...02603793328544677887"
      M33.q = 859433   : M33$ = "(2^859433)-1"   : Z33$ = "12949812560420764966...02414267243500142591"
      M34.q = 1257787  : M34$ = "(2^1257787)-1"  : Z34$ = "41224577362142867472...31257188976089366527"
      M35.q = 1398269  : M35$ = "(2^1398269)-1"  : Z35$ = "81471756441257307514...85532025868451315711"
      M36.q = 2976221  : M36$ = "(2^2976221)-1"  : Z36$ = "62334007624857864988...76506256743729201151"
      M37.q = 3021377  : M37$ = "(2^3021377)-1"  : Z37$ = "12741168303009336743...25422631973024694271"
      M38.q = 6972593  : M38$ = "(2^6972593)-1"  : Z38$ = "43707574412708137883...35366526142924193791"
      M39.q = 13466917 : M39$ = "(2^13466917)-1" : Z39$ = "92494773800670132224...30073855470256259071"
      M40.q = 20996011 : M40$ = "(2^20996011)-1" : Z40$ = "12597689545033010502...94714065762855682047"
      M41.q = 24036583 : M41$ = "(2^24036583)-1" : Z41$ = "29941042940415717208...67436921882733969407"
      M42.q = 25964951 : M42$ = "(2^25964951)-1" : Z42$ = "12216463006127794810...98933257280577077247"
      M43.q = 30402457 : M43$ = "(2^30402457)-1" : Z43$ = "31541647561884608093...11134297411652943871"
      M44.q = 32582657 : M44$ = "(2^32582657)-1" : Z44$ = "12457502601536945540...11752880154053967871"
      M45.q = 37156667 : M45$ = "(2^37156667)-1" : Z45$ = "20225440689097733553...21340265022308220927"
      M46.q = 42643801 : M46$ = "(2^42643801)-1" : Z46$ = "16987351645274162247...84101954765562314751"
      M47.q = 43112609 : M47$ = "(2^43112609)-1" : Z47$ = "31647026933025592314...80022181166697152511"
      M48.q = 57885161 : M48$ = "(2^57885161)-1" : Z48$ = "58188726623224644217...46141988071724285951"
      M49.q = 74207281 : M49$ = "(2^74207281)-1" : Z49$ = "30037641808460618205...87010073391086436351"
      M50.q = 77232917 : M50$ = "(2^77232917)-1" : Z50$ = "46733318335923109998...82730618069762179071"

      TextGadget(1, 20, 0, 300, 15, "Wähle eine Mersenne-Zahl:")

      OptionGadget(2, 20, 20, 780, 15, "Mersenne M25, " + M25$ + ", 6533 Dezimal-Ziffern, " + Z25$)
      OptionGadget(3, 20, 40, 780, 15, "Mersenne M26, " + M26$ + ", 6987 Dezimal-Ziffern, " + Z26$)
      OptionGadget(4, 20, 60, 780, 15, "Mersenne M27, " + M27$ + ", 13395 Dezimal-Ziffern, " + Z27$)
      OptionGadget(5, 20, 80, 780, 15, "Mersenne M28, " + M28$ + ", 25962 Dezimal-Ziffern, " + Z28$)
      OptionGadget(6, 20, 100, 780, 15, "Mersenne M29, " + M29$ + ", 33265 Dezimal-Ziffern, " + Z29$)
      OptionGadget(7, 20, 120, 780, 15, "Mersenne M30, " + M30$ + ", 39751 Dezimal-Ziffern, " + Z30$)
      OptionGadget(8, 20, 140, 780, 15, "Mersenne M31, " + M31$ + ", 65050 Dezimal-Ziffern, " + Z31$)
      OptionGadget(9, 20, 160, 780, 15, "Mersenne M32, " + M32$ + ", 227832 Dezimal-Ziffern, " + Z32$)
      OptionGadget(10, 20, 180, 780, 15, "Mersenne M33, " + M33$ + ", 258716 Dezimal-Ziffern, " + Z33$)
      OptionGadget(11, 20, 200, 780, 15, "Mersenne M34, " + M34$ + ", 378632 Dezimal-Ziffern, " + Z34$)
      OptionGadget(12, 20, 220, 780, 15, "Mersenne M35, " + M35$ + ", 420921 Dezimal-Ziffern, " + Z35$)
      OptionGadget(13, 20, 240, 780, 15, "Mersenne M36, " + M36$ + ", 895832 Dezimal-Ziffern, " + Z36$)
      OptionGadget(14, 20, 260, 780, 15, "Mersenne M37, " + M37$ + ", 909526 Dezimal-Ziffern, " + Z37$)
      OptionGadget(15, 20, 280, 780, 15, "Mersenne M38, " + M38$ + ", 2098960 Dezimal-Ziffern, " + Z38$)
      OptionGadget(16, 20, 300, 780, 15, "Mersenne M39, " + M39$ + ", 4053946 Dezimal-Ziffern, " + Z39$)
      OptionGadget(17, 20, 320, 780, 15, "Mersenne M40, " + M40$ + ", 6320430 Dezimal-Ziffern, " + Z40$)
      OptionGadget(18, 20, 340, 780, 15, "Mersenne M41, " + M41$ + ", 7235733 Dezimal-Ziffern, " + Z41$)
      OptionGadget(19, 20, 360, 780, 15, "Mersenne M42, " + M42$ + ", 7816230 Dezimal-Ziffern, " + Z42$)
      OptionGadget(20, 20, 380, 780, 15, "Mersenne M43, " + M43$ + ", 9152052 Dezimal-Ziffern, " + Z43$)
      OptionGadget(21, 20, 400, 780, 15, "Mersenne M44, " + M44$ + ", 9808358 Dezimal-Ziffern, " + Z44$)
      OptionGadget(22, 20, 420, 780, 15, "Mersenne M45, " + M45$ + ", 11185272 Dezimal-Ziffern, " + Z45$)
      OptionGadget(23, 20, 440, 780, 15, "Mersenne M46, " + M46$ + ", 12837064 Dezimal-Ziffern, " + Z46$)
      OptionGadget(24, 20, 460, 780, 15, "Mersenne M47, " + M47$ + ", 12978189 Dezimal-Ziffern, " + Z47$)
      OptionGadget(25, 20, 480, 780, 15, "Mersenne M48, " + M48$ + ", 17425170 Dezimal-Ziffern, " + Z48$)
      OptionGadget(26, 20, 500, 780, 15, "Mersenne M49, " + M49$ + ", 22338618 Dezimal-Ziffern, " + Z49$)
      OptionGadget(27, 20, 520, 780, 15, "Mersenne M50, " + M50$ + ", 23249425 Dezimal-Ziffern, " + Z50$)

      SetGadgetState(2, 1)

      ButtonGadget(28, 300, 560, 200, 30, "Berechnung starten")

      For i = 1 To 28
        SetGadgetFont(i, FontID(0))
      Next

      Basis$ = "2"
      Result$ = "2"
      Potenz = M25 : M$ = "M25=" + M25$ : Z$ = Z25$   ;Voreinstellung Mersenne, falls nichts ausgewählt wird

      Repeat
        Event = WaitWindowEvent()
        If Event = #PB_Event_CloseWindow
          CloseWindow(0)
          End
        EndIf
        If Event = #PB_Event_Gadget
          Select EventGadget()
            Case 2
              Potenz = M25 : M$ = "M25=" + M25$ : Z$ = Z25$
            Case 3
              Potenz = M26 : M$ = "M26=" + M26$ : Z$ = Z26$
            Case 4
              Potenz = M27 : M$ = "M27=" + M27$ : Z$ = Z27$
            Case 5
              Potenz = M28 : M$ = "M28=" + M28$ : Z$ = Z28$
            Case 6
              Potenz = M29 : M$ = "M29=" + M29$ : Z$ = Z29$
            Case 7
              Potenz = M30 : M$ = "M30=" + M30$ : Z$ = Z30$
            Case 8
              Potenz = M31 : M$ = "M31=" + M31$ : Z$ = Z31$
            Case 9
              Potenz = M32 : M$ = "M32=" + M32$ : Z$ = Z32$
            Case 10
              Potenz = M33 : M$ = "M33=" + M33$ : Z$ = Z33$
            Case 11
              Potenz = M34 : M$ = "M34=" + M34$ : Z$ = Z34$
            Case 12
              Potenz = M35 : M$ = "M35=" + M35$ : Z$ = Z35$
            Case 13
              Potenz = M36 : M$ = "M36=" + M36$ : Z$ = Z36$
            Case 14
              Potenz = M37 : M$ = "M37=" + M37$ : Z$ = Z37$
            Case 15
              Potenz = M38 : M$ = "M38=" + M38$ : Z$ = Z38$
            Case 16
              Potenz = M39 : M$ = "M39=" + M39$ : Z$ = Z39$
            Case 17
              Potenz = M40 : M$ = "M40=" + M40$ : Z$ = Z40$
            Case 18
              Potenz = M41 : M$ = "M41=" + M41$ : Z$ = Z41$
            Case 19
              Potenz = M42 : M$ = "M42=" + M42$ : Z$ = Z42$
            Case 20
              Potenz = M43 : M$ = "M43=" + M43$ : Z$ = Z43$
            Case 21
              Potenz = M44 : M$ = "M44=" + M44$ : Z$ = Z44$
            Case 22
              Potenz = M45 : M$ = "M45=" + M45$ : Z$ = Z45$
            Case 23
              Potenz = M46 : M$ = "M46=" + M46$ : Z$ = Z46$
            Case 24
              Potenz = M47 : M$ = "M47=" + M47$ : Z$ = Z47$
            Case 25
              Potenz = M48 : M$ = "M48=" + M48$ : Z$ = Z48$
            Case 26
              Potenz = M49 : M$ = "M49=" + M49$ : Z$ = Z49$
            Case 27
              Potenz = M50 : M$ = "M50=" + M50$ : Z$ = Z50$
            Case 28
              Info0$ + M$ + ": "
              SetWindowTitle(0, Info0$)
              Break          
          EndSelect
        EndIf
      ForEver

      For i = 1 To 28
        FreeGadget(i)
      Next

      Calc()

      ResLen.q = StringByteLength(Result$)
      PokeB(@Result$ + ResLen - Asc_Uni, PeekB(@Result$ + ResLen - Asc_Uni) - 1)    ;-1 für Mersenne
;=========================================================================
    Case 5                   ;Fakultät
      Info0$ = "Fakultät: "
      SetWindowTitle(0, Info0$)

      TextGadget(1, 50, 10, 300, 20, "Fakultät:")
      StringGadget(2, 50, 50, 300, 20, "2", #PB_String_Numeric)  
      ButtonGadget(3, 300, 560, 200, 30, "Berechnung starten")
      SetActiveGadget(2)

      For i = 1 To 3         ;das StringGadget einfach mal drin lassen
        SetGadgetFont(i, FontID(0))
      Next

      Repeat
        Repeat
          Event = WaitWindowEvent()
          If Event = #PB_Event_CloseWindow
            CloseWindow(0)
            End 
          EndIf
        Until EventGadget() = 3

        Fak = Val(GetGadgetText(2))    ;hier reicht der Wertebereich von Val() wohl aus ;-)   
        ;besser doch Test
        If (Len(GetGadgetText(2)) > 7) Or (Fak > 1E6) ;willkürlich gewählte Werte, sichere Seite
          MessageRequester("Abbruch!", "Wert zu groß! Bitte neue Eingabe!")
         Else
          Break
        EndIf     
      ForEver

      Info0$ + Str(Fak) + "! "
      SetWindowTitle(0, Info0$)      

      Result$ = "1"

      For i = 1 To 3
        FreeGadget(i)
      Next

      ListViewGadget(2, 10, 10, 250, 580)   ;Info-Fenster

      AddWindowTimer(0, 0, 5000)            ;für Anzeige im Info-Fenster
      BindEvent(#PB_Event_Timer, @FakInfo())

      GadgetZeile = 0
      AddGadgetItem (2, -1, "Berechnung gestartet... ")   
      SetGadgetState(2, GadgetZeile)

      ZeitA = ElapsedMilliseconds()
      For FakA = 1 To Fak
        Faktor1$ = Str(FakA)
        WindowEvent()                       ;nur für Anzeige im Info-Fenster
        Result$ = FFT_Mul(@Faktor1$, @Result$)
      Next
      ZeitE = ElapsedMilliseconds() - ZeitA

      RemoveWindowTimer(0, 0)
      UnbindEvent(#PB_Event_Timer, @FakInfo())   ;nötig!
      FreeGadget(2)
;=========================================================================
    Case 6                   ;Neuner-Test
      Info0$ = "Neuner-Test: "
      SetWindowTitle(0, Info0$)

      Test9 = 540000         ;540000 ist o.K. für Double

      Buffer9 = AllocateMemory(32 * Test9 * #Asc_Uni + 64)
      Buffer9A = Buffer9 + (64 - (Buffer9 & $3F))
      !mov rax,[v_Buffer9A]
      !mov byte[rax],9
      !mov rcx,[v_Test9]
      !cmp [v_Asc_Uni],1
     !je .ASCII
      !vpbroadcastw ymm0,[rax]    ;Unicode Word
      !shl rcx,1
     !jmp @f  
     !.ASCII:
      !vpbroadcastb ymm0,[rax]    ;ASCII Byte
     !@@:
      !vmovdqa [rax],ymm0
      !add rax,32
      !dec rcx
     !jnz @b

      Info0$ + "Berechne " + Str(32 * Test9) + "*" + Str(32 * Test9) + " Neuner-Bytes (Digits): "
      SetWindowTitle(0, Info0$)
      ZeitA = ElapsedMilliseconds()
      Result$ = FFT_Mul(Buffer9A, Buffer9A)
      ZeitE = ElapsedMilliseconds() - ZeitA

      Info9$ = ""
      L = StringByteLength(Result$)
      Buffer1 = @Result$

      If PeekB(Buffer1 + L - Asc_Uni) <> $31     ;hier könnte auch die Konstante genommen werden
        Info9$ = "Fehler 1! "
      EndIf

      If PeekB(Buffer1 + (L / 2) - Asc_Uni) <> $38
        Info9$ + "Fehler 8! "
      EndIf

      For i = 0 To (L / 2) - 2 * Asc_Uni Step #Asc_Uni
        If PeekB(Buffer1 + i) <> $39
          Info9$ + "Fehler 9! "
          Break
        EndIf
      Next

      For i = L / 2 To L - 2 * Asc_Uni Step #Asc_Uni
        If PeekB(Buffer1 + i) <> $30
          Info9$ + "Fehler 0! "
          Break
        EndIf
      Next

      If Info9$ = ""
        Info9$ = "Alles o.K.! "
      EndIf

      FreeMemory(Buffer9)
      Info0$ + Info9$
      SetWindowTitle(0, Info0$)
;=========================================================================
  EndSelect
;=========================================================================
  EditorGadget(1, 10, 10, 780, 540)

  Info1$ = Info0$ + "Benötigte Zeit: " + Str(ZeitE) + "ms  " + Z$
  SetWindowTitle(0, Info1$)

  Digits$ = ""
  If Asc(Mid(Result$, 1, 1)) < 58 ;also Ziffer (hoffentlich :-) !) und nicht z.B.Fehlermeldung
    Digits$ = "Digits:" + Str(Len(Result$)) + #LFCR$
  EndIf
  SetGadgetText(1, Digits$ + Result$)

  ButtonGadget(2, 100, 560, 250, 30, "In die Zwischenablage kopieren")
  SetGadgetFont(2, FontID(0))
  ButtonGadget(3, 450, 560, 250, 30, "Zurück zur Auswahl")
  SetGadgetFont(3, FontID(0))

  Repeat
    Event = WaitWindowEvent()
    If Event = #PB_Event_CloseWindow
      Break 2
    EndIf
    If Event = #PB_Event_Gadget
      Select EventGadget()
        Case 2
          SetClipboardText(Result$)
        Case 3
          SetWindowTitle(0, FensterTitel$)
          Info0$ = "" : Info1$ = "" : Z$ = ""
          Break          
      EndSelect
    EndIf
  ForEver

  For i = 1 To 3
    FreeGadget(i)
  Next

Until WaitWindowEvent() = #PB_Event_CloseWindow

CloseWindow(0)
End

Procedure.s FFT_Mul(PFac1, PFac2)
  EnableASM 

  Protected.q BufferCosSin, BufferCosSinA, BufferDez, BufferDezA, BufferFX, BufferFXA, BufferFaktor1, BufferFaktor1A
  Protected.q BufferFaktor2, BufferFaktor2A, BufferProdukt, BufferProduktA, BufferReverse, BufferReverseA, LG, N
  Protected.s sAusgabe

  ;Register-, Status- und Control- Speicher-Reservierung
  BufferFX = AllocateMemory(512 + 64)
  If BufferFX
    BufferFXA = BufferFX + (64 - (BufferFX & $3F))    ;64-er Alignment, klotzen, nicht kleckern!
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferFXA
  !pop [BufferFXA]
;!Align 16
  ;Save Registers
  !mov rax,[BufferFXA]
  !fxsave64 [rax]            ;sichert erstmal komplett FPU (mit MXCSR) und XMM (0-15). Lass ich mal so, obwohl FNSAVE auch reichen würde. Dann aber so keine Register-Sicherung
  !mov [rax+464],rdi         ;sichert callee-save registers, ohne RBX und RBP!
  !mov [rax+472],rsi
  !mov [rax+480],r12
  !mov [rax+488],r13
  !mov [rax+496],r14
  !mov [rax+504],r15

  ;übergebende Parameter
  !mov r8,[p.v_PFac1]
  !mov [PF1],r8
  !mov r9,[p.v_PFac2]
  !mov [PF2],r9

  !mov rax, -16
  !mov rdx, -16
  !pxor xmm1,xmm1            ;xmm1 null setzen, ist praktisch das Zero-Byte bzw. 2 Zero-Bytes für Unicode

  ;Länge der Strings ermitteln
  CompilerIf #PB_Compiler_Unicode
   !@@:
    !add rax,16
    !vpcmpistri xmm1,dqword[r8+rax],00001001b    ;Bit0 = 1 und Bit1 = 0, also String-Chars werden als unsigned Words interpretiert, Bit2 = 0 und Bit3 = 1, also Test auf equal each
   !jnz @b
    !shr rax,1                                   ;Pointer in String muss hier halbiert werden
    !add rax,rcx

    !pxor xmm1,xmm1                              ;xmm1 null setzen, ist praktisch das Zero-Byte bzw. 2 Zero-Bytes für Unicode
   !@@:
    !add rdx,16
    !vpcmpistri xmm1,dqword[r9+rdx],00001001b    ;Bit0 = 1 und Bit1 = 0, also String-Chars werden als unsigned Words interpretiert, Bit2 = 0 und Bit3 = 1, also Test auf equal each
   !jnz @b
    !shr rdx,1                                   ;Pointer in String muss hier halbiert werden
    !add rdx,rcx
   CompilerElse
    !@@:
    !add rax,16
    !vpcmpistri xmm1,dqword[r8+rax],00001000b    ;Bit0 = 0 und Bit1 = 0, also String-Chars werden als unsigned Bytes interpretiert, Bit2 = 0 und Bit3 = 1, also Test auf equal each
   !jnz @b
    !add rax,rcx

    !pxor xmm1,xmm1                              ;xmm1 null setzen, ist praktisch das Zero-Byte bzw. 2 Zero-Bytes für Unicode
   !@@:
    !add rdx,16
    !vpcmpistri xmm1,dqword[r9+rdx],00001000b    ;Bit0 = 0 und Bit1 = 0, also String-Chars werden als unsigned Bytes interpretiert, Bit2 = 0 und Bit3 = 1, also Test auf equal each
   !jnz @b
    !add rdx,rcx
  CompilerEndIf

  !mov [L1],rax
  !mov [L2],rdx
  !cmp rax,rdx               ;L1 soll der größere Faktor sein
 !jae @f
  !mov [L1],rdx
  !mov [L2],rax
  !mov [PF1],r9
  !mov [PF2],r8
 !@@:

  !add rax,rdx
  !push rax
  POP LG                     ;Gesamtlänge L1+L2 für "Verstringung" und Speicher-Management

  ;LenF = L1
  !mov rdx,[L1]              ;L1 ist größter Faktor
  !mov [LenF],rdx            ;LenF erstmal setzen, könnte ja schon 2-er Potenz sein 
  !bsr rcx,[L1]              ;höchstes gesetztes Bit
  !bsf rdx,[L2]              ;Test auf Länge Null
 !jz @f
  !mov [Exponent],rcx
  !bsf rdx,[L1]              ;niedrigstes gesetztes Bit
 !jmp .Faktoren_OK           ;local label
 !@@:
  sAusgabe = "Fehler! Mindestens einer der Faktoren hat Länge Null! " ;wegen Kompatibilität zu "normaler" Multiplikation gelassen
 ProcedureReturn sAusgabe

 !.Faktoren_OK:
  !cmp rcx,rdx
 !je @f                      ;also nur 1 Bit (identisch) gesetzt, ist 2-er Potenz
  !inc rcx
  !mov rax,1
  !shl rax,cl
  !mov [LenF],rax            ;LenF ist jetzt 2-er Potenz
  !add [Exponent],1
 !@@:
  !add [Exponent],1

  !mov rcx,5                 ;Produkt-Länge muss mindestens 32 sein, sonst an anderen Stellen Klimmzüge ("Verstringung")
  !sub rcx,[Exponent]
 !js @f                      ;Exponent ist > 5 (Hauptfall)
 !jz @f                      ;Exponent ist = 5
  !shl [LenF],cl
  !add [Exponent],rcx
 !@@:

  ;zur Sicherheit im MXCSR Bit 13 und 14 auf Null setzen (round to nearest)
  ;MXCSR[6]         : Denormals are zeros - 0
  ;MXCSR[7:12]      : Exception masks all 1's (all exceptions masked)
  ;MXCSR[13:14]   : Rounding  control - 0 (round to nearest)
  ;MXCSR[15]      : Flush to zero for masked underflow - 0 (off)
  ;R+	bit 14	Round Positive
  ;R-	bit 13	Round Negative
  ;RZ	bits 13 and 14	Round to Zero
  ;RN	bits 13 and 14 are 0	Round to Nearest
  ;load MXCSR
  !vstmxcsr [MXCSR]          ;load Control-and Status-Register in Memory
  !mov eax,[MXCSR] 
  !mov [MXCSR_Sicher],eax    ;trotz fxsave64 (s.o.) hier zur Sicherheit 
  ;MXCSR_Sicher = MXCSR
  !and [MXCSR],11111111111111111001111111111111b ;zur Sicherheit Bit 13 und 14 auf Null setzen (round to nearest)
  !vldmxcsr [MXCSR]          ;store MXCSR

  ;N = LenF << 1
  !mov rax,[LenF]
  !shl rax,1
  !mov [N],rax
  !push rax                  ;wird dann "PB-N"

  !vcvtsi2sd xmm1,xmm2,rax   ;Quad-N in Double-ND konvertieren für spätere Division
  !vmovsd [ND],xmm1

  ;Pointer_Produkt.q = N - 1
  !dec rax
  !mov [Pointer_Produkt],rax

  POP N                      ;wird u.a.für AllocateMemory benötigt!

  ;Speicher-Reservierung
  BufferCosSin = AllocateMemory((2 * N * 8) + 64)
  If BufferCosSin
    BufferCosSinA = BufferCosSin + (64 - (BufferCosSin & $3F))  ;64-er Alignment, klotzen, nicht kleckern!
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferCosSinA
  !pop [BufferCosSinA]

  BufferFaktor1 = AllocateMemory((2 * N * 8) + 64)
  If BufferFaktor1
    BufferFaktor1A = BufferFaktor1 + (64 - (BufferFaktor1 & $3F))    ;64-er Alignment
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferFaktor1A
  !pop [BufferFaktor1A]

  BufferFaktor2 = AllocateMemory((2 * N * 8) + 64)
  If BufferFaktor2
    BufferFaktor2A = BufferFaktor2 + (64 - (BufferFaktor2 & $3F))    ;64-er Alignment
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferFaktor2A
  !pop [BufferFaktor2A]
  ;========================================================
  ;Faktoren bearbeiten
  !tzcnt rax,[N]             ;BMI1, auch LZCNT möglich    
  !mov [Exponent],rax

  !mov rdx,32
  !sub rdx,rax
  !vmovq xmm9,rdx            ;32 - Exponent

  !mov rsi,[N]
  !mov r8,[PF1]
  !mov r10,[PF2]
  !mov r11,[BufferCosSinA]   ;Temp-Memory

  !mov r9,r11
  !add r9,rsi
  !add r9,rsi

  !mov r14,[L1]
  !mov rcx,[v_Asc_Uni]
  !dec rcx
  !shl r14,cl
  !sub r14,[v_Asc_Uni]

  !mov r15,[L2]
  !shl r15,cl
  !sub r15,[v_Asc_Uni]

  !mov rcx,2
  !sub rcx,[v_Asc_Uni]
  !shr rsi,cl   
  !sub rsi,[v_Asc_Uni]

  !mov rcx,7
 !@@:
  !mov [r11+rcx*4],ecx
  !dec rcx
 !jns @b 
  !vmovdqa ymm0,[r11]        ;Start-Values 0,1,2,3,4,5,6,7

  !vpbroadcastd ymm1,[V1]    ;AVX2
  !vpbroadcastd ymm2,[V2]
  !vpbroadcastd ymm3,[V3]
  !vpbroadcastd ymm4,[V4]
  !vpbroadcastd ymm11,[VAcht]

  !mov rdi,[BufferFaktor1A]  ;Temp-Memory
  !xor r12,r12
  ;------------------------------------
 !.Loop_Rev_Fak:
  !vpsrld ymm5,ymm0,1        ;>> 1
  !vpand ymm6,ymm5,ymm1      ;& V1 
  !vpand ymm7,ymm0,ymm1      ;& V1 
  !vpslld ymm8,ymm7,1        ;<< 1
  !vpor ymm10,ymm6,ymm8      ;| 

  !vpsrld ymm5,ymm10,2       ;>> 2
  !vpand ymm6,ymm5,ymm2      ;& V2 
  !vpand ymm7,ymm10,ymm2     ;& V2 
  !vpslld ymm8,ymm7,2        ;<< 2
  !vpor ymm10,ymm6,ymm8      ;| 

  !vpsrld ymm5,ymm10,4       ;>> 4
  !vpand ymm6,ymm5,ymm3      ;& V3 
  !vpand ymm7,ymm10,ymm3     ;& V3 
  !vpslld ymm8,ymm7,4        ;<< 4
  !vpor ymm10,ymm6,ymm8      ;| 

  !vpsrld ymm5,ymm10,8       ;>> 8
  !vpand ymm6,ymm5,ymm4      ;& V4 
  !vpand ymm7,ymm10,ymm4     ;& V4 
  !vpslld ymm8,ymm7,8        ;<< 8
  !vpor ymm10,ymm6,ymm8      ;|

  !vpsrld ymm5,ymm10,16      ;>> 16
  !vpslld ymm8,ymm10,16      ;<< 16
  !vpor ymm10,ymm5,ymm8      ;|

  !vpsrld ymm5,ymm10,xmm9    ;>> (32 - Exponent)
  !vmovdqa [rdi],ymm5

  !xor r13,r13
 !.Next_Byte8:
  !mov r12d,[rdi+r13*4]      ;8 reverse Positionen  
  !xor rax,rax
  !mov rdx,r14
  !mov rcx,[v_Asc_Uni]
  !dec rcx
  !shr rdx,cl
  !cmp [L1],rdx
 !jb @f 

  !movzx eax,byte[r8+r14]    ;next Byte Faktor1$
  !and eax,0fh               ;for String (=sub 30h)
 !@@:
  !mov [r9+r12*2],eax

  !xor rax,rax
  !mov rdx,r15
  !mov rcx,[v_Asc_Uni]
  !dec rcx
  !shr rdx,cl
  !cmp [L2],rdx
 !jb @f 
  !movzx eax,byte[r10+r15]   ;next Byte Faktor2$
  !and eax,0fh
 !@@:
  !mov [r11+r12*2],eax

  !sub rsi,[v_Asc_Uni]
  !sub r14,[v_Asc_Uni]
  !sub r15,[v_Asc_Uni]

  !inc r13
  !cmp r13,8
 !jne .Next_Byte8

  !vpaddd ymm0,ymm0,ymm11

  !or rsi,rsi
 !jns .Loop_Rev_Fak

  ;========================================================
  ;Integer in Double konvertieren, verdoppeln und abspeichern
  !xor rdx,rdx
  !xor r8,r8
  !mov rcx,[N]
  !shr rcx,3
  !mov r10,[BufferFaktor1A]
  !mov r11,[BufferCosSinA]
  !mov r12,[BufferFaktor2A]

  !mov r9,r11
  !mov rsi,[N]
  !add r9,rsi
  !add r9,rsi

 !@@:
  ;Faktor1
  !vcvtdq2pd ymm0,[r9+r8]    ;4 Ziffern (als DWords) von Faktor1 in 4 Doubles konvertieren 
  !vmovsd [r10+rdx],xmm0
  !vmovsd [r10+rdx+16],xmm0  ;+8 ist der I-Anteil, bleibt hier Null!
  !vpsrldq ymm1,ymm0,8
  !vmovsd [r10+rdx+32],xmm1
  !vmovsd [r10+rdx+48],xmm1
  !vperm2f128 ymm2,ymm0,ymm0,1
  !vmovsd [r10+rdx+64],xmm2
  !vmovsd [r10+rdx+80],xmm2
  !vpsrldq ymm1,ymm2,8
  !vmovsd [r10+rdx+96],xmm1
  !vmovsd [r10+rdx+112],xmm1
  ;Faktor2
  !vcvtdq2pd ymm0,[r11+r8]
  !vmovsd [r12+rdx],xmm0
  !vmovsd [r12+rdx+16],xmm0
  !vpsrldq ymm1,ymm0,8
  !vmovsd [r12+rdx+32],xmm1
  !vmovsd [r12+rdx+48],xmm1
  !vperm2f128 ymm2,ymm0,ymm0,1
  !vmovsd [r12+rdx+64],xmm2
  !vmovsd [r12+rdx+80],xmm2
  !vpsrldq ymm1,ymm2,8
  !vmovsd [r12+rdx+96],xmm1
  !vmovsd [r12+rdx+112],xmm1

  !add r8,16
  !add rdx,128
  !sub rcx,1
 !jnz @b
  ;========================================================
  ;Winkel erst hier, weil BufferCosSinA oben "missbraucht" wird
  !vcvtsi2sd xmm0,xmm0,qword[N]
  !vmovsd xmm1,[Pi2]
  !vdivsd xmm7,xmm1,xmm0     ;Schrittweite Konstante
  !vmovq xmm6,qword[Signum]  ;Konstante

  !mov rcx,[N]
  !shr rcx,3
  !xor r8,r8                 ;Schritt-Nr. Konstante
  !mov r9,[N]
  !shr r9,2                  ;N/4 für Adressierung unten
  !mov r10,[N]
  !shr r10,1                 ;N/2

  !mov r11,[BufferCosSinA]
  !inc rcx
 !@@:
  !vcvtsi2sd xmm15,xmm14,r8
  !vmulsd xmm0,xmm7,xmm15

  !vbroadcastsd ymm1,xmm0    ;YMM1: 1=x^1 2=x^1 3=x^1 4=x^1   AVX2, da Quelle Register ist

  !vmulpd ymm2,ymm1,ymm1     ;YMM2: 1=x^2 2=x^2 3=x^2 4=x^2
  !vmulsd xmm3,xmm2,xmm2     ;YMM3: 1=x^4 2=x^2 3=0   4=0
  !vmulpd ymm4,ymm2,ymm2     ;YMM4: 1=x^4 2=x^4 3=x^4 4=x^4
  !vmulpd ymm2,ymm4,ymm4     ;YMM2: 1=x^8 2=x^8 3=x^8 4=x^8  for the next 4 terms

  !vmovdqa xmm8,xmm3         ;save for Cos
  !vmulpd xmm9,xmm8,xmm4 
  !vperm2f128 ymm10,ymm9,ymm8,100000b

  !vmulpd xmm14,xmm3,xmm1              ;YMM14: 1=x^5 2=x^3 3=0 4=0
  !vmulpd xmm1,xmm14,xmm4              ;YMM1: 1=x^9 2=x^7 3=0 4=0
  !vperm2f128 ymm5,ymm1,ymm14,100000b  ;YMM5: 1=x^9 2=x^7 3=x^5 4=x^3  YMM5(0-127)=YMM1(0-127), YMM5(128-255)=YMM14(0-127)  

  !vmulpd ymm1,ymm5,yword[RezFakSin]   ;multiply the 4 values in YMM5 with RezFak -1/3! ... 1/9!, result in YMM1
  !vmulpd ymm11,ymm10,yword[RezFakCos]
  ;next 4 terms, without loop 
  !vmulpd ymm14,ymm5,ymm2              ;YMM14: 1=x^17 2=x^15 3=x^13 4=x^11
  !vmulpd ymm10,ymm10,ymm2

  !vfmadd231pd ymm1,ymm14,yword[RezFakSin+32]
  !vfmadd231pd ymm11,ymm10,yword[RezFakCos+32]
  ;next 4 terms
  !vmulpd ymm5,ymm14,ymm2              ;YMM5: 1=x^25 2=x^23 3=x^21 4=x^19
  !vmulpd ymm10,ymm10,ymm2

  !vfmadd231pd ymm1,ymm5,yword[RezFakSin+64]
  !vfmadd231pd ymm11,ymm10,yword[RezFakCos+64]
  ;result
  !vhaddpd ymm2,ymm1,ymm1              ;YMM2: 1=1+2 of YMM1, 3=3+4 of YMM1 
  !vhaddpd ymm12,ymm11,ymm11

  !vextractf128 xmm1,ymm2,1b           ;XMM1: 3+4 of YMM2
  !vextractf128 xmm11,ymm12,1b

  !vaddsd xmm3,xmm2,xmm1               ;XMM3: 1=sum of iterations 
  !vaddsd xmm13,xmm12,xmm11

  !vaddsd xmm12,xmm0,xmm3              ;Sine   XMM12: 1=sum of iterations plus start-value (1.term=x)
  !vaddsd xmm14,xmm13,[Eins]           ;Cosine
  ;------------------------------------
  ;store values 0-180°
  !mov rax,r8
  !shl rax,1
  !vmovsd qword[r11+rax*8+8],xmm12     ;Sin 0-45°
  !vorpd xmm1,xmm12,xmm6               ;neg Sin 
  !vmovsd qword[r11+rax*8],xmm14       ;Cos 0-45°

  !mov rax,r9
  !sub rax,r8
  !shl rax,1
  !vmovsd qword[r11+rax*8],xmm12       ;Cos 45-90°
  !mov rdx,r9
  !add rdx,r8
  !shl rdx,1
  !vmovsd qword[r11+rdx*8],xmm1        ;Cos 90-135°
  !mov r13,r10
  !sub r13,r8
  !shl r13,1
  !vorpd xmm2,xmm14,xmm6
  !vmovsd qword[r11+r13*8],xmm2        ;Cos 135-180°
  !vmovsd qword[r11+rax*8+8],xmm14     ;Sin 45-90° 
  !vmovsd qword[r11+rdx*8+8],xmm14     ;Sin 90-135°
  !vmovsd qword[r11+r13*8+8],xmm12     ;Sin 135-180°

  !inc r8
  !dec rcx
 !jnz @b
  ;========================================================
  ;jetzt die beiden Faktoren transformieren
  !mov rsi,[N]
  !mov r8,[BufferCosSinA]

  !mov rcx,1                           ;Exponent_FFT
  !mov r12,2                           ;Pointer0
  !mov r9,2                            ;Pointer2

  !mov r10,[BufferFaktor1A]
  !mov rdi,[BufferFaktor2A]

 !.OuterLoopF:
  !inc rcx                             ;Exponent_FFT
  !mov rax,rsi                         ;N
  !shr rax,cl
 !jz .EndF
  !mov r14,rax

   !.InnerLoopF:
    !xor r11,r11                       ;Pointer1
    !mov r15,r12                       ;Pointer0
  ;------------------------------------
     !@@:
      !mov rax,r14
      !mul r11                         ;Pointer1
      !shl rax,3                       ;EW*8, Pointer in CosSin

      !mov rdx,r9                      ;Pointer2
      !mov r13,rdx
      !shl rdx,4                       ;Pointer2*8
      !sub r13,r12                     ;[v_Pointer0]
      !shl r13,4                       ;(Pointer2 - Pointer0) * 8 

      !vperm2f128 ymm0,ymm9,yword[r8+rax],100010b     ;2x CosSinA in ymm0, ymm9 is Dummy

      !vmovapd xmm7,dqword[rdi+rdx]    ;BufferFaktor2A complete (R and I)
      !vperm2f128 ymm2,ymm7,yword[r10+rdx],10b

      !vmulpd ymm4,ymm2,ymm0
      !vshufpd ymm1,ymm0,ymm0,101b     ;change Sin and Cos
      !vhsubpd ymm3,ymm4,ymm4          ;WR
      
      !vmulpd ymm4,ymm2,ymm1
      !vhaddpd ymm5,ymm4,ymm4          ;WI
      !vshufpd ymm15,ymm3,ymm5,0b      ;WI in ymm15-high
      !vmovapd xmm1,dqword[r10+r13]
      !vmovapd xmm7,dqword[rdi+r13]
      !vperm2f128 ymm2,ymm7,ymm1,10b
      !vaddpd ymm14,ymm15,ymm2
      !vsubpd ymm13,ymm2,ymm15

      !vmovapd [r10+r13],xmm14

      !vperm2f128 ymm2,ymm7,ymm14,11b
      !vmovapd [rdi+r13],xmm2

      !vmovapd [r10+rdx],xmm13

      !vperm2f128 ymm2,ymm7,ymm13,11b
      !vmovapd [rdi+rdx],xmm2

      !add r11,2                       ;Pointer1
      !inc r9                          ;Pointer2
      !dec r15
     !jnz @b

    !add r9,r12                        ;Pointer0
    !cmp r9,rsi                        ;N
   !jb .InnerLoopF

  !shl r12,1
  !mov r9,r12

 !jmp .OuterLoopF

 !.EndF:
  ;========================================================
  ;Und jetzt die beiden transformierten Vektoren miteinander multiplizieren
  !mov r8,[BufferFaktor1A]        ;nimmt auch das Produkt auf
  !mov r9,[BufferFaktor2A]
  !xor r11,r11                    ;Pointer
  !mov rcx,[N]
  !shr rcx,1
  !vmovdqu ymm5,[Null1]           ;Null,Signum,Null,Signum
 !@@:
  !vmovapd ymm0,[r8+r11]
  !vmovapd ymm3,[r9+r11]
  !vmulpd ymm1,ymm0,ymm3          ;ymm1:R*R,I*I,R*R,I*I
  !vshufpd ymm2,ymm3,ymm3,101b    ;change R and I
  !vxorpd ymm4,ymm1,ymm5          ;Vorzeichen ändern, damit nicht sub sondern add für beide Terme verwendet werden kann
  
  !vmulpd ymm1,ymm0,ymm2          ;R*I, R*I
  !vhaddpd ymm6,ymm4,ymm1
  !vmovapd [r8+r11],ymm6

  !add r11,32
  !dec rcx
 !jnz @b
  ;========================================================
  ;Rück-Transformation
  ;Produkt bearbeiten
  BufferReverse = AllocateMemory(32 + 64)
  If BufferReverse
    BufferReverseA = BufferReverse + (64 - (BufferReverse & $3F))    ;64-er Alignment
   Else 
    sAusgabe = "-1"
   ProcedureReturn sAusgabe
  EndIf
  PUSH BufferReverseA
  !pop [BufferReverseA]

  !mov rsi,[N]               ;für Abbruch Schleife Loop_Rev_Prod; letzte Pos.steht ja fest (wie erste) 
  !shl rsi,1
  !sub rsi,2
  !mov r8,[BufferFaktor1A]   ;Ergebnis von Multiplikation oben
  !mov r9,[BufferFaktor2A]   ;hier neue Reihenfolge rein

  !mov rdx,32
  !sub rdx,[Exponent]
  !vmovq xmm9,rdx

  !mov rcx,7
 !@@:
  !mov [r9+rcx*4],ecx
  !dec rcx
 !jns @b 
  !vmovdqa ymm0,[r9]         ;Startwerte (8 DWords=32 Bytes) 0,1,2,3,4,5,6,7   r9 hier "missbraucht"

  !vpbroadcastd ymm1,[V1]    ;AVX2
  !vpbroadcastd ymm2,[V2]
  !vpbroadcastd ymm3,[V3]
  !vpbroadcastd ymm4,[V4]
  !vpbroadcastd ymm11,[VAcht]

  !mov rdi,[BufferReverseA]
  !xor r12,r12
  ;------------------------------------
 !.Loop_Rev_Prod:
  !vpsrld ymm5,ymm0,1        ;>> 1
  !vpand ymm6,ymm5,ymm1      ;& V1 
  !vpand ymm7,ymm0,ymm1      ;& V1 
  !vpslld ymm8,ymm7,1        ;<< 1
  !vpor ymm10,ymm6,ymm8      ;| 

  !vpsrld ymm5,ymm10,2       ;>> 2
  !vpand ymm6,ymm5,ymm2      ;& V2 
  !vpand ymm7,ymm10,ymm2     ;& V2 
  !vpslld ymm8,ymm7,2        ;<< 2
  !vpor ymm10,ymm6,ymm8      ;| 

  !vpsrld ymm5,ymm10,4       ;>> 4
  !vpand ymm6,ymm5,ymm3      ;& V3 
  !vpand ymm7,ymm10,ymm3     ;& V3 
  !vpslld ymm8,ymm7,4        ;<< 4
  !vpor ymm10,ymm6,ymm8      ;| 

  !vpsrld ymm5,ymm10,8       ;>> 8
  !vpand ymm6,ymm5,ymm4      ;& V4 
  !vpand ymm7,ymm10,ymm4     ;& V4 
  !vpslld ymm8,ymm7,8        ;<< 8
  !vpor ymm10,ymm6,ymm8      ;|

  !vpsrld ymm5,ymm10,16      ;>> 16
  !vpslld ymm8,ymm10,16      ;<< 16
  !vpor ymm10,ymm5,ymm8      ;|

  !vpsrld ymm5,ymm10,xmm9    ;>> (32 - Exponent)
  !vmovdqa [rdi],ymm5

  !xor rcx,rcx
 !@@:
  !mov r12d,[rdi+rcx*4]
  !shl r12,1
  !vmovapd xmm12,[r8]        ;16 Bytes (R.d und I.d) 
  !vmovapd [r9+r12*8],xmm12
  !add r8,16
  !inc rcx
  !cmp rcx,8
 !jne @b

  !vpaddd ymm0,ymm0,ymm11

  !cmp r12,rsi               ;Position letzter Wert?
 !jne .Loop_Rev_Prod
  ;========================================================
  !mov r8,[BufferFaktor2A]
  !mov rcx,[N]
  !vpbroadcastq ymm4,qword[Signum]     ;Signum kann auch hier für Gather verwendet werden
  !vmovdqu ymm6,[Adresse2]             ;0,8,32,40      auch add abspeichern
  !vmovdqu ymm7,[Adresse3]             ;16,24,48,56    auch sub abspeichern
  !vmovdqa ymm5,ymm4                   ;ymm4 sichern, wird auf Null gesetzt!
  !shr rcx,2
 !@@:
  !vgatherqpd ymm1,[r8+ymm6],ymm4      ;ymm4 s.o.
  !vmovdqa ymm4,ymm5
  !vgatherqpd ymm2,[r8+ymm7],ymm4

  !vaddpd ymm8,ymm1,ymm2
  !vmovapd [r8],xmm8
  !vextractf128 [r8+32],ymm8,1b

  !vsubpd ymm9,ymm1,ymm2
  !vmovapd [r8+16],xmm9
  !vextractf128 [r8+48],ymm9,1b

  !vmovdqa ymm4,ymm5

  !add r8,64
  !dec rcx
 !jnz @b
  ;------------------------------------

  ;------------------------------------
  !mov rsi,[N]
  !mov r8,[BufferCosSinA]
  !mov rcx,1                           ;Exponent_FFT
  !mov r12,2                           ;Pointer0 
  !mov r9,2                            ;Pointer2
  !mov r10,[BufferFaktor2A]

 !.OuterLoopP:
  !inc rcx                             ;Exponent_FFT
  !mov rax,rsi                         ;N
  !shr rax,cl
 !jz .EndP

  !mov r14,rax
  !shl r14,3                           ;EW*8
   !.InnerLoopP:
    !xor r11,r11                       ;Pointer1
    !mov r15,r12
  ;------------------------------------
     !@@:
      !mov rax,r14
      !mul r11                         ;Pointer1    oder imul rax,r11  ;dann wäre rdx frei

      !mov rdx,r9                      ;Pointer2
      !mov r13,rdx
      !shl rdx,4
      !sub r13,r12                     ;[v_Pointer0]
      !shl r13,4

      !vmovapd xmm7,dqword[r10+r13]

      !vmovapd xmm0,dqword[r10+rdx]    ;BufferFaktor2A
      !vmovapd xmm1,dqword[r8+rax]     ;CosSinA
      !vshufpd xmm2,xmm0,xmm0,1b       ;R und I vertauschen

      !vmulpd xmm3,xmm2,xmm1
      !vmulpd xmm4,xmm0,xmm1

      !vxorpd xmm5,xmm3,[Null1]        ;Signum, xmm3 negieren um vhaddpd verwenden zu können (sonst sub) wie xmm4

      !vhaddpd xmm6,xmm4,xmm5          ;WI in xmm15-high, WR in Low

      !vaddpd xmm14,xmm7,xmm6          ;dqword[r10+r13]
      !vsubpd xmm13,xmm7,xmm6

      !vmovapd [r10+r13],xmm14
      !vmovapd [r10+rdx],xmm13

      !add r11,2                       ;Pointer1
      !inc r9                          ;Pointer2

      !dec r15
     !jnz @b

;---------------------------------------
    !add r9,r12                        ;Pointer0
    !cmp r9,rsi                        ;N
   !jb .InnerLoopP

  !shl r12,1
  !mov r9,r12

 !jmp .OuterLoopP
 !.EndP:
  ;========================================================
  FreeMemory(BufferReverse)
  FreeMemory(BufferCosSin)
  FreeMemory(BufferFaktor1)

  BufferProdukt = AllocateMemory(N + 64)

  If BufferProdukt
    BufferProduktA = BufferProdukt + (64 - (BufferProdukt & $3F))    ;64-er Alignment
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferProduktA
  !pop [BufferProduktA]

  BufferDez = AllocateMemory(16 + 64)
  If BufferDez
    BufferDezA = BufferDez + (64 - (BufferDez & $3F)) ;64-er Alignment
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferDezA
  !pop [BufferDezA]
  ;========================================================
  ;Ergebnis für Ausgabe aufbereiten
  ;Konvertierung Double in Integer, Rundung gemäß Rounding Control Bits im MXCSR-Register, oben auf "nearest" gesetzt

  ;da BufferFaktor2A noch den I-Anteil enthält, wird hier "komprimiert", d.h. in die I-Bereiche werden die Produkte reinkopiert
  ;die Double-Werte in Long-Integer konvertieren, Division ist notwendig für höhere Genauigkeit; nicht später shiften!
  !vpbroadcastq ymm4,qword[Signum]     ;Signum kann auch hier für Gather verwendet werden
  !vmovdqu ymm6,[Adresse1]             ;Adressen der Real-Teile  0,16,32,48
  !vmovdqa ymm5,ymm4                   ;ymm4 wird auf Null gesetzt!
  !mov r15,[N] 
  !mov r14,[BufferFaktor2A]
  !shr r15,2
  !mov rcx,r14
  !vbroadcastsd ymm1,[ND]              ;"nur" AVX, da Quelle Memory ist
 !@@:
  !vgatherqpd ymm0,[r14+ymm6],ymm4     ;ymm4 s.o.
  !vdivpd ymm2,ymm0,ymm1
  !vcvtpd2dq xmm3,ymm2
  !vmovdqa [rcx],xmm3
  !vmovdqa ymm4,ymm5
  !add r14,64
  !add rcx,16
  !dec r15
 !jnz @b

  ;Hex2Dez und "Verstringung"
  ;Pointer_Produkt = N - 1
  ;For i = 0 To N - 1         ;vom Produkt alle Werte durch
  !mov r15,[N]
  !mov r14,[BufferFaktor2A]
  !mov rdx,[BufferProduktA]
  !mov rdi,r15
  !dec rdi
  !lea r10,[Dezi]
  !mov r11,[BufferDezA]

 !.DeziLoop: 
    !mov eax,[r14]           ;eax=Integer-Wert
    ;Hex2Dez
    ;Ziffer = 0
    !xor r12b,r12b 
    ;Pointer1 = 0
    !xor r8,r8
    ;Pointer2 = 0
    !xor r9,r9
 
    ;While PeekQ(?Dezi + Pointer1) <> 0 oder Zähler wie jetzt
    !mov rsi,10               ;Anzahl Dezi-Werte ohne Null, 9 reicht für M49, ansonsten auch Dezi anpassen an Faktoren(Produkt)-Länge!
   !.OuterLoopDez1:
      ;While (QuadWert - PeekQ(?Dezi + Pointer1)) >= 0
      !mov r13d,eax
     !.InnerLoopDez1:
        !sub r13d,[r10+r8]
       !js .InnerLoopDez1End
        ;Ziffer + 1
        !inc r12b
        ;QuadWert - PeekQ(?Dezi + Pointer1)
        ;!sub rax,[r10+r8]
        !mov eax,r13d         ;AL ist unten letzte Ziffer
      ;Wend     (QuadWert - PeekQ(?Dezi + Pointer1)) >= 0
      !jmp .InnerLoopDez1

      !.InnerLoopDez1End:
      ;PokeB(BufferDezA + Pointer2, Ziffer)
      !mov [r11+r9],r12b
      ;Pointer1 + 8
      !add r8,8
      ;Pointer2 + 1
      !inc r9
      ;Ziffer = 0
      !xor r12b,r12b 
      !dec rsi
    ;Wend      PeekQ(?Dezi + Pointer1) <> 0 oder Zähler wie jetzt
    !jnz .OuterLoopDez1

    ;PokeB(BufferDezA + Pointer2, QuadWert1)   ;letzte Dezimal-Ziffer abspeichern
    !mov [r11+r9],al         ;AL ist letzte Ziffer von oben
    ;Pointer_Produkt1 = Pointer_Produkt 
    !mov rsi,rdi
    !mov rax,10               ;Anzahl Dezi-Werte ohne Null, 9 reicht für M49, ansonsten auch Dezi anpassen an Faktoren(Produkt)-Länge!
    ;!mov rdx,[BufferProduktA]
    ;For k = 0 To 15
    !@@:
      ;Ziffer = PeekB(BufferProduktA + Pointer_Produkt1)
      !mov r12b,[rdx+rsi]
      ;Ziffer + PeekB(BufferDezA + Pointer2)
      !add r12b,[r11+r9]
      ;PokeB(BufferProduktA + Pointer_Produkt1, Ziffer)
      !mov [rdx+rsi],r12b
      ;If Ziffer > 9                           ;Übertrag
      !cmp r12b,9              ;Übertrag?
      !jbe .ZifferOK
        ;Ziffer - 10
        !sub r12b,10
        ;PokeB(BufferProduktA + Pointer_Produkt1, Ziffer)
        !mov [rdx+rsi],r12b
        ;Ziffer = PeekB(BufferProduktA + Pointer_Produkt1 - 1) + 1
        !mov r12b,[rdx+rsi-1]
        !inc r12b
        ;PokeB(BufferProduktA + Pointer_Produkt1 - 1, Ziffer)
        !mov [rdx+rsi-1],r12b
      ;EndIf
      !.ZifferOK:
      ;Pointer_Produkt1 - 1
      !dec rsi
      !js .Reicht
      ;Pointer2 - 1
      !dec r9
      !dec rax
    ;Next
    !jnz @b

    !.Reicht:
    ;Pointer_Produkt - 1
    !dec rdi
    !add r14,4
    !dec r15
  ;Next
  !jnz .DeziLoop

  ;------------------------
  ;Ergebnis (Produkt) als ASCII-Ziffern abspeichern
  !vpbroadcastq ymm0,[AsciiNull]  ;load 32x48 ("Null")   AVX2
  !mov r9,[BufferProduktA]
  !xor r10,r10
  !mov rcx,[N]
  !shr rcx,5                      ;/32, YMM-Length
 !@@:
  !vpor ymm1,ymm0,[r9+r10]        ;AVX2
  !vmovdqa [r9+r10],ymm1
  !add r10,32
  !dec rcx
 !jnz @b
  !mov byte[r9+r10],0             ;Zero-Byte für String, zur Sicherheit, evtl.word

  FreeMemory(BufferFaktor2)
  FreeMemory(BufferDez)

  sAusgabe = LTrim(PeekS(BufferProduktA + N - LG, LG, #PB_Ascii), "0")    ;evtl. führende Null(en) entfernen

  If sAusgabe = ""                ;mindestens einer der Faktoren war Null, Kompatibilität zu "normaler" Multiplikation
    sAusgabe = "0"
  EndIf
  FreeMemory(BufferProdukt)

  ;Restore Registers
  !mov rax,[BufferFXA]
  !fxrstor64 [rax]
  !mov rdi,[rax+464]         ;restore callee-save registers
  !mov rsi,[rax+472]
  !mov r12,[rax+480]
  !mov r13,[rax+488]
  !mov r14,[rax+496]
  !mov r15,[rax+504]
  FreeMemory(BufferFX)

  ;Restore MXCSR             ;nach fxrstor64!
  !vldmxcsr [MXCSR_Sicher]   ;schreibt den gesicherten Wert zurück ins MXCSR-Register

  !vzeroupper
  DisableASM
 ProcedureReturn sAusgabe

  ;========================================================
  DataSection
!RezFakSin:                       ;für Sinus
  !dq  2.755731922398589065e-6    ; 1/9!   4.Iteration
  !dq -1.984126984126984127e-4    ;-1/7!   3.Iteration
  !dq  8.333333333333333333e-3    ; 1/5!   2.Iteration
  !dq -1.666666666666666667e-1    ;-1/3!   1.Iteration

  !dq  2.811457254345520763e-15   ; 1/17!  8.Iteration
  !dq -7.647163731819816476e-13   ;-1/15!  7.Iteration
  !dq  1.605904383682161460e-10   ; 1/13!  6.Iteration
  !dq -2.505210838544171878e-8    ;-1/11!  5.Iteration

  !dq  6.446950284384473396e-26   ; 1/25!  12.Iteration
  !dq -3.868170170630684038e-23   ;-1/23!  11.Iteration
  !dq  1.957294106339126123e-20   ; 1/21!  10.Iteration
  !dq -8.220635246624329717e-18   ;-1/19!  9.Iteration

!RezFakCos:                       ;für Cosinus
  !dq  2.480158730158730159e-5    ; 1/8!   4.Iteration
  !dq -1.388888888888888889e-3    ;-1/6!   3.Iteration
  !dq  4.166666666666666667e-2    ; 1/4!   2.Iteration
  !dq -5.000000000000000000e-1    ;-1/2!   1.Iteration

  !dq  4.779477332387385297e-14   ; 1/16!  8.Iteration
  !dq -1.147074559772972471e-11   ;-1/14!  7.Iteration
  !dq  2.087675698786809898e-9    ; 1/12!  6.Iteration
  !dq -2.755731922398589065e-7    ;-1/10!  5.Iteration

  !dq  1.611737571096118349e-24   ; 1/24!  12.Iteration
  !dq -8.896791392450573287e-22   ;-1/22!  11.Iteration
  !dq  4.110317623312164858e-19   ; 1/20!  10.Iteration
  !dq -1.561920696858622646e-16   ;-1/18!  9.Iteration
!Eins:
  !dq 1.0
!Pi2:
  !dq 6.283185307179586476925
!Null1:
  !dq 0h 
!Signum:
  !dq 8000000000000000h           ;für neg Sin/Cos und Gather
  !dq 0h
  !dq 8000000000000000h
!AsciiNull:
  !dq 3030303030303030h           ;für String
!Dezi:
  !dq  10000000000
  !dq   1000000000                ;für M49 ist der größte Integer-Wert etwas über 200.000.000
  !dq    100000000
  !dq     10000000
  !dq      1000000
  !dq       100000
  !dq        10000
  !dq         1000
  !dq          100
  !dq           10
  !dq            0                ;blanke Sicherheit
!Adresse1:                         ;für Gather
  !dq 0
  !dq 16
  !dq 32
  !dq 48
!Adresse2:                         ;für Gather
  !dq 0
  !dq 8
  !dq 32
  !dq 40
!Adresse3:                         ;für Gather
  !dq 16
  !dq 24
  !dq 48
  !dq 56

  !BufferCosSinA   dq ?
  !BufferDezA      dq ?
  !BufferFaktor1A  dq ?
  !BufferFaktor2A  dq ?
  !BufferFXA       dq ?
  !BufferProduktA  dq ?
  !BufferReverseA  dq ?
  !Exponent        dq ?
  !L1              dq ?
  !L2              dq ?
  !LenF            dq ?
  !N               dq ?
  !ND              dq ?
  !PF1             dq ?
  !PF2             dq ?
  !Pointer_Produkt dq ?

  !MXCSR           dd ?
  !MXCSR_Sicher    dd ?

  !V1              dd 55555555h 
  !V2              dd 33333333h
  !V3              dd 0F0F0F0Fh
  !V4              dd 00FF00FFh
  !VAcht           dd 8
  EndDataSection
EndProcedure

Procedure Calc()             ;Aufruf durch Potenzen (Case 3) und Mersenne (Case 4)
  ZeitA = ElapsedMilliseconds()

  Stellen = Potenz * Len(Result$) + 1
  Stellen = Stellen + (64 - (Stellen & $3F))

  !mov rax,[v_Potenz]
  !bsr rdx,rax
  !mov [v_Bits_M],rdx

  CalcPot = Potenz           ;CalcPot wird geshiftet!

  BufferPot = AllocateMemory((Stellen * 2 * #Asc_Uni) + 64)
  If BufferPot
    BufferPotA = BufferPot + (64 - (BufferPot & $3F)) ;64-er Align
   Else
    MessageRequester("Fehler!", "Nicht genügend freier Arbeitsspeicher!", #MB_ICONERROR)
    End
  EndIf

  ListViewGadget(2, 10, 10, 250, 580)  ;Info-Fenster

  PokeS(BufferPotA + Stellen * #Asc_Uni, Result$)     ;hier noch Result$=Basis$, zum besseren Verständnis
  PokeS(BufferPotA, "1")     ;falls Bit0 von Potenz nicht gesetzt, X^0=1

  !shr [v_CalcPot],1
 !jnc @f
  PokeS(BufferPotA, Basis$)
  m = 1
 !@@:

  p = 1
  j = 2
  k = Bits_M - 1
  z = 0                      ;Zeile im ListViewGadget (Info)
  For i = 0 To k
    AddGadgetItem (2, -1, "Berechne " + Basis$ + "^" + Str(j))   
    SetGadgetState(2, z)
    z + 1
    Result$ = FFT_Mul(BufferPotA + Stellen * #Asc_Uni, BufferPotA + Stellen * #Asc_Uni)  ;Potenzen berechnen
    If Result$ = "-1"
      MessageRequester("Fehler!", "Nicht genügend freier Arbeitsspeicher!", #MB_ICONERROR)
      End
    EndIf
    PokeS(BufferPotA + Stellen * #Asc_Uni, Result$)
    j * 2
    !shr [v_CalcPot],1
   !jnc @f
    n + p
    AddGadgetItem (2, -1, "Berechne " + Basis$ + "^" + Str(n * 2 + m))   
    SetGadgetState(2, z)
    z + 1
    Result$ = FFT_Mul(BufferPotA, BufferPotA + Stellen * #Asc_Uni)   ;(Zwischen-) Ergebnis berechnen
    If Result$ = "-1"
      MessageRequester("Fehler!", "Nicht genügend freier Arbeitsspeicher!", #MB_ICONERROR)
      End
    EndIf
    If i < k
      PokeS(BufferPotA, Result$)
    EndIf
   !@@:
    p * 2
  Next

  FreeMemory(BufferPot)
  FreeGadget(2)
  ZeitE = ElapsedMilliseconds() - ZeitA
EndProcedure

Procedure.i RangeCharSSE42(StrPointer, ChrPointer)
  !mov rax,[p.v_ChrPointer]
  !movdqu xmm1,[rax]
  !mov rdx,[p.v_StrPointer]
  !mov rax, -16

 CompilerIf #PB_Compiler_Unicode                 ;Unicode
  !align 16                                      ;kann speed-mäßig was bringen
!@@:
  !add rax,16
  !pcmpistri xmm1,dqword[rdx+rax],00000101b      ;Bit0 = 1 und Bit1 = 0, also String-Chars werden als unsigned Words interpretiert, Bit2 = 1 und Bit3 = 0, also Test auf Ranges
  !ja @b
  !shr rax,1                                     ;Pointer in String muss hier halbiert werden (Word)
  !add rax,rcx
  !shl rcx,1                                     ;weil sonst nur 8 (Words), unten wird aber mit 16 verglichen
  CompilerElse                                   ;kein Unicode
  !align 16                                      ;kann speed-mäßig was bringen
!@@:
  !add rax,16
  !pcmpistri xmm1,dqword[rdx+rax],00000100b      ;Bit0 = 0 und Bit1 = 0, also String-Chars werden als unsigned Bytes interpretiert, Bit2 = 1 und Bit3 = 0, also Test auf Ranges
  !ja @b
  !add rax,rcx
 CompilerEndIf

  !cmp rcx,16                                    ;rcx=16 oder 8, wenn nichts gefunden! 8 schon erweitert
  !je @f
  !inc rax 
 ProcedureReturn                                 ;rax = Pos., beginnt bei 1
!@@:
  !xor rax,rax
 ProcedureReturn                                 ;rax = 0, Char nicht gefunden
EndProcedure 

Procedure FakInfo()
  GadgetZeile + 1
  AddGadgetItem (2, -1, "Berechne aktuell " + Str(FakA) + "!")
  SetGadgetState(2, GadgetZeile)
EndProcedure

Viel Spaß!
Helle

Edit 7.1.2018: Kosmetik
Edit 11.1.2018: Kosmetik und minimale Speicherplatz-Reduzierung
Zuletzt geändert von Helle am 11.01.2018 21:33, insgesamt 5-mal geändert.
ccode_new
Beiträge: 1214
Registriert: 27.11.2016 18:13
Wohnort: Erzgebirge

Re: FFT-Calculator

Beitrag von ccode_new »

Schade meine CPU unterstützt kein AVX2.

Trotzdem bist du hier der Assembler- und Mathe-Meister. :wink:

:allright:
Betriebssysteme: div. Windows, Linux, Unix - Systeme

no Keyboard, press any key
no mouse, you need a cat
Benutzeravatar
RSBasic
Admin
Beiträge: 8022
Registriert: 05.10.2006 18:55
Wohnort: Gernsbach
Kontaktdaten:

Re: FFT-Calculator

Beitrag von RSBasic »

Mein i5 2500K leider auch nicht.
Aus privaten Gründen habe ich leider nicht mehr so viel Zeit wie früher. Bitte habt Verständnis dafür.
Bild
Bild
Benutzeravatar
Kiffi
Beiträge: 10621
Registriert: 08.09.2004 08:21
Wohnort: Amphibios 9

Re: FFT-Calculator

Beitrag von Kiffi »

@Helle:

wenn Du #MB_ICONERROR und die Windows-Api Aufrufe (GetDeviceCaps_() & GetDC_()) weglässt, dann funktioniert es auch unter Linux ;-)

Das EditorGadget() unter Linux reagiert anscheinend allergisch auf zu viele Zeichen in einer Zeile. Ab einer gewissen Länge werden alle Zeichen ineinander geschoben.

Grüße ... Peter
Hygge
Benutzeravatar
NicTheQuick
Ein Admin
Beiträge: 8675
Registriert: 29.08.2004 20:20
Computerausstattung: Ryzen 7 5800X, 32 GB DDR4-3200
Ubuntu 22.04.3 LTS
GeForce RTX 3080 Ti
Wohnort: Saarbrücken
Kontaktdaten:

Re: FFT-Calculator

Beitrag von NicTheQuick »

Schade, ich hab auch kein AVX2, nur AVX. Aber der Laptop hier ist auch schon ein paar Jährchen alt. Dennoch leistet der Intel(R) Core(TM) i7-3820QM CPU @ 2.70GHz hier drin immer noch gute Dienste. Ist fast immer noch so leistungsstark wie die heutigen mobilen i7 laut den CPU-Benchs im Internet.
Bild
Benutzeravatar
Helle
Beiträge: 566
Registriert: 11.11.2004 16:13
Wohnort: Magdeburg

Re: FFT-Calculator

Beitrag von Helle »

Da wohl doch viele User (auch die Forums-Prominenz :mrgreen: ) noch CPUs ohne AVX2 haben, hier die Version mit AVX. Auch BMI und FMA habe ich ersetzt. Testen kann ich es z.Z. nicht, deshalb wäre eine Rückmeldung nett. Ansonsten gilt das oben geschriebene; es fehlt noch Feinschliff.

Code: Alles auswählen

;- Schnelle Unsigned Integer-Multiplication mit FFT
;- see https://en.wikipedia.org/wiki/Butterfly_diagram (radix-2 decimation-in-time FFT algorithm)
;- Windows 64-Bit, getestet mit PureBasic 5.46 LTS (x64) / PureBasic 5.61 (x64)
;- Ascii und Unicode
;- Benötigt CPU mit AVX-Unterstützung
;- "Helle" Klaus Helbing, 05.01.2018: Mersenne M50 zugefügt
;- Um in den Genuss der Info-Anzeige zu kommen, sollten Programme wie der Taskmanager nicht gestartet sein!
;- 6.1.2018 AVX-Variante, BMI und FMA ersetzt
;- Version vom 11.01.2018

Global.q Asc_Uni, Bits_M, CalcPot, FakA, GadgetZeile, Potenz, ZeitA, ZeitE
Global.s Basis$, Result$

Declare.s FFT_Mul(Para1, Para2) 
Declare Calc()
Declare RangeCharSSE42(Para1, Para2)
Declare FakInfo()

;Test auf AVX
!mov eax,1
!cpuid
!test ecx,10000000h     ;Bit28 für AVX
!jnz @f                 ;Yeah!
MessageRequester("Schade!", "Diese CPU unterstützt kein AVX!" + #LFCR$ + "Ende")
End 
!@@:

FensterTitel$ = "Helles FFT-Calculator AVX-Version"
OpenWindow(0, 0, 0, 800, 600, FensterTitel$, #PB_Window_SystemMenu | #PB_Window_ScreenCentered)

CompilerIf #PB_Compiler_Unicode
  Asc_Uni = 2 : #Asc_Uni = 2 ;mit Konstante für Step usw.
 CompilerElse
  Asc_Uni = 1 : #Asc_Uni = 1
CompilerEndIf

FontHigh = Int(9.0 / (GetDeviceCaps_(GetDC_(WindowID(0)), #LOGPIXELSY) / 96.0))     ;GetDeviceCaps kann bei der Bildschirmanpassung nützlich sein
;LoadFont(0, "Arial", FontHigh)
LoadFont(0, "Trebuchet MS Fett", FontHigh)       ;in Windows 7 verfügbar; evtl.testen

Repeat
  TextGadget(1, 20, 20, 300, 15, "Was soll berechnet werden:")
  OptionGadget(2, 20, 45, 780, 15, "Multiplikation zweier positiver Integer-Faktoren")
  OptionGadget(3, 20, 65, 780, 15, "Berechnung einer positiven ganzzahligen Potenz")
  OptionGadget(4, 20, 85, 780, 15, "Berechnung einer Mersenne-Prim-Zahl (M25...M50)")
  OptionGadget(5, 20, 105, 780, 15, "Fakultät")
  OptionGadget(6, 20, 125, 780, 15, "Neuner-Test")
  ButtonGadget(7, 300, 560, 200, 30, "W e i t e r")
  SetGadgetState(2, 1) : FFT = 2       ;Voreinstellung Multiplikation, falls nichts ausgewählt wird

  For i = 1 To 7
    SetGadgetFont(i, FontID(0))
  Next

  Repeat
    Event = WaitWindowEvent()
    If Event = #PB_Event_CloseWindow
      CloseWindow(0)
      End
    EndIf
    If Event = #PB_Event_Gadget
      Select EventGadget()
        Case 7
          Break          
        Default
          FFT = EventGadget()
      EndSelect
    EndIf
  ForEver

  For i = 1 To 7
    FreeGadget(i)
  Next

  Select FFT
;=========================================================================
    Case 2                   ;Multiplikation 
      Info0$ = "Multiplikation "
      SetWindowTitle(0, Info0$)

      CharRange$ = Chr(1) + Chr(47) + Chr(58) + Chr(255)   ;für Suche nach Chars, die keine Dezimal-Ziffer sind. Chr(0) geht nicht!

      TextGadget(1, 50, 0, 300, 20, "Faktor1 eintragen (z.B. einen String mit Strg+V)")
      TextGadget(2, 450, 0, 300, 20, "Faktor2 eintragen (z.B. einen String mit Strg+V)")
      ButtonGadget(3, 325, 560, 150, 30, "Berechnung starten")
      ButtonGadget(4, 125, 560, 150, 30, "Faktor1 überprüfen")
      ButtonGadget(5, 525, 560, 150, 30, "Faktor2 überprüfen")

      For i = 1 To 5
        SetGadgetFont(i, FontID(0))
      Next

      EditorGadget(6, 10, 20, 385, 530)
      EditorGadget(7, 405, 20, 385, 530)
      SetActiveGadget(6)

      Repeat
        Event = WaitWindowEvent()
        If Event = #PB_Event_CloseWindow
          CloseWindow(0)
          End 
        EndIf
        If Event = #PB_Event_Gadget
          Gadget = EventGadget()
          Select Gadget
            Case 3
              Break
            Case 4, 5
              If Gadget = 4
                Faktor$ = GetGadgetText(6)
               Else
                Faktor$ = GetGadgetText(7)
              EndIf
              RangeChar = RangeCharSSE42(@Faktor$, @CharRange$)
              If RangeChar
                Ascii = Asc(Mid(Faktor$, RangeChar, 1)) : Char$ = Mid(Faktor$, RangeChar, 1) 
                MessageRequester("Faktor" + Str(Gadget - 3) + ": Fehler!", "Digit an Stelle " + Str(RangeChar) + " ist keine Dezimal-Ziffer! (ASCII=" + Ascii + "), Char=" + Char$, #MB_ICONERROR)
               Else
                MessageRequester("Faktor" + Str(Gadget - 3) + ": O.K.!", "Alle Digits sind Dezimal-Ziffern!")
              EndIf
          EndSelect
        EndIf
      ForEver

      Faktor1$ = GetGadgetText(6)
      Faktor2$ = GetGadgetText(7)

      For i = 1 To 7
        FreeGadget(i)
      Next

      ZeitA = ElapsedMilliseconds()
      Result$ = FFT_Mul(@Faktor1$, @Faktor2$)
      ZeitE = ElapsedMilliseconds() - ZeitA
;=========================================================================
    Case 3                   ;Potenzen
      Info0$ = "Berechnung von X^Y "
      SetWindowTitle(0, Info0$)
      TextGadget(1, 50, 10, 300, 20, "Basis X eintragen für X^Y")
      TextGadget(2, 400, 10, 300, 20, "Exponent Y eintragen für X^Y")
      StringGadget(3, 50, 50, 300, 20, "2", #PB_String_Numeric) 
      StringGadget(4, 400, 50, 300, 20, "1", #PB_String_Numeric)
      ButtonGadget(5, 300, 560, 200, 30, "Berechnung starten")
      SetActiveGadget(2)

      For i = 1 To 5         ;das StringGadget einfach mal drin lassen
        SetGadgetFont(i, FontID(0))
      Next

      Repeat 
        Repeat
          Event = WaitWindowEvent()
          If Event = #PB_Event_CloseWindow
            CloseWindow(0)
            End 
          EndIf
        Until EventGadget() = 5

        Potenz = Val(GetGadgetText(4)) ;Exponent
        If (Len(GetGadgetText(4)) > 17) Or (Potenz > 1E16) ;willkürlich gewählte Werte, sichere Seite
          MessageRequester("Abbruch!", "Exponent zu groß! Bitte neue Eingabe!")
         Else
          Break
        EndIf     
      ForEver

      Basis$ = GetGadgetText(3)

      Basis = 0
      If Len(Basis$) < 2
        Basis = Len(Basis$)
      EndIf

      For i = 1 To 5
        FreeGadget(i)
      Next

      Info0$ = Basis$ + "^" + Str(Potenz) + ": "
      SetWindowTitle(0, Info0$)

      If Len(Basis$) > 0 Or Basis > 0
        If Potenz > 0
          Result$ = Basis$
          Calc()
         Else
          Result$ = "1"
        EndIf
       Else
        Result$ = "0"        ;0^Y
      EndIf
;=========================================================================
    Case 4                   ;Mersenne
      Info0$ = "Berechnung von Mersenne "
      SetWindowTitle(0, Info0$ + "Mxx")
      ;Werte für Mersenne M25-M49, die Strings sind nur zur Info
      M25.q = 21701    : M25$ = "(2^21701)-1"    : Z25$ = "44867916611904333479...57410828353511882751"
      M26.q = 23209    : M26$ = "(2^23209)-1"    : Z26$ = "40287411577898877818...36743355523779264511"
      M27.q = 44497    : M27$ = "(2^44497)-1"    : Z27$ = "85450982430363380319...44867686961011228671"
      M28.q = 86243    : M28$ = "(2^86243)-1"    : Z28$ = "53692799550275632152...99857021709433438207"
      M29.q = 110503   : M29$ = "(2^110503)-1"   : Z29$ = "52192831334175505976...69951621083465515007"
      M30.q = 132049   : M30$ = "(2^132049)-1"   : Z30$ = "51274027626932072381...52138578455730061311"
      M31.q = 216091   : M31$ = "(2^216091)-1"   : Z31$ = "74609310306466134368...91336204103815528447"
      M32.q = 756839   : M32$ = "(2^756839)-1"   : Z32$ = "17413590682008709732...02603793328544677887"
      M33.q = 859433   : M33$ = "(2^859433)-1"   : Z33$ = "12949812560420764966...02414267243500142591"
      M34.q = 1257787  : M34$ = "(2^1257787)-1"  : Z34$ = "41224577362142867472...31257188976089366527"
      M35.q = 1398269  : M35$ = "(2^1398269)-1"  : Z35$ = "81471756441257307514...85532025868451315711"
      M36.q = 2976221  : M36$ = "(2^2976221)-1"  : Z36$ = "62334007624857864988...76506256743729201151"
      M37.q = 3021377  : M37$ = "(2^3021377)-1"  : Z37$ = "12741168303009336743...25422631973024694271"
      M38.q = 6972593  : M38$ = "(2^6972593)-1"  : Z38$ = "43707574412708137883...35366526142924193791"
      M39.q = 13466917 : M39$ = "(2^13466917)-1" : Z39$ = "92494773800670132224...30073855470256259071"
      M40.q = 20996011 : M40$ = "(2^20996011)-1" : Z40$ = "12597689545033010502...94714065762855682047"
      M41.q = 24036583 : M41$ = "(2^24036583)-1" : Z41$ = "29941042940415717208...67436921882733969407"
      M42.q = 25964951 : M42$ = "(2^25964951)-1" : Z42$ = "12216463006127794810...98933257280577077247"
      M43.q = 30402457 : M43$ = "(2^30402457)-1" : Z43$ = "31541647561884608093...11134297411652943871"
      M44.q = 32582657 : M44$ = "(2^32582657)-1" : Z44$ = "12457502601536945540...11752880154053967871"
      M45.q = 37156667 : M45$ = "(2^37156667)-1" : Z45$ = "20225440689097733553...21340265022308220927"
      M46.q = 42643801 : M46$ = "(2^42643801)-1" : Z46$ = "16987351645274162247...84101954765562314751"
      M47.q = 43112609 : M47$ = "(2^43112609)-1" : Z47$ = "31647026933025592314...80022181166697152511"
      M48.q = 57885161 : M48$ = "(2^57885161)-1" : Z48$ = "58188726623224644217...46141988071724285951"
      M49.q = 74207281 : M49$ = "(2^74207281)-1" : Z49$ = "30037641808460618205...87010073391086436351"
      M50.q = 77232917 : M50$ = "(2^77232917)-1" : Z50$ = "46733318335923109998...82730618069762179071"

      TextGadget(1, 20, 0, 300, 15, "Wähle eine Mersenne-Zahl:")

      OptionGadget(2, 20, 20, 780, 15, "Mersenne M25, " + M25$ + ", 6533 Dezimal-Ziffern, " + Z25$)
      OptionGadget(3, 20, 40, 780, 15, "Mersenne M26, " + M26$ + ", 6987 Dezimal-Ziffern, " + Z26$)
      OptionGadget(4, 20, 60, 780, 15, "Mersenne M27, " + M27$ + ", 13395 Dezimal-Ziffern, " + Z27$)
      OptionGadget(5, 20, 80, 780, 15, "Mersenne M28, " + M28$ + ", 25962 Dezimal-Ziffern, " + Z28$)
      OptionGadget(6, 20, 100, 780, 15, "Mersenne M29, " + M29$ + ", 33265 Dezimal-Ziffern, " + Z29$)
      OptionGadget(7, 20, 120, 780, 15, "Mersenne M30, " + M30$ + ", 39751 Dezimal-Ziffern, " + Z30$)
      OptionGadget(8, 20, 140, 780, 15, "Mersenne M31, " + M31$ + ", 65050 Dezimal-Ziffern, " + Z31$)
      OptionGadget(9, 20, 160, 780, 15, "Mersenne M32, " + M32$ + ", 227832 Dezimal-Ziffern, " + Z32$)
      OptionGadget(10, 20, 180, 780, 15, "Mersenne M33, " + M33$ + ", 258716 Dezimal-Ziffern, " + Z33$)
      OptionGadget(11, 20, 200, 780, 15, "Mersenne M34, " + M34$ + ", 378632 Dezimal-Ziffern, " + Z34$)
      OptionGadget(12, 20, 220, 780, 15, "Mersenne M35, " + M35$ + ", 420921 Dezimal-Ziffern, " + Z35$)
      OptionGadget(13, 20, 240, 780, 15, "Mersenne M36, " + M36$ + ", 895832 Dezimal-Ziffern, " + Z36$)
      OptionGadget(14, 20, 260, 780, 15, "Mersenne M37, " + M37$ + ", 909526 Dezimal-Ziffern, " + Z37$)
      OptionGadget(15, 20, 280, 780, 15, "Mersenne M38, " + M38$ + ", 2098960 Dezimal-Ziffern, " + Z38$)
      OptionGadget(16, 20, 300, 780, 15, "Mersenne M39, " + M39$ + ", 4053946 Dezimal-Ziffern, " + Z39$)
      OptionGadget(17, 20, 320, 780, 15, "Mersenne M40, " + M40$ + ", 6320430 Dezimal-Ziffern, " + Z40$)
      OptionGadget(18, 20, 340, 780, 15, "Mersenne M41, " + M41$ + ", 7235733 Dezimal-Ziffern, " + Z41$)
      OptionGadget(19, 20, 360, 780, 15, "Mersenne M42, " + M42$ + ", 7816230 Dezimal-Ziffern, " + Z42$)
      OptionGadget(20, 20, 380, 780, 15, "Mersenne M43, " + M43$ + ", 9152052 Dezimal-Ziffern, " + Z43$)
      OptionGadget(21, 20, 400, 780, 15, "Mersenne M44, " + M44$ + ", 9808358 Dezimal-Ziffern, " + Z44$)
      OptionGadget(22, 20, 420, 780, 15, "Mersenne M45, " + M45$ + ", 11185272 Dezimal-Ziffern, " + Z45$)
      OptionGadget(23, 20, 440, 780, 15, "Mersenne M46, " + M46$ + ", 12837064 Dezimal-Ziffern, " + Z46$)
      OptionGadget(24, 20, 460, 780, 15, "Mersenne M47, " + M47$ + ", 12978189 Dezimal-Ziffern, " + Z47$)
      OptionGadget(25, 20, 480, 780, 15, "Mersenne M48, " + M48$ + ", 17425170 Dezimal-Ziffern, " + Z48$)
      OptionGadget(26, 20, 500, 780, 15, "Mersenne M49, " + M49$ + ", 22338618 Dezimal-Ziffern, " + Z49$)
      OptionGadget(27, 20, 520, 780, 15, "Mersenne M50, " + M50$ + ", 23249425 Dezimal-Ziffern, " + Z50$)

      SetGadgetState(2, 1)

      ButtonGadget(28, 300, 560, 200, 30, "Berechnung starten")

      For i = 1 To 28
        SetGadgetFont(i, FontID(0))
      Next

      Basis$ = "2"
      Result$ = "2"
      Potenz = M25 : M$ = "M25=" + M25$ : Z$ = Z25$   ;Voreinstellung Mersenne, falls nichts ausgewählt wird

      Repeat
        Event = WaitWindowEvent()
        If Event = #PB_Event_CloseWindow
          CloseWindow(0)
          End
        EndIf
        If Event = #PB_Event_Gadget
          Select EventGadget()
            Case 2
              Potenz = M25 : M$ = "M25=" + M25$ : Z$ = Z25$
            Case 3
              Potenz = M26 : M$ = "M26=" + M26$ : Z$ = Z26$
            Case 4
              Potenz = M27 : M$ = "M27=" + M27$ : Z$ = Z27$
            Case 5
              Potenz = M28 : M$ = "M28=" + M28$ : Z$ = Z28$
            Case 6
              Potenz = M29 : M$ = "M29=" + M29$ : Z$ = Z29$
            Case 7
              Potenz = M30 : M$ = "M30=" + M30$ : Z$ = Z30$
            Case 8
              Potenz = M31 : M$ = "M31=" + M31$ : Z$ = Z31$
            Case 9
              Potenz = M32 : M$ = "M32=" + M32$ : Z$ = Z32$
            Case 10
              Potenz = M33 : M$ = "M33=" + M33$ : Z$ = Z33$
            Case 11
              Potenz = M34 : M$ = "M34=" + M34$ : Z$ = Z34$
            Case 12
              Potenz = M35 : M$ = "M35=" + M35$ : Z$ = Z35$
            Case 13
              Potenz = M36 : M$ = "M36=" + M36$ : Z$ = Z36$
            Case 14
              Potenz = M37 : M$ = "M37=" + M37$ : Z$ = Z37$
            Case 15
              Potenz = M38 : M$ = "M38=" + M38$ : Z$ = Z38$
            Case 16
              Potenz = M39 : M$ = "M39=" + M39$ : Z$ = Z39$
            Case 17
              Potenz = M40 : M$ = "M40=" + M40$ : Z$ = Z40$
            Case 18
              Potenz = M41 : M$ = "M41=" + M41$ : Z$ = Z41$
            Case 19
              Potenz = M42 : M$ = "M42=" + M42$ : Z$ = Z42$
            Case 20
              Potenz = M43 : M$ = "M43=" + M43$ : Z$ = Z43$
            Case 21
              Potenz = M44 : M$ = "M44=" + M44$ : Z$ = Z44$
            Case 22
              Potenz = M45 : M$ = "M45=" + M45$ : Z$ = Z45$
            Case 23
              Potenz = M46 : M$ = "M46=" + M46$ : Z$ = Z46$
            Case 24
              Potenz = M47 : M$ = "M47=" + M47$ : Z$ = Z47$
            Case 25
              Potenz = M48 : M$ = "M48=" + M48$ : Z$ = Z48$
            Case 26
              Potenz = M49 : M$ = "M49=" + M49$ : Z$ = Z49$
            Case 27
              Potenz = M50 : M$ = "M50=" + M50$ : Z$ = Z50$
            Case 28
              Info0$ + M$ + ": "
              SetWindowTitle(0, Info0$)
              Break          
          EndSelect
        EndIf
      ForEver

      For i = 1 To 28
        FreeGadget(i)
      Next

      Calc()

      ResLen.q = StringByteLength(Result$)
      PokeB(@Result$ + ResLen - Asc_Uni, PeekB(@Result$ + ResLen - Asc_Uni) - 1)    ;-1 für Mersenne
;=========================================================================
    Case 5                   ;Fakultät
      Info0$ = "Fakultät: "
      SetWindowTitle(0, Info0$)

      TextGadget(1, 50, 10, 300, 20, "Fakultät:")
      StringGadget(2, 50, 50, 300, 20, "2", #PB_String_Numeric)  
      ButtonGadget(3, 300, 560, 200, 30, "Berechnung starten")
      SetActiveGadget(2)

      For i = 1 To 3         ;das StringGadget einfach mal drin lassen
        SetGadgetFont(i, FontID(0))
      Next

      Repeat
        Repeat
          Event = WaitWindowEvent()
          If Event = #PB_Event_CloseWindow
            CloseWindow(0)
            End 
          EndIf
        Until EventGadget() = 3

        Fak = Val(GetGadgetText(2))    ;hier reicht der Wertebereich von Val() wohl aus ;-)   
        ;besser doch Test
        If (Len(GetGadgetText(2)) > 7) Or (Fak > 1E6) ;willkürlich gewählte Werte, sichere Seite
          MessageRequester("Abbruch!", "Wert zu groß! Bitte neue Eingabe!")
         Else
          Break
        EndIf     
      ForEver

      Info0$ + Str(Fak) + "! "
      SetWindowTitle(0, Info0$)      

      Result$ = "1"

      For i = 1 To 3
        FreeGadget(i)
      Next

      ListViewGadget(2, 10, 10, 250, 580)   ;Info-Fenster

      AddWindowTimer(0, 0, 5000)            ;für Anzeige im Info-Fenster
      BindEvent(#PB_Event_Timer, @FakInfo())

      GadgetZeile = 0
      AddGadgetItem (2, -1, "Berechnung gestartet... ")   
      SetGadgetState(2, GadgetZeile)

      ZeitA = ElapsedMilliseconds()
      For FakA = 1 To Fak
        Faktor1$ = Str(FakA)
        WindowEvent()                       ;nur für Anzeige im Info-Fenster
        Result$ = FFT_Mul(@Faktor1$, @Result$)
      Next
      ZeitE = ElapsedMilliseconds() - ZeitA

      RemoveWindowTimer(0, 0)
      UnbindEvent(#PB_Event_Timer, @FakInfo())   ;nötig!
      FreeGadget(2)
;=========================================================================
    Case 6                   ;Neuner-Test
      Info0$ = "Neuner-Test: "
      SetWindowTitle(0, Info0$)

      Test9 = 540000 * 4     ;540000*4 ist o.K. für Double

      Buffer9 = AllocateMemory(8 * Test9 * #Asc_Uni + 64)
      Buffer9A = Buffer9 + (64 - (Buffer9 & $3F))
      !mov rax,[v_Buffer9A]
      !mov rcx,[v_Test9]
      !cmp [v_Asc_Uni],1
     !je .ASCII
      !mov rdx,009000900090009h   ;Unicode Word
      !shl rcx,1
     !jmp @f  
     !.ASCII:
      !mov rdx,909090909090909h   ;ASCII Byte
     !@@:
      !mov [rax],rdx
      !add rax,8
      !dec rcx
     !jnz @b

      Info0$ + "Berechne " + Str(8 * Test9) + "*" + Str(8 * Test9) + " Neuner-Bytes (Digits): "
      SetWindowTitle(0, Info0$)
      ZeitA = ElapsedMilliseconds()
      Result$ = FFT_Mul(Buffer9A, Buffer9A)
      ZeitE = ElapsedMilliseconds() - ZeitA

      Info9$ = ""
      L = StringByteLength(Result$)
      Buffer1 = @Result$

      If PeekB(Buffer1 + L - Asc_Uni) <> $31     ;hier könnte auch die Konstante genommen werden
        Info9$ = "Fehler 1! "
      EndIf

      If PeekB(Buffer1 + (L / 2) - Asc_Uni) <> $38
        Info9$ + "Fehler 8! "
      EndIf

      For i = 0 To (L / 2) - 2 * Asc_Uni Step #Asc_Uni
        If PeekB(Buffer1 + i) <> $39
          Info9$ + "Fehler 9! "
          Break
        EndIf
      Next

      For i = L / 2 To L - 2 * Asc_Uni Step #Asc_Uni
        If PeekB(Buffer1 + i) <> $30
          Info9$ + "Fehler 0! "
          Break
        EndIf
      Next

      If Info9$ = ""
        Info9$ = "Alles o.K.! "
      EndIf

      FreeMemory(Buffer9)
      Info0$ + Info9$
      SetWindowTitle(0, Info0$)
;=========================================================================
  EndSelect
;=========================================================================
  EditorGadget(1, 10, 10, 780, 540)

  Info1$ = Info0$ + "Benötigte Zeit: " + Str(ZeitE) + "ms  " + Z$
  SetWindowTitle(0, Info1$)

  Digits$ = ""
  If Asc(Mid(Result$, 1, 1)) < 58 ;also Ziffer (hoffentlich :-) !) und nicht z.B.Fehlermeldung
    Digits$ = "Digits:" + Str(Len(Result$)) + #LFCR$
  EndIf
  SetGadgetText(1, Digits$ + Result$)

  ButtonGadget(2, 100, 560, 250, 30, "In die Zwischenablage kopieren")
  SetGadgetFont(2, FontID(0))
  ButtonGadget(3, 450, 560, 250, 30, "Zurück zur Auswahl")
  SetGadgetFont(3, FontID(0))

  Repeat
    Event = WaitWindowEvent()
    If Event = #PB_Event_CloseWindow
      Break 2
    EndIf
    If Event = #PB_Event_Gadget
      Select EventGadget()
        Case 2
          SetClipboardText(Result$)
        Case 3
          SetWindowTitle(0, FensterTitel$)
          Info0$ = "" : Info1$ = "" : Z$ = ""
          Break          
      EndSelect
    EndIf
  ForEver

  For i = 1 To 3
    FreeGadget(i)
  Next

Until WaitWindowEvent() = #PB_Event_CloseWindow

CloseWindow(0)
End

Procedure.s FFT_Mul(PFac1, PFac2)
  EnableASM 

  Protected.q BufferCosSin, BufferCosSinA, BufferDez, BufferDezA, BufferFX, BufferFXA, BufferFaktor1, BufferFaktor1A
  Protected.q BufferFaktor2, BufferFaktor2A, BufferProdukt, BufferProduktA, BufferReverse, BufferReverseA, LG, N
  Protected.s sAusgabe

  ;Register-, Status- und Control- Speicher-Reservierung
  BufferFX = AllocateMemory(512 + 64)
  If BufferFX
    BufferFXA = BufferFX + (64 - (BufferFX & $3F))    ;64-er Alignment, klotzen, nicht kleckern!
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferFXA
  !pop [BufferFXA]
;!Align 16
  ;Save Registers
  !mov rax,[BufferFXA]
  !fxsave64 [rax]            ;sichert erstmal komplett FPU (mit MXCSR) und XMM (0-15). Lass ich mal so, obwohl FNSAVE auch reichen würde. Dann aber so keine Register-Sicherung
  !mov [rax+464],rdi         ;sichert callee-save registers, ohne RBX und RBP!
  !mov [rax+472],rsi
  !mov [rax+480],r12
  !mov [rax+488],r13
  !mov [rax+496],r14
  !mov [rax+504],r15

  ;übergebende Parameter
  !mov r8,[p.v_PFac1]
  !mov [PF1],r8
  !mov r9,[p.v_PFac2]
  !mov [PF2],r9

  !mov rax, -16
  !mov rdx, -16
  !pxor xmm1,xmm1            ;xmm1 null setzen, ist praktisch das Zero-Byte bzw. 2 Zero-Bytes für Unicode

  ;Länge der Strings ermitteln
  CompilerIf #PB_Compiler_Unicode
   !@@:
    !add rax,16
    !vpcmpistri xmm1,dqword[r8+rax],00001001b    ;Bit0 = 1 und Bit1 = 0, also String-Chars werden als unsigned Words interpretiert, Bit2 = 0 und Bit3 = 1, also Test auf equal each
   !jnz @b
    !shr rax,1                                   ;Pointer in String muss hier halbiert werden
    !add rax,rcx

    !pxor xmm1,xmm1                              ;xmm1 null setzen, ist praktisch das Zero-Byte bzw. 2 Zero-Bytes für Unicode
   !@@:
    !add rdx,16
    !vpcmpistri xmm1,dqword[r9+rdx],00001001b    ;Bit0 = 1 und Bit1 = 0, also String-Chars werden als unsigned Words interpretiert, Bit2 = 0 und Bit3 = 1, also Test auf equal each
   !jnz @b
    !shr rdx,1                                   ;Pointer in String muss hier halbiert werden
    !add rdx,rcx
   CompilerElse
   !@@:
    !add rax,16
    !vpcmpistri xmm1,dqword[r8+rax],00001000b    ;Bit0 = 0 und Bit1 = 0, also String-Chars werden als unsigned Bytes interpretiert, Bit2 = 0 und Bit3 = 1, also Test auf equal each
   !jnz @b
    !add rax,rcx

    !pxor xmm1,xmm1                              ;xmm1 null setzen, ist praktisch das Zero-Byte bzw. 2 Zero-Bytes für Unicode
   !@@:
    !add rdx,16
    !vpcmpistri xmm1,dqword[r9+rdx],00001000b    ;Bit0 = 0 und Bit1 = 0, also String-Chars werden als unsigned Bytes interpretiert, Bit2 = 0 und Bit3 = 1, also Test auf equal each
   !jnz @b
    !add rdx,rcx
  CompilerEndIf

  !mov [L1],rax
  !mov [L2],rdx
  !cmp rax,rdx               ;L1 soll der größere Faktor sein
 !jae @f
  !mov [L1],rdx
  !mov [L2],rax
  !mov [PF1],r9
  !mov [PF2],r8
 !@@:

  !add rax,rdx
  !push rax
  POP LG                     ;Gesamtlänge L1+L2 für "Verstringung" und Speicher-Management  

  ;LenF = L1
  !mov rdx,[L1]              ;L1 ist größter Faktor
  !mov [LenF],rdx            ;LenF erstmal setzen, könnte ja schon 2-er Potenz sein 
  !bsr rcx,[L1]              ;höchstes gesetztes Bit
  !bsf rdx,[L2]              ;Test auf Länge Null
 !jz @f
  !mov [Exponent],rcx
  !bsf rdx,[L1]              ;niedrigstes gesetztes Bit
 !jmp .Faktoren_OK           ;local label
 !@@:
  sAusgabe = "Fehler! Mindestens einer der Faktoren hat Länge Null! " ;wegen Kompatibilität zu "normaler" Multiplikation gelassen
 ProcedureReturn sAusgabe

 !.Faktoren_OK:
  !cmp rcx,rdx
 !je @f                      ;also nur 1 Bit (identisch) gesetzt, ist 2-er Potenz
  !inc rcx
  !mov rax,1
  !shl rax,cl
  !mov [LenF],rax            ;LenF ist jetzt 2-er Potenz
  !add [Exponent],1
 !@@:
  !add [Exponent],1

  !mov rcx,5                 ;Produkt-Länge muss mindestens 32 sein, sonst an anderen Stellen Klimmzüge ("Verstringung")
  !sub rcx,[Exponent]
 !js @f                      ;Exponent ist > 5 (Hauptfall)
 !jz @f                      ;Exponent ist = 5
  !shl [LenF],cl
  !add [Exponent],rcx
 !@@:

  ;zur Sicherheit im MXCSR Bit 13 und 14 auf Null setzen (round to nearest)
  ;MXCSR[6]         : Denormals are zeros - 0
  ;MXCSR[7:12]      : Exception masks all 1's (all exceptions masked)
  ;MXCSR[13:14]   : Rounding  control - 0 (round to nearest)
  ;MXCSR[15]      : Flush to zero for masked underflow - 0 (off)
  ;R+	bit 14	Round Positive
  ;R-	bit 13	Round Negative
  ;RZ	bits 13 and 14	Round to Zero
  ;RN	bits 13 and 14 are 0	Round to Nearest
  ;load MXCSR
  !vstmxcsr [MXCSR]          ;load Control-and Status-Register in Memory
  !mov eax,[MXCSR] 
  !mov [MXCSR_Sicher],eax    ;trotz fxsave64 (s.o.) hier zur Sicherheit 
  ;MXCSR_Sicher = MXCSR
  !and [MXCSR],11111111111111111001111111111111b ;zur Sicherheit Bit 13 und 14 auf Null setzen (round to nearest)
  !vldmxcsr [MXCSR]          ;store MXCSR

  ;N = LenF << 1
  !mov rax,[LenF]
  !shl rax,1
  !mov [N],rax
  !push rax                  ;wird dann "PB-N"

  !vcvtsi2sd xmm1,xmm2,rax   ;Quad-N in Double-ND konvertieren für spätere Division
  !vmovsd [ND],xmm1

  ;Pointer_Produkt.q = N - 1
  !dec rax
  !mov [Pointer_Produkt],rax

  POP N                      ;wird u.a.für AllocateMemory benötigt!

  ;Speicher-Reservierung
  BufferCosSin = AllocateMemory((2 * N * 8) + 64)
  If BufferCosSin
    BufferCosSinA = BufferCosSin + (64 - (BufferCosSin & $3F))  ;64-er Alignment, klotzen, nicht kleckern!
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferCosSinA
  !pop [BufferCosSinA]

  BufferFaktor1 = AllocateMemory((2 * N * 8) + 64)
  If BufferFaktor1
    BufferFaktor1A = BufferFaktor1 + (64 - (BufferFaktor1 & $3F))    ;64-er Alignment
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferFaktor1A
  !pop [BufferFaktor1A]

  BufferFaktor2 = AllocateMemory((2 * N * 8) + 64)
  If BufferFaktor2
    BufferFaktor2A = BufferFaktor2 + (64 - (BufferFaktor2 & $3F))    ;64-er Alignment
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferFaktor2A
  !pop [BufferFaktor2A]
  ;========================================================
  ;Faktoren bearbeiten
  !bsf rax,[N]               ;Da N hier nie Null sein kann tuts BSF 
  !mov [Exponent],rax

  !mov rdx,32
  !sub rdx,rax
  !vmovq xmm9,rdx            ;32 - Exponent

  !mov rsi,[N]
  !mov r8,[PF1]
  !mov r10,[PF2]
  !mov r11,[BufferCosSinA]   ;Temp-Memory

  !mov r9,r11
  !add r9,rsi
  !add r9,rsi

  !mov r14,[L1]
  !mov rcx,[v_Asc_Uni]
  !dec rcx
  !shl r14,cl
  !sub r14,[v_Asc_Uni]

  !mov r15,[L2]
  !shl r15,cl
  !sub r15,[v_Asc_Uni]

  !mov rcx,2
  !sub rcx,[v_Asc_Uni]
  !shr rsi,cl   
  !sub rsi,[v_Asc_Uni]

;------------------
  !mov rcx,3
 !@@:
  !mov [r11+rcx*4],ecx
  !dec rcx
 !jns @b 
  !vmovdqa xmm0,[r11]        ;Start-Values 0,1,2,3

  !vbroadcastss xmm1,[V1]    ;AVX  Missbrauch von Float, was solls... kostet aber Speed
  !vbroadcastss xmm2,[V2]
  !vbroadcastss xmm3,[V3]
  !vbroadcastss xmm4,[V4]
  !vbroadcastss xmm11,[VVier]

  !mov rdi,[BufferFaktor1A]  ;Temp-Memory
  !xor r12,r12
  ;------------------------------------

 !.Loop_Rev_Fak:
  !vpsrld xmm5,xmm0,1        ;>> 1
  !vpand xmm6,xmm5,xmm1      ;& V1 
  !vpand xmm7,xmm0,xmm1      ;& V1 
  !vpslld xmm8,xmm7,1        ;<< 1
  !vpor xmm10,xmm6,xmm8      ;| 

  !vpsrld xmm5,xmm10,2       ;>> 2
  !vpand xmm6,xmm5,xmm2      ;& V2 
  !vpand xmm7,xmm10,xmm2     ;& V2 
  !vpslld xmm8,xmm7,2        ;<< 2
  !vpor xmm10,xmm6,xmm8      ;| 

  !vpsrld xmm5,xmm10,4       ;>> 4
  !vpand xmm6,xmm5,xmm3      ;& V3 
  !vpand xmm7,xmm10,xmm3     ;& V3 
  !vpslld xmm8,xmm7,4        ;<< 4
  !vpor xmm10,xmm6,xmm8      ;| 

  !vpsrld xmm5,xmm10,8       ;>> 8
  !vpand xmm6,xmm5,xmm4      ;& V4 
  !vpand xmm7,xmm10,xmm4     ;& V4 
  !vpslld xmm8,xmm7,8        ;<< 8
  !vpor xmm10,xmm6,xmm8      ;|

  !vpsrld xmm5,xmm10,16      ;>> 16
  !vpslld xmm8,xmm10,16      ;<< 16
  !vpor xmm10,xmm5,xmm8      ;|

  !vpsrld xmm5,xmm10,xmm9    ;>> (32 - Exponent)
  !vmovdqa [rdi],xmm5

  !xor r13,r13
 !.Next_Byte4:
  !mov r12d,[rdi+r13*4]      ;4 reverse Positionen  
  !xor rax,rax
  !mov rdx,r14
  !mov rcx,[v_Asc_Uni]
  !dec rcx
  !shr rdx,cl
  !cmp [L1],rdx
 !jb @f 

  !movzx eax,byte[r8+r14]    ;next Byte Faktor1$
  !and eax,0fh               ;for String (=sub 30h)
 !@@:
  !mov [r9+r12*2],eax

  !xor rax,rax
  !mov rdx,r15
  !mov rcx,[v_Asc_Uni]
  !dec rcx
  !shr rdx,cl
  !cmp [L2],rdx
 !jb @f 
  !movzx eax,byte[r10+r15]   ;next Byte Faktor2$
  !and eax,0fh
 !@@:
  !mov [r11+r12*2],eax

  !sub rsi,[v_Asc_Uni]
  !sub r14,[v_Asc_Uni]
  !sub r15,[v_Asc_Uni]

  !inc r13
  !cmp r13,4
 !jne .Next_Byte4

  !vpaddd xmm0,xmm0,xmm11

  !or rsi,rsi
 !jns .Loop_Rev_Fak
  ;========================================================
  ;Integer in Double konvertieren, verdoppeln und abspeichern
  !xor rdx,rdx
  !xor r8,r8
  !mov rcx,[N]
  !shr rcx,3
  !mov r10,[BufferFaktor1A]
  !mov r11,[BufferCosSinA]
  !mov r12,[BufferFaktor2A]

  !mov r9,r11
  !mov rsi,[N]
  !add r9,rsi
  !add r9,rsi

 !@@:
  ;Faktor1
  !vcvtdq2pd ymm0,[r9+r8]    ;4 Ziffern (als DWords) von Faktor1 in 4 Doubles konvertieren 
  !vmovsd [r10+rdx],xmm0
  !vmovsd [r10+rdx+16],xmm0  ;+8 ist der I-Anteil, bleibt hier Null!
  !vunpckhpd ymm1,ymm0,ymm0  ;AVX 
  !vmovsd [r10+rdx+32],xmm1
  !vmovsd [r10+rdx+48],xmm1
  !vperm2f128 ymm2,ymm0,ymm0,1
  !vmovsd [r10+rdx+64],xmm2
  !vmovsd [r10+rdx+80],xmm2
  !vunpckhpd ymm1,ymm2,ymm0
  !vmovsd [r10+rdx+96],xmm1
  !vmovsd [r10+rdx+112],xmm1
  ;Faktor2
  !vcvtdq2pd ymm0,[r11+r8]
  !vmovsd [r12+rdx],xmm0
  !vmovsd [r12+rdx+16],xmm0
  !vunpckhpd ymm1,ymm0,ymm0
  !vmovsd [r12+rdx+32],xmm1
  !vmovsd [r12+rdx+48],xmm1
  !vperm2f128 ymm2,ymm0,ymm0,1
  !vmovsd [r12+rdx+64],xmm2
  !vmovsd [r12+rdx+80],xmm2
  !vunpckhpd ymm1,ymm2,ymm0
  !vmovsd [r12+rdx+96],xmm1
  !vmovsd [r12+rdx+112],xmm1

  !add r8,16
  !add rdx,128
  !sub rcx,1
 !jnz @b
  ;========================================================
  ;Winkel erst hier, weil BufferCosSinA oben "missbraucht" wird
  !vcvtsi2sd xmm0,xmm0,qword[N]
  !vmovsd xmm1,[Pi2]
  !vdivsd xmm7,xmm1,xmm0     ;Schrittweite Konstante
  !vmovq xmm6,qword[Signum]  ;Konstante

  !mov rcx,[N]
  !shr rcx,3
  !xor r8,r8                 ;Schritt-Nr. Konstante
  !mov r9,[N]
  !shr r9,2                  ;N/4 für Adressierung unten
  !mov r10,[N]
  !shr r10,1                 ;N/2

  !mov r11,[BufferCosSinA]
  !inc rcx
 !@@:
  !vcvtsi2sd xmm15,xmm14,r8
  !vmulsd xmm0,xmm7,xmm15

  !vmovddup xmm0,xmm0                  ;XMM0: 1=x^1 2=x^1  duplicate bits 0-63 in bits 64-127
  !vinsertf128 ymm1,ymm0,xmm0,1b       ;YMM1: 1=x^1 2=x^1 3=x^1 4=x^1  YMM1(0-127)=XMM0, YMM1(128-255)=XMM0

  !vmulpd ymm2,ymm1,ymm1     ;YMM2: 1=x^2 2=x^2 3=x^2 4=x^2
  !vmulsd xmm3,xmm2,xmm2     ;YMM3: 1=x^4 2=x^2 3=0   4=0
  !vmulpd ymm4,ymm2,ymm2     ;YMM4: 1=x^4 2=x^4 3=x^4 4=x^4
  !vmulpd ymm2,ymm4,ymm4     ;YMM2: 1=x^8 2=x^8 3=x^8 4=x^8  for the next 4 terms

  !vmovdqa xmm8,xmm3         ;save for Cos
  !vmulpd xmm9,xmm8,xmm4 
  !vperm2f128 ymm10,ymm9,ymm8,100000b

  !vmulpd xmm14,xmm3,xmm1              ;YMM14: 1=x^5 2=x^3 3=0 4=0
  !vmulpd xmm1,xmm14,xmm4              ;YMM1: 1=x^9 2=x^7 3=0 4=0
  !vperm2f128 ymm5,ymm1,ymm14,100000b  ;YMM5: 1=x^9 2=x^7 3=x^5 4=x^3  YMM5(0-127)=YMM1(0-127), YMM5(128-255)=YMM14(0-127)  

  !vmulpd ymm1,ymm5,yword[RezFakSin]   ;multiply the 4 values in YMM5 with RezFak -1/3! ... 1/9!, result in YMM1
  !vmulpd ymm11,ymm10,yword[RezFakCos]
  ;next 4 terms, without loop 
  !vmulpd ymm14,ymm5,ymm2              ;YMM14: 1=x^17 2=x^15 3=x^13 4=x^11
  !vmulpd ymm10,ymm10,ymm2

  !vmulpd ymm15,ymm14,yword[RezFakSin+32]
  !vaddpd ymm1,ymm15,ymm1
  !vmulpd ymm15,ymm10,yword[RezFakCos+32]
  !vaddpd ymm11,ymm15,ymm11
  ;next 4 terms
  !vmulpd ymm5,ymm14,ymm2              ;YMM5: 1=x^25 2=x^23 3=x^21 4=x^19
  !vmulpd ymm10,ymm10,ymm2
  !vmulpd ymm15,ymm5,yword[RezFakSin+64]
  !vaddpd ymm1,ymm15,ymm1
  !vmulpd ymm15,ymm10,yword[RezFakCos+64]
  !vaddpd ymm11,ymm15,ymm11
  ;result
  !vhaddpd ymm2,ymm1,ymm1              ;YMM2: 1=1+2 of YMM1, 3=3+4 of YMM1 
  !vhaddpd ymm12,ymm11,ymm11

  !vextractf128 xmm1,ymm2,1b           ;XMM1: 3+4 of YMM2
  !vextractf128 xmm11,ymm12,1b

  !vaddsd xmm3,xmm2,xmm1               ;XMM3: 1=sum of iterations 
  !vaddsd xmm13,xmm12,xmm11

  !vaddsd xmm12,xmm0,xmm3              ;Sine   XMM12: 1=sum of iterations plus start-value (1.term=x)
  !vaddsd xmm14,xmm13,[Eins]           ;Cosine
  ;------------------------------------
  ;store values 0-180°
  !mov rax,r8
  !shl rax,1
  !vmovsd qword[r11+rax*8+8],xmm12     ;Sin 0-45°
  !vorpd xmm1,xmm12,xmm6               ;neg Sin 
  !vmovsd qword[r11+rax*8],xmm14       ;Cos 0-45°

  !mov rax,r9
  !sub rax,r8
  !shl rax,1
  !vmovsd qword[r11+rax*8],xmm12       ;Cos 45-90°
  !mov rdx,r9
  !add rdx,r8
  !shl rdx,1
  !vmovsd qword[r11+rdx*8],xmm1        ;Cos 90-135°
  !mov r13,r10
  !sub r13,r8
  !shl r13,1
  !vorpd xmm2,xmm14,xmm6
  !vmovsd qword[r11+r13*8],xmm2        ;Cos 135-180°
  !vmovsd qword[r11+rax*8+8],xmm14     ;Sin 45-90° 
  !vmovsd qword[r11+rdx*8+8],xmm14     ;Sin 90-135°
  !vmovsd qword[r11+r13*8+8],xmm12     ;Sin 135-180°

  !inc r8
  !dec rcx
 !jnz @b
  ;========================================================
  ;jetzt die beiden Faktoren transformieren
  !mov rsi,[N]
  !mov r8,[BufferCosSinA]

  !mov rcx,1                           ;Exponent_FFT
  !mov r12,2                           ;Pointer0
  !mov r9,2                            ;Pointer2

  !mov r10,[BufferFaktor1A]
  !mov rdi,[BufferFaktor2A]

 !.OuterLoopF:
  !inc rcx                             ;Exponent_FFT
  !mov rax,rsi                         ;N
  !shr rax,cl
 !jz .EndF
  !mov r14,rax

   !.InnerLoopF:
    !xor r11,r11                       ;Pointer1
    !mov r15,r12                       ;Pointer0
  ;------------------------------------
     !@@:
      !mov rax,r14
      !mul r11                         ;Pointer1
      !shl rax,3                       ;EW*8, Pointer in CosSin

      !mov rdx,r9                      ;Pointer2
      !mov r13,rdx
      !shl rdx,4                       ;Pointer2*8
      !sub r13,r12                     ;[v_Pointer0]
      !shl r13,4                       ;(Pointer2 - Pointer0) * 8 

      !vperm2f128 ymm0,ymm9,yword[r8+rax],100010b     ;2x CosSinA in ymm0, ymm9 is Dummy

      !vmovapd xmm7,dqword[rdi+rdx]    ;BufferFaktor2A complete (R and I)
      !vperm2f128 ymm2,ymm7,yword[r10+rdx],10b

      !vmulpd ymm4,ymm2,ymm0
      !vshufpd ymm1,ymm0,ymm0,101b     ;change Sin and Cos
      !vhsubpd ymm3,ymm4,ymm4          ;WR
      
      !vmulpd ymm4,ymm2,ymm1
      !vhaddpd ymm5,ymm4,ymm4          ;WI
      !vshufpd ymm15,ymm3,ymm5,0b      ;WI in ymm15-high
      !vmovapd xmm1,dqword[r10+r13]
      !vmovapd xmm7,dqword[rdi+r13]
      !vperm2f128 ymm2,ymm7,ymm1,10b
      !vaddpd ymm14,ymm15,ymm2
      !vsubpd ymm13,ymm2,ymm15

      !vmovapd [r10+r13],xmm14

      !vperm2f128 ymm2,ymm7,ymm14,11b
      !vmovapd [rdi+r13],xmm2

      !vmovapd [r10+rdx],xmm13

      !vperm2f128 ymm2,ymm7,ymm13,11b
      !vmovapd [rdi+rdx],xmm2

      !add r11,2                       ;Pointer1
      !inc r9                          ;Pointer2
      !dec r15
     !jnz @b

    !add r9,r12                        ;Pointer0
    !cmp r9,rsi                        ;N
   !jb .InnerLoopF

  !shl r12,1
  !mov r9,r12

 !jmp .OuterLoopF

 !.EndF:
  ;========================================================
  ;Und jetzt die beiden transformierten Vektoren miteinander multiplizieren
  !mov r8,[BufferFaktor1A]        ;nimmt auch das Produkt auf
  !mov r9,[BufferFaktor2A]
  !xor r11,r11                    ;Pointer
  !mov rcx,[N]
  !shr rcx,1
  !vmovdqu ymm5,[Null1]           ;Null,Signum,Null,Signum
 !@@:
  !vmovapd ymm0,[r8+r11]
  !vmovapd ymm3,[r9+r11]
  !vmulpd ymm1,ymm0,ymm3          ;ymm1:R*R,I*I,R*R,I*I
  !vshufpd ymm2,ymm3,ymm3,101b    ;change R and I
  !vxorpd ymm4,ymm1,ymm5          ;Vorzeichen ändern, damit nicht sub sondern add für beide Terme verwendet werden kann
  
  !vmulpd ymm1,ymm0,ymm2          ;R*I, R*I
  !vhaddpd ymm6,ymm4,ymm1
  !vmovapd [r8+r11],ymm6

  !add r11,32
  !dec rcx
 !jnz @b
  ;========================================================
  ;Rück-Transformation
  ;Produkt bearbeiten
  BufferReverse = AllocateMemory(32 + 64)
  If BufferReverse
    BufferReverseA = BufferReverse + (64 - (BufferReverse & $3F))    ;64-er Alignment
   Else 
    sAusgabe = "-1"
   ProcedureReturn sAusgabe
  EndIf
  PUSH BufferReverseA
  !pop [BufferReverseA]

  !mov rsi,[N]               ;für Abbruch Schleife Loop_Rev_Prod; letzte Pos.steht ja fest (wie erste) 
  !shl rsi,1
  !sub rsi,2
  !mov r8,[BufferFaktor1A]   ;Ergebnis von Multiplikation oben
  !mov r9,[BufferFaktor2A]   ;hier neue Reihenfolge rein

  !mov rdx,32
  !sub rdx,[Exponent]
  !vmovq xmm9,rdx

  ;----------------------------------------------

  !mov rcx,3
 !@@:
  !mov [r9+rcx*4],ecx
  !dec rcx
 !jns @b 
  !vmovdqa xmm0,[r9]         ;Start-Values 0,1,2,3

  !vbroadcastss xmm1,[V1]    ;AVX  Missbrauch von Float, was solls... kostet aber Speed
  !vbroadcastss xmm2,[V2]
  !vbroadcastss xmm3,[V3]
  !vbroadcastss xmm4,[V4]
  !vbroadcastss xmm11,[VVier]

  !mov rdi,[BufferReverseA]
  !xor r12,r12
  ;------------------------------------

 !.Loop_Rev_Prod:
  !vpsrld xmm5,xmm0,1        ;>> 1
  !vpand xmm6,xmm5,xmm1      ;& V1 
  !vpand xmm7,xmm0,xmm1      ;& V1 
  !vpslld xmm8,xmm7,1        ;<< 1
  !vpor xmm10,xmm6,xmm8      ;| 

  !vpsrld xmm5,xmm10,2       ;>> 2
  !vpand xmm6,xmm5,xmm2      ;& V2 
  !vpand xmm7,xmm10,xmm2     ;& V2 
  !vpslld xmm8,xmm7,2        ;<< 2
  !vpor xmm10,xmm6,xmm8      ;| 

  !vpsrld xmm5,xmm10,4       ;>> 4
  !vpand xmm6,xmm5,xmm3      ;& V3 
  !vpand xmm7,xmm10,xmm3     ;& V3 
  !vpslld xmm8,xmm7,4        ;<< 4
  !vpor xmm10,xmm6,xmm8      ;| 

  !vpsrld xmm5,xmm10,8       ;>> 8
  !vpand xmm6,xmm5,xmm4      ;& V4 
  !vpand xmm7,xmm10,xmm4     ;& V4 
  !vpslld xmm8,xmm7,8        ;<< 8
  !vpor xmm10,xmm6,xmm8      ;|

  !vpsrld xmm5,xmm10,16      ;>> 16
  !vpslld xmm8,xmm10,16      ;<< 16
  !vpor xmm10,xmm5,xmm8      ;|

  !vpsrld xmm5,xmm10,xmm9    ;>> (32 - Exponent)
  !vmovdqa [rdi],xmm5

  !xor rcx,rcx
 !@@:
  !mov r12d,[rdi+rcx*4]
  !shl r12,1
  !vmovapd xmm12,[r8]        ;16 Bytes (R.d und I.d) 
  !vmovapd [r9+r12*8],xmm12
  !add r8,16
  !inc rcx
  !cmp rcx,4
 !jne @b

  !vpaddd xmm0,xmm0,xmm11

  !cmp r12,rsi               ;Position letzter Wert?
 !jne .Loop_Rev_Prod
  ;========================================================
  !mov r8,[BufferFaktor2A]
  !mov rcx,[N]
  !shr rcx,1
!@@:
  !vmovdqa xmm0,[r8]
  !vaddpd xmm1,xmm0,[r8+16]
  !vmovsd [r8],xmm1
  !vunpckhpd xmm2,xmm1,xmm1
  !vmovsd [r8+8],xmm2

  !vsubpd xmm1,xmm0,[r8+16]
  !vmovsd [r8+16],xmm1
  !vunpckhpd xmm2,xmm1,xmm1
  !vmovsd [r8+24],xmm2

  !add r8,32
  !dec rcx
 !jnz @b
  ;------------------------------------

  ;------------------------------------
  !mov rsi,[N]
  !mov r8,[BufferCosSinA]
  !mov rcx,1                           ;Exponent_FFT
  !mov r12,2                           ;Pointer0 
  !mov r9,2                            ;Pointer2
  !mov r10,[BufferFaktor2A]

 !.OuterLoopP:
  !inc rcx                             ;Exponent_FFT
  !mov rax,rsi                         ;N
  !shr rax,cl
 !jz .EndP

  !mov r14,rax
  !shl r14,3                           ;EW*8
   !.InnerLoopP:
    !xor r11,r11                       ;Pointer1
    !mov r15,r12
  ;------------------------------------
     !@@:
      !mov rax,r14
      !mul r11                         ;Pointer1    oder imul rax,r11  ;dann wäre rdx frei

      !mov rdx,r9                      ;Pointer2
      !mov r13,rdx
      !shl rdx,4
      !sub r13,r12                     ;[v_Pointer0]
      !shl r13,4

      !vmovapd xmm7,dqword[r10+r13]

      !vmovapd xmm0,dqword[r10+rdx]    ;BufferFaktor2A
      !vmovapd xmm1,dqword[r8+rax]     ;CosSinA
      !vshufpd xmm2,xmm0,xmm0,1b       ;R und I vertauschen

      !vmulpd xmm3,xmm2,xmm1
      !vmulpd xmm4,xmm0,xmm1

      !vxorpd xmm5,xmm3,[Null1]        ;Signum, xmm3 negieren um vhaddpd verwenden zu können (sonst sub) wie xmm4

      !vhaddpd xmm6,xmm4,xmm5          ;WI in xmm15-high, WR in Low

      !vaddpd xmm14,xmm7,xmm6          ;dqword[r10+r13]
      !vsubpd xmm13,xmm7,xmm6

      !vmovapd [r10+r13],xmm14
      !vmovapd [r10+rdx],xmm13

      !add r11,2                       ;Pointer1
      !inc r9                          ;Pointer2

      !dec r15
     !jnz @b

;---------------------------------------
    !add r9,r12                        ;Pointer0
    !cmp r9,rsi                        ;N
   !jb .InnerLoopP

  !shl r12,1
  !mov r9,r12

 !jmp .OuterLoopP
 !.EndP:
  ;========================================================
  FreeMemory(BufferReverse)
  FreeMemory(BufferCosSin)
  FreeMemory(BufferFaktor1)

  BufferProdukt = AllocateMemory(N + 64)
  If BufferProdukt
    BufferProduktA = BufferProdukt + (64 - (BufferProdukt & $3F))    ;64-er Alignment
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferProduktA
  !pop [BufferProduktA]

  BufferDez = AllocateMemory(16 + 64)
  If BufferDez
    BufferDezA = BufferDez + (64 - (BufferDez & $3F)) ;64-er Alignment
   Else 
    sAusgabe = "-1"
 ProcedureReturn sAusgabe
  EndIf
  PUSH BufferDezA
  !pop [BufferDezA]
  ;========================================================
  ;Ergebnis für Ausgabe aufbereiten
  ;Konvertierung Double in Integer, Rundung gemäß Rounding Control Bits im MXCSR-Register, oben auf "nearest" gesetzt

  ;da BufferFaktor2A noch den I-Anteil enthält, wird hier "komprimiert", d.h. in die I-Bereiche werden die Produkte reinkopiert
  ;die Double-Werte in Long-Integer konvertieren, Division ist notwendig für höhere Genauigkeit; nicht später shiften!
  !mov r15,[N] 
  !mov r14,[BufferFaktor2A]
  !mov r8,r14
  !add r14,16                ;der 1.Wert bleibt natürlich
  !add r8,8
 !@@:
  !mov rax,[r14]
  !mov [r8],rax
  !add r14,16
  !add r8,8
  !dec r15
 !jnz @b

  ;die Double-Werte in Long-Integer konvertieren, Division ist notwendig für Genauigkeit; nicht später shiften!
  !mov r15,[N] 
  !mov r14,[BufferFaktor2A]
  !shr r15,2
  !mov rcx,r14
 !@@:
  !vmovapd ymm0,[r14]
  !vbroadcastsd ymm1,[ND]
  !vdivpd ymm2,ymm0,ymm1
  !vcvtpd2dq xmm3,ymm2
  !vmovapd [rcx],xmm3
  !add r14,32
  !add rcx,16
  !dec r15
 !jnz @b

  ;Hex2Dez und "Verstringung"
  ;Pointer_Produkt = N - 1
  ;For i = 0 To N - 1         ;vom Produkt alle Werte durch
  !mov r15,[N]
  !mov r14,[BufferFaktor2A]
  !mov rdx,[BufferProduktA]
  !mov rdi,r15
  !dec rdi
  !lea r10,[Dezi]
  !mov r11,[BufferDezA]

 !.DeziLoop: 
    !mov eax,[r14]           ;eax=Integer-Wert
    ;Hex2Dez
    ;Ziffer = 0
    !xor r12b,r12b 
    ;Pointer1 = 0
    !xor r8,r8
    ;Pointer2 = 0
    !xor r9,r9
 
    ;While PeekQ(?Dezi + Pointer1) <> 0 oder Zähler wie jetzt
    !mov rsi,10               ;Anzahl Dezi-Werte ohne Null, 9 reicht für M49, ansonsten auch Dezi anpassen an Faktoren(Produkt)-Länge!
   !.OuterLoopDez1:
      ;While (QuadWert - PeekQ(?Dezi + Pointer1)) >= 0
      !mov r13d,eax
     !.InnerLoopDez1:
        !sub r13d,[r10+r8]
       !js .InnerLoopDez1End
        ;Ziffer + 1
        !inc r12b
        ;QuadWert - PeekQ(?Dezi + Pointer1)
        ;!sub rax,[r10+r8]
        !mov eax,r13d         ;AL ist unten letzte Ziffer
      ;Wend     (QuadWert - PeekQ(?Dezi + Pointer1)) >= 0
      !jmp .InnerLoopDez1

      !.InnerLoopDez1End:
      ;PokeB(BufferDezA + Pointer2, Ziffer)
      !mov [r11+r9],r12b
      ;Pointer1 + 8
      !add r8,8
      ;Pointer2 + 1
      !inc r9
      ;Ziffer = 0
      !xor r12b,r12b 
      !dec rsi
    ;Wend      PeekQ(?Dezi + Pointer1) <> 0 oder Zähler wie jetzt
    !jnz .OuterLoopDez1

    ;PokeB(BufferDezA + Pointer2, QuadWert1)   ;letzte Dezimal-Ziffer abspeichern
    !mov [r11+r9],al         ;AL ist letzte Ziffer von oben
    ;Pointer_Produkt1 = Pointer_Produkt 
    !mov rsi,rdi
    !mov rax,10               ;Anzahl Dezi-Werte ohne Null, 9 reicht für M49, ansonsten auch Dezi anpassen an Faktoren(Produkt)-Länge!
    ;!mov rdx,[BufferProduktA]
    ;For k = 0 To 15
    !@@:
      ;Ziffer = PeekB(BufferProduktA + Pointer_Produkt1)
      !mov r12b,[rdx+rsi]
      ;Ziffer + PeekB(BufferDezA + Pointer2)
      !add r12b,[r11+r9]
      ;PokeB(BufferProduktA + Pointer_Produkt1, Ziffer)
      !mov [rdx+rsi],r12b
      ;If Ziffer > 9                           ;Übertrag
      !cmp r12b,9              ;Übertrag?
      !jbe .ZifferOK
        ;Ziffer - 10
        !sub r12b,10
        ;PokeB(BufferProduktA + Pointer_Produkt1, Ziffer)
        !mov [rdx+rsi],r12b
        ;Ziffer = PeekB(BufferProduktA + Pointer_Produkt1 - 1) + 1
        !mov r12b,[rdx+rsi-1]
        !inc r12b
        ;PokeB(BufferProduktA + Pointer_Produkt1 - 1, Ziffer)
        !mov [rdx+rsi-1],r12b
      ;EndIf
      !.ZifferOK:
      ;Pointer_Produkt1 - 1
      !dec rsi
      !js .Reicht
      ;Pointer2 - 1
      !dec r9
      !dec rax
    ;Next
    !jnz @b

    !.Reicht:
    ;Pointer_Produkt - 1
    !dec rdi
    !add r14,4
    !dec r15
  ;Next
  !jnz .DeziLoop

  ;------------------------
  ;Ergebnis (Produkt) als ASCII-Ziffern abspeichern
  !vbroadcastsd ymm0,[AsciiNull]  ;AVX  Missbrauch von Float, was solls... kostet aber Speed
  !mov r9,[BufferProduktA]
  !xor r10,r10
  !mov rcx,[N]
  !shr rcx,4                      ;/16, XMM-Length
 !@@:
  !vpor xmm1,xmm0,[r9+r10]
  !vmovdqa [r9+r10],xmm1
  !add r10,16
  !dec rcx
 !jnz @b
  !mov byte[r9+r10],0             ;Zero-Byte für String, zur Sicherheit, evtl.word

  FreeMemory(BufferFaktor2)
  FreeMemory(BufferDez)

  sAusgabe = LTrim(PeekS(BufferProduktA + N - LG, LG, #PB_Ascii), "0")    ;evtl. führende Null(en) entfernen

  If sAusgabe = ""                ;mindestens einer der Faktoren war Null, Kompatibilität zu "normaler" Multiplikation
    sAusgabe = "0"
  EndIf
  FreeMemory(BufferProdukt)

  ;Restore Registers
  !mov rax,[BufferFXA]
  !fxrstor64 [rax]
  !mov rdi,[rax+464]         ;restore callee-save registers
  !mov rsi,[rax+472]
  !mov r12,[rax+480]
  !mov r13,[rax+488]
  !mov r14,[rax+496]
  !mov r15,[rax+504]
  FreeMemory(BufferFX)

  ;Restore MXCSR             ;nach fxrstor64!
  !vldmxcsr [MXCSR_Sicher]   ;schreibt den gesicherten Wert zurück ins MXCSR-Register

  !vzeroupper
  DisableASM
 ProcedureReturn sAusgabe

  ;========================================================
  DataSection
!RezFakSin:                       ;für Sinus
  !dq  2.755731922398589065e-6    ; 1/9!   4.Iteration
  !dq -1.984126984126984127e-4    ;-1/7!   3.Iteration
  !dq  8.333333333333333333e-3    ; 1/5!   2.Iteration
  !dq -1.666666666666666667e-1    ;-1/3!   1.Iteration

  !dq  2.811457254345520763e-15   ; 1/17!  8.Iteration
  !dq -7.647163731819816476e-13   ;-1/15!  7.Iteration
  !dq  1.605904383682161460e-10   ; 1/13!  6.Iteration
  !dq -2.505210838544171878e-8    ;-1/11!  5.Iteration

  !dq  6.446950284384473396e-26   ; 1/25!  12.Iteration
  !dq -3.868170170630684038e-23   ;-1/23!  11.Iteration
  !dq  1.957294106339126123e-20   ; 1/21!  10.Iteration
  !dq -8.220635246624329717e-18   ;-1/19!  9.Iteration

!RezFakCos:                       ;für Cosinus
  !dq  2.480158730158730159e-5    ; 1/8!   4.Iteration
  !dq -1.388888888888888889e-3    ;-1/6!   3.Iteration
  !dq  4.166666666666666667e-2    ; 1/4!   2.Iteration
  !dq -5.000000000000000000e-1    ;-1/2!   1.Iteration

  !dq  4.779477332387385297e-14   ; 1/16!  8.Iteration
  !dq -1.147074559772972471e-11   ;-1/14!  7.Iteration
  !dq  2.087675698786809898e-9    ; 1/12!  6.Iteration
  !dq -2.755731922398589065e-7    ;-1/10!  5.Iteration

  !dq  1.611737571096118349e-24   ; 1/24!  12.Iteration
  !dq -8.896791392450573287e-22   ;-1/22!  11.Iteration
  !dq  4.110317623312164858e-19   ; 1/20!  10.Iteration
  !dq -1.561920696858622646e-16   ;-1/18!  9.Iteration
!Eins:
  !dq 1.0
!Pi2:
  !dq 6.283185307179586476925
!Null1:
  !dq 0h 
!Signum:
  !dq 8000000000000000h           ;für neg Sin/Cos
  !dq 0h
  !dq 8000000000000000h
!AsciiNull:
  !dq 3030303030303030h           ;für String
!Dezi:
  !dq  10000000000
  !dq   1000000000                ;für M49 ist der größte Integer-Wert etwas über 200.000.000
  !dq    100000000
  !dq     10000000
  !dq      1000000
  !dq       100000
  !dq        10000
  !dq         1000
  !dq          100
  !dq           10
  !dq            0                ;blanke Sicherheit

  !BufferCosSinA   dq ?
  !BufferDezA      dq ?
  !BufferFaktor1A  dq ?
  !BufferFaktor2A  dq ?
  !BufferFXA       dq ?
  !BufferProduktA  dq ?
  !BufferReverseA  dq ?
  !Exponent        dq ?
  !L1              dq ?
  !L2              dq ?
  !LenF            dq ?
  !N               dq ?
  !ND              dq ?
  !PF1             dq ?
  !PF2             dq ?
  !Pointer_Produkt dq ?

  !MXCSR           dd ?
  !MXCSR_Sicher    dd ?

  !V1              dd 55555555h 
  !V2              dd 33333333h
  !V3              dd 0F0F0F0Fh
  !V4              dd 00FF00FFh
  !VVier           dd 4
  EndDataSection
EndProcedure

Procedure Calc()             ;Aufruf durch Potenzen (Case 3) und Mersenne (Case 4)
  ZeitA = ElapsedMilliseconds()

  Stellen = Potenz * Len(Result$) + 1
  Stellen = Stellen + (64 - (Stellen & $3F))

  !mov rax,[v_Potenz]
  !bsr rdx,rax
  !mov [v_Bits_M],rdx

  CalcPot = Potenz           ;CalcPot wird geshiftet!

  BufferPot = AllocateMemory((Stellen * 2 * #Asc_Uni) + 64)
  If BufferPot
    BufferPotA = BufferPot + (64 - (BufferPot & $3F)) ;64-er Align
   Else
    MessageRequester("Fehler!", "Nicht genügend freier Arbeitsspeicher!", #MB_ICONERROR)
    End
  EndIf

  ListViewGadget(2, 10, 10, 250, 580)  ;Info-Fenster

  PokeS(BufferPotA + Stellen * #Asc_Uni, Result$)     ;hier noch Result$=Basis$, zum besseren Verständnis
  PokeS(BufferPotA, "1")     ;falls Bit0 von Potenz nicht gesetzt, X^0=1

  !shr [v_CalcPot],1
 !jnc @f
  PokeS(BufferPotA, Basis$)
  m = 1
 !@@:

  p = 1
  j = 2
  k = Bits_M - 1
  z = 0                      ;Zeile im ListViewGadget (Info)
  For i = 0 To k
    AddGadgetItem (2, -1, "Berechne " + Basis$ + "^" + Str(j))   
    SetGadgetState(2, z)
    z + 1
    Result$ = FFT_Mul(BufferPotA + Stellen * #Asc_Uni, BufferPotA + Stellen * #Asc_Uni)  ;Potenzen berechnen
    If Result$ = "-1"
      MessageRequester("Fehler!", "Nicht genügend freier Arbeitsspeicher!", #MB_ICONERROR)
      End
    EndIf
    PokeS(BufferPotA + Stellen * #Asc_Uni, Result$)
    j * 2
    !shr [v_CalcPot],1
   !jnc @f
    n + p
    AddGadgetItem (2, -1, "Berechne " + Basis$ + "^" + Str(n * 2 + m))   
    SetGadgetState(2, z)
    z + 1
    Result$ = FFT_Mul(BufferPotA, BufferPotA + Stellen * #Asc_Uni)   ;(Zwischen-) Ergebnis berechnen
    If Result$ = "-1"
      MessageRequester("Fehler!", "Nicht genügend freier Arbeitsspeicher!", #MB_ICONERROR)
      End
    EndIf
    If i < k
      PokeS(BufferPotA, Result$)
    EndIf
   !@@:
    p * 2
  Next

  FreeMemory(BufferPot)
  FreeGadget(2)
  ZeitE = ElapsedMilliseconds() - ZeitA
EndProcedure

Procedure.i RangeCharSSE42(StrPointer, ChrPointer)
  !mov rax,[p.v_ChrPointer]
  !movdqu xmm1,[rax]
  !mov rdx,[p.v_StrPointer]
  !mov rax, -16

 CompilerIf #PB_Compiler_Unicode                 ;Unicode
  !align 16                                      ;kann speed-mäßig was bringen
!@@:
  !add rax,16
  !pcmpistri xmm1,dqword[rdx+rax],00000101b      ;Bit0 = 1 und Bit1 = 0, also String-Chars werden als unsigned Words interpretiert, Bit2 = 1 und Bit3 = 0, also Test auf Ranges
  !ja @b
  !shr rax,1                                     ;Pointer in String muss hier halbiert werden (Word)
  !add rax,rcx
  !shl rcx,1                                     ;weil sonst nur 8 (Words), unten wird aber mit 16 verglichen
  CompilerElse                                   ;kein Unicode
  !align 16                                      ;kann speed-mäßig was bringen
!@@:
  !add rax,16
  !pcmpistri xmm1,dqword[rdx+rax],00000100b      ;Bit0 = 0 und Bit1 = 0, also String-Chars werden als unsigned Bytes interpretiert, Bit2 = 1 und Bit3 = 0, also Test auf Ranges
  !ja @b
  !add rax,rcx
 CompilerEndIf

  !cmp rcx,16                                    ;rcx=16 oder 8, wenn nichts gefunden! 8 schon erweitert
  !je @f
  !inc rax 
 ProcedureReturn                                 ;rax = Pos., beginnt bei 1
!@@:
  !xor rax,rax
 ProcedureReturn                                 ;rax = 0, Char nicht gefunden
EndProcedure 

Procedure FakInfo()
  GadgetZeile + 1
  AddGadgetItem (2, -1, "Berechne aktuell " + Str(FakA) + "!")
  SetGadgetState(2, GadgetZeile)
EndProcedure
Viel Spaß!
Helle

Edit 7.1.2018: Kosmetik
Edit 11.1.2018: Kosmetik und minimale Speicherplatz-Reduzierung
Zuletzt geändert von Helle am 11.01.2018 21:36, insgesamt 4-mal geändert.
Benutzeravatar
RSBasic
Admin
Beiträge: 8022
Registriert: 05.10.2006 18:55
Wohnort: Gernsbach
Kontaktdaten:

Re: FFT-Calculator

Beitrag von RSBasic »

Damit geht es bei mir auch. :allright:

Mein i5 2500k spuckt die Ziffern von Mersenne M50 nach einer Minute und 25 Sekunden aus.
Aus privaten Gründen habe ich leider nicht mehr so viel Zeit wie früher. Bitte habt Verständnis dafür.
Bild
Bild
Benutzeravatar
Helle
Beiträge: 566
Registriert: 11.11.2004 16:13
Wohnort: Magdeburg

Re: FFT-Calculator, AVX2-und AVX-Version

Beitrag von Helle »

Danke für die Rückmeldung! Mein Problem ist, das ich nie eine AVX-Version erstellt hatte (gleich AVX2). Habe mir eben den Code noch mal angeschaut und nichts Verdächtiges (Umsetzungs-Fehler) gefunden.
Eigentlich auch eine gute Zeit für die doch schon etwas betagte CPU!
Gruß
Helle
ccode_new
Beiträge: 1214
Registriert: 27.11.2016 18:13
Wohnort: Erzgebirge

Re: FFT-Calculator, AVX2-und AVX-Version

Beitrag von ccode_new »

Cool, daran sehe ich wie schön langsam meine Laptop-CPU ist.

119,501 Sekunden (ohne Debugger) also ca. 2 Minuten.

Wenn man hierbei die CPU-Kerne(Multithreading) und die GPU (z.B. CUDA-API oder OpenCL) zur Berechnung mit nutzen könnte (ja ich weiß das wäre wahrscheinlich zu hoch) wäre das mal ganz interessant wie weit man das optimieren könnte (wenn man es könnte (vom Wissen her))

Das ganze läuft ja vollkommen ohne Parallisierung. (1 Kern läuft z.B immer zu 100%)

Aber ansonsten Top!
Zuletzt geändert von ccode_new am 06.01.2018 23:54, insgesamt 1-mal geändert.
Betriebssysteme: div. Windows, Linux, Unix - Systeme

no Keyboard, press any key
no mouse, you need a cat
Benutzeravatar
RSBasic
Admin
Beiträge: 8022
Registriert: 05.10.2006 18:55
Wohnort: Gernsbach
Kontaktdaten:

Re: FFT-Calculator, AVX2-und AVX-Version

Beitrag von RSBasic »

@Helle
Wie heißt dein Prozessor?
Aus privaten Gründen habe ich leider nicht mehr so viel Zeit wie früher. Bitte habt Verständnis dafür.
Bild
Bild
Antworten