Actual source code: iscomp.c


  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 comparison 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: .seealso: ISEqualUnsorted()
 31: @*/
 32: PetscErrorCode  ISEqual(IS is1,IS is2,PetscBool  *flg)
 33: {
 34:   PetscInt       sz1,sz2,*a1,*a2;
 35:   const PetscInt *ptr1,*ptr2;
 36:   PetscBool      flag;
 37:   MPI_Comm       comm;
 39:   PetscMPIInt    mflg;


 46:   if (is1 == is2) {
 47:     *flg = PETSC_TRUE;
 48:     return(0);
 49:   }

 51:   MPI_Comm_compare(PetscObjectComm((PetscObject)is1),PetscObjectComm((PetscObject)is2),&mflg);
 52:   if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
 53:     *flg = PETSC_FALSE;
 54:     return(0);
 55:   }

 57:   ISGetSize(is1,&sz1);
 58:   ISGetSize(is2,&sz2);
 59:   if (sz1 != sz2) *flg = PETSC_FALSE;
 60:   else {
 61:     ISGetLocalSize(is1,&sz1);
 62:     ISGetLocalSize(is2,&sz2);

 64:     if (sz1 != sz2) flag = PETSC_FALSE;
 65:     else {
 66:       ISGetIndices(is1,&ptr1);
 67:       ISGetIndices(is2,&ptr2);

 69:       PetscMalloc1(sz1,&a1);
 70:       PetscMalloc1(sz2,&a2);

 72:       PetscArraycpy(a1,ptr1,sz1);
 73:       PetscArraycpy(a2,ptr2,sz2);

 75:       PetscIntSortSemiOrdered(sz1,a1);
 76:       PetscIntSortSemiOrdered(sz2,a2);
 77:       PetscArraycmp(a1,a2,sz1,&flag);

 79:       ISRestoreIndices(is1,&ptr1);
 80:       ISRestoreIndices(is2,&ptr2);

 82:       PetscFree(a1);
 83:       PetscFree(a2);
 84:     }
 85:     PetscObjectGetComm((PetscObject)is1,&comm);
 86:     MPIU_Allreduce(&flag,flg,1,MPIU_BOOL,MPI_MIN,comm);
 87:   }
 88:   return(0);
 89: }

 91: /*@
 92:    ISEqualUnsorted  - Compares if two index sets have the same set of indices.

 94:    Collective on IS

 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:    This routine does NOT sort the contents of the index sets before
108:    the comparison is made.

110: .seealso: ISEqual()
111: @*/
112: PetscErrorCode  ISEqualUnsorted(IS is1,IS is2,PetscBool  *flg)
113: {
114:   PetscInt       sz1,sz2;
115:   const PetscInt *ptr1,*ptr2;
116:   PetscBool      flag;
117:   MPI_Comm       comm;
119:   PetscMPIInt    mflg;


126:   if (is1 == is2) {
127:     *flg = PETSC_TRUE;
128:     return(0);
129:   }

131:   MPI_Comm_compare(PetscObjectComm((PetscObject)is1),PetscObjectComm((PetscObject)is2),&mflg);
132:   if (mflg != MPI_CONGRUENT && mflg != MPI_IDENT) {
133:     *flg = PETSC_FALSE;
134:     return(0);
135:   }

137:   ISGetSize(is1,&sz1);
138:   ISGetSize(is2,&sz2);
139:   if (sz1 != sz2) *flg = PETSC_FALSE;
140:   else {
141:     ISGetLocalSize(is1,&sz1);
142:     ISGetLocalSize(is2,&sz2);

144:     if (sz1 != sz2) flag = PETSC_FALSE;
145:     else {
146:       ISGetIndices(is1,&ptr1);
147:       ISGetIndices(is2,&ptr2);

149:       PetscArraycmp(ptr1,ptr2,sz1,&flag);

151:       ISRestoreIndices(is1,&ptr1);
152:       ISRestoreIndices(is2,&ptr2);
153:     }
154:     PetscObjectGetComm((PetscObject)is1,&comm);
155:     MPIU_Allreduce(&flag,flg,1,MPIU_BOOL,MPI_MIN,comm);
156:   }
157:   return(0);
158: }