Actual source code: ex60.c
1: static char help[] = "Tests VecPlaceArray() and VecReciprocal().\n\n";
3: #include <petscvec.h>
5: int main(int argc, char **argv)
6: {
7: PetscInt n = 5, bs;
8: PetscBool iscuda, iskokkos, iship;
9: Vec x, x1, x2, x3;
10: PetscScalar *px;
12: PetscFunctionBeginUser;
13: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
15: /* create vector of length 2*n */
16: PetscCall(VecCreate(PETSC_COMM_SELF, &x));
17: PetscCall(VecSetSizes(x, 3 * n, PETSC_DECIDE));
18: PetscCall(VecSetFromOptions(x));
20: /* create two vectors of length n without array */
21: PetscCall(PetscObjectTypeCompare((PetscObject)x, VECSEQCUDA, &iscuda));
22: PetscCall(PetscObjectTypeCompare((PetscObject)x, VECSEQKOKKOS, &iskokkos));
23: PetscCall(PetscObjectTypeCompare((PetscObject)x, VECSEQHIP, &iship));
24: PetscCall(VecGetBlockSize(x, &bs));
25: if (iscuda) {
26: #if defined(PETSC_HAVE_CUDA)
27: PetscCall(VecCreateSeqCUDAWithArray(PETSC_COMM_SELF, bs, n, NULL, &x1));
28: PetscCall(VecCreateSeqCUDAWithArray(PETSC_COMM_SELF, bs, n, NULL, &x2));
29: PetscCall(VecCreateSeqCUDAWithArray(PETSC_COMM_SELF, bs, n, NULL, &x3));
30: #endif
31: } else if (iskokkos) {
32: #if defined(PETSC_HAVE_KOKKOS_KERNELS)
33: PetscCall(VecCreateSeqKokkosWithArray(PETSC_COMM_SELF, bs, n, NULL, &x1));
34: PetscCall(VecCreateSeqKokkosWithArray(PETSC_COMM_SELF, bs, n, NULL, &x2));
35: PetscCall(VecCreateSeqKokkosWithArray(PETSC_COMM_SELF, bs, n, NULL, &x3));
36: #endif
37: } else if (iship) {
38: #if defined(PETSC_HAVE_HIP)
39: PetscCall(VecCreateSeqHIPWithArray(PETSC_COMM_SELF, bs, n, NULL, &x1));
40: PetscCall(VecCreateSeqHIPWithArray(PETSC_COMM_SELF, bs, n, NULL, &x2));
41: PetscCall(VecCreateSeqHIPWithArray(PETSC_COMM_SELF, bs, n, NULL, &x3));
42: #endif
43: } else {
44: PetscCall(VecCreateSeqWithArray(PETSC_COMM_SELF, bs, n, NULL, &x1));
45: PetscCall(VecCreateSeqWithArray(PETSC_COMM_SELF, bs, n, NULL, &x2));
46: PetscCall(VecCreateSeqWithArray(PETSC_COMM_SELF, bs, n, NULL, &x3));
47: }
49: PetscCall(VecGetArrayWrite(x, &px));
50: PetscCall(VecPlaceArray(x1, px));
51: PetscCall(VecPlaceArray(x2, px + n));
52: PetscCall(VecPlaceArray(x3, px + 2 * n));
53: PetscCall(VecSet(x1, 1.0));
54: PetscCall(VecSet(x2, 0.0));
55: PetscCall(VecSet(x3, 2.0));
56: PetscCall(VecResetArray(x1));
57: PetscCall(VecResetArray(x2));
58: PetscCall(VecResetArray(x3));
59: PetscCall(VecRestoreArrayWrite(x, &px));
60: PetscCall(VecDestroy(&x1));
61: PetscCall(VecDestroy(&x2));
62: PetscCall(VecDestroy(&x3));
64: PetscCall(VecView(x, NULL));
65: PetscCall(VecReciprocal(x));
66: PetscCall(VecView(x, NULL));
68: PetscCall(VecDestroy(&x));
70: PetscCall(PetscFinalize());
71: return 0;
72: }
74: /*TEST
76: testset:
77: output_file: output/ex60_1.out
78: diff_args: -j
79: test:
80: suffix: 1
81: test:
82: suffix: 1_cuda
83: args: -vec_type cuda
84: filter: sed -e 's/seqcuda/seq/'
85: requires: cuda
86: test:
87: suffix: 1_kokkos
88: args: -vec_type kokkos
89: filter: sed -e 's/seqkokkos/seq/'
90: requires: kokkos_kernels
91: test:
92: suffix: 1_hip
93: args: -vec_type hip
94: filter: sed -e 's/seqhip/seq/'
95: requires: hip
97: TEST*/