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 on is1
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: $ Processor
26: $ 0 1
27: $ is1 = {0, 1} {2, 3}
28: $ is2 = {2, 3} {0, 1}
29: will return false.
31: .seealso: [](sec_scatter), `IS`, `ISEqualUnsorted()`
32: @*/
33: PetscErrorCode ISEqual(IS is1, IS is2, PetscBool *flg)
34: {
35: PetscInt sz1, sz2, *a1, *a2;
36: const PetscInt *ptr1, *ptr2;
37: PetscBool flag;
38: MPI_Comm comm;
39: PetscMPIInt mflg;
41: PetscFunctionBegin;
46: if (is1 == is2) {
47: *flg = PETSC_TRUE;
48: PetscFunctionReturn(PETSC_SUCCESS);
49: }
51: PetscCallMPI(MPI_Comm_compare(PetscObjectComm((PetscObject)is1), PetscObjectComm((PetscObject)is2), &mflg));
52: if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
53: *flg = PETSC_FALSE;
54: PetscFunctionReturn(PETSC_SUCCESS);
55: }
57: PetscCall(ISGetSize(is1, &sz1));
58: PetscCall(ISGetSize(is2, &sz2));
59: if (sz1 != sz2) *flg = PETSC_FALSE;
60: else {
61: PetscCall(ISGetLocalSize(is1, &sz1));
62: PetscCall(ISGetLocalSize(is2, &sz2));
64: if (sz1 != sz2) flag = PETSC_FALSE;
65: else {
66: PetscCall(ISGetIndices(is1, &ptr1));
67: PetscCall(ISGetIndices(is2, &ptr2));
69: PetscCall(PetscMalloc1(sz1, &a1));
70: PetscCall(PetscMalloc1(sz2, &a2));
72: PetscCall(PetscArraycpy(a1, ptr1, sz1));
73: PetscCall(PetscArraycpy(a2, ptr2, sz2));
75: PetscCall(PetscIntSortSemiOrdered(sz1, a1));
76: PetscCall(PetscIntSortSemiOrdered(sz2, a2));
77: PetscCall(PetscArraycmp(a1, a2, sz1, &flag));
79: PetscCall(ISRestoreIndices(is1, &ptr1));
80: PetscCall(ISRestoreIndices(is2, &ptr2));
82: PetscCall(PetscFree(a1));
83: PetscCall(PetscFree(a2));
84: }
85: PetscCall(PetscObjectGetComm((PetscObject)is1, &comm));
86: PetscCall(MPIU_Allreduce(&flag, flg, 1, MPIU_BOOL, MPI_MIN, comm));
87: }
88: PetscFunctionReturn(PETSC_SUCCESS);
89: }
91: /*@
92: ISEqualUnsorted - Compares if two index sets have the same indices.
94: Collective on is1
96: Input Parameters:
97: . is1, is2 - The index sets being compared
99: Output Parameters:
100: . flg - output flag, either `PETSC_TRUE` (if both index sets have the
101: same indices), or `PETSC_FALSE` if the index sets differ by size
102: or by the set of indices)
104: Level: intermediate
106: Note:
107: Unlike ISEqual(), this routine does NOT sort the contents of the index sets before
108: the comparison is made, i.e., the order of indices is important.
110: Each MPI rank must have the same indices.
112: .seealso: [](sec_scatter), `IS`, `ISEqual()`
113: @*/
114: PetscErrorCode ISEqualUnsorted(IS is1, IS is2, PetscBool *flg)
115: {
116: PetscInt sz1, sz2;
117: const PetscInt *ptr1, *ptr2;
118: PetscBool flag;
119: MPI_Comm comm;
120: PetscMPIInt mflg;
122: PetscFunctionBegin;
127: if (is1 == is2) {
128: *flg = PETSC_TRUE;
129: PetscFunctionReturn(PETSC_SUCCESS);
130: }
132: PetscCallMPI(MPI_Comm_compare(PetscObjectComm((PetscObject)is1), PetscObjectComm((PetscObject)is2), &mflg));
133: if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
134: *flg = PETSC_FALSE;
135: PetscFunctionReturn(PETSC_SUCCESS);
136: }
138: PetscCall(ISGetSize(is1, &sz1));
139: PetscCall(ISGetSize(is2, &sz2));
140: if (sz1 != sz2) *flg = PETSC_FALSE;
141: else {
142: PetscCall(ISGetLocalSize(is1, &sz1));
143: PetscCall(ISGetLocalSize(is2, &sz2));
145: if (sz1 != sz2) flag = PETSC_FALSE;
146: else {
147: PetscCall(ISGetIndices(is1, &ptr1));
148: PetscCall(ISGetIndices(is2, &ptr2));
150: PetscCall(PetscArraycmp(ptr1, ptr2, sz1, &flag));
152: PetscCall(ISRestoreIndices(is1, &ptr1));
153: PetscCall(ISRestoreIndices(is2, &ptr2));
154: }
155: PetscCall(PetscObjectGetComm((PetscObject)is1, &comm));
156: PetscCall(MPIU_Allreduce(&flag, flg, 1, MPIU_BOOL, MPI_MIN, comm));
157: }
158: PetscFunctionReturn(PETSC_SUCCESS);
159: }