133 void initIndexSet()
const override
137 const auto&
pis = this->m_cpuOwnerOverlapCopy.indexSet();
140 for (
const auto& index :
pis) {
141 if (index.local().attribute() == Dune::OwnerOverlapCopyAttributeSet::copy) {
147 if (index.local().attribute() == Dune::OwnerOverlapCopyAttributeSet::owner) {
155 this->m_indicesOwner = std::make_unique<GpuVector<int>>(
indicesOwnerCPU);
181 OPM_ERROR_IF(&source != &
dest,
"The provided GpuVectors' address did not match");
182 std::call_once(this->m_initializedIndices, [&]() { initIndexSet(); });
184 int rank = this->m_cpuOwnerOverlapCopy.communicator().rank();
185 dest.prepareSendBuf(*m_GPUSendBuf, *m_commpairIndicesOwner);
189 std::vector<MPI_Request>
sendRequests(m_messageInformation.size());
190 std::vector<MPI_Request>
recvRequests(m_messageInformation.size());
191 std::vector<int>
processMap(m_messageInformation.size());
199 for(
const_iterator info = m_messageInformation.begin(); info != end; ++info, ++i) {
201 if(info->second.second.m_size) {
202 MPI_Irecv(m_GPURecvBuf->data()+info->second.second.m_start,
207 this->m_cpuOwnerOverlapCopy.communicator(),
218 for(
const_iterator info = m_messageInformation.begin(); info != end; ++info, ++i) {
219 if(info->second.first.m_size) {
220 MPI_Issend(m_GPUSendBuf->data()+info->second.first.m_start,
225 this->m_cpuOwnerOverlapCopy.communicator(),
239 OPM_THROW(std::runtime_error, fmt::format(
"MPI_Error occurred while rank {} received a message from rank {}", rank,
processMap[finished]));
243 for(
size_t i = 0; i < m_messageInformation.size(); i++) {
245 OPM_THROW(std::runtime_error, fmt::format(
"MPI_Error occurred while rank {} sent a message from rank {}", rank,
processMap[finished]));
250 dest.syncFromRecvBuf(*m_GPURecvBuf, *m_commpairIndicesCopy);
254 mutable std::unique_ptr<GpuVector<int>> m_commpairIndicesCopy;
255 mutable std::unique_ptr<GpuVector<int>> m_commpairIndicesOwner;
256 mutable std::unique_ptr<GpuVector<field_type>> m_GPUSendBuf;
257 mutable std::unique_ptr<GpuVector<field_type>> m_GPURecvBuf;
259 struct MessageInformation
261 MessageInformation() : m_start(0), m_size(0) {}
262 MessageInformation(
size_t start,
size_t size) : m_start(start), m_size(size) {}
267 using InformationMap = std::map<int,std::pair<MessageInformation,MessageInformation> >;
268 mutable InformationMap m_messageInformation;
269 using IM = std::map<int,std::pair<std::vector<int>,std::vector<int> > >;
272 constexpr static int m_commTag = 0;
274 void buildCommPairIdxs()
const
276 auto &
ri = this->m_cpuOwnerOverlapCopy.remoteIndices();
280 for(
auto process :
ri) {
282 m_im[process.first] = std::pair(std::vector<int>(), std::vector<int>());
283 for(
int send = 0; send < 2; ++send) {
284 auto remoteEnd = send ? process.second.first->end()
285 : process.second.second->end();
286 auto remote = send ? process.second.first->begin()
287 : process.second.second->begin();
290 if (send ? (
remote->localIndexPair().local().attribute() == 1)
294 m_im[process.first].first.push_back(
remote->localIndexPair().local().local());
296 m_im[process.first].second.push_back(
remote->localIndexPair().local().local());
306 for (
auto it = m_im.begin(); it != m_im.end(); it++) {
307 int noSend = it->second.first.size();
308 int noRecv = it->second.second.size();
311 m_messageInformation.insert(
312 std::make_pair(it->first,
313 std::make_pair(MessageInformation(
315 noSend * block_size *
sizeof(field_type)),
318 noRecv * block_size *
sizeof(field_type)))));
320 for(
int x = 0; x <
noSend; x++) {
321 for(
int bs = 0;
bs < block_size;
bs++) {
325 for(
int x = 0; x <
noRecv; x++) {
326 for(
int bs = 0;
bs < block_size;
bs++) {
338 m_GPUSendBuf = std::make_unique<GpuVector<field_type>>(
sendBufIdx * block_size);
339 m_GPURecvBuf = std::make_unique<GpuVector<field_type>>(
recvBufIdx * block_size);
342 void initIndexSet()
const override
346 const auto&
pis = this->m_cpuOwnerOverlapCopy.indexSet();
349 for (
const auto& index :
pis) {
350 if (index.local().attribute() == Dune::OwnerOverlapCopyAttributeSet::copy) {
356 if (index.local().attribute() == Dune::OwnerOverlapCopyAttributeSet::owner) {
364 this->m_indicesOwner = std::make_unique<GpuVector<int>>(
indicesOwnerCPU);
void copyOwnerToAll(const X &source, X &dest) const override
copyOwnerToAll will copy source to the CPU, then call OwnerOverlapCopyCommunicationType::copyOwnerToA...
Definition GpuOwnerOverlapCopy.hpp:178
void copyOwnerToAll(const X &source, X &dest) const override
copyOwnerToAll will copy source to the CPU, then call OwnerOverlapCopyCommunicationType::copyOwnerToA...
Definition GpuOwnerOverlapCopy.hpp:123
virtual void copyOwnerToAll(const X &source, X &dest) const =0
copyOwnerToAll will copy source to the CPU, then call OwnerOverlapCopyCommunicationType::copyOwnerToA...
void dot(const X &x, const X &y, field_type &output) const
dot will carry out the dot product between x and y on the owned indices, then sum up the result acros...
Definition GpuOwnerOverlapCopy.hpp:72
int to_int(std::size_t s)
to_int converts a (on most relevant platforms) 64 bits unsigned size_t to a signed 32 bits signed int
Definition safe_conversion.hpp:52
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242