Actual source code: ex59.cxx

  1: static char help[] = "Test VecCreate{Seq|MPI}ViennaCLWithArrays.\n\n";

  3: #include "petsc.h"
  4: #include "petscviennacl.h"

  6: int main(int argc, char **argv)
  7: {
  8:   Vec         x, y;
  9:   PetscMPIInt size;
 10:   PetscInt    n        = 5;
 11:   PetscScalar xHost[5] = {0., 1., 2., 3., 4.};

 13:   PetscFunctionBeginUser;
 14:   PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
 15:   PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));

 17:   if (size == 1) {
 18:     PetscCall(VecCreateSeqViennaCLWithArrays(PETSC_COMM_WORLD, 1, n, xHost, NULL, &x));
 19:   } else {
 20:     PetscCall(VecCreateMPIViennaCLWithArrays(PETSC_COMM_WORLD, 1, n, PETSC_DECIDE, xHost, NULL, &x));
 21:   }
 22:   /* print x should be equivalent too xHost */
 23:   PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD));
 24:   PetscCall(VecSet(x, 42.0));
 25:   /* print x should be all 42 */
 26:   PetscCall(VecView(x, PETSC_VIEWER_STDOUT_WORLD));

 28:   if (size == 1) {
 29:     PetscCall(VecCreateSeqWithArray(PETSC_COMM_WORLD, 1, n, xHost, &y));
 30:   } else {
 31:     PetscCall(VecCreateMPIWithArray(PETSC_COMM_WORLD, 1, n, PETSC_DECIDE, xHost, &y));
 32:   }

 34:   /* print y should be all 42 */
 35:   PetscCall(VecView(y, PETSC_VIEWER_STDOUT_WORLD));

 37:   PetscCall(VecDestroy(&y));
 38:   PetscCall(VecDestroy(&x));
 39:   PetscCall(PetscFinalize());
 40:   return 0;
 41: }

 43: /*TEST

 45:    build:
 46:       requires: viennacl defined(PETSC_HAVE_VIENNACL_NO_CUDA)

 48:    test:
 49:       nsize: 1
 50:       suffix: 1
 51:       args: -viennacl_backend opencl

 53:    test:
 54:       nsize: 2
 55:       suffix: 2
 56:       args: -viennacl_backend opencl

 58: TEST*/