Prólogo
Necesitaba calcular el arco con mayor precisión en el procesador de la tarjeta de video en tiempo real.
El autor no estableció una meta para superar la función estándar System.Math.Sin () (C #) y no la alcanzó.
El resultado del trabajo y mi elección (para aquellos que no quieren leer):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); } }
Razones para su publicación.
- No existe una función Sin estándar para el doble en el lenguaje HLSL (pero esto no es exacto)
- Hay poca información disponible en Internet sobre este tema.
Enfoques considerados
Los parámetros analizados
- Precisión con Math.Sin
- Velocidad en relación con las matemáticas.
Además del análisis, mejoraremos su rendimiento.Taylor clasifica
Pros:
La más alta precisiónEsta función, utilizada para calcular el valor Sin, se puede utilizar para calcular el valor Sin infinitamente preciso . Cuantas más iteraciones sufra, más exactamente se obtendrá el valor en la salida (en una hipótesis). En la práctica de programación, vale la pena considerar los errores de redondeo de los cálculos dependiendo de los tipos de parámetros utilizados (doble, flotante, decimal y otros).
Calcula cualquier ánguloPuede ingresar cualquier valor como argumento de la función, por lo que no es necesario monitorear los parámetros de entrada.
IndependenciaNo requiere cálculos preliminares (como las funciones que se analizan a continuación), y a menudo es la base sobre la cual se ensamblan las funciones más rápidas.
Contras:
Muy baja velocidad (4-10%)Se necesitan muchas iteraciones para acercar la precisión a la precisión de Math.Sin, por lo que funciona 25 veces más lento que la función estándar.
Cuanto mayor es el ángulo, menor es la precisiónCuanto mayor sea el ángulo ingresado en la función, se necesitan más iteraciones para lograr la misma precisión que Math.Sin.
Aspecto original (velocidad: 4%):La función estándar implica calcular factoriales, así como elevar a una potencia cada iteración.
Modificado (velocidad: 10%):El cálculo de algunos grados se puede reducir en un ciclo (a * = aa;), y otros factores se pueden calcular previamente y poner en una matriz, mientras que el cambio de los signos (+, -, +, ...) no se puede elevar a una potencia y su cálculo también se puede reducir utilizando los valores anteriores.
El resultado es el siguiente código:
Sin (rad, pasos) // <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; } }
Vista superior (velocidad: 19%):Sabemos que cuanto menor es el ángulo, menos iteraciones son necesarias. El ángulo más pequeño que necesitamos es = 0.25 * PI, es decir 45 grados Considerando Sin y Cos en la región de 45 grados, podemos obtener todos los valores de -1 a 1 para Sin (en la región de 2 * PI). Para hacer esto, dividimos el círculo (2 * PI) en 8 partes y para cada parte indicamos nuestro propio método de cálculo del seno. Además, para acelerar el cálculo, nos negaremos a utilizar la función de obtener el resto (%) (para obtener la posición del ángulo dentro de la zona de 45 grados):
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; } }
Polinomios
Encontré este método en Internet, el autor necesitaba una función de búsqueda rápida de Sin para el doble de baja precisión (error <0.000 001) sin usar bibliotecas de valores precalculados.
Pros:
Alta velocidad (9-84%)Inicialmente, el polinomio arrojado sin cambios mostró una velocidad del 9% del Math.Sin original, que es 10 veces más lento. Gracias a pequeños cambios, la velocidad aumenta bruscamente al 84%, lo que no está mal si cierra los ojos con precisión.
No se requieren cálculos preliminares adicionales ni memoriaSi en la parte superior y más abajo necesitamos componer matrices de variables para acelerar los cálculos, entonces todos los coeficientes clave fueron calculados amablemente y colocados en la fórmula por el propio autor en forma de constantes.
Precisión superior a Mathf.Sin (flotante)Para comparar:
0.84147 1 - Mathf.Sin (1) (motor de Unity);
0.841470984807897 - Math.Sin (1) (función estándar de C #);
0.8414709 56802368 - sin (1) (GPU, lenguaje hlsl);
0.84147 1184637935 - Sin_0 (1) .
Contras:
No universalNo puede ajustar la precisión manualmente porque no se sabe qué herramientas utilizó el autor para calcular este polinomio.
Por qué¿Por qué el autor necesitaba una función que no requiera ninguna matriz y que tenga una precisión tan baja (en comparación con el doble)?
Vista originalSin_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); }
Vista superiorSin_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; }
Interpolación lineal
Este método se basa en la interpolación lineal entre los resultados de dos registros en una matriz.
Las entradas se dividen en mem_sin y mem_cos, contienen los resultados precalculados de la función estándar Math.Sin y Math.Cos en un segmento de parámetros de entrada de 0 a 0.25 * PI.
El principio de las manipulaciones con un ángulo de 0 a 45 grados no difiere de la versión mejorada de la serie Taylor, pero al mismo tiempo se llama una función que encuentra, entre los cuales dos registros hay un ángulo, y encuentra un valor entre ellos.
Pros:
Alta velocidad (65%)Debido a la simplicidad del algoritmo de interpolación, la velocidad alcanza el 65% de la velocidad de Math.Sin. Creo que la velocidad> 33% es satisfactoria.
La más alta precisiónUn ejemplo de un caso raro de rechazo:
0.255835595715180 - Math.Sin;
0.2558355957151 79 - Sin_3 .
Pie rápidoMe encanta esta función porque nació en agonía, escrita por mí y superó los requisitos: velocidad> 33%, precisión superior a 1e-14. Le daré un nombre orgulloso: Vēlōx Pes.
Contras:
Requiere un lugar en la memoriaPara trabajar, primero debe calcular dos matrices: para sin y para cos; cada matriz pesa ~ 16 mb (16 * 2 = 32 mb)
Vista originalSin_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: se corrigió el error al determinar el índice en Sin_Lerp (), Cos_Lerp (), Ini_Mem_Sin () e Ini_Mem_Cos ().