Actual source code: test13.c
slepc-3.10.1 2018-10-23
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2018, 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: }
29: int main(int argc,char **argv)
30: {
31: Mat A; /* problem matrices */
32: EPS eps; /* eigenproblem solver context */
33: PetscScalar seigr,seigi;
34: PetscReal tol=1000*PETSC_MACHINE_EPSILON;
35: Vec sxr,sxi;
36: PetscInt n=30,i,Istart,Iend,nconv;
39: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
41: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
42: PetscPrintf(PETSC_COMM_WORLD,"\nTridiagonal with zero diagonal, n=%D\n\n",n);
44: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
45: Create matrix tridiag([-1 0 -1])
46: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
47: MatCreate(PETSC_COMM_WORLD,&A);
48: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
49: MatSetFromOptions(A);
50: MatSetUp(A);
52: MatGetOwnershipRange(A,&Istart,&Iend);
53: for (i=Istart;i<Iend;i++) {
54: if (i>0) { MatSetValue(A,i,i-1,-1.0,INSERT_VALUES); }
55: if (i<n-1) { MatSetValue(A,i,i+1,-1.0,INSERT_VALUES); }
56: }
57: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
58: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
60: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
61: Create the eigensolver
62: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
63: EPSCreate(PETSC_COMM_WORLD,&eps);
64: EPSSetProblemType(eps,EPS_HEP);
65: EPSSetTolerances(eps,tol,PETSC_DEFAULT);
66: EPSSetOperators(eps,A,NULL);
67: EPSSetWhichEigenpairs(eps,EPS_SMALLEST_REAL);
68: EPSSetFromOptions(eps);
70: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
71: Solve eigenproblem and store some solution
72: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
73: EPSSolve(eps);
74: MatCreateVecs(A,&sxr,NULL);
75: MatCreateVecs(A,&sxi,NULL);
76: EPSGetConverged(eps,&nconv);
77: if (nconv>0) {
78: EPSGetEigenpair(eps,0,&seigr,&seigi,sxr,sxi);
79: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
81: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
82: Solve eigenproblem using an arbitrary selection
83: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
84: EPSSetArbitrarySelection(eps,MyArbitrarySelection,&sxr);
85: EPSSetWhichEigenpairs(eps,EPS_LARGEST_MAGNITUDE);
86: EPSSolve(eps);
87: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
88: } else {
89: PetscPrintf(PETSC_COMM_WORLD,"Problem: no eigenpairs converged.\n");
90: }
92: EPSDestroy(&eps);
93: VecDestroy(&sxr);
94: VecDestroy(&sxi);
95: MatDestroy(&A);
96: SlepcFinalize();
97: return ierr;
98: }
100: /*TEST
102: testset:
103: args: -eps_max_it 5000
104: output_file: output/test13_1.out
105: test:
106: suffix: 1
107: args: -eps_type {{krylovschur gd jd}}
108: test:
109: suffix: 1_gd2
110: args: -eps_type gd -eps_gd_double_expansion
111: test:
112: suffix: 2
113: args: -eps_non_hermitian -eps_type {{krylovschur gd jd}}
114: test:
115: suffix: 2_gd2
116: args: -eps_non_hermitian -eps_type gd -eps_gd_double_expansion
117: timeoutfactor: 2
119: TEST*/