Actual source code: iscomp.c

  1: /*$Id: iscomp.c,v 1.35 2001/03/23 23:21:16 balay Exp $*/

 3:  #include petscsys.h
 4:  #include petscis.h

  6: /*@C
  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: int ISEqual(IS is1,IS is2,PetscTruth *flg)
 37: {
 38:   int        sz1,sz2,ierr,*ptr1,*ptr2,*a1,*a2;
 39:   PetscTruth flag;
 40:   MPI_Comm   comm;

 47: 
 48:   ISGetSize(is1,&sz1);
 49:   ISGetSize(is2,&sz2);
 50:   if (sz1 != sz2) {
 51:     *flg = PETSC_FALSE;
 52:   } else {
 53:     ISGetLocalSize(is1,&sz1);
 54:     ISGetLocalSize(is2,&sz2);

 56:     if (sz1 != sz2) {
 57:       flag = PETSC_FALSE;
 58:     } else {
 59:       ISGetIndices(is1,&ptr1);
 60:       ISGetIndices(is2,&ptr2);
 61: 
 62:       PetscMalloc((sz1+1)*sizeof(int),&a1);
 63:       PetscMalloc((sz2+1)*sizeof(int),&a2);

 65:       PetscMemcpy(a1,ptr1,sz1*sizeof(int));
 66:       PetscMemcpy(a2,ptr2,sz2*sizeof(int));

 68:       PetscSortInt(sz1,a1);
 69:       PetscSortInt(sz2,a2);
 70:       PetscMemcmp(a1,a2,sz1*sizeof(int),&flag);

 72:       ISRestoreIndices(is1,&ptr1);
 73:       ISRestoreIndices(is2,&ptr2);
 74: 
 75:       PetscFree(a1);
 76:       PetscFree(a2);
 77:     }
 78:     PetscObjectGetComm((PetscObject)is1,&comm);
 79:     MPI_Allreduce(&flag,flg,1,MPI_INT,MPI_MIN,comm);
 80:   }
 81:   return(0);
 82: }
 83: