Actual source code: test13.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[] = "Test EPSSetArbitrarySelection.\n\n";
13: #include <slepceps.h>
15: PetscErrorCode MyArbitrarySelection(PetscScalar eigr,PetscScalar eigi,Vec xr,Vec xi,PetscScalar *rr,PetscScalar *ri,void *ctx)
16: {
17: PetscErrorCode ierr;
18: Vec xref = *(Vec*)ctx;
21: VecDot(xr,xref,rr);
22: *rr = PetscAbsScalar(*rr);
23: if (ri) *ri = 0.0;
24: return(0);
25: }
27: int main(int argc,char **argv)
28: {
29: Mat A; /* problem matrices */
30: EPS eps; /* eigenproblem solver context */
31: PetscScalar seigr,seigi;
32: PetscReal tol=1000*PETSC_MACHINE_EPSILON;
33: Vec sxr,sxi;
34: PetscInt n=30,i,Istart,Iend,nconv;
37: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
39: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
40: PetscPrintf(PETSC_COMM_WORLD,"\nTridiagonal with zero diagonal, n=%D\n\n",n);
42: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
43: Create matrix tridiag([-1 0 -1])
44: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
45: MatCreate(PETSC_COMM_WORLD,&A);
46: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
47: MatSetFromOptions(A);
48: MatSetUp(A);
50: MatGetOwnershipRange(A,&Istart,&Iend);
51: for (i=Istart;i<Iend;i++) {
52: if (i>0) { MatSetValue(A,i,i-1,-1.0,INSERT_VALUES); }
53: if (i<n-1) { MatSetValue(A,i,i+1,-1.0,INSERT_VALUES); }
54: }
55: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
56: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
58: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
59: Create the eigensolver
60: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
61: EPSCreate(PETSC_COMM_WORLD,&eps);
62: EPSSetProblemType(eps,EPS_HEP);
63: EPSSetTolerances(eps,tol,PETSC_DEFAULT);
64: EPSSetOperators(eps,A,NULL);
65: EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);
66: EPSSetFromOptions(eps);
68: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
69: Solve eigenproblem and store some solution
70: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
71: EPSSolve(eps);
72: MatCreateVecs(A,&sxr,NULL);
73: MatCreateVecs(A,&sxi,NULL);
74: EPSGetConverged(eps,&nconv);
75: if (nconv>0) {
76: EPSGetEigenpair(eps,0,&seigr,&seigi,sxr,sxi);
77: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
79: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
80: Solve eigenproblem using an arbitrary selection
81: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
82: EPSSetArbitrarySelection(eps,MyArbitrarySelection,&sxr);
83: EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE);
84: EPSSolve(eps);
85: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
86: } else {
87: PetscPrintf(PETSC_COMM_WORLD,"Problem: no eigenpairs converged.\n");
88: }
90: EPSDestroy(&eps);
91: VecDestroy(&sxr);
92: VecDestroy(&sxi);
93: MatDestroy(&A);
94: SlepcFinalize();
95: return ierr;
96: }
98: /*TEST
100: testset:
101: args: -eps_max_it 5000
102: output_file: output/test13_1.out
103: test:
104: suffix: 1
105: args: -eps_type {{krylovschur gd jd}}
106: test:
107: suffix: 1_gd2
108: args: -eps_type gd -eps_gd_double_expansion
109: test:
110: suffix: 2
111: args: -eps_non_hermitian -eps_type {{krylovschur gd jd}}
112: test:
113: suffix: 2_gd2
114: args: -eps_non_hermitian -eps_type gd -eps_gd_double_expansion
115: timeoutfactor: 2
117: TEST*/