Actual source code: sftype.c
petsc-3.11.4 2019-09-28
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: if (a == b) { /* this is common when using MPI builtin datatypes */
73: *match = PETSC_TRUE;
74: return(0);
75: }
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: }
127: /* Check whether a was created via MPI_Type_contiguous from b
128: *
129: */
130: PetscErrorCode MPIPetsc_Type_compare_contig(MPI_Datatype a,MPI_Datatype b,PetscInt *n)
131: {
133: MPI_Datatype atype,btype;
134: PetscMPIInt aintcount,aaddrcount,atypecount,acombiner;
135: PetscBool freeatype,freebtype;
137: MPIPetsc_Type_unwrap(a,&atype,&freeatype);
138: MPIPetsc_Type_unwrap(b,&btype,&freebtype);
139: *n = PETSC_FALSE;
140: if (atype == btype) {
141: *n = 1;
142: goto free_types;
143: }
144: MPI_Type_get_envelope(atype,&aintcount,&aaddrcount,&atypecount,&acombiner);
145: if (acombiner == MPI_COMBINER_CONTIGUOUS && aintcount >= 1) {
146: PetscMPIInt *aints;
147: MPI_Aint *aaddrs;
148: MPI_Datatype *atypes;
149: PetscInt i;
150: PetscBool same;
151: PetscMalloc3(aintcount,&aints,aaddrcount,&aaddrs,atypecount,&atypes);
152: MPI_Type_get_contents(atype,aintcount,aaddrcount,atypecount,aints,aaddrs,atypes);
153: /* Check for identity first. */
154: if (atypes[0] == btype) {
155: *n = aints[0];
156: } else {
157: /* atypes[0] merely has to be equivalent to the type used to create atype. */
158: MPIPetsc_Type_compare(atypes[0],btype,&same);
159: if (same) *n = aints[0];
160: }
161: for (i=0; i<atypecount; i++) {
162: MPIPetsc_Type_free(&(atypes[i]));
163: }
164: PetscFree3(aints,aaddrs,atypes);
165: }
166: free_types:
167: if (freeatype) {
168: MPIPetsc_Type_free(&atype);
169: }
170: if (freebtype) {
171: MPIPetsc_Type_free(&btype);
172: }
173: return(0);
174: }