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) MPI_Type_free(a);

 21:   *a = MPI_DATATYPE_NULL;
 22:   return 0;
 23: }

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

 28:    Input Parameter:
 29: .  a  - the datatype to be unwrapped

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

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

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

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

130: /* Check whether a was created via MPI_Type_contiguous from b
131:  *
132:  */
133: PetscErrorCode MPIPetsc_Type_compare_contig(MPI_Datatype a, MPI_Datatype b, PetscInt *n)
134: {
135:   MPI_Datatype atype, btype;
136:   PetscMPIInt  aintcount, aaddrcount, atypecount, acombiner;
137:   PetscBool    freeatype, freebtype;

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

167:   if (freeatype) MPIPetsc_Type_free(&atype);
168:   if (freebtype) MPIPetsc_Type_free(&btype);
169:   return 0;
170: }