Actual source code: ex28.c

petsc-3.6.4 2016-04-12
Report Typos and Errors
  2: static char help[] = "Tests repeated VecDotBegin()/VecDotEnd().\n\n";

  4: #include <petscvec.h>

  8: int main(int argc,char **argv)
  9: {
 11:   PetscInt       n   = 25,i,row0 = 0;
 12:   PetscScalar    one = 1.0,two = 2.0,result1,result2,results[40],value,ten = 10.0;
 13:   PetscScalar    result1a,result2a;
 14:   PetscReal      result3,result4,result[2],result3a,result4a,resulta[2];
 15:   Vec            x,y,vecs[40];

 17:   PetscInitialize(&argc,&argv,(char*)0,help);

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

 25:   VecSet(x,one);
 26:   VecSet(y,two);

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

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

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

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

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

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

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

 89:   VecDot(x,y,&result1a);
 90:   VecDot(y,x,&result2a);
 91:   VecNorm(x,NORM_MAX,&result3a);
 92:   VecNorm(x,NORM_1,&result4a);
 93:   if (result1 != result1a || result2 != result2a) {
 94:     PetscPrintf(PETSC_COMM_WORLD,"Error dot: result1 %g result2 %g\n",(double)PetscRealPart(result1),(double)PetscRealPart(result2));
 95:   }
 96:   if (result3 != result3a || result4 != result4a) {
 97:     PetscPrintf(PETSC_COMM_WORLD,"Error max 1 norms: result3 %g result4 %g\n",(double)result3,(double)result4);
 98:   }

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

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

120:   VecDestroy(&x);
121:   VecDestroy(&y);

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

147:   PetscFinalize();
148:   return 0;
149: }