Actual source code: ex44.c
2: static char help[] = "Test VecConcatenate both in serial and parallel.\n";
4: #include <petscvec.h>
6: int main(int argc, char **args)
7: {
8: Vec *x, x_test, y, y_test;
9: IS *x_is;
10: VecScatter y_to_x, x_to_y;
11: PetscInt i, j, nx, shift, x_size, y_size, *idx;
12: PetscScalar *vals;
13: PetscBool flg, x_equal, y_equal;
15: PetscFunctionBeginUser;
16: PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
17: PetscCall(PetscOptionsGetInt(NULL, NULL, "-nx", &nx, &flg));
18: if (!flg) nx = 3;
20: y_size = 0;
21: shift = 0;
22: PetscCall(PetscMalloc1(nx, &x));
23: for (i = 0; i < nx; i++) {
24: x_size = 2 * (i + 1);
25: y_size += x_size;
26: PetscCall(VecCreate(PETSC_COMM_WORLD, &x[i]));
27: PetscCall(VecSetSizes(x[i], PETSC_DECIDE, x_size));
28: PetscCall(VecSetFromOptions(x[i]));
29: PetscCall(VecSetUp(x[i]));
30: PetscCall(PetscMalloc2(x_size, &idx, x_size, &vals));
31: for (j = 0; j < x_size; j++) {
32: idx[j] = j;
33: vals[j] = (PetscScalar)(shift + j + 1);
34: }
35: shift += x_size;
36: PetscCall(VecSetValues(x[i], x_size, (const PetscInt *)idx, (const PetscScalar *)vals, INSERT_VALUES));
37: PetscCall(VecAssemblyBegin(x[i]));
38: PetscCall(VecAssemblyEnd(x[i]));
39: PetscCall(PetscFree2(idx, vals));
40: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Original X[%" PetscInt_FMT "] vector\n", i));
41: PetscCall(VecView(x[i], PETSC_VIEWER_STDOUT_WORLD));
42: }
43: PetscCall(VecCreate(PETSC_COMM_WORLD, &y));
44: PetscCall(VecSetSizes(y, PETSC_DECIDE, y_size));
45: PetscCall(VecSetFromOptions(y));
46: PetscCall(VecSetUp(y));
47: PetscCall(PetscMalloc2(y_size, &idx, y_size, &vals));
48: for (j = 0; j < y_size; j++) {
49: idx[j] = j;
50: vals[j] = (PetscScalar)(j + 1);
51: }
52: PetscCall(VecSetValues(y, y_size, (const PetscInt *)idx, (const PetscScalar *)vals, INSERT_VALUES));
53: PetscCall(VecAssemblyBegin(y));
54: PetscCall(VecAssemblyEnd(y));
55: PetscCall(PetscFree2(idx, vals));
56: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Expected Y vector\n"));
57: PetscCall(VecView(y, PETSC_VIEWER_STDOUT_WORLD));
58: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n"));
60: /* ---------- base VecConcatenate() test ----------- */
61: PetscCall(VecConcatenate(nx, (const Vec *)x, &y_test, &x_is));
62: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Testing VecConcatenate() for Y = [X[1], X[2], ...]\n"));
63: PetscCall(VecView(y_test, PETSC_VIEWER_STDOUT_WORLD));
64: y_equal = PETSC_FALSE;
65: PetscCall(VecEqual(y_test, y, &y_equal));
66: if (!y_equal) {
67: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " FAIL\n"));
68: } else {
69: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " PASS\n"));
70: }
71: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n"));
73: /* ---------- test VecConcatenate() without IS (checks for dangling memory from IS) ----------- */
74: PetscCall(VecDestroy(&y_test));
75: PetscCall(VecConcatenate(nx, (const Vec *)x, &y_test, NULL));
76: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Testing VecConcatenate() for Y = [X[1], X[2], ...] w/o IS\n"));
77: PetscCall(VecView(y_test, PETSC_VIEWER_STDOUT_WORLD));
78: y_equal = PETSC_FALSE;
79: PetscCall(VecEqual(y_test, y, &y_equal));
80: if (!y_equal) {
81: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " FAIL\n"));
82: } else {
83: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " PASS\n"));
84: }
85: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n"));
87: /* ---------- using index sets on expected Y instead of concatenated Y ----------- */
88: for (i = 0; i < nx; i++) {
89: PetscCall(VecGetSubVector(y, x_is[i], &x_test));
90: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Testing index set for X[%" PetscInt_FMT "] component\n", i));
91: PetscCall(VecView(x_test, PETSC_VIEWER_STDOUT_WORLD));
92: x_equal = PETSC_FALSE;
93: PetscCall(VecEqual(x_test, x[i], &x_equal));
94: if (!x_equal) {
95: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " FAIL\n"));
96: } else {
97: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " PASS\n"));
98: }
99: PetscCall(VecRestoreSubVector(y, x_is[i], &x_test));
100: }
101: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n"));
103: /* ---------- using VecScatter to communicate data from Y to X[i] for all i ----------- */
104: for (i = 0; i < nx; i++) {
105: PetscCall(VecDuplicate(x[i], &x_test));
106: PetscCall(VecZeroEntries(x_test));
107: PetscCall(VecScatterCreate(y, x_is[i], x[i], NULL, &y_to_x));
108: PetscCall(VecScatterBegin(y_to_x, y, x_test, INSERT_VALUES, SCATTER_FORWARD));
109: PetscCall(VecScatterEnd(y_to_x, y, x_test, INSERT_VALUES, SCATTER_FORWARD));
110: PetscCall(VecScatterDestroy(&y_to_x));
111: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Testing VecScatter for Y -> X[%" PetscInt_FMT "]\n", i));
112: PetscCall(VecView(x_test, PETSC_VIEWER_STDOUT_WORLD));
113: x_equal = PETSC_FALSE;
114: PetscCall(VecEqual(x_test, x[i], &x_equal));
115: if (!x_equal) {
116: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " FAIL\n"));
117: } else {
118: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " PASS\n"));
119: }
120: PetscCall(VecDestroy(&x_test));
121: }
122: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "------------------------------------------------------------\n"));
124: /* ---------- using VecScatter to communicate data from X[i] to Y for all i ----------- */
125: PetscCall(VecZeroEntries(y_test));
126: for (i = 0; i < nx; i++) {
127: PetscCall(VecScatterCreate(x[i], NULL, y, x_is[i], &x_to_y));
128: PetscCall(VecScatterBegin(x_to_y, x[i], y_test, INSERT_VALUES, SCATTER_FORWARD));
129: PetscCall(VecScatterEnd(x_to_y, x[i], y_test, INSERT_VALUES, SCATTER_FORWARD));
130: PetscCall(VecScatterDestroy(&x_to_y));
131: }
132: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Testing VecScatter for X[:] -> Y\n"));
133: PetscCall(VecView(y_test, PETSC_VIEWER_STDOUT_WORLD));
134: y_equal = PETSC_FALSE;
135: PetscCall(VecEqual(y_test, y, &y_equal));
136: if (!y_equal) {
137: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " FAIL\n"));
138: } else {
139: PetscCall(PetscPrintf(PETSC_COMM_WORLD, " PASS\n"));
140: }
141: PetscCall(VecDestroy(&y_test));
143: PetscCall(VecDestroy(&y));
144: for (i = 0; i < nx; i++) {
145: PetscCall(VecDestroy(&x[i]));
146: PetscCall(ISDestroy(&x_is[i]));
147: }
148: PetscCall(PetscFree(x));
149: PetscCall(PetscFree(x_is));
150: PetscCall(PetscFinalize());
151: return 0;
152: }
154: /*TEST
156: test:
157: suffix: serial
159: test:
160: suffix: parallel
161: nsize: 2
163: test:
164: suffix: cuda
165: nsize: 2
166: args: -vec_type cuda
167: requires: cuda
169: test:
170: suffix: uneven
171: nsize: 5
173: TEST*/