Actual source code: slp.c
slepc-3.14.2 2021-02-01
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, 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: Succesive 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),1,"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),1,"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,nep->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL/10.0:nep->tol/10.0,nep->max_it);
59: if (nep->twosided) {
60: nep->ops->solve = NEPSolve_SLP_Twosided;
61: nep->ops->computevectors = NULL;
62: if (!ctx->epsts) { NEPSLPGetEPSLeft(nep,&ctx->epsts); }
63: EPSGetST(ctx->epsts,&st);
64: PetscObjectTypeCompareAny((PetscObject)st,&flg,STSINVERT,STCAYLEY,"");
65: if (flg) SETERRQ(PetscObjectComm((PetscObject)nep),1,"SLP does not support spectral transformation");
66: EPSSetDimensions(ctx->epsts,1,PETSC_DECIDE,PETSC_DECIDE);
67: EPSSetWhichEigenpairs(ctx->epsts,EPS_LARGEST_MAGNITUDE);
68: EPSSetTolerances(ctx->epsts,nep->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL/10.0:nep->tol/10.0,nep->max_it);
69: } else {
70: nep->ops->solve = NEPSolve_SLP;
71: nep->ops->computevectors = NEPComputeVectors_Schur;
72: }
73: NEPAllocateSolution(nep,0);
74: return(0);
75: }
77: static PetscErrorCode MatMult_SLP(Mat M,Vec x,Vec y)
78: {
79: PetscErrorCode ierr;
80: NEP_SLP_MATSHELL *ctx;
83: MatShellGetContext(M,(void**)&ctx);
84: MatMult(ctx->extop->MJ,x,ctx->w);
85: NEPDeflationFunctionSolve(ctx->extop,ctx->w,y);
86: return(0);
87: }
89: static PetscErrorCode MatDestroy_SLP(Mat M)
90: {
91: PetscErrorCode ierr;
92: NEP_SLP_MATSHELL *ctx;
95: MatShellGetContext(M,(void**)&ctx);
96: VecDestroy(&ctx->w);
97: PetscFree(ctx);
98: return(0);
99: }
101: #if defined(PETSC_HAVE_CUDA)
102: static PetscErrorCode MatCreateVecs_SLP(Mat M,Vec *left,Vec *right)
103: {
104: PetscErrorCode ierr;
105: NEP_SLP_MATSHELL *ctx;
108: MatShellGetContext(M,(void**)&ctx);
109: if (right) {
110: VecDuplicate(ctx->w,right);
111: }
112: if (left) {
113: VecDuplicate(ctx->w,left);
114: }
115: return(0);
116: }
117: #endif
119: static PetscErrorCode NEPSLPSetUpLinearEP(NEP nep,NEP_EXT_OP extop,PetscScalar lambda,Vec u,PetscBool ini)
120: {
121: PetscErrorCode ierr;
122: NEP_SLP *slpctx = (NEP_SLP*)nep->data;
123: Mat Mshell;
124: PetscInt nloc,mloc;
125: NEP_SLP_MATSHELL *shellctx;
128: if (ini) {
129: /* Create mat shell */
130: PetscNew(&shellctx);
131: shellctx->extop = extop;
132: NEPDeflationCreateVec(extop,&shellctx->w);
133: MatGetLocalSize(nep->function,&mloc,&nloc);
134: nloc += extop->szd; mloc += extop->szd;
135: MatCreateShell(PetscObjectComm((PetscObject)nep),nloc,mloc,PETSC_DETERMINE,PETSC_DETERMINE,shellctx,&Mshell);
136: MatShellSetOperation(Mshell,MATOP_MULT,(void(*)(void))MatMult_SLP);
137: MatShellSetOperation(Mshell,MATOP_DESTROY,(void(*)(void))MatDestroy_SLP);
138: #if defined(PETSC_HAVE_CUDA)
139: MatShellSetOperation(Mshell,MATOP_CREATE_VECS,(void(*)(void))MatCreateVecs_SLP);
140: #endif
141: EPSSetOperators(slpctx->eps,Mshell,NULL);
142: MatDestroy(&Mshell);
143: }
144: NEPDeflationSolveSetUp(extop,lambda);
145: NEPDeflationComputeJacobian(extop,lambda,NULL);
146: EPSSetInitialSpace(slpctx->eps,1,&u);
147: return(0);
148: }
150: PetscErrorCode NEPSolve_SLP(NEP nep)
151: {
153: NEP_SLP *ctx = (NEP_SLP*)nep->data;
154: Mat F,H;
155: Vec uu,u,r;
156: PetscScalar sigma,lambda,mu,im,*Hp,*Ap;
157: PetscReal resnorm;
158: PetscInt nconv,ldh,ldds,i,j;
159: PetscBool skip=PETSC_FALSE,lock=PETSC_FALSE;
160: NEP_EXT_OP extop=NULL; /* Extended operator for deflation */
163: /* get initial approximation of eigenvalue and eigenvector */
164: NEPGetDefaultShift(nep,&sigma);
165: if (!nep->nini) {
166: BVSetRandomColumn(nep->V,0);
167: }
168: lambda = sigma;
169: if (!ctx->ksp) { NEPSLPGetKSP(nep,&ctx->ksp); }
170: NEPDeflationInitialize(nep,nep->V,ctx->ksp,PETSC_TRUE,nep->nev,&extop);
171: NEPDeflationCreateVec(extop,&u);
172: VecDuplicate(u,&r);
173: BVGetColumn(nep->V,0,&uu);
174: NEPDeflationCopyToExtendedVec(extop,uu,NULL,u,PETSC_FALSE);
175: BVRestoreColumn(nep->V,0,&uu);
177: /* Restart loop */
178: while (nep->reason == NEP_CONVERGED_ITERATING) {
179: nep->its++;
181: /* form residual, r = T(lambda)*u (used in convergence test only) */
182: NEPDeflationComputeFunction(extop,lambda,&F);
183: MatMult(F,u,r);
185: /* convergence test */
186: VecNorm(r,NORM_2,&resnorm);
187: (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
188: nep->eigr[nep->nconv] = lambda;
189: if (nep->errest[nep->nconv]<=nep->tol || nep->errest[nep->nconv]<=ctx->deftol) {
190: if (nep->errest[nep->nconv]<=ctx->deftol && !extop->ref && nep->nconv) {
191: NEPDeflationExtractEigenpair(extop,nep->nconv,u,lambda,nep->ds);
192: NEPDeflationSetRefine(extop,PETSC_TRUE);
193: MatMult(F,u,r);
194: VecNorm(r,NORM_2,&resnorm);
195: (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
196: if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
197: } else if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
198: }
200: if (lock) {
201: NEPDeflationSetRefine(extop,PETSC_FALSE);
202: nep->nconv = nep->nconv + 1;
203: skip = PETSC_TRUE;
204: lock = PETSC_FALSE;
205: NEPDeflationLocking(extop,u,lambda);
206: }
207: (*nep->stopping)(nep,nep->its,nep->max_it,nep->nconv,nep->nev,&nep->reason,nep->stoppingctx);
208: if (!skip || nep->reason>0) {
209: NEPMonitor(nep,nep->its,nep->nconv,nep->eigr,nep->eigi,nep->errest,(nep->reason>0)?nep->nconv:nep->nconv+1);
210: }
212: if (nep->reason == NEP_CONVERGED_ITERATING) {
213: if (!skip) {
214: /* evaluate T(lambda) and T'(lambda) */
215: NEPSLPSetUpLinearEP(nep,extop,lambda,u,nep->its==1?PETSC_TRUE:PETSC_FALSE);
216: /* compute new eigenvalue correction mu and eigenvector approximation u */
217: EPSSolve(ctx->eps);
218: EPSGetConverged(ctx->eps,&nconv);
219: if (!nconv) {
220: PetscInfo1(nep,"iter=%D, inner iteration failed, stopping solve\n",nep->its);
221: nep->reason = NEP_DIVERGED_LINEAR_SOLVE;
222: break;
223: }
224: EPSGetEigenpair(ctx->eps,0,&mu,&im,u,NULL);
225: mu = 1.0/mu;
226: if (PetscAbsScalar(im)>PETSC_MACHINE_EPSILON) SETERRQ(PetscObjectComm((PetscObject)nep),1,"Complex eigenvalue approximation - not implemented in real scalars");
227: } else {
228: nep->its--; /* do not count this as a full iteration */
229: /* use second eigenpair computed in previous iteration */
230: EPSGetConverged(ctx->eps,&nconv);
231: if (nconv>=2) {
232: EPSGetEigenpair(ctx->eps,1,&mu,&im,u,NULL);
233: mu = 1.0/mu;
234: } else {
235: NEPDeflationSetRandomVec(extop,u);
236: mu = lambda-sigma;
237: }
238: skip = PETSC_FALSE;
239: }
240: /* correct eigenvalue */
241: lambda = lambda - mu;
242: }
243: }
244: NEPDeflationGetInvariantPair(extop,NULL,&H);
245: MatGetSize(H,NULL,&ldh);
246: DSReset(nep->ds);
247: DSSetType(nep->ds,DSNHEP);
248: DSAllocate(nep->ds,PetscMax(nep->nconv,1));
249: DSGetLeadingDimension(nep->ds,&ldds);
250: MatDenseGetArray(H,&Hp);
251: DSGetArray(nep->ds,DS_MAT_A,&Ap);
252: for (j=0;j<nep->nconv;j++)
253: for (i=0;i<nep->nconv;i++) Ap[j*ldds+i] = Hp[j*ldh+i];
254: DSRestoreArray(nep->ds,DS_MAT_A,&Ap);
255: MatDenseRestoreArray(H,&Hp);
256: MatDestroy(&H);
257: DSSetDimensions(nep->ds,nep->nconv,0,0,nep->nconv);
258: DSSolve(nep->ds,nep->eigr,nep->eigi);
259: NEPDeflationReset(extop);
260: VecDestroy(&u);
261: VecDestroy(&r);
262: return(0);
263: }
265: PetscErrorCode NEPSetFromOptions_SLP(PetscOptionItems *PetscOptionsObject,NEP nep)
266: {
268: NEP_SLP *ctx = (NEP_SLP*)nep->data;
269: PetscBool flg;
270: PetscReal r;
273: PetscOptionsHead(PetscOptionsObject,"NEP SLP Options");
275: r = 0.0;
276: PetscOptionsReal("-nep_slp_deflation_threshold","Tolerance used as a threshold for including deflated eigenpairs","NEPSLPSetDeflationThreshold",ctx->deftol,&r,&flg);
277: if (flg) { NEPSLPSetDeflationThreshold(nep,r); }
279: PetscOptionsTail();
281: if (!ctx->eps) { NEPSLPGetEPS(nep,&ctx->eps); }
282: EPSSetFromOptions(ctx->eps);
283: if (nep->twosided) {
284: if (!ctx->epsts) { NEPSLPGetEPSLeft(nep,&ctx->epsts); }
285: EPSSetFromOptions(ctx->epsts);
286: }
287: if (!ctx->ksp) { NEPSLPGetKSP(nep,&ctx->ksp); }
288: KSPSetFromOptions(ctx->ksp);
289: return(0);
290: }
292: static PetscErrorCode NEPSLPSetDeflationThreshold_SLP(NEP nep,PetscReal deftol)
293: {
294: NEP_SLP *ctx = (NEP_SLP*)nep->data;
297: ctx->deftol = deftol;
298: return(0);
299: }
301: /*@
302: NEPSLPSetDeflationThreshold - Sets the threshold value used to switch between
303: deflated and non-deflated iteration.
305: Logically Collective on nep
307: Input Parameters:
308: + nep - nonlinear eigenvalue solver
309: - deftol - the threshold value
311: Options Database Keys:
312: . -nep_slp_deflation_threshold <deftol> - set the threshold
314: Notes:
315: Normally, the solver iterates on the extended problem in order to deflate
316: previously converged eigenpairs. If this threshold is set to a nonzero value,
317: then once the residual error is below this threshold the solver will
318: continue the iteration without deflation. The intention is to be able to
319: improve the current eigenpair further, despite having previous eigenpairs
320: with somewhat bad precision.
322: Level: advanced
324: .seealso: NEPSLPGetDeflationThreshold()
325: @*/
326: PetscErrorCode NEPSLPSetDeflationThreshold(NEP nep,PetscReal deftol)
327: {
333: PetscTryMethod(nep,"NEPSLPSetDeflationThreshold_C",(NEP,PetscReal),(nep,deftol));
334: return(0);
335: }
337: static PetscErrorCode NEPSLPGetDeflationThreshold_SLP(NEP nep,PetscReal *deftol)
338: {
339: NEP_SLP *ctx = (NEP_SLP*)nep->data;
342: *deftol = ctx->deftol;
343: return(0);
344: }
346: /*@
347: NEPSLPGetDeflationThreshold - Returns the threshold value that controls deflation.
349: Not Collective
351: Input Parameter:
352: . nep - nonlinear eigenvalue solver
354: Output Parameter:
355: . deftol - the threshold
357: Level: advanced
359: .seealso: NEPSLPSetDeflationThreshold()
360: @*/
361: PetscErrorCode NEPSLPGetDeflationThreshold(NEP nep,PetscReal *deftol)
362: {
368: PetscUseMethod(nep,"NEPSLPGetDeflationThreshold_C",(NEP,PetscReal*),(nep,deftol));
369: return(0);
370: }
372: static PetscErrorCode NEPSLPSetEPS_SLP(NEP nep,EPS eps)
373: {
375: NEP_SLP *ctx = (NEP_SLP*)nep->data;
378: PetscObjectReference((PetscObject)eps);
379: EPSDestroy(&ctx->eps);
380: ctx->eps = eps;
381: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->eps);
382: nep->state = NEP_STATE_INITIAL;
383: return(0);
384: }
386: /*@
387: NEPSLPSetEPS - Associate a linear eigensolver object (EPS) to the
388: nonlinear eigenvalue solver.
390: Collective on nep
392: Input Parameters:
393: + nep - nonlinear eigenvalue solver
394: - eps - the eigensolver object
396: Level: advanced
398: .seealso: NEPSLPGetEPS()
399: @*/
400: PetscErrorCode NEPSLPSetEPS(NEP nep,EPS eps)
401: {
408: PetscTryMethod(nep,"NEPSLPSetEPS_C",(NEP,EPS),(nep,eps));
409: return(0);
410: }
412: static PetscErrorCode NEPSLPGetEPS_SLP(NEP nep,EPS *eps)
413: {
415: NEP_SLP *ctx = (NEP_SLP*)nep->data;
418: if (!ctx->eps) {
419: EPSCreate(PetscObjectComm((PetscObject)nep),&ctx->eps);
420: PetscObjectIncrementTabLevel((PetscObject)ctx->eps,(PetscObject)nep,1);
421: EPSSetOptionsPrefix(ctx->eps,((PetscObject)nep)->prefix);
422: EPSAppendOptionsPrefix(ctx->eps,"nep_slp_");
423: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->eps);
424: PetscObjectSetOptions((PetscObject)ctx->eps,((PetscObject)nep)->options);
425: }
426: *eps = ctx->eps;
427: return(0);
428: }
430: /*@
431: NEPSLPGetEPS - Retrieve the linear eigensolver object (EPS) associated
432: to the nonlinear eigenvalue solver.
434: Not Collective
436: Input Parameter:
437: . nep - nonlinear eigenvalue solver
439: Output Parameter:
440: . eps - the eigensolver object
442: Level: advanced
444: .seealso: NEPSLPSetEPS()
445: @*/
446: PetscErrorCode NEPSLPGetEPS(NEP nep,EPS *eps)
447: {
453: PetscUseMethod(nep,"NEPSLPGetEPS_C",(NEP,EPS*),(nep,eps));
454: return(0);
455: }
457: static PetscErrorCode NEPSLPSetEPSLeft_SLP(NEP nep,EPS eps)
458: {
460: NEP_SLP *ctx = (NEP_SLP*)nep->data;
463: PetscObjectReference((PetscObject)eps);
464: EPSDestroy(&ctx->epsts);
465: ctx->epsts = eps;
466: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->epsts);
467: nep->state = NEP_STATE_INITIAL;
468: return(0);
469: }
471: /*@
472: NEPSLPSetEPSLeft - Associate a linear eigensolver object (EPS) to the
473: nonlinear eigenvalue solver, used to compute left eigenvectors in the
474: two-sided variant of SLP.
476: Collective on nep
478: Input Parameters:
479: + nep - nonlinear eigenvalue solver
480: - eps - the eigensolver object
482: Level: advanced
484: .seealso: NEPSLPGetEPSLeft(), NEPSetTwoSided()
485: @*/
486: PetscErrorCode NEPSLPSetEPSLeft(NEP nep,EPS eps)
487: {
494: PetscTryMethod(nep,"NEPSLPSetEPSLeft_C",(NEP,EPS),(nep,eps));
495: return(0);
496: }
498: static PetscErrorCode NEPSLPGetEPSLeft_SLP(NEP nep,EPS *eps)
499: {
501: NEP_SLP *ctx = (NEP_SLP*)nep->data;
504: if (!ctx->epsts) {
505: EPSCreate(PetscObjectComm((PetscObject)nep),&ctx->epsts);
506: PetscObjectIncrementTabLevel((PetscObject)ctx->epsts,(PetscObject)nep,1);
507: EPSSetOptionsPrefix(ctx->epsts,((PetscObject)nep)->prefix);
508: EPSAppendOptionsPrefix(ctx->epsts,"nep_slp_left_");
509: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->epsts);
510: PetscObjectSetOptions((PetscObject)ctx->epsts,((PetscObject)nep)->options);
511: }
512: *eps = ctx->epsts;
513: return(0);
514: }
516: /*@
517: NEPSLPGetEPSLeft - Retrieve the linear eigensolver object (EPS) associated
518: to the nonlinear eigenvalue solver, used to compute left eigenvectors in the
519: two-sided variant of SLP.
521: Not Collective
523: Input Parameter:
524: . nep - nonlinear eigenvalue solver
526: Output Parameter:
527: . eps - the eigensolver object
529: Level: advanced
531: .seealso: NEPSLPSetEPSLeft(), NEPSetTwoSided()
532: @*/
533: PetscErrorCode NEPSLPGetEPSLeft(NEP nep,EPS *eps)
534: {
540: PetscUseMethod(nep,"NEPSLPGetEPSLeft_C",(NEP,EPS*),(nep,eps));
541: return(0);
542: }
544: static PetscErrorCode NEPSLPSetKSP_SLP(NEP nep,KSP ksp)
545: {
547: NEP_SLP *ctx = (NEP_SLP*)nep->data;
550: PetscObjectReference((PetscObject)ksp);
551: KSPDestroy(&ctx->ksp);
552: ctx->ksp = ksp;
553: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);
554: nep->state = NEP_STATE_INITIAL;
555: return(0);
556: }
558: /*@
559: NEPSLPSetKSP - Associate a linear solver object (KSP) to the nonlinear
560: eigenvalue solver.
562: Collective on nep
564: Input Parameters:
565: + nep - eigenvalue solver
566: - ksp - the linear solver object
568: Level: advanced
570: .seealso: NEPSLPGetKSP()
571: @*/
572: PetscErrorCode NEPSLPSetKSP(NEP nep,KSP ksp)
573: {
580: PetscTryMethod(nep,"NEPSLPSetKSP_C",(NEP,KSP),(nep,ksp));
581: return(0);
582: }
584: static PetscErrorCode NEPSLPGetKSP_SLP(NEP nep,KSP *ksp)
585: {
587: NEP_SLP *ctx = (NEP_SLP*)nep->data;
590: if (!ctx->ksp) {
591: KSPCreate(PetscObjectComm((PetscObject)nep),&ctx->ksp);
592: PetscObjectIncrementTabLevel((PetscObject)ctx->ksp,(PetscObject)nep,1);
593: KSPSetOptionsPrefix(ctx->ksp,((PetscObject)nep)->prefix);
594: KSPAppendOptionsPrefix(ctx->ksp,"nep_slp_");
595: PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);
596: KSPSetErrorIfNotConverged(ctx->ksp,PETSC_TRUE);
597: KSPSetTolerances(ctx->ksp,SLEPC_DEFAULT_TOL,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
598: }
599: *ksp = ctx->ksp;
600: return(0);
601: }
603: /*@
604: NEPSLPGetKSP - Retrieve the linear solver object (KSP) associated with
605: the nonlinear eigenvalue solver.
607: Not Collective
609: Input Parameter:
610: . nep - nonlinear eigenvalue solver
612: Output Parameter:
613: . ksp - the linear solver object
615: Level: advanced
617: .seealso: NEPSLPSetKSP()
618: @*/
619: PetscErrorCode NEPSLPGetKSP(NEP nep,KSP *ksp)
620: {
626: PetscUseMethod(nep,"NEPSLPGetKSP_C",(NEP,KSP*),(nep,ksp));
627: return(0);
628: }
630: PetscErrorCode NEPView_SLP(NEP nep,PetscViewer viewer)
631: {
633: NEP_SLP *ctx = (NEP_SLP*)nep->data;
634: PetscBool isascii;
637: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
638: if (isascii) {
639: if (ctx->deftol) {
640: PetscViewerASCIIPrintf(viewer," deflation threshold: %g\n",(double)ctx->deftol);
641: }
642: if (!ctx->eps) { NEPSLPGetEPS(nep,&ctx->eps); }
643: PetscViewerASCIIPushTab(viewer);
644: EPSView(ctx->eps,viewer);
645: if (nep->twosided) {
646: if (!ctx->epsts) { NEPSLPGetEPSLeft(nep,&ctx->epsts); }
647: EPSView(ctx->epsts,viewer);
648: }
649: if (!ctx->ksp) { NEPSLPGetKSP(nep,&ctx->ksp); }
650: KSPView(ctx->ksp,viewer);
651: PetscViewerASCIIPopTab(viewer);
652: }
653: return(0);
654: }
656: PetscErrorCode NEPReset_SLP(NEP nep)
657: {
659: NEP_SLP *ctx = (NEP_SLP*)nep->data;
662: EPSReset(ctx->eps);
663: if (nep->twosided) { EPSReset(ctx->epsts); }
664: KSPReset(ctx->ksp);
665: return(0);
666: }
668: PetscErrorCode NEPDestroy_SLP(NEP nep)
669: {
671: NEP_SLP *ctx = (NEP_SLP*)nep->data;
674: KSPDestroy(&ctx->ksp);
675: EPSDestroy(&ctx->eps);
676: EPSDestroy(&ctx->epsts);
677: PetscFree(nep->data);
678: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetDeflationThreshold_C",NULL);
679: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetDeflationThreshold_C",NULL);
680: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPS_C",NULL);
681: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPS_C",NULL);
682: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPSLeft_C",NULL);
683: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPSLeft_C",NULL);
684: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetKSP_C",NULL);
685: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetKSP_C",NULL);
686: return(0);
687: }
689: SLEPC_EXTERN PetscErrorCode NEPCreate_SLP(NEP nep)
690: {
692: NEP_SLP *ctx;
695: PetscNewLog(nep,&ctx);
696: nep->data = (void*)ctx;
698: nep->useds = PETSC_TRUE;
699: ctx->deftol = 0.0;
701: nep->ops->solve = NEPSolve_SLP;
702: nep->ops->setup = NEPSetUp_SLP;
703: nep->ops->setfromoptions = NEPSetFromOptions_SLP;
704: nep->ops->reset = NEPReset_SLP;
705: nep->ops->destroy = NEPDestroy_SLP;
706: nep->ops->view = NEPView_SLP;
707: nep->ops->computevectors = NEPComputeVectors_Schur;
709: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetDeflationThreshold_C",NEPSLPSetDeflationThreshold_SLP);
710: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetDeflationThreshold_C",NEPSLPGetDeflationThreshold_SLP);
711: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPS_C",NEPSLPSetEPS_SLP);
712: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPS_C",NEPSLPGetEPS_SLP);
713: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPSLeft_C",NEPSLPSetEPSLeft_SLP);
714: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPSLeft_C",NEPSLPGetEPSLeft_SLP);
715: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetKSP_C",NEPSLPSetKSP_SLP);
716: PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetKSP_C",NEPSLPGetKSP_SLP);
717: return(0);
718: }