Actual source code: iscomp.c
petsc-3.7.5 2017-01-01
2: #include <petsc/private/isimpl.h> /*I "petscis.h" I*/
6: /*@
7: ISEqual - Compares if two index sets have the same set of indices.
9: Collective on IS
11: Input Parameters:
12: . is1, is2 - The index sets being compared
14: Output Parameters:
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: This routine sorts the contents of the index sets before
23: the comparision 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: $ Processor
27: $ 0 1
28: $ is1 = {0, 1} {2, 3}
29: $ is2 = {2, 3} {0, 1}
30: will return false.
32: Concepts: index sets^equal
33: Concepts: IS^equal
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;
43: PetscMPIInt mflg;
50: if (is1 == is2) {
51: *flg = PETSC_TRUE;
52: return(0);
53: }
55: MPI_Comm_compare(PetscObjectComm((PetscObject)is1),PetscObjectComm((PetscObject)is2),&mflg);
56: if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
57: *flg = PETSC_FALSE;
58: return(0);
59: }
61: ISGetSize(is1,&sz1);
62: ISGetSize(is2,&sz2);
63: if (sz1 != sz2) *flg = PETSC_FALSE;
64: else {
65: ISGetLocalSize(is1,&sz1);
66: ISGetLocalSize(is2,&sz2);
68: if (sz1 != sz2) flag = PETSC_FALSE;
69: else {
70: ISGetIndices(is1,&ptr1);
71: ISGetIndices(is2,&ptr2);
73: PetscMalloc1(sz1,&a1);
74: PetscMalloc1(sz2,&a2);
76: PetscMemcpy(a1,ptr1,sz1*sizeof(PetscInt));
77: PetscMemcpy(a2,ptr2,sz2*sizeof(PetscInt));
79: PetscSortInt(sz1,a1);
80: PetscSortInt(sz2,a2);
81: PetscMemcmp(a1,a2,sz1*sizeof(PetscInt),&flag);
83: ISRestoreIndices(is1,&ptr1);
84: ISRestoreIndices(is2,&ptr2);
86: PetscFree(a1);
87: PetscFree(a2);
88: }
89: PetscObjectGetComm((PetscObject)is1,&comm);
90: MPIU_Allreduce(&flag,flg,1,MPIU_BOOL,MPI_MIN,comm);
91: }
92: return(0);
93: }