1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365 366 367 368 369 370 371 372 373 374 375 376 377 378 379 380 381 382 383 384 385 386 387 388 389 390 391 392 393 394 395 396 397 398 399 400 401 402 403 404 405 406 407 408 409 410 411 412 413 414 415 416 417 418 419 420 421 422 423 424 425 426 427 428 429 430 431 432 433 434 435 436 437 438 439 440 441 442 443 444 445 446 447 448 449 450 451 452 453 454 455 456 457 458 459 460 461 462 463 464 465 466 467 468 469 470 471 472 473 474 475 476 477 478 479 480 481 482 483 484 485 486 487 488 489 490 491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506 507 508 509 510 511 512 513 514
|
// SPDX-FileCopyrightText: Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
// SPDX-License-Identifier: BSD-3-Clause
/**
* @class vtkMPIController
* @brief Process communication using MPI
*
* vtkMPIController is a concrete class which implements the
* abstract multi-process control methods defined in vtkMultiProcessController
* using MPI (Message Passing Interface)
* cf. Using MPI / Portable Parallel Programming with the Message-Passing
* Interface, Gropp et al, MIT Press.
* It also provide functionality specific to MPI and not present in
* vtkMultiProcessController.
* Before any MPI communication can occur Initialize() must be called
* by all processes. It is required to be called once, controllers
* created after this need not call Initialize().
* At the end of the program Finalize() must be called by all
* processes.
*
* The use of user-defined communicators are supported with the
* CreateSubController method. Note that a duplicate of the user defined
* communicator is used for internal communications (RMIs). This communicator
* has the same properties as the user one except that it has a new context
* which prevents the two communicators from interfering with each other.
*
* @sa
* vtkOutputPort vtkInputPort vtkMultiProcessController
* vtkMPICommunicator vtkProcessGroup
*/
#ifndef vtkMPIController_h
#define vtkMPIController_h
#include "vtkMultiProcessController.h"
#include "vtkParallelMPIModule.h" // For export macro
// Do not remove this header file. This class contains methods
// which take arguments defined in vtkMPICommunicator.h by
// reference.
#include "vtkMPICommunicator.h" // Needed for direct access to communicator
VTK_ABI_NAMESPACE_BEGIN
class vtkIntArray;
class VTKPARALLELMPI_EXPORT vtkMPIController : public vtkMultiProcessController
{
public:
static vtkMPIController* New();
vtkTypeMacro(vtkMPIController, vtkMultiProcessController);
void PrintSelf(ostream& os, vtkIndent indent) override;
/**
* This method is for setting up the processes.
* It needs to be called only once during program execution.
* Calling it more than once will have no effect. Controllers
* created after this call will be initialized automatically
* (i.e. they will have the proper LocalProcessId and NumberOfProcesses).
* The addresses of argc and argv should be passed to this method
* otherwise command line arguments will not be correct (because
* usually MPI implementations add their own arguments during
* startup).
*/
void Initialize(int* argc, char*** argv) override { this->Initialize(argc, argv, 0); }
void Initialize(
int* vtkNotUsed(argc), char*** vtkNotUsed(argv), int initializedExternally) override;
/**
* Same as Initialize(0, 0, 1). Mainly for calling from wrapped languages.
*/
virtual void Initialize();
/**
* This method is for cleaning up and has to be called before
* the end of the program if MPI was initialized with
* Initialize()
*/
void Finalize() override { this->Finalize(0); }
void Finalize(int finalizedExternally) override;
/**
* Execute the SingleMethod (as define by SetSingleMethod) using
* this->NumberOfProcesses processes.
*/
void SingleMethodExecute() override;
/**
* Execute the MultipleMethods (as define by calling SetMultipleMethod
* for each of the required this->NumberOfProcesses methods) using
* this->NumberOfProcesses processes.
*/
void MultipleMethodExecute() override;
/**
* This method can be used to tell the controller to create
* a special output window in which all messages are preceded
* by the process id.
*/
void CreateOutputWindow() override;
/**
* Given an MPI error code, return a string which contains
* an error message. This string has to be freed by the user.
*/
static char* ErrorString(int err);
/**
* MPIController uses this communicator in all sends and
* receives. By default, MPI_COMM_WORLD is used.
* THIS SHOULD ONLY BE CALLED ON THE PROCESSES INCLUDED
* IN THE COMMUNICATOR. FOR EXAMPLE, IF THE COMMUNICATOR
* CONTAINS PROCESSES 0 AND 1, INVOKING THIS METHOD ON
* ANY OTHER PROCESS WILL CAUSE AN MPI ERROR AND POSSIBLY
* LEAD TO A CRASH.
*/
void SetCommunicator(vtkMPICommunicator* comm);
vtkMPIController* CreateSubController(vtkProcessGroup* group) override;
vtkMPIController* PartitionController(int localColor, int localKey) override;
///@{
/**
* This method sends data to another process (non-blocking).
* Tag eliminates ambiguity when multiple sends or receives
* exist in the same process. The last argument,
* vtkMPICommunicator::Request& req can later be used (with
* req.Test() ) to test the success of the message. Return values
* are 1 for success and 0 otherwise.
* Note: These methods delegate to the communicator
*/
int NoBlockSend(
const int* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const unsigned long* data, int length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(
const char* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const unsigned char* data, int length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(
const float* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(
const double* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const vtkTypeInt64* data, int length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const int* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const unsigned long* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const char* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const unsigned char* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const float* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const double* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
int NoBlockSend(const vtkTypeInt64* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, remoteProcessId, tag, req);
}
///@}
/**
* Variant that permits dynamic type sends, like those create by MPI_Type_create_subarray
*/
int NoBlockSend(const void* data, vtkTypeInt64 length, MPI_Datatype mpiType, int remoteProcessId,
int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockSend(data, length, mpiType, remoteProcessId, tag, req);
}
///@{
/**
* This method receives data from a corresponding send (non-blocking).
* The last argument,
* vtkMPICommunicator::Request& req can later be used (with
* req.Test() ) to test the success of the message. Return values are
* 1 for success and 0 otherwise.
* Note: These methods delegate to the communicator
*/
int NoBlockReceive(
int* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(
unsigned long* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(
char* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(
unsigned char* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(
float* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(
double* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(
vtkTypeInt64* data, int length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(
int* data, vtkTypeInt64 length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(unsigned long* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(
char* data, vtkTypeInt64 length, int remoteProcessId, int tag, vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(unsigned char* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(float* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(double* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
int NoBlockReceive(vtkTypeInt64* data, vtkTypeInt64 length, int remoteProcessId, int tag,
vtkMPICommunicator::Request& req)
{
return ((vtkMPICommunicator*)this->Communicator)
->NoBlockReceive(data, length, remoteProcessId, tag, req);
}
///@}
///@{
/**
* Nonblocking test for a message. Inputs are: source -- the source rank
* or ANY_SOURCE; tag -- the tag value. Outputs are:
* flag -- True if a message matches; actualSource -- the rank
* sending the message (useful if ANY_SOURCE is used) if flag is True
* and actualSource isn't nullptr; size -- the length of the message in
* bytes if flag is true (only set if size isn't nullptr). The return
* value is 1 for success and 0 otherwise.
* Note: These methods delegate to the communicator
*/
int Iprobe(int source, int tag, int* flag, int* actualSource)
{
return ((vtkMPICommunicator*)this->Communicator)->Iprobe(source, tag, flag, actualSource);
}
int Iprobe(int source, int tag, int* flag, int* actualSource, int* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)
->Iprobe(source, tag, flag, actualSource, type, size);
}
int Iprobe(int source, int tag, int* flag, int* actualSource, unsigned long* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)
->Iprobe(source, tag, flag, actualSource, type, size);
}
int Iprobe(int source, int tag, int* flag, int* actualSource, const char* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)
->Iprobe(source, tag, flag, actualSource, type, size);
}
int Iprobe(int source, int tag, int* flag, int* actualSource, float* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)
->Iprobe(source, tag, flag, actualSource, type, size);
}
int Iprobe(int source, int tag, int* flag, int* actualSource, double* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)
->Iprobe(source, tag, flag, actualSource, type, size);
}
///@}
/**
* This controller does have probing capability
*/
bool CanProbe() override { return ((vtkMPICommunicator*)this->Communicator)->CanProbe(); }
///@{
/**
* Blocking test for a message. Inputs are: source -- the source rank
* or ANY_SOURCE; tag -- the tag value. Outputs are:
* actualSource -- the rank sending the message (useful if ANY_SOURCE is used)
* if actualSource isn't nullptr; size -- the length of the message in
* bytes if flag is true (only set if size isn't nullptr). The return
* value is 1 for success and 0 otherwise.
* Note: These methods delegate to the communicator
*/
int Probe(int source, int tag, int* actualSource) override
{
return ((vtkMPICommunicator*)this->Communicator)->Probe(source, tag, actualSource);
}
int Probe(int source, int tag, int* actualSource, int* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)->Probe(source, tag, actualSource, type, size);
}
int Probe(int source, int tag, int* actualSource, unsigned long* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)->Probe(source, tag, actualSource, type, size);
}
int Probe(int source, int tag, int* actualSource, const char* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)->Probe(source, tag, actualSource, type, size);
}
int Probe(int source, int tag, int* actualSource, float* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)->Probe(source, tag, actualSource, type, size);
}
int Probe(int source, int tag, int* actualSource, double* type, int* size)
{
return ((vtkMPICommunicator*)this->Communicator)->Probe(source, tag, actualSource, type, size);
}
///@}
/**
* Given the request objects of a set of non-blocking operations
* (send and/or receive) this method blocks until all requests are complete.
* Note: This method delegates to the communicator
*/
int WaitAll(int count, vtkMPICommunicator::Request requests[])
{
return ((vtkMPICommunicator*)this->Communicator)->WaitAll(count, requests);
}
/**
* Blocks until *one* of the specified requests in the given request array
* completes. Upon return, the index in the array of the completed request
* object is returned through the argument list.
* Note: this method delegates to the communicator
*/
int WaitAny(int count, vtkMPICommunicator::Request requests[], int& idx)
VTK_SIZEHINT(requests, count)
{
return ((vtkMPICommunicator*)this->Communicator)->WaitAny(count, requests, idx);
}
/**
* Blocks until *one or more* of the specified requests in the given request
* request array completes. Upon return, the list of handles that have
* completed is stored in the completed vtkIntArray.
*/
int WaitSome(int count, vtkMPICommunicator::Request requests[], vtkIntArray* completed)
VTK_SIZEHINT(requests, count);
/**
* Returns true iff *all* of the communication request objects are complete.
*/
bool TestAll(int count, vtkMPICommunicator::Request requests[]);
/**
* Returns true iff at least *one* of the communication request objects is
* complete. The index of the completed request, w.r.t. the requests array, is
* reflected in the out parameter idx. Otherwise, if none of the communication
* requests are complete false is returned.
*/
bool TestAny(int count, vtkMPICommunicator::Request requests[], int& idx)
VTK_SIZEHINT(requests, count);
/**
* Return true iff *one or more* of the communicator request objects is
* complete. The indices of the completed requests, w.r.t. the requests array,
* are given in the completed user-supplied vtkIntArray.
*/
bool TestSome(int count, vtkMPICommunicator::Request requests[], vtkIntArray* completed)
VTK_SIZEHINT(requests, count);
static const char* GetProcessorName();
/**
* When set to 1, TriggerRMI uses Ssend() instead of Send() calls.
* Off (0) by default.
*/
static void SetUseSsendForRMI(int use_send)
{
vtkMPIController::UseSsendForRMI = (use_send != 0) ? 1 : 0;
}
static int GetUseSsendForRMI() { return vtkMPIController::UseSsendForRMI; }
protected:
vtkMPIController();
~vtkMPIController() override;
// Set the communicator to comm and call InitializeNumberOfProcesses()
void InitializeCommunicator(vtkMPICommunicator* comm);
// Duplicate the current communicator, creating RMICommunicator
void InitializeRMICommunicator();
/**
* Implementation for TriggerRMI() provides subclasses an opportunity to
* modify the behaviour eg. MPIController provides ability to use Ssend
* instead of Send.
*/
void TriggerRMIInternal(
int remoteProcessId, void* arg, int argLength, int rmiTag, bool propagate) override;
// MPI communicator created when Initialize() called.
// This is a copy of MPI_COMM_WORLD but uses a new
// context, i.e. even if the tags are the same, the
// RMI messages will not interfere with user level
// messages.
static vtkMPICommunicator* WorldRMICommunicator;
friend class vtkMPIOutputWindow;
// Initialize only once.
static int Initialized;
static char ProcessorName[];
/**
* When set, TriggerRMI uses Ssend instead of Send.
*/
static int UseSsendForRMI;
private:
vtkMPIController(const vtkMPIController&) = delete;
void operator=(const vtkMPIController&) = delete;
};
VTK_ABI_NAMESPACE_END
#endif
|