Actual source code: ex59.cxx


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

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

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

 15:   PetscInitialize(&argc, &argv, (char *)0, help);
 16:   MPI_Comm_size(PETSC_COMM_WORLD, &size);

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

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

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

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

 44: /*TEST

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

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

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

 59: TEST*/