Actual source code: ex21.c
2: static char help[] = "Tests VecMax() with index.\n\
3: -n <length> : vector length\n\n";
5: #include <petscvec.h>
7: int main(int argc,char **argv)
8: {
10: PetscInt n = 5,idx;
11: PetscReal value,value2;
12: Vec x;
13: PetscScalar one = 1.0;
15: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
16: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
18: /* create vector */
19: VecCreate(PETSC_COMM_WORLD,&x);
20: VecSetSizes(x,PETSC_DECIDE,n);
21: VecSetFromOptions(x);
23: VecSet(x,one);
24: VecSetValue(x,0,0.0,INSERT_VALUES);
25: VecSetValue(x,n-1,2.0,INSERT_VALUES);
26: VecAssemblyBegin(x);
27: VecAssemblyEnd(x);
28: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
29: VecMax(x,&idx,&value);
30: VecMax(x,NULL,&value2);
31: PetscPrintf(PETSC_COMM_WORLD,"Maximum value %g index %D (no index %g)\n",(double)value,idx,(double)value2);
32: VecMin(x,&idx,&value);
33: VecMin(x,NULL,&value2);
34: PetscPrintf(PETSC_COMM_WORLD,"Minimum value %g index %D (no index %g)\n",(double)value,idx,(double)value2);
36: VecDestroy(&x);
38: PetscFinalize();
39: return ierr;
40: }
44: /*TEST
46: test:
47: diff_args: -j
48: filter: grep -v type
49: suffix: 1
51: test:
52: requires: cuda
53: diff_args: -j
54: filter: grep -v type
55: suffix: 1_cuda
56: args: -vec_type cuda
57: output_file: output/ex21_1.out
59: test:
60: requires: kokkos_kernels
61: diff_args: -j
62: filter: grep -v type
63: suffix: 1_kokkos
64: args: -vec_type kokkos
65: output_file: output/ex21_1.out
67: test:
68: diff_args: -j
69: filter: grep -v type
70: suffix: 2
71: nsize: 2
73: test:
74: requires: cuda
75: diff_args: -j
76: filter: grep -v type
77: suffix: 2_cuda
78: nsize: 2
79: args: -vec_type cuda
80: output_file: output/ex21_2.out
82: test:
83: requires: kokkos_kernels
84: diff_args: -j
85: filter: grep -v type
86: suffix: 2_kokkos
87: nsize: 2
88: args: -vec_type kokkos
89: output_file: output/ex21_2.out
91: TEST*/