Actual source code: ex58.c
2: static char help[] = "Test VecCreate{Seq|MPI}CUDAWithArrays.\n\n";
4: #include "petsc.h"
6: int main(int argc, char **argv)
7: {
8: Vec x, y, z;
9: PetscMPIInt size;
10: PetscInt n = 5;
11: PetscScalar xHost[5] = {0., 1., 2., 3., 4.};
12: PetscScalar zHost[5];
13: PetscBool equal;
15: PetscFunctionBeginUser;
16: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
17: PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
19: PetscCall(PetscArraycpy(zHost, xHost, n));
20: PetscCall(VecCreateMPIWithArray(PETSC_COMM_WORLD, 1, n, PETSC_DECIDE, zHost, &z)); /* build z for comparison */
22: if (size == 1) PetscCall(VecCreateSeqCUDAWithArrays(PETSC_COMM_WORLD, 1, n, xHost, NULL, &x));
23: else PetscCall(VecCreateMPICUDAWithArrays(PETSC_COMM_WORLD, 1, n, PETSC_DECIDE, xHost, NULL, &x));
25: PetscCall(VecEqual(z, x, &equal));
26: PetscCheck(equal, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "x, z are different");
28: PetscCall(VecSet(x, 42.0));
29: PetscCall(VecSet(z, 42.0));
30: PetscCall(VecEqual(z, x, &equal));
31: PetscCheck(equal, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "x, z are different");
33: PetscCall(VecCreateMPIWithArray(PETSC_COMM_WORLD, 1, n, PETSC_DECIDE, xHost, &y));
34: PetscCall(VecSetFromOptions(y)); /* changing y's type should not lose its value */
35: PetscCall(VecEqual(z, y, &equal));
36: PetscCheck(equal, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "y, z are different");
38: PetscCall(VecDestroy(&y));
39: PetscCall(VecDestroy(&x));
40: PetscCall(VecDestroy(&z));
41: PetscCall(PetscFinalize());
42: return 0;
43: }
45: /*TEST
47: build:
48: requires: cuda
50: testset:
51: output_file: output/empty.out
52: nsize: {{1 2}}
54: test:
55: suffix: y_host
57: test:
58: TODO: we need something like VecConvert()
59: requires: kokkos_kernels
60: suffix: y_dev
61: args: -vec_type {{standard mpi cuda kokkos}}
62: TEST*/