portelli/Grid
1
0
Fork 0
mirror of https://github.com/paboyle/Grid.git synced 2025-10-28 18:49:33 +00:00
Code Releases Activity

Works now with Clang-avx, Clang-sse and ICPC-avx, ICPC-sse

Browse Source
This commit is contained in:
Peter Boyle
2015-05-28 11:35:43 +01:00
parent bd81ac0f17
commit c0c1ebe757
3 changed files with 186 additions and 207 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
lib/Makefile.am
View File

@@ -96,7 +96,7 @@ nobase_include_HEADERS = algorithms/approx/bigfloat.h \
simd/Grid_vector_types.h \
simd/Grid_sse4.h \
simd/Grid_avx.h \
simd/Grid_knc.h
simd/Grid_avx512.h

0
lib/simd/Grid_knc.h → lib/simd/Grid_avx512.h
View File

391
lib/simd/Grid_vector_types.h
View File

@@ -14,7 +14,7 @@
#include "Grid_avx.h"
#endif
#if defined AVX512
#include "Grid_knc.h"
#include "Grid_avx512.h"
#endif
#if defined QPX
#include "Grid_qpx.h"
@@ -32,13 +32,24 @@ namespace Grid {
// type alias used to simplify the syntax of std::enable_if
template <typename T> using Invoke = typename T::type;
template <typename Condition, typename ReturnType> using EnableIf = Invoke<std::enable_if<Condition::value, ReturnType>>;
template <typename Condition, typename ReturnType> using NotEnableIf= Invoke<std::enable_if<!Condition::value, ReturnType>>;
template <typename Condition, typename ReturnType> using EnableIf = Invoke<std::enable_if<Condition::value, ReturnType> >;
template <typename Condition, typename ReturnType> using NotEnableIf= Invoke<std::enable_if<!Condition::value, ReturnType> >;
////////////////////////////////////////////////////////
// Check for complexity with type traits
template <typename T> struct is_complex : std::false_type {};
template < typename T > struct is_complex< std::complex<T> >: std::true_type {};
template <typename T> struct is_complex : public std::false_type {};
template <> struct is_complex<std::complex<double> >: public std::true_type {};
template <> struct is_complex<std::complex<float> > : public std::true_type {};
template <typename T> using IfReal = Invoke<std::enable_if<std::is_floating_point<T>::value,int> > ;
template <typename T> using IfComplex = Invoke<std::enable_if<is_complex<T>::value,int> > ;
template <typename T> using IfInteger = Invoke<std::enable_if<std::is_integral<T>::value,int> > ;
template <typename T> using IfNotReal = Invoke<std::enable_if<!std::is_floating_point<T>::value,int> > ;
template <typename T> using IfNotComplex = Invoke<std::enable_if<!is_complex<T>::value,int> > ;
template <typename T> using IfNotInteger = Invoke<std::enable_if<!std::is_integral<T>::value,int> > ;
////////////////////////////////////////////////////////
// Define the operation templates functors
// general forms to allow for vsplat syntax
@@ -54,11 +65,9 @@ namespace Grid {
Out unary(Input src, Operation op){
return op(src);
}
///////////////////////////////////////////////
/*
@brief Grid_simd class for the SIMD vector type operations
@@ -73,27 +82,28 @@ namespace Grid {
Vector_type v;
static inline int Nsimd(void) { return sizeof(Vector_type)/sizeof(Scalar_type);}
// Constructors
Grid_simd & operator = ( Zero & z){
vzero(*this);
return (*this);
}
Grid_simd& operator=(const Grid_simd&& rhs){v=rhs.v;return *this;};
Grid_simd& operator=(const Grid_simd& rhs){v=rhs.v;return *this;}; //faster than not declaring it and leaving to the compiler
Grid_simd()=default;
Grid_simd(const Grid_simd& rhs):v(rhs.v){}; //compiles in movaps
Grid_simd(const Grid_simd& rhs) :v(rhs.v){}; //compiles in movaps
Grid_simd(const Grid_simd&& rhs):v(rhs.v){};
/////////////////////////////
// Constructors
/////////////////////////////
Grid_simd & operator = ( Zero & z){
vzero(*this);
return (*this);
}
//Enable if complex type
template < class S = Scalar_type >
template < typename S = Scalar_type >
Grid_simd(const typename std::enable_if< is_complex < S >::value, S>::type a){
vsplat(*this,a);
};
Grid_simd(const Real a){
vsplat(*this,Scalar_type(a));
@@ -107,86 +117,16 @@ namespace Grid {
friend inline void sub (Grid_simd * __restrict__ y,const Grid_simd * __restrict__ l,const Grid_simd *__restrict__ r){ *y = (*l) - (*r); }
friend inline void add (Grid_simd * __restrict__ y,const Grid_simd * __restrict__ l,const Grid_simd *__restrict__ r){ *y = (*l) + (*r); }
friend inline void mac (Grid_simd *__restrict__ y,const Scalar_type *__restrict__ a,const Grid_simd *__restrict__ x){ *y = (*a)*(*x)+(*y); };
friend inline void mult(Grid_simd *__restrict__ y,const Scalar_type *__restrict__ l,const Grid_simd *__restrict__ r){ *y = (*l) * (*r); }
friend inline void sub (Grid_simd *__restrict__ y,const Scalar_type *__restrict__ l,const Grid_simd *__restrict__ r){ *y = (*l) - (*r); }
friend inline void add (Grid_simd *__restrict__ y,const Scalar_type *__restrict__ l,const Grid_simd *__restrict__ r){ *y = (*l) + (*r); }
friend inline void mac (Grid_simd *__restrict__ y,const Grid_simd *__restrict__ a,const Scalar_type *__restrict__ x){ *y = (*a)*(*x)+(*y); };
friend inline void mult(Grid_simd *__restrict__ y,const Grid_simd *__restrict__ l,const Scalar_type *__restrict__ r){ *y = (*l) * (*r); }
friend inline void sub (Grid_simd *__restrict__ y,const Grid_simd *__restrict__ l,const Scalar_type *__restrict__ r){ *y = (*l) - (*r); }
friend inline void add (Grid_simd *__restrict__ y,const Grid_simd *__restrict__ l,const Scalar_type *__restrict__ r){ *y = (*l) + (*r); }
//not for integer types...
template < class S = Scalar_type, NotEnableIf<std::is_integral < S >, int> = 0 >
friend inline Grid_simd adj(const Grid_simd &in){ return conjugate(in); }
///////////////////////////////////////////////
// Initialise to 1,0,i for the correct types
///////////////////////////////////////////////
// For complex types
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline void vone(Grid_simd &ret) { vsplat(ret,1.0,0.0); }
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline void vzero(Grid_simd &ret) { vsplat(ret,0.0,0.0); }// use xor?
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline void vcomplex_i(Grid_simd &ret){ vsplat(ret,0.0,1.0);}
// if not complex overload here
template < class S = Scalar_type, EnableIf<std::is_floating_point < S >,int> = 0 >
friend inline void vone(Grid_simd &ret) { vsplat(ret,1.0); }
template < class S = Scalar_type, EnableIf<std::is_floating_point < S >,int> = 0 >
friend inline void vzero(Grid_simd &ret) { vsplat(ret,0.0); }
// For integral types
template < class S = Scalar_type, EnableIf<std::is_integral < S >, int> = 0 >
friend inline void vone(Grid_simd &ret) { vsplat(ret,1); }
template < class S = Scalar_type, EnableIf<std::is_integral < S >, int> = 0 >
friend inline void vzero(Grid_simd &ret) { vsplat(ret,0); }
template < class S = Scalar_type, EnableIf<std::is_integral < S >, int> = 0 >
friend inline void vtrue (Grid_simd &ret){vsplat(ret,0xFFFFFFFF);}
template < class S = Scalar_type, EnableIf<std::is_integral < S >, int> = 0 >
friend inline void vfalse(Grid_simd &ret){vsplat(ret,0);}
////////////////////////////////////
// Arithmetic operator overloads +,-,*
////////////////////////////////////
friend inline Grid_simd operator + (Grid_simd a, Grid_simd b)
{
Grid_simd ret;
ret.v = binary<Vector_type>(a.v, b.v, SumSIMD());
return ret;
};
friend inline Grid_simd operator - (Grid_simd a, Grid_simd b)
{
Grid_simd ret;
ret.v = binary<Vector_type>(a.v, b.v, SubSIMD());
return ret;
};
// Distinguish between complex types and others
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline Grid_simd operator * (Grid_simd a, Grid_simd b)
{
Grid_simd ret;
ret.v = binary<Vector_type>(a.v,b.v, MultComplexSIMD());
return ret;
};
// Real/Integer types
template < class S = Scalar_type, NotEnableIf<is_complex < S >, int> = 0 >
friend inline Grid_simd operator * (Grid_simd a, Grid_simd b)
{
Grid_simd ret;
ret.v = binary<Vector_type>(a.v,b.v, MultSIMD());
return ret;
};
////////////////////////////////////////////////////////////////////////
// FIXME: gonna remove these load/store, get, set, prefetch
////////////////////////////////////////////////////////////////////////
@@ -194,28 +134,6 @@ namespace Grid {
ret.v = unary<Vector_type>(a, VsetSIMD());
}
///////////////////////
// Splat
///////////////////////
// overload if complex
template < class S = Scalar_type >
friend inline void vsplat(Grid_simd &ret, EnableIf<is_complex < S >, S> c){
Real a = real(c);
Real b = imag(c);
vsplat(ret,a,b);
}
// this is only for the complex version
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline void vsplat(Grid_simd &ret,Real a, Real b){
ret.v = binary<Vector_type>(a, b, VsplatSIMD());
}
//if real fill with a, if complex fill with a in the real part (first function above)
friend inline void vsplat(Grid_simd &ret,Real a){
ret.v = unary<Vector_type>(a, VsplatSIMD());
}
///////////////////////
// Vstore
///////////////////////
@@ -223,19 +141,6 @@ namespace Grid {
binary<void>(ret.v, (Real*)a, VstoreSIMD());
}
///////////////////////
// Vstream
///////////////////////
template < class S = Scalar_type, NotEnableIf<std::is_integral < S >, int> = 0 >
friend inline void vstream(Grid_simd &out,const Grid_simd &in){
binary<void>((Real*)&out.v, in.v, VstreamSIMD());
}
template < class S = Scalar_type, EnableIf<std::is_integral < S >, int> = 0 >
friend inline void vstream(Grid_simd &out,const Grid_simd &in){
out=in;
}
///////////////////////
// Vprefetch
///////////////////////
@@ -244,7 +149,6 @@ namespace Grid {
_mm_prefetch((const char*)&v.v,_MM_HINT_T0);
}
///////////////////////
// Reduce
///////////////////////
@@ -265,63 +169,6 @@ namespace Grid {
return a*b;
}
///////////////////////
// Conjugate
///////////////////////
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline Grid_simd conjugate(const Grid_simd &in){
Grid_simd ret ;
ret.v = unary<Vector_type>(in.v, ConjSIMD());
return ret;
}
template < class S = Scalar_type, NotEnableIf<is_complex < S >, int> = 0 >
friend inline Grid_simd conjugate(const Grid_simd &in){
return in; // for real objects
}
///////////////////////
// timesMinusI
///////////////////////
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline void timesMinusI( Grid_simd &ret,const Grid_simd &in){
ret.v = binary<Vector_type>(in.v, ret.v, TimesMinusISIMD());
}
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline Grid_simd timesMinusI(const Grid_simd &in){
Grid_simd ret;
timesMinusI(ret,in);
return ret;
}
template < class S = Scalar_type, NotEnableIf<is_complex < S >, int> = 0 >
friend inline Grid_simd timesMinusI(const Grid_simd &in){
return in;
}
///////////////////////
// timesI
///////////////////////
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline void timesI(Grid_simd &ret,const Grid_simd &in){
ret.v = binary<Vector_type>(in.v, ret.v, TimesISIMD());
}
template < class S = Scalar_type, EnableIf<is_complex < S >, int> = 0 >
friend inline Grid_simd timesI(const Grid_simd &in){
Grid_simd ret;
timesI(ret,in);
return ret;
}
template < class S = Scalar_type, NotEnableIf<is_complex < S >, int> = 0 >
friend inline Grid_simd timesI(const Grid_simd &in){
return in;
}
///////////////////////
// Unary negation
///////////////////////
@@ -346,9 +193,6 @@ namespace Grid {
return *this;
}
////////////////////////////////////////////////////////////////////
// General permute; assumes vector length is same across
// all subtypes; may not be a good assumption, but could
@@ -359,48 +203,183 @@ namespace Grid {
Gpermute<Grid_simd>(y,b,perm);
}
};// end of Grid_simd class definition
///////////////////////
// Splat
///////////////////////
// this is only for the complex version
template <class S, class V, IfComplex<S> =0, class ABtype>
inline void vsplat(Grid_simd<S,V> &ret,ABtype a, ABtype b){
ret.v = binary<V>(a, b, VsplatSIMD());
}
template<class scalar_type, class vector_type >
inline Grid_simd< scalar_type, vector_type> innerProduct(const Grid_simd< scalar_type, vector_type> & l, const Grid_simd< scalar_type, vector_type> & r)
// overload if complex
template <class S,class V> inline void vsplat(Grid_simd<S,V> &ret, EnableIf<is_complex < S >, S> c) {
Real a = real(c);
Real b = imag(c);
vsplat(ret,a,b);
}
//if real fill with a, if complex fill with a in the real part (first function above)
template <class S,class V>
inline void vsplat(Grid_simd<S,V> &ret,NotEnableIf<is_complex< S>,S> a){
ret.v = unary<V>(a, VsplatSIMD());
}
//////////////////////////
///////////////////////////////////////////////
// Initialise to 1,0,i for the correct types
///////////////////////////////////////////////
// For complex types
template <class S,class V, IfComplex<S> = 0 > inline void vone(Grid_simd<S,V> &ret) { vsplat(ret,S(1.0,0.0)); }
template <class S,class V, IfComplex<S> = 0 > inline void vzero(Grid_simd<S,V> &ret) { vsplat(ret,S(0.0,0.0)); }// use xor?
template <class S,class V, IfComplex<S> = 0 > inline void vcomplex_i(Grid_simd<S,V> &ret){ vsplat(ret,S(0.0,1.0));}
// if not complex overload here
template <class S,class V, IfReal<S> = 0 > inline void vone (Grid_simd<S,V> &ret){ vsplat(ret,1.0); }
template <class S,class V, IfReal<S> = 0 > inline void vzero(Grid_simd<S,V> &ret) { vsplat(ret,0.0); }
// For integral types
template <class S,class V,IfInteger<S> = 0 > inline void vone(Grid_simd<S,V> &ret) {vsplat(ret,1); }
template <class S,class V,IfInteger<S> = 0 > inline void vzero(Grid_simd<S,V> &ret) {vsplat(ret,0); }
template <class S,class V,IfInteger<S> = 0 > inline void vtrue (Grid_simd<S,V> &ret){vsplat(ret,0xFFFFFFFF);}
template <class S,class V,IfInteger<S> = 0 > inline void vfalse(Grid_simd<S,V> &ret){vsplat(ret,0);}
template<class S,class V> inline void zeroit(Grid_simd<S,V> &z){ vzero(z);}
///////////////////////
// Vstream
///////////////////////
template <class S,class V, IfNotInteger<S> = 0 >
inline void vstream(Grid_simd<S,V> &out,const Grid_simd<S,V> &in){
binary<void>((Real*)&out.v, in.v, VstreamSIMD());
}
template <class S,class V, IfInteger<S> = 0 >
inline void vstream(Grid_simd<S,V> &out,const Grid_simd<S,V> &in){
out=in;
}
////////////////////////////////////
// Arithmetic operator overloads +,-,*
////////////////////////////////////
template<class S,class V> inline Grid_simd<S,V> operator + (Grid_simd<S,V> a, Grid_simd<S,V> b) {
Grid_simd<S,V> ret;
ret.v = binary<V>(a.v, b.v, SumSIMD());
return ret;
};
template<class S,class V> inline Grid_simd<S,V> operator - (Grid_simd<S,V> a, Grid_simd<S,V> b) {
Grid_simd<S,V> ret;
ret.v = binary<V>(a.v, b.v, SubSIMD());
return ret;
};
// Distinguish between complex types and others
template<class S,class V, IfComplex<S> = 0 > inline Grid_simd<S,V> operator * (Grid_simd<S,V> a, Grid_simd<S,V> b) {
Grid_simd<S,V> ret;
ret.v = binary<V>(a.v,b.v, MultComplexSIMD());
return ret;
};
// Real/Integer types
template<class S,class V, IfNotComplex<S> = 0 > inline Grid_simd<S,V> operator * (Grid_simd<S,V> a, Grid_simd<S,V> b) {
Grid_simd<S,V> ret;
ret.v = binary<V>(a.v,b.v, MultSIMD());
return ret;
};
///////////////////////
// Conjugate
///////////////////////
template <class S,class V, IfComplex<S> = 0 >
inline Grid_simd<S,V> conjugate(const Grid_simd<S,V> &in){
Grid_simd<S,V> ret ;
ret.v = unary<V>(in.v, ConjSIMD());
return ret;
}
template <class S,class V, IfNotComplex<S> = 0 > inline Grid_simd<S,V> conjugate(const Grid_simd<S,V> &in){
return in; // for real objects
}
//Suppress adj for integer types... // odd; why conjugate above but not adj??
template < class S, class V, IfNotInteger<S> = 0 >
inline Grid_simd<S,V> adj(const Grid_simd<S,V> &in){ return conjugate(in); }
///////////////////////
// timesMinusI
///////////////////////
template<class S,class V,IfComplex<S> = 0 >
inline void timesMinusI( Grid_simd<S,V> &ret,const Grid_simd<S,V> &in){
ret.v = binary<V>(in.v, ret.v, TimesMinusISIMD());
}
template<class S,class V,IfComplex<S> = 0 >
inline Grid_simd<S,V> timesMinusI(const Grid_simd<S,V> &in){
Grid_simd<S,V> ret;
timesMinusI(ret,in);
return ret;
}
template<class S,class V,IfNotComplex<S> = 0 >
inline Grid_simd<S,V> timesMinusI(const Grid_simd<S,V> &in){
return in;
}
///////////////////////
// timesI
///////////////////////
template<class S,class V,IfComplex<S> = 0 >
inline void timesI(Grid_simd<S,V> &ret,const Grid_simd<S,V> &in){
ret.v = binary<V>(in.v, ret.v, TimesISIMD());
}
template<class S,class V,IfComplex<S> = 0 >
inline Grid_simd<S,V> timesI(const Grid_simd<S,V> &in){
Grid_simd<S,V> ret;
timesI(ret,in);
return ret;
}
template<class S,class V,IfNotComplex<S> = 0 >
inline Grid_simd<S,V> timesI(const Grid_simd<S,V> &in){
return in;
}
/////////////////////
// Inner, outer
/////////////////////
template<class S, class V >
inline Grid_simd< S, V> innerProduct(const Grid_simd< S, V> & l, const Grid_simd< S, V> & r)
{
return conjugate(l)*r;
}
template<class scalar_type, class vector_type >
inline void zeroit(Grid_simd< scalar_type, vector_type> &z){ vzero(z);}
template<class scalar_type, class vector_type >
inline Grid_simd< scalar_type, vector_type> outerProduct(const Grid_simd< scalar_type, vector_type> &l, const Grid_simd< scalar_type, vector_type>& r)
template<class S, class V >
inline Grid_simd< S, V> outerProduct(const Grid_simd< S, V> &l, const Grid_simd< S, V> & r)
{
return l*r;
}
template<class scalar_type, class vector_type >
inline Grid_simd< scalar_type, vector_type> trace(const Grid_simd< scalar_type, vector_type> &arg){
template<class S, class V >
inline Grid_simd< S, V> trace(const Grid_simd< S, V> &arg){
return arg;
}
///////////////////////////////
// Define available types
///////////////////////////////
typedef Grid_simd< float , SIMD_Ftype > vRealF;
typedef Grid_simd< double , SIMD_Dtype > vRealD;
typedef Grid_simd< std::complex< float > , SIMD_Ftype > vComplexF;
typedef Grid_simd< std::complex< double >, SIMD_Dtype > vComplexD;
typedef Grid_simd< Integer , SIMD_Itype > vInteger;
}
#endif