Actual source code: petscis.h
1: /*
2: An index set is a generalization of a subset of integers. Index sets
3: are used for defining scatters and gathers.
4: */
7: #include petsc.h
12: EXTERN PetscErrorCode ISInitializePackage(const char[]);
14: /*S
15: IS - Abstract PETSc object that indexing.
17: Level: beginner
19: Concepts: indexing, stride
21: .seealso: ISCreateGeneral(), ISCreateBlock(), ISCreateStride(), ISGetIndices(), ISDestroy()
22: S*/
23: typedef struct _p_IS* IS;
25: /*
26: Default index set data structures that PETSc provides.
27: */
28: typedef enum {IS_GENERAL=0,IS_STRIDE=1,IS_BLOCK = 2} ISType;
29: EXTERN PetscErrorCode ISCreateGeneral(MPI_Comm,PetscInt,const PetscInt[],IS *);
30: EXTERN PetscErrorCode ISCreateGeneralNC(MPI_Comm,PetscInt,const PetscInt[],IS *);
31: EXTERN PetscErrorCode ISCreateGeneralWithArray(MPI_Comm,PetscInt,PetscInt[],IS *);
32: EXTERN PetscErrorCode ISCreateBlock(MPI_Comm,PetscInt,PetscInt,const PetscInt[],IS *);
33: EXTERN PetscErrorCode ISCreateStride(MPI_Comm,PetscInt,PetscInt,PetscInt,IS *);
35: EXTERN PetscErrorCode ISDestroy(IS);
37: EXTERN PetscErrorCode ISSetPermutation(IS);
38: EXTERN PetscErrorCode ISPermutation(IS,PetscTruth*);
39: EXTERN PetscErrorCode ISSetIdentity(IS);
40: EXTERN PetscErrorCode ISIdentity(IS,PetscTruth*);
42: EXTERN PetscErrorCode ISGetIndices(IS,PetscInt *[]);
43: EXTERN PetscErrorCode ISRestoreIndices(IS,PetscInt *[]);
44: EXTERN PetscErrorCode ISGetSize(IS,PetscInt *);
45: EXTERN PetscErrorCode ISGetLocalSize(IS,PetscInt *);
46: EXTERN PetscErrorCode ISInvertPermutation(IS,PetscInt,IS*);
47: EXTERN PetscErrorCode ISView(IS,PetscViewer);
48: EXTERN PetscErrorCode ISEqual(IS,IS,PetscTruth *);
49: EXTERN PetscErrorCode ISSort(IS);
50: EXTERN PetscErrorCode ISSorted(IS,PetscTruth *);
51: EXTERN PetscErrorCode ISDifference(IS,IS,IS*);
52: EXTERN PetscErrorCode ISSum(IS,IS,IS*);
53: EXTERN PetscErrorCode ISExpand(IS,IS,IS*);
55: EXTERN PetscErrorCode ISBlock(IS,PetscTruth*);
56: EXTERN PetscErrorCode ISBlockGetIndices(IS,PetscInt *[]);
57: EXTERN PetscErrorCode ISBlockRestoreIndices(IS,PetscInt *[]);
58: EXTERN PetscErrorCode ISBlockGetSize(IS,PetscInt *);
59: EXTERN PetscErrorCode ISBlockGetBlockSize(IS,PetscInt *);
61: EXTERN PetscErrorCode ISStride(IS,PetscTruth*);
62: EXTERN PetscErrorCode ISStrideGetInfo(IS,PetscInt *,PetscInt*);
64: EXTERN PetscErrorCode ISStrideToGeneral(IS);
66: EXTERN PetscErrorCode ISDuplicate(IS,IS*);
67: EXTERN PetscErrorCode ISAllGather(IS,IS*);
68: EXTERN PetscErrorCode ISAllGatherIndices(MPI_Comm,PetscInt,const PetscInt[],PetscInt*,PetscInt*[]);
70: /* --------------------------------------------------------------------------*/
73: /*S
74: ISLocalToGlobalMapping - mappings from an arbitrary
75: local ordering from 0 to n-1 to a global PETSc ordering
76: used by a vector or matrix.
78: Level: intermediate
80: Note: mapping from Local to Global is scalable; but Global
81: to Local may not be if the range of global values represented locally
82: is very large.
84: Note: the ISLocalToGlobalMapping is actually a private object; it is included
85: here for the MACRO ISLocalToGlobalMappingApply() to allow it to be inlined since
86: it is used so often.
88: .seealso: ISLocalToGlobalMappingCreate()
89: S*/
90: struct _p_ISLocalToGlobalMapping{
91: PETSCHEADER(int);
92: PetscInt n; /* number of local indices */
93: PetscInt *indices; /* global index of each local index */
94: PetscInt globalstart; /* first global referenced in indices */
95: PetscInt globalend; /* last + 1 global referenced in indices */
96: PetscInt *globals; /* local index for each global index between start and end */
97: };
98: typedef struct _p_ISLocalToGlobalMapping* ISLocalToGlobalMapping;
100: /*E
101: ISGlobalToLocalMappingType - Indicates if missing global indices are
103: IS_GTOLM_MASK - missing global indices are replaced with -1
104: IS_GTOLM_DROP - missing global indices are dropped
106: Level: beginner
108: .seealso: ISGlobalToLocalMappingApply()
110: E*/
111: typedef enum {IS_GTOLM_MASK,IS_GTOLM_DROP} ISGlobalToLocalMappingType;
113: EXTERN PetscErrorCode ISLocalToGlobalMappingCreate(MPI_Comm,PetscInt,const PetscInt[],ISLocalToGlobalMapping*);
114: EXTERN PetscErrorCode ISLocalToGlobalMappingCreateNC(MPI_Comm,PetscInt,const PetscInt[],ISLocalToGlobalMapping*);
115: EXTERN PetscErrorCode ISLocalToGlobalMappingCreateIS(IS,ISLocalToGlobalMapping *);
116: EXTERN PetscErrorCode ISLocalToGlobalMappingView(ISLocalToGlobalMapping,PetscViewer);
117: EXTERN PetscErrorCode ISLocalToGlobalMappingDestroy(ISLocalToGlobalMapping);
118: EXTERN PetscErrorCode ISLocalToGlobalMappingApplyIS(ISLocalToGlobalMapping,IS,IS*);
119: EXTERN PetscErrorCode ISGlobalToLocalMappingApply(ISLocalToGlobalMapping,ISGlobalToLocalMappingType,PetscInt,const PetscInt[],PetscInt*,PetscInt[]);
120: EXTERN PetscErrorCode ISLocalToGlobalMappingGetSize(ISLocalToGlobalMapping,PetscInt*);
121: EXTERN PetscErrorCode ISLocalToGlobalMappingGetInfo(ISLocalToGlobalMapping,PetscInt*,PetscInt*[],PetscInt*[],PetscInt**[]);
122: EXTERN PetscErrorCode ISLocalToGlobalMappingRestoreInfo(ISLocalToGlobalMapping,PetscInt*,PetscInt*[],PetscInt*[],PetscInt**[]);
123: EXTERN PetscErrorCode ISLocalToGlobalMappingBlock(ISLocalToGlobalMapping,PetscInt,ISLocalToGlobalMapping*);
125: PETSC_STATIC_INLINE PetscErrorCode ISLocalToGlobalMappingApply(ISLocalToGlobalMapping mapping,PetscInt N,const PetscInt in[],PetscInt out[]){
126: PetscInt i,*idx = mapping->indices,Nmax = mapping->n;
127: for (i=0; i<N; i++) {
128: if (in[i] < 0) {out[i] = in[i]; continue;}
129: if (in[i] >= Nmax) SETERRQ3(PETSC_ERR_ARG_OUTOFRANGE,"Local index %D too large %D (max) at %D",in[i],Nmax,i);
130: out[i] = idx[in[i]];
131: }
132: return(0);
133: }
135: /* --------------------------------------------------------------------------*/
136: /*E
137: ISColoringType - determines if the coloring is for the entire parallel grid/graph/matrix
138: or for just the local ghosted portion
140: Level: beginner
142: $ IS_COLORING_GLOBAL - does not include the colors for ghost points, this is used when the function
143: $ is called synchronously in parallel. This requires generating a "parallel coloring".
144: $ IS_COLORING_GHOSTED - includes colors for ghost points, this is used when the function can be called
145: $ seperately on individual processes with the ghost points already filled in. Does not
146: $ require a "parallel coloring", rather each process colors its local + ghost part.
147: $ Using this can result in much less parallel communication. In the paradigm of
148: $ DAGetLocalVector() and DAGetGlobalVector() this could be called IS_COLORING_LOCAL
150: .seealso: DAGetColoring()
151: E*/
152: typedef enum {IS_COLORING_GLOBAL,IS_COLORING_GHOSTED} ISColoringType;
154: typedef unsigned PETSC_IS_COLOR_VALUE_TYPE ISColoringValue;
155: EXTERN PetscErrorCode ISAllGatherColors(MPI_Comm,PetscInt,ISColoringValue*,PetscInt*,ISColoringValue*[]);
157: /*S
158: ISColoring - sets of IS's that define a coloring
159: of the underlying indices
161: Level: intermediate
163: Notes:
164: One should not access the *is records below directly because they may not yet
165: have been created. One should use ISColoringGetIS() to make sure they are
166: created when needed.
168: .seealso: ISColoringCreate(), ISColoringGetIS(), ISColoringView(), ISColoringGetIS()
169: S*/
170: struct _n_ISColoring {
171: PetscInt refct;
172: PetscInt n; /* number of colors */
173: IS *is; /* for each color indicates columns */
174: MPI_Comm comm;
175: ISColoringValue *colors; /* for each column indicates color */
176: PetscInt N; /* number of columns */
177: ISColoringType ctype;
178: };
179: typedef struct _n_ISColoring* ISColoring;
181: EXTERN PetscErrorCode ISColoringCreate(MPI_Comm,PetscInt,PetscInt,const ISColoringValue[],ISColoring*);
182: EXTERN PetscErrorCode ISColoringDestroy(ISColoring);
183: EXTERN PetscErrorCode ISColoringView(ISColoring,PetscViewer);
184: EXTERN PetscErrorCode ISColoringGetIS(ISColoring,PetscInt*,IS*[]);
185: EXTERN PetscErrorCode ISColoringRestoreIS(ISColoring,IS*[]);
186: #define ISColoringReference(coloring) ((coloring)->refct++,0)
187: #define ISColoringSetType(coloring,type) ((coloring)->ctype = type,0)
189: /* --------------------------------------------------------------------------*/
191: EXTERN PetscErrorCode ISPartitioningToNumbering(IS,IS*);
192: EXTERN PetscErrorCode ISPartitioningCount(IS,PetscInt[]);
194: EXTERN PetscErrorCode ISCompressIndicesGeneral(PetscInt,PetscInt,PetscInt,const IS[],IS[]);
195: EXTERN PetscErrorCode ISCompressIndicesSorted(PetscInt,PetscInt,PetscInt,const IS[],IS[]);
196: EXTERN PetscErrorCode ISExpandIndicesGeneral(PetscInt,PetscInt,PetscInt,const IS[],IS[]);
199: #endif