Actual source code: ex28.c
1: /*$Id: ex28.c,v 1.19 2001/08/07 21:29:25 bsmith Exp $*/
3: static char help[] = "Tests repeated VecDotBegin()/VecDotEnd().nn";
5: #include petscvec.h
6: #include petscsys.h
8: int main(int argc,char **argv)
9: {
10: int ierr,n = 25,i,row0 = 0;
11: PetscScalar one = 1.0,two = 2.0,result1,result2,results[40],value,ten = 10.0;
12: PetscScalar result1a,result2a;
13: PetscReal result3,result4,result[2],result3a,result4a,resulta[2];
14: Vec x,y,vecs[40];
16: PetscInitialize(&argc,&argv,(char*)0,help);
18: /* create vector */
19: VecCreate(PETSC_COMM_WORLD,&x);
20: VecSetSizes(x,n,PETSC_DECIDE);
21: VecSetFromOptions(x);
22: VecDuplicate(x,&y);
24: VecSet(&one,x);
25: VecSet(&two,y);
27: /*
28: Test mixing dot products and norms that require sums
29: */
30: result1 = result2 = 0.0;
31: result3 = result4 = 0.0;
32: VecDotBegin(x,y,&result1);
33: VecDotBegin(y,x,&result2);
34: VecNormBegin(y,NORM_2,&result3);
35: VecNormBegin(x,NORM_1,&result4);
36: VecDotEnd(x,y,&result1);
37: VecDotEnd(y,x,&result2);
38: VecNormEnd(y,NORM_2,&result3);
39: VecNormEnd(x,NORM_1,&result4);
40:
41: VecDot(x,y,&result1a);
42: VecDot(y,x,&result2a);
43: VecNorm(y,NORM_2,&result3a);
44: VecNorm(x,NORM_1,&result4a);
45:
46: if (result1 != result1a || result2 != result2a) {
47: PetscPrintf(PETSC_COMM_WORLD,"Error dot: result1 %g result2 %gn",PetscRealPart(result1),PetscRealPart(result2));
48: }
49: if (result3 != result3a || result4 != result4a) {
50: PetscPrintf(PETSC_COMM_WORLD,"Error 1,2 norms: result3 %g result4 %gn",result3,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 %gn",result3,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);
91: if (result1 != result1a || result2 != result2a) {
92: PetscPrintf(PETSC_COMM_WORLD,"Error dot: result1 %g result2 %gn",PetscRealPart(result1),PetscRealPart(result2));
93: }
94: if (result3 != result3a || result4 != result4a) {
95: PetscPrintf(PETSC_COMM_WORLD,"Error max 1 norms: result3 %g result4 %gn",result3,result4);
96: }
98: /*
99: tests 1_and_2 norm
100: */
101: VecNormBegin(x,NORM_MAX,&result3);
102: VecNormBegin(x,NORM_1_AND_2,result);
103: VecNormBegin(y,NORM_MAX,&result4);
104: VecNormEnd(x,NORM_MAX,&result3);
105: VecNormEnd(x,NORM_1_AND_2,result);
106: VecNormEnd(y,NORM_MAX,&result4);
108: VecNorm(x,NORM_MAX,&result3a);
109: VecNorm(x,NORM_1_AND_2,resulta);
110: VecNorm(y,NORM_MAX,&result4a);
111: if (result3 != result3a || result4 != result4a) {
112: PetscPrintf(PETSC_COMM_WORLD,"Error max: result1 %g result2 %gn",result3,result4);
113: }
114: if (PetscAbsReal(result[0]-resulta[0]) > .01 || PetscAbsReal(result[1]-resulta[1]) > .01) {
115: PetscPrintf(PETSC_COMM_WORLD,"Error 1 and 2 norms: result[0] %g result[1] %gn",result[0],result[1]);
116: }
118: VecDestroy(x);
119: VecDestroy(y);
121: /*
122: Tests computing a large number of operations that require
123: allocating a larger data structure internally
124: */
125: for (i=0; i<40; i++) {
126: ierr = VecCreate(PETSC_COMM_WORLD,vecs+i);
127: ierr = VecSetSizes(vecs[i],PETSC_DECIDE,n);
128: ierr = VecSetFromOptions(vecs[i]);
129: value = (PetscReal)i;
130: ierr = VecSet(&value,vecs[i]);
131: }
132: for (i=0; i<39; i++) {
133: VecDotBegin(vecs[i],vecs[i+1],results+i);
134: }
135: for (i=0; i<39; i++) {
136: VecDotEnd(vecs[i],vecs[i+1],results+i);
137: if (results[i] != 25.0*i*(i+1)) {
138: PetscPrintf(PETSC_COMM_WORLD,"i %d expected %g got %gn",i,25.0*i*(i+1),PetscRealPart(results[i]));
139: }
140: }
141: for (i=0; i<40; i++) {
142: VecDestroy(vecs[i]);
143: }
145: PetscFinalize();
146: return 0;
147: }
148: