Actual source code: ex28.c

petsc-3.12.2 2019-11-22
Report Typos and Errors

  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];
 14:   PetscRandom    rctx;

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

 18:   /* create vector */
 19:   VecCreate(PETSC_COMM_WORLD,&x);
 20:   VecSetSizes(x,n,PETSC_DECIDE);
 21:   VecSetFromOptions(x);
 22:   VecDuplicate(x,&y);

 24:   PetscRandomCreate(PETSC_COMM_WORLD,&rctx);
 25:   PetscRandomSetFromOptions(rctx);
 26:   VecSetRandom(x,rctx);
 27:   PetscRandomDestroy(&rctx);
 28:   VecSet(y,two);

 30:   /*
 31:         Test mixing dot products and norms that require sums
 32:   */
 33:   result1 = result2 = 0.0;
 34:   result3 = result4 = 0.0;
 35:   VecDotBegin(x,y,&result1);
 36:   VecDotBegin(y,x,&result2);
 37:   VecNormBegin(y,NORM_2,&result3);
 38:   VecNormBegin(x,NORM_1,&result4);
 39:   PetscCommSplitReductionBegin(PetscObjectComm((PetscObject)x));
 40:   VecDotEnd(x,y,&result1);
 41:   VecDotEnd(y,x,&result2);
 42:   VecNormEnd(y,NORM_2,&result3);
 43:   VecNormEnd(x,NORM_1,&result4);

 45:   VecDot(x,y,&result1a);
 46:   VecDot(y,x,&result2a);
 47:   VecNorm(y,NORM_2,&result3a);
 48:   VecNorm(x,NORM_1,&result4a);

 50:   if (result1 != result1a || result2 != result2a) {
 51:     PetscPrintf(PETSC_COMM_WORLD,"Error dot: result1 %g result2 %g\n",(double)PetscRealPart(result1),(double)PetscRealPart(result2));
 52:   }
 53:   if (result3 != result3a || result4 != result4a) {
 54:     PetscPrintf(PETSC_COMM_WORLD,"Error 1,2 norms: result3 %g result4 %g\n",(double)result3,(double)result4);
 55:   }

 57:   /*
 58:         Test norms that only require abs
 59:   */
 60:   result1 = result2 = 0.0;
 61:   result3 = result4 = 0.0;
 62:   VecNormBegin(y,NORM_MAX,&result3);
 63:   VecNormBegin(x,NORM_MAX,&result4);
 64:   VecNormEnd(y,NORM_MAX,&result3);
 65:   VecNormEnd(x,NORM_MAX,&result4);

 67:   VecNorm(x,NORM_MAX,&result4a);
 68:   VecNorm(y,NORM_MAX,&result3a);
 69:   if (result3 != result3a || result4 != result4a) {
 70:     PetscPrintf(PETSC_COMM_WORLD,"Error max norm: result3 %g result4 %g\n",(double)result3,(double)result4);
 71:   }

 73:   /*
 74:         Tests dot,  max, 1, norm
 75:   */
 76:   result1 = result2 = 0.0;
 77:   result3 = result4 = 0.0;
 78:   VecSetValues(x,1,&row0,&ten,INSERT_VALUES);
 79:   VecAssemblyBegin(x);
 80:   VecAssemblyEnd(x);

 82:   VecDotBegin(x,y,&result1);
 83:   VecDotBegin(y,x,&result2);
 84:   VecNormBegin(x,NORM_MAX,&result3);
 85:   VecNormBegin(x,NORM_1,&result4);
 86:   VecDotEnd(x,y,&result1);
 87:   VecDotEnd(y,x,&result2);
 88:   VecNormEnd(x,NORM_MAX,&result3);
 89:   VecNormEnd(x,NORM_1,&result4);

 91:   VecDot(x,y,&result1a);
 92:   VecDot(y,x,&result2a);
 93:   VecNorm(x,NORM_MAX,&result3a);
 94:   VecNorm(x,NORM_1,&result4a);

 96:   if (result1 != result1a || result2 != result2a) {
 97:     PetscPrintf(PETSC_COMM_WORLD,"Error dot: result1 %g result2 %g\n",(double)PetscRealPart(result1),(double)PetscRealPart(result2));
 98:   }
 99:   if (result3 != result3a || result4 != result4a) {
100:     PetscPrintf(PETSC_COMM_WORLD,"Error max 1 norms: result3 %g result4 %g\n",(double)result3,(double)result4);
101:   }

103:   /*
104:        tests 1_and_2 norm
105:   */
106:   VecNormBegin(x,NORM_MAX,&result3);
107:   VecNormBegin(x,NORM_1_AND_2,result);
108:   VecNormBegin(y,NORM_MAX,&result4);
109:   VecNormEnd(x,NORM_MAX,&result3);
110:   VecNormEnd(x,NORM_1_AND_2,result);
111:   VecNormEnd(y,NORM_MAX,&result4);

113:   VecNorm(x,NORM_MAX,&result3a);
114:   VecNorm(x,NORM_1_AND_2,resulta);
115:   VecNorm(y,NORM_MAX,&result4a);
116:   if (result3 != result3a || result4 != result4a) {
117:     PetscPrintf(PETSC_COMM_WORLD,"Error max: result1 %g result2 %g\n",(double)result3,(double)result4);
118:   }
119:   if (PetscAbsReal(result[0]-resulta[0]) > .01 || PetscAbsReal(result[1]-resulta[1]) > .01) {
120:     PetscPrintf(PETSC_COMM_WORLD,"Error 1 and 2 norms: result[0] %g result[1] %g\n",(double)result[0],(double)result[1]);
121:   }

123:   VecDestroy(&x);
124:   VecDestroy(&y);

126:   /*
127:        Tests computing a large number of operations that require
128:     allocating a larger data structure internally
129:   */
130:   for (i=0; i<40; i++) {
131:     VecCreate(PETSC_COMM_WORLD,vecs+i);
132:     VecSetSizes(vecs[i],PETSC_DECIDE,n);
133:     VecSetFromOptions(vecs[i]);
134:     value = (PetscReal)i;
135:     VecSet(vecs[i],value);
136:   }
137:   for (i=0; i<39; i++) {
138:     VecDotBegin(vecs[i],vecs[i+1],results+i);
139:   }
140:   for (i=0; i<39; i++) {
141:     VecDotEnd(vecs[i],vecs[i+1],results+i);
142:     if (results[i] != 25.0*i*(i+1)) {
143:       PetscPrintf(PETSC_COMM_WORLD,"i %D expected %g got %g\n",i,25.0*i*(i+1),(double)PetscRealPart(results[i]));
144:     }
145:   }
146:   for (i=0; i<40; i++) {
147:     VecDestroy(&vecs[i]);
148:   }

150:   PetscFinalize();
151:   return ierr;
152: }






159: /*TEST

161:    test:
162:       nsize: 3

164:    test:
165:       suffix: 2
166:       nsize: 3
167:       args: -splitreduction_async
168:       output_file: output/ex28_1.out

170:    test:
171:       suffix: 2_cuda
172:       nsize: 3
173:       args: -splitreduction_async -vec_type cuda
174:       requires: cuda
175:       output_file: output/ex28_1.out

177:    test:
178:       suffix: cuda
179:       nsize: 3
180:       args: -vec_type cuda
181:       requires: cuda
182:       output_file: output/ex28_1.out
183:  
184: TEST*/