1 /* 2 A star forest (SF) describes a communication pattern 3 */ 4 #if !defined(__PETSCSF_H) 5 #define __PETSCSF_H 6 #include "petscsys.h" 7 PETSC_EXTERN_CXX_BEGIN 8 9 extern PetscClassId PETSCSF_CLASSID; 10 11 /*S 12 PetscSF - PETSc object for communication using star forests 13 14 Level: intermediate 15 16 Concepts: star forest 17 18 .seealso: PetscSFCreate() 19 S*/ 20 typedef struct _p_PetscSF* PetscSF; 21 22 /*S 23 PetscSFNode - specifier of owner and index 24 25 Level: beginner 26 27 Concepts: indexing, stride, distribution 28 29 .seealso: PetscSFSetGraph() 30 S*/ 31 typedef struct { 32 PetscInt rank; /* Rank of owner */ 33 PetscInt index; /* Index of node on rank */ 34 } PetscSFNode; 35 36 /*E 37 PetscSFSynchronizationType - Type of synchronization for PetscSF 38 39 $ PETSCSF_SYNCHRONIZATION_FENCE - simplest model, synchronizing across communicator 40 $ PETSCSF_SYNCHRONIZATION_LOCK - passive model, less synchronous, requires less setup than PETSCSF_SYNCHRONIZATION_ACTIVE, but may require more handshakes 41 $ PETSCSF_SYNCHRONIZATION_ACTIVE - active model, provides most information to MPI implementation, needs to construct 2-way process groups (more setup than PETSCSF_SYNCHRONIZATION_LOCK) 42 43 Level: beginner 44 45 .seealso: PetscSFSetSynchronizationType() 46 E*/ 47 typedef enum {PETSCSF_SYNCHRONIZATION_FENCE,PETSCSF_SYNCHRONIZATION_LOCK,PETSCSF_SYNCHRONIZATION_ACTIVE} PetscSFSynchronizationType; 48 extern const char *const PetscSFSynchronizationTypes[]; 49 50 extern PetscErrorCode PetscSFInitializePackage(const char*); 51 extern PetscErrorCode PetscSFFinalizePackage(void); 52 extern PetscErrorCode PetscSFCreate(MPI_Comm comm,PetscSF*); 53 extern PetscErrorCode PetscSFDestroy(PetscSF*); 54 extern PetscErrorCode PetscSFView(PetscSF,PetscViewer); 55 extern PetscErrorCode PetscSFSetFromOptions(PetscSF); 56 extern PetscErrorCode PetscSFSetSynchronizationType(PetscSF,PetscSFSynchronizationType); 57 extern PetscErrorCode PetscSFSetRankOrder(PetscSF,PetscBool); 58 extern PetscErrorCode PetscSFSetGraph(PetscSF,PetscInt nroots,PetscInt nleaves,const PetscInt *ilocal,PetscCopyMode modelocal,const PetscSFNode *remote,PetscCopyMode moderemote); 59 extern PetscErrorCode PetscSFCreateArray(PetscSF,MPI_Datatype,void*,void*); 60 extern PetscErrorCode PetscSFDestroyArray(PetscSF,MPI_Datatype,void*,void*); 61 extern PetscErrorCode PetscSFReset(PetscSF); 62 extern PetscErrorCode PetscSFGetRanks(PetscSF,PetscInt*,const PetscInt**,const PetscInt**,const PetscMPIInt**,const PetscMPIInt**); 63 extern PetscErrorCode PetscSFGetDataTypes(PetscSF,MPI_Datatype,const MPI_Datatype**,const MPI_Datatype**); 64 extern PetscErrorCode PetscSFGetWindow(PetscSF,MPI_Datatype,void*,PetscBool,PetscMPIInt,PetscMPIInt,PetscMPIInt,MPI_Win*); 65 extern PetscErrorCode PetscSFFindWindow(PetscSF,MPI_Datatype,const void*,MPI_Win*); 66 extern PetscErrorCode PetscSFRestoreWindow(PetscSF,MPI_Datatype,const void*,PetscBool,PetscMPIInt,MPI_Win*); 67 extern PetscErrorCode PetscSFGetGroups(PetscSF,MPI_Group*,MPI_Group*); 68 extern PetscErrorCode PetscSFGetMultiSF(PetscSF,PetscSF*); 69 70 /* broadcasts rootdata to leafdata */ 71 extern PetscErrorCode PetscSFBcastBegin(PetscSF,MPI_Datatype,const void *rootdata,void *leafdata); 72 extern PetscErrorCode PetscSFBcastEnd(PetscSF,MPI_Datatype,const void *rootdata,void *leafdata); 73 /* Reduce leafdata into rootdata using provided operation */ 74 extern PetscErrorCode PetscSFReduceBegin(PetscSF,MPI_Datatype,const void *leafdata,void *rootdata,MPI_Op); 75 extern PetscErrorCode PetscSFReduceEnd(PetscSF,MPI_Datatype,const void *leafdata,void *rootdata,MPI_Op); 76 /* Atomically modifies (using provided operation) rootdata using leafdata from each leaf, value at root at time of modification is returned in leafupdate. */ 77 extern PetscErrorCode PetscSFFetchAndOpBegin(PetscSF,MPI_Datatype,void *rootdata,const void *leafdata,void *leafupdate,MPI_Op); 78 extern PetscErrorCode PetscSFFetchAndOpEnd(PetscSF,MPI_Datatype,void *rootdata,const void *leafdata,void *leafupdate,MPI_Op); 79 /* Compute the degree of every root vertex (number of leaves in its star) */ 80 extern PetscErrorCode PetscSFComputeDegreeBegin(PetscSF,const PetscInt **degree); 81 extern PetscErrorCode PetscSFComputeDegreeEnd(PetscSF,const PetscInt **degree); 82 /* Concatenate data from all leaves into roots */ 83 extern PetscErrorCode PetscSFGatherBegin(PetscSF,MPI_Datatype,const void *leafdata,void *multirootdata); 84 extern PetscErrorCode PetscSFGatherEnd(PetscSF,MPI_Datatype,const void *leafdata,void *multirootdata); 85 /* Distribute distinct values to each leaf from roots */ 86 extern PetscErrorCode PetscSFScatterBegin(PetscSF,MPI_Datatype,const void *multirootdata,void *leafdata); 87 extern PetscErrorCode PetscSFScatterEnd(PetscSF,MPI_Datatype,const void *multirootdata,void *leafdata); 88 89 PETSC_EXTERN_CXX_END 90 91 #endif 92