Actual source code: test3.c
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2013, 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: PetscScalar value[3];
33: PetscReal tol=1000*PETSC_MACHINE_EPSILON,v;
34: Vec d;
35: PetscInt n=30,i,Istart,Iend,col[3];
36: PetscBool FirstBlock=PETSC_FALSE,LastBlock=PETSC_FALSE;
37: PetscRandom myrand;
40: SlepcInitialize(&argc,&argv,(char*)0,help);
42: PetscOptionsGetInt(NULL,"-n",&n,NULL);
43: PetscPrintf(PETSC_COMM_WORLD,"\nTridiagonal with random diagonal, n=%D\n\n",n);
45: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
46: Create matrix tridiag([-1 0 -1])
47: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
48: MatCreate(PETSC_COMM_WORLD,&A1);
49: MatSetSizes(A1,PETSC_DECIDE,PETSC_DECIDE,n,n);
50: MatSetFromOptions(A1);
51: MatSetUp(A1);
53: MatGetOwnershipRange(A1,&Istart,&Iend);
54: if (Istart==0) FirstBlock=PETSC_TRUE;
55: if (Iend==n) LastBlock=PETSC_TRUE;
56: value[0]=-1.0; value[1]=0.0; value[2]=-1.0;
57: for (i=(FirstBlock? Istart+1: Istart); i<(LastBlock? Iend-1: Iend); i++) {
58: col[0]=i-1; col[1]=i; col[2]=i+1;
59: MatSetValues(A1,1,&i,3,col,value,INSERT_VALUES);
60: }
61: if (LastBlock) {
62: i=n-1; col[0]=n-2; col[1]=n-1;
63: MatSetValues(A1,1,&i,2,col,value,INSERT_VALUES);
64: }
65: if (FirstBlock) {
66: i=0; col[0]=0; col[1]=1; value[0]=0.0; value[1]=-1.0;
67: MatSetValues(A1,1,&i,2,col,value,INSERT_VALUES);
68: }
70: MatAssemblyBegin(A1,MAT_FINAL_ASSEMBLY);
71: MatAssemblyEnd(A1,MAT_FINAL_ASSEMBLY);
73: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
74: Create two matrices by filling the diagonal with rand values
75: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
76: MatDuplicate(A1,MAT_COPY_VALUES,&A2);
77: MatGetVecs(A1,NULL,&d);
78: PetscRandomCreate(PETSC_COMM_WORLD,&myrand);
79: PetscRandomSetFromOptions(myrand);
80: PetscRandomSetInterval(myrand,0.0,1.0);
81: for (i=0; i<n; i++) {
82: PetscRandomGetValueReal(myrand,&v);
83: VecSetValue(d,i,v,INSERT_VALUES);
84: }
85: VecAssemblyBegin(d);
86: VecAssemblyEnd(d);
87: MatDiagonalSet(A1,d,INSERT_VALUES);
88: for (i=0; i<n; i++) {
89: PetscRandomGetValueReal(myrand,&v);
90: VecSetValue(d,i,v,INSERT_VALUES);
91: }
92: VecAssemblyBegin(d);
93: VecAssemblyEnd(d);
94: MatDiagonalSet(A2,d,INSERT_VALUES);
95: VecDestroy(&d);
96: PetscRandomDestroy(&myrand);
98: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
99: Create the eigensolver
100: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
101: EPSCreate(PETSC_COMM_WORLD,&eps);
102: EPSSetProblemType(eps,EPS_HEP);
103: EPSSetTolerances(eps,tol,PETSC_DECIDE);
104: EPSSetOperators(eps,A1,NULL);
105: EPSSetFromOptions(eps);
107: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
108: Solve first eigenproblem
109: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
110: EPSSolve(eps);
111: PetscPrintf(PETSC_COMM_WORLD," - - - First matrix - - -\n");
112: EPSPrintSolution(eps,NULL);
114: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
115: Solve second eigenproblem
116: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
117: EPSSetOperators(eps,A2,NULL);
118: EPSSolve(eps);
119: PetscPrintf(PETSC_COMM_WORLD," - - - Second matrix - - -\n");
120: EPSPrintSolution(eps,NULL);
122: EPSDestroy(&eps);
123: MatDestroy(&A1);
124: MatDestroy(&A2);
125: SlepcFinalize();
126: return 0;
127: }