Actual source code: test16.c
slepc-3.5.2 2014-10-10
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2014, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
8: SLEPc is free software: you can redistribute it and/or modify it under the
9: terms of version 3 of the GNU Lesser General Public License as published by
10: the Free Software Foundation.
12: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
13: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
15: more details.
17: You should have received a copy of the GNU Lesser General Public License
18: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
19: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: */
22: static char help[] = "Tests a user-defined convergence test.\n\n";
24: #include <slepceps.h>
29: /*
30: MyConvergedAbsolute - Bizarre convergence test that requires more accuracy
31: to positive eigenvalues compared to negative ones.
32: */
33: PetscErrorCode MyConvergedAbsolute(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
34: {
36: *errest = (PetscRealPart(eigr)<0.0)?res:100*res;
37: return(0);
38: }
42: int main(int argc,char **argv)
43: {
44: Mat A; /* problem matrix */
45: EPS eps; /* eigenproblem solver context */
46: PetscScalar value[2];
47: PetscInt n=30,i,Istart,Iend,col[2];
48: PetscBool FirstBlock=PETSC_FALSE,LastBlock=PETSC_FALSE;
51: SlepcInitialize(&argc,&argv,(char*)0,help);
52: PetscOptionsGetInt(NULL,"-n",&n,NULL);
53: PetscPrintf(PETSC_COMM_WORLD,"\nTridiagonal Eigenproblem, n=%D\n\n",n);
55: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
56: Compute the operator matrix that defines the eigensystem, Ax=kx
57: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
59: MatCreate(PETSC_COMM_WORLD,&A);
60: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
61: MatSetFromOptions(A);
62: MatSetUp(A);
64: MatGetOwnershipRange(A,&Istart,&Iend);
65: if (Istart==0) FirstBlock=PETSC_TRUE;
66: if (Iend==n) LastBlock=PETSC_TRUE;
67: value[0]=-1.0; value[1]=-1.0;
68: for (i=(FirstBlock? Istart+1: Istart); i<(LastBlock? Iend-1: Iend); i++) {
69: col[0]=i-1; col[1]=i+1;
70: MatSetValues(A,1,&i,2,col,value,INSERT_VALUES);
71: }
72: if (LastBlock) {
73: i=n-1; col[0]=n-2;
74: MatSetValues(A,1,&i,1,col,value,INSERT_VALUES);
75: }
76: if (FirstBlock) {
77: i=0; col[0]=1;
78: MatSetValues(A,1,&i,1,col,value,INSERT_VALUES);
79: }
81: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
82: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
84: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
85: Create the eigensolver
86: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
87: EPSCreate(PETSC_COMM_WORLD,&eps);
88: EPSSetOperators(eps,A,NULL);
89: EPSSetProblemType(eps,EPS_HEP);
90: /* set user-defined convergence test */
91: EPSSetConvergenceTestFunction(eps,MyConvergedAbsolute,NULL,NULL);
92: EPSSetFromOptions(eps);
94: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
95: Solve the problem
96: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
97: EPSSolve(eps);
98: EPSPrintSolution(eps,NULL);
100: EPSDestroy(&eps);
101: MatDestroy(&A);
102: SlepcFinalize();
103: return 0;
104: }