6.1.1.4.2 : Le fichier source

Écrivons le fichier sgemm_intrinsics.cpp :



Il nous faut tout d'abord inclure un header qui permet de définir l'alignement des données en float avec PLIB_VECTOR_SIZE_FLOAT. Ainsi que des fonctions intrisèques génériques telles que :
  • plib_load_ps : pour lire des registres vectoriels (par exemple 8 floats en AVX
  • plib_mul_ps : pour multiplier des registres vectoriels
  • plib_add_ps : pour additionner des registres vectoriel
  • plib_store_ps : pour sauvegarder des registres vectoriels 


1

#include "phoenix_intrinsics.h"


Et string.h pour avoir la fonction memset :

1

#include <string.h>


Ensuite nous incluons notre header :

1

#include "sgemm_intrinsics.h"


Et nous implémentons le produit de matrices :

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30

///Compute the Matrix-Matrix product of the x,y matrices
/**	@param[out] matOut : result
 * 	@param matX : left matrix
 * 	@param matY : right matrix
 * 	@param size : size of the square matrices
*/
void sgemm_intrinsics(float* matOut, const float* matX, const float* matY, long unsigned int size){
	memset(matOut, 0, sizeof(float)*size*size);
	long unsigned int vecSize(PLIB_VECTOR_SIZE_FLOAT);
	long unsigned int nbVec(size/vecSize);
	
	for(long unsigned int i(0lu); i < size; ++i){			//Iterate over X rows
		const float * rowX = matX + i*size;				//Get current X row
		float * rowOut = matOut + i*size;					//Get current out row
		for(long unsigned int k(0lu); k < size; ++k){					//Part of dot product
			PRegVecf regX = plib_broadcast_ss(rowX[k]);
			
			const float* rowY = matY + k*size;				//Get current Y row
			for(long unsigned int j(0lu); j < nbVec; ++j){
				PRegVecf regY = plib_load_ps(rowY + vecSize*j);
				
				PRegVecf regRes = plib_load_ps(rowOut + vecSize*j);
				
				regRes = plib_add_ps(regRes, plib_mul_ps(regX, regY));
				
				plib_store_ps(rowOut + vecSize*j, regRes);
			}
		}
	}
}


Le fichier sgemm_intrinsics.cpp complet :

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39

/***************************************
	Auteur : Pierre Aubert
	Mail : pierre.aubert@lapp.in2p3.fr
	Licence : CeCILL-C
****************************************/

#include "phoenix_intrinsics.h"
#include <string.h>
#include "sgemm_intrinsics.h"
///Compute the Matrix-Matrix product of the x,y matrices
/**	@param[out] matOut : result
 * 	@param matX : left matrix
 * 	@param matY : right matrix
 * 	@param size : size of the square matrices
*/
void sgemm_intrinsics(float* matOut, const float* matX, const float* matY, long unsigned int size){
	memset(matOut, 0, sizeof(float)*size*size);
	long unsigned int vecSize(PLIB_VECTOR_SIZE_FLOAT);
	long unsigned int nbVec(size/vecSize);
	
	for(long unsigned int i(0lu); i < size; ++i){			//Iterate over X rows
		const float * rowX = matX + i*size;				//Get current X row
		float * rowOut = matOut + i*size;					//Get current out row
		for(long unsigned int k(0lu); k < size; ++k){					//Part of dot product
			PRegVecf regX = plib_broadcast_ss(rowX[k]);
			
			const float* rowY = matY + k*size;				//Get current Y row
			for(long unsigned int j(0lu); j < nbVec; ++j){
				PRegVecf regY = plib_load_ps(rowY + vecSize*j);
				
				PRegVecf regRes = plib_load_ps(rowOut + vecSize*j);
				
				regRes = plib_add_ps(regRes, plib_mul_ps(regX, regY));
				
				plib_store_ps(rowOut + vecSize*j, regRes);
			}
		}
	}
}


Vous pouvez le télécharger ici.