Actual source code: test31.c
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 STFILTER interface functions.\n\n"
12: "Based on ex2.\n"
13: "The command line options are:\n"
14: " -n <n>, where <n> = number of grid subdivisions in x dimension.\n"
15: " -m <m>, where <m> = number of grid subdivisions in y dimension.\n\n";
17: #include <slepceps.h>
19: int main(int argc,char **argv)
20: {
21: Mat A;
22: EPS eps;
23: ST st;
24: PetscInt N,n=10,m,Istart,Iend,II,i,j,degree;
25: PetscBool flag,modify=PETSC_FALSE,terse;
26: PetscReal inta,intb,rleft,rright;
27: STFilterType ftype;
28: STFilterDamping damp;
30: PetscFunctionBeginUser;
31: PetscCall(SlepcInitialize(&argc,&argv,NULL,help));
32: PetscCall(PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL));
33: PetscCall(PetscOptionsGetInt(NULL,NULL,"-m",&m,&flag));
34: if (!flag) m=n;
35: N = n*m;
36: PetscCall(PetscPrintf(PETSC_COMM_WORLD,"\n2-D Laplacian Eigenproblem, N=%" PetscInt_FMT " (%" PetscInt_FMT "x%" PetscInt_FMT " grid)\n\n",N,n,m));
37: PetscCall(PetscOptionsGetBool(NULL,NULL,"-modify",&modify,&flag));
39: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
40: Create the 2-D Laplacian
41: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
43: PetscCall(MatCreate(PETSC_COMM_WORLD,&A));
44: PetscCall(MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,N,N));
45: PetscCall(MatSetFromOptions(A));
46: PetscCall(MatGetOwnershipRange(A,&Istart,&Iend));
47: for (II=Istart;II<Iend;II++) {
48: i = II/n; j = II-i*n;
49: if (i>0) PetscCall(MatSetValue(A,II,II-n,-1.0,INSERT_VALUES));
50: if (i<m-1) PetscCall(MatSetValue(A,II,II+n,-1.0,INSERT_VALUES));
51: if (j>0) PetscCall(MatSetValue(A,II,II-1,-1.0,INSERT_VALUES));
52: if (j<n-1) PetscCall(MatSetValue(A,II,II+1,-1.0,INSERT_VALUES));
53: PetscCall(MatSetValue(A,II,II,4.0,INSERT_VALUES));
54: }
55: PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
56: PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
58: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
59: Create the eigensolver and set various options
60: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
62: PetscCall(EPSCreate(PETSC_COMM_WORLD,&eps));
63: PetscCall(EPSSetOperators(eps,A,NULL));
64: PetscCall(EPSSetProblemType(eps,EPS_HEP));
65: PetscCall(EPSSetType(eps,EPSKRYLOVSCHUR));
66: PetscCall(EPSSetWhichEigenpairs(eps,EPS_ALL));
67: PetscCall(EPSSetInterval(eps,0.5,1.3));
68: PetscCall(EPSGetST(eps,&st));
69: PetscCall(STSetType(st,STFILTER));
70: PetscCall(EPSSetFromOptions(eps));
72: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
73: Solve the problem and display the solution
74: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
76: PetscCall(EPSSolve(eps));
78: /* print filter information */
79: PetscCall(PetscObjectTypeCompare((PetscObject)st,STFILTER,&flag));
80: if (flag) {
81: PetscCall(STFilterGetType(st,&ftype));
82: PetscCall(STFilterGetDegree(st,°ree));
83: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Filter type: %s, degree: %" PetscInt_FMT "\n",STFilterTypes[ftype],degree));
84: PetscCall(STFilterGetInterval(st,&inta,&intb));
85: PetscCall(STFilterGetRange(st,&rleft,&rright));
86: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Requested interval: [%g,%g], range: [%g,%g]\n\n",(double)inta,(double)intb,(double)rleft,(double)rright));
87: if (ftype==ST_FILTER_CHEBYSHEV) {
88: PetscCall(STFilterGetDamping(st,&damp));
89: PetscCall(PetscPrintf(PETSC_COMM_WORLD," Damping: %s\n",STFilterDampings[damp]));
90: }
91: }
93: PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse));
94: if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL));
95: else {
96: PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL));
97: PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD));
98: PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD));
99: PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
100: }
102: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
103: Solve the problem again after changing the matrix
104: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
105: if (modify) {
106: PetscCall(MatSetValue(A,0,0,0.3,INSERT_VALUES));
107: PetscCall(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY));
108: PetscCall(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY));
109: PetscCall(EPSSetOperators(eps,A,NULL));
110: PetscCall(EPSSolve(eps));
111: PetscCall(PetscOptionsHasName(NULL,NULL,"-terse",&terse));
112: if (terse) PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL));
113: else {
114: PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_ASCII_INFO_DETAIL));
115: PetscCall(EPSConvergedReasonView(eps,PETSC_VIEWER_STDOUT_WORLD));
116: PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,PETSC_VIEWER_STDOUT_WORLD));
117: PetscCall(PetscViewerPopFormat(PETSC_VIEWER_STDOUT_WORLD));
118: }
119: }
121: PetscCall(EPSDestroy(&eps));
122: PetscCall(MatDestroy(&A));
123: PetscCall(SlepcFinalize());
124: return 0;
125: }
127: /*TEST
129: test:
130: suffix: 1
131: args: -terse
132: filter: sed -e "s/0.161982,7.83797/0.162007,7.83897/"
133: requires: !single
135: test:
136: suffix: 2
137: args: -modify -st_filter_range -0.3,8 -terse
138: requires: !single
140: test:
141: suffix: 3
142: args: -modify -terse
143: filter: sed -e "s/0.161982,7.83797/0.162007,7.83897/"
144: requires: !single
146: test:
147: suffix: 4
148: args: -terse -eps_ncv 18 -st_filter_type chebyshev -st_filter_damping {{fejer jackson}}
149: filter: sed -e "s/0.161982,7.83797/0.162007,7.83897/" | sed -e "s/CHEBYSHEV/FILTLAN/" | grep -v Damping
150: output_file: output/test31_1.out
151: requires: !single
153: test:
154: suffix: 5
155: args: -modify -terse -st_filter_range -0.3,8 -st_filter_type chebyshev -st_filter_damping {{lanczos none}}
156: filter: sed -e "s/CHEBYSHEV/FILTLAN/" | grep -v Damping
157: output_file: output/test31_2.out
158: requires: !single
160: TEST*/