Actual source code: test3.c
slepc-3.9.1 2018-05-02
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[] = "Tests multiple calls to EPSSolve with different matrix.\n\n";
13: #include <slepceps.h>
15: int main(int argc,char **argv)
16: {
17: Mat A1,A2; /* problem matrices */
18: EPS eps; /* eigenproblem solver context */
19: PetscReal tol=1000*PETSC_MACHINE_EPSILON,v;
20: Vec d;
21: PetscInt n=30,i,Istart,Iend;
22: PetscRandom myrand;
25: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
27: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
28: PetscPrintf(PETSC_COMM_WORLD,"\nTridiagonal with random diagonal, n=%D\n\n",n);
30: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
31: Create matrix tridiag([-1 0 -1])
32: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
33: MatCreate(PETSC_COMM_WORLD,&A1);
34: MatSetSizes(A1,PETSC_DECIDE,PETSC_DECIDE,n,n);
35: MatSetFromOptions(A1);
36: MatSetUp(A1);
38: MatGetOwnershipRange(A1,&Istart,&Iend);
39: for (i=Istart;i<Iend;i++) {
40: if (i>0) { MatSetValue(A1,i,i-1,-1.0,INSERT_VALUES); }
41: if (i<n-1) { MatSetValue(A1,i,i+1,-1.0,INSERT_VALUES); }
42: }
43: MatAssemblyBegin(A1,MAT_FINAL_ASSEMBLY);
44: MatAssemblyEnd(A1,MAT_FINAL_ASSEMBLY);
46: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
47: Create two matrices by filling the diagonal with rand values
48: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
49: MatDuplicate(A1,MAT_COPY_VALUES,&A2);
50: MatCreateVecs(A1,NULL,&d);
51: PetscRandomCreate(PETSC_COMM_WORLD,&myrand);
52: PetscRandomSetFromOptions(myrand);
53: PetscRandomSetInterval(myrand,0.0,1.0);
54: for (i=Istart;i<Iend;i++) {
55: PetscRandomGetValueReal(myrand,&v);
56: VecSetValue(d,i,v,INSERT_VALUES);
57: }
58: VecAssemblyBegin(d);
59: VecAssemblyEnd(d);
60: MatDiagonalSet(A1,d,INSERT_VALUES);
61: for (i=Istart;i<Iend;i++) {
62: PetscRandomGetValueReal(myrand,&v);
63: VecSetValue(d,i,v,INSERT_VALUES);
64: }
65: VecAssemblyBegin(d);
66: VecAssemblyEnd(d);
67: MatDiagonalSet(A2,d,INSERT_VALUES);
68: VecDestroy(&d);
69: PetscRandomDestroy(&myrand);
71: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
72: Create the eigensolver
73: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
74: EPSCreate(PETSC_COMM_WORLD,&eps);
75: EPSSetProblemType(eps,EPS_HEP);
76: EPSSetTolerances(eps,tol,PETSC_DEFAULT);
77: EPSSetOperators(eps,A1,NULL);
78: EPSSetFromOptions(eps);
80: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
81: Solve first eigenproblem
82: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
83: EPSSolve(eps);
84: PetscPrintf(PETSC_COMM_WORLD," - - - First matrix - - -\n");
85: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
87: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
88: Solve second eigenproblem
89: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
90: EPSSetOperators(eps,A2,NULL);
91: EPSSolve(eps);
92: PetscPrintf(PETSC_COMM_WORLD," - - - Second matrix - - -\n");
93: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
95: EPSDestroy(&eps);
96: MatDestroy(&A1);
97: MatDestroy(&A2);
98: SlepcFinalize();
99: return ierr;
100: }