Actual source code: sftype.c

  1: #include <petsc/private/sfimpl.h>

  3: #if !defined(PETSC_HAVE_MPI_COMBINER_DUP) && !defined(MPI_COMBINER_DUP)  /* We have no way to interpret output of MPI_Type_get_envelope without this. */
  4: #  define MPI_COMBINER_DUP   0
  5: #endif
  6: #if !defined(PETSC_HAVE_MPI_COMBINER_NAMED) && !defined(MPI_COMBINER_NAMED)
  7: #define MPI_COMBINER_NAMED -2
  8: #endif
  9: #if !defined(PETSC_HAVE_MPI_COMBINER_CONTIGUOUS) && !defined(MPI_COMBINER_CONTIGUOUS) && MPI_VERSION < 2
 10: #  define MPI_COMBINER_CONTIGUOUS -1
 11: #endif

 13: static PetscErrorCode MPIPetsc_Type_free(MPI_Datatype *a)
 14: {
 15:   PetscMPIInt    nints,naddrs,ntypes,combiner;

 17:   MPI_Type_get_envelope(*a,&nints,&naddrs,&ntypes,&combiner);

 19:   if (combiner != MPI_COMBINER_NAMED) {
 20:     MPI_Type_free(a);
 21:   }

 23:   *a = MPI_DATATYPE_NULL;
 24:   return 0;
 25: }

 27: /*
 28:   Unwrap an MPI datatype recursively in case it is dupped or MPI_Type_contiguous(1,...)'ed from another type.

 30:    Input Parameter:
 31: .  a  - the datatype to be unwrapped

 33:    Output Parameters:
 34: + atype - the unwrapped datatype, which is either equal(=) to a or equivalent to a.
 35: - flg   - true if atype != a, which implies caller should MPIPetsc_Type_free(atype) after use. Note atype might be MPI builtin.
 36: */
 37: PetscErrorCode MPIPetsc_Type_unwrap(MPI_Datatype a,MPI_Datatype *atype,PetscBool *flg)
 38: {
 39:   PetscMPIInt    nints,naddrs,ntypes,combiner,ints[1];
 40:   MPI_Aint       addrs[1];
 41:   MPI_Datatype   types[1];

 43:   *flg = PETSC_FALSE;
 44:   *atype = a;
 45:   if (a == MPIU_INT || a == MPIU_REAL || a == MPIU_SCALAR) return 0;
 46:   MPI_Type_get_envelope(a,&nints,&naddrs,&ntypes,&combiner);
 47:   if (combiner == MPI_COMBINER_DUP) {
 49:     MPI_Type_get_contents(a,0,0,1,ints,addrs,types);
 50:     /* Recursively unwrap dupped types. */
 51:     MPIPetsc_Type_unwrap(types[0],atype,flg);
 52:     if (*flg) {
 53:       /* If the recursive call returns a new type, then that means that atype[0] != types[0] and we're on the hook to
 54:        * free types[0].  Note that this case occurs if combiner(types[0]) is MPI_COMBINER_DUP, so we're safe to
 55:        * directly call MPI_Type_free rather than MPIPetsc_Type_free here. */
 56:       MPI_Type_free(&(types[0]));
 57:     }
 58:     /* In any case, it's up to the caller to free the returned type in this case. */
 59:     *flg = PETSC_TRUE;
 60:   } else if (combiner == MPI_COMBINER_CONTIGUOUS) {
 62:     MPI_Type_get_contents(a,1,0,1,ints,addrs,types);
 63:     if (ints[0] == 1) { /* If a is created by MPI_Type_contiguous(1,..) */
 64:       MPIPetsc_Type_unwrap(types[0],atype,flg);
 65:       if (*flg) MPIPetsc_Type_free(&(types[0]));
 66:       *flg = PETSC_TRUE;
 67:     } else {
 68:       MPIPetsc_Type_free(&(types[0]));
 69:     }
 70:   }
 71:   return 0;
 72: }

 74: PetscErrorCode MPIPetsc_Type_compare(MPI_Datatype a,MPI_Datatype b,PetscBool *match)
 75: {
 76:   MPI_Datatype   atype,btype;
 77:   PetscMPIInt    aintcount,aaddrcount,atypecount,acombiner;
 78:   PetscMPIInt    bintcount,baddrcount,btypecount,bcombiner;
 79:   PetscBool      freeatype, freebtype;

 81:   if (a == b) { /* this is common when using MPI builtin datatypes */
 82:     *match = PETSC_TRUE;
 83:     return 0;
 84:   }
 85:   MPIPetsc_Type_unwrap(a,&atype,&freeatype);
 86:   MPIPetsc_Type_unwrap(b,&btype,&freebtype);
 87:   *match = PETSC_FALSE;
 88:   if (atype == btype) {
 89:     *match = PETSC_TRUE;
 90:     goto free_types;
 91:   }
 92:   MPI_Type_get_envelope(atype,&aintcount,&aaddrcount,&atypecount,&acombiner);
 93:   MPI_Type_get_envelope(btype,&bintcount,&baddrcount,&btypecount,&bcombiner);
 94:   if (acombiner == bcombiner && aintcount == bintcount && aaddrcount == baddrcount && atypecount == btypecount && (aintcount > 0 || aaddrcount > 0 || atypecount > 0)) {
 95:     PetscMPIInt  *aints,*bints;
 96:     MPI_Aint     *aaddrs,*baddrs;
 97:     MPI_Datatype *atypes,*btypes;
 98:     PetscInt     i;
 99:     PetscBool    same;
100:     PetscMalloc6(aintcount,&aints,bintcount,&bints,aaddrcount,&aaddrs,baddrcount,&baddrs,atypecount,&atypes,btypecount,&btypes);
101:     MPI_Type_get_contents(atype,aintcount,aaddrcount,atypecount,aints,aaddrs,atypes);
102:     MPI_Type_get_contents(btype,bintcount,baddrcount,btypecount,bints,baddrs,btypes);
103:     PetscArraycmp(aints,bints,aintcount,&same);
104:     if (same) {
105:       PetscArraycmp(aaddrs,baddrs,aaddrcount,&same);
106:       if (same) {
107:         /* Check for identity first */
108:         PetscArraycmp(atypes,btypes,atypecount,&same);
109:         if (!same) {
110:           /* If the atype or btype were not predefined data types, then the types returned from MPI_Type_get_contents
111:            * will merely be equivalent to the types used in the construction, so we must recursively compare. */
112:           for (i=0; i<atypecount; i++) {
113:             MPIPetsc_Type_compare(atypes[i],btypes[i],&same);
114:             if (!same) break;
115:           }
116:         }
117:       }
118:     }
119:     for (i=0; i<atypecount; i++) {
120:       MPIPetsc_Type_free(&(atypes[i]));
121:       MPIPetsc_Type_free(&(btypes[i]));
122:     }
123:     PetscFree6(aints,bints,aaddrs,baddrs,atypes,btypes);
124:     if (same) *match = PETSC_TRUE;
125:   }
126: free_types:
127:   if (freeatype) {
128:     MPIPetsc_Type_free(&atype);
129:   }
130:   if (freebtype) {
131:     MPIPetsc_Type_free(&btype);
132:   }
133:   return 0;
134: }

136: /* Check whether a was created via MPI_Type_contiguous from b
137:  *
138:  */
139: PetscErrorCode MPIPetsc_Type_compare_contig(MPI_Datatype a,MPI_Datatype b,PetscInt *n)
140: {
141:   MPI_Datatype   atype,btype;
142:   PetscMPIInt    aintcount,aaddrcount,atypecount,acombiner;
143:   PetscBool      freeatype,freebtype;

145:   if (a == b) {*n = 1; return 0;}
146:   *n = 0;
147:   MPIPetsc_Type_unwrap(a,&atype,&freeatype);
148:   MPIPetsc_Type_unwrap(b,&btype,&freebtype);
149:   MPI_Type_get_envelope(atype,&aintcount,&aaddrcount,&atypecount,&acombiner);
150:   if (acombiner == MPI_COMBINER_CONTIGUOUS && aintcount >= 1) {
151:     PetscMPIInt  *aints;
152:     MPI_Aint     *aaddrs;
153:     MPI_Datatype *atypes;
154:     PetscInt      i;
155:     PetscBool     same;
156:     PetscMalloc3(aintcount,&aints,aaddrcount,&aaddrs,atypecount,&atypes);
157:     MPI_Type_get_contents(atype,aintcount,aaddrcount,atypecount,aints,aaddrs,atypes);
158:     /* Check for identity first. */
159:     if (atypes[0] == btype) {
160:       *n = aints[0];
161:     } else {
162:       /* atypes[0] merely has to be equivalent to the type used to create atype. */
163:       MPIPetsc_Type_compare(atypes[0],btype,&same);
164:       if (same) *n = aints[0];
165:     }
166:     for (i=0; i<atypecount; i++) {
167:       MPIPetsc_Type_free(&(atypes[i]));
168:     }
169:     PetscFree3(aints,aaddrs,atypes);
170:   }

172:   if (freeatype) {
173:     MPIPetsc_Type_free(&atype);
174:   }
175:   if (freebtype) {
176:     MPIPetsc_Type_free(&btype);
177:   }
178:   return 0;
179: }