Actual source code: iscomp.c
2: #include <petsc/private/isimpl.h>
3: #include <petscviewer.h>
5: /*@
6: ISEqual - Compares if two index sets have the same set of indices.
8: Collective
10: Input Parameters:
11: . is1, is2 - The index sets being compared
13: Output Parameters:
14: . flg - output flag, either `PETSC_TRUE` (if both index sets have the
15: same indices), or `PETSC_FALSE` if the index sets differ by size
16: or by the set of indices)
18: Level: intermediate
20: Note:
21: Unlike `ISEqualUnsorted()`, this routine sorts the contents of the index sets (only within each MPI rank) before
22: the comparison is made, so the order of the indices on a processor is immaterial.
24: Each processor has to have the same indices in the two sets, for example,
25: .vb
26: Processor
27: 0 1
28: is1 = {0, 1} {2, 3}
29: is2 = {2, 3} {0, 1}
30: .ve
31: will return false.
33: .seealso: [](sec_scatter), `IS`, `ISEqualUnsorted()`
34: @*/
35: PetscErrorCode ISEqual(IS is1, IS is2, PetscBool *flg)
36: {
37: PetscInt sz1, sz2, *a1, *a2;
38: const PetscInt *ptr1, *ptr2;
39: PetscBool flag;
40: MPI_Comm comm;
41: PetscMPIInt mflg;
43: PetscFunctionBegin;
48: if (is1 == is2) {
49: *flg = PETSC_TRUE;
50: PetscFunctionReturn(PETSC_SUCCESS);
51: }
53: PetscCallMPI(MPI_Comm_compare(PetscObjectComm((PetscObject)is1), PetscObjectComm((PetscObject)is2), &mflg));
54: if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
55: *flg = PETSC_FALSE;
56: PetscFunctionReturn(PETSC_SUCCESS);
57: }
59: PetscCall(ISGetSize(is1, &sz1));
60: PetscCall(ISGetSize(is2, &sz2));
61: if (sz1 != sz2) *flg = PETSC_FALSE;
62: else {
63: PetscCall(ISGetLocalSize(is1, &sz1));
64: PetscCall(ISGetLocalSize(is2, &sz2));
66: if (sz1 != sz2) flag = PETSC_FALSE;
67: else {
68: PetscCall(ISGetIndices(is1, &ptr1));
69: PetscCall(ISGetIndices(is2, &ptr2));
71: PetscCall(PetscMalloc1(sz1, &a1));
72: PetscCall(PetscMalloc1(sz2, &a2));
74: PetscCall(PetscArraycpy(a1, ptr1, sz1));
75: PetscCall(PetscArraycpy(a2, ptr2, sz2));
77: PetscCall(PetscIntSortSemiOrdered(sz1, a1));
78: PetscCall(PetscIntSortSemiOrdered(sz2, a2));
79: PetscCall(PetscArraycmp(a1, a2, sz1, &flag));
81: PetscCall(ISRestoreIndices(is1, &ptr1));
82: PetscCall(ISRestoreIndices(is2, &ptr2));
84: PetscCall(PetscFree(a1));
85: PetscCall(PetscFree(a2));
86: }
87: PetscCall(PetscObjectGetComm((PetscObject)is1, &comm));
88: PetscCall(MPIU_Allreduce(&flag, flg, 1, MPIU_BOOL, MPI_MIN, comm));
89: }
90: PetscFunctionReturn(PETSC_SUCCESS);
91: }
93: /*@
94: ISEqualUnsorted - Compares if two index sets have the same indices.
96: Collective
98: Input Parameters:
99: . is1, is2 - The index sets being compared
101: Output Parameters:
102: . flg - output flag, either `PETSC_TRUE` (if both index sets have the
103: same indices), or `PETSC_FALSE` if the index sets differ by size
104: or by the set of indices)
106: Level: intermediate
108: Note:
109: Unlike `ISEqual()`, this routine does NOT sort the contents of the index sets before
110: the comparison is made, i.e., the order of indices is important.
112: Each MPI rank must have the same indices.
114: .seealso: [](sec_scatter), `IS`, `ISEqual()`
115: @*/
116: PetscErrorCode ISEqualUnsorted(IS is1, IS is2, PetscBool *flg)
117: {
118: PetscInt sz1, sz2;
119: const PetscInt *ptr1, *ptr2;
120: PetscBool flag;
121: MPI_Comm comm;
122: PetscMPIInt mflg;
124: PetscFunctionBegin;
129: if (is1 == is2) {
130: *flg = PETSC_TRUE;
131: PetscFunctionReturn(PETSC_SUCCESS);
132: }
134: PetscCallMPI(MPI_Comm_compare(PetscObjectComm((PetscObject)is1), PetscObjectComm((PetscObject)is2), &mflg));
135: if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
136: *flg = PETSC_FALSE;
137: PetscFunctionReturn(PETSC_SUCCESS);
138: }
140: PetscCall(ISGetSize(is1, &sz1));
141: PetscCall(ISGetSize(is2, &sz2));
142: if (sz1 != sz2) *flg = PETSC_FALSE;
143: else {
144: PetscCall(ISGetLocalSize(is1, &sz1));
145: PetscCall(ISGetLocalSize(is2, &sz2));
147: if (sz1 != sz2) flag = PETSC_FALSE;
148: else {
149: PetscCall(ISGetIndices(is1, &ptr1));
150: PetscCall(ISGetIndices(is2, &ptr2));
152: PetscCall(PetscArraycmp(ptr1, ptr2, sz1, &flag));
154: PetscCall(ISRestoreIndices(is1, &ptr1));
155: PetscCall(ISRestoreIndices(is2, &ptr2));
156: }
157: PetscCall(PetscObjectGetComm((PetscObject)is1, &comm));
158: PetscCall(MPIU_Allreduce(&flag, flg, 1, MPIU_BOOL, MPI_MIN, comm));
159: }
160: PetscFunctionReturn(PETSC_SUCCESS);
161: }