Actual source code: slp-twosided.c

slepc-3.17.2 2022-08-09
Report Typos and Errors
  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:    Two-sided variant of the NEPSLP solver.
 12: */

 14: #include <slepc/private/nepimpl.h>
 15: #include "slp.h"

 17: typedef struct _n_nep_def_ctx *NEP_NEDEF_CTX;

 19: struct _n_nep_def_ctx {
 20:   PetscInt    n;
 21:   PetscBool   ref;
 22:   PetscScalar *eig;
 23:   BV          V,W;
 24: };

 26: typedef struct {   /* context for two-sided solver */
 27:   Mat         Ft;
 28:   Mat         Jt;
 29:   Vec         w;
 30: } NEP_SLPTS_MATSHELL;

 32: typedef struct {   /* context for non-equivalence deflation */
 33:   NEP_NEDEF_CTX defctx;
 34:   Mat           F;
 35:   Mat           P;
 36:   Mat           J;
 37:   KSP           ksp;
 38:   PetscBool     isJ;
 39:   PetscScalar   lambda;
 40:   Vec           w[2];
 41: } NEP_NEDEF_MATSHELL;

 43: static PetscErrorCode MatMult_SLPTS_Right(Mat M,Vec x,Vec y)
 44: {
 45:   NEP_SLPTS_MATSHELL *ctx;

 47:   MatShellGetContext(M,&ctx);
 48:   MatMult(ctx->Jt,x,ctx->w);
 49:   MatSolve(ctx->Ft,ctx->w,y);
 50:   PetscFunctionReturn(0);
 51: }

 53: static PetscErrorCode MatMult_SLPTS_Left(Mat M,Vec x,Vec y)
 54: {
 55:   NEP_SLPTS_MATSHELL *ctx;

 57:   MatShellGetContext(M,&ctx);
 58:   MatMultTranspose(ctx->Jt,x,ctx->w);
 59:   MatSolveTranspose(ctx->Ft,ctx->w,y);
 60:   PetscFunctionReturn(0);
 61: }

 63: static PetscErrorCode MatDestroy_SLPTS(Mat M)
 64: {
 65:   NEP_SLPTS_MATSHELL *ctx;

 67:   MatShellGetContext(M,&ctx);
 68:   VecDestroy(&ctx->w);
 69:   PetscFree(ctx);
 70:   PetscFunctionReturn(0);
 71: }

 73: #if defined(PETSC_HAVE_CUDA)
 74: static PetscErrorCode MatCreateVecs_SLPTS(Mat M,Vec *left,Vec *right)
 75: {
 76:   NEP_SLPTS_MATSHELL *ctx;

 78:   MatShellGetContext(M,&ctx);
 79:   if (right) VecDuplicate(ctx->w,right);
 80:   if (left) VecDuplicate(ctx->w,left);
 81:   PetscFunctionReturn(0);
 82: }
 83: #endif

 85: static PetscErrorCode NEPSLPSetUpEPSMat(NEP nep,Mat F,Mat J,PetscBool left,Mat *M)
 86: {
 87:   Mat                Mshell;
 88:   PetscInt           nloc,mloc;
 89:   NEP_SLPTS_MATSHELL *shellctx;

 91:   /* Create mat shell */
 92:   PetscNew(&shellctx);
 93:   shellctx->Ft = F;
 94:   shellctx->Jt = J;
 95:   MatGetLocalSize(nep->function,&mloc,&nloc);
 96:   MatCreateShell(PetscObjectComm((PetscObject)nep),nloc,mloc,PETSC_DETERMINE,PETSC_DETERMINE,shellctx,&Mshell);
 97:   if (left) MatShellSetOperation(Mshell,MATOP_MULT,(void(*)(void))MatMult_SLPTS_Left);
 98:   else MatShellSetOperation(Mshell,MATOP_MULT,(void(*)(void))MatMult_SLPTS_Right);
 99:   MatShellSetOperation(Mshell,MATOP_DESTROY,(void(*)(void))MatDestroy_SLPTS);
100: #if defined(PETSC_HAVE_CUDA)
101:   MatShellSetOperation(Mshell,MATOP_CREATE_VECS,(void(*)(void))MatCreateVecs_SLPTS);
102: #endif
103:   *M = Mshell;
104:   MatCreateVecs(nep->function,&shellctx->w,NULL);
105:   PetscFunctionReturn(0);
106: }

108: /* Functions for deflation */
109: static PetscErrorCode NEPDeflationNEDestroy(NEP_NEDEF_CTX defctx)
110: {
111:   if (!defctx) PetscFunctionReturn(0);
112:   PetscFree(defctx->eig);
113:   PetscFree(defctx);
114:   PetscFunctionReturn(0);
115: }

117: static PetscErrorCode NEPDeflationNECreate(NEP nep,BV V,BV W,PetscInt sz,NEP_NEDEF_CTX *defctx)
118: {
119:   NEP_NEDEF_CTX  op;

121:   PetscNew(&op);
122:   *defctx = op;
123:   op->n   = 0;
124:   op->ref = PETSC_FALSE;
125:   PetscCalloc1(sz,&op->eig);
126:   PetscObjectStateIncrease((PetscObject)V);
127:   PetscObjectStateIncrease((PetscObject)W);
128:   op->V = V;
129:   op->W = W;
130:   PetscFunctionReturn(0);
131: }

133: static PetscErrorCode NEPDeflationNEComputeFunction(NEP nep,Mat M,PetscScalar lambda)
134: {
135:   NEP_NEDEF_MATSHELL *matctx;

137:   MatShellGetContext(M,&matctx);
138:   if (lambda==matctx->lambda) PetscFunctionReturn(0);
139:   NEPComputeFunction(nep,lambda,matctx->F,matctx->P);
140:   if (matctx->isJ) NEPComputeJacobian(nep,lambda,matctx->J);
141:   if (matctx->ksp) NEP_KSPSetOperators(matctx->ksp,matctx->F,matctx->P);
142:   matctx->lambda = lambda;
143:   PetscFunctionReturn(0);
144: }

146: static PetscErrorCode MatMult_NEPDeflationNE(Mat M,Vec x,Vec r)
147: {
148:   Vec                t,tt;
149:   PetscScalar        *h,*alpha,lambda,*eig;
150:   PetscInt           i,k;
151:   NEP_NEDEF_MATSHELL *matctx;

153:   MatShellGetContext(M,&matctx);
154:   if (matctx->defctx->n && !matctx->defctx->ref) {
155:     k = matctx->defctx->n;
156:     lambda = matctx->lambda;
157:     eig = matctx->defctx->eig;
158:     t = matctx->w[0];
159:     VecCopy(x,t);
160:     PetscMalloc2(k,&h,k,&alpha);
161:     for (i=0;i<k;i++) alpha[i] = (lambda-eig[i]-1.0)/(lambda-eig[i]);
162:     BVDotVec(matctx->defctx->V,t,h);
163:     for (i=0;i<k;i++) h[i] *= alpha[i];
164:     BVMultVec(matctx->defctx->W,-1.0,1.0,t,h);
165:     MatMult(matctx->isJ?matctx->J:matctx->F,t,r);
166:     if (matctx->isJ) {
167:       for (i=0;i<k;i++) h[i] *= (1.0/((lambda-eig[i])*(lambda-eig[i])))/alpha[i];
168:       tt = matctx->w[1];
169:       BVMultVec(matctx->defctx->W,-1.0,0.0,tt,h);
170:       MatMult(matctx->F,tt,t);
171:       VecAXPY(r,1.0,t);
172:     }
173:     PetscFree2(h,alpha);
174:   } else MatMult(matctx->isJ?matctx->J:matctx->F,x,r);
175:   PetscFunctionReturn(0);
176: }

178: static PetscErrorCode MatMultTranspose_NEPDeflationNE(Mat M,Vec x,Vec r)
179: {
180:   Vec                t,tt;
181:   PetscScalar        *h,*alphaC,lambda,*eig;
182:   PetscInt           i,k;
183:   NEP_NEDEF_MATSHELL *matctx;

185:   MatShellGetContext(M,&matctx);
186:   t    = matctx->w[0];
187:   VecCopy(x,t);
188:   if (matctx->defctx->n && !matctx->defctx->ref) {
189:     VecConjugate(t);
190:     k = matctx->defctx->n;
191:     lambda = matctx->lambda;
192:     eig = matctx->defctx->eig;
193:     PetscMalloc2(k,&h,k,&alphaC);
194:     for (i=0;i<k;i++) alphaC[i] = PetscConj((lambda-eig[i]-1.0)/(lambda-eig[i]));
195:     BVDotVec(matctx->defctx->W,t,h);
196:     for (i=0;i<k;i++) h[i] *= alphaC[i];
197:     BVMultVec(matctx->defctx->V,-1.0,1.0,t,h);
198:     VecConjugate(t);
199:     MatMultTranspose(matctx->isJ?matctx->J:matctx->F,t,r);
200:     if (matctx->isJ) {
201:       for (i=0;i<k;i++) h[i] *= PetscConj(1.0/((lambda-eig[i])*(lambda-eig[i])))/alphaC[i];
202:       tt = matctx->w[1];
203:       BVMultVec(matctx->defctx->V,-1.0,0.0,tt,h);
204:       VecConjugate(tt);
205:       MatMultTranspose(matctx->F,tt,t);
206:       VecAXPY(r,1.0,t);
207:     }
208:     PetscFree2(h,alphaC);
209:   } else MatMultTranspose(matctx->isJ?matctx->J:matctx->F,t,r);
210:   PetscFunctionReturn(0);
211: }

213: static PetscErrorCode MatSolve_NEPDeflationNE(Mat M,Vec b,Vec x)
214: {
215:   PetscScalar        *h,*alpha,lambda,*eig;
216:   PetscInt           i,k;
217:   NEP_NEDEF_MATSHELL *matctx;

219:   MatShellGetContext(M,&matctx);
220:   if (!matctx->ksp) {
221:     VecCopy(b,x);
222:     PetscFunctionReturn(0);
223:   }
224:   KSPSolve(matctx->ksp,b,x);
225:   if (matctx->defctx->n && !matctx->defctx->ref) {
226:     k = matctx->defctx->n;
227:     lambda = matctx->lambda;
228:     eig = matctx->defctx->eig;
229:     PetscMalloc2(k,&h,k,&alpha);
230:     BVDotVec(matctx->defctx->V,x,h);
231:     for (i=0;i<k;i++) alpha[i] = (lambda-eig[i]-1.0)/(lambda-eig[i]);
232:     for (i=0;i<k;i++) h[i] *= alpha[i]/(1.0-alpha[i]);
233:     BVMultVec(matctx->defctx->W,1.0,1.0,x,h);
234:     PetscFree2(h,alpha);
235:   }
236:   PetscFunctionReturn(0);
237: }

239: static PetscErrorCode MatSolveTranspose_NEPDeflationNE(Mat M,Vec b,Vec x)
240: {
241:   PetscScalar        *h,*alphaC,lambda,*eig;
242:   PetscInt           i,k;
243:   NEP_NEDEF_MATSHELL *matctx;

245:   MatShellGetContext(M,&matctx);
246:   if (!matctx->ksp) {
247:     VecCopy(b,x);
248:     PetscFunctionReturn(0);
249:   }
250:   KSPSolveTranspose(matctx->ksp,b,x);
251:   if (matctx->defctx->n && !matctx->defctx->ref) {
252:     VecConjugate(x);
253:     k = matctx->defctx->n;
254:     lambda = matctx->lambda;
255:     eig = matctx->defctx->eig;
256:     PetscMalloc2(k,&h,k,&alphaC);
257:     BVDotVec(matctx->defctx->W,x,h);
258:     for (i=0;i<k;i++) alphaC[i] = PetscConj((lambda-eig[i]-1.0)/(lambda-eig[i]));
259:     for (i=0;i<k;i++) h[i] *= alphaC[i]/(1.0-alphaC[i]);
260:     BVMultVec(matctx->defctx->V,1.0,1.0,x,h);
261:     PetscFree2(h,alphaC);
262:     VecConjugate(x);
263:   }
264:   PetscFunctionReturn(0);
265: }

267: static PetscErrorCode MatDestroy_NEPDeflationNE(Mat M)
268: {
269:   NEP_NEDEF_MATSHELL *matctx;

271:   MatShellGetContext(M,&matctx);
272:   VecDestroy(&matctx->w[0]);
273:   VecDestroy(&matctx->w[1]);
274:   PetscFree(matctx);
275:   PetscFunctionReturn(0);
276: }

278: static PetscErrorCode MatCreateVecs_NEPDeflationNE(Mat M,Vec *right,Vec *left)
279: {
280:   NEP_NEDEF_MATSHELL *matctx;

282:   MatShellGetContext(M,&matctx);
283:   MatCreateVecs(matctx->F,right,left);
284:   PetscFunctionReturn(0);
285: }

287: static PetscErrorCode NEPDeflationNEFunctionCreate(NEP_NEDEF_CTX defctx,NEP nep,Mat F,Mat P,Mat J,KSP ksp,PetscBool isJ,Mat *Mshell)
288: {
289:   NEP_NEDEF_MATSHELL *matctx;
290:   PetscInt           nloc,mloc;

292:   /* Create mat shell */
293:   PetscNew(&matctx);
294:   MatGetLocalSize(nep->function,&mloc,&nloc);
295:   MatCreateShell(PetscObjectComm((PetscObject)nep),nloc,mloc,PETSC_DETERMINE,PETSC_DETERMINE,matctx,Mshell);
296:   matctx->F   = F;
297:   matctx->P   = P;
298:   matctx->J   = J;
299:   matctx->isJ = isJ;
300:   matctx->ksp = ksp;
301:   matctx->defctx = defctx;
302:   matctx->lambda = PETSC_MAX_REAL;
303:   MatCreateVecs(F,&matctx->w[0],NULL);
304:   VecDuplicate(matctx->w[0],&matctx->w[1]);
305:   MatShellSetOperation(*Mshell,MATOP_MULT,(void(*)(void))MatMult_NEPDeflationNE);
306:   MatShellSetOperation(*Mshell,MATOP_MULT_TRANSPOSE,(void(*)(void))MatMultTranspose_NEPDeflationNE);
307:   MatShellSetOperation(*Mshell,MATOP_SOLVE,(void(*)(void))MatSolve_NEPDeflationNE);
308:   MatShellSetOperation(*Mshell,MATOP_SOLVE_TRANSPOSE,(void(*)(void))MatSolveTranspose_NEPDeflationNE);
309:   MatShellSetOperation(*Mshell,MATOP_DESTROY,(void(*)(void))MatDestroy_NEPDeflationNE);
310:   MatShellSetOperation(*Mshell,MATOP_CREATE_VECS,(void(*)(void))MatCreateVecs_NEPDeflationNE);
311:   PetscFunctionReturn(0);
312: }

314: static PetscErrorCode NEPDeflationNERecoverEigenvectors(NEP_NEDEF_CTX defctx,Vec u,Vec w,PetscScalar lambda)
315: {
316:   PetscScalar    *h,*alpha,*eig;
317:   PetscInt       i,k;

319:   if (w) VecConjugate(w);
320:   if (defctx->n && !defctx->ref) {
321:     eig = defctx->eig;
322:     k = defctx->n;
323:     PetscMalloc2(k,&h,k,&alpha);
324:     for (i=0;i<k;i++) alpha[i] = (lambda-eig[i]-1.0)/(lambda-eig[i]);
325:     BVDotVec(defctx->V,u,h);
326:     for (i=0;i<k;i++) h[i] *= alpha[i];
327:     BVMultVec(defctx->W,-1.0,1.0,u,h);
328:     VecNormalize(u,NULL);
329:     if (w) {
330:       BVDotVec(defctx->W,w,h);
331:       for (i=0;i<k;i++) alpha[i] = PetscConj((lambda-eig[i]-1.0)/(lambda-eig[i]));
332:       for (i=0;i<k;i++) h[i] *= alpha[i];
333:       BVMultVec(defctx->V,-1.0,1.0,w,h);
334:       VecNormalize(w,NULL);
335:     }
336:     PetscFree2(h,alpha);
337:   }
338:   PetscFunctionReturn(0);
339: }

341: static PetscErrorCode NEPDeflationNELocking(NEP_NEDEF_CTX defctx,Vec u,Vec w,PetscScalar lambda)
342: {
343:   PetscInt       n;

345:   n = defctx->n++;
346:   defctx->eig[n] = lambda;
347:   BVInsertVec(defctx->V,n,u);
348:   BVInsertVec(defctx->W,n,w);
349:   BVSetActiveColumns(defctx->V,0,defctx->n);
350:   BVSetActiveColumns(defctx->W,0,defctx->n);
351:   BVBiorthonormalizeColumn(defctx->V,defctx->W,n,NULL);
352:   PetscFunctionReturn(0);
353: }

355: static PetscErrorCode NEPDeflationNESetRefine(NEP_NEDEF_CTX defctx,PetscBool ref)
356: {
357:   defctx->ref = ref;
358:   PetscFunctionReturn(0);
359: }

361: PetscErrorCode NEPSolve_SLP_Twosided(NEP nep)
362: {
363:   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
364:   Mat            mF,mJ,M,Mt;
365:   Vec            u,r,t,w;
366:   BV             X,Y;
367:   PetscScalar    sigma,lambda,mu,im=0.0,mu2,im2;
368:   PetscReal      resnorm,resl;
369:   PetscInt       nconv,nconv2,i;
370:   PetscBool      skip=PETSC_FALSE,lock=PETSC_FALSE;
371:   NEP_NEDEF_CTX  defctx=NULL;    /* Extended operator for deflation */

373:   /* get initial approximation of eigenvalue and eigenvector */
374:   NEPGetDefaultShift(nep,&sigma);
375:   if (!nep->nini) BVSetRandomColumn(nep->V,0);
376:   BVSetRandomColumn(nep->W,0);
377:   lambda = sigma;
378:   if (!ctx->ksp) NEPSLPGetKSP(nep,&ctx->ksp);
379:   BVDuplicate(nep->V,&X);
380:   BVDuplicate(nep->W,&Y);
381:   NEPDeflationNECreate(nep,X,Y,nep->nev,&defctx);
382:   BVGetColumn(nep->V,0,&t);
383:   VecDuplicate(t,&u);
384:   VecDuplicate(t,&w);
385:   BVRestoreColumn(nep->V,0,&t);
386:   BVCopyVec(nep->V,0,u);
387:   BVCopyVec(nep->W,0,w);
388:   VecDuplicate(u,&r);
389:   NEPDeflationNEFunctionCreate(defctx,nep,nep->function,nep->function_pre?nep->function_pre:nep->function,NULL,ctx->ksp,PETSC_FALSE,&mF);
390:   NEPDeflationNEFunctionCreate(defctx,nep,nep->function,nep->function,nep->jacobian,NULL,PETSC_TRUE,&mJ);
391:   NEPSLPSetUpEPSMat(nep,mF,mJ,PETSC_FALSE,&M);
392:   NEPSLPSetUpEPSMat(nep,mF,mJ,PETSC_TRUE,&Mt);
393:   EPSSetOperators(ctx->eps,M,NULL);
394:   MatDestroy(&M);
395:   EPSSetOperators(ctx->epsts,Mt,NULL);
396:   MatDestroy(&Mt);

398:   /* Restart loop */
399:   while (nep->reason == NEP_CONVERGED_ITERATING) {
400:     nep->its++;

402:     /* form residual,  r = T(lambda)*u (used in convergence test only) */
403:     NEPDeflationNEComputeFunction(nep,mF,lambda);
404:     MatMultTranspose(mF,w,r);
405:     VecNorm(r,NORM_2,&resl);
406:     MatMult(mF,u,r);

408:     /* convergence test */
409:     VecNorm(r,NORM_2,&resnorm);
410:     resnorm = PetscMax(resnorm,resl);
411:     (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
412:     nep->eigr[nep->nconv] = lambda;
413:     if (nep->errest[nep->nconv]<=nep->tol || nep->errest[nep->nconv]<=ctx->deftol) {
414:       if (nep->errest[nep->nconv]<=ctx->deftol && !defctx->ref && nep->nconv) {
415:         NEPDeflationNERecoverEigenvectors(defctx,u,w,lambda);
416:         VecConjugate(w);
417:         NEPDeflationNESetRefine(defctx,PETSC_TRUE);
418:         MatMultTranspose(mF,w,r);
419:         VecNorm(r,NORM_2,&resl);
420:         MatMult(mF,u,r);
421:         VecNorm(r,NORM_2,&resnorm);
422:         resnorm = PetscMax(resnorm,resl);
423:         (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);
424:         if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
425:       } else if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
426:     }
427:     if (lock) {
428:       lock = PETSC_FALSE;
429:       skip = PETSC_TRUE;
430:       NEPDeflationNERecoverEigenvectors(defctx,u,w,lambda);
431:       NEPDeflationNELocking(defctx,u,w,lambda);
432:       NEPDeflationNESetRefine(defctx,PETSC_FALSE);
433:       BVInsertVec(nep->V,nep->nconv,u);
434:       BVInsertVec(nep->W,nep->nconv,w);
435:       VecConjugate(w);
436:       nep->nconv = nep->nconv + 1;
437:     }
438:     (*nep->stopping)(nep,nep->its,nep->max_it,nep->nconv,nep->nev,&nep->reason,nep->stoppingctx);
439:     if (!skip || nep->reason>0) NEPMonitor(nep,nep->its,nep->nconv,nep->eigr,nep->eigi,nep->errest,(nep->reason>0)?nep->nconv:nep->nconv+1);

441:     if (nep->reason == NEP_CONVERGED_ITERATING) {
442:       if (!skip) {
443:         /* evaluate T(lambda) and T'(lambda) */
444:         NEPDeflationNEComputeFunction(nep,mF,lambda);
445:         NEPDeflationNEComputeFunction(nep,mJ,lambda);
446:         EPSSetInitialSpace(ctx->eps,1,&u);
447:         EPSSetInitialSpace(ctx->epsts,1,&w);

449:         /* compute new eigenvalue correction mu and eigenvector approximation u */
450:         EPSSolve(ctx->eps);
451:         EPSSolve(ctx->epsts);
452:         EPSGetConverged(ctx->eps,&nconv);
453:         EPSGetConverged(ctx->epsts,&nconv2);
454:         if (!nconv||!nconv2) {
455:           PetscInfo(nep,"iter=%" PetscInt_FMT ", inner iteration failed, stopping solve\n",nep->its);
456:           nep->reason = NEP_DIVERGED_LINEAR_SOLVE;
457:           break;
458:         }
459:         EPSGetEigenpair(ctx->eps,0,&mu,&im,u,NULL);
460:         for (i=0;i<nconv2;i++) {
461:           EPSGetEigenpair(ctx->epsts,i,&mu2,&im2,w,NULL);
462:           if (SlepcAbsEigenvalue(mu-mu2,im-im2)/SlepcAbsEigenvalue(mu,im)<nep->tol*1000) break;
463:         }
464:         if (i==nconv2) {
465:           PetscInfo(nep,"iter=%" PetscInt_FMT ", inner iteration failed, stopping solve\n",nep->its);
466:           nep->reason = NEP_DIVERGED_LINEAR_SOLVE;
467:           break;
468:         }

470:         mu = 1.0/mu;
472:       } else {
473:         nep->its--;  /* do not count this as a full iteration */
474:         /* use second eigenpair computed in previous iteration */
475:         EPSGetConverged(ctx->eps,&nconv);
476:         if (nconv>=2 && nconv2>=2) {
477:           EPSGetEigenpair(ctx->eps,1,&mu,&im,u,NULL);
478:           EPSGetEigenpair(ctx->epsts,1,&mu2,&im2,w,NULL);
479:           mu = 1.0/mu;
480:         } else {
481:           BVSetRandomColumn(nep->V,nep->nconv);
482:           BVSetRandomColumn(nep->W,nep->nconv);
483:           BVCopyVec(nep->V,nep->nconv,u);
484:           BVCopyVec(nep->W,nep->nconv,w);
485:           mu = lambda-sigma;
486:         }
487:         skip = PETSC_FALSE;
488:       }
489:       /* correct eigenvalue */
490:       lambda = lambda - mu;
491:     }
492:   }
493:   VecDestroy(&u);
494:   VecDestroy(&w);
495:   VecDestroy(&r);
496:   MatDestroy(&mF);
497:   MatDestroy(&mJ);
498:   BVDestroy(&X);
499:   BVDestroy(&Y);
500:   NEPDeflationNEDestroy(defctx);
501:   PetscFunctionReturn(0);
502: }