Actual source code: rii.c
slepc-3.18.0 2022-10-01
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: SLEPc nonlinear eigensolver: "rii"
13: Method: Residual inverse iteration
15: Algorithm:
17: Simple residual inverse iteration with varying shift.
19: References:
21: [1] A. Neumaier, "Residual inverse iteration for the nonlinear
22: eigenvalue problem", SIAM J. Numer. Anal. 22(5):914-923, 1985.
23: */
25: #include <slepc/private/nepimpl.h>
26: #include <../src/nep/impls/nepdefl.h>
28: typedef struct {
29: PetscInt max_inner_it; /* maximum number of Newton iterations */
30: PetscInt lag; /* interval to rebuild preconditioner */
31: PetscBool cctol; /* constant correction tolerance */
32: PetscBool herm; /* whether the Hermitian version of the scalar equation must be used */
33: PetscReal deftol; /* tolerance for the deflation (threshold) */
34: KSP ksp; /* linear solver object */
35: } NEP_RII;
37: PetscErrorCode NEPSetUp_RII(NEP nep)
38: {
39: if (nep->ncv!=PETSC_DEFAULT) PetscInfo(nep,"Setting ncv = nev, ignoring user-provided value\n");
40: nep->ncv = nep->nev;
41: if (nep->mpd!=PETSC_DEFAULT) PetscInfo(nep,"Setting mpd = nev, ignoring user-provided value\n");
42: nep->mpd = nep->nev;
43: if (nep->max_it==PETSC_DEFAULT) nep->max_it = PetscMax(5000,2*nep->n/nep->ncv);
44: if (!nep->which) nep->which = NEP_TARGET_MAGNITUDE;
46: NEPCheckUnsupported(nep,NEP_FEATURE_REGION | NEP_FEATURE_TWOSIDED);
47: NEPAllocateSolution(nep,0);
48: NEPSetWorkVecs(nep,2);
49: return 0;
50: }
52: PetscErrorCode NEPSolve_RII(NEP nep)
53: {
54: NEP_RII *ctx = (NEP_RII*)nep->data;
55: Mat T,Tp,H,A;
56: Vec uu,u,r,delta,t;
57: PetscScalar lambda,lambda2,sigma,a1,a2,corr;
58: PetscReal nrm,resnorm=1.0,ktol=0.1,perr,rtol;
59: PetscBool skip=PETSC_FALSE,lock=PETSC_FALSE;
60: PetscInt inner_its,its=0;
61: NEP_EXT_OP extop=NULL;
62: KSPConvergedReason kspreason;
64: /* get initial approximation of eigenvalue and eigenvector */
65: NEPGetDefaultShift(nep,&sigma);
66: lambda = sigma;
67: if (!nep->nini) {
68: BVSetRandomColumn(nep->V,0);
69: BVNormColumn(nep->V,0,NORM_2,&nrm);
70: BVScaleColumn(nep->V,0,1.0/nrm);
71: }
72: if (!ctx->ksp) NEPRIIGetKSP(nep,&ctx->ksp);
73: NEPDeflationInitialize(nep,nep->V,ctx->ksp,PETSC_FALSE,nep->nev,&extop);
74: NEPDeflationCreateVec(extop,&u);
75: VecDuplicate(u,&r);
76: VecDuplicate(u,&delta);
77: BVGetColumn(nep->V,0,&uu);
78: NEPDeflationCopyToExtendedVec(extop,uu,NULL,u,PETSC_FALSE);
79: BVRestoreColumn(nep->V,0,&uu);
81: /* prepare linear solver */
82: NEPDeflationSolveSetUp(extop,sigma);
83: KSPGetTolerances(ctx->ksp,&rtol,NULL,NULL,NULL);
85: VecCopy(u,r);
86: NEPDeflationFunctionSolve(extop,r,u);
87: VecNorm(u,NORM_2,&nrm);
88: VecScale(u,1.0/nrm);
90: /* Restart loop */
91: while (nep->reason == NEP_CONVERGED_ITERATING) {
92: its++;
94: /* Use Newton's method to compute nonlinear Rayleigh functional. Current eigenvalue
95: estimate as starting value. */
96: inner_its=0;
97: lambda2 = lambda;
98: do {
99: NEPDeflationComputeFunction(extop,lambda,&T);
100: MatMult(T,u,r);
101: if (!ctx->herm) {
102: NEPDeflationFunctionSolve(extop,r,delta);
103: KSPGetConvergedReason(ctx->ksp,&kspreason);
104: if (kspreason<0) PetscInfo(nep,"iter=%" PetscInt_FMT ", linear solve failed\n",nep->its);
105: t = delta;
106: } else t = r;
107: VecDot(t,u,&a1);
108: NEPDeflationComputeJacobian(extop,lambda,&Tp);
109: MatMult(Tp,u,r);
110: if (!ctx->herm) {
111: NEPDeflationFunctionSolve(extop,r,delta);
112: KSPGetConvergedReason(ctx->ksp,&kspreason);
113: if (kspreason<0) PetscInfo(nep,"iter=%" PetscInt_FMT ", linear solve failed\n",nep->its);
114: t = delta;
115: } else t = r;
116: VecDot(t,u,&a2);
117: corr = a1/a2;
118: lambda = lambda - corr;
119: inner_its++;
120: } while (PetscAbsScalar(corr)/PetscAbsScalar(lambda)>PETSC_SQRT_MACHINE_EPSILON && inner_its<ctx->max_inner_it);
122: /* form residual, r = T(lambda)*u */
123: NEPDeflationComputeFunction(extop,lambda,&T);
124: MatMult(T,u,r);
126: /* convergence test */
127: perr = nep->errest[nep->nconv];
128: VecNorm(r,NORM_2,&resnorm);
129: (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
130: nep->eigr[nep->nconv] = lambda;
131: if (its>1 && (nep->errest[nep->nconv]<=nep->tol || nep->errest[nep->nconv]<=ctx->deftol)) {
132: if (nep->errest[nep->nconv]<=ctx->deftol && !extop->ref && nep->nconv) {
133: NEPDeflationExtractEigenpair(extop,nep->nconv,u,lambda,nep->ds);
134: NEPDeflationSetRefine(extop,PETSC_TRUE);
135: MatMult(T,u,r);
136: VecNorm(r,NORM_2,&resnorm);
137: (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
138: if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
139: } else if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
140: }
141: if (lock) {
142: NEPDeflationSetRefine(extop,PETSC_FALSE);
143: nep->nconv = nep->nconv + 1;
144: NEPDeflationLocking(extop,u,lambda);
145: nep->its += its;
146: skip = PETSC_TRUE;
147: lock = PETSC_FALSE;
148: }
149: (*nep->stopping)(nep,nep->its+its,nep->max_it,nep->nconv,nep->nev,&nep->reason,nep->stoppingctx);
150: if (!skip || nep->reason>0) NEPMonitor(nep,nep->its+its,nep->nconv,nep->eigr,nep->eigi,nep->errest,(nep->reason>0)?nep->nconv:nep->nconv+1);
152: if (nep->reason == NEP_CONVERGED_ITERATING) {
153: if (!skip) {
154: /* update preconditioner and set adaptive tolerance */
155: if (ctx->lag && !(its%ctx->lag) && its>=2*ctx->lag && perr && nep->errest[nep->nconv]>.5*perr) NEPDeflationSolveSetUp(extop,lambda2);
156: if (!ctx->cctol) {
157: ktol = PetscMax(ktol/2.0,rtol);
158: KSPSetTolerances(ctx->ksp,ktol,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
159: }
161: /* eigenvector correction: delta = T(sigma)\r */
162: NEPDeflationFunctionSolve(extop,r,delta);
163: KSPGetConvergedReason(ctx->ksp,&kspreason);
164: if (kspreason<0) {
165: PetscInfo(nep,"iter=%" PetscInt_FMT ", linear solve failed, stopping solve\n",nep->its);
166: nep->reason = NEP_DIVERGED_LINEAR_SOLVE;
167: break;
168: }
170: /* update eigenvector: u = u - delta */
171: VecAXPY(u,-1.0,delta);
173: /* normalize eigenvector */
174: VecNormalize(u,NULL);
175: } else {
176: its = -1;
177: NEPDeflationSetRandomVec(extop,u);
178: NEPDeflationSolveSetUp(extop,sigma);
179: VecCopy(u,r);
180: NEPDeflationFunctionSolve(extop,r,u);
181: VecNorm(u,NORM_2,&nrm);
182: VecScale(u,1.0/nrm);
183: lambda = sigma;
184: skip = PETSC_FALSE;
185: }
186: }
187: }
188: NEPDeflationGetInvariantPair(extop,NULL,&H);
189: DSSetType(nep->ds,DSNHEP);
190: DSAllocate(nep->ds,PetscMax(nep->nconv,1));
191: DSSetDimensions(nep->ds,nep->nconv,0,nep->nconv);
192: DSGetMat(nep->ds,DS_MAT_A,&A);
193: MatCopy(H,A,SAME_NONZERO_PATTERN);
194: DSRestoreMat(nep->ds,DS_MAT_A,&A);
195: MatDestroy(&H);
196: DSSolve(nep->ds,nep->eigr,nep->eigi);
197: NEPDeflationReset(extop);
198: VecDestroy(&u);
199: VecDestroy(&r);
200: VecDestroy(&delta);
201: return 0;
202: }
204: PetscErrorCode NEPSetFromOptions_RII(NEP nep,PetscOptionItems *PetscOptionsObject)
205: {
206: NEP_RII *ctx = (NEP_RII*)nep->data;
207: PetscBool flg;
208: PetscInt i;
209: PetscReal r;
211: PetscOptionsHeadBegin(PetscOptionsObject,"NEP RII Options");
213: i = 0;
214: PetscOptionsInt("-nep_rii_max_it","Maximum number of Newton iterations for updating Rayleigh functional","NEPRIISetMaximumIterations",ctx->max_inner_it,&i,&flg);
215: if (flg) NEPRIISetMaximumIterations(nep,i);
217: PetscOptionsBool("-nep_rii_const_correction_tol","Constant correction tolerance for the linear solver","NEPRIISetConstCorrectionTol",ctx->cctol,&ctx->cctol,NULL);
219: PetscOptionsBool("-nep_rii_hermitian","Use Hermitian version of the scalar nonlinear equation","NEPRIISetHermitian",ctx->herm,&ctx->herm,NULL);
221: i = 0;
222: PetscOptionsInt("-nep_rii_lag_preconditioner","Interval to rebuild preconditioner","NEPRIISetLagPreconditioner",ctx->lag,&i,&flg);
223: if (flg) NEPRIISetLagPreconditioner(nep,i);
225: r = 0.0;
226: PetscOptionsReal("-nep_rii_deflation_threshold","Tolerance used as a threshold for including deflated eigenpairs","NEPRIISetDeflationThreshold",ctx->deftol,&r,&flg);
227: if (flg) NEPRIISetDeflationThreshold(nep,r);
229: PetscOptionsHeadEnd();
231: if (!ctx->ksp) NEPRIIGetKSP(nep,&ctx->ksp);
232: KSPSetFromOptions(ctx->ksp);
233: return 0;
234: }
236: static PetscErrorCode NEPRIISetMaximumIterations_RII(NEP nep,PetscInt its)
237: {
238: NEP_RII *ctx = (NEP_RII*)nep->data;
240: if (its==PETSC_DEFAULT) ctx->max_inner_it = 10;
241: else {
243: ctx->max_inner_it = its;
244: }
245: return 0;
246: }
248: /*@
249: NEPRIISetMaximumIterations - Sets the maximum number of inner iterations to be
250: used in the RII solver. These are the Newton iterations related to the computation
251: of the nonlinear Rayleigh functional.
253: Logically Collective on nep
255: Input Parameters:
256: + nep - nonlinear eigenvalue solver
257: - its - maximum inner iterations
259: Level: advanced
261: .seealso: NEPRIIGetMaximumIterations()
262: @*/
263: PetscErrorCode NEPRIISetMaximumIterations(NEP nep,PetscInt its)
264: {
267: PetscTryMethod(nep,"NEPRIISetMaximumIterations_C",(NEP,PetscInt),(nep,its));
268: return 0;
269: }
271: static PetscErrorCode NEPRIIGetMaximumIterations_RII(NEP nep,PetscInt *its)
272: {
273: NEP_RII *ctx = (NEP_RII*)nep->data;
275: *its = ctx->max_inner_it;
276: return 0;
277: }
279: /*@
280: NEPRIIGetMaximumIterations - Gets the maximum number of inner iterations of RII.
282: Not Collective
284: Input Parameter:
285: . nep - nonlinear eigenvalue solver
287: Output Parameter:
288: . its - maximum inner iterations
290: Level: advanced
292: .seealso: NEPRIISetMaximumIterations()
293: @*/
294: PetscErrorCode NEPRIIGetMaximumIterations(NEP nep,PetscInt *its)
295: {
298: PetscUseMethod(nep,"NEPRIIGetMaximumIterations_C",(NEP,PetscInt*),(nep,its));
299: return 0;
300: }
302: static PetscErrorCode NEPRIISetLagPreconditioner_RII(NEP nep,PetscInt lag)
303: {
304: NEP_RII *ctx = (NEP_RII*)nep->data;
307: ctx->lag = lag;
308: return 0;
309: }
311: /*@
312: NEPRIISetLagPreconditioner - Determines when the preconditioner is rebuilt in the
313: nonlinear solve.
315: Logically Collective on nep
317: Input Parameters:
318: + nep - nonlinear eigenvalue solver
319: - lag - 0 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is
320: computed within the nonlinear iteration, 2 means every second time
321: the Jacobian is built, etc.
323: Options Database Keys:
324: . -nep_rii_lag_preconditioner <lag> - the lag value
326: Notes:
327: The default is 1.
328: The preconditioner is ALWAYS built in the first iteration of a nonlinear solve.
330: Level: intermediate
332: .seealso: NEPRIIGetLagPreconditioner()
333: @*/
334: PetscErrorCode NEPRIISetLagPreconditioner(NEP nep,PetscInt lag)
335: {
338: PetscTryMethod(nep,"NEPRIISetLagPreconditioner_C",(NEP,PetscInt),(nep,lag));
339: return 0;
340: }
342: static PetscErrorCode NEPRIIGetLagPreconditioner_RII(NEP nep,PetscInt *lag)
343: {
344: NEP_RII *ctx = (NEP_RII*)nep->data;
346: *lag = ctx->lag;
347: return 0;
348: }
350: /*@
351: NEPRIIGetLagPreconditioner - Indicates how often the preconditioner is rebuilt.
353: Not Collective
355: Input Parameter:
356: . nep - nonlinear eigenvalue solver
358: Output Parameter:
359: . lag - the lag parameter
361: Level: intermediate
363: .seealso: NEPRIISetLagPreconditioner()
364: @*/
365: PetscErrorCode NEPRIIGetLagPreconditioner(NEP nep,PetscInt *lag)
366: {
369: PetscUseMethod(nep,"NEPRIIGetLagPreconditioner_C",(NEP,PetscInt*),(nep,lag));
370: return 0;
371: }
373: static PetscErrorCode NEPRIISetConstCorrectionTol_RII(NEP nep,PetscBool cct)
374: {
375: NEP_RII *ctx = (NEP_RII*)nep->data;
377: ctx->cctol = cct;
378: return 0;
379: }
381: /*@
382: NEPRIISetConstCorrectionTol - Sets a flag to keep the tolerance used
383: in the linear solver constant.
385: Logically Collective on nep
387: Input Parameters:
388: + nep - nonlinear eigenvalue solver
389: - cct - a boolean value
391: Options Database Keys:
392: . -nep_rii_const_correction_tol <bool> - set the boolean flag
394: Notes:
395: By default, an exponentially decreasing tolerance is set in the KSP used
396: within the nonlinear iteration, so that each Newton iteration requests
397: better accuracy than the previous one. The constant correction tolerance
398: flag stops this behaviour.
400: Level: intermediate
402: .seealso: NEPRIIGetConstCorrectionTol()
403: @*/
404: PetscErrorCode NEPRIISetConstCorrectionTol(NEP nep,PetscBool cct)
405: {
408: PetscTryMethod(nep,"NEPRIISetConstCorrectionTol_C",(NEP,PetscBool),(nep,cct));
409: return 0;
410: }
412: static PetscErrorCode NEPRIIGetConstCorrectionTol_RII(NEP nep,PetscBool *cct)
413: {
414: NEP_RII *ctx = (NEP_RII*)nep->data;
416: *cct = ctx->cctol;
417: return 0;
418: }
420: /*@
421: NEPRIIGetConstCorrectionTol - Returns the constant tolerance flag.
423: Not Collective
425: Input Parameter:
426: . nep - nonlinear eigenvalue solver
428: Output Parameter:
429: . cct - the value of the constant tolerance flag
431: Level: intermediate
433: .seealso: NEPRIISetConstCorrectionTol()
434: @*/
435: PetscErrorCode NEPRIIGetConstCorrectionTol(NEP nep,PetscBool *cct)
436: {
439: PetscUseMethod(nep,"NEPRIIGetConstCorrectionTol_C",(NEP,PetscBool*),(nep,cct));
440: return 0;
441: }
443: static PetscErrorCode NEPRIISetHermitian_RII(NEP nep,PetscBool herm)
444: {
445: NEP_RII *ctx = (NEP_RII*)nep->data;
447: ctx->herm = herm;
448: return 0;
449: }
451: /*@
452: NEPRIISetHermitian - Sets a flag to indicate if the Hermitian version of the
453: scalar nonlinear equation must be used by the solver.
455: Logically Collective on nep
457: Input Parameters:
458: + nep - nonlinear eigenvalue solver
459: - herm - a boolean value
461: Options Database Keys:
462: . -nep_rii_hermitian <bool> - set the boolean flag
464: Notes:
465: By default, the scalar nonlinear equation x'*inv(T(sigma))*T(z)*x=0 is solved
466: at each step of the nonlinear iteration. When this flag is set the simpler
467: form x'*T(z)*x=0 is used, which is supposed to be valid only for Hermitian
468: problems.
470: Level: intermediate
472: .seealso: NEPRIIGetHermitian()
473: @*/
474: PetscErrorCode NEPRIISetHermitian(NEP nep,PetscBool herm)
475: {
478: PetscTryMethod(nep,"NEPRIISetHermitian_C",(NEP,PetscBool),(nep,herm));
479: return 0;
480: }
482: static PetscErrorCode NEPRIIGetHermitian_RII(NEP nep,PetscBool *herm)
483: {
484: NEP_RII *ctx = (NEP_RII*)nep->data;
486: *herm = ctx->herm;
487: return 0;
488: }
490: /*@
491: NEPRIIGetHermitian - Returns the flag about using the Hermitian version of
492: the scalar nonlinear equation.
494: Not Collective
496: Input Parameter:
497: . nep - nonlinear eigenvalue solver
499: Output Parameter:
500: . herm - the value of the hermitian flag
502: Level: intermediate
504: .seealso: NEPRIISetHermitian()
505: @*/
506: PetscErrorCode NEPRIIGetHermitian(NEP nep,PetscBool *herm)
507: {
510: PetscUseMethod(nep,"NEPRIIGetHermitian_C",(NEP,PetscBool*),(nep,herm));
511: return 0;
512: }
514: static PetscErrorCode NEPRIISetDeflationThreshold_RII(NEP nep,PetscReal deftol)
515: {
516: NEP_RII *ctx = (NEP_RII*)nep->data;
518: ctx->deftol = deftol;
519: return 0;
520: }
522: /*@
523: NEPRIISetDeflationThreshold - Sets the threshold value used to switch between
524: deflated and non-deflated iteration.
526: Logically Collective on nep
528: Input Parameters:
529: + nep - nonlinear eigenvalue solver
530: - deftol - the threshold value
532: Options Database Keys:
533: . -nep_rii_deflation_threshold <deftol> - set the threshold
535: Notes:
536: Normally, the solver iterates on the extended problem in order to deflate
537: previously converged eigenpairs. If this threshold is set to a nonzero value,
538: then once the residual error is below this threshold the solver will
539: continue the iteration without deflation. The intention is to be able to
540: improve the current eigenpair further, despite having previous eigenpairs
541: with somewhat bad precision.
543: Level: advanced
545: .seealso: NEPRIIGetDeflationThreshold()
546: @*/
547: PetscErrorCode NEPRIISetDeflationThreshold(NEP nep,PetscReal deftol)
548: {
551: PetscTryMethod(nep,"NEPRIISetDeflationThreshold_C",(NEP,PetscReal),(nep,deftol));
552: return 0;
553: }
555: static PetscErrorCode NEPRIIGetDeflationThreshold_RII(NEP nep,PetscReal *deftol)
556: {
557: NEP_RII *ctx = (NEP_RII*)nep->data;
559: *deftol = ctx->deftol;
560: return 0;
561: }
563: /*@
564: NEPRIIGetDeflationThreshold - Returns the threshold value that controls deflation.
566: Not Collective
568: Input Parameter:
569: . nep - nonlinear eigenvalue solver
571: Output Parameter:
572: . deftol - the threshold
574: Level: advanced
576: .seealso: NEPRIISetDeflationThreshold()
577: @*/
578: PetscErrorCode NEPRIIGetDeflationThreshold(NEP nep,PetscReal *deftol)
579: {
582: PetscUseMethod(nep,"NEPRIIGetDeflationThreshold_C",(NEP,PetscReal*),(nep,deftol));
583: return 0;
584: }
586: static PetscErrorCode NEPRIISetKSP_RII(NEP nep,KSP ksp)
587: {
588: NEP_RII *ctx = (NEP_RII*)nep->data;
590: PetscObjectReference((PetscObject)ksp);
591: KSPDestroy(&ctx->ksp);
592: ctx->ksp = ksp;
593: nep->state = NEP_STATE_INITIAL;
594: return 0;
595: }
597: /*@
598: NEPRIISetKSP - Associate a linear solver object (KSP) to the nonlinear
599: eigenvalue solver.
601: Collective on nep
603: Input Parameters:
604: + nep - eigenvalue solver
605: - ksp - the linear solver object
607: Level: advanced
609: .seealso: NEPRIIGetKSP()
610: @*/
611: PetscErrorCode NEPRIISetKSP(NEP nep,KSP ksp)
612: {
616: PetscTryMethod(nep,"NEPRIISetKSP_C",(NEP,KSP),(nep,ksp));
617: return 0;
618: }
620: static PetscErrorCode NEPRIIGetKSP_RII(NEP nep,KSP *ksp)
621: {
622: NEP_RII *ctx = (NEP_RII*)nep->data;
624: if (!ctx->ksp) {
625: KSPCreate(PetscObjectComm((PetscObject)nep),&ctx->ksp);
626: PetscObjectIncrementTabLevel((PetscObject)ctx->ksp,(PetscObject)nep,1);
627: KSPSetOptionsPrefix(ctx->ksp,((PetscObject)nep)->prefix);
628: KSPAppendOptionsPrefix(ctx->ksp,"nep_rii_");
629: PetscObjectSetOptions((PetscObject)ctx->ksp,((PetscObject)nep)->options);
630: KSPSetErrorIfNotConverged(ctx->ksp,PETSC_TRUE);
631: KSPSetTolerances(ctx->ksp,SlepcDefaultTol(nep->tol),PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
632: }
633: *ksp = ctx->ksp;
634: return 0;
635: }
637: /*@
638: NEPRIIGetKSP - Retrieve the linear solver object (KSP) associated with
639: the nonlinear eigenvalue solver.
641: Not Collective
643: Input Parameter:
644: . nep - nonlinear eigenvalue solver
646: Output Parameter:
647: . ksp - the linear solver object
649: Level: advanced
651: .seealso: NEPRIISetKSP()
652: @*/
653: PetscErrorCode NEPRIIGetKSP(NEP nep,KSP *ksp)
654: {
657: PetscUseMethod(nep,"NEPRIIGetKSP_C",(NEP,KSP*),(nep,ksp));
658: return 0;
659: }
661: PetscErrorCode NEPView_RII(NEP nep,PetscViewer viewer)
662: {
663: NEP_RII *ctx = (NEP_RII*)nep->data;
664: PetscBool isascii;
666: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
667: if (isascii) {
668: PetscViewerASCIIPrintf(viewer," maximum number of inner iterations: %" PetscInt_FMT "\n",ctx->max_inner_it);
669: if (ctx->cctol) PetscViewerASCIIPrintf(viewer," using a constant tolerance for the linear solver\n");
670: if (ctx->herm) PetscViewerASCIIPrintf(viewer," using the Hermitian version of the scalar nonlinear equation\n");
671: if (ctx->lag) PetscViewerASCIIPrintf(viewer," updating the preconditioner every %" PetscInt_FMT " iterations\n",ctx->lag);
672: if (ctx->deftol) PetscViewerASCIIPrintf(viewer," deflation threshold: %g\n",(double)ctx->deftol);
673: if (!ctx->ksp) NEPRIIGetKSP(nep,&ctx->ksp);
674: PetscViewerASCIIPushTab(viewer);
675: KSPView(ctx->ksp,viewer);
676: PetscViewerASCIIPopTab(viewer);
677: }
678: return 0;
679: }
681: PetscErrorCode NEPReset_RII(NEP nep)
682: {
683: NEP_RII *ctx = (NEP_RII*)nep->data;
685: KSPReset(ctx->ksp);
686: return 0;
687: }
689: PetscErrorCode NEPDestroy_RII(NEP nep)
690: {
691: NEP_RII *ctx = (NEP_RII*)nep->data;
693: KSPDestroy(&ctx->ksp);
694: PetscFree(nep->data);
695: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetMaximumIterations_C",NULL);
696: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetMaximumIterations_C",NULL);
697: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetLagPreconditioner_C",NULL);
698: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetLagPreconditioner_C",NULL);
699: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetConstCorrectionTol_C",NULL);
700: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetConstCorrectionTol_C",NULL);
701: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetHermitian_C",NULL);
702: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetHermitian_C",NULL);
703: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetDeflationThreshold_C",NULL);
704: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetDeflationThreshold_C",NULL);
705: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetKSP_C",NULL);
706: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetKSP_C",NULL);
707: return 0;
708: }
710: SLEPC_EXTERN PetscErrorCode NEPCreate_RII(NEP nep)
711: {
712: NEP_RII *ctx;
714: PetscNew(&ctx);
715: nep->data = (void*)ctx;
716: ctx->max_inner_it = 10;
717: ctx->lag = 1;
718: ctx->cctol = PETSC_FALSE;
719: ctx->herm = PETSC_FALSE;
720: ctx->deftol = 0.0;
722: nep->useds = PETSC_TRUE;
724: nep->ops->solve = NEPSolve_RII;
725: nep->ops->setup = NEPSetUp_RII;
726: nep->ops->setfromoptions = NEPSetFromOptions_RII;
727: nep->ops->reset = NEPReset_RII;
728: nep->ops->destroy = NEPDestroy_RII;
729: nep->ops->view = NEPView_RII;
730: nep->ops->computevectors = NEPComputeVectors_Schur;
732: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetMaximumIterations_C",NEPRIISetMaximumIterations_RII);
733: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetMaximumIterations_C",NEPRIIGetMaximumIterations_RII);
734: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetLagPreconditioner_C",NEPRIISetLagPreconditioner_RII);
735: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetLagPreconditioner_C",NEPRIIGetLagPreconditioner_RII);
736: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetConstCorrectionTol_C",NEPRIISetConstCorrectionTol_RII);
737: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetConstCorrectionTol_C",NEPRIIGetConstCorrectionTol_RII);
738: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetHermitian_C",NEPRIISetHermitian_RII);
739: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetHermitian_C",NEPRIIGetHermitian_RII);
740: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetDeflationThreshold_C",NEPRIISetDeflationThreshold_RII);
741: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetDeflationThreshold_C",NEPRIIGetDeflationThreshold_RII);
742: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetKSP_C",NEPRIISetKSP_RII);
743: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetKSP_C",NEPRIIGetKSP_RII);
744: return 0;
745: }