Actual source code: test3.c
slepc-3.7.4 2017-05-17
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2016, 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 multiple calls to EPSSolve with different matrix.\n\n";
24: #include <slepceps.h>
28: int main(int argc,char **argv)
29: {
30: Mat A1,A2; /* problem matrices */
31: EPS eps; /* eigenproblem solver context */
32: PetscReal tol=1000*PETSC_MACHINE_EPSILON,v;
33: Vec d;
34: PetscInt n=30,i,Istart,Iend;
35: PetscRandom myrand;
38: SlepcInitialize(&argc,&argv,(char*)0,help);
40: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
41: PetscPrintf(PETSC_COMM_WORLD,"\nTridiagonal with random diagonal, n=%D\n\n",n);
43: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
44: Create matrix tridiag([-1 0 -1])
45: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
46: MatCreate(PETSC_COMM_WORLD,&A1);
47: MatSetSizes(A1,PETSC_DECIDE,PETSC_DECIDE,n,n);
48: MatSetFromOptions(A1);
49: MatSetUp(A1);
51: MatGetOwnershipRange(A1,&Istart,&Iend);
52: for (i=Istart;i<Iend;i++) {
53: if (i>0) { MatSetValue(A1,i,i-1,-1.0,INSERT_VALUES); }
54: if (i<n-1) { MatSetValue(A1,i,i+1,-1.0,INSERT_VALUES); }
55: }
56: MatAssemblyBegin(A1,MAT_FINAL_ASSEMBLY);
57: MatAssemblyEnd(A1,MAT_FINAL_ASSEMBLY);
59: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
60: Create two matrices by filling the diagonal with rand values
61: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
62: MatDuplicate(A1,MAT_COPY_VALUES,&A2);
63: MatCreateVecs(A1,NULL,&d);
64: PetscRandomCreate(PETSC_COMM_WORLD,&myrand);
65: PetscRandomSetFromOptions(myrand);
66: PetscRandomSetInterval(myrand,0.0,1.0);
67: for (i=Istart;i<Iend;i++) {
68: PetscRandomGetValueReal(myrand,&v);
69: VecSetValue(d,i,v,INSERT_VALUES);
70: }
71: VecAssemblyBegin(d);
72: VecAssemblyEnd(d);
73: MatDiagonalSet(A1,d,INSERT_VALUES);
74: for (i=Istart;i<Iend;i++) {
75: PetscRandomGetValueReal(myrand,&v);
76: VecSetValue(d,i,v,INSERT_VALUES);
77: }
78: VecAssemblyBegin(d);
79: VecAssemblyEnd(d);
80: MatDiagonalSet(A2,d,INSERT_VALUES);
81: VecDestroy(&d);
82: PetscRandomDestroy(&myrand);
84: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
85: Create the eigensolver
86: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
87: EPSCreate(PETSC_COMM_WORLD,&eps);
88: EPSSetProblemType(eps,EPS_HEP);
89: EPSSetTolerances(eps,tol,PETSC_DEFAULT);
90: EPSSetOperators(eps,A1,NULL);
91: EPSSetFromOptions(eps);
93: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
94: Solve first eigenproblem
95: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
96: EPSSolve(eps);
97: PetscPrintf(PETSC_COMM_WORLD," - - - First matrix - - -\n");
98: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
100: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
101: Solve second eigenproblem
102: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
103: EPSSetOperators(eps,A2,NULL);
104: EPSSolve(eps);
105: PetscPrintf(PETSC_COMM_WORLD," - - - Second matrix - - -\n");
106: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
108: EPSDestroy(&eps);
109: MatDestroy(&A1);
110: MatDestroy(&A2);
111: SlepcFinalize();
112: return ierr;
113: }