Actual source code: sftype.c
petsc-3.7.4 2016-10-02
1: #include <petsc/private/sfimpl.h>
3: #if !defined(PETSC_HAVE_MPI_TYPE_GET_ENVELOPE)
4: #define MPI_Type_get_envelope(datatype,num_ints,num_addrs,num_dtypes,combiner) (*(num_ints)=0,*(num_addrs)=0,*(num_dtypes)=0,*(combiner)=0,1);SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP_SYS,"Need an MPI-2 implementation")
5: #define MPI_Type_get_contents(datatype,num_ints,num_addrs,num_dtypes,ints,addrs,dtypes) (*(ints)=0,*(addrs)=0,*(dtypes)=0,1);SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP_SYS,"Need an MPI-2 implementation")
6: #endif
7: #if !defined(PETSC_HAVE_MPI_COMBINER_DUP) /* We have no way to interpret output of MPI_Type_get_envelope without this. */
8: # define MPI_COMBINER_DUP 0
9: #endif
10: #if !defined(PETSC_HAVE_MPI_COMBINER_NAMED)
11: #define MPI_COMBINER_NAMED -2
12: #endif
13: #if !defined(PETSC_HAVE_MPI_COMBINER_CONTIGUOUS) && MPI_VERSION < 2
14: # define MPI_COMBINER_CONTIGUOUS -1
15: #endif
19: static PetscErrorCode MPIPetsc_Type_free(MPI_Datatype *a)
20: {
21: PetscMPIInt nints,naddrs,ntypes,combiner;
25: MPI_Type_get_envelope(*a,&nints,&naddrs,&ntypes,&combiner);
27: if (combiner != MPI_COMBINER_NAMED) {
28: MPI_Type_free(a);
29: }
31: *a = MPI_DATATYPE_NULL;
32: return(0);
33: }
37: PetscErrorCode MPIPetsc_Type_unwrap(MPI_Datatype a,MPI_Datatype *atype,PetscBool *flg)
38: {
39: PetscMPIInt nints,naddrs,ntypes,combiner;
43: *flg = PETSC_FALSE;
44: MPI_Type_get_envelope(a,&nints,&naddrs,&ntypes,&combiner);
45: if (combiner == MPI_COMBINER_DUP) {
46: PetscMPIInt ints[1];
47: MPI_Aint addrs[1];
48: MPI_Datatype types[1];
49: if (nints != 0 || naddrs != 0 || ntypes != 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Unexpected returns from MPI_Type_get_envelope()");
50: MPI_Type_get_contents(a,0,0,1,ints,addrs,types);
51: /* Recursively unwrap dupped types. */
52: MPIPetsc_Type_unwrap(types[0],atype,flg);
53: if (*flg) {
54: /* If the recursive call returns a new type, then that means that atype[0] != types[0] and we're on the hook to
55: * free types[0]. Note that this case occurs if combiner(types[0]) is MPI_COMBINER_DUP, so we're safe to
56: * directly call MPI_Type_free rather than MPIPetsc_Type_free here. */
57: MPI_Type_free(&(types[0]));
58: }
59: /* In any case, it's up to the caller to free the returned type in this case. */
60: *flg = PETSC_TRUE;
61: } else *atype = a;
62: return(0);
63: }
67: PetscErrorCode MPIPetsc_Type_compare(MPI_Datatype a,MPI_Datatype b,PetscBool *match)
68: {
70: MPI_Datatype atype,btype;
71: PetscMPIInt aintcount,aaddrcount,atypecount,acombiner;
72: PetscMPIInt bintcount,baddrcount,btypecount,bcombiner;
73: PetscBool freeatype, freebtype;
76: MPIPetsc_Type_unwrap(a,&atype,&freeatype);
77: MPIPetsc_Type_unwrap(b,&btype,&freebtype);
78: *match = PETSC_FALSE;
79: if (atype == btype) {
80: *match = PETSC_TRUE;
81: goto free_types;
82: }
83: MPI_Type_get_envelope(atype,&aintcount,&aaddrcount,&atypecount,&acombiner);
84: MPI_Type_get_envelope(btype,&bintcount,&baddrcount,&btypecount,&bcombiner);
85: if (acombiner == bcombiner && aintcount == bintcount && aaddrcount == baddrcount && atypecount == btypecount && (aintcount > 0 || aaddrcount > 0 || atypecount > 0)) {
86: PetscMPIInt *aints,*bints;
87: MPI_Aint *aaddrs,*baddrs;
88: MPI_Datatype *atypes,*btypes;
89: PetscInt i;
90: PetscBool same;
91: PetscMalloc6(aintcount,&aints,bintcount,&bints,aaddrcount,&aaddrs,baddrcount,&baddrs,atypecount,&atypes,btypecount,&btypes);
92: MPI_Type_get_contents(atype,aintcount,aaddrcount,atypecount,aints,aaddrs,atypes);
93: MPI_Type_get_contents(btype,bintcount,baddrcount,btypecount,bints,baddrs,btypes);
94: PetscMemcmp(aints,bints,aintcount*sizeof(aints[0]),&same);
95: if (same) {
96: PetscMemcmp(aaddrs,baddrs,aaddrcount*sizeof(aaddrs[0]),&same);
97: if (same) {
98: /* Check for identity first */
99: PetscMemcmp(atypes,btypes,atypecount*sizeof(atypes[0]),&same);
100: if (!same) {
101: /* If the atype or btype were not predefined data types, then the types returned from MPI_Type_get_contents
102: * will merely be equivalent to the types used in the construction, so we must recursively compare. */
103: for (i=0; i<atypecount; i++) {
104: MPIPetsc_Type_compare(atypes[i],btypes[i],&same);
105: if (!same) break;
106: }
107: }
108: }
109: }
110: for (i=0; i<atypecount; i++) {
111: MPIPetsc_Type_free(&(atypes[i]));
112: MPIPetsc_Type_free(&(btypes[i]));
113: }
114: PetscFree6(aints,bints,aaddrs,baddrs,atypes,btypes);
115: if (same) *match = PETSC_TRUE;
116: }
117: free_types:
118: if (freeatype) {
119: MPIPetsc_Type_free(&atype);
120: }
121: if (freebtype) {
122: MPIPetsc_Type_free(&btype);
123: }
124: return(0);
125: }
129: /* Check whether a was created via MPI_Type_contiguous from b
130: *
131: */
132: PetscErrorCode MPIPetsc_Type_compare_contig(MPI_Datatype a,MPI_Datatype b,PetscInt *n)
133: {
135: MPI_Datatype atype,btype;
136: PetscMPIInt aintcount,aaddrcount,atypecount,acombiner;
137: PetscBool freeatype,freebtype;
139: MPIPetsc_Type_unwrap(a,&atype,&freeatype);
140: MPIPetsc_Type_unwrap(b,&btype,&freebtype);
141: *n = PETSC_FALSE;
142: if (atype == btype) {
143: *n = 1;
144: goto free_types;
145: }
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++) {
164: MPIPetsc_Type_free(&(atypes[i]));
165: }
166: PetscFree3(aints,aaddrs,atypes);
167: }
168: free_types:
169: if (freeatype) {
170: MPIPetsc_Type_free(&atype);
171: }
172: if (freebtype) {
173: MPIPetsc_Type_free(&btype);
174: }
175: return(0);
176: }