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: {
9: PetscInt n = 5,idx;
10: PetscReal value,value2;
11: Vec x;
12: PetscScalar one = 1.0;
14: PetscInitialize(&argc,&argv,(char*)0,help);
15: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
17: /* create vector */
18: VecCreate(PETSC_COMM_WORLD,&x);
19: VecSetSizes(x,PETSC_DECIDE,n);
20: VecSetFromOptions(x);
22: VecSet(x,one);
23: VecSetValue(x,0,0.0,INSERT_VALUES);
24: VecSetValue(x,n-1,2.0,INSERT_VALUES);
25: VecAssemblyBegin(x);
26: VecAssemblyEnd(x);
27: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
28: VecMax(x,&idx,&value);
29: VecMax(x,NULL,&value2);
30: PetscPrintf(PETSC_COMM_WORLD,"Maximum value %g index %" PetscInt_FMT " (no index %g)\n",(double)value,idx,(double)value2);
31: VecMin(x,&idx,&value);
32: VecMin(x,NULL,&value2);
33: PetscPrintf(PETSC_COMM_WORLD,"Minimum value %g index %" PetscInt_FMT " (no index %g)\n",(double)value,idx,(double)value2);
35: VecDestroy(&x);
37: PetscFinalize();
38: return 0;
39: }
41: /*TEST
43: testset:
44: diff_args: -j
45: filter: grep -v type | grep -v "MPI processes" | grep -v Process
46: output_file: output/ex21_1.out
48: test:
49: suffix: 1
50: args: -vec_type {{seq mpi}}
52: test:
53: requires: cuda
54: suffix: 1_cuda
55: args: -vec_type {{cuda mpicuda}}
57: test:
58: requires: kokkos_kernels
59: suffix: 1_kokkos
60: args: -vec_type {{kokkos mpikokkos}}
62: test:
63: requires: hip
64: suffix: 1_hip
65: args: -vec_type {{hip mpihip}}
67: testset:
68: diff_args: -j
69: filter: grep -v type
70: output_file: output/ex21_2.out
71: nsize: 2
73: test:
74: suffix: 2
76: test:
77: requires: cuda
78: suffix: 2_cuda
79: args: -vec_type cuda
81: test:
82: requires: kokkos_kernels
83: suffix: 2_kokkos
84: args: -vec_type kokkos
86: test:
87: requires: hip
88: suffix: 2_hip
89: args: -vec_type hip
91: TEST*/