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