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__,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__,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__,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: {
 23:   PetscInt       n = 25,i,row0 = 0;
 24:   PetscScalar    two = 2.0,result1,result2,results[40],value,ten = 10.0;
 25:   PetscScalar    result1a,result2a;
 26:   PetscReal      result3,result4,result[2],result3a,result4a,resulta[2];
 27:   Vec            x,y,vecs[40];
 28:   PetscReal      tol = PETSC_SMALL;

 30:   PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;

 32:   /* create vectors */
 33:   VecCreate(PETSC_COMM_WORLD,&x);
 34:   VecSetSizes(x,n,PETSC_DECIDE);
 35:   VecSetFromOptions(x);
 36:   VecDuplicate(x,&y);
 37:   VecSetRandom(x,NULL);
 38:   VecViewFromOptions(x,NULL,"-x_view");
 39:   VecSet(y,two);

 41:   /*
 42:         Test mixing dot products and norms that require sums
 43:   */
 44:   result1 = result2 = 0.0;
 45:   result3 = result4 = 0.0;
 46:   VecDotBegin(x,y,&result1);
 47:   VecDotBegin(y,x,&result2);
 48:   VecNormBegin(y,NORM_2,&result3);
 49:   VecNormBegin(x,NORM_1,&result4);
 50:   PetscCommSplitReductionBegin(PetscObjectComm((PetscObject)x));
 51:   VecDotEnd(x,y,&result1);
 52:   VecDotEnd(y,x,&result2);
 53:   VecNormEnd(y,NORM_2,&result3);
 54:   VecNormEnd(x,NORM_1,&result4);

 56:   VecDot(x,y,&result1a);
 57:   VecDot(y,x,&result2a);
 58:   VecNorm(y,NORM_2,&result3a);
 59:   VecNorm(x,NORM_1,&result4a);

 61:   CheckErrorScalar(result1,result1a,tol);
 62:   CheckErrorScalar(result2,result2a,tol);
 63:   CheckError(result3,result3a,tol);
 64:   CheckError(result4,result4a,tol);

 66:   /*
 67:         Test norms that only require abs
 68:   */
 69:   result1 = result2 = 0.0;
 70:   result3 = result4 = 0.0;
 71:   VecNormBegin(y,NORM_MAX,&result3);
 72:   VecNormBegin(x,NORM_MAX,&result4);
 73:   VecNormEnd(y,NORM_MAX,&result3);
 74:   VecNormEnd(x,NORM_MAX,&result4);

 76:   VecNorm(x,NORM_MAX,&result4a);
 77:   VecNorm(y,NORM_MAX,&result3a);
 78:   CheckError(result3,result3a,tol);
 79:   CheckError(result4,result4a,tol);

 81:   /*
 82:         Tests dot,  max, 1, norm
 83:   */
 84:   result1 = result2 = 0.0;
 85:   result3 = result4 = 0.0;
 86:   VecSetValues(x,1,&row0,&ten,INSERT_VALUES);
 87:   VecAssemblyBegin(x);
 88:   VecAssemblyEnd(x);

 90:   VecDotBegin(x,y,&result1);
 91:   VecDotBegin(y,x,&result2);
 92:   VecNormBegin(x,NORM_MAX,&result3);
 93:   VecNormBegin(x,NORM_1,&result4);
 94:   VecDotEnd(x,y,&result1);
 95:   VecDotEnd(y,x,&result2);
 96:   VecNormEnd(x,NORM_MAX,&result3);
 97:   VecNormEnd(x,NORM_1,&result4);

 99:   VecDot(x,y,&result1a);
100:   VecDot(y,x,&result2a);
101:   VecNorm(x,NORM_MAX,&result3a);
102:   VecNorm(x,NORM_1,&result4a);

104:   CheckErrorScalar(result1,result1a,tol);
105:   CheckErrorScalar(result2,result2a,tol);
106:   CheckError(result3,result3a,tol);
107:   CheckError(result4,result4a,tol);

109:   /*
110:        tests 1_and_2 norm
111:   */
112:   VecNormBegin(x,NORM_MAX,&result3);
113:   VecNormBegin(x,NORM_1_AND_2,result);
114:   VecNormBegin(y,NORM_MAX,&result4);
115:   VecNormEnd(x,NORM_MAX,&result3);
116:   VecNormEnd(x,NORM_1_AND_2,result);
117:   VecNormEnd(y,NORM_MAX,&result4);

119:   VecNorm(x,NORM_MAX,&result3a);
120:   VecNorm(x,NORM_1_AND_2,resulta);
121:   VecNorm(y,NORM_MAX,&result4a);

123:   CheckError(result3,result3a,tol);
124:   CheckError(result4,result4a,tol);
125:   CheckError(result[0],resulta[0],tol);
126:   CheckError(result[1],resulta[1],tol);

128:   VecDestroy(&x);
129:   VecDestroy(&y);

131:   /*
132:        Tests computing a large number of operations that require
133:     allocating a larger data structure internally
134:   */
135:   for (i=0; i<40; i++) {
136:     VecCreate(PETSC_COMM_WORLD,vecs+i);
137:     VecSetSizes(vecs[i],PETSC_DECIDE,n);
138:     VecSetFromOptions(vecs[i]);
139:     value = (PetscReal)i;
140:     VecSet(vecs[i],value);
141:   }
142:   for (i=0; i<39; i++) {
143:     VecDotBegin(vecs[i],vecs[i+1],results+i);
144:   }
145:   for (i=0; i<39; i++) {
146:     PetscScalar expected = 25.0*i*(i+1);
147:     VecDotEnd(vecs[i],vecs[i+1],results+i);
148:     CheckErrorScalar(results[i],expected,tol);
149:   }
150:   for (i=0; i<40; i++) {
151:     VecDestroy(&vecs[i]);
152:   }

154:   PetscFinalize();
155:   return ierr;
156: }

158: /*TEST

160:    test:
161:       nsize: 3

163:    testset:
164:      nsize: 3
165:      output_file: output/ex28_1.out

167:      test:
168:         suffix: 2
169:         args: -splitreduction_async

171:      test:
172:         suffix: 2_cuda
173:         args: -splitreduction_async -vec_type cuda
174:         requires: cuda

176:      test:
177:         suffix: cuda
178:         args: -vec_type cuda
179:         requires: cuda

181:      test:
182:         suffix: 2_hip
183:         args: -splitreduction_async -vec_type hip
184:         requires: hip

186:      test:
187:         suffix: hip
188:         args: -vec_type hip
189:         requires: hip

191:      test:
192:         suffix: 2_kokkos
193:         args: -splitreduction_async -vec_type kokkos
194:         requires: kokkos_kernels

196:      test:
197:         suffix: kokkos
198:         args: -vec_type kokkos
199:         requires: kokkos_kernels
200: TEST*/