Actual source code: test6.c
2: static char help[] = "Diagonal eigenproblem.\n\n"
3: "The command line options are:\n"
4: " -n <n>, where <n> = number of grid subdivisions = matrix dimension.\n"
5: " -seed <s>, where <s> = seed for random number generation.\n\n";
7: #include <slepceps.h>
11: int main(int argc,char **argv)
12: {
13: Mat A; /* problem matrix */
14: EPS eps; /* eigenproblem solver context */
15: Vec v0; /* initial vector */
16: PetscRandom rand;
17: PetscReal tol=1000*PETSC_MACHINE_EPSILON;
18: PetscInt n=30,i,Istart,Iend,seed=0x12345678;
21: SlepcInitialize(&argc,&argv,(char*)0,help);
23: PetscOptionsGetInt(PETSC_NULL,"-n",&n,PETSC_NULL);
24: PetscPrintf(PETSC_COMM_WORLD,"\nDiagonal Eigenproblem, n=%d\n\n",n);
26: MatCreate(PETSC_COMM_WORLD,&A);
27: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
28: MatSetFromOptions(A);
29: MatGetOwnershipRange(A,&Istart,&Iend);
30: for (i=Istart;i<Iend;i++) {
31: MatSetValue(A,i,i,i+1,INSERT_VALUES);
32: }
33: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
34: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
36: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
37: Solve the eigensystem
38: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
39: EPSCreate(PETSC_COMM_WORLD,&eps);
40: EPSSetOperators(eps,A,PETSC_NULL);
41: EPSSetProblemType(eps,EPS_HEP);
42: EPSSetTolerances(eps,tol,PETSC_DECIDE);
43: EPSSetFromOptions(eps);
44: /* set random initial vector */
45: MatGetVecs(A,&v0,PETSC_NULL);
46: PetscRandomCreate(PETSC_COMM_WORLD,&rand);
47: PetscRandomSetFromOptions(rand);
48: PetscOptionsGetInt(PETSC_NULL,"-seed",&seed,PETSC_NULL);
49: PetscRandomSetSeed(rand,seed);
50: PetscRandomSeed(rand);
51: VecSetRandom(v0,rand);
52: EPSSetInitialSpace(eps,1,&v0);
53: /* call the solver */
54: EPSSolve(eps);
56: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
57: Display solution and clean up
58: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
59: EPSPrintSolution(eps,PETSC_NULL);
60: EPSDestroy(&eps);
61: MatDestroy(&A);
62: VecDestroy(&v0);
63: PetscRandomDestroy(&rand);
64: SlepcFinalize();
65: return 0;
66: }