1 
2 #include <petsc/private/isimpl.h>    /*I "petscis.h"  I*/
3 
4 /*@
5    ISEqual  - Compares if two index sets have the same set of indices.
6 
7    Collective on IS
8 
9    Input Parameters:
10 .  is1, is2 - The index sets being compared
11 
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)
16 
17    Level: intermediate
18 
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.
22 
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.
29 
30 
31 .seealso: ISEqualUnsorted()
32 @*/
ISEqual(IS is1,IS is2,PetscBool * flg)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   PetscErrorCode ierr;
40   PetscMPIInt    mflg;
41 
42   PetscFunctionBegin;
43   PetscValidHeaderSpecific(is1,IS_CLASSID,1);
44   PetscValidHeaderSpecific(is2,IS_CLASSID,2);
45   PetscValidBoolPointer(flg,3);
46 
47   if (is1 == is2) {
48     *flg = PETSC_TRUE;
49     PetscFunctionReturn(0);
50   }
51 
52   ierr = MPI_Comm_compare(PetscObjectComm((PetscObject)is1),PetscObjectComm((PetscObject)is2),&mflg);CHKERRQ(ierr);
53   if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
54     *flg = PETSC_FALSE;
55     PetscFunctionReturn(0);
56   }
57 
58   ierr = ISGetSize(is1,&sz1);CHKERRQ(ierr);
59   ierr = ISGetSize(is2,&sz2);CHKERRQ(ierr);
60   if (sz1 != sz2) *flg = PETSC_FALSE;
61   else {
62     ierr = ISGetLocalSize(is1,&sz1);CHKERRQ(ierr);
63     ierr = ISGetLocalSize(is2,&sz2);CHKERRQ(ierr);
64 
65     if (sz1 != sz2) flag = PETSC_FALSE;
66     else {
67       ierr = ISGetIndices(is1,&ptr1);CHKERRQ(ierr);
68       ierr = ISGetIndices(is2,&ptr2);CHKERRQ(ierr);
69 
70       ierr = PetscMalloc1(sz1,&a1);CHKERRQ(ierr);
71       ierr = PetscMalloc1(sz2,&a2);CHKERRQ(ierr);
72 
73       ierr = PetscArraycpy(a1,ptr1,sz1);CHKERRQ(ierr);
74       ierr = PetscArraycpy(a2,ptr2,sz2);CHKERRQ(ierr);
75 
76       ierr = PetscIntSortSemiOrdered(sz1,a1);CHKERRQ(ierr);
77       ierr = PetscIntSortSemiOrdered(sz2,a2);CHKERRQ(ierr);
78       ierr = PetscArraycmp(a1,a2,sz1,&flag);CHKERRQ(ierr);
79 
80       ierr = ISRestoreIndices(is1,&ptr1);CHKERRQ(ierr);
81       ierr = ISRestoreIndices(is2,&ptr2);CHKERRQ(ierr);
82 
83       ierr = PetscFree(a1);CHKERRQ(ierr);
84       ierr = PetscFree(a2);CHKERRQ(ierr);
85     }
86     ierr = PetscObjectGetComm((PetscObject)is1,&comm);CHKERRQ(ierr);
87     ierr = MPIU_Allreduce(&flag,flg,1,MPIU_BOOL,MPI_MIN,comm);CHKERRQ(ierr);
88   }
89   PetscFunctionReturn(0);
90 }
91 
92 /*@
93    ISEqualUnsorted  - Compares if two index sets have the same set of indices.
94 
95    Collective on IS
96 
97    Input Parameters:
98 .  is1, is2 - The index sets being compared
99 
100    Output Parameters:
101 .  flg - output flag, either PETSC_TRUE (if both index sets have the
102          same indices), or PETSC_FALSE if the index sets differ by size
103          or by the set of indices)
104 
105    Level: intermediate
106 
107    Note:
108    This routine does NOT sort the contents of the index sets before
109    the comparision is made.
110 
111 
112 .seealso: ISEqual()
113 @*/
ISEqualUnsorted(IS is1,IS is2,PetscBool * flg)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   PetscErrorCode ierr;
121   PetscMPIInt    mflg;
122 
123   PetscFunctionBegin;
124   PetscValidHeaderSpecific(is1,IS_CLASSID,1);
125   PetscValidHeaderSpecific(is2,IS_CLASSID,2);
126   PetscValidBoolPointer(flg,3);
127 
128   if (is1 == is2) {
129     *flg = PETSC_TRUE;
130     PetscFunctionReturn(0);
131   }
132 
133   ierr = MPI_Comm_compare(PetscObjectComm((PetscObject)is1),PetscObjectComm((PetscObject)is2),&mflg);CHKERRQ(ierr);
134   if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
135     *flg = PETSC_FALSE;
136     PetscFunctionReturn(0);
137   }
138 
139   ierr = ISGetSize(is1,&sz1);CHKERRQ(ierr);
140   ierr = ISGetSize(is2,&sz2);CHKERRQ(ierr);
141   if (sz1 != sz2) *flg = PETSC_FALSE;
142   else {
143     ierr = ISGetLocalSize(is1,&sz1);CHKERRQ(ierr);
144     ierr = ISGetLocalSize(is2,&sz2);CHKERRQ(ierr);
145 
146     if (sz1 != sz2) flag = PETSC_FALSE;
147     else {
148       ierr = ISGetIndices(is1,&ptr1);CHKERRQ(ierr);
149       ierr = ISGetIndices(is2,&ptr2);CHKERRQ(ierr);
150 
151       ierr = PetscArraycmp(ptr1,ptr2,sz1,&flag);CHKERRQ(ierr);
152 
153       ierr = ISRestoreIndices(is1,&ptr1);CHKERRQ(ierr);
154       ierr = ISRestoreIndices(is2,&ptr2);CHKERRQ(ierr);
155     }
156     ierr = PetscObjectGetComm((PetscObject)is1,&comm);CHKERRQ(ierr);
157     ierr = MPIU_Allreduce(&flag,flg,1,MPIU_BOOL,MPI_MIN,comm);CHKERRQ(ierr);
158   }
159   PetscFunctionReturn(0);
160 }
161 
162 
163