Actual source code: test10.c
slepc-3.18.2 2023-01-26
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;
24: SlepcInitialize(&argc,&argv,(char*)0,help);
25: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
26: PetscOptionsGetInt(NULL,NULL,"-k",&k,NULL);
28: PetscPrintf(PETSC_COMM_WORLD,"BV split ops (%" PetscInt_FMT " columns of dimension %" PetscInt_FMT ").\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) VecSetValue(w,i+j,(PetscScalar)(3*i+j-2),INSERT_VALUES);
49: }
50: VecAssemblyBegin(w);
51: VecAssemblyEnd(w);
52: BVRestoreColumn(X,j,&w);
53: }
55: /* Use regular operations */
56: VecCreateSeq(PETSC_COMM_SELF,k+6,&z);
57: VecGetArray(z,&zarray);
58: BVGetColumn(X,0,&w);
59: VecDot(w,v,zarray);
60: BVRestoreColumn(X,0,&w);
61: BVDotVec(X,v,zarray+1);
62: BVDotColumn(X,2,zarray+1+k);
64: BVGetColumn(X,1,&y);
65: VecNorm(y,NORM_2,&nrm);
66: BVRestoreColumn(X,1,&y);
67: zarray[k+3] = nrm;
68: BVNormVec(X,v,NORM_2,&nrm);
69: zarray[k+4] = nrm;
70: BVNormColumn(X,0,NORM_2,&nrm);
71: zarray[k+5] = nrm;
72: VecRestoreArray(z,&zarray);
74: /* Use split operations */
75: VecCreateSeq(PETSC_COMM_SELF,k+6,&zsplit);
76: VecGetArray(zsplit,&zarray);
77: BVGetColumn(X,0,&w);
78: VecDotBegin(w,v,zarray);
79: BVDotVecBegin(X,v,zarray+1);
80: BVDotColumnBegin(X,2,zarray+1+k);
81: VecDotEnd(w,v,zarray);
82: BVRestoreColumn(X,0,&w);
83: BVDotVecEnd(X,v,zarray+1);
84: BVDotColumnEnd(X,2,zarray+1+k);
86: BVGetColumn(X,1,&y);
87: VecNormBegin(y,NORM_2,&nrm);
88: BVNormVecBegin(X,v,NORM_2,&nrm);
89: BVNormColumnBegin(X,0,NORM_2,&nrm);
90: VecNormEnd(y,NORM_2,&nrm);
91: BVRestoreColumn(X,1,&y);
92: zarray[k+3] = nrm;
93: BVNormVecEnd(X,v,NORM_2,&nrm);
94: zarray[k+4] = nrm;
95: BVNormColumnEnd(X,0,NORM_2,&nrm);
96: zarray[k+5] = nrm;
97: VecRestoreArray(zsplit,&zarray);
99: /* Show difference */
100: VecAXPY(z,-1.0,zsplit);
101: VecNorm(z,NORM_1,&nrm);
102: PetscSynchronizedPrintf(PETSC_COMM_WORLD,"%g\n",(double)nrm);
103: PetscSynchronizedFlush(PETSC_COMM_WORLD,NULL);
105: BVDestroy(&X);
106: VecDestroy(&t);
107: VecDestroy(&v);
108: VecDestroy(&z);
109: VecDestroy(&zsplit);
110: SlepcFinalize();
111: return 0;
112: }
114: /*TEST
116: testset:
117: nsize: 2
118: output_file: output/test10_1.out
119: test:
120: suffix: 1
121: args: -bv_type {{vecs contiguous svec mat}shared output}
122: test:
123: suffix: 1_cuda
124: args: -bv_type svec -vec_type cuda
125: requires: cuda
127: TEST*/