1 static char help[] = "This example shows 1) how to transfer vectors from a parent communicator to vectors on a child communicator and vice versa;\n\
2 2) how to transfer vectors from a subcommunicator to vectors on another subcommunicator. The two subcommunicators are not\n\
3 required to cover all processes in PETSC_COMM_WORLD; 3) how to copy a vector from a parent communicator to vectors on its child communicators.\n\
4 To run any example with VECCUDA vectors, add -vectype cuda to the argument list\n\n";
5
6 #include <petscvec.h>
main(int argc,char ** argv)7 int main(int argc, char **argv)
8 {
9 PetscMPIInt nproc, grank, mycolor;
10 PetscInt i, n, N = 20, low, high;
11 MPI_Comm subcomm;
12 Vec x = NULL; /* global vectors on PETSC_COMM_WORLD */
13 Vec yg = NULL; /* global vectors on PETSC_COMM_WORLD */
14 VecScatter vscat;
15 IS ix, iy;
16 PetscBool iscuda = PETSC_FALSE; /* Option to use VECCUDA vectors */
17 PetscBool optionflag, compareflag;
18 char vectypename[PETSC_MAX_PATH_LEN];
19 PetscBool world2sub = PETSC_FALSE; /* Copy a vector from WORLD to a subcomm? */
20 PetscBool sub2sub = PETSC_FALSE; /* Copy a vector from a subcomm to another subcomm? */
21 PetscBool world2subs = PETSC_FALSE; /* Copy a vector from WORLD to multiple subcomms? */
22
23 PetscFunctionBeginUser;
24 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
25 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &nproc));
26 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &grank));
27
28 PetscCheck(nproc >= 2, PETSC_COMM_WORLD, PETSC_ERR_ARG_SIZ, "This test must have at least two processes to run");
29
30 PetscCall(PetscOptionsGetBool(NULL, 0, "-world2sub", &world2sub, NULL));
31 PetscCall(PetscOptionsGetBool(NULL, 0, "-sub2sub", &sub2sub, NULL));
32 PetscCall(PetscOptionsGetBool(NULL, 0, "-world2subs", &world2subs, NULL));
33 PetscCall(PetscOptionsGetString(NULL, NULL, "-vectype", vectypename, sizeof(vectypename), &optionflag));
34 if (optionflag) {
35 PetscCall(PetscStrncmp(vectypename, "cuda", (size_t)4, &compareflag));
36 if (compareflag) iscuda = PETSC_TRUE;
37 }
38
39 /* Split PETSC_COMM_WORLD into three subcomms. Each process can only see the subcomm it belongs to */
40 mycolor = grank % 3;
41 PetscCallMPI(MPI_Comm_split(PETSC_COMM_WORLD, mycolor, grank, &subcomm));
42
43 /*===========================================================================
44 * Transfer a vector x defined on PETSC_COMM_WORLD to a vector y defined on
45 * a subcommunicator of PETSC_COMM_WORLD and vice versa.
46 *===========================================================================*/
47 if (world2sub) {
48 PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
49 PetscCall(VecSetSizes(x, PETSC_DECIDE, N));
50 if (iscuda) {
51 PetscCall(VecSetType(x, VECCUDA));
52 } else {
53 PetscCall(VecSetType(x, VECSTANDARD));
54 }
55 PetscCall(VecSetUp(x));
56 PetscCall(PetscObjectSetName((PetscObject)x, "x_commworld")); /* Give a name to view x clearly */
57
58 /* Initialize x to [-0.0, -1.0, -2.0, ..., -19.0] */
59 PetscCall(VecGetOwnershipRange(x, &low, &high));
60 for (i = low; i < high; i++) {
61 PetscScalar val = -i;
62 PetscCall(VecSetValue(x, i, val, INSERT_VALUES));
63 }
64 PetscCall(VecAssemblyBegin(x));
65 PetscCall(VecAssemblyEnd(x));
66
67 /* Transfer x to a vector y only defined on subcomm0 and vice versa */
68 if (mycolor == 0) { /* subcomm0 contains ranks 0, 3, 6, ... in PETSC_COMM_WORLD */
69 Vec y;
70 PetscScalar *yvalue;
71 PetscCall(VecCreate(subcomm, &y));
72 PetscCall(VecSetSizes(y, PETSC_DECIDE, N));
73 if (iscuda) {
74 PetscCall(VecSetType(y, VECCUDA));
75 } else {
76 PetscCall(VecSetType(y, VECSTANDARD));
77 }
78 PetscCall(VecSetUp(y));
79 PetscCall(PetscObjectSetName((PetscObject)y, "y_subcomm_0")); /* Give a name to view y clearly */
80 PetscCall(VecGetLocalSize(y, &n));
81 if (iscuda) {
82 #if defined(PETSC_HAVE_CUDA)
83 PetscCall(VecCUDAGetArray(y, &yvalue));
84 #endif
85 } else {
86 PetscCall(VecGetArray(y, &yvalue));
87 }
88 /* Create yg on PETSC_COMM_WORLD and alias yg with y. They share the memory pointed by yvalue.
89 Note this is a collective call. All processes have to call it and supply consistent N.
90 */
91 if (iscuda) {
92 #if defined(PETSC_HAVE_CUDA)
93 PetscCall(VecCreateMPICUDAWithArray(PETSC_COMM_WORLD, 1, n, N, yvalue, &yg));
94 #endif
95 } else {
96 PetscCall(VecCreateMPIWithArray(PETSC_COMM_WORLD, 1, n, N, yvalue, &yg));
97 }
98
99 /* Create an identity map that makes yg[i] = x[i], i=0..N-1 */
100 PetscCall(VecGetOwnershipRange(yg, &low, &high)); /* low, high are global indices */
101 PetscCall(ISCreateStride(PETSC_COMM_SELF, high - low, low, 1, &ix));
102 PetscCall(ISDuplicate(ix, &iy));
103
104 /* Union of ix's on subcomm0 covers the full range of [0,N) */
105 PetscCall(VecScatterCreate(x, ix, yg, iy, &vscat));
106 PetscCall(VecScatterBegin(vscat, x, yg, INSERT_VALUES, SCATTER_FORWARD));
107 PetscCall(VecScatterEnd(vscat, x, yg, INSERT_VALUES, SCATTER_FORWARD));
108
109 /* Once yg got the data from x, we return yvalue to y so that we can use y in other operations.
110 VecGetArray must be paired with VecRestoreArray.
111 */
112 if (iscuda) {
113 #if defined(PETSC_HAVE_CUDA)
114 PetscCall(VecCUDARestoreArray(y, &yvalue));
115 #endif
116 } else {
117 PetscCall(VecRestoreArray(y, &yvalue));
118 }
119
120 /* Libraries on subcomm0 can safely use y now, for example, view and scale it */
121 PetscCall(VecView(y, PETSC_VIEWER_STDOUT_(subcomm)));
122 PetscCall(VecScale(y, 2.0));
123
124 /* Send the new y back to x */
125 PetscCall(VecGetArray(y, &yvalue)); /* If VecScale is done on GPU, PETSc will prepare a valid yvalue for access */
126 /* Supply new yvalue to yg without memory copying */
127 PetscCall(VecPlaceArray(yg, yvalue));
128 PetscCall(VecScatterBegin(vscat, yg, x, INSERT_VALUES, SCATTER_REVERSE));
129 PetscCall(VecScatterEnd(vscat, yg, x, INSERT_VALUES, SCATTER_REVERSE));
130 PetscCall(VecResetArray(yg));
131 PetscCall(VecRestoreArray(y, &yvalue));
132 PetscCall(VecDestroy(&y));
133 } else {
134 /* Ranks outside of subcomm0 do not supply values to yg */
135 if (iscuda) {
136 #if defined(PETSC_HAVE_CUDA)
137 PetscCall(VecCreateMPICUDAWithArray(PETSC_COMM_WORLD, 1, 0 /*n*/, N, NULL, &yg));
138 #endif
139 } else {
140 PetscCall(VecCreateMPIWithArray(PETSC_COMM_WORLD, 1, 0 /*n*/, N, NULL, &yg));
141 }
142
143 /* Ranks in subcomm0 already specified the full range of the identity map. The remaining
144 ranks just need to create empty ISes to cheat VecScatterCreate.
145 */
146 PetscCall(ISCreateStride(PETSC_COMM_SELF, 0, 0, 1, &ix));
147 PetscCall(ISDuplicate(ix, &iy));
148
149 PetscCall(VecScatterCreate(x, ix, yg, iy, &vscat));
150 PetscCall(VecScatterBegin(vscat, x, yg, INSERT_VALUES, SCATTER_FORWARD));
151 PetscCall(VecScatterEnd(vscat, x, yg, INSERT_VALUES, SCATTER_FORWARD));
152
153 /* Send the new y back to x. Ranks outside of subcomm0 actually have nothing to send.
154 But they have to call VecScatterBegin/End since these routines are collective.
155 */
156 PetscCall(VecScatterBegin(vscat, yg, x, INSERT_VALUES, SCATTER_REVERSE));
157 PetscCall(VecScatterEnd(vscat, yg, x, INSERT_VALUES, SCATTER_REVERSE));
158 }
159
160 PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD));
161 PetscCall(ISDestroy(&ix));
162 PetscCall(ISDestroy(&iy));
163 PetscCall(VecDestroy(&x));
164 PetscCall(VecDestroy(&yg));
165 PetscCall(VecScatterDestroy(&vscat));
166 } /* world2sub */
167
168 /*===========================================================================
169 * Transfer a vector x defined on subcomm0 to a vector y defined on
170 * subcomm1. The two subcomms are not overlapping and their union is
171 * not necessarily equal to PETSC_COMM_WORLD.
172 *===========================================================================*/
173 if (sub2sub) {
174 if (mycolor == 0) {
175 /* Intentionally declare N as a local variable so that processes in subcomm1 do not know its value */
176 PetscInt n, N = 22;
177 Vec x, xg, yg;
178 IS ix, iy;
179 VecScatter vscat;
180 const PetscScalar *xvalue;
181 MPI_Comm intercomm, parentcomm;
182 PetscMPIInt lrank;
183
184 PetscCallMPI(MPI_Comm_rank(subcomm, &lrank));
185 /* x is on subcomm */
186 PetscCall(VecCreate(subcomm, &x));
187 PetscCall(VecSetSizes(x, PETSC_DECIDE, N));
188 if (iscuda) {
189 PetscCall(VecSetType(x, VECCUDA));
190 } else {
191 PetscCall(VecSetType(x, VECSTANDARD));
192 }
193 PetscCall(VecSetUp(x));
194 PetscCall(VecGetOwnershipRange(x, &low, &high));
195
196 /* initialize x = [0.0, 1.0, 2.0, ..., 21.0] */
197 for (i = low; i < high; i++) {
198 PetscScalar val = i;
199 PetscCall(VecSetValue(x, i, val, INSERT_VALUES));
200 }
201 PetscCall(VecAssemblyBegin(x));
202 PetscCall(VecAssemblyEnd(x));
203
204 PetscCallMPI(MPI_Intercomm_create(subcomm, 0, PETSC_COMM_WORLD /*peer_comm*/, 1, 100 /*tag*/, &intercomm));
205
206 /* Tell rank 0 of subcomm1 the global size of x */
207 if (!lrank) PetscCallMPI(MPI_Send(&N, 1, MPIU_INT, 0 /*receiver's rank in remote comm, i.e., subcomm1*/, 200 /*tag*/, intercomm));
208
209 /* Create an intracomm PETSc can work on. Ranks in subcomm0 are ordered before ranks in subcomm1 in parentcomm.
210 But this order actually does not matter, since what we care is vector y, which is defined on subcomm1.
211 */
212 PetscCallMPI(MPI_Intercomm_merge(intercomm, 0 /*low*/, &parentcomm));
213
214 /* Create a vector xg on parentcomm, which shares memory with x */
215 PetscCall(VecGetLocalSize(x, &n));
216 if (iscuda) {
217 #if defined(PETSC_HAVE_CUDA)
218 PetscCall(VecCUDAGetArrayRead(x, &xvalue));
219 PetscCall(VecCreateMPICUDAWithArray(parentcomm, 1, n, N, xvalue, &xg));
220 #endif
221 } else {
222 PetscCall(VecGetArrayRead(x, &xvalue));
223 PetscCall(VecCreateMPIWithArray(parentcomm, 1, n, N, xvalue, &xg));
224 }
225
226 /* Ranks in subcomm 0 have nothing on yg, so they simply have n=0, array=NULL */
227 if (iscuda) {
228 #if defined(PETSC_HAVE_CUDA)
229 PetscCall(VecCreateMPICUDAWithArray(parentcomm, 1, 0 /*n*/, N, NULL /*array*/, &yg));
230 #endif
231 } else {
232 PetscCall(VecCreateMPIWithArray(parentcomm, 1, 0 /*n*/, N, NULL /*array*/, &yg));
233 }
234
235 /* Create the vecscatter, which does identity map by setting yg[i] = xg[i], i=0..N-1. */
236 PetscCall(VecGetOwnershipRange(xg, &low, &high)); /* low, high are global indices of xg */
237 PetscCall(ISCreateStride(PETSC_COMM_SELF, high - low, low, 1, &ix));
238 PetscCall(ISDuplicate(ix, &iy));
239 PetscCall(VecScatterCreate(xg, ix, yg, iy, &vscat));
240
241 /* Scatter values from xg to yg */
242 PetscCall(VecScatterBegin(vscat, xg, yg, INSERT_VALUES, SCATTER_FORWARD));
243 PetscCall(VecScatterEnd(vscat, xg, yg, INSERT_VALUES, SCATTER_FORWARD));
244
245 /* After the VecScatter is done, xg is idle so we can safely return xvalue to x */
246 if (iscuda) {
247 #if defined(PETSC_HAVE_CUDA)
248 PetscCall(VecCUDARestoreArrayRead(x, &xvalue));
249 #endif
250 } else {
251 PetscCall(VecRestoreArrayRead(x, &xvalue));
252 }
253 PetscCall(VecDestroy(&x));
254 PetscCall(ISDestroy(&ix));
255 PetscCall(ISDestroy(&iy));
256 PetscCall(VecDestroy(&xg));
257 PetscCall(VecDestroy(&yg));
258 PetscCall(VecScatterDestroy(&vscat));
259 PetscCallMPI(MPI_Comm_free(&intercomm));
260 PetscCallMPI(MPI_Comm_free(&parentcomm));
261 } else if (mycolor == 1) { /* subcomm 1, containing ranks 1, 4, 7, ... in PETSC_COMM_WORLD */
262 PetscInt n, N;
263 Vec y, xg, yg;
264 IS ix, iy;
265 VecScatter vscat;
266 PetscScalar *yvalue;
267 MPI_Comm intercomm, parentcomm;
268 PetscMPIInt lrank;
269
270 PetscCallMPI(MPI_Comm_rank(subcomm, &lrank));
271 PetscCallMPI(MPI_Intercomm_create(subcomm, 0, PETSC_COMM_WORLD /*peer_comm*/, 0 /*remote_leader*/, 100 /*tag*/, &intercomm));
272
273 /* Two rank-0 are talking */
274 if (!lrank) PetscCallMPI(MPI_Recv(&N, 1, MPIU_INT, 0 /*sender's rank in remote comm, i.e. subcomm0*/, 200 /*tag*/, intercomm, MPI_STATUS_IGNORE));
275 /* Rank 0 of subcomm1 bcasts N to its members */
276 PetscCallMPI(MPI_Bcast(&N, 1, MPIU_INT, 0 /*local root*/, subcomm));
277
278 /* Create a intracomm PETSc can work on */
279 PetscCallMPI(MPI_Intercomm_merge(intercomm, 1 /*high*/, &parentcomm));
280
281 /* Ranks in subcomm1 have nothing on xg, so they simply have n=0, array=NULL.*/
282 if (iscuda) {
283 #if defined(PETSC_HAVE_CUDA)
284 PetscCall(VecCreateMPICUDAWithArray(parentcomm, 1 /*bs*/, 0 /*n*/, N, NULL /*array*/, &xg));
285 #endif
286 } else {
287 PetscCall(VecCreateMPIWithArray(parentcomm, 1 /*bs*/, 0 /*n*/, N, NULL /*array*/, &xg));
288 }
289
290 PetscCall(VecCreate(subcomm, &y));
291 PetscCall(VecSetSizes(y, PETSC_DECIDE, N));
292 if (iscuda) {
293 PetscCall(VecSetType(y, VECCUDA));
294 } else {
295 PetscCall(VecSetType(y, VECSTANDARD));
296 }
297 PetscCall(VecSetUp(y));
298
299 PetscCall(PetscObjectSetName((PetscObject)y, "y_subcomm_1")); /* Give a name to view y clearly */
300 PetscCall(VecGetLocalSize(y, &n));
301 if (iscuda) {
302 #if defined(PETSC_HAVE_CUDA)
303 PetscCall(VecCUDAGetArray(y, &yvalue));
304 #endif
305 } else {
306 PetscCall(VecGetArray(y, &yvalue));
307 }
308 /* Create a vector yg on parentcomm, which shares memory with y. xg and yg must be
309 created in the same order in subcomm0/1. For example, we can not reverse the order of
310 creating xg and yg in subcomm1.
311 */
312 if (iscuda) {
313 #if defined(PETSC_HAVE_CUDA)
314 PetscCall(VecCreateMPICUDAWithArray(parentcomm, 1 /*bs*/, n, N, yvalue, &yg));
315 #endif
316 } else {
317 PetscCall(VecCreateMPIWithArray(parentcomm, 1 /*bs*/, n, N, yvalue, &yg));
318 }
319
320 /* Ranks in subcomm0 already specified the full range of the identity map.
321 ranks in subcomm1 just need to create empty ISes to cheat VecScatterCreate.
322 */
323 PetscCall(ISCreateStride(PETSC_COMM_SELF, 0, 0, 1, &ix));
324 PetscCall(ISDuplicate(ix, &iy));
325 PetscCall(VecScatterCreate(xg, ix, yg, iy, &vscat));
326
327 /* Scatter values from xg to yg */
328 PetscCall(VecScatterBegin(vscat, xg, yg, INSERT_VALUES, SCATTER_FORWARD));
329 PetscCall(VecScatterEnd(vscat, xg, yg, INSERT_VALUES, SCATTER_FORWARD));
330
331 /* After the VecScatter is done, values in yg are available. y is our interest, so we return yvalue to y */
332 if (iscuda) {
333 #if defined(PETSC_HAVE_CUDA)
334 PetscCall(VecCUDARestoreArray(y, &yvalue));
335 #endif
336 } else {
337 PetscCall(VecRestoreArray(y, &yvalue));
338 }
339
340 /* Libraries on subcomm1 can safely use y now, for example, view it */
341 PetscCall(VecView(y, PETSC_VIEWER_STDOUT_(subcomm)));
342
343 PetscCall(VecDestroy(&y));
344 PetscCall(ISDestroy(&ix));
345 PetscCall(ISDestroy(&iy));
346 PetscCall(VecDestroy(&xg));
347 PetscCall(VecDestroy(&yg));
348 PetscCall(VecScatterDestroy(&vscat));
349 PetscCallMPI(MPI_Comm_free(&intercomm));
350 PetscCallMPI(MPI_Comm_free(&parentcomm));
351 } else if (mycolor == 2) { /* subcomm2 */
352 /* Processes in subcomm2 do not participate in the VecScatter. They can freely do unrelated things on subcomm2 */
353 }
354 } /* sub2sub */
355
356 /*===========================================================================
357 * Copy a vector x defined on PETSC_COMM_WORLD to vectors y defined on
358 * every subcommunicator of PETSC_COMM_WORLD. We could use multiple transfers
359 * as we did in case 1, but that is not efficient. Instead, we use one vecscatter
360 * to achieve that.
361 *===========================================================================*/
362 if (world2subs) {
363 Vec y;
364 PetscInt n, N = 15, xstart, ystart, low, high;
365 PetscScalar *yvalue;
366
367 /* Initialize x to [0, 1, 2, 3, ..., N-1] */
368 PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
369 PetscCall(VecSetSizes(x, PETSC_DECIDE, N));
370 if (iscuda) {
371 PetscCall(VecSetType(x, VECCUDA));
372 } else {
373 PetscCall(VecSetType(x, VECSTANDARD));
374 }
375 PetscCall(VecSetUp(x));
376 PetscCall(VecGetOwnershipRange(x, &low, &high));
377 for (i = low; i < high; i++) PetscCall(VecSetValue(x, i, (PetscScalar)i, INSERT_VALUES));
378 PetscCall(VecAssemblyBegin(x));
379 PetscCall(VecAssemblyEnd(x));
380
381 /* Every subcomm has a y as long as x */
382 PetscCall(VecCreate(subcomm, &y));
383 PetscCall(VecSetSizes(y, PETSC_DECIDE, N));
384 if (iscuda) {
385 PetscCall(VecSetType(y, VECCUDA));
386 } else {
387 PetscCall(VecSetType(y, VECSTANDARD));
388 }
389 PetscCall(VecSetUp(y));
390 PetscCall(VecGetLocalSize(y, &n));
391
392 /* Create a global vector yg on PETSC_COMM_WORLD using y's memory. yg's global size = N*(number of subcommunicators).
393 Eeach rank in subcomms contributes a piece to construct the global yg. Keep in mind that pieces from a subcomm are not
394 necessarily consecutive in yg. That depends on how PETSC_COMM_WORLD is split. In our case, subcomm0 is made of rank
395 0, 3, 6 etc from PETSC_COMM_WORLD. So subcomm0's pieces are interleaved with pieces from other subcomms in yg.
396 */
397 if (iscuda) {
398 #if defined(PETSC_HAVE_CUDA)
399 PetscCall(VecCUDAGetArray(y, &yvalue));
400 PetscCall(VecCreateMPICUDAWithArray(PETSC_COMM_WORLD, 1, n, PETSC_DECIDE, yvalue, &yg));
401 #endif
402 } else {
403 PetscCall(VecGetArray(y, &yvalue));
404 PetscCall(VecCreateMPIWithArray(PETSC_COMM_WORLD, 1, n, PETSC_DECIDE, yvalue, &yg));
405 }
406 PetscCall(PetscObjectSetName((PetscObject)yg, "yg_on_subcomms")); /* Give a name to view yg clearly */
407
408 /* The following two lines are key. From xstart, we know where to pull entries from x. Note that we get xstart from y,
409 since first entry of y on this rank is from x[xstart]. From ystart, we know where ot put entries to yg.
410 */
411 PetscCall(VecGetOwnershipRange(y, &xstart, NULL));
412 PetscCall(VecGetOwnershipRange(yg, &ystart, NULL));
413
414 PetscCall(ISCreateStride(PETSC_COMM_SELF, n, xstart, 1, &ix));
415 PetscCall(ISCreateStride(PETSC_COMM_SELF, n, ystart, 1, &iy));
416 PetscCall(VecScatterCreate(x, ix, yg, iy, &vscat));
417 PetscCall(VecScatterBegin(vscat, x, yg, INSERT_VALUES, SCATTER_FORWARD));
418 PetscCall(VecScatterEnd(vscat, x, yg, INSERT_VALUES, SCATTER_FORWARD));
419
420 /* View yg on PETSC_COMM_WORLD before destroying it. We shall see the interleaving effect in output. */
421 PetscCall(VecView(yg, PETSC_VIEWER_STDOUT_WORLD));
422 PetscCall(VecDestroy(&yg));
423
424 /* Restory yvalue so that processes in subcomm can use y from now on. */
425 if (iscuda) {
426 #if defined(PETSC_HAVE_CUDA)
427 PetscCall(VecCUDARestoreArray(y, &yvalue));
428 #endif
429 } else {
430 PetscCall(VecRestoreArray(y, &yvalue));
431 }
432 PetscCall(VecScale(y, 3.0));
433
434 PetscCall(ISDestroy(&ix)); /* One can also destroy ix, iy immediately after VecScatterCreate() */
435 PetscCall(ISDestroy(&iy));
436 PetscCall(VecDestroy(&x));
437 PetscCall(VecDestroy(&y));
438 PetscCall(VecScatterDestroy(&vscat));
439 } /* world2subs */
440
441 PetscCallMPI(MPI_Comm_free(&subcomm));
442 PetscCall(PetscFinalize());
443 return 0;
444 }
445
446 /*TEST
447
448 build:
449 requires: !defined(PETSC_HAVE_MPIUNI)
450
451 testset:
452 nsize: 7
453
454 test:
455 suffix: 1
456 args: -world2sub
457
458 test:
459 suffix: 2
460 args: -sub2sub
461 # deadlocks with NECMPI and INTELMPI (20210400300)
462 requires: !defined(PETSC_HAVE_NECMPI) !defined(PETSC_HAVE_I_MPI)
463
464 test:
465 suffix: 3
466 args: -world2subs
467
468 test:
469 suffix: 4
470 args: -world2sub -vectype cuda
471 requires: cuda
472
473 test:
474 suffix: 5
475 args: -sub2sub -vectype cuda
476 requires: cuda
477
478 test:
479 suffix: 6
480 args: -world2subs -vectype cuda
481 requires: cuda
482
483 testset:
484 nsize: 7
485 args: -world2sub -sf_type neighbor
486 output_file: output/ex9_1.out
487 # segfaults with NECMPI
488 requires: defined(PETSC_HAVE_MPI_NEIGHBORHOOD_COLLECTIVES) !defined(PETSC_HAVE_NECMPI)
489
490 test:
491 suffix: 71
492
493 test:
494 suffix: 72
495 requires: defined(PETSC_HAVE_MPI_PERSISTENT_NEIGHBORHOOD_COLLECTIVES)
496 args: -sf_neighbor_persistent
497
498 TEST*/
499