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 - first index set to compare
12: - is2 - second index set to compare
14: Output Parameter:
15: . flg - output flag, either `PETSC_TRUE` (if both index sets have the
16: same indices), or `PETSC_FALSE` if the index sets differ by size
17: or by the set of indices)
19: Level: intermediate
21: Note:
22: Unlike `ISEqualUnsorted()`, this routine sorts the contents of the index sets (only within each MPI rank) before
23: the comparison is made, so the order of the indices on a processor is immaterial.
25: Each processor has to have the same indices in the two sets, for example,
26: .vb
27: Processor
28: 0 1
29: is1 = {0, 1} {2, 3}
30: is2 = {2, 3} {0, 1}
31: .ve
32: will return false.
34: .seealso: [](sec_scatter), `IS`, `ISEqualUnsorted()`
35: @*/
36: PetscErrorCode ISEqual(IS is1, IS is2, PetscBool *flg)
37: {
38: PetscInt sz1, sz2, *a1, *a2;
39: const PetscInt *ptr1, *ptr2;
40: PetscBool flag;
41: MPI_Comm comm;
42: PetscMPIInt mflg;
44: PetscFunctionBegin;
49: if (is1 == is2) {
50: *flg = PETSC_TRUE;
51: PetscFunctionReturn(PETSC_SUCCESS);
52: }
54: PetscCallMPI(MPI_Comm_compare(PetscObjectComm((PetscObject)is1), PetscObjectComm((PetscObject)is2), &mflg));
55: if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
56: *flg = PETSC_FALSE;
57: PetscFunctionReturn(PETSC_SUCCESS);
58: }
60: PetscCall(ISGetSize(is1, &sz1));
61: PetscCall(ISGetSize(is2, &sz2));
62: if (sz1 != sz2) *flg = PETSC_FALSE;
63: else {
64: PetscCall(ISGetLocalSize(is1, &sz1));
65: PetscCall(ISGetLocalSize(is2, &sz2));
67: if (sz1 != sz2) flag = PETSC_FALSE;
68: else {
69: PetscCall(ISGetIndices(is1, &ptr1));
70: PetscCall(ISGetIndices(is2, &ptr2));
72: PetscCall(PetscMalloc1(sz1, &a1));
73: PetscCall(PetscMalloc1(sz2, &a2));
75: PetscCall(PetscArraycpy(a1, ptr1, sz1));
76: PetscCall(PetscArraycpy(a2, ptr2, sz2));
78: PetscCall(PetscIntSortSemiOrdered(sz1, a1));
79: PetscCall(PetscIntSortSemiOrdered(sz2, a2));
80: PetscCall(PetscArraycmp(a1, a2, sz1, &flag));
82: PetscCall(ISRestoreIndices(is1, &ptr1));
83: PetscCall(ISRestoreIndices(is2, &ptr2));
85: PetscCall(PetscFree(a1));
86: PetscCall(PetscFree(a2));
87: }
88: PetscCall(PetscObjectGetComm((PetscObject)is1, &comm));
89: PetscCall(MPIU_Allreduce(&flag, flg, 1, MPIU_BOOL, MPI_MIN, comm));
90: }
91: PetscFunctionReturn(PETSC_SUCCESS);
92: }
94: /*@
95: ISEqualUnsorted - Compares if two index sets have the same indices.
97: Collective
99: Input Parameters:
100: + is1 - first index set to compare
101: - is2 - second index set to compare
103: Output Parameter:
104: . flg - output flag, either `PETSC_TRUE` (if both index sets have the
105: same indices), or `PETSC_FALSE` if the index sets differ by size
106: or by the set of indices)
108: Level: intermediate
110: Note:
111: Unlike `ISEqual()`, this routine does NOT sort the contents of the index sets before
112: the comparison is made, i.e., the order of indices is important.
114: Each MPI rank must have the same indices.
116: .seealso: [](sec_scatter), `IS`, `ISEqual()`
117: @*/
118: PetscErrorCode ISEqualUnsorted(IS is1, IS is2, PetscBool *flg)
119: {
120: PetscInt sz1, sz2;
121: const PetscInt *ptr1, *ptr2;
122: PetscBool flag;
123: MPI_Comm comm;
124: PetscMPIInt mflg;
126: PetscFunctionBegin;
131: if (is1 == is2) {
132: *flg = PETSC_TRUE;
133: PetscFunctionReturn(PETSC_SUCCESS);
134: }
136: PetscCallMPI(MPI_Comm_compare(PetscObjectComm((PetscObject)is1), PetscObjectComm((PetscObject)is2), &mflg));
137: if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
138: *flg = PETSC_FALSE;
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: PetscCall(ISGetSize(is1, &sz1));
143: PetscCall(ISGetSize(is2, &sz2));
144: if (sz1 != sz2) *flg = PETSC_FALSE;
145: else {
146: PetscCall(ISGetLocalSize(is1, &sz1));
147: PetscCall(ISGetLocalSize(is2, &sz2));
149: if (sz1 != sz2) flag = PETSC_FALSE;
150: else {
151: PetscCall(ISGetIndices(is1, &ptr1));
152: PetscCall(ISGetIndices(is2, &ptr2));
154: PetscCall(PetscArraycmp(ptr1, ptr2, sz1, &flag));
156: PetscCall(ISRestoreIndices(is1, &ptr1));
157: PetscCall(ISRestoreIndices(is2, &ptr2));
158: }
159: PetscCall(PetscObjectGetComm((PetscObject)is1, &comm));
160: PetscCall(MPIU_Allreduce(&flag, flg, 1, MPIU_BOOL, MPI_MIN, comm));
161: }
162: PetscFunctionReturn(PETSC_SUCCESS);
163: }