Actual source code: ex28.c
1: /*$Id: ex28.c,v 1.17 2001/04/10 19:35:02 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: Scalar one = 1.0,two = 2.0,result1,result2,results[40],value,ten = 10.0;
12: Scalar result1a,result2a;
13: double 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,n,PETSC_DECIDE,&x);
20: VecSetFromOptions(x);
21: VecDuplicate(x,&y);
23: VecSet(&one,x);
24: VecSet(&two,y);
26: /*
27: Test mixing dot products and norms that require sums
28: */
29: result1 = result2 = 0.0;
30: result3 = result4 = 0.0;
31: VecDotBegin(x,y,&result1);
32: VecDotBegin(y,x,&result2);
33: VecNormBegin(y,NORM_2,&result3);
34: VecNormBegin(x,NORM_1,&result4);
35: VecDotEnd(x,y,&result1);
36: VecDotEnd(y,x,&result2);
37: VecNormEnd(y,NORM_2,&result3);
38: VecNormEnd(x,NORM_1,&result4);
39:
40: VecDot(x,y,&result1a);
41: VecDot(y,x,&result2a);
42: VecNorm(y,NORM_2,&result3a);
43: VecNorm(x,NORM_1,&result4a);
44:
45: if (result1 != result1a || result2 != result2a) {
46: PetscPrintf(PETSC_COMM_WORLD,"Error dot: result1 %g result2 %gn",PetscRealPart(result1),PetscRealPart(result2));
47: }
48: if (result3 != result3a || result4 != result4a) {
49: PetscPrintf(PETSC_COMM_WORLD,"Error 1,2 norms: result3 %g result4 %gn",result3,result4);
50: }
52: /*
53: Test norms that only require abs
54: */
55: result1 = result2 = 0.0;
56: result3 = result4 = 0.0;
57: VecNormBegin(y,NORM_MAX,&result3);
58: VecNormBegin(x,NORM_MAX,&result4);
59: VecNormEnd(y,NORM_MAX,&result3);
60: VecNormEnd(x,NORM_MAX,&result4);
62: VecNorm(x,NORM_MAX,&result4a);
63: VecNorm(y,NORM_MAX,&result3a);
64: if (result3 != result3a || result4 != result4a) {
65: PetscPrintf(PETSC_COMM_WORLD,"Error max norm: result3 %g result4 %gn",result3,result4);
66: }
68: /*
69: Tests dot, max, 1, norm
70: */
71: result1 = result2 = 0.0;
72: result3 = result4 = 0.0;
73: VecSetValues(x,1,&row0,&ten,INSERT_VALUES);
74: VecAssemblyBegin(x);
75: VecAssemblyEnd(x);
77: VecDotBegin(x,y,&result1);
78: VecDotBegin(y,x,&result2);
79: VecNormBegin(x,NORM_MAX,&result3);
80: VecNormBegin(x,NORM_1,&result4);
81: VecDotEnd(x,y,&result1);
82: VecDotEnd(y,x,&result2);
83: VecNormEnd(x,NORM_MAX,&result3);
84: VecNormEnd(x,NORM_1,&result4);
86: VecDot(x,y,&result1a);
87: VecDot(y,x,&result2a);
88: VecNorm(x,NORM_MAX,&result3a);
89: VecNorm(x,NORM_1,&result4a);
90: if (result1 != result1a || result2 != result2a) {
91: PetscPrintf(PETSC_COMM_WORLD,"Error dot: result1 %g result2 %gn",PetscRealPart(result1),PetscRealPart(result2));
92: }
93: if (result3 != result3a || result4 != result4a) {
94: PetscPrintf(PETSC_COMM_WORLD,"Error max 1 norms: result3 %g result4 %gn",result3,result4);
95: }
97: /*
98: tests 1_and_2 norm
99: */
100: VecNormBegin(x,NORM_MAX,&result3);
101: VecNormBegin(x,NORM_1_AND_2,result);
102: VecNormBegin(y,NORM_MAX,&result4);
103: VecNormEnd(x,NORM_MAX,&result3);
104: VecNormEnd(x,NORM_1_AND_2,result);
105: VecNormEnd(y,NORM_MAX,&result4);
107: VecNorm(x,NORM_MAX,&result3a);
108: VecNorm(x,NORM_1_AND_2,resulta);
109: VecNorm(y,NORM_MAX,&result4a);
110: if (result3 != result3a || result4 != result4a) {
111: PetscPrintf(PETSC_COMM_WORLD,"Error max: result1 %g result2 %gn",result3,result4);
112: }
113: if (PetscAbsDouble(result[0]-resulta[0]) > .01 || PetscAbsDouble(result[1]-resulta[1]) > .01) {
114: PetscPrintf(PETSC_COMM_WORLD,"Error 1 and 2 norms: result[0] %g result[1] %gn",result[0],result[1]);
115: }
117: VecDestroy(x);
118: VecDestroy(y);
120: /*
121: Tests computing a large number of operations that require
122: allocating a larger data structure internally
123: */
124: for (i=0; i<40; i++) {
125: ierr = VecCreate(PETSC_COMM_WORLD,PETSC_DECIDE,n,vecs+i);
126: ierr = VecSetFromOptions(vecs[i]);
127: value = (double)i;
128: ierr = VecSet(&value,vecs[i]);
129: }
130: for (i=0; i<39; i++) {
131: VecDotBegin(vecs[i],vecs[i+1],results+i);
132: }
133: for (i=0; i<39; i++) {
134: VecDotEnd(vecs[i],vecs[i+1],results+i);
135: if (results[i] != 25.0*i*(i+1)) {
136: PetscPrintf(PETSC_COMM_WORLD,"i %d expected %g got %gn",i,25.0*i*(i+1),PetscRealPart(results[i]));
137: }
138: }
139: for (i=0; i<40; i++) {
140: VecDestroy(vecs[i]);
141: }
143: PetscFinalize();
144: return 0;
145: }
146: