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;


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

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

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

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

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

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

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

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

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

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

 93:    Collective on is1

 95:    Input Parameters:
 96: .  is1, is2 - The index sets being compared

 98:    Output Parameters:
 99: .  flg - output flag, either `PETSC_TRUE` (if both index sets have the
100:          same indices), or `PETSC_FALSE` if the index sets differ by size
101:          or by the set of indices)

103:    Level: intermediate

105:    Note:
106:    Unlike ISEqual(), this routine does NOT sort the contents of the index sets before
107:    the comparison is made, i.e., the order of indices is important.

109:    Each MPI rank must have the same indices.

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


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

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

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

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

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

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