ASM-Tricks

    Teil 1: extrem schnelle Sinus Funktion
    Original: http://music.columbia.edu/pipermail/musi…ber/046673.html

    Die Funktion f=X^2 ergibt eine Parabel; das kann man nützen, um näherungsweise eine Sinuskurve zu berechnen.
    Der ASM-Code dazu ist sehr kurz und sieht so aus:


    Aufgerufen wird der Code mit folgender Funktion:

    AutoIt
    Func _Sin($fX)
	Local $aResult = DllCallAddress("int", $pASM_Sin, "uint", $fX * 2 ^ 31 / $cPI)
	Return $aResult[0] / -2 ^ 31
EndFunc   ;==>_Sin


    Das TestScript dazu ergibt dieses Bild:
    SinCos_1.au3
    SinCos_1.png


    Sieht schon ganz OK aus! Allerdings nicht ganz...
    Wenn man damit einen Kreis zeichnen möchte sieht man besser, dass eine Parabel doch keine Sinuskurve ist:
    SinCos_4.au3
    SinCos_4.png


    Es geht natürlich noch besser, um nicht zu sagen: Fast perfekt!
    SinCos_2.au3
    SinCos_2.png


    Nun wollen wir testen, wie schnell die Funktionen im Vergleich zu der Sinusberechnung der FPU ist.
    Dieses Script erzeugt eine Wav-Datei von mehreren Sekunden in den 3 Varianten FPU, Parabel, und "Fast Perfekt":
    SinCos_3.au3


    Wir sehen, dass die ASM-Versionen 10x schneller sind, als die Berechnung der FPU!
    Weiters kann man hören, dass die Parabel-Version harmonische Obertöne erzeugt, siehe FFT-Analyse im nächsten Post!


    Fazit: Die ParabelVersion ist am schnellsten, jedoch sehr ungenau. Das reicht etwa für Animationen, aber nicht für Sound. (Außer man kann mit den harmonischen Obertönen leben)
    Die Version "Fast Perfekt" ist wirklich nahezu perfekt! Nur ein einziger Oberton mit -70dB, das ist zu vernachlässigen. Der Geschwindigkeitsvorteil ist wirklich enorm!
    Nachteil bei beiden Versionen ist, dass man etwas tricksen muss, um am Ende das gewünschte Ergebnis zu erhalten. Das kann aber auch ein Vorteil sein: Siehe SinCos_3.au3: Bei der Variable $fPhaseInc kann man PI komplett wegkürzen.


    Hier die FFT-Analyse der 3 Wav-Dateien: FFT_1.png FFT_2.png FFT_3.png


    Wie funktioniert eigentlich die Parabel-Version?
    Ich hab die FUnktion in AutoIt nachgebaut:


    Bei den einzelnen Schritten kommt man nun zu diesen Outputs:
    Sin_1.png Sin_2.png Sin_3.png Sin_4.png Sin_5.png

    Erklärung kommt evtl später ;)

    Integer Division

    Zuerst noch ganz kurz Multiplikation:
    Man sollte prüfen, ob sich eine Multiplikation in mehrere Shifts zerlegen lässt.

    Shift_Left_X entspricht 2^x
    shl 1 = *2
    shl 2 = *4
    shl 3 = *8
    usw...

    z.B.: X*10 = (X*2) + (X*8)


    Wie Minx sehr schön gezeigt hat, kann man eine Division durch eine Multiplikation mit einer Konstanten erreichen.
    Ich habe dazu mein Helper-Script etwas grafisch aufgemotzt:

    Hi,
    wie hier bereits beschrieben, existiert für die Berechnung einer Folge von Werten und dessen Sinus/Cosinus ein mathematischer "Trick". Dabei wird lediglich für den Startwert der Sinus/Cosinus ermittelt und dann benötigt man zur Ermittlung der folgenden Werte lediglich Additionen und Multiplikationen.
    Als "Goodie" ist die DEG-Berechnung statt RAD bereits enthalten!

    Ggf. gibt es ja mittlerweile das Pendant zu PMADDWD (Multiply and Add Packed Integers) für Floats/Double bei den SSE-Befehlen^^. Damit wäre imho das geschwindigkeitsmäßige Optimum erreicht ohne jegliche Interpolation...
    Du kannst ja mal schauen, inwieweit du das integrieren kannst/willst!

    ciao
    Andy


    "Schlechtes Benehmen halten die Leute doch nur deswegen für eine Art Vorrecht, weil keiner ihnen aufs Maul haut." Klaus Kinski
    "Hint: Write comments after each line. So you can (better) see what your program does and what it not does. And we can see what you're thinking what your program does and we can point to the missunderstandings." A-Jay

    Wie man Fragen richtig stellt... Tutorial: Wie man Script-Fehler findet und beseitigt...X-Y-Problem

    2 Mal editiert, zuletzt von Andy (22. November 2015 um 15:39)

    Ein paar Notizen zu verschiedenen Geschwindigkeitsoptimierungen die ich mal aus einem älteren Buch abgeschrieben hatte.

    1. LEA

    LEA kann in vielen Fällen einfache Berechnungen ersetzten bzw. zusammenfassen. Z.b. kann man folgenden Code

    Code
    mov eax, ecx
shl eax, 3
add eax, ebx
sub eax, 1000

    zu

    Code
    lea eax, [ebx+8*ecx-1000]

    optimieren.

    2. INC / DEC

    INC und DEC verbrauchen auf so ziemlich allen modernen CPUs eine zusätzliche Mikroinstruktion. Daher lieber ADD / SUB benutzen, damit bleiben nämlich auch alle Flags intakt.

    Besonders interessant wirds hier, wenn ein Loop über eine Struct iteriert. Sofern der Step immer konstant ist, sollte man immer zu ADD greifen um zum Ende der Daten zu springen, statt INC bzw. MUL.

    3. XCHG

    Sollte man generell niemals verwenden. Denn XCHG impliziert ein vorhergehendes LOCK damit alle Threads synchronisiert werden. Das sorgt für kleinere, unnötige Verzögerungen. Wenn es um Dateigröße statt Speed geht, spart man mit XCHG hier und da auch ein ein Byte.

    4. Bit Tests

    Diese sind bis heute auf AMD CPUs nicht optimal umgesetzt und verursachen immer 2 zusätzliche Mikroinstruktionen (im Folgenden µOps). Hier lieber selbst mit TEST und bitwise Logik Hand anlegen. Wenn nur wenig getestet wird ists aber halb so schlimm.

    5. L/SAHF

    Wieder langsam auf AMD CPUs. Dort statt LAHF einfach mit SETcc arbeiten. Auf jeder CPU kann man SAHF durch TEST ersetzen, wenn nur auf ein Bit in AH geprüft werden soll.

    Nützlich ist auch FCOMI, dann dies hier:

    Code
    FCOM
FNSTSW AX
SAHF

    kann einfach durch ein FCOMI ersetzt werden. FCOMI kann man an äußerst vielen anderen Stellen benutzen um FPU Berechnungen zu optimieren, aber das erklären andere Tutorials besser.

    6. Integer Multiplikation

    MUL (und IMUL) sollte generell komplett vermieden werden und ist in jedem Fall durch schnelleren Code ersetzbar, so lange ein Faktor konstant ist. So ist

    Code
    IMUL EAX, 5

    gleich

    Code
    LEA EAX, [EAX+4*EAX]


    7. Division

    Ah ja, mein Lieblingsthema.

    7.1 Integer Division mit 2^N

    Das ist ein Spezialfall, der sich gut optimieren lässt. Hier kann man sogar gleichzeitig für Größe optimieren wenn man sich auf unsigned beschränkt:

    Code
    shr EAX, N

    Geht auch für signed, ist eben bloß ein bissl länger:

    Code
    ; Divide signed ints by 2^N
cdq
and EDX, (1 shl N) - 1
add EAX, EDX
sar EAX, N


    7.2 Integer Division mit beliebigen Konstanten

    Viele Wege führen nach Rom. Ich beschreibe hier mal den Mathisen Algorithmus. Dieser wird verwendet um ASM Code zu konstruieren, der nicht nur für die jeweilige Konstante optimiert ist, sondern auch auf die Breite (ie. wie viel Bit) der Zahl.

    Unsere Konstante ist d. Daraus erechnen wir mal folgende Zahlen:

    b = (signifikante bits in d) - 1
    r = w + b
    f = 2^r / d

    Dann ergeben sich drei Fälle, die letztendlich den Code vorgeben:

    Fall A, f ist ein Integer (frac(f) = 0)

    • Ergebnis = x SHR b


    Fall B, frac(f) < 0.5

    • f abrunden
    • Ergebnis = ((x + 1) * f) SHR r


    Fall C, frac(f) > 0.5

    • f aufrunden
    • Ergebnis = (x * f) SHR r


    Eigentlich ganz einfach. Aufgepasst: Das ist nur der erste Schritt. Nachdem wir hier fertig sind, müssen wir trotzdem noch alle oben genannten Optimierungen anwenden (z.B. das MUL loswerden).

    Exerzieren wir das mal durch. Wir wollen (EAX) durch 5 teilen:

    5 = 101b
    w = 32 (ie. die Breite in bits)
    b = 1 (da 5 drei signifikante bits hat - 1)
    r = 32 + 2 = 34
    f = 2^34 / 5

    f ist also 3435973836.8. Da frac(f) = 0.8 > 0.5 müssen wir Fall C anwenden, also runden wir f auf 0CCCCCCCDh auf.

    Also folgender simpler Code:

    Code
    ; x liegt als Wert in EAX
mov EDX, 0CCCCCCCDh           ;
mul EDX                       ; x * f
shr EDX, 2

    Warum nur SHR EDX, 2? Weil das Ergebnis in EDX schon um 32bits geshifted vorliegt. Somit brauchen wir nur noch um b zu shiften. Das Ergebnis liegt jetzt in EDX.

    Für den Fall B ergibt sich folgender Code:

    Code
    add EAX, 1
mov EDX, f
mul EDX
shr EDX, b

    Jetzt steht da aber noch ein hässliches MUL, das werden wir aber gleich los. Zunächst ein paar Hinweise.




    Optimierung I

    Obiger Code schließt leider x = 0FFFFFFFFh aus, da durch das ADD ein Overflow ensteht. Das lässt sich wie folgt umgehen (sollte man NUR einbauen, wenn man einen Overflow erwartet!):

    Code
    mov EDX, f
	add EAX, 1
	jc overflow
	mul EDX
overflow:
	shr EDX, b

    Optimierung II

    Es ergibt Sinn r kleiner zu wählen. Wenn nämlich r = w = 32, dann kann der letzte shift einfach weg :)

    Optimierung III

    Man kann auch mit IMUL was erreichen. Das hat den Vorteil, dass es vieles vom obigen Code einfach in einem Schritt erledigt, aber mit niedrigerer Präzision. Wenn also r = 16 + b, dann geht es ganz simpel:

    Code
    ; Dividieren durch 5
imul eax, 0CCCDh
shr EAX, 18

    Optimierung IV

    Vielleicht ist es schon aufgefallen, aber wenn man zwischen B und C wählen kann, sollte man immer so rechnen, dass man C verwendet. Warum? Weil man sich ein ADD spart. Dafür muss wieder r intelligent gewählt werden. So. Und nun optimieren wir noch das MUL weg und shiften stattdessen.

    Dazu ein Beispiel. Wir wollen 1.) durch 10 teilen mit 2.) niedriger Präzision für mehr Speed und 3.) das Ergebnis in EAX statt EDX:


    Optimierung V

    Das waren die Grundlagen, jetzt wollen wir mal mehrere Zahlen auf einmal teilen. Als Beispiel werden hier 8 unsigned 16bit Integer durch 100 geteilt:

    Code
    .data
align 16
HDIV dw 8 dup (0CCCDh)


.code
pmulhuw xmm0, HDIV
psrlw xmm0, 3

    Optimierung VI

    Logischerweise kann man damit jede Division beschleunigen, die wiederholt mit einer Konstante ausgeführt wird. In Loops z.B.. Dann einfach einmal obige Zwischenergebnisse ausrechnen und den Code anpassen. Wie man das im Algorithmen Design umsetzt ist hier erklärt.

    Optimierung VII

    Fast geschafft. Nur noch ein kurzer Exkurs mit SSE. Wenn man obige Methode auf SSE anwendet hat man normalerweise nur 12bit zur Verfügung. Intel hat dazu aber eine Empfehlung ausgesprochen (siehe AP-803), wie man mit dem Newton-Raphson Verfahren nochmal was rausholt um auf 23bits zu kommen. Die Berechnung dazu sieht so aus:

    x0 = rcpss(d)
    x1 = x0 * (2 - d * x0) = 2 * x0 - d * x0 * x0

    xN ist hier die N-beste Annäherung an das Reziproke von d.

    Umgesetzt wirds so:

    Code
    ; Ergebnisse werden in xmm0 stehen
movaps xmm1, [divisors]
rcpps xmm0, xmm1
mulps xmm1, xmm0
mulps xmm1, xmm0
addps xmm0, xmm0
subps xmm0, xmm1
mulps xmm0, [dividends]

    Das sind 4 23bit Divisionen in 18 Takten.


    8. Exponenten in Kommazahlen

    BTW. FSCALE immer vermeiden, einfach direkt in die Exponentenbits schreiben.

    Für alle 2^N (N ist signed) gilt folgendes:

    a) |N| < 2^7 - 1 (single precision)

    Code
    mov EAX, [N]
shl EAX, 23
add EAX, 3f800000h
mov dowrd ptr [TEMP], eax
fld dword ptr [TEMP]

    b) |N| < 2^10 - 1 (double precision)

    Code
    mov EAX, [N]
shl EAX, 20
add EAX, 3ff00000h
mov dword ptr [TEMP], 0
mod dword ptr [TEMP+4], EAX
fld qword ptr [TEMP]

    Und schon können wir eine allgemeine, schnelle und (so ziemlich) bug-freie Exponentialfunktion mit double-Genauigkeit ableiten:

    Zitat von Andy

    Hi,
    wie hier bereits beschrieben, existiert für die Berechnung einer Folge von Werten und dessen Sinus/Cosinus ein mathematischer "Trick". Dabei wird lediglich für den Startwert der Sinus/Cosinus ermittelt und dann benötigt man zur Ermittlung der folgenden Werte lediglich Additionen und Multiplikationen.
    Als "Goodie" ist die DEG-Berechnung statt RAD bereits enthalten!

    Ggf. gibt es ja mittlerweile das Pendant zu PMADDWD (Multiply and Add Packed Integers) für Floats/Double bei den SSE-Befehlen^^. Damit wäre imho das geschwindigkeitsmäßige Optimum erreicht ohne jegliche Interpolation...
    Du kannst ja mal schauen, inwieweit du das integrieren kannst/willst!


    Sehr interessant!
    Wenn ich das aber richtig verstehe, dann sind mit dieser Methode keine Sweeps möglich bzw. nur umständlich zu erreichen...
    Also sowas wie das hier:

    Zitat von eukalyptus

    Wenn ich das aber richtig verstehe, dann sind mit dieser Methode keine Sweeps möglich bzw. nur umständlich zu erreichen...

    Ja, da hast du Recht! Denn der neue "Startpunkt" für die "gleichmäßige Folge" ist ja immer wieder ein neuer....
    Der Trick funktioniert also nur bei gleichmäßigen Folgen :D
    Aber besser als nix, damit wären gefühlte 90% abgedeckt 8) . In Assembler(Compilersprache) heisst das Faktor 10 für reines Float/Double, bei Festkommazahlen rechnet man in "Integern" da ist dann wesentlich mehr drin, je nach Anwendung!