Newton Raphson con SSE2 - ¿puede alguien explicarme estas 3 líneas

Estoy leyendo este documento: http://software.intel.com/en-us/articles/interactive-ray-tracing

Y me topé con estas tres líneas de código:

La versión SIMD ya es bastante más rápida, pero podemos hacerlo mejor. Intel ha añadido una función rápida 1 / sqrt (x) al conjunto de instrucciones SSE2. El único inconveniente es que su precisión es limitada. Necesitamos el precisión, así que lo refinamos usando Newton-Rhapson:

 __m128 nr = _mm_rsqrt_ps( x ); 
 __m128 muls = _mm_mul_ps( _mm_mul_ps( x, nr ), nr ); 
 result = _mm_mul_ps( _mm_mul_ps( half, nr ), _mm_sub_ps( three, muls ) );

Este código asume la existencia de una variable _ _ m128 llamada 'half' (cuatro veces 0.5 f) y una variable 'tres' (cuatro veces 3.0 f).

Sé cómo usar Newton Raphson para calcular el cero de una función y sé cómo usarlo para calcular la raíz cuadrada de un número, pero no puedo ver cómo lo realiza este código.

¿Puede alguien explicármelo, por favor?

Author: Marco A., 2013-02-07
Source

2 answers

Dada la iteración de Newton y_n + 1=y_n(3-x (y_n)^2) / 2, debería ser bastante sencillo ver esto en el código fuente.

 __m128 nr   = _mm_rsqrt_ps( x );                  // The initial approximation y_0
 __m128 muls = _mm_mul_ps( _mm_mul_ps( x, nr ), nr ); // muls = x*nr*nr == x(y_n)^2
 result = _mm_mul_ps(
               _mm_sub_ps( three, muls )    // this is 3.0 - mul;
   /*multiplied by */ __mm_mul_ps(half,nr)  // y_0 / 2 or y_0 * 0.5
 );

Y para ser precisos, este algoritmo es para la raíz cuadrada inversa.

Tenga en cuenta que este todavía no da un resultado completamente exacto. rsqrtps con una iteración NR da casi 23 bits de precisión, frente a los 24 bits de sqrtps con el redondeo correcto para el último bit.

La precisión limitada es un problema si desea truncar el resultado a integer . (int)4.99999 es 4. Además, tenga cuidado con el caso x == 0.0 si usa sqrt(x) ~= x * sqrt(x), porque 0 * +Inf = NaN.

 35
Author: Aki Suihkonen,
Warning: date(): Invalid date.timezone value 'Europe/Kyiv', we selected the timezone 'UTC' for now. in /var/www/agent_stack/data/www/ajaxhispano.com/template/agent.layouts/content.php on line 61
2017-05-23 12:01:39

Para calcular la raíz cuadrada inversa de a, el método de Newton se aplica a la ecuación 0=f(x)=a-x^(-2) con derivada f'(x)=2*x^(-3) y por lo tanto el paso de iteración 

N(x) = x - f(x)/f'(x) = x - (a*x^3-x)/2 
     = x/2 * (3 - a*x^2)

Este método sin división tiene contrast en contraste con el método de Heron que converge globalmente a una región limitada de convergencia, por lo que necesita una buena aproximación de la raíz cuadrada inversa para obtener una mejor aproximación.

 3
Author: LutzL,
Warning: date(): Invalid date.timezone value 'Europe/Kyiv', we selected the timezone 'UTC' for now. in /var/www/agent_stack/data/www/ajaxhispano.com/template/agent.layouts/content.php on line 61
2017-06-04 13:53:23