Actual source code: test30.c
slepc-3.18.1 2022-11-02
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;
24: SlepcInitialize(&argc,&argv,(char*)0,help);
25: PetscPrintf(PETSC_COMM_WORLD,"\nDiagonal Eigenproblem, n=%" PetscInt_FMT "\n\n",n);
27: MatCreate(PETSC_COMM_WORLD,&A);
28: MatSetSizes(A,PETSC_DECIDE,PETSC_DECIDE,n,n);
29: MatSetFromOptions(A);
30: MatSetUp(A);
31: MatGetOwnershipRange(A,&Istart,&Iend);
32: for (i=Istart;i<Iend;i++) MatSetValue(A,i,i,i+1,INSERT_VALUES);
33: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
34: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
36: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
37: Create eigensolver and view options
38: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
39: EPSCreate(PETSC_COMM_WORLD,&eps);
40: EPSSetOperators(eps,A,NULL);
41: EPSSetProblemType(eps,EPS_HEP);
42: EPSSetTarget(eps,4.8);
43: EPSSetWhichEigenpairs(eps,EPS_TARGET_MAGNITUDE);
44: EPSSetDimensions(eps,4,PETSC_DEFAULT,PETSC_DEFAULT);
46: EPSSetType(eps,EPSRQCG);
47: EPSRQCGSetReset(eps,10);
48: EPSSetFromOptions(eps);
49: EPSRQCGGetReset(eps,&nrest);
50: PetscPrintf(PETSC_COMM_WORLD,"RQCG reset parameter = %" PetscInt_FMT "\n\n",nrest);
51: EPSView(eps,NULL);
53: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
54: Change eigensolver type and solve problem
55: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
56: EPSSetType(eps,EPSGD);
57: EPSGetST(eps,&st);
58: STGetKSP(st,&ksp);
59: KSPSetType(ksp,KSPPREONLY);
60: EPSSolve(eps);
61: EPSView(eps,NULL);
63: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
64: Display solution and clean up
65: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
66: EPSErrorView(eps,EPS_ERROR_RELATIVE,NULL);
67: EPSDestroy(&eps);
68: MatDestroy(&A);
69: SlepcFinalize();
70: return 0;
71: }
73: /*TEST
75: test:
76: suffix: 1
77: args: -eps_harmonic -eps_conv_norm
78: filter: grep -v tolerance | sed -e "s/symmetric/hermitian/"
80: TEST*/