Actual source code: ex52.c
1: static const char help[] = "Tests VecSum()\n\n";
3: #include <petscvec.h>
5: static PetscErrorCode CheckVecSumReturn(Vec v, PetscScalar expected)
6: {
7: PetscScalar actual;
9: PetscFunctionBegin;
10: PetscCall(VecSum(v, &actual));
11: PetscCheck(PetscIsCloseAtTolScalar(actual, expected, 1e-12, 0.0), PETSC_COMM_SELF, PETSC_ERR_PLIB, "VecSum() returned %g + %gi, expected %g + %gi", (double)PetscRealPart(actual), (double)PetscImaginaryPart(actual), (double)PetscRealPart(expected), (double)PetscImaginaryPart(expected));
12: PetscFunctionReturn(PETSC_SUCCESS);
13: }
15: int main(int argc, char **argv)
16: {
17: PetscInt n, lo, N = 10;
18: PetscScalar sum = 0.0;
19: Vec x;
20: PetscScalar *array;
22: PetscFunctionBeginUser;
23: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
25: PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
26: PetscCall(VecSetSizes(x, PETSC_DECIDE, N));
27: PetscCall(VecSetFromOptions(x));
29: PetscCall(VecZeroEntries(x));
30: PetscCall(CheckVecSumReturn(x, 0.0));
32: PetscCall(VecSet(x, 1.0));
33: PetscCall(CheckVecSumReturn(x, 10.0));
35: PetscCall(VecGetLocalSize(x, &n));
36: PetscCall(VecGetOwnershipRange(x, &lo, NULL));
37: PetscCall(VecGetArrayWrite(x, &array));
38: for (PetscInt i = 0; i < n; ++i) array[i] = (PetscScalar)(lo + i);
39: PetscCall(VecRestoreArrayWrite(x, &array));
41: PetscCall(VecGetSize(x, &N));
42: for (PetscInt i = 0; i < N; ++i) sum += (PetscScalar)i;
43: PetscCall(CheckVecSumReturn(x, sum));
45: PetscCall(VecDestroy(&x));
46: PetscCall(PetscFinalize());
47: return 0;
48: }
50: /*TEST
52: testset:
53: output_file: ./output/empty.out
54: nsize: {{1 2}}
55: test:
56: suffix: standard
57: test:
58: requires: defined(PETSC_USE_SHARED_MEMORY)
59: args: -vec_type shared
60: suffix: shared
61: test:
62: requires: viennacl
63: args: -vec_type viennacl
64: suffix: viennacl
65: test:
66: requires: kokkos_kernels
67: args: -vec_type kokkos
68: suffix: kokkos
69: test:
70: requires: cuda
71: args: -vec_type cuda
72: suffix: cuda
73: test:
74: requires: hip
75: args: -vec_type hip
76: suffix: hip
78: TEST*/