Actual source code: test8.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 NEP view and monitor functionality.\n\n";
13: #include <slepcnep.h>
15: int main(int argc,char **argv)
16: {
17: Mat A[3];
18: FN f[3];
19: NEP nep;
20: Vec xr,xi;
21: PetscScalar kr,ki,coeffs[3];
22: PetscInt n=6,i,Istart,Iend,nconv,its;
23: PetscErrorCode ierr;
25: SlepcInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr;
26: PetscPrintf(PETSC_COMM_WORLD,"\nDiagonal Nonlinear Eigenproblem, n=%D\n\n",n);
28: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
29: Generate the matrices
30: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
32: MatCreate(PETSC_COMM_WORLD,&A[0]);
33: MatSetSizes(A[0],PETSC_DECIDE,PETSC_DECIDE,n,n);
34: MatSetFromOptions(A[0]);
35: MatSetUp(A[0]);
36: MatGetOwnershipRange(A[0],&Istart,&Iend);
37: for (i=Istart;i<Iend;i++) {
38: MatSetValue(A[0],i,i,i+1,INSERT_VALUES);
39: }
40: MatAssemblyBegin(A[0],MAT_FINAL_ASSEMBLY);
41: MatAssemblyEnd(A[0],MAT_FINAL_ASSEMBLY);
43: MatCreate(PETSC_COMM_WORLD,&A[1]);
44: MatSetSizes(A[1],PETSC_DECIDE,PETSC_DECIDE,n,n);
45: MatSetFromOptions(A[1]);
46: MatSetUp(A[1]);
47: for (i=Istart;i<Iend;i++) {
48: MatSetValue(A[1],i,i,-1.5,INSERT_VALUES);
49: }
50: MatAssemblyBegin(A[1],MAT_FINAL_ASSEMBLY);
51: MatAssemblyEnd(A[1],MAT_FINAL_ASSEMBLY);
53: MatCreate(PETSC_COMM_WORLD,&A[2]);
54: MatSetSizes(A[2],PETSC_DECIDE,PETSC_DECIDE,n,n);
55: MatSetFromOptions(A[2]);
56: MatSetUp(A[2]);
57: for (i=Istart;i<Iend;i++) {
58: MatSetValue(A[2],i,i,-1.0/(i+1),INSERT_VALUES);
59: }
60: MatAssemblyBegin(A[2],MAT_FINAL_ASSEMBLY);
61: MatAssemblyEnd(A[2],MAT_FINAL_ASSEMBLY);
63: /*
64: Functions: f0=1.0, f1=lambda, f2=lambda^2
65: */
66: FNCreate(PETSC_COMM_WORLD,&f[0]);
67: FNSetType(f[0],FNRATIONAL);
68: coeffs[0] = 1.0;
69: FNRationalSetNumerator(f[0],1,coeffs);
71: FNCreate(PETSC_COMM_WORLD,&f[1]);
72: FNSetType(f[1],FNRATIONAL);
73: coeffs[0] = 1.0; coeffs[1] = 0.0;
74: FNRationalSetNumerator(f[1],2,coeffs);
76: FNCreate(PETSC_COMM_WORLD,&f[2]);
77: FNSetType(f[2],FNRATIONAL);
78: coeffs[0] = 1.0; coeffs[1] = 0.0; coeffs[2] = 0.0;
79: FNRationalSetNumerator(f[2],3,coeffs);
81: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
82: Create the NEP solver
83: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
84: NEPCreate(PETSC_COMM_WORLD,&nep);
85: PetscObjectSetName((PetscObject)nep,"nep");
86: NEPSetSplitOperator(nep,3,A,f,SAME_NONZERO_PATTERN);
87: NEPSetTarget(nep,1.1);
88: NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE);
89: NEPSetFromOptions(nep);
91: /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
92: Solve the eigensystem and display solution
93: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
94: NEPSolve(nep);
95: NEPGetConverged(nep,&nconv);
96: NEPGetIterationNumber(nep,&its);
97: PetscPrintf(PETSC_COMM_WORLD," %D converged eigenpairs after %D iterations\n",nconv,its);
98: if (nconv>0) {
99: MatCreateVecs(A[0],&xr,&xi);
100: NEPGetEigenpair(nep,0,&kr,&ki,xr,xi);
101: VecDestroy(&xr);
102: VecDestroy(&xi);
103: }
105: NEPErrorView(nep,NEP_ERROR_BACKWARD,NULL);
106: NEPDestroy(&nep);
107: MatDestroy(&A[0]);
108: MatDestroy(&A[1]);
109: MatDestroy(&A[2]);
110: FNDestroy(&f[0]);
111: FNDestroy(&f[1]);
112: FNDestroy(&f[2]);
113: SlepcFinalize();
114: return ierr;
115: }
117: /*TEST
119: test:
120: suffix: 1
121: args: -nep_type slp -nep_target -.5 -nep_error_backward ::ascii_info_detail -nep_view_values -nep_error_absolute ::ascii_matlab -nep_monitor_all -nep_converged_reason -nep_view
122: filter: grep -v "tolerance" | grep -v "problem type" | sed -e "s/[0-9]\.[0-9]*e[+-]\([0-9]*\)/removed/g"
123: requires: double !complex !define(PETSC_USE_64BIT_INDICES)
125: TEST*/