¿Cuál es la forma más rápida de calcular sin y cos juntos?

100

Me gustaría calcular el seno y el coseno de un valor juntos (por ejemplo, para crear una matriz de rotación). Por supuesto, podría calcularlos por separado uno tras otro a = cos(x); b = sin(x);, pero me pregunto si hay una forma más rápida cuando se necesitan ambos valores.

Editar: para resumir las respuestas hasta ahora:

  • Vlad dijo que existe el comando asm que losFSINCOScalcula a ambos (casi al mismo tiempo que una llamada aFSINsolo)

  • Como advirtió Chi , esta optimización a veces ya la realiza el compilador (cuando se utilizan indicadores de optimización).

  • caf señaló, que las funcionessincosysincosfprobablemente están disponibles y se pueden llamar directamente con solo incluirmath.h

  • El enfoque de tanascius del uso de una tabla de consulta es controvertido. (Sin embargo, en mi computadora y en un escenario de referencia, se ejecuta 3 veces más rápido quesincoscon casi la misma precisión para puntos flotantes de 32 bits).

  • Joel Goodwin se vinculó a un enfoque interesante de una técnica de aproximación extremadamente rápida con una precisión bastante buena (para mí, esto es incluso más rápido que la búsqueda de tablas)

Danvil
fuente
1
Consulte también esta pregunta sobre la implementación nativa de sin / cos: stackoverflow.com/questions/1640595
Joel Goodwin
1
pruebe sinx ~ x-x^3/6y cosx~1-x^2/4como aproximaciones si le importa la velocidad más que la precisión. Puede agregar términos en cualquiera de las series a medida que le da más peso a la precisión ( en.wikipedia.org/wiki/Taylor_series, desplácese hacia abajo hasta la serie trig taylor). Tenga en cuenta que esta es una forma general de aproximar cualquier función que desee que sea de ntiempos diferenciables . Entonces, si tiene una función más grande a la que pertenecen el seno y el coseno, obtendrá una velocidad mucho mayor si la aproxima en lugar del pecado, el cos de forma independiente.
ldog
Esta es una técnica pobre con muy poca precisión. Ver publicación de Joel Goodwin. Las series de Taylor se han publicado a continuación. Por favor publíquelo como respuesta.
Danvil
1
Bueno, depende de sus requisitos, si desea precisión, la serie de Taylor será una buena aproximación solo si necesita valores xcercanos a algún punto x_0, luego expanda su serie de Taylor en x_0lugar de 0. Esto le dará una excelente precisión cerca, x_0pero cuanto más lejos esté vaya peor los resultados. Probablemente pensó que la precisión es una porquería cuando miró la respuesta dada y la probó con valores alejados de 0. Esa respuesta es sin, cos expandido alrededor de 0.
ldog

Respuestas:

52

Los procesadores Intel / AMD modernos tienen instrucciones FSINCOSpara calcular funciones de seno y coseno simultáneamente. Si necesita una optimización sólida, tal vez debería usarla.

Aquí hay un pequeño ejemplo: http://home.broadpark.no/~alein/fsincos.html

Aquí hay otro ejemplo (para MSVC): http://www.codeguru.com/forum/showthread.php?t=328669

Aquí hay otro ejemplo (con gcc): http://www.allegro.cc/forums/thread/588470

Espero que alguno de ellos ayude. (No utilicé esta instrucción, lo siento).

Como son compatibles a nivel de procesador, espero que sean mucho más rápidos que las búsquedas de tablas.

Editar:
Wikipedia sugiere que FSINCOSse agregó en 387 procesadores, por lo que difícilmente puede encontrar un procesador que no lo admita.

Editar:
la documentación de Intel indica que FSINCOSes aproximadamente 5 veces más lento que FDIV(es decir, división de punto flotante).

Editar:
tenga en cuenta que no todos los compiladores modernos optimizan el cálculo de seno y coseno en una llamada a FSINCOS. En particular, mi VS 2008 no lo hizo de esa manera.

Editar:
El primer enlace de ejemplo está muerto, pero todavía hay una versión en Wayback Machine .

Vlad
fuente
1
@phkahler: Eso sería genial. No sé si los compiladores modernos utilizan dicha optimización.
Vlad
12
La fsincosinstrucción no es "bastante rápida". El propio manual de optimización de Intel lo cita como que requiere entre 119 y 250 ciclos en microarquitecturas recientes. La biblioteca matemática de Intel (distribuida con ICC), en comparación, puede calcular por separadosin y cosen menos de 100 ciclos, utilizando una implementación de software que usa SSE en lugar de la unidad x87. Una implementación de software similar que calculó ambos simultáneamente podría ser aún más rápida.
Stephen Canon
2
@Vlad: Las bibliotecas matemáticas ICC no son de código abierto y no tengo una licencia para redistribuirlas, por lo que no puedo publicar el ensamblaje. sinSin embargo, puedo decirles que no hay ningún cálculo integrado que puedan aprovechar; utilizan las mismas instrucciones SSE que todos los demás. Para su segundo comentario, la velocidad relativa a fdives irrelevante; si hay dos formas de hacer algo y una es dos veces más rápida que la otra, no tiene sentido llamar "rápida" a la más lenta, independientemente del tiempo que tome en relación con una tarea completamente no relacionada.
Stephen Canon
1
La sinfunción de software de su biblioteca ofrece una precisión total de doble precisión. La fsincosinstrucción ofrece algo más de precisión (doble extendida), pero esa precisión adicional se descarta en la mayoría de los programas que llaman a la sinfunción, ya que su resultado generalmente se redondea a doble precisión mediante operaciones aritméticas posteriores o un almacenamiento en la memoria. En la mayoría de las situaciones, ofrecen la misma precisión para un uso práctico.
Stephen Canon
4
Tenga en cuenta también que fsincosno es una implementación completa por sí sola; necesita un paso de reducción de rango adicional para poner el argumento en el rango de entrada válido para la fsincosinstrucción. La biblioteca siny las cosfunciones incluyen esta reducción, así como el cálculo del núcleo, por lo que son incluso más rápidos (en comparación) de lo que podrían indicar los tiempos de ciclo que enumeré.
Stephen Canon
39

Los procesadores x86 modernos tienen una instrucción fsincos que hará exactamente lo que está pidiendo: calcular sin y cos al mismo tiempo. Un buen compilador de optimización debería detectar el código que calcula sin y cos para el mismo valor y usar el comando fsincos para ejecutarlo.

Se necesitaron algunos juegos de banderas del compilador para que esto funcionara, pero:

$ gcc --version
i686-apple-darwin9-gcc-4.0.1 (GCC) 4.0.1 (Apple Inc. build 5488)
Copyright (C) 2005 Free Software Foundation, Inc.
This is free software; see the source for copying conditions.  There is NO
warranty; not even for MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.

$ cat main.c
#include <math.h> 

struct Sin_cos {double sin; double cos;};

struct Sin_cos fsincos(double val) {
  struct Sin_cos r;
  r.sin = sin(val);
  r.cos = cos(val);
  return r;
}

$ gcc -c -S -O3 -ffast-math -mfpmath=387 main.c -o main.s

$ cat main.s
    .text
    .align 4,0x90
.globl _fsincos
_fsincos:
    pushl   %ebp
    movl    %esp, %ebp
    fldl    12(%ebp)
    fsincos
    movl    8(%ebp), %eax
    fstpl   8(%eax)
    fstpl   (%eax)
    leave
    ret $4
    .subsections_via_symbols

¡Tada, usa la instrucción fsincos!

Chi
fuente
¡Esto es genial! ¿Podría explicar qué está haciendo -mfpmath = 387? ¿Y también funciona con MSVC?
Danvil
1
Tenga en cuenta eso -ffast-mathy -mfpmathlleve a resultados diferentes en algunos casos.
Debilski
3
mfpmath = 387 forzará a gcc a usar instrucciones x87 en lugar de instrucciones SSE. Sospecho que MSVC tiene optimizaciones y banderas similares, pero no tengo MSVC a mano para estar seguro. Sin embargo, es probable que el uso de instrucciones x87 sea un detrimento del rendimiento en otro código, también debe mirar mi otra respuesta, para usar MKL de Intel.
Chi
Mi viejo gcc 3.4.4 de cygwin produce 2 llamadas separadas a fsiny fcos. :-(
Vlad
Probado con Visual Studio 2008 con las optimizaciones más altas habilitadas. Llama a 2 funciones de biblioteca __CIsiny __CIcos.
Vlad
13

Cuando necesite rendimiento, puede usar una tabla sin / cos precalculada (una tabla servirá, almacenada como un diccionario). Bueno, depende de la precisión que necesites (quizás la mesa sea demasiado grande), pero debería ser muy rápido.

tanascius
fuente
Luego, el valor de entrada debe asignarse a [0,2 * pi] (o menos con verificaciones adicionales) y esta llamada a fmod destruye el rendimiento. En mi implementación (probablemente subóptima) no pude obtener rendimiento con la tabla de búsqueda. ¿Tiene algún consejo aquí?
Danvil
11
Es casi seguro que una tabla precalculada sea más lenta que simplemente llamar sinporque la tabla precalculada destruirá la caché.
Andreas Brinck
1
Depende del tamaño de la mesa. Una tabla de 256 entradas suele ser lo suficientemente precisa y usa solo 1 Kb ... si la usa mucho, ¿no se quedaría atascada en la caché sin afectar negativamente al resto del rendimiento de la aplicación?
Mr. Boy
@Danvil: Aquí hay un ejemplo de una tabla de búsqueda sinusoidal en.wikipedia.org/wiki/Lookup_table#Computing_sines . Sin embargo, asume que ya asignó su entrada a [0; 2pi] también.
tanascius
@AndreasBrinck No iría tan lejos. Depende (TM). Los cachés modernos son enormes y las tablas de búsqueda son pequeñas. Muy a menudo, si tiene un poco de cuidado en el diseño de la memoria, su tabla de búsqueda no necesita hacer ninguna diferencia en la utilización de la memoria caché del resto de su cálculo. El hecho de que la tabla de búsqueda encaje dentro de la caché es una de las razones por las que es tan rápida. Incluso en Java, donde es difícil controlar el diseño de mem con precisión, he tenido ganancias de rendimiento masivas con tablas de búsqueda.
Jarrod Smith
13

Técnicamente, lograrías esto usando números complejos y la fórmula de Euler . Por lo tanto, algo como (C ++)

complex<double> res = exp(complex<double>(0, x));
// or equivalent
complex<double> res = polar<double>(1, x);
double sin_x = res.imag();
double cos_x = res.real();

debería darte seno y coseno en un solo paso. Cómo se hace esto internamente depende del compilador y la biblioteca que se utilicen. Podría (y podría) llevar más tiempo hacerlo de esta manera (solo porque la Fórmula de Euler se usa principalmente para calcular el complejo expusando siny cos, y no al revés) pero podría haber alguna optimización teórica posible.


Editar

Los encabezados <complex>de GNU C ++ 4.2 utilizan cálculos explícitos de siny cosadentro polar, por lo que no se ve demasiado bien para las optimizaciones allí a menos que el compilador haga algo de magia (vea los interruptores -ffast-mathy -mfpmathcomo está escrito en la respuesta de Chi ).

Debilski
fuente
lo siento, pero la Fórmula de Euler en realidad no te dice cómo calcular algo, es solo una identidad (aunque muy útil) que relaciona exponenciales complejas con funciones trigonométricas reales. Hay beneficios de calcular el seno y el coseno juntos, pero involucran subexpresiones comunes y su respuesta no discute esto.
Jason S
12

Puede calcular cualquiera y luego usar la identidad:

cos (x) 2 = 1 - sin (x) 2

pero como dice @tanascius, una tabla precalculada es el camino a seguir.

Trigo Mitch
fuente
8
Y tenga en cuenta que el uso de este método implica calcular una potencia y una raíz cuadrada, por lo que si el rendimiento es importante, asegúrese de verificar que sea más rápido que calcular la otra función trigonométrica directamente.
Tyler McHenry
4
sqrt()a menudo está optimizado en hardware, por lo que puede ser más rápido entonces sin()o cos(). El poder es solo una auto multiplicación, así que no lo uses pow(). Existen algunos trucos para obtener raíces cuadradas razonablemente precisas muy rápidamente sin soporte de hardware. Por último, asegúrese de crear un perfil antes de hacer nada de esto.
deft_code
12
Tenga en cuenta que √ (1 - cos ^ 2 x) es menos preciso que calcular sen x directamente, en particular cuando x ~ 0.
kennytm
1
Para x pequeña, la serie de Taylor para y = sqrt (1-x * x) es muy buena. Puede obtener una buena precisión con los primeros 3 términos y solo requiere unas pocas multiplicaciones y un turno. Lo he usado en código de punto fijo.
phkahler
1
@phkahler: su serie de Taylor no se aplica porque cuando x ~ 0, cos x ~ 1.
kennytm
10

Si usa la biblioteca GNU C, entonces puede hacer:

#define _GNU_SOURCE
#include <math.h>

y obtendrá declaraciones de las funciones sincos(), sincosf()y sincosl()que calculan ambos valores juntos, presumiblemente de la manera más rápida para su arquitectura de destino.

coste y flete
fuente
8

Hay cosas muy interesantes en esta página del foro, que se centra en encontrar buenas aproximaciones que sean rápidas: http://www.devmaster.net/forums/showthread.php?t=5784

Descargo de responsabilidad: no utilicé nada de esto yo mismo.

Actualización 22 de febrero de 2018: Wayback Machine es la única forma de visitar la página original ahora: https://web.archive.org/web/20130927121234/http://devmaster.net/posts/9648/fast-and-accurate- seno-coseno

Joel Goodwin
fuente
También probé este y me dio un rendimiento bastante bueno. Pero sin y cos se calculan de forma independiente.
Danvil
Mi sensación es que este cálculo de seno / coseno será más rápido que obtener el seno y usar una aproximación de raíz cuadrada para obtener el coseno, pero una prueba lo verificará. La relación principal entre seno y coseno es de fase; ¿Es posible codificar para poder reutilizar los valores sinusoidales que calcula para las llamadas de coseno con desplazamiento de fase teniendo esto en cuenta? (Esto puede ser exagerado, pero tuve que preguntar)
Joel Goodwin
No directamente (a pesar de que la pregunta hace exactamente esto). Necesito sin y cos de un valor x y no hay forma de saber si en algún otro lugar
calculé
Lo usé en mi juego para dibujar un círculo de partículas. Dado que es solo un efecto visual, el resultado es lo suficientemente cercano y el rendimiento es realmente impresionante.
Maxim Kamalov
No me impresiona; Las aproximaciones de Chebyshev generalmente le brindan la mayor precisión para una interpretación determinada.
Jason S
7

Muchas bibliotecas matemáticas de C, como indica caf, ya tienen sincos (). La excepción notable es MSVC.

  • Sun ha tenido sincos () desde al menos 1987 (veintitrés años; tengo una página de manual impresa)
  • HPUX 11 lo tenía en 1997 (pero no está en HPUX 10.20)
  • Agregado a glibc en la versión 2.1 (febrero de 1999)
  • Se convirtió en una función integrada en gcc 3.4 (2004), __builtin_sincos ().

Y con respecto a la búsqueda, Eric S. Raymond en Art of Unix Programming (2004) (Capítulo 12) dice explícitamente que esta es una mala idea (en el momento actual):

"Otro ejemplo es el cálculo previo de tablas pequeñas; por ejemplo, una tabla de sin (x) por grado para optimizar las rotaciones en un motor de gráficos 3D ocuparía 365 × 4 bytes en una máquina moderna. Antes de que los procesadores fueran lo suficientemente más rápidos que la memoria para exigir el almacenamiento en caché , esta fue una optimización de velocidad obvia. Hoy en día, puede ser más rápido volver a calcular cada vez que pagar el porcentaje de pérdidas de caché adicionales causadas por la tabla.

"Pero en el futuro, esto podría cambiar de nuevo a medida que los cachés crezcan. De manera más general, muchas optimizaciones son temporales y pueden convertirse fácilmente en pesimizaciones a medida que cambian los índices de costos. La única forma de saberlo es medir y ver". (del arte de la programación Unix )

Pero, a juzgar por la discusión anterior, no todos están de acuerdo.

Joseph Quinsey
fuente
10
"365 x 4 bytes". Debe tener en cuenta los años bisiestos, por lo que en realidad debería ser 365,25 x 4 bytes. O tal vez pretendía usar el número de grados en un círculo en lugar del número de días en un año terrestre.
Ponkadoodle
@Wallacoloo: Buena observación. Me lo perdi. Pero el error está en el original .
Joseph Quinsey
LOL. Además, ignora el hecho de que en muchos de los juegos de computadora de esa área, solo necesitarás un número finito de ángulos. Entonces no hay pérdidas de caché, si conoce los posibles ángulos. Usaría tablas exactamente en este caso, y fsincosprobaría (¡instrucción de CPU!) Para los demás. A menudo es tan rápido como interpolar pecado y cos de una tabla grande.
Erich Schubert
5

No creo que las tablas de búsqueda sean necesariamente una buena idea para este problema. A menos que sus requisitos de precisión sean muy bajos, la mesa debe ser muy grande. Y las CPU modernas pueden realizar muchos cálculos mientras se obtiene un valor de la memoria principal. Esta no es una de esas preguntas que pueden responderse adecuadamente con un argumento (ni siquiera el mío), probar y medir y considerar los datos.

Pero buscaría las implementaciones rápidas de SinCos que se encuentran en bibliotecas como ACML de AMD y MKL de Intel.

Marca de alto rendimiento
fuente
3

Si está dispuesto a utilizar un producto comercial y está calculando una serie de cálculos sin / cos al mismo tiempo (para poder utilizar funciones vectoriales), debería consultar la biblioteca de kernel matemática de Intel.

Tiene una función sincos

De acuerdo con esa documentación, promedia 13.08 relojes / elemento en core 2 duo en modo de alta precisión, que creo que será incluso más rápido que fsincos.

Chi
fuente
1
De manera similar, en OSX se puede usar vvsincoso vvsincosfdesde Accelerate.framework. Creo que AMD también tiene funciones similares en su biblioteca de vectores.
Stephen Canon
2

Cuando el rendimiento es fundamental para este tipo de cosas, no es inusual introducir una tabla de búsqueda.

Tom Cabanski
fuente
2

Para un enfoque creativo, ¿qué tal expandir la serie Taylor? Dado que tienen términos similares, podría hacer algo como el siguiente pseudo:

numerator = x
denominator = 1
sine = x
cosine = 1
op = -1
fact = 1

while (not enough precision) {
    fact++
    denominator *= fact
    numerator *= x

    cosine += op * numerator / denominator

    fact++
    denominator *= fact
    numerator *= x

    sine += op * numerator / denominator

    op *= -1
}

Esto significa que haces algo como esto: comenzando en xy 1 para el pecado y el coseno, sigue el patrón: reste x ^ 2/2. del coseno, reste x ^ 3/3! desde el seno, agregue x ^ 4/4! al coseno, suma x ^ 5/5! a seno ...

No tengo idea de si esto funcionaría bien. Si necesita menos precisión que la incorporada sin () y cos (), puede ser una opción.

Tesserex
fuente
En realidad, el i-el factor de extensión del seno es x / i multiplicado por el i-el factor de extensión del coseno. Pero dudaría que usar la serie Taylor sea realmente rápido ...
Danvil
1
Chebyshev es mucho mejor que Taylor para la aproximación de funciones polinomiales. No utilice la aproximación de Taylor.
Timmmm
Hay un montón de pasos en falso numéricos aquí; tanto el numerador como el denominador se vuelven grandes rápidamente y eso conduce a errores de punto flotante. Sin mencionar cómo se decide qué es "precisión insuficiente" y cómo calcularla. La aproximación de Taylor es buena en la vecindad alrededor de un solo punto; lejos de ese punto, rápidamente se vuelven inexactos y requieren una gran cantidad de términos, razón por la cual la sugerencia de Timmmm sobre la aproximación de Chebyshev (que crea buenas aproximaciones en un intervalo dado) es buena.
Jason S
2

Hay una buena solución en la biblioteca CEPHES que puede ser bastante rápida y puede agregar / eliminar precisión de manera bastante flexible por un poco más / menos de tiempo de CPU.

Recuerde que cos (x) y sin (x) son las partes real e imaginaria de exp (ix). Entonces queremos calcular exp (ix) para obtener ambos. Calculamos previamente exp (iy) para algunos valores discretos de y entre 0 y 2pi. Cambiamos x al intervalo [0, 2pi). Luego seleccionamos la y que está más cerca de x y escribimos
exp (ix) = exp (iy + (ix-iy)) = exp (iy) exp (i (xy)).

Obtenemos exp (iy) de la tabla de búsqueda. Y desde | xy | es pequeña (como mucho la mitad de la distancia entre los valores de y), la serie de Taylor convergerá muy bien en unos pocos términos, por lo que la usamos para exp (i (xy)). Y luego solo necesitamos una multiplicación compleja para obtener exp (ix).

Otra buena propiedad de esto es que puedes vectorizarlo usando SSE.

Jsl
fuente
2

Es posible que desee echar un vistazo a http://gruntthepeon.free.fr/ssemath/ , que ofrece una implementación vectorizada SSE inspirada en la biblioteca CEPHES. Tiene buena precisión (desviación máxima de sin / cos en el orden de 5e-8) y velocidad (supera ligeramente a fsincos en una sola llamada y un claro ganador sobre múltiples valores).

SleuthEye
fuente
1

Se puede encontrar una aproximación precisa pero rápida de la función sin y cos simultáneamente, en javascript, aquí: http://danisraelmalta.github.io/Fmath/ (se importa fácilmente a c / c ++)

user2781980
fuente
0

¿Ha pensado en declarar tablas de búsqueda para las dos funciones? Aún tendría que "calcular" sin (x) y cos (x), pero sería decididamente más rápido, si no necesita un alto grado de precisión.

Frank Shearar
fuente
0

El compilador de MSVC puede utilizar las funciones SSE2 (internas)

 ___libm_sse2_sincos_ (for x86)
 __libm_sse2_sincos_  (for x64)

en compilaciones optimizadas si se especifican los indicadores del compilador apropiados (al mínimo / O2 / arch: SSE2 / fp: fast). Los nombres de estas funciones parecen implicar que no calculan sen y cos por separado, sino que ambos "en un solo paso".

Por ejemplo:

void sincos(double const x, double & s, double & c)
{
  s = std::sin(x);
  c = std::cos(x);
}

Ensamblado (para x86) con / fp: rápido:

movsd   xmm0, QWORD PTR _x$[esp-4]
call    ___libm_sse2_sincos_
mov     eax, DWORD PTR _s$[esp-4]
movsd   QWORD PTR [eax], xmm0
mov     eax, DWORD PTR _c$[esp-4]
shufpd  xmm0, xmm0, 1
movsd   QWORD PTR [eax], xmm0
ret     0

Ensamblaje (para x86) sin / fp: rápido pero con / fp: preciso en su lugar (que es el predeterminado) llama sin y cos por separado:

movsd   xmm0, QWORD PTR _x$[esp-4]
call    __libm_sse2_sin_precise
mov     eax, DWORD PTR _s$[esp-4]
movsd   QWORD PTR [eax], xmm0
movsd   xmm0, QWORD PTR _x$[esp-4]
call    __libm_sse2_cos_precise
mov     eax, DWORD PTR _c$[esp-4]
movsd   QWORD PTR [eax], xmm0
ret     0

Entonces / fp: fast es obligatorio para la optimización de sincos.

Pero tenga en cuenta que

___libm_sse2_sincos_

tal vez no sea tan preciso como

__libm_sse2_sin_precise
__libm_sse2_cos_precise

debido a la falta "precisa" al final de su nombre.

En mi sistema "ligeramente" más antiguo (Intel Core 2 Duo E6750) con el último compilador MSVC 2019 y las optimizaciones apropiadas, mi punto de referencia muestra que la llamada sincos es aproximadamente 2,4 veces más rápida que las llamadas sin y cos por separado.

xy
fuente