Actual source code: test10.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: static char help[] = "Test split reductions in BV.\n\n";
13: #include <slepcbv.h>
15: int main(int argc,char **argv)
16: {
18: Vec t,v,w,y,z,zsplit;
19: BV X;
20: PetscInt i,j,n=10,k=5;
21: PetscScalar *zarray;
22: PetscReal nrm;
24: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
25: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
26: PetscOptionsGetInt(NULL,NULL,"-k",&k,NULL);
27: if (k<3) SETERRQ(PETSC_COMM_WORLD,1,"Should specify at least k=3 columns");
28: PetscPrintf(PETSC_COMM_WORLD,"BV split ops (%D columns of dimension %D).\n",k,n);
30: /* Create template vector */
31: VecCreate(PETSC_COMM_WORLD,&t);
32: VecSetSizes(t,PETSC_DECIDE,n);
33: VecSetFromOptions(t);
34: VecDuplicate(t,&v);
35: VecSet(v,1.0);
37: /* Create BV object X */
38: BVCreate(PETSC_COMM_WORLD,&X);
39: PetscObjectSetName((PetscObject)X,"X");
40: BVSetSizesFromVec(X,t,k);
41: BVSetFromOptions(X);
43: /* Fill X entries */
44: for (j=0;j<k;j++) {
45: BVGetColumn(X,j,&w);
46: VecSet(w,0.0);
47: for (i=0;i<4;i++) {
48: if (i+j<n) {
49: VecSetValue(w,i+j,(PetscScalar)(3*i+j-2),INSERT_VALUES);
50: }
51: }
52: VecAssemblyBegin(w);
53: VecAssemblyEnd(w);
54: BVRestoreColumn(X,j,&w);
55: }
57: /* Use regular operations */
58: VecCreateSeq(PETSC_COMM_SELF,k+6,&z);
59: VecGetArray(z,&zarray);
60: BVGetColumn(X,0,&w);
61: VecDot(w,v,zarray);
62: BVRestoreColumn(X,0,&w);
63: BVDotVec(X,v,zarray+1);
64: BVDotColumn(X,2,zarray+1+k);
66: BVGetColumn(X,1,&y);
67: VecNorm(y,NORM_2,&nrm);
68: BVRestoreColumn(X,1,&y);
69: zarray[k+3] = nrm;
70: BVNormVec(X,v,NORM_2,&nrm);
71: zarray[k+4] = nrm;
72: BVNormColumn(X,0,NORM_2,&nrm);
73: zarray[k+5] = nrm;
74: VecRestoreArray(z,&zarray);
76: /* Use split operations */
77: VecCreateSeq(PETSC_COMM_SELF,k+6,&zsplit);
78: VecGetArray(zsplit,&zarray);
79: BVGetColumn(X,0,&w);
80: VecDotBegin(w,v,zarray);
81: BVDotVecBegin(X,v,zarray+1);
82: BVDotColumnBegin(X,2,zarray+1+k);
83: VecDotEnd(w,v,zarray);
84: BVRestoreColumn(X,0,&w);
85: BVDotVecEnd(X,v,zarray+1);
86: BVDotColumnEnd(X,2,zarray+1+k);
88: BVGetColumn(X,1,&y);
89: VecNormBegin(y,NORM_2,&nrm);
90: BVNormVecBegin(X,v,NORM_2,&nrm);
91: BVNormColumnBegin(X,0,NORM_2,&nrm);
92: VecNormEnd(y,NORM_2,&nrm);
93: BVRestoreColumn(X,1,&y);
94: zarray[k+3] = nrm;
95: BVNormVecEnd(X,v,NORM_2,&nrm);
96: zarray[k+4] = nrm;
97: BVNormColumnEnd(X,0,NORM_2,&nrm);
98: zarray[k+5] = nrm;
99: VecRestoreArray(zsplit,&zarray);
101: /* Show difference */
102: VecAXPY(z,-1.0,zsplit);
103: VecNorm(z,NORM_1,&nrm);
104: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"%g\n",(double)nrm);
105: PetscSynchronizedFlush(PETSC_COMM_WORLD,NULL);
107: BVDestroy(&X);
108: VecDestroy(&t);
109: VecDestroy(&v);
110: VecDestroy(&z);
111: VecDestroy(&zsplit);
112: SlepcFinalize();
113: return ierr;
114: }
116: /*TEST
118: test:
119: suffix: 1
120: nsize: 2
121: args: -bv_type {{vecs contiguous svec mat}shared output}
122: output_file: output/test10_1.out
124: test:
125: suffix: 1_cuda
126: nsize: 2
127: args: -bv_type svec -vec_type cuda
128: requires: cuda
129: output_file: output/test10_1.out
131: TEST*/