KratosMultiphysics
KRATOS Multiphysics (Kratos) is a framework for building parallel, multi-disciplinary simulation software, aiming at modularity, extensibility, and high performance. Kratos is written in C++, and counts with an extensive Python interface.
data_communicator.h
Go to the documentation of this file.
1 // | / |
2 // ' / __| _` | __| _ \ __|
3 // . \ | ( | | ( |\__ `
4 // _|\_\_| \__,_|\__|\___/ ____/
5 // Multi-Physics
6 //
7 // License: BSD License
8 // Kratos default license: kratos/license.txt
9 //
10 // Main author: Jordi Cotela
11 //
12 
13 #pragma once
14 
15 // System includes
16 #include <string>
17 #include <iostream>
18 #include <type_traits>
19 
20 // External includes
21 
22 // Project includes
23 #include "containers/array_1d.h"
24 #include "containers/flags.h"
25 #include "includes/define.h"
27 
28 // Using a macro instead of a function to get the correct line in the error message.
29 #ifndef KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK
30 #define KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(Size1, Size2, CheckedFunction) \
31  KRATOS_DEBUG_ERROR_IF(Size1 != Size2) \
32  << "Input error in call to DataCommunicator::" << CheckedFunction \
33  << ": The sizes of the local and distributed buffers do not match." << std::endl;
34 #endif
35 
36 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SYNC_SHAPE_INTERFACE_FOR_TYPE
37 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SYNC_SHAPE_INTERFACE_FOR_TYPE(...) \
38 virtual bool SynchronizeShape(__VA_ARGS__&) const { return false; } \
39 virtual bool SynchronizeShape( \
40  const __VA_ARGS__& rSendValue, const int SendDestination, const int SendTag, \
41  __VA_ARGS__& rRecvValue, const int RecvSource, const int RecvTag) const { return false; } \
42 
43 #endif
44 
45 // Methods based on MPI_Reduce, supporting sum, max or min operations.
46 /* Variants for each method are provided, either returning the reduced value or filling a provided vector buffer.
47  * The returned value is only meaningful on the Root rank.
48  */
49 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_REDUCE_INTERFACE_FOR_TYPE
50 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_REDUCE_INTERFACE_FOR_TYPE(...) \
51 virtual __VA_ARGS__ Sum(const __VA_ARGS__& rLocalValue, const int Root) const { return rLocalValue; } \
52 virtual std::vector<__VA_ARGS__> Sum(const std::vector<__VA_ARGS__>& rLocalValues, const int Root) const { \
53  return rLocalValues; \
54 } \
55 virtual void Sum(const std::vector<__VA_ARGS__>& rLocalValues, std::vector<__VA_ARGS__>& rGlobalValues, const int Root) const { \
56  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rLocalValues.size(), rGlobalValues.size(), "Sum"); \
57  rGlobalValues = Sum(rLocalValues, Root); \
58 } \
59 virtual __VA_ARGS__ Min(const __VA_ARGS__& rLocalValue, const int Root) const { return rLocalValue; } \
60 virtual std::vector<__VA_ARGS__> Min(const std::vector<__VA_ARGS__>& rLocalValues, const int Root) const { \
61  return rLocalValues; \
62 } \
63 virtual void Min(const std::vector<__VA_ARGS__>& rLocalValues, std::vector<__VA_ARGS__>& rGlobalValues, const int Root) const { \
64  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rLocalValues.size(), rGlobalValues.size(), "Min"); \
65  rGlobalValues = Min(rLocalValues, Root); \
66 } \
67 virtual __VA_ARGS__ Max(const __VA_ARGS__& rLocalValue, const int Root) const { return rLocalValue; } \
68 virtual std::vector<__VA_ARGS__> Max(const std::vector<__VA_ARGS__>& rLocalValues, const int Root) const { \
69  return rLocalValues; \
70 } \
71 virtual void Max(const std::vector<__VA_ARGS__>& rLocalValues, std::vector<__VA_ARGS__>& rGlobalValues, const int Root) const { \
72  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rLocalValues.size(), rGlobalValues.size(), "Max"); \
73  rGlobalValues = Max(rLocalValues, Root); \
74 } \
75 
76 #endif
77 
78 // Methods based on MPI_Allreduce, supporting sum, max or min operations.
79 /* Variants for each method are provided, either returning the reduced value or filling a provided vector buffer.
80  * The returned value is defined on all ranks.
81  */
82 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_ALLREDUCE_INTERFACE_FOR_TYPE
83 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_ALLREDUCE_INTERFACE_FOR_TYPE(...) \
84 virtual __VA_ARGS__ SumAll(const __VA_ARGS__& rLocalValue) const { return rLocalValue; } \
85 virtual std::vector<__VA_ARGS__> SumAll(const std::vector<__VA_ARGS__>& rLocalValues) const { \
86  return rLocalValues; \
87 } \
88 virtual void SumAll(const std::vector<__VA_ARGS__>& rLocalValues, std::vector<__VA_ARGS__>& rGlobalValues) const { \
89  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rLocalValues.size(), rGlobalValues.size(), "SumAll"); \
90  rGlobalValues = SumAll(rLocalValues); \
91 } \
92 virtual __VA_ARGS__ MinAll(const __VA_ARGS__& rLocalValue) const { return rLocalValue; } \
93 virtual std::vector<__VA_ARGS__> MinAll(const std::vector<__VA_ARGS__>& rLocalValues) const { \
94  return rLocalValues; \
95 } \
96 virtual void MinAll(const std::vector<__VA_ARGS__>& rLocalValues, std::vector<__VA_ARGS__>& rGlobalValues) const { \
97  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rLocalValues.size(), rGlobalValues.size(), "MinAll"); \
98  rGlobalValues = MinAll(rLocalValues); \
99 } \
100 virtual __VA_ARGS__ MaxAll(const __VA_ARGS__& rLocalValue) const { return rLocalValue; } \
101 virtual std::vector<__VA_ARGS__> MaxAll(const std::vector<__VA_ARGS__>& rLocalValues) const { \
102  return rLocalValues; \
103 } \
104 virtual void MaxAll(const std::vector<__VA_ARGS__>& rLocalValues, std::vector<__VA_ARGS__>& rGlobalValues) const { \
105  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rLocalValues.size(), rGlobalValues.size(), "MaxAll"); \
106  rGlobalValues = MaxAll(rLocalValues); \
107 }
108 #endif
109 
110 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_ALLREDUCE_LOC_INTERFACE_FOR_TYPE
111 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_ALLREDUCE_LOC_INTERFACE_FOR_TYPE(...) \
112 virtual std::pair<__VA_ARGS__, int> MinLocAll(const __VA_ARGS__& rLocalValue) const { return std::pair<__VA_ARGS__, int>(rLocalValue, 0); } \
113 virtual std::pair<__VA_ARGS__, int> MaxLocAll(const __VA_ARGS__& rLocalValue) const { return std::pair<__VA_ARGS__, int>(rLocalValue, 0); }
114 #endif
115 
116 // Compute the partial sum of the given quantity from rank 0 to the current rank (included).
117 /* This is a wrapper to MPI_Scan.
118  * Variants for each method are provided, either returning the reduced value or filling a provided vector buffer.
119  * The returned value is defined on all ranks.
120  */
121 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SCANSUM_INTERFACE_FOR_TYPE
122 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SCANSUM_INTERFACE_FOR_TYPE(...) \
123 virtual __VA_ARGS__ ScanSum(const __VA_ARGS__& rLocalValue) const { return rLocalValue; } \
124 virtual std::vector<__VA_ARGS__> ScanSum(const std::vector<__VA_ARGS__>& rLocalValues) const { \
125  return rLocalValues; \
126 } \
127 virtual void ScanSum(const std::vector<__VA_ARGS__>& rLocalValues, std::vector<__VA_ARGS__>& rPartialSums) const { \
128  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rLocalValues.size(), rPartialSums.size(), "ScanSum"); \
129  rPartialSums = ScanSum(rLocalValues); \
130 } \
131 
132 #endif
133 
134 // Exchange data with other ranks. This is a wrapper for MPI_Sendrecv, MPI_Send and MPI_Recv.
135 /* Versions which outputting the result as a return argument or by filling an output buffer argument are provided.
136  * The return version has a performance overhead, since the dimensions of the receiving buffer have to be
137  * communicated. If the dimensions of the receiving buffer are known at the destination rank, the output buffer
138  * variant should be preferred.
139  */
140 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SENDRECV_INTERFACE_FOR_TYPE
141 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SENDRECV_INTERFACE_FOR_TYPE(...) \
142 virtual __VA_ARGS__ SendRecvImpl( \
143  const __VA_ARGS__& rSendValues, const int SendDestination, const int SendTag, \
144  const int RecvSource, const int RecvTag) const { \
145  KRATOS_ERROR_IF( (Rank() != SendDestination) || (Rank() != RecvSource)) \
146  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
147  return rSendValues; \
148 } \
149 virtual std::vector<__VA_ARGS__> SendRecvImpl( \
150  const std::vector<__VA_ARGS__>& rSendValues, const int SendDestination, const int SendTag, \
151  const int RecvSource, const int RecvTag) const { \
152  KRATOS_ERROR_IF( (Rank() != SendDestination) || (Rank() != RecvSource)) \
153  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
154  return rSendValues; \
155 } \
156 virtual void SendRecvImpl( \
157  const __VA_ARGS__& rSendValues, const int SendDestination, const int SendTag, \
158  __VA_ARGS__& rRecvValues, const int RecvSource, const int RecvTag) const { \
159  rRecvValues = SendRecvImpl(rSendValues, SendDestination, SendTag, RecvSource, RecvTag); \
160 } \
161 virtual void SendRecvImpl( \
162  const std::vector<__VA_ARGS__>& rSendValues, const int SendDestination, const int SendTag, \
163  std::vector<__VA_ARGS__>& rRecvValues, const int RecvSource, const int RecvTag) const { \
164  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rSendValues.size(), rRecvValues.size(), "SendRecv"); \
165  rRecvValues = SendRecvImpl(rSendValues, SendDestination, SendTag, RecvSource, RecvTag); \
166 } \
167 virtual void SendImpl( \
168  const __VA_ARGS__& rSendValues, const int SendDestination, const int SendTag = 0) const { \
169  KRATOS_ERROR_IF(Rank() != SendDestination) \
170  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
171 } \
172 virtual void SendImpl( \
173  const std::vector<__VA_ARGS__>& rSendValues, const int SendDestination, const int SendTag = 0) const { \
174  KRATOS_ERROR_IF(Rank() != SendDestination) \
175  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
176 } \
177 virtual void RecvImpl(__VA_ARGS__& rRecvValues, const int RecvSource, const int RecvTag = 0) const { \
178  KRATOS_ERROR << "Calling serial DataCommunicator::Recv, which has no meaningful return." << std::endl; \
179 } \
180 virtual void RecvImpl(std::vector<__VA_ARGS__>& rRecvValues, const int RecvSource, const int RecvTag = 0) const { \
181  KRATOS_ERROR << "Calling serial DataCommunicator::Recv, which has no meaningful return." << std::endl; \
182 } \
183 
184 #endif
185 
186 // Synchronize a buffer to the value held by the broadcasting rank.
187 /* This is a wrapper for MPI_Bcast.
188  * @param[in/out] The broadcast value (input on SourceRank, output on all other ranks).
189  * @param[in] SourceRank The rank transmitting the value.
190  */
191 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_BROADCAST_INTERFACE_FOR_TYPE
192 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_BROADCAST_INTERFACE_FOR_TYPE(...) \
193 virtual void BroadcastImpl(__VA_ARGS__& rBuffer, const int SourceRank) const {} \
194 virtual void BroadcastImpl(std::vector<__VA_ARGS__>& rBuffer, const int SourceRank) const {} \
195 
196 #endif
197 
199 /* Versions which outputting the result as a return argument or by filling an output buffer argument are provided.
200  * The return version has a performance overhead, since the dimensions of the receiving buffer have to be
201  * communicated. If the dimensions of the receiving buffers are known at the destination rank, the output buffer
202  * variant should be preferred.
203  */
204 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SCATTER_INTERFACE_FOR_TYPE
205 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SCATTER_INTERFACE_FOR_TYPE(...) \
206 virtual std::vector<__VA_ARGS__> Scatter(const std::vector<__VA_ARGS__>& rSendValues, const int SourceRank) const { \
207  KRATOS_ERROR_IF( Rank() != SourceRank ) \
208  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
209  return rSendValues; \
210 } \
211 virtual void Scatter( \
212  const std::vector<__VA_ARGS__>& rSendValues, std::vector<__VA_ARGS__>& rRecvValues, const int SourceRank) const { \
213  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rSendValues.size(),rRecvValues.size(),"Scatter"); \
214  rRecvValues = Scatter(rSendValues, SourceRank); \
215 } \
216 virtual std::vector<__VA_ARGS__> Scatterv(const std::vector<std::vector<__VA_ARGS__>>& rSendValues, const int SourceRank) const { \
217  KRATOS_ERROR_IF( Rank() != SourceRank ) \
218  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
219  KRATOS_ERROR_IF( static_cast<unsigned int>(Size()) != rSendValues.size() ) \
220  << "Unexpected number of sends in DataCommuncatior::Scatterv (serial DataCommunicator always assumes a single process)." << std::endl; \
221  return rSendValues[0]; \
222 } \
223 virtual void Scatterv( \
224  const std::vector<__VA_ARGS__>& rSendValues, const std::vector<int>& rSendCounts, const std::vector<int>& rSendOffsets, \
225  std::vector<__VA_ARGS__>& rRecvValues, const int SourceRank) const { \
226  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rRecvValues.size(), rSendValues.size(), "Scatterv (values check)"); \
227  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rSendCounts.size(), 1, "Scatterv (counts check)"); \
228  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rSendOffsets.size(), 1, "Scatterv (offsets check)"); \
229  KRATOS_ERROR_IF( Rank() != SourceRank ) \
230  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
231  rRecvValues = rSendValues; \
232 } \
233 
234 #endif
235 
237 /* Versions which outputting the result as a return argument or by filling an output buffer argument are provided.
238  * The return version has a performance overhead, since the dimensions of the receiving buffer have to be
239  * communicated. If the dimensions of the receiving buffers are known at the destination rank, the output buffer
240  * variant should be preferred.
241  */
242 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_GATHER_INTERFACE_FOR_TYPE
243 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_GATHER_INTERFACE_FOR_TYPE(...) \
244 virtual std::vector<__VA_ARGS__> Gather(const std::vector<__VA_ARGS__>& rSendValues, const int DestinationRank) const { \
245  KRATOS_ERROR_IF( Rank() != DestinationRank ) \
246  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
247  return rSendValues; \
248 } \
249 virtual void Gather( \
250  const std::vector<__VA_ARGS__>& rSendValues, std::vector<__VA_ARGS__>& rRecvValues, const int DestinationRank) const { \
251  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rSendValues.size(),rRecvValues.size(),"Gather"); \
252  rRecvValues = Gather(rSendValues, DestinationRank); \
253 } \
254 virtual std::vector<std::vector<__VA_ARGS__>> Gatherv( \
255  const std::vector<__VA_ARGS__>& rSendValues, const int DestinationRank) const { \
256  KRATOS_ERROR_IF( Rank() != DestinationRank ) \
257  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
258  return std::vector<std::vector<__VA_ARGS__>>{rSendValues}; \
259 } \
260 virtual void Gatherv( \
261  const std::vector<__VA_ARGS__>& rSendValues, std::vector<__VA_ARGS__>& rRecvValues, \
262  const std::vector<int>& rRecvCounts, const std::vector<int>& rRecvOffsets, const int DestinationRank) const { \
263  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rRecvValues.size(), rSendValues.size(), "Gatherv (values check)"); \
264  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rRecvCounts.size(), 1, "Gatherv (counts check)"); \
265  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rRecvOffsets.size(), 1, "Gatherv (offset check)"); \
266  KRATOS_ERROR_IF( Rank() != DestinationRank ) \
267  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl; \
268  rRecvValues = rSendValues; \
269 } \
270 virtual std::vector<__VA_ARGS__> AllGather(const std::vector<__VA_ARGS__>& rSendValues) const { return rSendValues; } \
271 virtual void AllGather(const std::vector<__VA_ARGS__>& rSendValues, std::vector<__VA_ARGS__>& rRecvValues) const { \
272  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rSendValues.size(),rRecvValues.size(),"AllGather"); \
273  rRecvValues = AllGather(rSendValues); \
274 } \
275 virtual std::vector<std::vector<__VA_ARGS__>> AllGatherv(const std::vector<__VA_ARGS__>& rSendValues) const { \
276  return std::vector<std::vector<__VA_ARGS__>>{rSendValues}; \
277 } \
278 virtual void AllGatherv(const std::vector<__VA_ARGS__>& rSendValues, std::vector<__VA_ARGS__>& rRecvValues, \
279  const std::vector<int>& rRecvCounts, const std::vector<int>& rRecvOffsets) const { \
280  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rRecvValues.size(), rSendValues.size(), "AllGatherv (values check)"); \
281  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rRecvCounts.size(), 1, "AllGatherv (counts check)"); \
282  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rRecvOffsets.size(), 1, "AllGatherv (offset check)"); \
283  rRecvValues = rSendValues; \
284 }
285 #endif
286 
287 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_PUBLIC_INTERFACE_FOR_TYPE
288 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_PUBLIC_INTERFACE_FOR_TYPE(...) \
289 KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_REDUCE_INTERFACE_FOR_TYPE(__VA_ARGS__) \
290 KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_ALLREDUCE_INTERFACE_FOR_TYPE(__VA_ARGS__) \
291 KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SCANSUM_INTERFACE_FOR_TYPE(__VA_ARGS__) \
292 KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SCATTER_INTERFACE_FOR_TYPE(__VA_ARGS__) \
293 KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_GATHER_INTERFACE_FOR_TYPE(__VA_ARGS__) \
294 KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SYNC_SHAPE_INTERFACE_FOR_TYPE(__VA_ARGS__) \
295 
296 #endif
297 
298 #ifndef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_IMPLEMENTATION_FOR_TYPE
299 #define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_IMPLEMENTATION_FOR_TYPE(...) \
300 KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SENDRECV_INTERFACE_FOR_TYPE(__VA_ARGS__) \
301 KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_BROADCAST_INTERFACE_FOR_TYPE(__VA_ARGS__) \
302 
303 #endif
304 
305 
306 namespace Kratos
307 {
310 
313 
315 
317 class KRATOS_API(KRATOS_CORE) DataCommunicator
318 {
319  private:
320 
321  template<typename T> class serialization_is_required {
322  private:
323 
324  template<typename U> struct serialization_traits {
325  constexpr static bool is_std_vector = false;
326  constexpr static bool value_type_is_compound = false;
327  constexpr static bool value_type_is_bool = false;
328  };
329 
330  template<typename U> struct serialization_traits<std::vector<U>> {
331  constexpr static bool is_std_vector = true;
332  constexpr static bool value_type_is_compound = std::is_compound<U>::value;
333  constexpr static bool value_type_is_bool = std::is_same<U, bool>::value;
334  };
335 
336  constexpr static bool is_vector_of_simple_types = serialization_traits<T>::is_std_vector && !serialization_traits<T>::value_type_is_compound;
337  constexpr static bool is_vector_of_bools = serialization_traits<T>::is_std_vector && serialization_traits<T>::value_type_is_bool;
338 
339  constexpr static bool is_vector_of_directly_communicable_type = is_vector_of_simple_types && !is_vector_of_bools;
340 
341  public:
342  constexpr static bool value = std::is_compound<T>::value && !is_vector_of_directly_communicable_type;
343  };
344 
345  template<bool value> struct TypeFromBool {};
346 
347  template<typename T> void CheckSerializationForSimpleType(const T& rSerializedType, TypeFromBool<true>) const {}
348 
349  template<typename T>
350  KRATOS_DEPRECATED_MESSAGE("Calling serialization-based communication for a simple type. Please implement direct communication support for this type.")
351  void CheckSerializationForSimpleType(const T& rSerializedType, TypeFromBool<false>) const {}
352 
353  public:
356 
359 
363 
366 
368  virtual ~DataCommunicator() {}
369 
373 
375 
379  static DataCommunicator::UniquePointer Create()
380  {
381  return Kratos::make_unique<DataCommunicator>();
382  }
383 
385 
386  virtual void Barrier() const {}
387 
388  // Complete interface for basic types
389 
401 
402  // MinLoc and MaxLoc AllReduce operations
408 
409  // Reduce operations
410 
411  virtual bool AndReduce(
412  const bool Value,
413  const int Root) const
414  {
415  return Value;
416  }
417 
419  const Kratos::Flags Values,
420  const Kratos::Flags Mask,
421  const int Root) const
422  {
423  return Values;
424  }
425 
426  virtual bool OrReduce(
427  const bool Value,
428  const int Root) const
429  {
430  return Value;
431  }
432 
434  const Kratos::Flags Values,
435  const Kratos::Flags Mask,
436  const int Root) const
437  {
438  return Values;
439  }
440 
441  // Allreduce operations
442 
443  virtual bool AndReduceAll(const bool Value) const
444  {
445  return Value;
446  }
447 
448  virtual Kratos::Flags AndReduceAll(const Kratos::Flags Values, const Kratos::Flags Mask) const
449  {
450  return Values;
451  }
452 
453  virtual bool OrReduceAll(const bool Value) const
454  {
455  return Value;
456  }
457 
458  virtual Kratos::Flags OrReduceAll(const Kratos::Flags Values, const Kratos::Flags Mask) const
459  {
460  return Values;
461  }
462 
463  // Broadcast operations
464 
466 
472  template<typename TObject>
473  void Broadcast(TObject& rBroadcastObject, const int SourceRank) const
474  {
475  this->BroadcastImpl(rBroadcastObject, SourceRank);
476  }
477 
478  // Sendrecv operations
479 
481 
491  template<typename TObject>
492  TObject SendRecv(
493  const TObject& rSendObject, const int SendDestination, const int SendTag,
494  const int RecvSource, const int RecvTag) const
495  {
496  return this->SendRecvImpl(rSendObject, SendDestination, SendTag, RecvSource, RecvTag);
497  }
498 
500 
508  template<class TObject>
509  TObject SendRecv(
510  const TObject& rSendObject, const int SendDestination, const int RecvSource) const
511  {
512  return this->SendRecvImpl(rSendObject, SendDestination, 0, RecvSource, 0);
513  }
514 
516 
526  template<class TObject>
527  void SendRecv(
528  const TObject& rSendObject, const int SendDestination, const int SendTag,
529  TObject& rRecvObject, const int RecvSource, const int RecvTag) const
530  {
531  this->SendRecvImpl(rSendObject, SendDestination, SendTag, rRecvObject, RecvSource, RecvTag);
532  }
533 
535 
543  template<class TObject>
544  void SendRecv(
545  const TObject& rSendObject, const int SendDestination, TObject& rRecvObject, const int RecvSource) const
546  {
547  this->SendRecvImpl(rSendObject, SendDestination, 0, rRecvObject, RecvSource, 0);
548  }
549 
551 
558  template<typename TObject>
559  void Send(const TObject& rSendValues, const int SendDestination, const int SendTag = 0) const
560  {
561  this->SendImpl(rSendValues, SendDestination, SendTag);
562  }
563 
565 
572  template<typename TObject>
573  void Recv(TObject& rRecvObject, const int RecvSource, const int RecvTag = 0) const
574  {
575  this->RecvImpl(rRecvObject, RecvSource, RecvTag);
576  }
577 
581 
587  virtual int Rank() const
588  {
589  return 0;
590  }
591 
597  virtual int Size() const
598  {
599  return 1;
600  }
601 
606  virtual bool IsDistributed() const
607  {
608  return false;
609  }
610 
616  virtual bool IsDefinedOnThisRank() const
617  {
618  return true;
619  }
620 
626  virtual bool IsNullOnThisRank() const
627  {
628  return false;
629  }
630 
639  const std::vector<int>& rRanks,
640  const std::string& rNewCommunicatorName
641  ) const
642  {
643  return *this;
644  }
645 
649 
651 
653  KRATOS_DEPRECATED_MESSAGE("This function is deprecated, please retrieve the DataCommunicator through the ModelPart (or by name in special cases)")
654  static DataCommunicator& GetDefault();
655 
659 
661 
680  virtual bool BroadcastErrorIfTrue(bool Condition, const int SourceRank) const
681  {
682  return Condition;
683  }
684 
686 
705  virtual bool BroadcastErrorIfFalse(bool Condition, const int SourceRank) const
706  {
707  return Condition;
708  }
709 
711 
729  virtual bool ErrorIfTrueOnAnyRank(bool Condition) const
730  {
731  return Condition;
732  }
733 
735 
753  virtual bool ErrorIfFalseOnAnyRank(bool Condition) const
754  {
755  return Condition;
756  }
757 
761 
763  virtual std::string Info() const
764  {
765  std::stringstream buffer;
766  PrintInfo(buffer);
767  return buffer.str();
768  }
769 
771  virtual void PrintInfo(std::ostream &rOStream) const
772  {
773  rOStream << "DataCommunicator";
774  }
775 
777  virtual void PrintData(std::ostream &rOStream) const
778  {
779  rOStream
780  << "Serial do-nothing version of the Kratos wrapper for MPI communication.\n"
781  << "Rank 0 of 1 assumed." << std::endl;
782  }
783 
785 
786  protected:
787 
790 
802 
803 
808  virtual void BroadcastImpl(std::string& rBuffer, const int SourceRank) const {};
809 
811 
815  template<class TObject>
816  void BroadcastImpl(TObject& rBroadcastObject, const int SourceRank) const
817  {
818  CheckSerializationForSimpleType(rBroadcastObject, TypeFromBool<serialization_is_required<TObject>::value>());
819  if (this->IsDistributed())
820  {
821  unsigned int message_size;
822  std::string broadcast_message;
823  int rank = this->Rank();
824  if (rank == SourceRank)
825  {
826  MpiSerializer send_serializer;
827  send_serializer.save("data", rBroadcastObject);
828  broadcast_message = send_serializer.GetStringRepresentation();
829 
830  message_size = broadcast_message.size();
831  }
832 
833  this->Broadcast(message_size, SourceRank);
834 
835  if (rank != SourceRank)
836  {
837  broadcast_message.resize(message_size);
838  }
839 
840  this->Broadcast(broadcast_message, SourceRank);
841 
842  if (rank != SourceRank)
843  {
844  MpiSerializer recv_serializer(broadcast_message);
845  recv_serializer.load("data", rBroadcastObject);
846  }
847  }
848  }
849 
851 
859  virtual void SendRecvImpl(
860  const std::string& rSendValues, const int SendDestination, const int SendTag,
861  std::string& rRecvValues, const int RecvSource, const int RecvTag) const
862  {
863  KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(rSendValues.size(), rRecvValues.size(), "SendRecv");
864  rRecvValues = SendRecvImpl(rSendValues, SendDestination, SendTag, RecvSource, RecvTag);
865  }
866 
868 
876  virtual std::string SendRecvImpl(
877  const std::string& rSendValues, const int SendDestination, const int SendTag,
878  const int RecvSource, const int RecvTag) const
879  {
880  KRATOS_ERROR_IF( (Rank() != SendDestination) || (Rank() != RecvSource))
881  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl;
882  return rSendValues;
883  }
884 
886 
895  template<class TObject> TObject SendRecvImpl(
896  const TObject& rSendObject,
897  const int SendDestination, const int SendTag,
898  const int RecvSource, const int RecvTag) const
899  {
900  CheckSerializationForSimpleType(rSendObject, TypeFromBool<serialization_is_required<TObject>::value>());
901  if (this->IsDistributed())
902  {
903  MpiSerializer send_serializer;
904  send_serializer.save("data", rSendObject);
905  std::string send_message = send_serializer.GetStringRepresentation();
906 
907  std::string recv_message = this->SendRecv(send_message, SendDestination, RecvSource);
908 
909  MpiSerializer recv_serializer(recv_message);
910  TObject recv_object;
911  recv_serializer.load("data", recv_object);
912  return recv_object;
913  }
914  else
915  {
916  KRATOS_ERROR_IF( (Rank() != SendDestination) || (Rank() != RecvSource))
917  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl;
918 
919  return rSendObject;
920  }
921  }
922 
924 
929  virtual void SendImpl(const std::string& rSendValues, const int SendDestination, const int SendTag) const
930  {
931  KRATOS_ERROR_IF(Rank() != SendDestination)
932  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl;
933  }
934 
936 
942  template<class TObject> void SendImpl(
943  const TObject& rSendObject, const int SendDestination, const int SendTag) const
944  {
945  CheckSerializationForSimpleType(rSendObject, TypeFromBool<serialization_is_required<TObject>::value>());
946  if (this->IsDistributed())
947  {
948  MpiSerializer send_serializer;
949  send_serializer.save("data", rSendObject);
950  std::string send_message = send_serializer.GetStringRepresentation();
951 
952  this->SendImpl(send_message, SendDestination, SendTag);
953  }
954  else
955  {
956  KRATOS_ERROR_IF(Rank() != SendDestination)
957  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl;
958  }
959  }
960 
962 
967  virtual void RecvImpl(std::string& rRecvValues, const int RecvSource, const int RecvTag = 0) const
968  {
969  KRATOS_ERROR << "Calling serial DataCommunicator::Recv, which has no meaningful return." << std::endl;
970  }
971 
973 
979  template<class TObject> void RecvImpl(
980  TObject& rRecvObject, const int RecvSource, const int RecvTag = 0) const
981  {
982  CheckSerializationForSimpleType(rRecvObject, TypeFromBool<serialization_is_required<TObject>::value>());
983  if (this->IsDistributed())
984  {
985  std::string recv_message;
986 
987  this->Recv(recv_message, RecvSource, RecvTag);
988 
989  MpiSerializer recv_serializer(recv_message);
990  recv_serializer.load("data", rRecvObject);
991  }
992  else
993  {
994  KRATOS_ERROR_IF(Rank() != RecvSource)
995  << "Communication between different ranks is not possible with a serial DataCommunicator." << std::endl;
996  }
997  }
998 
1000 
1001  private:
1002 
1005 
1007  DataCommunicator(DataCommunicator const &rOther) = delete;
1008 
1010  DataCommunicator &operator=(DataCommunicator const &rOther) = delete;
1011 
1013 
1014 }; // Class DataCommunicator
1015 
1017 
1020 
1024 
1026 inline std::istream &operator>>(std::istream &rIStream,
1027  DataCommunicator &rThis)
1028 {
1029  return rIStream;
1030 }
1031 
1033 inline std::ostream &operator<<(std::ostream &rOStream,
1034  const DataCommunicator &rThis)
1035 {
1036  rThis.PrintInfo(rOStream);
1037  rOStream << std::endl;
1038  rThis.PrintData(rOStream);
1039 
1040  return rOStream;
1041 }
1043 
1045 
1046 } // namespace Kratos.
1047 
1048 #undef KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK
1049 
1050 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SYNC_SHAPE_INTERFACE_FOR_TYPE
1051 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_REDUCE_INTERFACE_FOR_TYPE
1052 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_ALLREDUCE_INTERFACE_FOR_TYPE
1053 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SCANSUM_INTERFACE_FOR_TYPE
1054 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SENDRECV_INTERFACE_FOR_TYPE
1055 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_BROADCAST_INTERFACE_FOR_TYPE
1056 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_SCATTER_INTERFACE_FOR_TYPE
1057 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_GATHER_INTERFACE_FOR_TYPE
1058 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_PUBLIC_INTERFACE_FOR_TYPE
1059 #undef KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_IMPLEMENTATION_FOR_TYPE
PeriodicInterfaceProcess & operator=(const PeriodicInterfaceProcess &)=delete
Base class for all Conditions.
Definition: condition.h:59
Serial (do-nothing) version of a wrapper class for MPI communication.
Definition: data_communicator.h:318
static DataCommunicator::UniquePointer Create()
Create a new DataCommunicator as a copy of this one.
Definition: data_communicator.h:379
virtual Kratos::Flags AndReduceAll(const Kratos::Flags Values, const Kratos::Flags Mask) const
Definition: data_communicator.h:448
virtual bool AndReduceAll(const bool Value) const
Definition: data_communicator.h:443
virtual void PrintData(std::ostream &rOStream) const
Print object's data.
Definition: data_communicator.h:777
virtual bool OrReduceAll(const bool Value) const
Definition: data_communicator.h:453
virtual void RecvImpl(std::string &rRecvValues, const int RecvSource, const int RecvTag=0) const
Receive data from other ranks (string version).
Definition: data_communicator.h:967
virtual bool ErrorIfTrueOnAnyRank(bool Condition) const
This function throws an error on ranks where Condition evaluates to false, if it evaluated to true on...
Definition: data_communicator.h:729
TObject SendRecvImpl(const TObject &rSendObject, const int SendDestination, const int SendTag, const int RecvSource, const int RecvTag) const
Exchange data with other ranks (generic version).
Definition: data_communicator.h:895
virtual void PrintInfo(std::ostream &rOStream) const
Print information about this object.
Definition: data_communicator.h:771
virtual bool IsDefinedOnThisRank() const
Check whether this DataCommunicator involves the current rank.
Definition: data_communicator.h:616
virtual void Barrier() const
Pause program execution until all threads reach this call.
Definition: data_communicator.h:386
virtual bool OrReduce(const bool Value, const int Root) const
Definition: data_communicator.h:426
void SendRecv(const TObject &rSendObject, const int SendDestination, TObject &rRecvObject, const int RecvSource) const
Exchange data with other ranks.
Definition: data_communicator.h:544
virtual Kratos::Flags OrReduce(const Kratos::Flags Values, const Kratos::Flags Mask, const int Root) const
Definition: data_communicator.h:433
virtual Kratos::Flags OrReduceAll(const Kratos::Flags Values, const Kratos::Flags Mask) const
Definition: data_communicator.h:458
void RecvImpl(TObject &rRecvObject, const int RecvSource, const int RecvTag=0) const
Exchange data with other ranks (generic version).
Definition: data_communicator.h:979
void SendImpl(const TObject &rSendObject, const int SendDestination, const int SendTag) const
Exchange data with other ranks (generic version).
Definition: data_communicator.h:942
virtual void SendRecvImpl(const std::string &rSendValues, const int SendDestination, const int SendTag, std::string &rRecvValues, const int RecvSource, const int RecvTag) const
Exchange data with other ranks (string version).
Definition: data_communicator.h:859
virtual std::string Info() const
Turn back information as a string.
Definition: data_communicator.h:763
virtual const DataCommunicator & GetSubDataCommunicator(const std::vector< int > &rRanks, const std::string &rNewCommunicatorName) const
Get a sub-data communicator.
Definition: data_communicator.h:638
virtual bool IsNullOnThisRank() const
Check whether this DataCommunicator is MPI_COMM_NULL for the current rank.
Definition: data_communicator.h:626
void Recv(TObject &rRecvObject, const int RecvSource, const int RecvTag=0) const
Exchange data with other ranks.
Definition: data_communicator.h:573
virtual int Size() const
Get the parallel size of this DataCommunicator.
Definition: data_communicator.h:597
DataCommunicator()
Default constructor.
Definition: data_communicator.h:365
void SendRecv(const TObject &rSendObject, const int SendDestination, const int SendTag, TObject &rRecvObject, const int RecvSource, const int RecvTag) const
Exchange data with other ranks.
Definition: data_communicator.h:527
virtual bool BroadcastErrorIfFalse(bool Condition, const int SourceRank) const
This function throws an error on ranks != Sourcerank if Condition evaluates to false.
Definition: data_communicator.h:705
virtual Kratos::Flags AndReduce(const Kratos::Flags Values, const Kratos::Flags Mask, const int Root) const
Definition: data_communicator.h:418
TObject SendRecv(const TObject &rSendObject, const int SendDestination, const int RecvSource) const
Exchange data with other ranks.
Definition: data_communicator.h:509
KRATOS_CLASS_POINTER_DEFINITION(DataCommunicator)
Pointer definition of DataCommunicator.
TObject SendRecv(const TObject &rSendObject, const int SendDestination, const int SendTag, const int RecvSource, const int RecvTag) const
Exchange data with other ranks.
Definition: data_communicator.h:492
virtual void SendImpl(const std::string &rSendValues, const int SendDestination, const int SendTag) const
Send data to other ranks (string version).
Definition: data_communicator.h:929
virtual bool ErrorIfFalseOnAnyRank(bool Condition) const
This function throws an error on ranks where Condition evaluates to true, if it evaluated to false on...
Definition: data_communicator.h:753
virtual int Rank() const
Get the parallel rank for this DataCommunicator.
Definition: data_communicator.h:587
void Broadcast(TObject &rBroadcastObject, const int SourceRank) const
Synchronize a buffer to the value held by the broadcasting rank.
Definition: data_communicator.h:473
virtual std::string SendRecvImpl(const std::string &rSendValues, const int SendDestination, const int SendTag, const int RecvSource, const int RecvTag) const
Exchange data with other ranks (string version).
Definition: data_communicator.h:876
virtual ~DataCommunicator()
Destructor.
Definition: data_communicator.h:368
void BroadcastImpl(TObject &rBroadcastObject, const int SourceRank) const
Synchronize a buffer to the value held by the broadcasting rank (generic version).
Definition: data_communicator.h:816
void Send(const TObject &rSendValues, const int SendDestination, const int SendTag=0) const
Exchange data with other ranks.
Definition: data_communicator.h:559
virtual bool IsDistributed() const
Check whether this DataCommunicator is aware of parallelism.
Definition: data_communicator.h:606
Definition: flags.h:58
Definition: mpi_serializer.h:33
void load(std::string const &rTag, TDataType &rObject)
Definition: serializer.h:207
void save(std::string const &rTag, std::array< TDataType, TDataSize > const &rObject)
Definition: serializer.h:545
std::string GetStringRepresentation()
Definition: stream_serializer.h:53
#define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_IMPLEMENTATION_FOR_TYPE(...)
Definition: data_communicator.h:299
#define KRATOS_DATA_COMMUNICATOR_DEBUG_SIZE_CHECK(Size1, Size2, CheckedFunction)
Definition: data_communicator.h:30
#define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_PUBLIC_INTERFACE_FOR_TYPE(...)
Definition: data_communicator.h:288
#define KRATOS_BASE_DATA_COMMUNICATOR_DECLARE_ALLREDUCE_LOC_INTERFACE_FOR_TYPE(...)
Definition: data_communicator.h:111
#define KRATOS_ERROR
Definition: exception.h:161
#define KRATOS_ERROR_IF(conditional)
Definition: exception.h:162
REF: G. R. Cowper, GAUSSIAN QUADRATURE FORMULAS FOR TRIANGLES.
Definition: mesh_condition.cpp:21
std::istream & operator>>(std::istream &rIStream, LinearMasterSlaveConstraint &rThis)
input stream function
namespace KRATOS_DEPRECATED_MESSAGE("Please use std::filesystem directly") filesystem
Definition: kratos_filesystem.h:33
std::ostream & operator<<(std::ostream &rOStream, const LinearMasterSlaveConstraint &rThis)
output stream function
Definition: linear_master_slave_constraint.h:432
tuple const
Definition: ode_solve.py:403
namespace
Definition: array_1d.h:793