Actual source code: ex57.cxx

  1: #include "petsc.h"
  2: #include "petscviennacl.h"
  3: #include <viennacl/vector.hpp>
  4: typedef viennacl::vector<PetscScalar> ViennaclVector;

  6: int main(int argc,char *argv[])
  7: {
  8:   Vec            x,y;
  9:   PetscInt       n = 5;
 10:   ViennaclVector *x_vcl;

 12:   PetscInitialize(&argc,&argv,(char*)0,NULL);
 13:   VecCreate(PETSC_COMM_WORLD,&x);
 14:   VecSetSizes(x,n,PETSC_DECIDE);
 15:   VecSetType(x,VECVIENNACL);
 16:   VecSet(x,42.0);

 18:   VecViennaCLGetArray(x,&x_vcl);

 20:   VecCreateSeqViennaCLWithArray(PETSC_COMM_WORLD,1,n,(const ViennaclVector *)x_vcl,&y);

 22:   // Operated on 'y', but 'x' would also be changed since both
 23:   // 'x' and 'y' share the same viennacl vector.
 24:   VecScale(y,2.0);

 26:   VecViennaCLRestoreArray(x,&x_vcl);

 28:   // Expected output: 'x' is a 5-vector with all entries as '84'.
 29:   VecView(x,PETSC_VIEWER_STDOUT_WORLD);
 30:   VecDestroy(&y);
 31:   VecDestroy(&x);

 33:   PetscFinalize();
 34:   return 0;
 35: }

 37: /*TEST

 39:    build:
 40:       requires: viennacl defined(PETSC_HAVE_VIENNACL_NO_CUDA)

 42:    test:
 43:       nsize: 1
 44:       suffix: 1
 45:       args: -viennacl_backend opencl -viennacl_opencl_device_type gpu

 47: TEST*/