Actual source code: rii.c
slepc-3.16.0 2021-09-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2021, 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: {
42: if (nep->ncv!=PETSC_DEFAULT) { PetscInfo(nep,"Setting ncv = nev, ignoring user-provided value\n"); }
43: nep->ncv = nep->nev;
44: if (nep->mpd!=PETSC_DEFAULT) { PetscInfo(nep,"Setting mpd = nev, ignoring user-provided value\n"); }
45: nep->mpd = nep->nev;
46: if (nep->max_it==PETSC_DEFAULT) nep->max_it = PetscMax(5000,2*nep->n/nep->ncv);
47: if (!nep->which) nep->which = NEP_TARGET_MAGNITUDE;
48: if (nep->which!=NEP_TARGET_MAGNITUDE) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"This solver supports only target magnitude eigenvalues");
49: NEPCheckUnsupported(nep,NEP_FEATURE_REGION | NEP_FEATURE_TWOSIDED);
50: NEPAllocateSolution(nep,0);
51: NEPSetWorkVecs(nep,2);
52: return(0);
53: }
55: PetscErrorCode NEPSolve_RII(NEP nep)
56: {
57: PetscErrorCode ierr;
58: NEP_RII *ctx = (NEP_RII*)nep->data;
59: Mat T,Tp,H;
60: Vec uu,u,r,delta,t;
61: PetscScalar lambda,lambda2,sigma,a1,a2,corr,*Ap;
62: const PetscScalar *Hp;
63: PetscReal nrm,resnorm=1.0,ktol=0.1,perr,rtol;
64: PetscBool skip=PETSC_FALSE,lock=PETSC_FALSE;
65: PetscInt inner_its,its=0,ldh,ldds,i,j;
66: NEP_EXT_OP extop=NULL;
67: KSPConvergedReason kspreason;
70: /* get initial approximation of eigenvalue and eigenvector */
71: NEPGetDefaultShift(nep,&sigma);
72: lambda = sigma;
73: if (!nep->nini) {
74: BVSetRandomColumn(nep->V,0);
75: BVNormColumn(nep->V,0,NORM_2,&nrm);
76: BVScaleColumn(nep->V,0,1.0/nrm);
77: }
78: if (!ctx->ksp) { NEPRIIGetKSP(nep,&ctx->ksp); }
79: NEPDeflationInitialize(nep,nep->V,ctx->ksp,PETSC_FALSE,nep->nev,&extop);
80: NEPDeflationCreateVec(extop,&u);
81: VecDuplicate(u,&r);
82: VecDuplicate(u,&delta);
83: BVGetColumn(nep->V,0,&uu);
84: NEPDeflationCopyToExtendedVec(extop,uu,NULL,u,PETSC_FALSE);
85: BVRestoreColumn(nep->V,0,&uu);
87: /* prepare linear solver */
88: NEPDeflationSolveSetUp(extop,sigma);
89: KSPGetTolerances(ctx->ksp,&rtol,NULL,NULL,NULL);
91: VecCopy(u,r);
92: NEPDeflationFunctionSolve(extop,r,u);
93: VecNorm(u,NORM_2,&nrm);
94: VecScale(u,1.0/nrm);
96: /* Restart loop */
97: while (nep->reason == NEP_CONVERGED_ITERATING) {
98: its++;
100: /* Use Newton's method to compute nonlinear Rayleigh functional. Current eigenvalue
101: estimate as starting value. */
102: inner_its=0;
103: lambda2 = lambda;
104: do {
105: NEPDeflationComputeFunction(extop,lambda,&T);
106: MatMult(T,u,r);
107: if (!ctx->herm) {
108: NEPDeflationFunctionSolve(extop,r,delta);
109: KSPGetConvergedReason(ctx->ksp,&kspreason);
110: if (kspreason<0) {
111: PetscInfo1(nep,"iter=%D, linear solve failed\n",nep->its);
112: }
113: t = delta;
114: } else t = r;
115: VecDot(t,u,&a1);
116: NEPDeflationComputeJacobian(extop,lambda,&Tp);
117: MatMult(Tp,u,r);
118: if (!ctx->herm) {
119: NEPDeflationFunctionSolve(extop,r,delta);
120: KSPGetConvergedReason(ctx->ksp,&kspreason);
121: if (kspreason<0) {
122: PetscInfo1(nep,"iter=%D, linear solve failed\n",nep->its);
123: }
124: t = delta;
125: } else t = r;
126: VecDot(t,u,&a2);
127: corr = a1/a2;
128: lambda = lambda - corr;
129: inner_its++;
130: } while (PetscAbsScalar(corr)/PetscAbsScalar(lambda)>PETSC_SQRT_MACHINE_EPSILON && inner_its<ctx->max_inner_it);
132: /* form residual, r = T(lambda)*u */
133: NEPDeflationComputeFunction(extop,lambda,&T);
134: MatMult(T,u,r);
136: /* convergence test */
137: perr = nep->errest[nep->nconv];
138: VecNorm(r,NORM_2,&resnorm);
139: (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
140: nep->eigr[nep->nconv] = lambda;
141: if (its>1 && (nep->errest[nep->nconv]<=nep->tol || nep->errest[nep->nconv]<=ctx->deftol)) {
142: if (nep->errest[nep->nconv]<=ctx->deftol && !extop->ref && nep->nconv) {
143: NEPDeflationExtractEigenpair(extop,nep->nconv,u,lambda,nep->ds);
144: NEPDeflationSetRefine(extop,PETSC_TRUE);
145: MatMult(T,u,r);
146: VecNorm(r,NORM_2,&resnorm);
147: (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
148: if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
149: } else if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
150: }
151: if (lock) {
152: NEPDeflationSetRefine(extop,PETSC_FALSE);
153: nep->nconv = nep->nconv + 1;
154: NEPDeflationLocking(extop,u,lambda);
155: nep->its += its;
156: skip = PETSC_TRUE;
157: lock = PETSC_FALSE;
158: }
159: (*nep->stopping)(nep,nep->its+its,nep->max_it,nep->nconv,nep->nev,&nep->reason,nep->stoppingctx);
160: if (!skip || nep->reason>0) {
161: NEPMonitor(nep,nep->its+its,nep->nconv,nep->eigr,nep->eigi,nep->errest,(nep->reason>0)?nep->nconv:nep->nconv+1);
162: }
164: if (nep->reason == NEP_CONVERGED_ITERATING) {
165: if (!skip) {
166: /* update preconditioner and set adaptive tolerance */
167: if (ctx->lag && !(its%ctx->lag) && its>=2*ctx->lag && perr && nep->errest[nep->nconv]>.5*perr) {
168: NEPDeflationSolveSetUp(extop,lambda2);
169: }
170: if (!ctx->cctol) {
171: ktol = PetscMax(ktol/2.0,rtol);
172: KSPSetTolerances(ctx->ksp,ktol,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
173: }
175: /* eigenvector correction: delta = T(sigma)\r */
176: NEPDeflationFunctionSolve(extop,r,delta);
177: KSPGetConvergedReason(ctx->ksp,&kspreason);
178: if (kspreason<0) {
179: PetscInfo1(nep,"iter=%D, linear solve failed, stopping solve\n",nep->its);
180: nep->reason = NEP_DIVERGED_LINEAR_SOLVE;
181: break;
182: }
184: /* update eigenvector: u = u - delta */
185: VecAXPY(u,-1.0,delta);
187: /* normalize eigenvector */
188: VecNormalize(u,NULL);
189: } else {
190: its = -1;
191: NEPDeflationSetRandomVec(extop,u);
192: NEPDeflationSolveSetUp(extop,sigma);
193: VecCopy(u,r);
194: NEPDeflationFunctionSolve(extop,r,u);
195: VecNorm(u,NORM_2,&nrm);
196: VecScale(u,1.0/nrm);
197: lambda = sigma;
198: skip = PETSC_FALSE;
199: }
200: }
201: }
202: NEPDeflationGetInvariantPair(extop,NULL,&H);
203: MatGetSize(H,NULL,&ldh);
204: DSSetType(nep->ds,DSNHEP);
205: DSReset(nep->ds);
206: DSAllocate(nep->ds,PetscMax(nep->nconv,1));
207: DSGetLeadingDimension(nep->ds,&ldds);
208: MatDenseGetArrayRead(H,&Hp);
209: DSGetArray(nep->ds,DS_MAT_A,&Ap);
210: for (j=0;j<nep->nconv;j++)
211: for (i=0;i<nep->nconv;i++) Ap[j*ldds+i] = Hp[j*ldh+i];
212: DSRestoreArray(nep->ds,DS_MAT_A,&Ap);
213: MatDenseRestoreArrayRead(H,&Hp);
214: MatDestroy(&H);
215: DSSetDimensions(nep->ds,nep->nconv,0,nep->nconv);
216: DSSolve(nep->ds,nep->eigr,nep->eigi);
217: NEPDeflationReset(extop);
218: VecDestroy(&u);
219: VecDestroy(&r);
220: VecDestroy(&delta);
221: return(0);
222: }
224: PetscErrorCode NEPSetFromOptions_RII(PetscOptionItems *PetscOptionsObject,NEP nep)
225: {
227: NEP_RII *ctx = (NEP_RII*)nep->data;
228: PetscBool flg;
229: PetscInt i;
230: PetscReal r;
233: PetscOptionsHead(PetscOptionsObject,"NEP RII Options");
235: i = 0;
236: PetscOptionsInt("-nep_rii_max_it","Maximum number of Newton iterations for updating Rayleigh functional","NEPRIISetMaximumIterations",ctx->max_inner_it,&i,&flg);
237: if (flg) { NEPRIISetMaximumIterations(nep,i); }
239: PetscOptionsBool("-nep_rii_const_correction_tol","Constant correction tolerance for the linear solver","NEPRIISetConstCorrectionTol",ctx->cctol,&ctx->cctol,NULL);
241: PetscOptionsBool("-nep_rii_hermitian","Use Hermitian version of the scalar nonlinear equation","NEPRIISetHermitian",ctx->herm,&ctx->herm,NULL);
243: i = 0;
244: PetscOptionsInt("-nep_rii_lag_preconditioner","Interval to rebuild preconditioner","NEPRIISetLagPreconditioner",ctx->lag,&i,&flg);
245: if (flg) { NEPRIISetLagPreconditioner(nep,i); }
247: r = 0.0;
248: PetscOptionsReal("-nep_rii_deflation_threshold","Tolerance used as a threshold for including deflated eigenpairs","NEPRIISetDeflationThreshold",ctx->deftol,&r,&flg);
249: if (flg) { NEPRIISetDeflationThreshold(nep,r); }
251: PetscOptionsTail();
253: if (!ctx->ksp) { NEPRIIGetKSP(nep,&ctx->ksp); }
254: KSPSetFromOptions(ctx->ksp);
255: return(0);
256: }
258: static PetscErrorCode NEPRIISetMaximumIterations_RII(NEP nep,PetscInt its)
259: {
260: NEP_RII *ctx = (NEP_RII*)nep->data;
263: if (its==PETSC_DEFAULT) ctx->max_inner_it = 10;
264: else {
265: if (its<=0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Number of iterations must be >0");
266: ctx->max_inner_it = its;
267: }
268: return(0);
269: }
271: /*@
272: NEPRIISetMaximumIterations - Sets the maximum number of inner iterations to be
273: used in the RII solver. These are the Newton iterations related to the computation
274: of the nonlinear Rayleigh functional.
276: Logically Collective on nep
278: Input Parameters:
279: + nep - nonlinear eigenvalue solver
280: - its - maximum inner iterations
282: Level: advanced
284: .seealso: NEPRIIGetMaximumIterations()
285: @*/
286: PetscErrorCode NEPRIISetMaximumIterations(NEP nep,PetscInt its)
287: {
293: PetscTryMethod(nep,"NEPRIISetMaximumIterations_C",(NEP,PetscInt),(nep,its));
294: return(0);
295: }
297: static PetscErrorCode NEPRIIGetMaximumIterations_RII(NEP nep,PetscInt *its)
298: {
299: NEP_RII *ctx = (NEP_RII*)nep->data;
302: *its = ctx->max_inner_it;
303: return(0);
304: }
306: /*@
307: NEPRIIGetMaximumIterations - Gets the maximum number of inner iterations of RII.
309: Not Collective
311: Input Parameter:
312: . nep - nonlinear eigenvalue solver
314: Output Parameter:
315: . its - maximum inner iterations
317: Level: advanced
319: .seealso: NEPRIISetMaximumIterations()
320: @*/
321: PetscErrorCode NEPRIIGetMaximumIterations(NEP nep,PetscInt *its)
322: {
328: PetscUseMethod(nep,"NEPRIIGetMaximumIterations_C",(NEP,PetscInt*),(nep,its));
329: return(0);
330: }
332: static PetscErrorCode NEPRIISetLagPreconditioner_RII(NEP nep,PetscInt lag)
333: {
334: NEP_RII *ctx = (NEP_RII*)nep->data;
337: if (lag<0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Lag must be non-negative");
338: ctx->lag = lag;
339: return(0);
340: }
342: /*@
343: NEPRIISetLagPreconditioner - Determines when the preconditioner is rebuilt in the
344: nonlinear solve.
346: Logically Collective on nep
348: Input Parameters:
349: + nep - nonlinear eigenvalue solver
350: - lag - 0 indicates NEVER rebuild, 1 means rebuild every time the Jacobian is
351: computed within the nonlinear iteration, 2 means every second time
352: the Jacobian is built, etc.
354: Options Database Keys:
355: . -nep_rii_lag_preconditioner <lag> - the lag value
357: Notes:
358: The default is 1.
359: The preconditioner is ALWAYS built in the first iteration of a nonlinear solve.
361: Level: intermediate
363: .seealso: NEPRIIGetLagPreconditioner()
364: @*/
365: PetscErrorCode NEPRIISetLagPreconditioner(NEP nep,PetscInt lag)
366: {
372: PetscTryMethod(nep,"NEPRIISetLagPreconditioner_C",(NEP,PetscInt),(nep,lag));
373: return(0);
374: }
376: static PetscErrorCode NEPRIIGetLagPreconditioner_RII(NEP nep,PetscInt *lag)
377: {
378: NEP_RII *ctx = (NEP_RII*)nep->data;
381: *lag = ctx->lag;
382: return(0);
383: }
385: /*@
386: NEPRIIGetLagPreconditioner - Indicates how often the preconditioner is rebuilt.
388: Not Collective
390: Input Parameter:
391: . nep - nonlinear eigenvalue solver
393: Output Parameter:
394: . lag - the lag parameter
396: Level: intermediate
398: .seealso: NEPRIISetLagPreconditioner()
399: @*/
400: PetscErrorCode NEPRIIGetLagPreconditioner(NEP nep,PetscInt *lag)
401: {
407: PetscUseMethod(nep,"NEPRIIGetLagPreconditioner_C",(NEP,PetscInt*),(nep,lag));
408: return(0);
409: }
411: static PetscErrorCode NEPRIISetConstCorrectionTol_RII(NEP nep,PetscBool cct)
412: {
413: NEP_RII *ctx = (NEP_RII*)nep->data;
416: ctx->cctol = cct;
417: return(0);
418: }
420: /*@
421: NEPRIISetConstCorrectionTol - Sets a flag to keep the tolerance used
422: in the linear solver constant.
424: Logically Collective on nep
426: Input Parameters:
427: + nep - nonlinear eigenvalue solver
428: - cct - a boolean value
430: Options Database Keys:
431: . -nep_rii_const_correction_tol <bool> - set the boolean flag
433: Notes:
434: By default, an exponentially decreasing tolerance is set in the KSP used
435: within the nonlinear iteration, so that each Newton iteration requests
436: better accuracy than the previous one. The constant correction tolerance
437: flag stops this behaviour.
439: Level: intermediate
441: .seealso: NEPRIIGetConstCorrectionTol()
442: @*/
443: PetscErrorCode NEPRIISetConstCorrectionTol(NEP nep,PetscBool cct)
444: {
450: PetscTryMethod(nep,"NEPRIISetConstCorrectionTol_C",(NEP,PetscBool),(nep,cct));
451: return(0);
452: }
454: static PetscErrorCode NEPRIIGetConstCorrectionTol_RII(NEP nep,PetscBool *cct)
455: {
456: NEP_RII *ctx = (NEP_RII*)nep->data;
459: *cct = ctx->cctol;
460: return(0);
461: }
463: /*@
464: NEPRIIGetConstCorrectionTol - Returns the constant tolerance flag.
466: Not Collective
468: Input Parameter:
469: . nep - nonlinear eigenvalue solver
471: Output Parameter:
472: . cct - the value of the constant tolerance flag
474: Level: intermediate
476: .seealso: NEPRIISetConstCorrectionTol()
477: @*/
478: PetscErrorCode NEPRIIGetConstCorrectionTol(NEP nep,PetscBool *cct)
479: {
485: PetscUseMethod(nep,"NEPRIIGetConstCorrectionTol_C",(NEP,PetscBool*),(nep,cct));
486: return(0);
487: }
489: static PetscErrorCode NEPRIISetHermitian_RII(NEP nep,PetscBool herm)
490: {
491: NEP_RII *ctx = (NEP_RII*)nep->data;
494: ctx->herm = herm;
495: return(0);
496: }
498: /*@
499: NEPRIISetHermitian - Sets a flag to indicate if the Hermitian version of the
500: scalar nonlinear equation must be used by the solver.
502: Logically Collective on nep
504: Input Parameters:
505: + nep - nonlinear eigenvalue solver
506: - herm - a boolean value
508: Options Database Keys:
509: . -nep_rii_hermitian <bool> - set the boolean flag
511: Notes:
512: By default, the scalar nonlinear equation x'*inv(T(sigma))*T(z)*x=0 is solved
513: at each step of the nonlinear iteration. When this flag is set the simpler
514: form x'*T(z)*x=0 is used, which is supposed to be valid only for Hermitian
515: problems.
517: Level: intermediate
519: .seealso: NEPRIIGetHermitian()
520: @*/
521: PetscErrorCode NEPRIISetHermitian(NEP nep,PetscBool herm)
522: {
528: PetscTryMethod(nep,"NEPRIISetHermitian_C",(NEP,PetscBool),(nep,herm));
529: return(0);
530: }
532: static PetscErrorCode NEPRIIGetHermitian_RII(NEP nep,PetscBool *herm)
533: {
534: NEP_RII *ctx = (NEP_RII*)nep->data;
537: *herm = ctx->herm;
538: return(0);
539: }
541: /*@
542: NEPRIIGetHermitian - Returns the flag about using the Hermitian version of
543: the scalar nonlinear equation.
545: Not Collective
547: Input Parameter:
548: . nep - nonlinear eigenvalue solver
550: Output Parameter:
551: . herm - the value of the hermitian flag
553: Level: intermediate
555: .seealso: NEPRIISetHermitian()
556: @*/
557: PetscErrorCode NEPRIIGetHermitian(NEP nep,PetscBool *herm)
558: {
564: PetscUseMethod(nep,"NEPRIIGetHermitian_C",(NEP,PetscBool*),(nep,herm));
565: return(0);
566: }
568: static PetscErrorCode NEPRIISetDeflationThreshold_RII(NEP nep,PetscReal deftol)
569: {
570: NEP_RII *ctx = (NEP_RII*)nep->data;
573: ctx->deftol = deftol;
574: return(0);
575: }
577: /*@
578: NEPRIISetDeflationThreshold - Sets the threshold value used to switch between
579: deflated and non-deflated iteration.
581: Logically Collective on nep
583: Input Parameters:
584: + nep - nonlinear eigenvalue solver
585: - deftol - the threshold value
587: Options Database Keys:
588: . -nep_rii_deflation_threshold <deftol> - set the threshold
590: Notes:
591: Normally, the solver iterates on the extended problem in order to deflate
592: previously converged eigenpairs. If this threshold is set to a nonzero value,
593: then once the residual error is below this threshold the solver will
594: continue the iteration without deflation. The intention is to be able to
595: improve the current eigenpair further, despite having previous eigenpairs
596: with somewhat bad precision.
598: Level: advanced
600: .seealso: NEPRIIGetDeflationThreshold()
601: @*/
602: PetscErrorCode NEPRIISetDeflationThreshold(NEP nep,PetscReal deftol)
603: {
609: PetscTryMethod(nep,"NEPRIISetDeflationThreshold_C",(NEP,PetscReal),(nep,deftol));
610: return(0);
611: }
613: static PetscErrorCode NEPRIIGetDeflationThreshold_RII(NEP nep,PetscReal *deftol)
614: {
615: NEP_RII *ctx = (NEP_RII*)nep->data;
618: *deftol = ctx->deftol;
619: return(0);
620: }
622: /*@
623: NEPRIIGetDeflationThreshold - Returns the threshold value that controls deflation.
625: Not Collective
627: Input Parameter:
628: . nep - nonlinear eigenvalue solver
630: Output Parameter:
631: . deftol - the threshold
633: Level: advanced
635: .seealso: NEPRIISetDeflationThreshold()
636: @*/
637: PetscErrorCode NEPRIIGetDeflationThreshold(NEP nep,PetscReal *deftol)
638: {
644: PetscUseMethod(nep,"NEPRIIGetDeflationThreshold_C",(NEP,PetscReal*),(nep,deftol));
645: return(0);
646: }
648: static PetscErrorCode NEPRIISetKSP_RII(NEP nep,KSP ksp)
649: {
651: NEP_RII *ctx = (NEP_RII*)nep->data;
654: PetscObjectReference((PetscObject)ksp);
655: KSPDestroy(&ctx->ksp);
656: ctx->ksp = ksp;
657: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);
658: nep->state = NEP_STATE_INITIAL;
659: return(0);
660: }
662: /*@
663: NEPRIISetKSP - Associate a linear solver object (KSP) to the nonlinear
664: eigenvalue solver.
666: Collective on nep
668: Input Parameters:
669: + nep - eigenvalue solver
670: - ksp - the linear solver object
672: Level: advanced
674: .seealso: NEPRIIGetKSP()
675: @*/
676: PetscErrorCode NEPRIISetKSP(NEP nep,KSP ksp)
677: {
684: PetscTryMethod(nep,"NEPRIISetKSP_C",(NEP,KSP),(nep,ksp));
685: return(0);
686: }
688: static PetscErrorCode NEPRIIGetKSP_RII(NEP nep,KSP *ksp)
689: {
691: NEP_RII *ctx = (NEP_RII*)nep->data;
694: if (!ctx->ksp) {
695: KSPCreate(PetscObjectComm((PetscObject)nep),&ctx->ksp);
696: PetscObjectIncrementTabLevel((PetscObject)ctx->ksp,(PetscObject)nep,1);
697: KSPSetOptionsPrefix(ctx->ksp,((PetscObject)nep)->prefix);
698: KSPAppendOptionsPrefix(ctx->ksp,"nep_rii_");
699: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);
700: PetscObjectSetOptions((PetscObject)ctx->ksp,((PetscObject)nep)->options);
701: KSPSetErrorIfNotConverged(ctx->ksp,PETSC_TRUE);
702: KSPSetTolerances(ctx->ksp,SlepcDefaultTol(nep->tol),PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
703: }
704: *ksp = ctx->ksp;
705: return(0);
706: }
708: /*@
709: NEPRIIGetKSP - Retrieve the linear solver object (KSP) associated with
710: the nonlinear eigenvalue solver.
712: Not Collective
714: Input Parameter:
715: . nep - nonlinear eigenvalue solver
717: Output Parameter:
718: . ksp - the linear solver object
720: Level: advanced
722: .seealso: NEPRIISetKSP()
723: @*/
724: PetscErrorCode NEPRIIGetKSP(NEP nep,KSP *ksp)
725: {
731: PetscUseMethod(nep,"NEPRIIGetKSP_C",(NEP,KSP*),(nep,ksp));
732: return(0);
733: }
735: PetscErrorCode NEPView_RII(NEP nep,PetscViewer viewer)
736: {
738: NEP_RII *ctx = (NEP_RII*)nep->data;
739: PetscBool isascii;
742: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
743: if (isascii) {
744: PetscViewerASCIIPrintf(viewer," maximum number of inner iterations: %D\n",ctx->max_inner_it);
745: if (ctx->cctol) {
746: PetscViewerASCIIPrintf(viewer," using a constant tolerance for the linear solver\n");
747: }
748: if (ctx->herm) {
749: PetscViewerASCIIPrintf(viewer," using the Hermitian version of the scalar nonlinear equation\n");
750: }
751: if (ctx->lag) {
752: PetscViewerASCIIPrintf(viewer," updating the preconditioner every %D iterations\n",ctx->lag);
753: }
754: if (ctx->deftol) {
755: PetscViewerASCIIPrintf(viewer," deflation threshold: %g\n",(double)ctx->deftol);
756: }
757: if (!ctx->ksp) { NEPRIIGetKSP(nep,&ctx->ksp); }
758: PetscViewerASCIIPushTab(viewer);
759: KSPView(ctx->ksp,viewer);
760: PetscViewerASCIIPopTab(viewer);
761: }
762: return(0);
763: }
765: PetscErrorCode NEPReset_RII(NEP nep)
766: {
768: NEP_RII *ctx = (NEP_RII*)nep->data;
771: KSPReset(ctx->ksp);
772: return(0);
773: }
775: PetscErrorCode NEPDestroy_RII(NEP nep)
776: {
778: NEP_RII *ctx = (NEP_RII*)nep->data;
781: KSPDestroy(&ctx->ksp);
782: PetscFree(nep->data);
783: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetMaximumIterations_C",NULL);
784: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetMaximumIterations_C",NULL);
785: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetLagPreconditioner_C",NULL);
786: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetLagPreconditioner_C",NULL);
787: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetConstCorrectionTol_C",NULL);
788: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetConstCorrectionTol_C",NULL);
789: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetHermitian_C",NULL);
790: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetHermitian_C",NULL);
791: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetDeflationThreshold_C",NULL);
792: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetDeflationThreshold_C",NULL);
793: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetKSP_C",NULL);
794: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetKSP_C",NULL);
795: return(0);
796: }
798: SLEPC_EXTERN PetscErrorCode NEPCreate_RII(NEP nep)
799: {
801: NEP_RII *ctx;
804: PetscNewLog(nep,&ctx);
805: nep->data = (void*)ctx;
806: ctx->max_inner_it = 10;
807: ctx->lag = 1;
808: ctx->cctol = PETSC_FALSE;
809: ctx->herm = PETSC_FALSE;
810: ctx->deftol = 0.0;
812: nep->useds = PETSC_TRUE;
814: nep->ops->solve = NEPSolve_RII;
815: nep->ops->setup = NEPSetUp_RII;
816: nep->ops->setfromoptions = NEPSetFromOptions_RII;
817: nep->ops->reset = NEPReset_RII;
818: nep->ops->destroy = NEPDestroy_RII;
819: nep->ops->view = NEPView_RII;
820: nep->ops->computevectors = NEPComputeVectors_Schur;
822: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetMaximumIterations_C",NEPRIISetMaximumIterations_RII);
823: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetMaximumIterations_C",NEPRIIGetMaximumIterations_RII);
824: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetLagPreconditioner_C",NEPRIISetLagPreconditioner_RII);
825: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetLagPreconditioner_C",NEPRIIGetLagPreconditioner_RII);
826: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetConstCorrectionTol_C",NEPRIISetConstCorrectionTol_RII);
827: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetConstCorrectionTol_C",NEPRIIGetConstCorrectionTol_RII);
828: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetHermitian_C",NEPRIISetHermitian_RII);
829: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetHermitian_C",NEPRIIGetHermitian_RII);
830: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetDeflationThreshold_C",NEPRIISetDeflationThreshold_RII);
831: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetDeflationThreshold_C",NEPRIIGetDeflationThreshold_RII);
832: PetscObjectComposeFunction((PetscObject)nep,"NEPRIISetKSP_C",NEPRIISetKSP_RII);
833: PetscObjectComposeFunction((PetscObject)nep,"NEPRIIGetKSP_C",NEPRIIGetKSP_RII);
834: return(0);
835: }