Math.Sin (Doppel) -Funktion für GPU

Vorwort


Ich musste den Lichtbogen mit erhöhter Genauigkeit auf dem Grafikkartenprozessor in Echtzeit berechnen.

Der Autor hat sich kein Ziel gesetzt, die Standardfunktion System.Math.Sin () (C #) zu übertreffen, und hat es nicht erreicht.

Das Ergebnis der Arbeit und meiner Wahl (für diejenigen, die nicht lesen wollen):

Sin_3 (rad)
using System; class Math_d { const double PI025 = Math.PI / 4; /// <summary> 2^17 = 131072 (1 ),   10000 ( ),  2^21 = 22097152 (16 )   +-1 ( ) (  ) </summary> const int length_mem = 22097152; const int length_mem_M1 = length_mem - 1; /// <summary>    sin,    . </summary> static double[] mem_sin; /// <summary>    cos,    . </summary> static double[] mem_cos; /// <summary>  ,   sin,    . </summary> public static void Initialise() { Ini_Mem_Sin(); Ini_Mem_Cos(); } /// <summary>       Cos,   . </summary> /// <param name="rad"></param> public static double Sin_3(double rad) { double rad_025; int i; //    //if (rad < 0) { rad = -rad + Math.PI; } i = (int)(rad / PI025); //   rad_025 = rad - PI025 * i; //     ( ) i = i & 7; //      8 //    switch (i) { case 0: return Sin_Lerp(rad_025); case 1: return Cos_Lerp(PI025 - rad_025); case 2: return Cos_Lerp(rad_025); case 3: return Sin_Lerp(PI025 - rad_025); case 4: return -Sin_Lerp(rad_025); case 5: return -Cos_Lerp(PI025 - rad_025); case 6: return -Cos_Lerp(rad_025); case 7: return -Sin_Lerp(PI025 - rad_025); } return 0; } /// <summary>   sin    </summary> static void Ini_Mem_Sin() { double rad; mem_sin = new double[length_mem]; for (int i = 0; i < length_mem; i++) { rad = (i * PI025) / length_mem_M1; mem_sin[i] = Math.Sin(rad); } } /// <summary>   cos    </summary> static void Ini_Mem_Cos() { double rad; mem_cos = new double[length_mem]; for (int i = 0; i < length_mem; i++) { rad = (i * PI025) / length_mem_M1; mem_cos[i] = Math.Cos(rad); } } /// <summary>      sin  0  pi/4. </summary> /// <param name="rad">     0  pi/4. </param> static double Sin_Lerp(double rad) { int i_0; int i_1; double i_0d; double percent; double a; double b; double s; percent = rad / PI025; i_0d = percent * length_mem_M1; i_0 = (int)i_0d; i_1 = i_0 + 1; a = mem_sin[i_0]; b = mem_sin[i_1]; s = i_0d - i_0; return Lerp(a, b, s); } /// <summary>      cos  0  pi/4. </summary> /// <param name="rad">     0  pi/4. </param> static double Cos_Lerp(double rad) { int i_0; int i_1; double i_0d; double percent; double a; double b; double s; percent = rad / PI025; i_0d = percent * length_mem_M1; i_0 = (int)i_0d; i_1 = i_0 + 1; a = mem_cos[i_0]; b = mem_cos[i_1]; s = i_0d - i_0; return Lerp(a, b, s); } /// <summary>      . (return a + s * (b - a)) </summary> /// <param name="a">  . </param> /// <param name="b">  . </param> /// <param name="s">  . 0 = a, 1 = b, 0.5 =   a  b. </param> public static double Lerp(double a, double b, double s) { return a + s * (b - a); } } 


Gründe für die Veröffentlichung


  • In der HLSL-Sprache gibt es keine Standard-Sin-Funktion für double (dies ist jedoch nicht korrekt).
  • Im Internet sind zu diesem Thema nur wenige Informationen verfügbar.

Überlegte Ansätze



Die analysierten Parameter


  • Genauigkeit mit Math.Sin
  • Geschwindigkeit in Bezug auf Math.Sin

Zusätzlich zur Analyse werden wir ihre Leistung verbessern.

Taylor zählt


Vorteile:

  • Höchste Präzision
    Mit dieser Funktion zur Berechnung des Sin-Werts kann der unendlich genaue Sin- Wert berechnet werden . Je mehr Iterationen durchlaufen werden, desto genauer wird der Wert am Ausgang erhalten (in einer Hypothese). In der Programmierpraxis lohnt es sich, Rundungsfehler von Berechnungen in Abhängigkeit von den verwendeten Parametertypen (double, float, decimal und andere) zu berücksichtigen.
  • Berechnet einen beliebigen Winkel
    Sie können einen beliebigen Wert als Argument für die Funktion eingeben, sodass die Eingabeparameter nicht überwacht werden müssen.
  • Unabhängigkeit
    Es sind keine vorläufigen Berechnungen erforderlich (wie nachstehend erläutert), und es ist häufig die Basis, auf der schnellere Funktionen zusammengesetzt werden.

Nachteile:

  • Sehr niedrige Geschwindigkeit (4-10%)
    Es sind viele Iterationen erforderlich, um die Genauigkeit näher an die Genauigkeit von Math.Sin heranzuführen. Dadurch arbeitet es 25-mal langsamer als die Standardfunktion.
  • Je größer der Winkel, desto geringer die Genauigkeit
    Je größer der in die Funktion eingegebene Winkel ist, desto mehr Iterationen sind erforderlich, um die gleiche Genauigkeit wie bei Math.Sin zu erzielen.

Originaler Look (Geschwindigkeit: 4%):

Die Standardfunktion beinhaltet das Berechnen von Fakultäten sowie das Erhöhen jeder Iteration auf eine Potenz.

Bild

Geändert (Geschwindigkeit: 10%):

Die Berechnung einiger Grade kann in einem Zyklus reduziert werden (a * = aa;), und andere Fakultäten können vorberechnet und in ein Array eingefügt werden, während das Ändern der Vorzeichen (+, -, +, ...) nicht auf eine Potenz angehoben werden kann und ihre Berechnung auch reduziert werden kann unter Verwendung der vorherigen Werte.

Das Ergebnis ist der folgende Code:

Sünde (rad, Schritte)
  // <summary>  ,    Fact </summary> static double[] fact; /// <summary>            . ///   rad,   . ///  ( Math): 4% (fps)  steps = 17 </summary> /// <param name="rad">   .      pi/4. </param> /// <param name="steps">  :  ,   .  pi/4   E-15  8. </param> public static double Sin(double rad, int steps) { double ret; double a; //,     double aa; // *  int i_f; //  int sign; // (  -  +,     = +) ret = 0; sign = -1; aa = rad * rad; a = rad; i_f = 1; //      for (int n = 0; n < steps; n++) { sign *= -1; ret += sign * a / Fact(i_f); a *= aa; i_f += 2; } return ret; } /// <summary>   (n!).  n > fact.Length,  -1. </summary> /// <param name="n"> ,     . </param> public static double Fact(int n) { if (n >= 0 && n < fact.Length) { return fact[n]; } else { Debug.Log("    . n = " + n + ",   = " + fact.Length); return -1; } } /// <summary>  . </summary> static void Init_Fact() { int steps; steps = 46; fact = new double[steps]; fact[0] = 1; for (int n = 1; n < steps; n++) { fact[n] = fact[n - 1] * n; } } 


Überlegene Ansicht (Geschwindigkeit: 19%):

Wir wissen, dass je kleiner der Winkel ist, desto weniger Iterationen erforderlich sind. Der kleinste Winkel, den wir brauchen, ist = 0,25 * PI, d.h. 45 Grad. Unter Berücksichtigung von Sin und Cos im Bereich von 45 Grad können wir alle Werte von -1 bis 1 für Sin (im Bereich von 2 * PI) erhalten. Dazu teilen wir den Kreis (2 * PI) in 8 Teile und geben für jeden Teil unsere eigene Methode zur Berechnung des Sinus an. Um die Berechnung zu beschleunigen, werden wir uns außerdem weigern, die Funktion zum Erhalten des Restes (%) zu verwenden (um die Position des Winkels innerhalb der Zone von 45 Grad zu erhalten):

Sin_2 (rad)
  // <summary>  ,    Fact </summary> static double[] fact; /// <summary>   Sin </summary> /// <param name="rad"></param> public static double Sin_2(double rad) { double rad_025; int i; //rad = rad % PI2; //% -   .  , fps   90  150 (  100 000   ) //rad_025 = rad % PI025; i = (int)(rad / PI025); rad_025 = rad - PI025 * i; i = i & 7; //     8  //    switch (i) { case 0: return Sin(rad_025, 8); case 1: return Cos(PI025 - rad_025, 8); case 2: return Cos(rad_025, 8); case 3: return Sin(PI025 - rad_025, 8); case 4: return -Sin(rad_025, 8); case 5: return -Cos(PI025 - rad_025, 8); case 6: return -Cos(rad_025, 8); case 7: return -Sin(PI025 - rad_025, 8); } return 0; } /// <summary>            . ///   rad,   . ///  ( Math): 10% (fps)  steps = 17 </summary> /// <param name="rad">   .      pi/4. </param> /// <param name="steps">  :  ,   .  pi/4   E-15  8. </param> public static double Sin(double rad, int steps) { double ret; double a; //,     double aa; // *  int i_f; //  int sign; // (  -  +,     = +) ret = 0; sign = -1; aa = rad * rad; a = rad; i_f = 1; //      for (int n = 0; n < steps; n++) { sign *= -1; ret += sign * a / Fact(i_f); a *= aa; i_f += 2; } return ret; } /// <summary>            . ///   rad,   . ///  ( Math): 10% (fps), 26% (test)  steps = 17 </summary> /// <param name="rad">   .      pi/4. </param> /// <param name="steps">  :  ,   .  pi/4   E-15  8. </param> public static double Cos(double rad, int steps) { double ret; double a; double aa; // *  int i_f; //  int sign; // (  -  +,     = +) ret = 0; sign = -1; aa = rad * rad; a = 1; i_f = 0; //      for (int n = 0; n < steps; n++) { sign *= -1; ret += sign * a / Fact(i_f); a *= aa; i_f += 2; } return ret; } /// <summary>   (n!).  n > fact.Length,  -1. </summary> /// <param name="n"> ,     . </param> public static double Fact(int n) { if (n >= 0 && n < fact.Length) { return fact[n]; } else { Debug.Log("    . n = " + n + ",   = " + fact.Length); return -1; } } /// <summary>  . </summary> static void Init_Fact() { int steps; steps = 46; fact = new double[steps]; fact[0] = 1; for (int n = 1; n < steps; n++) { fact[n] = fact[n - 1] * n; } } 


Polynome


Ich bin im Internet auf diese Methode gestoßen. Der Autor benötigte eine schnelle Sin-Suchfunktion für eine doppelte Genauigkeit (Fehler <0,000 001), ohne Bibliotheken mit vorberechneten Werten zu verwenden.

Vorteile:

  • Hohe Geschwindigkeit (9-84%)
    Anfänglich zeigte das geworfene Polynom ohne Änderungen eine Geschwindigkeit von 9% des ursprünglichen Math.Sin, was 10-mal langsamer ist. Dank kleiner Änderungen steigt die Geschwindigkeit stark auf 84%, was nicht schlecht ist, wenn Sie Ihre Augen genau schließen.
  • Keine zusätzlichen vorläufigen Berechnungen und Speicher erforderlich
    Wenn wir oben und unten Arrays von Variablen zusammenstellen müssen, um die Berechnungen zu beschleunigen, wurden hier alle Schlüsselkoeffizienten freundlicherweise berechnet und vom Autor selbst in Form von Konstanten in die Formel eingefügt.
  • Genauigkeit höher als Mathf.Sin (float)
    Zum Vergleich:

    0,84147 1 - Mathf.Sin (1) (Unity-Engine);
    0,841470984807897 - Math.Sin (1) (Standard-C # -Funktion);
    0,8414709 56802368 - sin (1) (GPU, hlsl-Sprache);
    0,84147 1184637935 - Sin_0 (1) .

Nachteile:

  • Nicht universell
    Sie können die Genauigkeit nicht manuell anpassen, da nicht bekannt ist, mit welchen Werkzeugen der Autor dieses Polynom berechnet hat.
  • Warum?
    Warum brauchte der Autor eine Funktion, die keine Arrays benötigt und eine so geringe (im Vergleich zur doppelten) Genauigkeit aufweist?

Originalansicht:

Sin_1 (x)
 /// <summary>  ( Math): 9% (fps)</summary> /// <param name="x">     -2*Pi  2*Pi </param> public static double Sin_1(double x) { return 0.9999997192673006 * x - 0.1666657564532464 * Math.Pow(x, 3) + 0.008332803647181511 * Math.Pow(x, 5) - 0.00019830197237204295 * Math.Pow(x, 7) + 2.7444305061093514e-6 * Math.Pow(x, 9) - 2.442176561869478e-8 * Math.Pow(x, 11) + 1.407555708887347e-10 * Math.Pow(x, 13) - 4.240664814288337e-13 * Math.Pow(x, 15); } 


Überlegene Aussicht:

Sin_0 (rad)
 /// <summary>  ( Math): 83% (fps)</summary> /// <param name="rad">     -2*Pi  2*Pi </param> public static double Sin_0(double rad) { double x; double xx; double ret; xx = rad * rad; x = rad; //1 ret = 0.9999997192673006 * x; x *= xx; //3 ret -= 0.1666657564532464 * x; x *= xx; //5 ret += 0.008332803647181511 * x; x *= xx; //7 ret -= 0.00019830197237204295 * x; x *= xx; //9 ret += 2.7444305061093514e-6 * x; x *= xx; //11 ret -= 2.442176561869478e-8 * x; x *= xx; //13 ret += 1.407555708887347e-10 * x; x *= xx; //15 ret -= 4.240664814288337e-13 * x; return ret; } 


Lineare Interpolation


Diese Methode basiert auf einer linearen Interpolation zwischen den Ergebnissen zweier Datensätze in einem Array.
Die Einträge sind in mem_sin und mem_cos unterteilt und enthalten die vorberechneten Ergebnisse der Standardfunktion Math.Sin und Math.Cos für ein Segment von Eingabeparametern von 0 bis 0,25 * PI.

Das Prinzip der Manipulationen mit einem Winkel von 0 bis 45 Grad unterscheidet sich nicht von der verbesserten Version der Taylor-Reihe, aber gleichzeitig wird eine Funktion aufgerufen, die findet - zwischen zwei Datensätzen befindet sich ein Winkel - und einen Wert zwischen ihnen findet.

Vorteile:

  • Hohe Geschwindigkeit (65%)
    Aufgrund der Einfachheit des Interpolationsalgorithmus erreicht die Geschwindigkeit 65% der Geschwindigkeit von Math.Sin. Ich denke die Geschwindigkeit> 33% ist zufriedenstellend.
  • Höchste Präzision
    Ein Beispiel für einen seltenen Fall der Ablehnung:
    0,255835595715180 - Math.Sin;
    0,2558355957151 79 - Sin_3 .
  • Schneller Fuß
    Ich liebe diese Funktion, weil sie in Qual geboren wurde, von mir geschrieben wurde und die Anforderungen übertraf: Geschwindigkeit> 33%, Genauigkeit über 1e-14. Ich werde ihr einen stolzen Namen geben - Vēlōx Pes.

Nachteile:

  • Es erfordert einen Platz im Gedächtnis
    Um zu arbeiten, müssen Sie zuerst zwei Arrays berechnen: für sin und für cos; Jedes Array wiegt ~ 16 MB (16 * 2 = 32 MB).

Originalansicht:

Sin_3 (rad)
 class Math_d { const double PI025 = Math.PI / 4; /// <summary> 2^17 = 131072 (1 ),   10000 ( ),  2^21 = 22097152 (16 )   +-1 ( ) (  ) </summary> const int length_mem = 22097152; const int length_mem_M1 = length_mem - 1; /// <summary>    sin,    . </summary> static double[] mem_sin; /// <summary>    cos,    . </summary> static double[] mem_cos; /// <summary>  ,   sin,    . </summary> public static void Initialise() { Ini_Mem_Sin(); Ini_Mem_Cos(); } /// <summary>       Cos,   . </summary> /// <param name="rad"></param> public static double Sin_3(double rad) { double rad_025; int i; //    //if (rad < 0) { rad = -rad + Math.PI; } i = (int)(rad / PI025); //   rad_025 = rad - PI025 * i; //     ( ) i = i & 7; //      8 //    switch (i) { case 0: return Sin_Lerp(rad_025); case 1: return Cos_Lerp(PI025 - rad_025); case 2: return Cos_Lerp(rad_025); case 3: return Sin_Lerp(PI025 - rad_025); case 4: return -Sin_Lerp(rad_025); case 5: return -Cos_Lerp(PI025 - rad_025); case 6: return -Cos_Lerp(rad_025); case 7: return -Sin_Lerp(PI025 - rad_025); } return 0; } /// <summary>   sin    </summary> static void Ini_Mem_Sin() { double rad; mem_sin = new double[length_mem]; for (int i = 0; i < length_mem; i++) { rad = (i * PI025) / length_mem_M1; mem_sin[i] = Math.Sin(rad); } } /// <summary>   cos    </summary> static void Ini_Mem_Cos() { double rad; mem_cos = new double[length_mem]; for (int i = 0; i < length_mem; i++) { rad = (i * PI025) / length_mem_M1; mem_cos[i] = Math.Cos(rad); } } /// <summary>      sin  0  pi/4. </summary> /// <param name="rad">     0  pi/4. </param> static double Sin_Lerp(double rad) { int i_0; int i_1; double i_0d; double percent; double a; double b; double s; percent = rad / PI025; i_0d = percent * length_mem_M1; i_0 = (int)i_0d; i_1 = i_0 + 1; a = mem_sin[i_0]; b = mem_sin[i_1]; s = i_0d - i_0; return Lerp(a, b, s); } /// <summary>      cos  0  pi/4. </summary> /// <param name="rad">     0  pi/4. </param> static double Cos_Lerp(double rad) { int i_0; int i_1; double i_0d; double percent; double a; double b; double s; percent = rad / PI025; i_0d = percent * length_mem_M1; i_0 = (int)i_0d; i_1 = i_0 + 1; a = mem_cos[i_0]; b = mem_cos[i_1]; s = i_0d - i_0; return Lerp(a, b, s); } /// <summary>      . (return a + s * (b - a)) </summary> /// <param name="a">  . </param> /// <param name="b">  . </param> /// <param name="s">  . 0 = a, 1 = b, 0.5 =   a  b. </param> public static double Lerp(double a, double b, double s) { return a + s * (b - a); } } 



UPD: Fehler beim Ermitteln des Index in Sin_Lerp (), Cos_Lerp (), Ini_Mem_Sin () und Ini_Mem_Cos () behoben.

Source: https://habr.com/ru/post/de426355/


All Articles