Actual source code: nepdefault.c
slepc-3.18.3 2023-03-24
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: */
10: /*
11: Simple default routines for common NEP operations
12: */
14: #include <slepc/private/nepimpl.h>
16: /*@
17: NEPSetWorkVecs - Sets a number of work vectors into a NEP object
19: Collective on nep
21: Input Parameters:
22: + nep - nonlinear eigensolver context
23: - nw - number of work vectors to allocate
25: Developer Notes:
26: This is SLEPC_EXTERN because it may be required by user plugin NEP
27: implementations.
29: Level: developer
31: .seealso: NEPSetUp()
32: @*/
33: PetscErrorCode NEPSetWorkVecs(NEP nep,PetscInt nw)
34: {
35: Vec t;
40: if (nep->nwork < nw) {
41: VecDestroyVecs(nep->nwork,&nep->work);
42: nep->nwork = nw;
43: BVGetColumn(nep->V,0,&t);
44: VecDuplicateVecs(t,nw,&nep->work);
45: BVRestoreColumn(nep->V,0,&t);
46: }
47: return 0;
48: }
50: /*
51: NEPGetDefaultShift - Return the value of sigma to start the nonlinear iteration.
52: */
53: PetscErrorCode NEPGetDefaultShift(NEP nep,PetscScalar *sigma)
54: {
56: switch (nep->which) {
57: case NEP_LARGEST_MAGNITUDE:
58: case NEP_LARGEST_IMAGINARY:
59: case NEP_ALL:
60: case NEP_WHICH_USER:
61: *sigma = 1.0; /* arbitrary value */
62: break;
63: case NEP_SMALLEST_MAGNITUDE:
64: case NEP_SMALLEST_IMAGINARY:
65: *sigma = 0.0;
66: break;
67: case NEP_LARGEST_REAL:
68: *sigma = PETSC_MAX_REAL;
69: break;
70: case NEP_SMALLEST_REAL:
71: *sigma = PETSC_MIN_REAL;
72: break;
73: case NEP_TARGET_MAGNITUDE:
74: case NEP_TARGET_REAL:
75: case NEP_TARGET_IMAGINARY:
76: *sigma = nep->target;
77: break;
78: }
79: return 0;
80: }
82: /*
83: NEPConvergedRelative - Checks convergence relative to the eigenvalue.
84: */
85: PetscErrorCode NEPConvergedRelative(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
86: {
87: PetscReal w;
89: w = SlepcAbsEigenvalue(eigr,eigi);
90: *errest = res/w;
91: return 0;
92: }
94: /*
95: NEPConvergedAbsolute - Checks convergence absolutely.
96: */
97: PetscErrorCode NEPConvergedAbsolute(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
98: {
99: *errest = res;
100: return 0;
101: }
103: /*
104: NEPConvergedNorm - Checks convergence relative to the matrix norms.
105: */
106: PetscErrorCode NEPConvergedNorm(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
107: {
108: PetscScalar s;
109: PetscReal w=0.0;
110: PetscInt j;
111: PetscBool flg;
113: if (nep->fui!=NEP_USER_INTERFACE_SPLIT) {
114: NEPComputeFunction(nep,eigr,nep->function,nep->function);
115: MatHasOperation(nep->function,MATOP_NORM,&flg);
117: MatNorm(nep->function,NORM_INFINITY,&w);
118: } else {
119: /* initialization of matrix norms */
120: if (!nep->nrma[0]) {
121: for (j=0;j<nep->nt;j++) {
122: MatHasOperation(nep->A[j],MATOP_NORM,&flg);
124: MatNorm(nep->A[j],NORM_INFINITY,&nep->nrma[j]);
125: }
126: }
127: for (j=0;j<nep->nt;j++) {
128: FNEvaluateFunction(nep->f[j],eigr,&s);
129: w = w + nep->nrma[j]*PetscAbsScalar(s);
130: }
131: }
132: *errest = res/w;
133: return 0;
134: }
136: /*@C
137: NEPStoppingBasic - Default routine to determine whether the outer eigensolver
138: iteration must be stopped.
140: Collective on nep
142: Input Parameters:
143: + nep - nonlinear eigensolver context obtained from NEPCreate()
144: . its - current number of iterations
145: . max_it - maximum number of iterations
146: . nconv - number of currently converged eigenpairs
147: . nev - number of requested eigenpairs
148: - ctx - context (not used here)
150: Output Parameter:
151: . reason - result of the stopping test
153: Notes:
154: A positive value of reason indicates that the iteration has finished successfully
155: (converged), and a negative value indicates an error condition (diverged). If
156: the iteration needs to be continued, reason must be set to NEP_CONVERGED_ITERATING
157: (zero).
159: NEPStoppingBasic() will stop if all requested eigenvalues are converged, or if
160: the maximum number of iterations has been reached.
162: Use NEPSetStoppingTest() to provide your own test instead of using this one.
164: Level: advanced
166: .seealso: NEPSetStoppingTest(), NEPConvergedReason, NEPGetConvergedReason()
167: @*/
168: PetscErrorCode NEPStoppingBasic(NEP nep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,NEPConvergedReason *reason,void *ctx)
169: {
170: *reason = NEP_CONVERGED_ITERATING;
171: if (nconv >= nev) {
172: PetscInfo(nep,"Nonlinear eigensolver finished successfully: %" PetscInt_FMT " eigenpairs converged at iteration %" PetscInt_FMT "\n",nconv,its);
173: *reason = NEP_CONVERGED_TOL;
174: } else if (its >= max_it) {
175: *reason = NEP_DIVERGED_ITS;
176: PetscInfo(nep,"Nonlinear eigensolver iteration reached maximum number of iterations (%" PetscInt_FMT ")\n",its);
177: }
178: return 0;
179: }
181: PetscErrorCode NEPComputeVectors_Schur(NEP nep)
182: {
183: Mat Z;
185: DSVectors(nep->ds,DS_MAT_X,NULL,NULL);
186: DSGetMat(nep->ds,DS_MAT_X,&Z);
187: BVMultInPlace(nep->V,Z,0,nep->nconv);
188: DSRestoreMat(nep->ds,DS_MAT_X,&Z);
189: BVNormalize(nep->V,nep->eigi);
190: return 0;
191: }