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;
13: PetscInitialize(&argc, &argv, (char *)0, help);
15: /* create vector of length 2*n */
16: VecCreate(PETSC_COMM_SELF, &x);
17: VecSetSizes(x, 3 * n, PETSC_DECIDE);
18: VecSetFromOptions(x);
20: /* create two vectors of length n without array */
21: PetscObjectTypeCompare((PetscObject)x, VECSEQCUDA, &iscuda);
22: PetscObjectTypeCompare((PetscObject)x, VECSEQKOKKOS, &iskokkos);
23: PetscObjectTypeCompare((PetscObject)x, VECSEQHIP, &iship);
24: VecGetBlockSize(x, &bs);
25: if (iscuda) {
26: #if defined(PETSC_HAVE_CUDA)
27: VecCreateSeqCUDAWithArray(PETSC_COMM_SELF, bs, n, NULL, &x1);
28: VecCreateSeqCUDAWithArray(PETSC_COMM_SELF, bs, n, NULL, &x2);
29: VecCreateSeqCUDAWithArray(PETSC_COMM_SELF, bs, n, NULL, &x3);
30: #endif
31: } else if (iskokkos) {
32: #if defined(PETSC_HAVE_KOKKOS_KERNELS)
33: VecCreateSeqKokkosWithArray(PETSC_COMM_SELF, bs, n, NULL, &x1);
34: VecCreateSeqKokkosWithArray(PETSC_COMM_SELF, bs, n, NULL, &x2);
35: VecCreateSeqKokkosWithArray(PETSC_COMM_SELF, bs, n, NULL, &x3);
36: #endif
37: } else if (iship) {
38: #if defined(PETSC_HAVE_HIP)
39: VecCreateSeqHIPWithArray(PETSC_COMM_SELF, bs, n, NULL, &x1);
40: VecCreateSeqHIPWithArray(PETSC_COMM_SELF, bs, n, NULL, &x2);
41: VecCreateSeqHIPWithArray(PETSC_COMM_SELF, bs, n, NULL, &x3);
42: #endif
43: } else {
44: VecCreateSeqWithArray(PETSC_COMM_SELF, bs, n, NULL, &x1);
45: VecCreateSeqWithArray(PETSC_COMM_SELF, bs, n, NULL, &x2);
46: VecCreateSeqWithArray(PETSC_COMM_SELF, bs, n, NULL, &x3);
47: }
49: VecGetArrayWrite(x, &px);
50: VecPlaceArray(x1, px);
51: VecPlaceArray(x2, px + n);
52: VecPlaceArray(x3, px + 2 * n);
53: VecSet(x1, 1.0);
54: VecSet(x2, 0.0);
55: VecSet(x3, 2.0);
56: VecResetArray(x1);
57: VecResetArray(x2);
58: VecResetArray(x3);
59: VecRestoreArrayWrite(x, &px);
60: VecDestroy(&x1);
61: VecDestroy(&x2);
62: VecDestroy(&x3);
64: VecView(x, NULL);
65: VecReciprocal(x);
66: VecView(x, NULL);
68: VecDestroy(&x);
70: 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*/