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.};
14: PetscInitialize(&argc, &argv, (char*)0, help);
15: MPI_Comm_size(PETSC_COMM_WORLD,&size);
17: if (size == 1) {
18: VecCreateSeqViennaCLWithArrays(PETSC_COMM_WORLD,1,n,xHost,NULL,&x);
19: } else {
20: VecCreateMPIViennaCLWithArrays(PETSC_COMM_WORLD,1,n,PETSC_DECIDE,xHost,NULL,&x);
21: }
22: /* print x should be equivalent too xHost */
23: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
24: VecSet(x,42.0);
25: /* print x should be all 42 */
26: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
28: if (size == 1) {
29: VecCreateSeqWithArray(PETSC_COMM_WORLD,1,n,xHost,&y);
30: } else {
31: VecCreateMPIWithArray(PETSC_COMM_WORLD,1,n,PETSC_DECIDE,xHost,&y);
32: }
34: /* print y should be all 42 */
35: VecView(y, PETSC_VIEWER_STDOUT_WORLD);
37: VecDestroy(&y);
38: VecDestroy(&x);
39: 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*/