Actual source code: test30.c
slepc-3.17.2 2022-08-09
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, 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[] = "Test changing EPS type.\n\n";
13: #include <slepceps.h>
15: int main(int argc,char **argv)
16: {
17: Mat A; /* problem matrix */
18: EPS eps; /* eigenproblem solver context */
19: ST st;
20: KSP ksp;
21: PetscInt n=20,i,Istart,Iend,nrest;
23: SlepcInitialize(&argc,&argv,(char*)0,help);
24: PetscPrintf(PETSC_COMM_WORLD,"\nDiagonal Eigenproblem, n=%" PetscInt_FMT "\n\n",n);
26: MatCreate(PETSC_COMM_WORLD,&A);
27: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
28: MatSetFromOptions(A);
29: MatSetUp(A);
30: MatGetOwnershipRange(A,&Istart,&Iend);
31: for (i=Istart;i<Iend;i++) MatSetValue(A,i,i,i+1,INSERT_VALUES);
32: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
33: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
35: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
36: Create eigensolver and view options
37: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
38: EPSCreate(PETSC_COMM_WORLD,&eps);
39: EPSSetOperators(eps,A,NULL);
40: EPSSetProblemType(eps,EPS_HEP);
41: EPSSetTarget(eps,4.8);
42: EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);
43: EPSSetDimensions(eps,4,PETSC_DEFAULT,PETSC_DEFAULT);
45: EPSSetType(eps,EPSRQCG);
46: EPSRQCGSetReset(eps,10);
47: EPSSetFromOptions(eps);
48: EPSRQCGGetReset(eps,&nrest);
49: PetscPrintf(PETSC_COMM_WORLD,"RQCG reset parameter = %" PetscInt_FMT "\n\n",nrest);
50: EPSView(eps,NULL);
52: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
53: Change eigensolver type and solve problem
54: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
55: EPSSetType(eps,EPSGD);
56: EPSGetST(eps,&st);
57: STGetKSP(st,&ksp);
58: KSPSetType(ksp,KSPPREONLY);
59: EPSSolve(eps);
60: EPSView(eps,NULL);
62: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
63: Display solution and clean up
64: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
65: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
66: EPSDestroy(&eps);
67: MatDestroy(&A);
68: SlepcFinalize();
69: return 0;
70: }
72: /*TEST
74: test:
75: suffix: 1
76: args: -eps_harmonic -eps_conv_norm
77: filter: grep -v tolerance | sed -e "s/symmetric/hermitian/"
79: TEST*/