Actual source code: ex28.c
petsc-3.15.0 2021-03-30
2: static char help[] = "Tests repeated VecDotBegin()/VecDotEnd().\n\n";
4: #include <petscvec.h>
6: int main(int argc,char **argv)
7: {
9: PetscInt n = 25,i,row0 = 0;
10: PetscScalar two = 2.0,result1,result2,results[40],value,ten = 10.0;
11: PetscScalar result1a,result2a;
12: PetscReal result3,result4,result[2],result3a,result4a,resulta[2];
13: Vec x,y,vecs[40];
15: PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
17: /* create vectors */
18: VecCreate(PETSC_COMM_WORLD,&x);
19: VecSetSizes(x,n,PETSC_DECIDE);
20: VecSetFromOptions(x);
21: VecDuplicate(x,&y);
22: VecSetRandom(x,NULL);
23: VecViewFromOptions(x,NULL,"-x_view");
24: VecSet(y,two);
26: /*
27: Test mixing dot products and norms that require sums
28: */
29: result1 = result2 = 0.0;
30: result3 = result4 = 0.0;
31: VecDotBegin(x,y,&result1);
32: VecDotBegin(y,x,&result2);
33: VecNormBegin(y,NORM_2,&result3);
34: VecNormBegin(x,NORM_1,&result4);
35: PetscCommSplitReductionBegin(PetscObjectComm((PetscObject)x));
36: VecDotEnd(x,y,&result1);
37: VecDotEnd(y,x,&result2);
38: VecNormEnd(y,NORM_2,&result3);
39: VecNormEnd(x,NORM_1,&result4);
41: VecDot(x,y,&result1a);
42: VecDot(y,x,&result2a);
43: VecNorm(y,NORM_2,&result3a);
44: VecNorm(x,NORM_1,&result4a);
46: if (result1 != result1a || result2 != result2a) {
47: PetscPrintf(PETSC_COMM_WORLD,"Error dot: result1 %g result2 %g\n",(double)PetscRealPart(result1),(double)PetscRealPart(result2));
48: }
49: if (result3 != result3a || result4 != result4a) {
50: PetscPrintf(PETSC_COMM_WORLD,"Error 1,2 norms: result3 %g result4 %g\n",(double)result3,(double)result4);
51: }
53: /*
54: Test norms that only require abs
55: */
56: result1 = result2 = 0.0;
57: result3 = result4 = 0.0;
58: VecNormBegin(y,NORM_MAX,&result3);
59: VecNormBegin(x,NORM_MAX,&result4);
60: VecNormEnd(y,NORM_MAX,&result3);
61: VecNormEnd(x,NORM_MAX,&result4);
63: VecNorm(x,NORM_MAX,&result4a);
64: VecNorm(y,NORM_MAX,&result3a);
65: if (result3 != result3a || result4 != result4a) {
66: PetscPrintf(PETSC_COMM_WORLD,"Error max norm: result3 %g result4 %g\n",(double)result3,(double)result4);
67: }
69: /*
70: Tests dot, max, 1, norm
71: */
72: result1 = result2 = 0.0;
73: result3 = result4 = 0.0;
74: VecSetValues(x,1,&row0,&ten,INSERT_VALUES);
75: VecAssemblyBegin(x);
76: VecAssemblyEnd(x);
78: VecDotBegin(x,y,&result1);
79: VecDotBegin(y,x,&result2);
80: VecNormBegin(x,NORM_MAX,&result3);
81: VecNormBegin(x,NORM_1,&result4);
82: VecDotEnd(x,y,&result1);
83: VecDotEnd(y,x,&result2);
84: VecNormEnd(x,NORM_MAX,&result3);
85: VecNormEnd(x,NORM_1,&result4);
87: VecDot(x,y,&result1a);
88: VecDot(y,x,&result2a);
89: VecNorm(x,NORM_MAX,&result3a);
90: VecNorm(x,NORM_1,&result4a);
92: if (result1 != result1a || result2 != result2a) {
93: PetscPrintf(PETSC_COMM_WORLD,"Error dot: result1 %g result2 %g\n",(double)PetscRealPart(result1),(double)PetscRealPart(result2));
94: }
95: if (result3 != result3a || result4 != result4a) {
96: PetscPrintf(PETSC_COMM_WORLD,"Error max 1 norms: result3 %g result4 %g\n",(double)result3,(double)result4);
97: }
99: /*
100: tests 1_and_2 norm
101: */
102: VecNormBegin(x,NORM_MAX,&result3);
103: VecNormBegin(x,NORM_1_AND_2,result);
104: VecNormBegin(y,NORM_MAX,&result4);
105: VecNormEnd(x,NORM_MAX,&result3);
106: VecNormEnd(x,NORM_1_AND_2,result);
107: VecNormEnd(y,NORM_MAX,&result4);
109: VecNorm(x,NORM_MAX,&result3a);
110: VecNorm(x,NORM_1_AND_2,resulta);
111: VecNorm(y,NORM_MAX,&result4a);
112: if (result3 != result3a || result4 != result4a) {
113: PetscPrintf(PETSC_COMM_WORLD,"Error max: result1 %g result2 %g\n",(double)result3,(double)result4);
114: }
115: if (PetscAbsReal(result[0]-resulta[0]) > .01 || PetscAbsReal(result[1]-resulta[1]) > .01) {
116: PetscPrintf(PETSC_COMM_WORLD,"Error 1 and 2 norms: result[0] %g result[1] %g\n",(double)result[0],(double)result[1]);
117: }
119: VecDestroy(&x);
120: VecDestroy(&y);
122: /*
123: Tests computing a large number of operations that require
124: allocating a larger data structure internally
125: */
126: for (i=0; i<40; i++) {
127: VecCreate(PETSC_COMM_WORLD,vecs+i);
128: VecSetSizes(vecs[i],PETSC_DECIDE,n);
129: VecSetFromOptions(vecs[i]);
130: value = (PetscReal)i;
131: VecSet(vecs[i],value);
132: }
133: for (i=0; i<39; i++) {
134: VecDotBegin(vecs[i],vecs[i+1],results+i);
135: }
136: for (i=0; i<39; i++) {
137: VecDotEnd(vecs[i],vecs[i+1],results+i);
138: if (results[i] != 25.0*i*(i+1)) {
139: PetscPrintf(PETSC_COMM_WORLD,"i %D expected %g got %g\n",i,25.0*i*(i+1),(double)PetscRealPart(results[i]));
140: }
141: }
142: for (i=0; i<40; i++) {
143: VecDestroy(&vecs[i]);
144: }
146: PetscFinalize();
147: return ierr;
148: }
150: /*TEST
152: test:
153: nsize: 3
155: testset:
156: nsize: 3
157: output_file: output/ex28_1.out
159: test:
160: suffix: 2
161: args: -splitreduction_async
163: test:
164: suffix: 2_cuda
165: args: -splitreduction_async -vec_type cuda
166: requires: cuda
168: test:
169: suffix: cuda
170: args: -vec_type cuda
171: requires: cuda
173: test:
174: suffix: 2_kokkos
175: args: -splitreduction_async -vec_type kokkos
176: requires: kokkos_kernels
178: test:
179: suffix: kokkos
180: args: -vec_type kokkos
181: requires: kokkos_kernels
182: TEST*/