Actual source code: test16.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[] = "Tests a user-defined convergence test.\n\n";
13: #include <slepceps.h>
15: /*
16: MyConvergedAbsolute - Bizarre convergence test that requires more accuracy
17: to positive eigenvalues compared to negative ones.
18: */
19: PetscErrorCode MyConvergedAbsolute(EPS eps,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
20: {
22: *errest = (PetscRealPart(eigr)<0.0)?res:100*res;
23: return(0);
24: }
26: int main(int argc,char **argv)
27: {
28: Mat A; /* problem matrix */
29: EPS eps; /* eigenproblem solver context */
30: PetscInt n=30,i,Istart,Iend;
33: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
34: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
35: PetscPrintf(PETSC_COMM_WORLD,"\nTridiagonal Eigenproblem, n=%D\n\n",n);
37: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
38: Compute the operator matrix that defines the eigensystem, Ax=kx
39: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
41: MatCreate(PETSC_COMM_WORLD,&A);
42: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
43: MatSetFromOptions(A);
44: MatSetUp(A);
46: MatGetOwnershipRange(A,&Istart,&Iend);
47: for (i=Istart;i<Iend;i++) {
48: if (i>0) { MatSetValue(A,i,i-1,-1.0,INSERT_VALUES); }
49: if (i<n-1) { MatSetValue(A,i,i+1,-1.0,INSERT_VALUES); }
50: }
51: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
52: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
53: MatShift(A,-1e-3);
55: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
56: Create the eigensolver
57: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
58: EPSCreate(PETSC_COMM_WORLD,&eps);
59: EPSSetOperators(eps,A,NULL);
60: EPSSetProblemType(eps,EPS_HEP);
61: /* set user-defined convergence test */
62: EPSSetConvergenceTestFunction(eps,MyConvergedAbsolute,NULL,NULL);
63: EPSSetFromOptions(eps);
65: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
66: Solve the problem
67: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
68: EPSSolve(eps);
69: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
71: EPSDestroy(&eps);
72: MatDestroy(&A);
73: SlepcFinalize();
74: return ierr;
75: }
77: /*TEST
79: test:
80: suffix: 1
81: args: -n 200 -eps_nev 6 -eps_ncv 24 -eps_smallest_magnitude
82: requires: !single
84: TEST*/