Comment normaliser 4 vecteurs rapidement via les instructions SIMD

Une grande partie des performances dégagées par notre moteur sont dues au concept de SIMD; Simple Instruction Multiple Data. Ce sont des instructions processeurs (SSE, SSE2, SSE3, SSE4a, 3dNow!) qui permettent d'accélérer des calculs bas niveau.


L'exemple suivant montre comment faire une addition de vecteur classique

 struct vector  
 {  
   float x,y,z;  
 };  
 vector add(const vector &a, const vector &b)  
 {  
   vector c;  
   c.x = a.x + b.x;  
   c.y = a.y + b.y;  
   c.z = a.z + b.z;  
   return c;  
 }  

Et maintenant la version SIMD:

 # include <xmmintrin.h>  
 struct vector  
 {  
   __m128 m;  
   // Pour acceder aux composantes x, y, z  
   // vous pouvez simplement faire un cast ((float*)&m)[ 0 1 2 ou 3]  
 };  
 vector add(const vector &a, const vector &b)  
 {  
   vector c;  
   c.m = _mm_add_ps(a.m, b.m);  
   return c;  
 }  

Entre les deux versions, le facteur de speedup se situe entre 3 et 4. Si vous développez des applications nécessitants de nombreux calculs 3D, utilisez absolument ces instructions!


Voici maintenant une astuce pour gagner du temps sur les normalisations.


Il est plus performant de multiplier par la sqrt inverse (_mm_rsqrt_ps) que de diviser par la sqrt. En revanche la fonction _mm_rsqrt_ps a une précision moindre. Pour combler ce manque de précision, on fait une nouvelle itération dans le calcul de la racine, appelée "Iteration de Newton-Raphson".

 __m128    fastrsqrt( const __m128 v )  
 {  
       static const __m128 _half4 = _mm_set_ps1( 0.5f );  
       static const __m128 _three = _mm_set_ps1( 3.0f );  
       const __m128 nr = _mm_rsqrt_ps( v ); // Newton-Raphson iteration  
       const __m128 muls = _mm_mul_ps( _mm_mul_ps( v, nr ), nr );  
       return _mm_mul_ps( _mm_mul_ps( _half4, nr ), _mm_sub_ps( _three, muls ) );  
 }  
 void Normalize4(vector p[4])  
 {  
     // Store x, y, and z separately in 3 __m128  
     __m128 x4 = _mm_set_ps(p[3][0], p[2][0], p[1][0], p[0][0]);  
     __m128 y4 = _mm_set_ps(p[3][1], p[2][1], p[1][1], p[0][1]);  
     __m128 z4 = _mm_set_ps(p[3][2], p[2][2], p[1][2], p[0][2]);  
     x4 = _mm_mul_ps(x4, x4); // x * x  
     y4 = _mm_mul_ps(y4, y4); // y * y  
     z4 = _mm_mul_ps(z4, z4); // z * z  
     const vector l = fastrsqrt(_mm_add_ps(z4, _mm_add_ps(y4, x4))); // sqrt(x*x+y*y+z*z)  
     p[0].m = _mm_mul_ps(p[0].m, _mm_set_ps1(l[0])); // load results in packet  
     p[1].m = _mm_mul_ps(p[1].m, _mm_set_ps1(l[1]));  
     p[2].m = _mm_mul_ps(p[2].m, _mm_set_ps1(l[2]));  
     p[3].m = _mm_mul_ps(p[3].m, _mm_set_ps1(l[3]));  
 }  

Si vous avez aimé l'article n'hésitez pas à laisser un commentaire, ou me suivre sur twitter

Rencontre avec AMD

La semaine dernière, nous avons eu le plaisir de rencontrer à l'ACSEL un employé de chez AMD! Un entretien très intéressant autour de l'optimisation matérielle, et de l'utilisation du trio CPU / GPU / APU pour lancer des rayons.

Les 5 points à retenir absolument avant de se lancer dans le GPGPU:

1. Un GPU est très bon pour les calculs parallèles sur de gros blocs de données, si toutefois les calculs ne diffèrent pas trop. Ce qui veut dire qu'il est très bon pour les rayons primaires, beaucoup moins pour les secondaires.
2. Les algorithmes avec beaucoup de branchements (type parcours d'arbre) sont à BANNIR sur GPU. Le CPU est beaucoup plus à l'aise dessus.
3. Un APU à environ 1/4 de la puissance d'un GPU. Mais du fait de sa localisation, il sera plus performant pour certains calculs que le GPU.
4. La taille de la VRAM ne doit pas être un soucis. Des solutions arriveront.
5. Le temps de transfert ne doit pas être un soucis. Des solutions arriveront.

Nous voilà donc partis dans nos premiers TPs sur OpenCL. Les premiers résultats sont concluants, mais il faudra un gros travail de conception avant de voir le GPGPU d'imposer dans le moteur.

Si vous avez aimé l'article n'hésitez pas à laisser un commentaire, ou me suivre sur twitter