Compare commits
No commits in common. "933d556c8496dfb1e033876c986e4879e7073a61" and "be96e4bf8c14de0c554136cce76298f810af9703" have entirely different histories.
933d556c84
...
be96e4bf8c
@ -38,8 +38,6 @@ namespace acc {
|
||||
__MAYBE_DEVICE__ __MAYBE_HOST__ __INLINE__
|
||||
F maybeConjugateScalar(const F &a) { return a; }
|
||||
|
||||
// TODO: instantiate for std::complex<double>
|
||||
|
||||
#if defined(HAVE_CUDA)
|
||||
template <>
|
||||
__MAYBE_DEVICE__ __MAYBE_HOST__ __INLINE__
|
||||
|
||||
@ -200,7 +200,7 @@ template <typename F=double>
|
||||
: Slice<F>::Fetch
|
||||
;
|
||||
if (blank.info.state == Slice<F>::SelfSufficient) {
|
||||
#if defined(HAVE_CUDA) && !defined(ATRIP_SOURCES_IN_GPU)
|
||||
#if defined(HAVE_CUDA)
|
||||
const size_t _size = sizeof(F) * sliceSize;
|
||||
// TODO: this is code duplication with downstairs
|
||||
if (freePointers.size() == 0) {
|
||||
@ -221,6 +221,7 @@ template <typename F=double>
|
||||
(void*)SOURCES_DATA(sources[from.source]),
|
||||
sizeof(F) * sliceSize));
|
||||
))
|
||||
|
||||
#else
|
||||
blank.data = SOURCES_DATA(sources[from.source]);
|
||||
#endif
|
||||
@ -464,30 +465,30 @@ template <typename F=double>
|
||||
|
||||
#if defined(HAVE_CUDA)
|
||||
LOG(1,"Atrip") << "warming communication up " << slices.size() << "\n";
|
||||
WITH_CHRONO("cuda:warmup",
|
||||
int nRanks=Atrip::np, requestCount=0;
|
||||
int nSends=sliceBuffers.size()*nRanks;
|
||||
MPI_Request *requests = (MPI_Request*) malloc(nSends*2 * sizeof(MPI_Request));
|
||||
MPI_Status *statuses = (MPI_Status*) malloc(nSends*2 * sizeof(MPI_Status));
|
||||
for (int sliceId=0; sliceId<sliceBuffers.size(); sliceId++){
|
||||
for (int rankId=0; rankId<nRanks; rankId++){
|
||||
MPI_Isend((void*)SOURCES_DATA(sources[0]),
|
||||
sliceSize,
|
||||
traits::mpi::datatypeOf<F>(),
|
||||
rankId,
|
||||
100,
|
||||
universe,
|
||||
&requests[requestCount++]);
|
||||
MPI_Irecv((void*)sliceBuffers[sliceId],
|
||||
sliceSize,
|
||||
traits::mpi::datatypeOf<F>(),
|
||||
rankId,
|
||||
100,
|
||||
universe,
|
||||
&requests[requestCount++]);
|
||||
}
|
||||
}
|
||||
MPI_Waitall(nSends*2, requests, statuses);
|
||||
 WITH_CHRONO("cuda:warmup",
|
||||
int nRanks=Atrip::np, requestCount=0;
|
||||
int nSends=sliceBuffers.size()*nRanks;
|
||||
MPI_Request *requests = (MPI_Request*) malloc(nSends*2 * sizeof(MPI_Request));
|
||||
MPI_Status *statuses = (MPI_Status*) malloc(nSends*2 * sizeof(MPI_Status));
|
||||
for (int sliceId=0; sliceId<sliceBuffers.size(); sliceId++){
|
||||
for (int rankId=0; rankId<nRanks; rankId++){
|
||||
MPI_Isend((void*)SOURCES_DATA(sources[0]),
|
||||
sliceSize,
|
||||
traits::mpi::datatypeOf<F>(),
|
||||
rankId,
|
||||
100,
|
||||
universe,
|
||||
&requests[requestCount++]);
|
||||
MPI_Irecv((void*)sliceBuffers[sliceId],
|
||||
sliceSize,
|
||||
traits::mpi::datatypeOf<F>(),
|
||||
rankId,
|
||||
100,
|
||||
universe,
|
||||
&requests[requestCount++]);
|
||||
}
|
||||
}
|
||||
MPI_Waitall(nSends*2, requests, statuses);
|
||||
)
|
||||
#endif
|
||||
|
||||
@ -573,15 +574,12 @@ template <typename F=double>
|
||||
// TODO: do it through the slice class
|
||||
slice.info.state = Slice<F>::Dispatched;
|
||||
#if defined(HAVE_CUDA) && defined(ATRIP_SOURCES_IN_GPU)
|
||||
# if !defined(ATRIP_CUDA_AWARE_MPI)
|
||||
# if !defined(ATRIP_CUDA_AWARE_MPI)
|
||||
# error "You need CUDA aware MPI to have slices on the GPU"
|
||||
# endif
|
||||
MPI_Irecv((void*)slice.data,
|
||||
#elif defined(HAVE_CUDA) && !defined(ATRIP_SOURCES_IN_GPU)
|
||||
slice.mpi_data = (F*)malloc(sizeof(F) * slice.size);
|
||||
MPI_Irecv(slice.mpi_data,
|
||||
#else
|
||||
MPI_Irecv((void*)slice.data,
|
||||
MPI_Irecv(slice.data,
|
||||
#endif
|
||||
slice.size,
|
||||
traits::mpi::datatypeOf<F>(),
|
||||
|
||||
@ -202,7 +202,7 @@ Atrip::Output Atrip::run(Atrip::Input<F> const& in) {
|
||||
_CHECK_CUDA_SUCCESS("Zijk",
|
||||
cuMemAlloc(&Zijk, sizeof(F) * No * No * No));
|
||||
#else
|
||||
DataPtr<F> Tai = _Tai.data(), epsi = _epsi.data(), epsa = _epsa.data();
|
||||
std::vector<F> &Tai = _Tai, &epsi = _epsi, &epsa = _epsa;
|
||||
Zijk = (DataFieldType<F>*)malloc(No*No*No * sizeof(DataFieldType<F>));
|
||||
Tijk = (DataFieldType<F>*)malloc(No*No*No * sizeof(DataFieldType<F>));
|
||||
#endif
|
||||
@ -693,7 +693,7 @@ Atrip::Output Atrip::run(Atrip::Input<F> const& in) {
|
||||
(DataFieldType<F>*)Tai,
|
||||
#else
|
||||
singlesContribution<F>(No, Nv, abc[0], abc[1], abc[2],
|
||||
Tai,
|
||||
Tai.data(),
|
||||
#endif
|
||||
(DataFieldType<F>*)abhh.unwrapSlice(Slice<F>::AB,
|
||||
abc),
|
||||
@ -733,30 +733,18 @@ Atrip::Output Atrip::run(Atrip::Input<F> const& in) {
|
||||
1, 1, // for cuda
|
||||
_epsabc,
|
||||
No,
|
||||
#if defined(HAVE_CUDA)
|
||||
(DataFieldType<F>*)epsi,
|
||||
(DataFieldType<F>*)Tijk,
|
||||
(DataFieldType<F>*)Zijk,
|
||||
#else
|
||||
epsi,
|
||||
Tijk,
|
||||
Zijk,
|
||||
#endif
|
||||
tupleEnergy);
|
||||
} else {
|
||||
ACC_FUNCALL(getEnergySame<DataFieldType<F>>,
|
||||
1, 1, // for cuda
|
||||
_epsabc,
|
||||
No,
|
||||
#if defined(HAVE_CUDA)
|
||||
(DataFieldType<F>*)epsi,
|
||||
(DataFieldType<F>*)Tijk,
|
||||
(DataFieldType<F>*)Zijk,
|
||||
#else
|
||||
epsi,
|
||||
Tijk,
|
||||
Zijk,
|
||||
#endif
|
||||
tupleEnergy);
|
||||
})
|
||||
|
||||
@ -903,5 +891,5 @@ Atrip::Output Atrip::run(Atrip::Input<F> const& in) {
|
||||
}
|
||||
// instantiate
|
||||
template Atrip::Output Atrip::run(Atrip::Input<double> const& in);
|
||||
// template Atrip::Output Atrip::run(Atrip::Input<Complex> const& in);
|
||||
template Atrip::Output Atrip::run(Atrip::Input<Complex> const& in);
|
||||
// Main:1 ends here
|
||||
|
||||
@ -25,8 +25,11 @@ namespace atrip {
|
||||
|
||||
#if defined(HAVE_CUDA)
|
||||
#define FOR_K() \
|
||||
const size_t k = blockIdx.x * blockDim.x + threadIdx.x; \
|
||||
size_t idx = 0;
|
||||
for (size_t kmin = blockIdx.x * blockDim.x + threadIdx.x, \
|
||||
k = kmin, \
|
||||
idx = kmin * size * size * size; \
|
||||
k < (kmin < size) ? kmin + 1 : size; \
|
||||
k++)
|
||||
#else
|
||||
#define FOR_K() for (size_t k=0, idx=0; k < size; k++)
|
||||
#endif
|
||||
@ -99,7 +102,6 @@ namespace atrip {
|
||||
# define MIN(a, b) std::min((a), (b))
|
||||
#endif
|
||||
|
||||
#if defined(ATRIP_NEW_ENERGY)
|
||||
|
||||
// [[file:~/cuda/atrip/atrip.org::*Energy][Energy:2]]
|
||||
template <typename F>
|
||||
@ -248,131 +250,6 @@ void getEnergySame
|
||||
}
|
||||
// Energy:2 ends here
|
||||
|
||||
#else
|
||||
|
||||
// [[file:~/cuda/atrip/atrip.org::*Energy][Energy:2]]
|
||||
template <typename F>
|
||||
__MAYBE_GLOBAL__
|
||||
void getEnergyDistinct
|
||||
( F const epsabc
|
||||
, size_t const No
|
||||
, F* const epsi
|
||||
, F* const Tijk
|
||||
, F* const Zijk
|
||||
, double* _energy
|
||||
) {
|
||||
constexpr size_t blockSize=16;
|
||||
F energy(0.);
|
||||
for (size_t kk=0; kk<No; kk+=blockSize){
|
||||
const size_t kend( MIN(No, kk+blockSize) );
|
||||
for (size_t jj(kk); jj<No; jj+=blockSize){
|
||||
const size_t jend( MIN( No, jj+blockSize) );
|
||||
for (size_t ii(jj); ii<No; ii+=blockSize){
|
||||
const size_t iend( MIN( No, ii+blockSize) );
|
||||
for (size_t k(kk); k < kend; k++){
|
||||
const F ek(epsi[k]);
|
||||
const size_t jstart = jj > k ? jj : k;
|
||||
for (size_t j(jstart); j < jend; j++){
|
||||
F const ej(epsi[j]);
|
||||
F const facjk = j == k ? F(0.5) : F(1.0);
|
||||
size_t istart = ii > j ? ii : j;
|
||||
for (size_t i(istart); i < iend; i++){
|
||||
const F
|
||||
ei(epsi[i])
|
||||
, facij = i == j ? F(0.5) : F(1.0)
|
||||
, denominator(epsabc - ei - ej - ek)
|
||||
, U(Zijk[i + No*j + No*No*k])
|
||||
, V(Zijk[i + No*k + No*No*j])
|
||||
, W(Zijk[j + No*i + No*No*k])
|
||||
, X(Zijk[j + No*k + No*No*i])
|
||||
, Y(Zijk[k + No*i + No*No*j])
|
||||
, Z(Zijk[k + No*j + No*No*i])
|
||||
, A(acc::maybeConjugateScalar<F>(Tijk[i + No*j + No*No*k]))
|
||||
, B(acc::maybeConjugateScalar<F>(Tijk[i + No*k + No*No*j]))
|
||||
, C(acc::maybeConjugateScalar<F>(Tijk[j + No*i + No*No*k]))
|
||||
, D(acc::maybeConjugateScalar<F>(Tijk[j + No*k + No*No*i]))
|
||||
, E(acc::maybeConjugateScalar<F>(Tijk[k + No*i + No*No*j]))
|
||||
, _F(acc::maybeConjugateScalar<F>(Tijk[k + No*j + No*No*i]))
|
||||
, value
|
||||
= 3.0 * ( A * U
|
||||
+ B * V
|
||||
+ C * W
|
||||
+ D * X
|
||||
+ E * Y
|
||||
+ _F * Z )
|
||||
+ ( ( U + X + Y )
|
||||
- 2.0 * ( V + W + Z )
|
||||
) * ( A + D + E )
|
||||
+ ( ( V + W + Z )
|
||||
- 2.0 * ( U + X + Y )
|
||||
) * ( B + C + _F )
|
||||
;
|
||||
energy += 2.0 * value / denominator * facjk * facij;
|
||||
} // i
|
||||
} // j
|
||||
} // k
|
||||
} // ii
|
||||
} // jj
|
||||
} // kk
|
||||
*_energy = acc::real(energy);
|
||||
}
|
||||
|
||||
|
||||
template <typename F>
|
||||
__MAYBE_GLOBAL__
|
||||
void getEnergySame
|
||||
( F const epsabc
|
||||
, size_t const No
|
||||
, F* const epsi
|
||||
, F* const Tijk
|
||||
, F* const Zijk
|
||||
, double* _energy
|
||||
) {
|
||||
constexpr size_t blockSize = 16;
|
||||
F energy = F(0.);
|
||||
for (size_t kk=0; kk<No; kk+=blockSize){
|
||||
const size_t kend( MIN( kk+blockSize, No) );
|
||||
for (size_t jj(kk); jj<No; jj+=blockSize){
|
||||
const size_t jend( MIN( jj+blockSize, No) );
|
||||
for (size_t ii(jj); ii<No; ii+=blockSize){
|
||||
const size_t iend( MIN( ii+blockSize, No) );
|
||||
for (size_t k(kk); k < kend; k++){
|
||||
const F ek(epsi[k]);
|
||||
const size_t jstart = jj > k ? jj : k;
|
||||
for(size_t j(jstart); j < jend; j++){
|
||||
const F facjk( j == k ? F(0.5) : F(1.0));
|
||||
const F ej(epsi[j]);
|
||||
const size_t istart = ii > j ? ii : j;
|
||||
for(size_t i(istart); i < iend; i++){
|
||||
const F
|
||||
ei(epsi[i])
|
||||
, facij ( i==j ? F(0.5) : F(1.0))
|
||||
, denominator(epsabc - ei - ej - ek)
|
||||
, U(Zijk[i + No*j + No*No*k])
|
||||
, V(Zijk[j + No*k + No*No*i])
|
||||
, W(Zijk[k + No*i + No*No*j])
|
||||
, A(acc::maybeConjugateScalar<F>(Tijk[i + No*j + No*No*k]))
|
||||
, B(acc::maybeConjugateScalar<F>(Tijk[j + No*k + No*No*i]))
|
||||
, C(acc::maybeConjugateScalar<F>(Tijk[k + No*i + No*No*j]))
|
||||
, value
|
||||
= F(3.0) * ( A * U
|
||||
+ B * V
|
||||
+ C * W
|
||||
)
|
||||
- ( A + B + C ) * ( U + V + W )
|
||||
;
|
||||
energy += F(2.0) * value / denominator * facjk * facij;
|
||||
} // i
|
||||
} // j
|
||||
} // k
|
||||
} // ii
|
||||
} // jj
|
||||
} // kk
|
||||
*_energy = acc::real(energy);
|
||||
}
|
||||
// Energy:2 ends here
|
||||
#endif /* defined(ATRIP_NEW_ENERGY) */
|
||||
|
||||
// [[file:~/cuda/atrip/atrip.org::*Energy][Energy:3]]
|
||||
// instantiate double
|
||||
template
|
||||
@ -397,8 +274,6 @@ void getEnergySame
|
||||
, DataFieldType<double>* energy
|
||||
);
|
||||
|
||||
// TODO: put this back in
|
||||
#if defined(ATRIP_WITH_COMPLEX)
|
||||
// instantiate Complex
|
||||
template
|
||||
__MAYBE_GLOBAL__
|
||||
@ -422,7 +297,6 @@ void getEnergySame
|
||||
, DataFieldType<double>* energy
|
||||
);
|
||||
// Energy:3 ends here
|
||||
#endif
|
||||
|
||||
// [[file:~/cuda/atrip/atrip.org::*Singles%20contribution][Singles contribution:2]]
|
||||
template <typename F> __MAYBE_GLOBAL__
|
||||
@ -542,7 +416,7 @@ void getEnergySame
|
||||
#if defined(ATRIP_USE_DGEMM)
|
||||
#if defined(HAVE_CUDA)
|
||||
#define REORDER(__II, __JJ, __KK) \
|
||||
reorder<<<1, No>>>(reorder_proxy< \
|
||||
reorder<<<bs, ths>>>(reorder_proxy< \
|
||||
DataFieldType<F>, \
|
||||
__II ## __JJ ## __KK \
|
||||
>{}, \
|
||||
@ -659,15 +533,11 @@ void getEnergySame
|
||||
// cuMemAlloc((CUdeviceptr*)&_vhhh,
|
||||
// NoNoNo * sizeof(DataFieldType<F>)));
|
||||
// )
|
||||
#if !defined(ATRIP_ONLY_DGEMM)
|
||||
// we still have to zero this
|
||||
const size_t
|
||||
bs = Atrip::kernelDimensions.ooo.blocks,
|
||||
ths = Atrip::kernelDimensions.ooo.threads;
|
||||
acc::zeroing<<<bs, ths>>>((DataFieldType<F>*)_t_buffer, NoNoNo);
|
||||
acc::zeroing<<<bs, ths>>>((DataFieldType<F>*)_vhhh, NoNoNo);
|
||||
#endif
|
||||
|
||||
// const size_t
|
||||
// bs = Atrip::kernelDimensions.ooo.blocks,
|
||||
// ths = Atrip::kernelDimensions.ooo.threads;
|
||||
//cuda::zeroing<<<bs, ths>>>((DataFieldType<F>*)_t_buffer, NoNoNo);
|
||||
//cuda::zeroing<<<bs, ths>>>((DataFieldType<F>*)_vhhh, NoNoNo);
|
||||
#else
|
||||
DataFieldType<F>* _t_buffer = (DataFieldType<F>*)malloc(NoNoNo * sizeof(F));
|
||||
DataFieldType<F>* _vhhh = (DataFieldType<F>*)malloc(NoNoNo * sizeof(F));
|
||||
@ -684,14 +554,12 @@ void getEnergySame
|
||||
acc::zeroing<<<bs, ths>>>((DataFieldType<F>*)Tijk,
|
||||
NoNoNo);
|
||||
)
|
||||
#endif
|
||||
|
||||
#if !defined(HAVE_CUDA)
|
||||
#else
|
||||
WITH_CHRONO("double:reorder",
|
||||
for (size_t k = 0; k < NoNoNo; k++) {
|
||||
Tijk[k] = DataFieldType<F>{0.0};
|
||||
})
|
||||
#endif /* !defined(HAVE_CUDA) */
|
||||
#endif /* defined(HAVE_CUDA) && !defined(ATRIP_ONLY_DGEMM) */
|
||||
|
||||
|
||||
#if defined(ATRIP_ONLY_DGEMM)
|
||||
|
||||
Loading…
Reference in New Issue
Block a user