Actual source code: sftype.c
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: /*
34: Unwrap an MPI datatype recursively in case it is dupped or MPI_Type_contiguous(1,...)'ed from another type.
36: Input Parameter:
37: . a - the datatype to be unwrapped
39: Output Parameters:
40: + atype - the unwrapped datatype, which is either equal(=) to a or equivalent to a.
41: - flg - true if atype != a, which implies caller should MPIPetsc_Type_free(atype) after use. Note atype might be MPI builtin.
42: */
43: PetscErrorCode MPIPetsc_Type_unwrap(MPI_Datatype a,MPI_Datatype *atype,PetscBool *flg)
44: {
46: PetscMPIInt nints,naddrs,ntypes,combiner,ints[1];
47: MPI_Aint addrs[1];
48: MPI_Datatype types[1];
51: *flg = PETSC_FALSE;
52: *atype = a;
53: if (a == MPIU_INT || a == MPIU_REAL || a == MPIU_SCALAR) return(0);
54: MPI_Type_get_envelope(a,&nints,&naddrs,&ntypes,&combiner);
55: if (combiner == MPI_COMBINER_DUP) {
56: if (nints != 0 || naddrs != 0 || ntypes != 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Unexpected returns from MPI_Type_get_envelope()");
57: MPI_Type_get_contents(a,0,0,1,ints,addrs,types);
58: /* Recursively unwrap dupped types. */
59: MPIPetsc_Type_unwrap(types[0],atype,flg);
60: if (*flg) {
61: /* If the recursive call returns a new type, then that means that atype[0] != types[0] and we're on the hook to
62: * free types[0]. Note that this case occurs if combiner(types[0]) is MPI_COMBINER_DUP, so we're safe to
63: * directly call MPI_Type_free rather than MPIPetsc_Type_free here. */
64: MPI_Type_free(&(types[0]));
65: }
66: /* In any case, it's up to the caller to free the returned type in this case. */
67: *flg = PETSC_TRUE;
68: } else if (combiner == MPI_COMBINER_CONTIGUOUS) {
69: if (nints != 1 || naddrs != 0 || ntypes != 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Unexpected returns from MPI_Type_get_envelope()");
70: MPI_Type_get_contents(a,1,0,1,ints,addrs,types);
71: if (ints[0] == 1) { /* If a is created by MPI_Type_contiguous(1,..) */
72: MPIPetsc_Type_unwrap(types[0],atype,flg);
73: if (*flg) {MPIPetsc_Type_free(&(types[0]));}
74: *flg = PETSC_TRUE;
75: } else {
76: MPIPetsc_Type_free(&(types[0]));
77: }
78: }
79: return(0);
80: }
82: PetscErrorCode MPIPetsc_Type_compare(MPI_Datatype a,MPI_Datatype b,PetscBool *match)
83: {
85: MPI_Datatype atype,btype;
86: PetscMPIInt aintcount,aaddrcount,atypecount,acombiner;
87: PetscMPIInt bintcount,baddrcount,btypecount,bcombiner;
88: PetscBool freeatype, freebtype;
91: if (a == b) { /* this is common when using MPI builtin datatypes */
92: *match = PETSC_TRUE;
93: return(0);
94: }
95: MPIPetsc_Type_unwrap(a,&atype,&freeatype);
96: MPIPetsc_Type_unwrap(b,&btype,&freebtype);
97: *match = PETSC_FALSE;
98: if (atype == btype) {
99: *match = PETSC_TRUE;
100: goto free_types;
101: }
102: MPI_Type_get_envelope(atype,&aintcount,&aaddrcount,&atypecount,&acombiner);
103: MPI_Type_get_envelope(btype,&bintcount,&baddrcount,&btypecount,&bcombiner);
104: if (acombiner == bcombiner && aintcount == bintcount && aaddrcount == baddrcount && atypecount == btypecount && (aintcount > 0 || aaddrcount > 0 || atypecount > 0)) {
105: PetscMPIInt *aints,*bints;
106: MPI_Aint *aaddrs,*baddrs;
107: MPI_Datatype *atypes,*btypes;
108: PetscInt i;
109: PetscBool same;
110: PetscMalloc6(aintcount,&aints,bintcount,&bints,aaddrcount,&aaddrs,baddrcount,&baddrs,atypecount,&atypes,btypecount,&btypes);
111: MPI_Type_get_contents(atype,aintcount,aaddrcount,atypecount,aints,aaddrs,atypes);
112: MPI_Type_get_contents(btype,bintcount,baddrcount,btypecount,bints,baddrs,btypes);
113: PetscArraycmp(aints,bints,aintcount,&same);
114: if (same) {
115: PetscArraycmp(aaddrs,baddrs,aaddrcount,&same);
116: if (same) {
117: /* Check for identity first */
118: PetscArraycmp(atypes,btypes,atypecount,&same);
119: if (!same) {
120: /* If the atype or btype were not predefined data types, then the types returned from MPI_Type_get_contents
121: * will merely be equivalent to the types used in the construction, so we must recursively compare. */
122: for (i=0; i<atypecount; i++) {
123: MPIPetsc_Type_compare(atypes[i],btypes[i],&same);
124: if (!same) break;
125: }
126: }
127: }
128: }
129: for (i=0; i<atypecount; i++) {
130: MPIPetsc_Type_free(&(atypes[i]));
131: MPIPetsc_Type_free(&(btypes[i]));
132: }
133: PetscFree6(aints,bints,aaddrs,baddrs,atypes,btypes);
134: if (same) *match = PETSC_TRUE;
135: }
136: free_types:
137: if (freeatype) {
138: MPIPetsc_Type_free(&atype);
139: }
140: if (freebtype) {
141: MPIPetsc_Type_free(&btype);
142: }
143: return(0);
144: }
146: /* Check whether a was created via MPI_Type_contiguous from b
147: *
148: */
149: PetscErrorCode MPIPetsc_Type_compare_contig(MPI_Datatype a,MPI_Datatype b,PetscInt *n)
150: {
152: MPI_Datatype atype,btype;
153: PetscMPIInt aintcount,aaddrcount,atypecount,acombiner;
154: PetscBool freeatype,freebtype;
157: if (a == b) {*n = 1; return(0);}
158: *n = 0;
159: MPIPetsc_Type_unwrap(a,&atype,&freeatype);
160: MPIPetsc_Type_unwrap(b,&btype,&freebtype);
161: MPI_Type_get_envelope(atype,&aintcount,&aaddrcount,&atypecount,&acombiner);
162: if (acombiner == MPI_COMBINER_CONTIGUOUS && aintcount >= 1) {
163: PetscMPIInt *aints;
164: MPI_Aint *aaddrs;
165: MPI_Datatype *atypes;
166: PetscInt i;
167: PetscBool same;
168: PetscMalloc3(aintcount,&aints,aaddrcount,&aaddrs,atypecount,&atypes);
169: MPI_Type_get_contents(atype,aintcount,aaddrcount,atypecount,aints,aaddrs,atypes);
170: /* Check for identity first. */
171: if (atypes[0] == btype) {
172: *n = aints[0];
173: } else {
174: /* atypes[0] merely has to be equivalent to the type used to create atype. */
175: MPIPetsc_Type_compare(atypes[0],btype,&same);
176: if (same) *n = aints[0];
177: }
178: for (i=0; i<atypecount; i++) {
179: MPIPetsc_Type_free(&(atypes[i]));
180: }
181: PetscFree3(aints,aaddrs,atypes);
182: }
184: if (freeatype) {
185: MPIPetsc_Type_free(&atype);
186: }
187: if (freebtype) {
188: MPIPetsc_Type_free(&btype);
189: }
190: return(0);
191: }