Actual source code: iscomp.c
petsc-3.11.4 2019-09-28
2: #include <petsc/private/isimpl.h>
4: /*@
5: ISEqual - Compares if two index sets have the same set of indices.
7: Collective on IS
9: Input Parameters:
10: . is1, is2 - The index sets being compared
12: Output Parameters:
13: . flg - output flag, either PETSC_TRUE (if both index sets have the
14: same indices), or PETSC_FALSE if the index sets differ by size
15: or by the set of indices)
17: Level: intermediate
19: Note:
20: This routine sorts the contents of the index sets before
21: the comparision is made, so the order of the indices on a processor is immaterial.
23: Each processor has to have the same indices in the two sets, for example,
24: $ Processor
25: $ 0 1
26: $ is1 = {0, 1} {2, 3}
27: $ is2 = {2, 3} {0, 1}
28: will return false.
30: Concepts: index sets^equal
31: Concepts: IS^equal
33: .seealso: 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;
42: PetscMPIInt mflg;
49: if (is1 == is2) {
50: *flg = PETSC_TRUE;
51: return(0);
52: }
54: MPI_Comm_compare(PetscObjectComm((PetscObject)is1),PetscObjectComm((PetscObject)is2),&mflg);
55: if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
56: *flg = PETSC_FALSE;
57: return(0);
58: }
60: ISGetSize(is1,&sz1);
61: ISGetSize(is2,&sz2);
62: if (sz1 != sz2) *flg = PETSC_FALSE;
63: else {
64: ISGetLocalSize(is1,&sz1);
65: ISGetLocalSize(is2,&sz2);
67: if (sz1 != sz2) flag = PETSC_FALSE;
68: else {
69: ISGetIndices(is1,&ptr1);
70: ISGetIndices(is2,&ptr2);
72: PetscMalloc1(sz1,&a1);
73: PetscMalloc1(sz2,&a2);
75: PetscMemcpy(a1,ptr1,sz1*sizeof(PetscInt));
76: PetscMemcpy(a2,ptr2,sz2*sizeof(PetscInt));
78: PetscSortInt(sz1,a1);
79: PetscSortInt(sz2,a2);
80: PetscMemcmp(a1,a2,sz1*sizeof(PetscInt),&flag);
82: ISRestoreIndices(is1,&ptr1);
83: ISRestoreIndices(is2,&ptr2);
85: PetscFree(a1);
86: PetscFree(a2);
87: }
88: PetscObjectGetComm((PetscObject)is1,&comm);
89: MPIU_Allreduce(&flag,flg,1,MPIU_BOOL,MPI_MIN,comm);
90: }
91: return(0);
92: }
94: /*@
95: ISEqualUnsorted - Compares if two index sets have the same set of indices.
97: Collective on IS
99: Input Parameters:
100: . is1, is2 - The index sets being compared
102: Output Parameters:
103: . flg - output flag, either PETSC_TRUE (if both index sets have the
104: same indices), or PETSC_FALSE if the index sets differ by size
105: or by the set of indices)
107: Level: intermediate
109: Note:
110: This routine does NOT sort the contents of the index sets before
111: the comparision is made.
113: Concepts: index sets^equal
114: Concepts: IS^equal
116: .seealso: 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;
125: PetscMPIInt mflg;
132: if (is1 == is2) {
133: *flg = PETSC_TRUE;
134: return(0);
135: }
137: MPI_Comm_compare(PetscObjectComm((PetscObject)is1),PetscObjectComm((PetscObject)is2),&mflg);
138: if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
139: *flg = PETSC_FALSE;
140: return(0);
141: }
143: ISGetSize(is1,&sz1);
144: ISGetSize(is2,&sz2);
145: if (sz1 != sz2) *flg = PETSC_FALSE;
146: else {
147: ISGetLocalSize(is1,&sz1);
148: ISGetLocalSize(is2,&sz2);
150: if (sz1 != sz2) flag = PETSC_FALSE;
151: else {
152: ISGetIndices(is1,&ptr1);
153: ISGetIndices(is2,&ptr2);
155: PetscMemcmp(ptr1,ptr2,sz1*sizeof(PetscInt),&flag);
157: ISRestoreIndices(is1,&ptr1);
158: ISRestoreIndices(is2,&ptr2);
159: }
160: PetscObjectGetComm((PetscObject)is1,&comm);
161: MPIU_Allreduce(&flag,flg,1,MPIU_BOOL,MPI_MIN,comm);
162: }
163: return(0);
164: }