141 void initIndexSet()
const override
145 const auto&
pis = this->m_cpuOwnerOverlapCopy.indexSet();
148 for (
const auto& index :
pis) {
149 if (index.local().attribute() == Dune::OwnerOverlapCopyAttributeSet::copy) {
155 if (index.local().attribute() == Dune::OwnerOverlapCopyAttributeSet::owner) {
163 this->m_indicesOwner = std::make_unique<GpuVector<int>>(
indicesOwnerCPU);
188 OPM_ERROR_IF(&source != &
dest,
"The provided GpuVectors' address did not match");
189 std::call_once(this->m_initializedIndices, [&]() { initIndexSet(); });
191 int rank = this->m_cpuOwnerOverlapCopy.communicator().rank();
192 dest.prepareSendBuf(*m_GPUSendBuf, *m_commpairIndicesOwner);
196 std::vector<MPI_Request>
sendRequests(m_messageInformation.size());
197 std::vector<MPI_Request>
recvRequests(m_messageInformation.size());
198 std::vector<int>
processMap(m_messageInformation.size());
206 for (
const_iterator info = m_messageInformation.begin(); info != end; ++info, ++i) {
208 if (info->second.second.m_size) {
209 MPI_Irecv(m_GPURecvBuf->data()+info->second.second.m_start,
214 this->m_cpuOwnerOverlapCopy.communicator(),
226 for (
const_iterator info = m_messageInformation.begin(); info != end; ++info, ++i) {
227 if (info->second.first.m_size) {
228 MPI_Issend(m_GPUSendBuf->data()+info->second.first.m_start,
233 this->m_cpuOwnerOverlapCopy.communicator(),
248 fmt::format(
"MPI_Error occurred while rank {} received a message from rank {}",
253 for (
size_t i = 0; i < m_messageInformation.size(); i++) {
256 fmt::format(
"MPI_Error occurred while rank {} sent a message from rank {}",
262 dest.syncFromRecvBuf(*m_GPURecvBuf, *m_commpairIndicesCopy);
266 mutable std::unique_ptr<GpuVector<int>> m_commpairIndicesCopy;
267 mutable std::unique_ptr<GpuVector<int>> m_commpairIndicesOwner;
268 mutable std::unique_ptr<GpuVector<field_type>> m_GPUSendBuf;
269 mutable std::unique_ptr<GpuVector<field_type>> m_GPURecvBuf;
271 struct MessageInformation
273 MessageInformation() : m_start(0), m_size(0) {}
274 MessageInformation(
size_t start,
size_t size) : m_start(start), m_size(size) {}
279 using InformationMap = std::map<int,std::pair<MessageInformation,MessageInformation> >;
280 mutable InformationMap m_messageInformation;
281 using IM = std::map<int,std::pair<std::vector<int>,std::vector<int> > >;
284 constexpr static int m_commTag = 0;
286 void buildCommPairIdxs()
const
288 const auto&
ri = this->m_cpuOwnerOverlapCopy.remoteIndices();
292 for (
const auto& process :
ri) {
293 m_im[process.first] = std::pair(std::vector<int>(), std::vector<int>());
294 for (
int send = 0; send < 2; ++send) {
295 auto remoteEnd = send ? process.second.first->end()
296 : process.second.second->end();
297 auto remote = send ? process.second.first->begin()
298 : process.second.second->begin();
301 if (send ? (
remote->localIndexPair().local().attribute() == 1)
304 m_im[process.first].first.push_back(
remote->localIndexPair().local().local());
307 m_im[process.first].second.push_back(
remote->localIndexPair().local().local());
317 for (
auto it = m_im.begin(); it != m_im.end(); it++) {
318 int noSend = it->second.first.size();
319 int noRecv = it->second.second.size();
322 m_messageInformation.insert(
323 std::make_pair(it->first,
324 std::make_pair(MessageInformation(
326 noSend * block_size *
sizeof(field_type)),
329 noRecv * block_size *
sizeof(field_type)))));
331 for (
int x = 0; x <
noSend; x++) {
332 for (
int bs = 0;
bs < block_size;
bs++) {
336 for (
int x = 0; x <
noRecv; x++) {
337 for (
int bs = 0;
bs < block_size;
bs++) {
349 m_GPUSendBuf = std::make_unique<GpuVector<field_type>>(
sendBufIdx * block_size);
350 m_GPURecvBuf = std::make_unique<GpuVector<field_type>>(
recvBufIdx * block_size);
353 void initIndexSet()
const override
357 const auto&
pis = this->m_cpuOwnerOverlapCopy.indexSet();
360 for (
const auto& index :
pis) {
361 if (index.local().attribute() == Dune::OwnerOverlapCopyAttributeSet::copy) {
367 if (index.local().attribute() == Dune::OwnerOverlapCopyAttributeSet::owner) {
375 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:186
void copyOwnerToAll(const X &source, X &dest) const override
copyOwnerToAll will copy source to the CPU, then call OwnerOverlapCopyCommunicationType::copyOwnerToA...
Definition GpuOwnerOverlapCopy.hpp:130
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:79
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