Actual source code: test3.c
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2011, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7:
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(PETSC_NULL,"-n",&n,PETSC_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:
52: MatGetOwnershipRange(A1,&Istart,&Iend);
53: if (Istart==0) FirstBlock=PETSC_TRUE;
54: if (Iend==n) LastBlock=PETSC_TRUE;
55: value[0]=-1.0; value[1]=0.0; value[2]=-1.0;
56: for (i=(FirstBlock? Istart+1: Istart); i<(LastBlock? Iend-1: Iend); i++) {
57: col[0]=i-1; col[1]=i; col[2]=i+1;
58: MatSetValues(A1,1,&i,3,col,value,INSERT_VALUES);
59: }
60: if (LastBlock) {
61: i=n-1; col[0]=n-2; col[1]=n-1;
62: MatSetValues(A1,1,&i,2,col,value,INSERT_VALUES);
63: }
64: if (FirstBlock) {
65: i=0; col[0]=0; col[1]=1; value[0]=0.0; value[1]=-1.0;
66: MatSetValues(A1,1,&i,2,col,value,INSERT_VALUES);
67: }
69: MatAssemblyBegin(A1,MAT_FINAL_ASSEMBLY);
70: MatAssemblyEnd(A1,MAT_FINAL_ASSEMBLY);
72: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
73: Create two matrices by filling the diagonal with rand values
74: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
75: MatDuplicate(A1,MAT_COPY_VALUES,&A2);
76: MatGetVecs(A1,PETSC_NULL,&d);
77: PetscRandomCreate(PETSC_COMM_WORLD,&myrand);
78: PetscRandomSetFromOptions(myrand);
79: PetscRandomSetInterval(myrand,0.0,1.0);
80: for (i=0; i<n; i++) {
81: PetscRandomGetValueReal(myrand,&v);
82: VecSetValue(d,i,v,INSERT_VALUES);
83: }
84: VecAssemblyBegin(d);
85: VecAssemblyEnd(d);
86: MatDiagonalSet(A1,d,INSERT_VALUES);
87: for (i=0; i<n; i++) {
88: PetscRandomGetValueReal(myrand,&v);
89: VecSetValue(d,i,v,INSERT_VALUES);
90: }
91: VecAssemblyBegin(d);
92: VecAssemblyEnd(d);
93: MatDiagonalSet(A2,d,INSERT_VALUES);
94: VecDestroy(&d);
95: PetscRandomDestroy(&myrand);
97: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
98: Create the eigensolver
99: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
100: EPSCreate(PETSC_COMM_WORLD,&eps);
101: EPSSetProblemType(eps,EPS_HEP);
102: EPSSetTolerances(eps,tol,PETSC_DECIDE);
103: EPSSetOperators(eps,A1,PETSC_NULL);
104: EPSSetFromOptions(eps);
106: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
107: Solve first eigenproblem
108: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
109: EPSSolve(eps);
110: PetscPrintf(PETSC_COMM_WORLD," - - - First matrix - - -\n");
111: EPSPrintSolution(eps,PETSC_NULL);
113: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
114: Solve second eigenproblem
115: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
116: EPSSetOperators(eps,A2,PETSC_NULL);
117: EPSSolve(eps);
118: PetscPrintf(PETSC_COMM_WORLD," - - - Second matrix - - -\n");
119: EPSPrintSolution(eps,PETSC_NULL);
120:
121: EPSDestroy(&eps);
122: MatDestroy(&A1);
123: MatDestroy(&A2);
124: SlepcFinalize();
125: return 0;
126: }