mirror of
https://github.com/paboyle/Grid.git
synced 2026-08-24 19:39:35 +01:00
Accumulate in sobjD throughout rather than accumulating in sobj and converting the final sum. For float fields this matters: summing N floats then casting loses O(N*eps_float) relative precision vs accumulating in double from the start. Co-Authored-By: Claude Sonnet 4.6 <noreply@anthropic.com>
111 lines
3.3 KiB
C++
111 lines
3.3 KiB
C++
NAMESPACE_BEGIN(Grid);
|
|
|
|
/////////////////////////////////////////////////////////////////////////////////////////////////////////
|
|
// Possibly promote to double and sum
|
|
/////////////////////////////////////////////////////////////////////////////////////////////////////////
|
|
|
|
|
|
template <class vobj>
|
|
inline typename vobj::scalar_objectD sumD_gpu_tensor(const vobj *lat, Integer osites)
|
|
{
|
|
typedef typename vobj::scalar_object sobj;
|
|
typedef typename vobj::scalar_objectD sobjD;
|
|
|
|
sobjD identity; zeroit(identity);
|
|
sobjD ret; zeroit(ret);
|
|
{
|
|
sycl::buffer<sobjD, 1> abuff(&ret, {1});
|
|
theGridAccelerator->submit([&](sycl::handler &cgh) {
|
|
auto Reduction = sycl::reduction(abuff, cgh, identity, std::plus<>());
|
|
cgh.parallel_for(sycl::range<1>{(size_t)osites},
|
|
Reduction,
|
|
[=](sycl::id<1> item, auto &sum) {
|
|
sobj s = Reduce(lat[item[0]]);
|
|
sobjD sd; sd = s;
|
|
sum += sd;
|
|
});
|
|
});
|
|
}
|
|
return ret;
|
|
}
|
|
|
|
template <class vobj>
|
|
inline typename vobj::scalar_objectD sumD_gpu_large(const vobj *lat, Integer osites)
|
|
{
|
|
return sumD_gpu_tensor(lat,osites);
|
|
}
|
|
template <class vobj>
|
|
inline typename vobj::scalar_objectD sumD_gpu_small(const vobj *lat, Integer osites)
|
|
{
|
|
return sumD_gpu_large(lat,osites);
|
|
}
|
|
|
|
template <class vobj>
|
|
inline typename vobj::scalar_objectD sumD_gpu(const vobj *lat, Integer osites)
|
|
{
|
|
return sumD_gpu_large(lat,osites);
|
|
}
|
|
|
|
/////////////////////////////////////////////////////////////////////////////////////////////////////////
|
|
// Return as same precision as input performing reduction in double precision though
|
|
/////////////////////////////////////////////////////////////////////////////////////////////////////////
|
|
template <class vobj>
|
|
inline typename vobj::scalar_object sum_gpu(const vobj *lat, Integer osites)
|
|
{
|
|
typedef typename vobj::scalar_object sobj;
|
|
sobj result;
|
|
result = sumD_gpu(lat,osites);
|
|
return result;
|
|
}
|
|
|
|
template <class vobj>
|
|
inline typename vobj::scalar_object sum_gpu_large(const vobj *lat, Integer osites)
|
|
{
|
|
typedef typename vobj::scalar_object sobj;
|
|
sobj result;
|
|
result = sumD_gpu_large(lat,osites);
|
|
return result;
|
|
}
|
|
|
|
|
|
template<class Word> Word svm_xor(Word *vec,uint64_t L)
|
|
{
|
|
Word identity; identity=0;
|
|
Word ret = 0;
|
|
{
|
|
sycl::buffer<Word, 1> abuff(&ret, {1});
|
|
theGridAccelerator->submit([&](sycl::handler &cgh) {
|
|
auto Reduction = sycl::reduction(abuff,cgh,identity,std::bit_xor<>());
|
|
cgh.parallel_for(sycl::range<1>{L},
|
|
Reduction,
|
|
[=] (sycl::id<1> index, auto &sum) {
|
|
sum ^=vec[index];
|
|
});
|
|
});
|
|
}
|
|
theGridAccelerator->wait();
|
|
return ret;
|
|
}
|
|
template<class Word> Word checksum_gpu(Word *vec,uint64_t L)
|
|
{
|
|
Word identity; identity=0;
|
|
Word ret = 0;
|
|
{
|
|
sycl::buffer<Word, 1> abuff(&ret, {1});
|
|
theGridAccelerator->submit([&](sycl::handler &cgh) {
|
|
auto Reduction = sycl::reduction(abuff,cgh,identity,std::bit_xor<>());
|
|
cgh.parallel_for(sycl::range<1>{L},
|
|
Reduction,
|
|
[=] (sycl::id<1> index, auto &sum) {
|
|
auto l = index % 61;
|
|
sum ^= vec[index]<<l | vec[index]>>(64-l);
|
|
});
|
|
});
|
|
}
|
|
theGridAccelerator->wait();
|
|
return ret;
|
|
}
|
|
|
|
NAMESPACE_END(Grid);
|
|
|