Actual source code: ex28.c
2: static char help[] = "Tests repeated VecDotBegin()/VecDotEnd().\n\n";
4: #include <petscvec.h>
5: #define CheckError(a,b,tol) do {\
6: if (!PetscIsCloseAtTol(a,b,0,tol)) {\
7: PetscPrintf(PETSC_COMM_WORLD,"Real error at line %d, tol %g: %s %g %s %g diff %g\n",__LINE__,(double)tol,#a,(double)(a),#b,(double)(b),(double)((a)-(b))); \
8: }\
9: } while (0)
11: #define CheckErrorScalar(a,b,tol) do {\
12: if (!PetscIsCloseAtTol(PetscRealPart(a),PetscRealPart(b),0,tol)) {\
13: PetscPrintf(PETSC_COMM_WORLD,"Real error at line %d, tol %g: %s %g %s %g diff %g\n",__LINE__,(double)tol,#a,(double)PetscRealPart(a),#b,(double)PetscRealPart(b),(double)PetscRealPart((a)-(b))); \
14: }\
15: if (!PetscIsCloseAtTol(PetscImaginaryPart(a),PetscImaginaryPart(b),0,PETSC_SMALL)) {\
16: PetscPrintf(PETSC_COMM_WORLD,"Imag error at line %d, tol %g: %s %g %s %g diff %g\n",__LINE__,(double)tol,#a,(double)PetscImaginaryPart(a),#b,(double)PetscImaginaryPart(b),(double)PetscImaginaryPart((a)-(b))); \
17: }\
18: } while (0)
20: int main(int argc,char **argv)
21: {
22: PetscInt n = 25,i,row0 = 0;
23: PetscScalar two = 2.0,result1,result2,results[40],value,ten = 10.0;
24: PetscScalar result1a,result2a;
25: PetscReal result3,result4,result[2],result3a,result4a,resulta[2];
26: Vec x,y,vecs[40];
27: PetscReal tol = PETSC_SMALL;
29: PetscInitialize(&argc,&argv,(char*)0,help);
31: /* create vectors */
32: VecCreate(PETSC_COMM_WORLD,&x);
33: VecSetSizes(x,n,PETSC_DECIDE);
34: VecSetFromOptions(x);
35: VecDuplicate(x,&y);
36: VecSetRandom(x,NULL);
37: VecViewFromOptions(x,NULL,"-x_view");
38: VecSet(y,two);
40: /*
41: Test mixing dot products and norms that require sums
42: */
43: result1 = result2 = 0.0;
44: result3 = result4 = 0.0;
45: VecDotBegin(x,y,&result1);
46: VecDotBegin(y,x,&result2);
47: VecNormBegin(y,NORM_2,&result3);
48: VecNormBegin(x,NORM_1,&result4);
49: PetscCommSplitReductionBegin(PetscObjectComm((PetscObject)x));
50: VecDotEnd(x,y,&result1);
51: VecDotEnd(y,x,&result2);
52: VecNormEnd(y,NORM_2,&result3);
53: VecNormEnd(x,NORM_1,&result4);
55: VecDot(x,y,&result1a);
56: VecDot(y,x,&result2a);
57: VecNorm(y,NORM_2,&result3a);
58: VecNorm(x,NORM_1,&result4a);
60: CheckErrorScalar(result1,result1a,tol);
61: CheckErrorScalar(result2,result2a,tol);
62: CheckError(result3,result3a,tol);
63: CheckError(result4,result4a,tol);
65: /*
66: Test norms that only require abs
67: */
68: result1 = result2 = 0.0;
69: result3 = result4 = 0.0;
70: VecNormBegin(y,NORM_MAX,&result3);
71: VecNormBegin(x,NORM_MAX,&result4);
72: VecNormEnd(y,NORM_MAX,&result3);
73: VecNormEnd(x,NORM_MAX,&result4);
75: VecNorm(x,NORM_MAX,&result4a);
76: VecNorm(y,NORM_MAX,&result3a);
77: CheckError(result3,result3a,tol);
78: CheckError(result4,result4a,tol);
80: /*
81: Tests dot, max, 1, norm
82: */
83: result1 = result2 = 0.0;
84: result3 = result4 = 0.0;
85: VecSetValues(x,1,&row0,&ten,INSERT_VALUES);
86: VecAssemblyBegin(x);
87: VecAssemblyEnd(x);
89: VecDotBegin(x,y,&result1);
90: VecDotBegin(y,x,&result2);
91: VecNormBegin(x,NORM_MAX,&result3);
92: VecNormBegin(x,NORM_1,&result4);
93: VecDotEnd(x,y,&result1);
94: VecDotEnd(y,x,&result2);
95: VecNormEnd(x,NORM_MAX,&result3);
96: VecNormEnd(x,NORM_1,&result4);
98: VecDot(x,y,&result1a);
99: VecDot(y,x,&result2a);
100: VecNorm(x,NORM_MAX,&result3a);
101: VecNorm(x,NORM_1,&result4a);
103: CheckErrorScalar(result1,result1a,tol);
104: CheckErrorScalar(result2,result2a,tol);
105: CheckError(result3,result3a,tol);
106: CheckError(result4,result4a,tol);
108: /*
109: tests 1_and_2 norm
110: */
111: VecNormBegin(x,NORM_MAX,&result3);
112: VecNormBegin(x,NORM_1_AND_2,result);
113: VecNormBegin(y,NORM_MAX,&result4);
114: VecNormEnd(x,NORM_MAX,&result3);
115: VecNormEnd(x,NORM_1_AND_2,result);
116: VecNormEnd(y,NORM_MAX,&result4);
118: VecNorm(x,NORM_MAX,&result3a);
119: VecNorm(x,NORM_1_AND_2,resulta);
120: VecNorm(y,NORM_MAX,&result4a);
122: CheckError(result3,result3a,tol);
123: CheckError(result4,result4a,tol);
124: CheckError(result[0],resulta[0],tol);
125: CheckError(result[1],resulta[1],tol);
127: VecDestroy(&x);
128: VecDestroy(&y);
130: /*
131: Tests computing a large number of operations that require
132: allocating a larger data structure internally
133: */
134: for (i=0; i<40; i++) {
135: VecCreate(PETSC_COMM_WORLD,vecs+i);
136: VecSetSizes(vecs[i],PETSC_DECIDE,n);
137: VecSetFromOptions(vecs[i]);
138: value = (PetscReal)i;
139: VecSet(vecs[i],value);
140: }
141: for (i=0; i<39; i++) {
142: VecDotBegin(vecs[i],vecs[i+1],results+i);
143: }
144: for (i=0; i<39; i++) {
145: PetscScalar expected = 25.0*i*(i+1);
146: VecDotEnd(vecs[i],vecs[i+1],results+i);
147: CheckErrorScalar(results[i],expected,tol);
148: }
149: for (i=0; i<40; i++) {
150: VecDestroy(&vecs[i]);
151: }
153: PetscFinalize();
154: return 0;
155: }
157: /*TEST
159: test:
160: nsize: 3
162: testset:
163: nsize: 3
164: output_file: output/ex28_1.out
166: test:
167: suffix: 2
168: args: -splitreduction_async
170: test:
171: suffix: 2_cuda
172: args: -splitreduction_async -vec_type cuda
173: requires: cuda
175: test:
176: suffix: cuda
177: args: -vec_type cuda
178: requires: cuda
180: test:
181: suffix: 2_hip
182: args: -splitreduction_async -vec_type hip
183: requires: hip
185: test:
186: suffix: hip
187: args: -vec_type hip
188: requires: hip
190: test:
191: suffix: 2_kokkos
192: args: -splitreduction_async -vec_type kokkos
193: requires: kokkos_kernels
195: test:
196: suffix: kokkos
197: args: -vec_type kokkos
198: requires: kokkos_kernels
199: TEST*/