Actual source code: test3.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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 PEP interface functions.\n\n";
13: #include <slepcpep.h>
15: int main(int argc,char **argv)
16: {
17: Mat A[3],B; /* problem matrices */
18: PEP pep; /* eigenproblem solver context */
19: ST st;
20: KSP ksp;
21: DS ds;
22: PetscReal tol,alpha;
23: PetscScalar target;
24: PetscInt n=20,i,its,nev,ncv,mpd,Istart,Iend,nmat;
25: PEPWhich which;
26: PEPConvergedReason reason;
27: PEPType type;
28: PEPExtract extr;
29: PEPBasis basis;
30: PEPScale scale;
31: PEPRefine refine;
32: PEPRefineScheme rscheme;
33: PEPConv conv;
34: PEPStop stop;
35: PEPProblemType ptype;
36: PetscErrorCode ierr;
37: PetscViewerAndFormat *vf;
39: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
40: PetscPrintf(PETSC_COMM_WORLD,"\nDiagonal Quadratic Eigenproblem, n=%D\n\n",n);
42: MatCreate(PETSC_COMM_WORLD,&A[0]);
43: MatSetSizes(A[0],PETSC_DECIDE,PETSC_DECIDE,n,n);
44: MatSetFromOptions(A[0]);
45: MatSetUp(A[0]);
46: MatGetOwnershipRange(A[0],&Istart,&Iend);
47: for (i=Istart;i<Iend;i++) {
48: MatSetValue(A[0],i,i,i+1,INSERT_VALUES);
49: }
50: MatAssemblyBegin(A[0],MAT_FINAL_ASSEMBLY);
51: MatAssemblyEnd(A[0],MAT_FINAL_ASSEMBLY);
53: MatCreate(PETSC_COMM_WORLD,&A[1]);
54: MatSetSizes(A[1],PETSC_DECIDE,PETSC_DECIDE,n,n);
55: MatSetFromOptions(A[1]);
56: MatSetUp(A[1]);
57: MatGetOwnershipRange(A[1],&Istart,&Iend);
58: for (i=Istart;i<Iend;i++) {
59: MatSetValue(A[1],i,i,1.0,INSERT_VALUES);
60: }
61: MatAssemblyBegin(A[1],MAT_FINAL_ASSEMBLY);
62: MatAssemblyEnd(A[1],MAT_FINAL_ASSEMBLY);
64: MatCreate(PETSC_COMM_WORLD,&A[2]);
65: MatSetSizes(A[2],PETSC_DECIDE,PETSC_DECIDE,n,n);
66: MatSetFromOptions(A[2]);
67: MatSetUp(A[2]);
68: MatGetOwnershipRange(A[1],&Istart,&Iend);
69: for (i=Istart;i<Iend;i++) {
70: MatSetValue(A[2],i,i,n/(PetscReal)(i+1),INSERT_VALUES);
71: }
72: MatAssemblyBegin(A[2],MAT_FINAL_ASSEMBLY);
73: MatAssemblyEnd(A[2],MAT_FINAL_ASSEMBLY);
75: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
76: Create eigensolver and test interface functions
77: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
78: PEPCreate(PETSC_COMM_WORLD,&pep);
79: PEPSetOperators(pep,3,A);
80: PEPGetNumMatrices(pep,&nmat);
81: PetscPrintf(PETSC_COMM_WORLD," Polynomial of degree %d\n",(int)nmat-1);
82: PEPGetOperators(pep,0,&B);
83: MatView(B,NULL);
85: PEPSetType(pep,PEPTOAR);
86: PEPGetType(pep,&type);
87: PetscPrintf(PETSC_COMM_WORLD," Type set to %s\n",type);
89: PEPGetProblemType(pep,&ptype);
90: PetscPrintf(PETSC_COMM_WORLD," Problem type before changing = %d",(int)ptype);
91: PEPSetProblemType(pep,PEP_HERMITIAN);
92: PEPGetProblemType(pep,&ptype);
93: PetscPrintf(PETSC_COMM_WORLD," ... changed to %d.\n",(int)ptype);
95: PEPGetExtract(pep,&extr);
96: PetscPrintf(PETSC_COMM_WORLD," Extraction before changing = %d",(int)extr);
97: PEPSetExtract(pep,PEP_EXTRACT_STRUCTURED);
98: PEPGetExtract(pep,&extr);
99: PetscPrintf(PETSC_COMM_WORLD," ... changed to %d\n",(int)extr);
101: PEPSetScale(pep,PEP_SCALE_SCALAR,.1,NULL,NULL,5,1.0);
102: PEPGetScale(pep,&scale,&alpha,NULL,NULL,&its,NULL);
103: PetscPrintf(PETSC_COMM_WORLD," Scaling: %s, alpha=%g, its=%D\n",PEPScaleTypes[scale],(double)alpha,its);
105: PEPSetBasis(pep,PEP_BASIS_CHEBYSHEV1);
106: PEPGetBasis(pep,&basis);
107: PetscPrintf(PETSC_COMM_WORLD," Polynomial basis: %s\n",PEPBasisTypes[basis]);
109: PEPSetRefine(pep,PEP_REFINE_SIMPLE,1,1e-9,2,PEP_REFINE_SCHEME_SCHUR);
110: PEPGetRefine(pep,&refine,NULL,&tol,&its,&rscheme);
111: PetscPrintf(PETSC_COMM_WORLD," Refinement: %s, tol=%g, its=%D, scheme=%s\n",PEPRefineTypes[refine],(double)tol,its,PEPRefineSchemes[rscheme]);
113: PEPSetTarget(pep,4.8);
114: PEPGetTarget(pep,&target);
115: PEPSetWhichEigenpairs(pep,PEP_TARGET_MAGNITUDE);
116: PEPGetWhichEigenpairs(pep,&which);
117: PetscPrintf(PETSC_COMM_WORLD," Which = %d, target = %g\n",(int)which,(double)PetscRealPart(target));
119: PEPSetDimensions(pep,4,PETSC_DEFAULT,PETSC_DEFAULT);
120: PEPGetDimensions(pep,&nev,&ncv,&mpd);
121: PetscPrintf(PETSC_COMM_WORLD," Dimensions: nev=%D, ncv=%D, mpd=%D\n",nev,ncv,mpd);
123: PEPSetTolerances(pep,2.2e-4,200);
124: PEPGetTolerances(pep,&tol,&its);
125: PetscPrintf(PETSC_COMM_WORLD," Tolerance = %.5f, max_its = %D\n",(double)tol,its);
127: PEPSetConvergenceTest(pep,PEP_CONV_ABS);
128: PEPGetConvergenceTest(pep,&conv);
129: PEPSetStoppingTest(pep,PEP_STOP_BASIC);
130: PEPGetStoppingTest(pep,&stop);
131: PetscPrintf(PETSC_COMM_WORLD," Convergence test = %d, stopping test = %d\n",(int)conv,(int)stop);
133: PetscViewerAndFormatCreate(PETSC_VIEWER_STDOUT_WORLD,PETSC_VIEWER_DEFAULT,&vf);
134: PEPMonitorSet(pep,(PetscErrorCode (*)(PEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*))PEPMonitorFirst,vf,(PetscErrorCode (*)(void**))PetscViewerAndFormatDestroy);
135: PEPMonitorCancel(pep);
137: PEPGetST(pep,&st);
138: STGetKSP(st,&ksp);
139: KSPSetTolerances(ksp,1e-8,1e-35,PETSC_DEFAULT,PETSC_DEFAULT);
140: STView(st,NULL);
141: PEPGetDS(pep,&ds);
142: DSView(ds,NULL);
144: PEPSetFromOptions(pep);
145: PEPSolve(pep);
146: PEPGetConvergedReason(pep,&reason);
147: PetscPrintf(PETSC_COMM_WORLD," Finished - converged reason = %d\n",(int)reason);
149: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
150: Display solution and clean up
151: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
152: PEPErrorView(pep,PEP_ERROR_RELATIVE,NULL);
153: PEPDestroy(&pep);
154: MatDestroy(&A[0]);
155: MatDestroy(&A[1]);
156: MatDestroy(&A[2]);
157: SlepcFinalize();
158: return ierr;
159: }
161: /*TEST
163: test:
164: suffix: 1
165: args: -pep_tol 1e-6 -pep_ncv 22
167: TEST*/