Actual source code: test16.c
slepc-3.17.0 2022-03-31
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[] = "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: {
21: *errest = (PetscRealPart(eigr)<0.0)?res:100*res;
22: PetscFunctionReturn(0);
23: }
25: int main(int argc,char **argv)
26: {
27: Mat A; /* problem matrix */
28: EPS eps; /* eigenproblem solver context */
29: PetscInt n=30,i,Istart,Iend;
31: SlepcInitialize(&argc,&argv,(char*)0,help);
32: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
33: PetscPrintf(PETSC_COMM_WORLD,"\nTridiagonal Eigenproblem, n=%" PetscInt_FMT "\n\n",n);
35: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
36: Compute the operator matrix that defines the eigensystem, Ax=kx
37: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
39: MatCreate(PETSC_COMM_WORLD,&A);
40: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
41: MatSetFromOptions(A);
42: MatSetUp(A);
44: MatGetOwnershipRange(A,&Istart,&Iend);
45: for (i=Istart;i<Iend;i++) {
46: if (i>0) MatSetValue(A,i,i-1,-1.0,INSERT_VALUES);
47: if (i<n-1) MatSetValue(A,i,i+1,-1.0,INSERT_VALUES);
48: }
49: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
50: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
51: MatShift(A,-1e-3);
53: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
54: Create the eigensolver
55: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
56: EPSCreate(PETSC_COMM_WORLD,&eps);
57: EPSSetOperators(eps,A,NULL);
58: EPSSetProblemType(eps,EPS_HEP);
59: /* set user-defined convergence test */
60: EPSSetConvergenceTestFunction(eps,MyConvergedAbsolute,NULL,NULL);
61: EPSSetFromOptions(eps);
63: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
64: Solve the problem
65: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
66: EPSSolve(eps);
67: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
69: EPSDestroy(&eps);
70: MatDestroy(&A);
71: SlepcFinalize();
72: return 0;
73: }
75: /*TEST
77: test:
78: suffix: 1
79: args: -n 200 -eps_nev 6 -eps_ncv 24 -eps_smallest_magnitude
80: requires: !single
82: TEST*/