Actual source code: sftype.c
petsc-3.9.3 2018-07-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) && !defined(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) && !defined(MPI_COMBINER_NAMED)
11: #define MPI_COMBINER_NAMED -2
12: #endif
13: #if !defined(PETSC_HAVE_MPI_COMBINER_CONTIGUOUS) && !defined(MPI_COMBINER_CONTIGUOUS) && MPI_VERSION < 2
14: # define MPI_COMBINER_CONTIGUOUS -1
15: #endif
17: static PetscErrorCode MPIPetsc_Type_free(MPI_Datatype *a)
18: {
19: PetscMPIInt nints,naddrs,ntypes,combiner;
23: MPI_Type_get_envelope(*a,&nints,&naddrs,&ntypes,&combiner);
25: if (combiner != MPI_COMBINER_NAMED) {
26: MPI_Type_free(a);
27: }
29: *a = MPI_DATATYPE_NULL;
30: return(0);
31: }
33: PetscErrorCode MPIPetsc_Type_unwrap(MPI_Datatype a,MPI_Datatype *atype,PetscBool *flg)
34: {
35: PetscMPIInt nints,naddrs,ntypes,combiner;
39: *flg = PETSC_FALSE;
40: *atype = a;
41: if (a == MPIU_INT || a == MPIU_REAL || a == MPIU_SCALAR) return(0);
42: MPI_Type_get_envelope(a,&nints,&naddrs,&ntypes,&combiner);
43: if (combiner == MPI_COMBINER_DUP) {
44: PetscMPIInt ints[1];
45: MPI_Aint addrs[1];
46: MPI_Datatype types[1];
47: if (nints != 0 || naddrs != 0 || ntypes != 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Unexpected returns from MPI_Type_get_envelope()");
48: MPI_Type_get_contents(a,0,0,1,ints,addrs,types);
49: /* Recursively unwrap dupped types. */
50: MPIPetsc_Type_unwrap(types[0],atype,flg);
51: if (*flg) {
52: /* If the recursive call returns a new type, then that means that atype[0] != types[0] and we're on the hook to
53: * free types[0]. Note that this case occurs if combiner(types[0]) is MPI_COMBINER_DUP, so we're safe to
54: * directly call MPI_Type_free rather than MPIPetsc_Type_free here. */
55: MPI_Type_free(&(types[0]));
56: }
57: /* In any case, it's up to the caller to free the returned type in this case. */
58: *flg = PETSC_TRUE;
59: }
60: return(0);
61: }
63: PetscErrorCode MPIPetsc_Type_compare(MPI_Datatype a,MPI_Datatype b,PetscBool *match)
64: {
66: MPI_Datatype atype,btype;
67: PetscMPIInt aintcount,aaddrcount,atypecount,acombiner;
68: PetscMPIInt bintcount,baddrcount,btypecount,bcombiner;
69: PetscBool freeatype, freebtype;
72: MPIPetsc_Type_unwrap(a,&atype,&freeatype);
73: MPIPetsc_Type_unwrap(b,&btype,&freebtype);
74: *match = PETSC_FALSE;
75: if (atype == btype) {
76: *match = PETSC_TRUE;
77: goto free_types;
78: }
79: MPI_Type_get_envelope(atype,&aintcount,&aaddrcount,&atypecount,&acombiner);
80: MPI_Type_get_envelope(btype,&bintcount,&baddrcount,&btypecount,&bcombiner);
81: if (acombiner == bcombiner && aintcount == bintcount && aaddrcount == baddrcount && atypecount == btypecount && (aintcount > 0 || aaddrcount > 0 || atypecount > 0)) {
82: PetscMPIInt *aints,*bints;
83: MPI_Aint *aaddrs,*baddrs;
84: MPI_Datatype *atypes,*btypes;
85: PetscInt i;
86: PetscBool same;
87: PetscMalloc6(aintcount,&aints,bintcount,&bints,aaddrcount,&aaddrs,baddrcount,&baddrs,atypecount,&atypes,btypecount,&btypes);
88: MPI_Type_get_contents(atype,aintcount,aaddrcount,atypecount,aints,aaddrs,atypes);
89: MPI_Type_get_contents(btype,bintcount,baddrcount,btypecount,bints,baddrs,btypes);
90: PetscMemcmp(aints,bints,aintcount*sizeof(aints[0]),&same);
91: if (same) {
92: PetscMemcmp(aaddrs,baddrs,aaddrcount*sizeof(aaddrs[0]),&same);
93: if (same) {
94: /* Check for identity first */
95: PetscMemcmp(atypes,btypes,atypecount*sizeof(atypes[0]),&same);
96: if (!same) {
97: /* If the atype or btype were not predefined data types, then the types returned from MPI_Type_get_contents
98: * will merely be equivalent to the types used in the construction, so we must recursively compare. */
99: for (i=0; i<atypecount; i++) {
100: MPIPetsc_Type_compare(atypes[i],btypes[i],&same);
101: if (!same) break;
102: }
103: }
104: }
105: }
106: for (i=0; i<atypecount; i++) {
107: MPIPetsc_Type_free(&(atypes[i]));
108: MPIPetsc_Type_free(&(btypes[i]));
109: }
110: PetscFree6(aints,bints,aaddrs,baddrs,atypes,btypes);
111: if (same) *match = PETSC_TRUE;
112: }
113: free_types:
114: if (freeatype) {
115: MPIPetsc_Type_free(&atype);
116: }
117: if (freebtype) {
118: MPIPetsc_Type_free(&btype);
119: }
120: return(0);
121: }
123: /* Check whether a was created via MPI_Type_contiguous from b
124: *
125: */
126: PetscErrorCode MPIPetsc_Type_compare_contig(MPI_Datatype a,MPI_Datatype b,PetscInt *n)
127: {
129: MPI_Datatype atype,btype;
130: PetscMPIInt aintcount,aaddrcount,atypecount,acombiner;
131: PetscBool freeatype,freebtype;
133: MPIPetsc_Type_unwrap(a,&atype,&freeatype);
134: MPIPetsc_Type_unwrap(b,&btype,&freebtype);
135: *n = PETSC_FALSE;
136: if (atype == btype) {
137: *n = 1;
138: goto free_types;
139: }
140: MPI_Type_get_envelope(atype,&aintcount,&aaddrcount,&atypecount,&acombiner);
141: if (acombiner == MPI_COMBINER_CONTIGUOUS && aintcount >= 1) {
142: PetscMPIInt *aints;
143: MPI_Aint *aaddrs;
144: MPI_Datatype *atypes;
145: PetscInt i;
146: PetscBool same;
147: PetscMalloc3(aintcount,&aints,aaddrcount,&aaddrs,atypecount,&atypes);
148: MPI_Type_get_contents(atype,aintcount,aaddrcount,atypecount,aints,aaddrs,atypes);
149: /* Check for identity first. */
150: if (atypes[0] == btype) {
151: *n = aints[0];
152: } else {
153: /* atypes[0] merely has to be equivalent to the type used to create atype. */
154: MPIPetsc_Type_compare(atypes[0],btype,&same);
155: if (same) *n = aints[0];
156: }
157: for (i=0; i<atypecount; i++) {
158: MPIPetsc_Type_free(&(atypes[i]));
159: }
160: PetscFree3(aints,aaddrs,atypes);
161: }
162: free_types:
163: if (freeatype) {
164: MPIPetsc_Type_free(&atype);
165: }
166: if (freebtype) {
167: MPIPetsc_Type_free(&btype);
168: }
169: return(0);
170: }