Actual source code: slp.c
slepc-3.16.1 2021-11-17
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: "slp"
13: Method: Successive linear problems
15: Algorithm:
17: Newton-type iteration based on first order Taylor approximation.
19: References:
21: [1] A. Ruhe, "Algorithms for the nonlinear eigenvalue problem", SIAM J.
22: Numer. Anal. 10(4):674-689, 1973.
23: */
25: #include <slepc/private/nepimpl.h>
26: #include <../src/nep/impls/nepdefl.h>
27: #include "slp.h"
29: typedef struct {
30: NEP_EXT_OP extop;
31: Vec w;
32: } NEP_SLP_MATSHELL;
34: PetscErrorCode NEPSetUp_SLP(NEP nep)
35: {
37: NEP_SLP *ctx = (NEP_SLP*)nep->data;
38: PetscBool flg;
39: ST st;
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->ncv>nep->nev+nep->mpd) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_USER_INPUT,"The value of ncv must not be larger than nev+mpd");
47: if (nep->max_it==PETSC_DEFAULT) nep->max_it = PetscMax(5000,2*nep->n/nep->ncv);
48: if (!nep->which) nep->which = NEP_TARGET_MAGNITUDE;
49: if (nep->which!=NEP_TARGET_MAGNITUDE) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"This solver supports only target magnitude eigenvalues");
50: NEPCheckUnsupported(nep,NEP_FEATURE_REGION);
52: if (!ctx->eps) { NEPSLPGetEPS(nep,&ctx->eps); }
53: EPSGetST(ctx->eps,&st);
54: PetscObjectTypeCompareAny((PetscObject)st,&flg,STSINVERT,STCAYLEY,"");
55: if (flg) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"SLP does not support spectral transformation");
56: EPSSetDimensions(ctx->eps,1,PETSC_DECIDE,PETSC_DECIDE);
57: EPSSetWhichEigenpairs(ctx->eps,EPS_LARGEST_MAGNITUDE);
58: EPSSetTolerances(ctx->eps,SlepcDefaultTol(nep->tol)/10.0,nep->max_it);
59: if (nep->tol==PETSC_DEFAULT) nep->tol = SLEPC_DEFAULT_TOL;
60: if (ctx->deftol==PETSC_DEFAULT) ctx->deftol = nep->tol;
62: if (nep->twosided) {
63: nep->ops->solve = NEPSolve_SLP_Twosided;
64: nep->ops->computevectors = NULL;
65: if (!ctx->epsts) { NEPSLPGetEPSLeft(nep,&ctx->epsts); }
66: EPSGetST(ctx->epsts,&st);
67: PetscObjectTypeCompareAny((PetscObject)st,&flg,STSINVERT,STCAYLEY,"");
68: if (flg) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"SLP does not support spectral transformation");
69: EPSSetDimensions(ctx->epsts,1,PETSC_DECIDE,PETSC_DECIDE);
70: EPSSetWhichEigenpairs(ctx->epsts,EPS_LARGEST_MAGNITUDE);
71: EPSSetTolerances(ctx->epsts,SlepcDefaultTol(nep->tol)/10.0,nep->max_it);
72: } else {
73: nep->ops->solve = NEPSolve_SLP;
74: nep->ops->computevectors = NEPComputeVectors_Schur;
75: }
76: NEPAllocateSolution(nep,0);
77: return(0);
78: }
80: static PetscErrorCode MatMult_SLP(Mat M,Vec x,Vec y)
81: {
82: PetscErrorCode ierr;
83: NEP_SLP_MATSHELL *ctx;
86: MatShellGetContext(M,&ctx);
87: MatMult(ctx->extop->MJ,x,ctx->w);
88: NEPDeflationFunctionSolve(ctx->extop,ctx->w,y);
89: return(0);
90: }
92: static PetscErrorCode MatDestroy_SLP(Mat M)
93: {
94: PetscErrorCode ierr;
95: NEP_SLP_MATSHELL *ctx;
98: MatShellGetContext(M,&ctx);
99: VecDestroy(&ctx->w);
100: PetscFree(ctx);
101: return(0);
102: }
104: #if defined(PETSC_HAVE_CUDA)
105: static PetscErrorCode MatCreateVecs_SLP(Mat M,Vec *left,Vec *right)
106: {
107: PetscErrorCode ierr;
108: NEP_SLP_MATSHELL *ctx;
111: MatShellGetContext(M,&ctx);
112: if (right) {
113: VecDuplicate(ctx->w,right);
114: }
115: if (left) {
116: VecDuplicate(ctx->w,left);
117: }
118: return(0);
119: }
120: #endif
122: static PetscErrorCode NEPSLPSetUpLinearEP(NEP nep,NEP_EXT_OP extop,PetscScalar lambda,Vec u,PetscBool ini)
123: {
124: PetscErrorCode ierr;
125: NEP_SLP *slpctx = (NEP_SLP*)nep->data;
126: Mat Mshell;
127: PetscInt nloc,mloc;
128: NEP_SLP_MATSHELL *shellctx;
131: if (ini) {
132: /* Create mat shell */
133: PetscNew(&shellctx);
134: shellctx->extop = extop;
135: NEPDeflationCreateVec(extop,&shellctx->w);
136: MatGetLocalSize(nep->function,&mloc,&nloc);
137: nloc += extop->szd; mloc += extop->szd;
138: MatCreateShell(PetscObjectComm((PetscObject)nep),nloc,mloc,PETSC_DETERMINE,PETSC_DETERMINE,shellctx,&Mshell);
139: MatShellSetOperation(Mshell,MATOP_MULT,(void(*)(void))MatMult_SLP);
140: MatShellSetOperation(Mshell,MATOP_DESTROY,(void(*)(void))MatDestroy_SLP);
141: #if defined(PETSC_HAVE_CUDA)
142: MatShellSetOperation(Mshell,MATOP_CREATE_VECS,(void(*)(void))MatCreateVecs_SLP);
143: #endif
144: EPSSetOperators(slpctx->eps,Mshell,NULL);
145: MatDestroy(&Mshell);
146: }
147: NEPDeflationSolveSetUp(extop,lambda);
148: NEPDeflationComputeJacobian(extop,lambda,NULL);
149: EPSSetInitialSpace(slpctx->eps,1,&u);
150: return(0);
151: }
153: PetscErrorCode NEPSolve_SLP(NEP nep)
154: {
155: PetscErrorCode ierr;
156: NEP_SLP *ctx = (NEP_SLP*)nep->data;
157: Mat F,H;
158: Vec uu,u,r;
159: PetscScalar sigma,lambda,mu,im,*Ap;
160: const PetscScalar *Hp;
161: PetscReal resnorm;
162: PetscInt nconv,ldh,ldds,i,j;
163: PetscBool skip=PETSC_FALSE,lock=PETSC_FALSE;
164: NEP_EXT_OP extop=NULL; /* Extended operator for deflation */
167: /* get initial approximation of eigenvalue and eigenvector */
168: NEPGetDefaultShift(nep,&sigma);
169: if (!nep->nini) {
170: BVSetRandomColumn(nep->V,0);
171: }
172: lambda = sigma;
173: if (!ctx->ksp) { NEPSLPGetKSP(nep,&ctx->ksp); }
174: NEPDeflationInitialize(nep,nep->V,ctx->ksp,PETSC_TRUE,nep->nev,&extop);
175: NEPDeflationCreateVec(extop,&u);
176: VecDuplicate(u,&r);
177: BVGetColumn(nep->V,0,&uu);
178: NEPDeflationCopyToExtendedVec(extop,uu,NULL,u,PETSC_FALSE);
179: BVRestoreColumn(nep->V,0,&uu);
181: /* Restart loop */
182: while (nep->reason == NEP_CONVERGED_ITERATING) {
183: nep->its++;
185: /* form residual, r = T(lambda)*u (used in convergence test only) */
186: NEPDeflationComputeFunction(extop,lambda,&F);
187: MatMult(F,u,r);
189: /* convergence test */
190: VecNorm(r,NORM_2,&resnorm);
191: (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
192: nep->eigr[nep->nconv] = lambda;
193: if (nep->errest[nep->nconv]<=nep->tol || nep->errest[nep->nconv]<=ctx->deftol) {
194: if (nep->errest[nep->nconv]<=ctx->deftol && !extop->ref && nep->nconv) {
195: NEPDeflationExtractEigenpair(extop,nep->nconv,u,lambda,nep->ds);
196: NEPDeflationSetRefine(extop,PETSC_TRUE);
197: MatMult(F,u,r);
198: VecNorm(r,NORM_2,&resnorm);
199: (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
200: if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
201: } else if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
202: }
204: if (lock) {
205: NEPDeflationSetRefine(extop,PETSC_FALSE);
206: nep->nconv = nep->nconv + 1;
207: skip = PETSC_TRUE;
208: lock = PETSC_FALSE;
209: NEPDeflationLocking(extop,u,lambda);
210: }
211: (*nep->stopping)(nep,nep->its,nep->max_it,nep->nconv,nep->nev,&nep->reason,nep->stoppingctx);
212: if (!skip || nep->reason>0) {
213: NEPMonitor(nep,nep->its,nep->nconv,nep->eigr,nep->eigi,nep->errest,(nep->reason>0)?nep->nconv:nep->nconv+1);
214: }
216: if (nep->reason == NEP_CONVERGED_ITERATING) {
217: if (!skip) {
218: /* evaluate T(lambda) and T'(lambda) */
219: NEPSLPSetUpLinearEP(nep,extop,lambda,u,nep->its==1?PETSC_TRUE:PETSC_FALSE);
220: /* compute new eigenvalue correction mu and eigenvector approximation u */
221: EPSSolve(ctx->eps);
222: EPSGetConverged(ctx->eps,&nconv);
223: if (!nconv) {
224: PetscInfo1(nep,"iter=%D, inner iteration failed, stopping solve\n",nep->its);
225: nep->reason = NEP_DIVERGED_LINEAR_SOLVE;
226: break;
227: }
228: EPSGetEigenpair(ctx->eps,0,&mu,&im,u,NULL);
229: mu = 1.0/mu;
230: if (PetscAbsScalar(im)>PETSC_MACHINE_EPSILON) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"Complex eigenvalue approximation - not implemented in real scalars");
231: } else {
232: nep->its--; /* do not count this as a full iteration */
233: /* use second eigenpair computed in previous iteration */
234: EPSGetConverged(ctx->eps,&nconv);
235: if (nconv>=2) {
236: EPSGetEigenpair(ctx->eps,1,&mu,&im,u,NULL);
237: mu = 1.0/mu;
238: } else {
239: NEPDeflationSetRandomVec(extop,u);
240: mu = lambda-sigma;
241: }
242: skip = PETSC_FALSE;
243: }
244: /* correct eigenvalue */
245: lambda = lambda - mu;
246: }
247: }
248: NEPDeflationGetInvariantPair(extop,NULL,&H);
249: MatGetSize(H,NULL,&ldh);
250: DSReset(nep->ds);
251: DSSetType(nep->ds,DSNHEP);
252: DSAllocate(nep->ds,PetscMax(nep->nconv,1));
253: DSGetLeadingDimension(nep->ds,&ldds);
254: MatDenseGetArrayRead(H,&Hp);
255: DSGetArray(nep->ds,DS_MAT_A,&Ap);
256: for (j=0;j<nep->nconv;j++)
257: for (i=0;i<nep->nconv;i++) Ap[j*ldds+i] = Hp[j*ldh+i];
258: DSRestoreArray(nep->ds,DS_MAT_A,&Ap);
259: MatDenseRestoreArrayRead(H,&Hp);
260: MatDestroy(&H);
261: DSSetDimensions(nep->ds,nep->nconv,0,nep->nconv);
262: DSSolve(nep->ds,nep->eigr,nep->eigi);
263: NEPDeflationReset(extop);
264: VecDestroy(&u);
265: VecDestroy(&r);
266: return(0);
267: }
269: PetscErrorCode NEPSetFromOptions_SLP(PetscOptionItems *PetscOptionsObject,NEP nep)
270: {
272: NEP_SLP *ctx = (NEP_SLP*)nep->data;
273: PetscBool flg;
274: PetscReal r;
277: PetscOptionsHead(PetscOptionsObject,"NEP SLP Options");
279: r = 0.0;
280: PetscOptionsReal("-nep_slp_deflation_threshold","Tolerance used as a threshold for including deflated eigenpairs","NEPSLPSetDeflationThreshold",ctx->deftol,&r,&flg);
281: if (flg) { NEPSLPSetDeflationThreshold(nep,r); }
283: PetscOptionsTail();
285: if (!ctx->eps) { NEPSLPGetEPS(nep,&ctx->eps); }
286: EPSSetFromOptions(ctx->eps);
287: if (nep->twosided) {
288: if (!ctx->epsts) { NEPSLPGetEPSLeft(nep,&ctx->epsts); }
289: EPSSetFromOptions(ctx->epsts);
290: }
291: if (!ctx->ksp) { NEPSLPGetKSP(nep,&ctx->ksp); }
292: KSPSetFromOptions(ctx->ksp);
293: return(0);
294: }
296: static PetscErrorCode NEPSLPSetDeflationThreshold_SLP(NEP nep,PetscReal deftol)
297: {
298: NEP_SLP *ctx = (NEP_SLP*)nep->data;
301: if (deftol == PETSC_DEFAULT) {
302: ctx->deftol = PETSC_DEFAULT;
303: nep->state = NEP_STATE_INITIAL;
304: } else {
305: if (deftol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of deftol. Must be > 0");
306: ctx->deftol = deftol;
307: }
308: return(0);
309: }
311: /*@
312: NEPSLPSetDeflationThreshold - Sets the threshold value used to switch between
313: deflated and non-deflated iteration.
315: Logically Collective on nep
317: Input Parameters:
318: + nep - nonlinear eigenvalue solver
319: - deftol - the threshold value
321: Options Database Keys:
322: . -nep_slp_deflation_threshold <deftol> - set the threshold
324: Notes:
325: Normally, the solver iterates on the extended problem in order to deflate
326: previously converged eigenpairs. If this threshold is set to a nonzero value,
327: then once the residual error is below this threshold the solver will
328: continue the iteration without deflation. The intention is to be able to
329: improve the current eigenpair further, despite having previous eigenpairs
330: with somewhat bad precision.
332: Level: advanced
334: .seealso: NEPSLPGetDeflationThreshold()
335: @*/
336: PetscErrorCode NEPSLPSetDeflationThreshold(NEP nep,PetscReal deftol)
337: {
343: PetscTryMethod(nep,"NEPSLPSetDeflationThreshold_C",(NEP,PetscReal),(nep,deftol));
344: return(0);
345: }
347: static PetscErrorCode NEPSLPGetDeflationThreshold_SLP(NEP nep,PetscReal *deftol)
348: {
349: NEP_SLP *ctx = (NEP_SLP*)nep->data;
352: *deftol = ctx->deftol;
353: return(0);
354: }
356: /*@
357: NEPSLPGetDeflationThreshold - Returns the threshold value that controls deflation.
359: Not Collective
361: Input Parameter:
362: . nep - nonlinear eigenvalue solver
364: Output Parameter:
365: . deftol - the threshold
367: Level: advanced
369: .seealso: NEPSLPSetDeflationThreshold()
370: @*/
371: PetscErrorCode NEPSLPGetDeflationThreshold(NEP nep,PetscReal *deftol)
372: {
378: PetscUseMethod(nep,"NEPSLPGetDeflationThreshold_C",(NEP,PetscReal*),(nep,deftol));
379: return(0);
380: }
382: static PetscErrorCode NEPSLPSetEPS_SLP(NEP nep,EPS eps)
383: {
385: NEP_SLP *ctx = (NEP_SLP*)nep->data;
388: PetscObjectReference((PetscObject)eps);
389: EPSDestroy(&ctx->eps);
390: ctx->eps = eps;
391: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->eps);
392: nep->state = NEP_STATE_INITIAL;
393: return(0);
394: }
396: /*@
397: NEPSLPSetEPS - Associate a linear eigensolver object (EPS) to the
398: nonlinear eigenvalue solver.
400: Collective on nep
402: Input Parameters:
403: + nep - nonlinear eigenvalue solver
404: - eps - the eigensolver object
406: Level: advanced
408: .seealso: NEPSLPGetEPS()
409: @*/
410: PetscErrorCode NEPSLPSetEPS(NEP nep,EPS eps)
411: {
418: PetscTryMethod(nep,"NEPSLPSetEPS_C",(NEP,EPS),(nep,eps));
419: return(0);
420: }
422: static PetscErrorCode NEPSLPGetEPS_SLP(NEP nep,EPS *eps)
423: {
425: NEP_SLP *ctx = (NEP_SLP*)nep->data;
428: if (!ctx->eps) {
429: EPSCreate(PetscObjectComm((PetscObject)nep),&ctx->eps);
430: PetscObjectIncrementTabLevel((PetscObject)ctx->eps,(PetscObject)nep,1);
431: EPSSetOptionsPrefix(ctx->eps,((PetscObject)nep)->prefix);
432: EPSAppendOptionsPrefix(ctx->eps,"nep_slp_");
433: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->eps);
434: PetscObjectSetOptions((PetscObject)ctx->eps,((PetscObject)nep)->options);
435: }
436: *eps = ctx->eps;
437: return(0);
438: }
440: /*@
441: NEPSLPGetEPS - Retrieve the linear eigensolver object (EPS) associated
442: to the nonlinear eigenvalue solver.
444: Not Collective
446: Input Parameter:
447: . nep - nonlinear eigenvalue solver
449: Output Parameter:
450: . eps - the eigensolver object
452: Level: advanced
454: .seealso: NEPSLPSetEPS()
455: @*/
456: PetscErrorCode NEPSLPGetEPS(NEP nep,EPS *eps)
457: {
463: PetscUseMethod(nep,"NEPSLPGetEPS_C",(NEP,EPS*),(nep,eps));
464: return(0);
465: }
467: static PetscErrorCode NEPSLPSetEPSLeft_SLP(NEP nep,EPS eps)
468: {
470: NEP_SLP *ctx = (NEP_SLP*)nep->data;
473: PetscObjectReference((PetscObject)eps);
474: EPSDestroy(&ctx->epsts);
475: ctx->epsts = eps;
476: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->epsts);
477: nep->state = NEP_STATE_INITIAL;
478: return(0);
479: }
481: /*@
482: NEPSLPSetEPSLeft - Associate a linear eigensolver object (EPS) to the
483: nonlinear eigenvalue solver, used to compute left eigenvectors in the
484: two-sided variant of SLP.
486: Collective on nep
488: Input Parameters:
489: + nep - nonlinear eigenvalue solver
490: - eps - the eigensolver object
492: Level: advanced
494: .seealso: NEPSLPGetEPSLeft(), NEPSetTwoSided()
495: @*/
496: PetscErrorCode NEPSLPSetEPSLeft(NEP nep,EPS eps)
497: {
504: PetscTryMethod(nep,"NEPSLPSetEPSLeft_C",(NEP,EPS),(nep,eps));
505: return(0);
506: }
508: static PetscErrorCode NEPSLPGetEPSLeft_SLP(NEP nep,EPS *eps)
509: {
511: NEP_SLP *ctx = (NEP_SLP*)nep->data;
514: if (!ctx->epsts) {
515: EPSCreate(PetscObjectComm((PetscObject)nep),&ctx->epsts);
516: PetscObjectIncrementTabLevel((PetscObject)ctx->epsts,(PetscObject)nep,1);
517: EPSSetOptionsPrefix(ctx->epsts,((PetscObject)nep)->prefix);
518: EPSAppendOptionsPrefix(ctx->epsts,"nep_slp_left_");
519: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->epsts);
520: PetscObjectSetOptions((PetscObject)ctx->epsts,((PetscObject)nep)->options);
521: }
522: *eps = ctx->epsts;
523: return(0);
524: }
526: /*@
527: NEPSLPGetEPSLeft - Retrieve the linear eigensolver object (EPS) associated
528: to the nonlinear eigenvalue solver, used to compute left eigenvectors in the
529: two-sided variant of SLP.
531: Not Collective
533: Input Parameter:
534: . nep - nonlinear eigenvalue solver
536: Output Parameter:
537: . eps - the eigensolver object
539: Level: advanced
541: .seealso: NEPSLPSetEPSLeft(), NEPSetTwoSided()
542: @*/
543: PetscErrorCode NEPSLPGetEPSLeft(NEP nep,EPS *eps)
544: {
550: PetscUseMethod(nep,"NEPSLPGetEPSLeft_C",(NEP,EPS*),(nep,eps));
551: return(0);
552: }
554: static PetscErrorCode NEPSLPSetKSP_SLP(NEP nep,KSP ksp)
555: {
557: NEP_SLP *ctx = (NEP_SLP*)nep->data;
560: PetscObjectReference((PetscObject)ksp);
561: KSPDestroy(&ctx->ksp);
562: ctx->ksp = ksp;
563: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);
564: nep->state = NEP_STATE_INITIAL;
565: return(0);
566: }
568: /*@
569: NEPSLPSetKSP - Associate a linear solver object (KSP) to the nonlinear
570: eigenvalue solver.
572: Collective on nep
574: Input Parameters:
575: + nep - eigenvalue solver
576: - ksp - the linear solver object
578: Level: advanced
580: .seealso: NEPSLPGetKSP()
581: @*/
582: PetscErrorCode NEPSLPSetKSP(NEP nep,KSP ksp)
583: {
590: PetscTryMethod(nep,"NEPSLPSetKSP_C",(NEP,KSP),(nep,ksp));
591: return(0);
592: }
594: static PetscErrorCode NEPSLPGetKSP_SLP(NEP nep,KSP *ksp)
595: {
597: NEP_SLP *ctx = (NEP_SLP*)nep->data;
600: if (!ctx->ksp) {
601: KSPCreate(PetscObjectComm((PetscObject)nep),&ctx->ksp);
602: PetscObjectIncrementTabLevel((PetscObject)ctx->ksp,(PetscObject)nep,1);
603: KSPSetOptionsPrefix(ctx->ksp,((PetscObject)nep)->prefix);
604: KSPAppendOptionsPrefix(ctx->ksp,"nep_slp_");
605: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);
606: PetscObjectSetOptions((PetscObject)ctx->ksp,((PetscObject)nep)->options);
607: KSPSetErrorIfNotConverged(ctx->ksp,PETSC_TRUE);
608: KSPSetTolerances(ctx->ksp,SlepcDefaultTol(nep->tol),PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
609: }
610: *ksp = ctx->ksp;
611: return(0);
612: }
614: /*@
615: NEPSLPGetKSP - Retrieve the linear solver object (KSP) associated with
616: the nonlinear eigenvalue solver.
618: Not Collective
620: Input Parameter:
621: . nep - nonlinear eigenvalue solver
623: Output Parameter:
624: . ksp - the linear solver object
626: Level: advanced
628: .seealso: NEPSLPSetKSP()
629: @*/
630: PetscErrorCode NEPSLPGetKSP(NEP nep,KSP *ksp)
631: {
637: PetscUseMethod(nep,"NEPSLPGetKSP_C",(NEP,KSP*),(nep,ksp));
638: return(0);
639: }
641: PetscErrorCode NEPView_SLP(NEP nep,PetscViewer viewer)
642: {
644: NEP_SLP *ctx = (NEP_SLP*)nep->data;
645: PetscBool isascii;
648: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
649: if (isascii) {
650: if (ctx->deftol) {
651: PetscViewerASCIIPrintf(viewer," deflation threshold: %g\n",(double)ctx->deftol);
652: }
653: if (!ctx->eps) { NEPSLPGetEPS(nep,&ctx->eps); }
654: PetscViewerASCIIPushTab(viewer);
655: EPSView(ctx->eps,viewer);
656: if (nep->twosided) {
657: if (!ctx->epsts) { NEPSLPGetEPSLeft(nep,&ctx->epsts); }
658: EPSView(ctx->epsts,viewer);
659: }
660: if (!ctx->ksp) { NEPSLPGetKSP(nep,&ctx->ksp); }
661: KSPView(ctx->ksp,viewer);
662: PetscViewerASCIIPopTab(viewer);
663: }
664: return(0);
665: }
667: PetscErrorCode NEPReset_SLP(NEP nep)
668: {
670: NEP_SLP *ctx = (NEP_SLP*)nep->data;
673: EPSReset(ctx->eps);
674: if (nep->twosided) { EPSReset(ctx->epsts); }
675: KSPReset(ctx->ksp);
676: return(0);
677: }
679: PetscErrorCode NEPDestroy_SLP(NEP nep)
680: {
682: NEP_SLP *ctx = (NEP_SLP*)nep->data;
685: KSPDestroy(&ctx->ksp);
686: EPSDestroy(&ctx->eps);
687: EPSDestroy(&ctx->epsts);
688: PetscFree(nep->data);
689: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetDeflationThreshold_C",NULL);
690: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetDeflationThreshold_C",NULL);
691: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPS_C",NULL);
692: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPS_C",NULL);
693: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPSLeft_C",NULL);
694: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPSLeft_C",NULL);
695: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetKSP_C",NULL);
696: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetKSP_C",NULL);
697: return(0);
698: }
700: SLEPC_EXTERN PetscErrorCode NEPCreate_SLP(NEP nep)
701: {
703: NEP_SLP *ctx;
706: PetscNewLog(nep,&ctx);
707: nep->data = (void*)ctx;
709: nep->useds = PETSC_TRUE;
710: ctx->deftol = PETSC_DEFAULT;
712: nep->ops->solve = NEPSolve_SLP;
713: nep->ops->setup = NEPSetUp_SLP;
714: nep->ops->setfromoptions = NEPSetFromOptions_SLP;
715: nep->ops->reset = NEPReset_SLP;
716: nep->ops->destroy = NEPDestroy_SLP;
717: nep->ops->view = NEPView_SLP;
718: nep->ops->computevectors = NEPComputeVectors_Schur;
720: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetDeflationThreshold_C",NEPSLPSetDeflationThreshold_SLP);
721: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetDeflationThreshold_C",NEPSLPGetDeflationThreshold_SLP);
722: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPS_C",NEPSLPSetEPS_SLP);
723: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPS_C",NEPSLPGetEPS_SLP);
724: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPSLeft_C",NEPSLPSetEPSLeft_SLP);
725: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPSLeft_C",NEPSLPGetEPSLeft_SLP);
726: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetKSP_C",NEPSLPSetKSP_SLP);
727: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetKSP_C",NEPSLPGetKSP_SLP);
728: return(0);
729: }