Actual source code: sfimpl.h
petsc-3.11.1 2019-04-12
1: #if !defined(_PETSCSFIMPL_H)
2: #define _PETSCSFIMPL_H
4: #include <petscsf.h>
5: #include <petsc/private/petscimpl.h>
6: #include <petscviewer.h>
8: PETSC_EXTERN PetscLogEvent PETSCSF_SetGraph;
9: PETSC_EXTERN PetscLogEvent PETSCSF_SetUp;
10: PETSC_EXTERN PetscLogEvent PETSCSF_BcastBegin;
11: PETSC_EXTERN PetscLogEvent PETSCSF_BcastEnd;
12: PETSC_EXTERN PetscLogEvent PETSCSF_BcastAndOpBegin;
13: PETSC_EXTERN PetscLogEvent PETSCSF_BcastAndOpEnd;
14: PETSC_EXTERN PetscLogEvent PETSCSF_ReduceBegin;
15: PETSC_EXTERN PetscLogEvent PETSCSF_ReduceEnd;
16: PETSC_EXTERN PetscLogEvent PETSCSF_FetchAndOpBegin;
17: PETSC_EXTERN PetscLogEvent PETSCSF_FetchAndOpEnd;
19: struct _PetscSFOps {
20: PetscErrorCode (*Reset)(PetscSF);
21: PetscErrorCode (*Destroy)(PetscSF);
22: PetscErrorCode (*SetUp)(PetscSF);
23: PetscErrorCode (*SetFromOptions)(PetscOptionItems*,PetscSF);
24: PetscErrorCode (*View)(PetscSF,PetscViewer);
25: PetscErrorCode (*Duplicate)(PetscSF,PetscSFDuplicateOption,PetscSF);
26: PetscErrorCode (*BcastBegin)(PetscSF,MPI_Datatype,const void*,void*);
27: PetscErrorCode (*BcastEnd)(PetscSF,MPI_Datatype,const void*,void*);
28: PetscErrorCode (*BcastAndOpBegin)(PetscSF,MPI_Datatype,const void*,void*,MPI_Op);
29: PetscErrorCode (*BcastAndOpEnd)(PetscSF,MPI_Datatype,const void*,void*,MPI_Op);
30: PetscErrorCode (*ReduceBegin)(PetscSF,MPI_Datatype,const void*,void*,MPI_Op);
31: PetscErrorCode (*ReduceEnd)(PetscSF,MPI_Datatype,const void*,void*,MPI_Op);
32: PetscErrorCode (*FetchAndOpBegin)(PetscSF,MPI_Datatype,void*,const void*,void*,MPI_Op);
33: PetscErrorCode (*FetchAndOpEnd)(PetscSF,MPI_Datatype,void*,const void *,void *,MPI_Op);
34: PetscErrorCode (*GetLeafRanks)(PetscSF,PetscInt*,const PetscMPIInt**,const PetscInt**,const PetscInt **);
35: };
37: struct _p_PetscSF {
38: PETSCHEADER(struct _PetscSFOps);
39: PetscInt nroots; /* Number of root vertices on current process (candidates for incoming edges) */
40: PetscInt nleaves; /* Number of leaf vertices on current process (this process specifies a root for each leaf) */
41: PetscInt *mine; /* Location of leaves in leafdata arrays provided to the communication routines */
42: PetscInt *mine_alloc;
43: PetscInt minleaf,maxleaf;
44: PetscSFNode *remote; /* Remote references to roots for each local leaf */
45: PetscSFNode *remote_alloc;
46: PetscInt nranks; /* Number of ranks owning roots connected to my leaves */
47: PetscInt ndranks; /* Number of ranks in distinguished group holding roots connected to my leaves */
48: PetscMPIInt *ranks; /* List of ranks referenced by "remote" */
49: PetscInt *roffset; /* Array of length nranks+1, offset in rmine/rremote for each rank */
50: PetscInt *rmine; /* Concatenated array holding local indices referencing each remote rank */
51: PetscInt *rremote; /* Concatenated array holding remote indices referenced for each remote rank */
52: PetscBool degreeknown; /* The degree is currently known, do not have to recompute */
53: PetscInt *degree; /* Degree of each of my root vertices */
54: PetscInt *degreetmp; /* Temporary local array for computing degree */
55: PetscBool rankorder; /* Sort ranks for gather and scatter operations */
56: MPI_Group ingroup; /* Group of processes connected to my roots */
57: MPI_Group outgroup; /* Group of processes connected to my leaves */
58: PetscSF multi; /* Internal graph used to implement gather and scatter operations */
59: PetscBool graphset; /* Flag indicating that the graph has been set, required before calling communication routines */
60: PetscBool setupcalled; /* Type and communication structures have been set up */
62: void *data; /* Pointer to implementation */
63: };
65: PETSC_EXTERN PetscBool PetscSFRegisterAllCalled;
66: PETSC_EXTERN PetscErrorCode PetscSFRegisterAll(void);
68: PETSC_EXTERN PetscErrorCode MPIPetsc_Type_unwrap(MPI_Datatype,MPI_Datatype*,PetscBool*);
69: PETSC_EXTERN PetscErrorCode MPIPetsc_Type_compare(MPI_Datatype,MPI_Datatype,PetscBool*);
70: PETSC_EXTERN PetscErrorCode MPIPetsc_Type_compare_contig(MPI_Datatype,MPI_Datatype,PetscInt*);
72: #endif