Actual source code: test10.c
slepc-3.17.1 2022-04-11
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, 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: {
17: Vec t,v,w,y,z,zsplit;
18: BV X;
19: PetscInt i,j,n=10,k=5;
20: PetscScalar *zarray;
21: PetscReal nrm;
23: SlepcInitialize(&argc,&argv,(char*)0,help);
24: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
25: PetscOptionsGetInt(NULL,NULL,"-k",&k,NULL);
27: PetscPrintf(PETSC_COMM_WORLD,"BV split ops (%" PetscInt_FMT " columns of dimension %" PetscInt_FMT ").\n",k,n);
29: /* Create template vector */
30: VecCreate(PETSC_COMM_WORLD,&t);
31: VecSetSizes(t,PETSC_DECIDE,n);
32: VecSetFromOptions(t);
33: VecDuplicate(t,&v);
34: VecSet(v,1.0);
36: /* Create BV object X */
37: BVCreate(PETSC_COMM_WORLD,&X);
38: PetscObjectSetName((PetscObject)X,"X");
39: BVSetSizesFromVec(X,t,k);
40: BVSetFromOptions(X);
42: /* Fill X entries */
43: for (j=0;j<k;j++) {
44: BVGetColumn(X,j,&w);
45: VecSet(w,0.0);
46: for (i=0;i<4;i++) {
47: if (i+j<n) VecSetValue(w,i+j,(PetscScalar)(3*i+j-2),INSERT_VALUES);
48: }
49: VecAssemblyBegin(w);
50: VecAssemblyEnd(w);
51: BVRestoreColumn(X,j,&w);
52: }
54: /* Use regular operations */
55: VecCreateSeq(PETSC_COMM_SELF,k+6,&z);
56: VecGetArray(z,&zarray);
57: BVGetColumn(X,0,&w);
58: VecDot(w,v,zarray);
59: BVRestoreColumn(X,0,&w);
60: BVDotVec(X,v,zarray+1);
61: BVDotColumn(X,2,zarray+1+k);
63: BVGetColumn(X,1,&y);
64: VecNorm(y,NORM_2,&nrm);
65: BVRestoreColumn(X,1,&y);
66: zarray[k+3] = nrm;
67: BVNormVec(X,v,NORM_2,&nrm);
68: zarray[k+4] = nrm;
69: BVNormColumn(X,0,NORM_2,&nrm);
70: zarray[k+5] = nrm;
71: VecRestoreArray(z,&zarray);
73: /* Use split operations */
74: VecCreateSeq(PETSC_COMM_SELF,k+6,&zsplit);
75: VecGetArray(zsplit,&zarray);
76: BVGetColumn(X,0,&w);
77: VecDotBegin(w,v,zarray);
78: BVDotVecBegin(X,v,zarray+1);
79: BVDotColumnBegin(X,2,zarray+1+k);
80: VecDotEnd(w,v,zarray);
81: BVRestoreColumn(X,0,&w);
82: BVDotVecEnd(X,v,zarray+1);
83: BVDotColumnEnd(X,2,zarray+1+k);
85: BVGetColumn(X,1,&y);
86: VecNormBegin(y,NORM_2,&nrm);
87: BVNormVecBegin(X,v,NORM_2,&nrm);
88: BVNormColumnBegin(X,0,NORM_2,&nrm);
89: VecNormEnd(y,NORM_2,&nrm);
90: BVRestoreColumn(X,1,&y);
91: zarray[k+3] = nrm;
92: BVNormVecEnd(X,v,NORM_2,&nrm);
93: zarray[k+4] = nrm;
94: BVNormColumnEnd(X,0,NORM_2,&nrm);
95: zarray[k+5] = nrm;
96: VecRestoreArray(zsplit,&zarray);
98: /* Show difference */
99: VecAXPY(z,-1.0,zsplit);
100: VecNorm(z,NORM_1,&nrm);
101: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"%g\n",(double)nrm);
102: PetscSynchronizedFlush(PETSC_COMM_WORLD,NULL);
104: BVDestroy(&X);
105: VecDestroy(&t);
106: VecDestroy(&v);
107: VecDestroy(&z);
108: VecDestroy(&zsplit);
109: SlepcFinalize();
110: return 0;
111: }
113: /*TEST
115: testset:
116: nsize: 2
117: output_file: output/test10_1.out
118: test:
119: suffix: 1
120: args: -bv_type {{vecs contiguous svec mat}shared output}
121: test:
122: suffix: 1_cuda
123: args: -bv_type svec -vec_type cuda
124: requires: cuda
126: TEST*/