Templatize doubles

This commit is contained in:
Alejandro Gallo 2022-01-27 20:45:38 +01:00
parent 9d684b6624
commit c7c6db77dc

View File

@ -1608,11 +1608,11 @@ namespace atrip {
( size_t No ( size_t No
, size_t Nv , size_t Nv
, const ABCTuple &abc , const ABCTuple &abc
, double const* Tph , F const* Tph
, double const* VABij , F const* VABij
, double const* VACij , F const* VACij
, double const* VBCij , F const* VBCij
, double *Zijk , F *Zijk
) { ) {
const size_t a(abc[0]), b(abc[1]), c(abc[2]); const size_t a(abc[0]), b(abc[1]), c(abc[2]);
for (size_t k=0; k < No; k++) for (size_t k=0; k < No; k++)
@ -1627,31 +1627,32 @@ namespace atrip {
} }
} }
template <typename F=double>
void doublesContribution void doublesContribution
( const ABCTuple &abc ( const ABCTuple &abc
, size_t const No , size_t const No
, size_t const Nv , size_t const Nv
// -- VABCI // -- VABCI
, double const* VABph , F const* VABph
, double const* VACph , F const* VACph
, double const* VBCph , F const* VBCph
, double const* VBAph , F const* VBAph
, double const* VCAph , F const* VCAph
, double const* VCBph , F const* VCBph
// -- VHHHA // -- VHHHA
, double const* VhhhA , F const* VhhhA
, double const* VhhhB , F const* VhhhB
, double const* VhhhC , F const* VhhhC
// -- TA // -- TA
, double const* TAphh , F const* TAphh
, double const* TBphh , F const* TBphh
, double const* TCphh , F const* TCphh
// -- TABIJ // -- TABIJ
, double const* TABhh , F const* TABhh
, double const* TAChh , F const* TAChh
, double const* TBChh , F const* TBChh
// -- TIJK // -- TIJK
, double *Tijk , F *Tijk
, atrip::Timings& chrono , atrip::Timings& chrono
) { ) {
@ -1671,7 +1672,7 @@ namespace atrip {
} \ } \
t_reorder.stop(); t_reorder.stop();
#define DGEMM_PARTICLES(__A, __B) \ #define DGEMM_PARTICLES(__A, __B) \
atrip::dgemm_( "T" \ atrip::xgemm<F>( "T" \
, "N" \ , "N" \
, (int const*)&NoNo \ , (int const*)&NoNo \
, (int const*)&No \ , (int const*)&No \
@ -1686,7 +1687,7 @@ namespace atrip {
, (int const*)&NoNo \ , (int const*)&NoNo \
); );
#define DGEMM_HOLES(__A, __B, __TRANSB) \ #define DGEMM_HOLES(__A, __B, __TRANSB) \
atrip::dgemm_( "N" \ atrip::xgemm<F>( "N" \
, __TRANSB \ , __TRANSB \
, (int const*)&NoNo \ , (int const*)&NoNo \
, (int const*)&No \ , (int const*)&No \
@ -1700,10 +1701,17 @@ namespace atrip {
, _t_buffer.data() \ , _t_buffer.data() \
, (int const*)&NoNo \ , (int const*)&NoNo \
); );
#define MAYBE_CONJ(_conj, _buffer) \
if (traits::isComplex<F>()) { \
for (size_t __i = 0; __i < NoNoNo; ++__i) \
_conj[__i] = std::conj(_buffer[__i]); \
} else { \
for (size_t __i = 0; __i < NoNoNo; ++__i) \
_conj[__i] = _buffer[__i]; \
}
using F = double;
const size_t NoNoNo = No*NoNo; const size_t NoNoNo = No*NoNo;
std::vector<double> _t_buffer; std::vector<F> _t_buffer;
_t_buffer.reserve(NoNoNo); _t_buffer.reserve(NoNoNo);
F one{1.0}, m_one{-1.0}, zero{0.0}; F one{1.0}, m_one{-1.0}, zero{0.0};
@ -1716,38 +1724,48 @@ namespace atrip {
chrono["doubles:holes"].start(); chrono["doubles:holes"].start();
{ // Holes part ============================================================ { // Holes part ============================================================
std::vector<F> _vhhh(NoNoNo);
// VhhhC[i + k*No + L*NoNo] * TABhh[L + j*No]; H1 // VhhhC[i + k*No + L*NoNo] * TABhh[L + j*No]; H1
MAYBE_CONJ(_vhhh, VhhhC)
chrono["doubles:holes:1"].start(); chrono["doubles:holes:1"].start();
DGEMM_HOLES(VhhhC, TABhh, "N") DGEMM_HOLES(_vhhh.data(), TABhh, "N")
REORDER(i, k, j) REORDER(i, k, j)
chrono["doubles:holes:1"].stop(); chrono["doubles:holes:1"].stop();
// VhhhC[j + k*No + L*NoNo] * TABhh[i + L*No]; H0 // VhhhC[j + k*No + L*NoNo] * TABhh[i + L*No]; H0
chrono["doubles:holes:2"].start(); chrono["doubles:holes:2"].start();
DGEMM_HOLES(VhhhC, TABhh, "T") DGEMM_HOLES(_vhhh.data(), TABhh, "T")
REORDER(j, k, i) REORDER(j, k, i)
chrono["doubles:holes:2"].stop(); chrono["doubles:holes:2"].stop();
// VhhhB[i + j*No + L*NoNo] * TAChh[L + k*No]; H5 // VhhhB[i + j*No + L*NoNo] * TAChh[L + k*No]; H5
MAYBE_CONJ(_vhhh, VhhhB)
chrono["doubles:holes:3"].start(); chrono["doubles:holes:3"].start();
DGEMM_HOLES(VhhhB, TAChh, "N") DGEMM_HOLES(_vhhh.data(), TAChh, "N")
REORDER(i, j, k) REORDER(i, j, k)
chrono["doubles:holes:3"].stop(); chrono["doubles:holes:3"].stop();
// VhhhB[k + j*No + L*NoNo] * TAChh[i + L*No]; H3 // VhhhB[k + j*No + L*NoNo] * TAChh[i + L*No]; H3
chrono["doubles:holes:4"].start(); chrono["doubles:holes:4"].start();
DGEMM_HOLES(VhhhB, TAChh, "T") DGEMM_HOLES(_vhhh.data(), TAChh, "T")
REORDER(k, j, i) REORDER(k, j, i)
chrono["doubles:holes:4"].stop(); chrono["doubles:holes:4"].stop();
// VhhhA[j + i*No + L*NoNo] * TBChh[L + k*No]; H1 // VhhhA[j + i*No + L*NoNo] * TBChh[L + k*No]; H1
MAYBE_CONJ(_vhhh, VhhhA)
chrono["doubles:holes:5"].start(); chrono["doubles:holes:5"].start();
DGEMM_HOLES(VhhhA, TBChh, "N") DGEMM_HOLES(_vhhh.data(), TBChh, "N")
REORDER(j, i, k) REORDER(j, i, k)
chrono["doubles:holes:5"].stop(); chrono["doubles:holes:5"].stop();
// VhhhA[k + i*No + L*NoNo] * TBChh[j + L*No]; H4 // VhhhA[k + i*No + L*NoNo] * TBChh[j + L*No]; H4
chrono["doubles:holes:6"].start(); chrono["doubles:holes:6"].start();
DGEMM_HOLES(VhhhA, TBChh, "T") DGEMM_HOLES(_vhhh.data(), TBChh, "T")
REORDER(k, i, j) REORDER(k, i, j)
chrono["doubles:holes:6"].stop(); chrono["doubles:holes:6"].stop();
} }
chrono["doubles:holes"].stop(); chrono["doubles:holes"].stop();
#undef MAYBE_CONJ
chrono["doubles:particles"].start(); chrono["doubles:particles"].start();
{ // Particle part ========================================================= { // Particle part =========================================================