Actual source code: krylovschur.c

slepc-3.6.1 2015-09-03
Report Typos and Errors
  1: /*

  3:    SLEPc eigensolver: "krylovschur"

  5:    Method: Krylov-Schur

  7:    Algorithm:

  9:        Single-vector Krylov-Schur method for non-symmetric problems,
 10:        including harmonic extraction.

 12:    References:

 14:        [1] "Krylov-Schur Methods in SLEPc", SLEPc Technical Report STR-7,
 15:            available at http://slepc.upv.es.

 17:        [2] G.W. Stewart, "A Krylov-Schur Algorithm for Large Eigenproblems",
 18:            SIAM J. Matrix Anal. App. 23(3):601-614, 2001.

 20:        [3] "Practical Implementation of Harmonic Krylov-Schur", SLEPc Technical
 21:             Report STR-9, available at http://slepc.upv.es.

 23:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 24:    SLEPc - Scalable Library for Eigenvalue Problem Computations
 25:    Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain

 27:    This file is part of SLEPc.

 29:    SLEPc is free software: you can redistribute it and/or modify it under  the
 30:    terms of version 3 of the GNU Lesser General Public License as published by
 31:    the Free Software Foundation.

 33:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY
 34:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS
 35:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for
 36:    more details.

 38:    You  should have received a copy of the GNU Lesser General  Public  License
 39:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 40:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 41: */

 43: #include <slepc/private/epsimpl.h>                /*I "slepceps.h" I*/
 44:  #include krylovschur.h

 48: PetscErrorCode EPSGetArbitraryValues(EPS eps,PetscScalar *rr,PetscScalar *ri)
 49: {
 51:   PetscInt       i,newi,ld,n,l;
 52:   Vec            xr=eps->work[0],xi=eps->work[1];
 53:   PetscScalar    re,im,*Zr,*Zi,*X;

 56:   DSGetLeadingDimension(eps->ds,&ld);
 57:   DSGetDimensions(eps->ds,&n,NULL,&l,NULL,NULL);
 58:   for (i=l;i<n;i++) {
 59:     re = eps->eigr[i];
 60:     im = eps->eigi[i];
 61:     STBackTransform(eps->st,1,&re,&im);
 62:     newi = i;
 63:     DSVectors(eps->ds,DS_MAT_X,&newi,NULL);
 64:     DSGetArray(eps->ds,DS_MAT_X,&X);
 65:     Zr = X+i*ld;
 66:     if (newi==i+1) Zi = X+newi*ld;
 67:     else Zi = NULL;
 68:     EPSComputeRitzVector(eps,Zr,Zi,eps->V,xr,xi);
 69:     DSRestoreArray(eps->ds,DS_MAT_X,&X);
 70:     (*eps->arbitrary)(re,im,xr,xi,rr+i,ri+i,eps->arbitraryctx);
 71:   }
 72:   return(0);
 73: }

 77: PetscErrorCode EPSSetUp_KrylovSchur(EPS eps)
 78: {
 79:   PetscErrorCode    ierr;
 80:   PetscReal         eta;
 81:   BVOrthogType      otype;
 82:   BVOrthogBlockType obtype;
 83:   EPS_KRYLOVSCHUR   *ctx = (EPS_KRYLOVSCHUR*)eps->data;
 84:   enum { EPS_KS_DEFAULT,EPS_KS_SYMM,EPS_KS_SLICE,EPS_KS_INDEF } variant;

 87:   /* spectrum slicing requires special treatment of default values */
 88:   if (eps->which==EPS_ALL) {
 89:     EPSSetUp_KrylovSchur_Slice(eps);
 90:   } else {
 91:     EPSSetDimensions_Default(eps,eps->nev,&eps->ncv,&eps->mpd);
 92:     if (eps->ncv>eps->nev+eps->mpd) SETERRQ(PetscObjectComm((PetscObject)eps),1,"The value of ncv must not be larger than nev+mpd");
 93:     if (!eps->max_it) eps->max_it = PetscMax(100,2*eps->n/eps->ncv);
 94:     if (!eps->which) { EPSSetWhichEigenpairs_Default(eps); }
 95:   }
 96:   if (!ctx->lock && eps->mpd<eps->ncv) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Should not use mpd parameter in non-locking variant");

 98:   if (eps->isgeneralized && eps->ishermitian && !eps->ispositive && eps->arbitrary) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Arbitrary selection of eigenpairs not implemented for indefinite problems");
 99:   if (eps->ishermitian && eps->ispositive && (eps->which==EPS_LARGEST_IMAGINARY || eps->which==EPS_SMALLEST_IMAGINARY)) SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->which");

101:   if (!eps->extraction) {
102:     EPSSetExtraction(eps,EPS_RITZ);
103:   } else if (eps->extraction!=EPS_RITZ && eps->extraction!=EPS_HARMONIC)
104:     SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Unsupported extraction type");

106:   if (!ctx->keep) ctx->keep = 0.5;

108:   EPSAllocateSolution(eps,1);
109:   EPS_SetInnerProduct(eps);
110:   if (eps->arbitrary) {
111:     EPSSetWorkVecs(eps,2);
112:   } else if (eps->ishermitian && !eps->ispositive){
113:     EPSSetWorkVecs(eps,1);
114:   }

116:   /* dispatch solve method */
117:   if (eps->ishermitian) {
118:     if (eps->which==EPS_ALL) {
119:       if (eps->isgeneralized && !eps->ispositive) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Spectrum slicing not implemented for indefinite problems");
120:       else variant = EPS_KS_SLICE;
121:     } else if (eps->isgeneralized && !eps->ispositive) {
122:       variant = EPS_KS_INDEF;
123:     } else {
124:       switch (eps->extraction) {
125:         case EPS_RITZ:     variant = EPS_KS_SYMM; break;
126:         case EPS_HARMONIC: variant = EPS_KS_DEFAULT; break;
127:         default: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Unsupported extraction type");
128:       }
129:     }
130:   } else {
131:     switch (eps->extraction) {
132:       case EPS_RITZ:     variant = EPS_KS_DEFAULT; break;
133:       case EPS_HARMONIC: variant = EPS_KS_DEFAULT; break;
134:       default: SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Unsupported extraction type");
135:     }
136:   }
137:   switch (variant) {
138:     case EPS_KS_DEFAULT:
139:       eps->ops->solve = EPSSolve_KrylovSchur_Default;
140:       eps->ops->computevectors = EPSComputeVectors_Schur;
141:       DSSetType(eps->ds,DSNHEP);
142:       DSAllocate(eps->ds,eps->ncv+1);
143:       break;
144:     case EPS_KS_SYMM:
145:       eps->ops->solve = EPSSolve_KrylovSchur_Symm;
146:       eps->ops->computevectors = EPSComputeVectors_Hermitian;
147:       DSSetType(eps->ds,DSHEP);
148:       DSSetCompact(eps->ds,PETSC_TRUE);
149:       DSSetExtraRow(eps->ds,PETSC_TRUE);
150:       DSAllocate(eps->ds,eps->ncv+1);
151:       break;
152:     case EPS_KS_SLICE:
153:       eps->ops->solve = EPSSolve_KrylovSchur_Slice;
154:       eps->ops->computevectors = EPSComputeVectors_Slice;
155:       break;
156:     case EPS_KS_INDEF:
157:       eps->ops->solve = EPSSolve_KrylovSchur_Indefinite;
158:       eps->ops->computevectors = EPSComputeVectors_Indefinite;
159:       DSSetType(eps->ds,DSGHIEP);
160:       DSSetCompact(eps->ds,PETSC_TRUE);
161:       DSAllocate(eps->ds,eps->ncv+1);
162:       /* force reorthogonalization for pseudo-Lanczos */
163:       BVGetOrthogonalization(eps->V,&otype,NULL,&eta,&obtype);
164:       BVSetOrthogonalization(eps->V,otype,BV_ORTHOG_REFINE_ALWAYS,eta,obtype);
165:       break;
166:     default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Unexpected error");
167:   }
168:   return(0);
169: }

173: PetscErrorCode EPSSolve_KrylovSchur_Default(EPS eps)
174: {
175:   PetscErrorCode  ierr;
176:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
177:   PetscInt        i,j,*pj,k,l,nv,ld,nconv;
178:   Mat             U;
179:   PetscScalar     *S,*Q,*g;
180:   PetscReal       beta,gamma=1.0;
181:   PetscBool       breakdown,harmonic;

184:   DSGetLeadingDimension(eps->ds,&ld);
185:   harmonic = (eps->extraction==EPS_HARMONIC || eps->extraction==EPS_REFINED_HARMONIC)?PETSC_TRUE:PETSC_FALSE;
186:   if (harmonic) { PetscMalloc1(ld,&g); }
187:   if (eps->arbitrary) pj = &j;
188:   else pj = NULL;

190:   /* Get the starting Arnoldi vector */
191:   EPSGetStartVector(eps,0,NULL);
192:   l = 0;

194:   /* Restart loop */
195:   while (eps->reason == EPS_CONVERGED_ITERATING) {
196:     eps->its++;

198:     /* Compute an nv-step Arnoldi factorization */
199:     nv = PetscMin(eps->nconv+eps->mpd,eps->ncv);
200:     DSGetArray(eps->ds,DS_MAT_A,&S);
201:     EPSBasicArnoldi(eps,PETSC_FALSE,S,ld,eps->nconv+l,&nv,&beta,&breakdown);
202:     DSRestoreArray(eps->ds,DS_MAT_A,&S);
203:     DSSetDimensions(eps->ds,nv,0,eps->nconv,eps->nconv+l);
204:     if (l==0) {
205:       DSSetState(eps->ds,DS_STATE_INTERMEDIATE);
206:     } else {
207:       DSSetState(eps->ds,DS_STATE_RAW);
208:     }
209:     BVSetActiveColumns(eps->V,eps->nconv,nv);

211:     /* Compute translation of Krylov decomposition if harmonic extraction used */
212:     if (harmonic) {
213:       DSTranslateHarmonic(eps->ds,eps->target,beta,PETSC_FALSE,g,&gamma);
214:     }

216:     /* Solve projected problem */
217:     DSSolve(eps->ds,eps->eigr,eps->eigi);
218:     if (eps->arbitrary) {
219:       EPSGetArbitraryValues(eps,eps->rr,eps->ri);
220:       j=1;
221:     }
222:     DSSort(eps->ds,eps->eigr,eps->eigi,eps->rr,eps->ri,pj);

224:     /* Check convergence */
225:     EPSKrylovConvergence(eps,PETSC_FALSE,eps->nconv,nv-eps->nconv,beta,gamma,&k);
226:     if (eps->its >= eps->max_it) eps->reason = EPS_DIVERGED_ITS;
227:     if (k >= eps->nev) eps->reason = EPS_CONVERGED_TOL;
228:     nconv = k;

230:     /* Update l */
231:     if (eps->reason != EPS_CONVERGED_ITERATING || breakdown) l = 0;
232:     else {
233:       l = PetscMax(1,(PetscInt)((nv-k)*ctx->keep));
234: #if !defined(PETSC_USE_COMPLEX)
235:       DSGetArray(eps->ds,DS_MAT_A,&S);
236:       if (S[k+l+(k+l-1)*ld] != 0.0) {
237:         if (k+l<nv-1) l = l+1;
238:         else l = l-1;
239:       }
240:       DSRestoreArray(eps->ds,DS_MAT_A,&S);
241: #endif
242:     }
243:     if (!ctx->lock && l>0) { l += k; k = 0; } /* non-locking variant: reset no. of converged pairs */

245:     if (eps->reason == EPS_CONVERGED_ITERATING) {
246:       if (breakdown) {
247:         /* Start a new Arnoldi factorization */
248:         PetscInfo2(eps,"Breakdown in Krylov-Schur method (it=%D norm=%g)\n",eps->its,(double)beta);
249:         if (k<eps->nev) {
250:           EPSGetStartVector(eps,k,&breakdown);
251:           if (breakdown) {
252:             eps->reason = EPS_DIVERGED_BREAKDOWN;
253:             PetscInfo(eps,"Unable to generate more start vectors\n");
254:           }
255:         }
256:       } else {
257:         /* Undo translation of Krylov decomposition */
258:         if (harmonic) {
259:           DSSetDimensions(eps->ds,nv,0,k,l);
260:           DSTranslateHarmonic(eps->ds,0.0,beta,PETSC_TRUE,g,&gamma);
261:           /* gamma u^ = u - U*g~ */
262:           BVMultColumn(eps->V,-1.0,1.0,nv,g);
263:           BVScaleColumn(eps->V,nv,1.0/gamma);
264:         }
265:         /* Prepare the Rayleigh quotient for restart */
266:         DSGetArray(eps->ds,DS_MAT_A,&S);
267:         DSGetArray(eps->ds,DS_MAT_Q,&Q);
268:         for (i=k;i<k+l;i++) {
269:           S[k+l+i*ld] = Q[nv-1+i*ld]*beta*gamma;
270:         }
271:         DSRestoreArray(eps->ds,DS_MAT_A,&S);
272:         DSRestoreArray(eps->ds,DS_MAT_Q,&Q);
273:       }
274:     }
275:     /* Update the corresponding vectors V(:,idx) = V*Q(:,idx) */
276:     DSGetMat(eps->ds,DS_MAT_Q,&U);
277:     BVMultInPlace(eps->V,U,eps->nconv,k+l);
278:     MatDestroy(&U);

280:     if (eps->reason == EPS_CONVERGED_ITERATING && !breakdown) {
281:       BVCopyColumn(eps->V,nv,k+l);
282:     }
283:     eps->nconv = k;
284:     EPSMonitor(eps,eps->its,nconv,eps->eigr,eps->eigi,eps->errest,nv);
285:   }

287:   if (harmonic) { PetscFree(g); }
288:   /* truncate Schur decomposition and change the state to raw so that
289:      DSVectors() computes eigenvectors from scratch */
290:   DSSetDimensions(eps->ds,eps->nconv,0,0,0);
291:   DSSetState(eps->ds,DS_STATE_RAW);
292:   return(0);
293: }

297: static PetscErrorCode EPSKrylovSchurSetRestart_KrylovSchur(EPS eps,PetscReal keep)
298: {
299:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;

302:   if (keep==PETSC_DEFAULT) ctx->keep = 0.5;
303:   else {
304:     if (keep<0.1 || keep>0.9) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"The keep argument must be in the range [0.1,0.9]");
305:     ctx->keep = keep;
306:   }
307:   return(0);
308: }

312: /*@
313:    EPSKrylovSchurSetRestart - Sets the restart parameter for the Krylov-Schur
314:    method, in particular the proportion of basis vectors that must be kept
315:    after restart.

317:    Logically Collective on EPS

319:    Input Parameters:
320: +  eps - the eigenproblem solver context
321: -  keep - the number of vectors to be kept at restart

323:    Options Database Key:
324: .  -eps_krylovschur_restart - Sets the restart parameter

326:    Notes:
327:    Allowed values are in the range [0.1,0.9]. The default is 0.5.

329:    Level: advanced

331: .seealso: EPSKrylovSchurGetRestart()
332: @*/
333: PetscErrorCode EPSKrylovSchurSetRestart(EPS eps,PetscReal keep)
334: {

340:   PetscTryMethod(eps,"EPSKrylovSchurSetRestart_C",(EPS,PetscReal),(eps,keep));
341:   return(0);
342: }

346: static PetscErrorCode EPSKrylovSchurGetRestart_KrylovSchur(EPS eps,PetscReal *keep)
347: {
348:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;

351:   *keep = ctx->keep;
352:   return(0);
353: }

357: /*@
358:    EPSKrylovSchurGetRestart - Gets the restart parameter used in the
359:    Krylov-Schur method.

361:    Not Collective

363:    Input Parameter:
364: .  eps - the eigenproblem solver context

366:    Output Parameter:
367: .  keep - the restart parameter

369:    Level: advanced

371: .seealso: EPSKrylovSchurSetRestart()
372: @*/
373: PetscErrorCode EPSKrylovSchurGetRestart(EPS eps,PetscReal *keep)
374: {

380:   PetscTryMethod(eps,"EPSKrylovSchurGetRestart_C",(EPS,PetscReal*),(eps,keep));
381:   return(0);
382: }

386: static PetscErrorCode EPSKrylovSchurSetLocking_KrylovSchur(EPS eps,PetscBool lock)
387: {
388:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;

391:   ctx->lock = lock;
392:   return(0);
393: }

397: /*@
398:    EPSKrylovSchurSetLocking - Choose between locking and non-locking variants of
399:    the Krylov-Schur method.

401:    Logically Collective on EPS

403:    Input Parameters:
404: +  eps  - the eigenproblem solver context
405: -  lock - true if the locking variant must be selected

407:    Options Database Key:
408: .  -eps_krylovschur_locking - Sets the locking flag

410:    Notes:
411:    The default is to lock converged eigenpairs when the method restarts.
412:    This behaviour can be changed so that all directions are kept in the
413:    working subspace even if already converged to working accuracy (the
414:    non-locking variant).

416:    Level: advanced

418: .seealso: EPSKrylovSchurGetLocking()
419: @*/
420: PetscErrorCode EPSKrylovSchurSetLocking(EPS eps,PetscBool lock)
421: {

427:   PetscTryMethod(eps,"EPSKrylovSchurSetLocking_C",(EPS,PetscBool),(eps,lock));
428:   return(0);
429: }

433: static PetscErrorCode EPSKrylovSchurGetLocking_KrylovSchur(EPS eps,PetscBool *lock)
434: {
435:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;

438:   *lock = ctx->lock;
439:   return(0);
440: }

444: /*@
445:    EPSKrylovSchurGetLocking - Gets the locking flag used in the Krylov-Schur
446:    method.

448:    Not Collective

450:    Input Parameter:
451: .  eps - the eigenproblem solver context

453:    Output Parameter:
454: .  lock - the locking flag

456:    Level: advanced

458: .seealso: EPSKrylovSchurSetLocking()
459: @*/
460: PetscErrorCode EPSKrylovSchurGetLocking(EPS eps,PetscBool *lock)
461: {

467:   PetscTryMethod(eps,"EPSKrylovSchurGetLocking_C",(EPS,PetscBool*),(eps,lock));
468:   return(0);
469: }

473: static PetscErrorCode EPSKrylovSchurSetPartitions_KrylovSchur(EPS eps,PetscInt npart)
474: {
475:   PetscErrorCode  ierr;
476:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
477:   PetscMPIInt     size;

480:   if (ctx->npart!=npart) {
481:     if (ctx->commset) { PetscSubcommDestroy(&ctx->subc); }
482:     if (ctx->eps) { EPSDestroy(&ctx->eps); }
483:   } 
484:   if (npart == PETSC_DEFAULT || npart == PETSC_DECIDE) {
485:     ctx->npart = 1;
486:   } else {
487:     MPI_Comm_size(PetscObjectComm((PetscObject)eps),&size);
488:     if (npart<1 || npart>size) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of npart");
489:     ctx->npart = npart;
490:   }
491:   eps->state = EPS_STATE_INITIAL;
492:   return(0);
493: }

497: /*@
498:    EPSKrylovSchurSetPartitions - Sets the number of partitions for the
499:    case of doing spectrum slicing for a computational interval with the
500:    communicator split in several sub-communicators.

502:    Logically Collective on EPS

504:    Input Parameters:
505: +  eps   - the eigenproblem solver context
506: -  npart - number of partitions

508:    Options Database Key:
509: .  -eps_krylovschur_partitions <npart> - Sets the number of partitions

511:    Notes:
512:    By default, npart=1 so all processes in the communicator participate in
513:    the processing of the whole interval. If npart>1 then the interval is
514:    divided into npart subintervals, each of them being processed by a
515:    subset of processes.

517:    The interval is split proportionally unless the separation points are
518:    specified with EPSKrylovSchurSetSubintervals().

520:    Level: advanced

522: .seealso: EPSKrylovSchurSetSubintervals(), EPSSetInterval()
523: @*/
524: PetscErrorCode EPSKrylovSchurSetPartitions(EPS eps,PetscInt npart)
525: {

531:   PetscTryMethod(eps,"EPSKrylovSchurSetPartitions_C",(EPS,PetscInt),(eps,npart));
532:   return(0);
533: }

537: static PetscErrorCode EPSKrylovSchurGetPartitions_KrylovSchur(EPS eps,PetscInt *npart)
538: {
539:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;

542:   *npart  = ctx->npart;
543:   return(0);
544: }

548: /*@
549:    EPSKrylovSchurGetPartitions - Gets the number of partitions of the
550:    communicator in case of spectrum slicing.

552:    Not Collective

554:    Input Parameter:
555: .  eps - the eigenproblem solver context

557:    Output Parameter:
558: .  npart - number of partitions

560:    Level: advanced

562: .seealso: EPSKrylovSchurSetPartitions()
563: @*/
564: PetscErrorCode EPSKrylovSchurGetPartitions(EPS eps,PetscInt *npart)
565: {

571:   PetscTryMethod(eps,"EPSKrylovSchurGetPartitions_C",(EPS,PetscInt*),(eps,npart));
572:   return(0);
573: }

577: static PetscErrorCode EPSKrylovSchurSetDetectZeros_KrylovSchur(EPS eps,PetscBool detect)
578: {
579:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;

582:   ctx->detect = detect;
583:   eps->state  = EPS_STATE_INITIAL;
584:   return(0);
585: }

589: /*@
590:    EPSKrylovSchurSetDetectZeros - Sets a flag to enforce detection of
591:    zeros during the factorizations throughout the spectrum slicing computation.

593:    Logically Collective on EPS

595:    Input Parameters:
596: +  eps    - the eigenproblem solver context
597: -  detect - check for zeros

599:    Options Database Key:
600: .  -eps_krylovschur_detect_zeros - Check for zeros; this takes an optional
601:    bool value (0/1/no/yes/true/false)

603:    Notes:
604:    A zero in the factorization indicates that a shift coincides with an eigenvalue.

606:    This flag is turned off by default, and may be necessary in some cases,
607:    especially when several partitions are being used. This feature currently
608:    requires an external package for factorizations with support for zero
609:    detection, e.g. MUMPS.

611:    Level: advanced

613: .seealso: EPSKrylovSchurSetPartitions(), EPSSetInterval()
614: @*/
615: PetscErrorCode EPSKrylovSchurSetDetectZeros(EPS eps,PetscBool detect)
616: {

622:   PetscTryMethod(eps,"EPSKrylovSchurSetDetectZeros_C",(EPS,PetscBool),(eps,detect));
623:   return(0);
624: }

628: static PetscErrorCode EPSKrylovSchurGetDetectZeros_KrylovSchur(EPS eps,PetscBool *detect)
629: {
630:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;

633:   *detect = ctx->detect;
634:   return(0);
635: }

639: /*@
640:    EPSKrylovSchurGetDetectZeros - Gets the flag that enforces zero detection
641:    in spectrum slicing.

643:    Not Collective

645:    Input Parameter:
646: .  eps - the eigenproblem solver context

648:    Output Parameter:
649: .  detect - whether zeros detection is enforced during factorizations

651:    Level: advanced

653: .seealso: EPSKrylovSchurSetDetectZeros()
654: @*/
655: PetscErrorCode EPSKrylovSchurGetDetectZeros(EPS eps,PetscBool *detect)
656: {

662:   PetscTryMethod(eps,"EPSKrylovSchurGetDetectZeros_C",(EPS,PetscBool*),(eps,detect));
663:   return(0);
664: }

668: static PetscErrorCode EPSKrylovSchurSetDimensions_KrylovSchur(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)
669: {
670:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;

673:   if (nev<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
674:   ctx->nev = nev;
675:   if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
676:     ctx->ncv = 0;
677:   } else {
678:     if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
679:     ctx->ncv = ncv;
680:   }
681:   if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
682:     ctx->mpd = 0;
683:   } else {
684:     if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
685:     ctx->mpd = mpd;
686:   }
687:   eps->state = EPS_STATE_INITIAL;
688:   return(0);
689: }

693: /*@
694:    EPSKrylovSchurSetDimensions - Sets the dimensions used for each subsolve
695:    step in case of doing spectrum slicing for a computational interval.
696:    The meaning of the parameters is the same as in EPSSetDimensions().

698:    Logically Collective on EPS

700:    Input Parameters:
701: +  eps - the eigenproblem solver context
702: .  nev - number of eigenvalues to compute
703: .  ncv - the maximum dimension of the subspace to be used by the subsolve
704: -  mpd - the maximum dimension allowed for the projected problem

706:    Options Database Key:
707: +  -eps_krylovschur_nev <nev> - Sets the number of eigenvalues
708: .  -eps_krylovschur_ncv <ncv> - Sets the dimension of the subspace
709: -  -eps_krylovschur_mpd <mpd> - Sets the maximum projected dimension

711:    Level: advanced

713: .seealso: EPSKrylovSchurGetDimensions(), EPSSetDimensions(), EPSSetInterval()
714: @*/
715: PetscErrorCode EPSKrylovSchurSetDimensions(EPS eps,PetscInt nev,PetscInt ncv,PetscInt mpd)
716: {

724:   PetscTryMethod(eps,"EPSKrylovSchurSetDimensions_C",(EPS,PetscInt,PetscInt,PetscInt),(eps,nev,ncv,mpd));
725:   return(0);
726: }

730: static PetscErrorCode EPSKrylovSchurGetDimensions_KrylovSchur(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
731: {
732:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;

735:   if (nev) *nev = ctx->nev;
736:   if (ncv) *ncv = ctx->ncv;
737:   if (mpd) *mpd = ctx->mpd;
738:   return(0);
739: }

743: /*@
744:    EPSKrylovSchurGetDimensions - Gets the dimensions used for each subsolve
745:    step in case of doing spectrum slicing for a computational interval.

747:    Not Collective

749:    Input Parameter:
750: .  eps - the eigenproblem solver context

752:    Output Parameters:
753: +  nev - number of eigenvalues to compute
754: .  ncv - the maximum dimension of the subspace to be used by the subsolve
755: -  mpd - the maximum dimension allowed for the projected problem

757:    Level: advanced

759: .seealso: EPSKrylovSchurSetDimensions()
760: @*/
761: PetscErrorCode EPSKrylovSchurGetDimensions(EPS eps,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
762: {

767:   PetscTryMethod(eps,"EPSKrylovSchurGetDimensions_C",(EPS,PetscInt*,PetscInt*,PetscInt*),(eps,nev,ncv,mpd));
768:   return(0);
769: }

773: static PetscErrorCode EPSKrylovSchurSetSubintervals_KrylovSchur(EPS eps,PetscReal* subint)
774: {
775:   PetscErrorCode  ierr;
776:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
777:   PetscInt        i;

780:   if (subint[0]!=eps->inta || subint[ctx->npart]!=eps->intb) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"First and last values must match the endpoints of EPSSetInterval()");
781:   for (i=0;i<ctx->npart;i++) if (subint[i]>=subint[i+1]) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Array must contain values in ascending order");
782:   if (ctx->subintervals) { PetscFree(ctx->subintervals); }
783:   PetscMalloc1(ctx->npart+1,&ctx->subintervals);
784:   for (i=0;i<ctx->npart+1;i++) ctx->subintervals[i] = subint[i];
785:   ctx->subintset = PETSC_TRUE;
786:   eps->state = EPS_STATE_INITIAL;
787:   return(0);
788: }

792: /*@C
793:    EPSKrylovSchurSetSubintervals - Sets the points that delimit the
794:    subintervals to be used in spectrum slicing with several partitions.

796:    Logically Collective on EPS

798:    Input Parameters:
799: +  eps    - the eigenproblem solver context
800: -  subint - array of real values specifying subintervals

802:    Notes:
803:    This function must be called after EPSKrylovSchurSetPartitions(). For npart
804:    partitions, the argument subint must contain npart+1 real values sorted in
805:    ascending order: subint_0, subint_1, ..., subint_npart, where the first
806:    and last values must coincide with the interval endpoints set with
807:    EPSSetInterval().

809:    The subintervals are then defined by two consecutive points: [subint_0,subint_1],
810:    [subint_1,subint_2], and so on.

812:    Level: advanced

814: .seealso: EPSKrylovSchurSetPartitions(), EPSKrylovSchurGetSubintervals(), EPSSetInterval()
815: @*/
816: PetscErrorCode EPSKrylovSchurSetSubintervals(EPS eps,PetscReal *subint)
817: {

822:   PetscTryMethod(eps,"EPSKrylovSchurSetSubintervals_C",(EPS,PetscReal*),(eps,subint));
823:   return(0);
824: }

828: static PetscErrorCode EPSKrylovSchurGetSubintervals_KrylovSchur(EPS eps,PetscReal **subint)
829: {
830:   PetscErrorCode  ierr;
831:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
832:   PetscInt        i;

835:   if (!ctx->subintset) {
836:     if (!eps->state) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONGSTATE,"Must call EPSSetUp() first");
837:     if (!ctx->sr) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONGSTATE,"Only available in interval computations, see EPSSetInterval()");
838:   }
839:   PetscMalloc1(ctx->npart+1,subint);
840:   for (i=0;i<=ctx->npart;i++) (*subint)[i] = ctx->subintervals[i];
841:   return(0);
842: }

846: /*@C
847:    EPSKrylovSchurGetSubintervals - Returns the points that delimit the
848:    subintervals used in spectrum slicing with several partitions.

850:    Logically Collective on EPS

852:    Input Parameter:
853: .  eps    - the eigenproblem solver context

855:    Output Parameter:
856: .  subint - array of real values specifying subintervals

858:    Notes:
859:    If the user passed values with EPSKrylovSchurSetSubintervals(), then the
860:    same values are returned. Otherwise, the values computed internally are
861:    obtained.

863:    This function is only available for spectrum slicing runs.

865:    The returned array has length npart+1 (see EPSKrylovSchurGetPartitions())
866:    and should be freed by the user.

868:    Level: advanced

870: .seealso: EPSKrylovSchurSetSubintervals(), EPSKrylovSchurGetPartitions(), EPSSetInterval()
871: @*/
872: PetscErrorCode EPSKrylovSchurGetSubintervals(EPS eps,PetscReal** subint)
873: {

878:   PetscTryMethod(eps,"EPSKrylovSchurGetSubintervals_C",(EPS,PetscReal**),(eps,subint));
879:   return(0);
880: }

884: static PetscErrorCode EPSKrylovSchurGetInertias_KrylovSchur(EPS eps,PetscInt *n,PetscReal **shifts,PetscInt **inertias)
885: {
886:   PetscErrorCode  ierr;
887:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
888:   PetscInt        i;
889:   EPS_SR          sr = ctx->sr;

892:   if (!eps->state) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONGSTATE,"Must call EPSSetUp() first");
893:   if (!ctx->sr) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONGSTATE,"Only available in interval computations, see EPSSetInterval()");
894:   switch (eps->state) {
895:   case EPS_STATE_INITIAL:
896:     break;
897:   case EPS_STATE_SETUP:
898:     *n = ctx->npart+1;
899:     PetscMalloc1(*n,shifts);
900:     PetscMalloc1(*n,inertias);
901:     (*shifts)[0] = eps->inta;
902:     (*inertias)[0] = (sr->dir==1)?sr->inertia0:sr->inertia1;
903:     if (ctx->npart==1) {
904:       (*shifts)[1] = eps->intb;
905:       (*inertias)[1] = (sr->dir==1)?sr->inertia1:sr->inertia0;
906:     } else {
907:       for (i=1;i<*n;i++) {
908:         (*shifts)[i] = ctx->subintervals[i];
909:         (*inertias)[i] = (*inertias)[i-1]+ctx->nconv_loc[i-1];
910:       }
911:     }
912:     break;
913:   case EPS_STATE_SOLVED:
914:   case EPS_STATE_EIGENVECTORS:
915:     *n = ctx->nshifts;
916:     PetscMalloc1(*n,shifts);
917:     PetscMalloc1(*n,inertias);
918:     for (i=0;i<*n;i++) {
919:       (*shifts)[i] = ctx->shifts[i];
920:       (*inertias)[i] = ctx->inertias[i];
921:     }
922:     break;
923:   }
924:   return(0);
925: }

929: /*@C
930:    EPSKrylovSchurGetInertias - Gets the values of the shifts and their
931:    corresponding inertias in case of doing spectrum slicing for a
932:    computational interval.

934:    Not Collective

936:    Input Parameter:
937: .  eps - the eigenproblem solver context

939:    Output Parameters:
940: +  n        - number of shifts, including the endpoints of the interval
941: .  shifts   - the values of the shifts used internally in the solver
942: -  inertias - the values of the inertia in each shift

944:    Notes:
945:    If called after EPSSolve(), all shifts used internally by the solver are
946:    returned (including both endpoints and any intermediate ones). If called
947:    before EPSSolve() and after EPSSetUp() then only the information of the
948:    endpoints of subintervals is available.

950:    This function is only available for spectrum slicing runs.

952:    The returned arrays should be freed by the user.

954:    Level: advanced

956: .seealso: EPSSetInterval(), EPSKrylovSchurSetSubintervals()
957: @*/
958: PetscErrorCode EPSKrylovSchurGetInertias(EPS eps,PetscInt *n,PetscReal **shifts,PetscInt **inertias)
959: {

965:   PetscTryMethod(eps,"EPSKrylovSchurGetInertias_C",(EPS,PetscInt*,PetscReal**,PetscInt**),(eps,n,shifts,inertias));
966:   return(0);
967: }

971: static PetscErrorCode EPSKrylovSchurGetSubcommInfo_KrylovSchur(EPS eps,PetscInt *k,PetscInt *n,Vec *v)
972: {
973:   PetscErrorCode  ierr;
974:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
975:   EPS_SR          sr = ((EPS_KRYLOVSCHUR*)ctx->eps->data)->sr;

978:   if (!eps->state) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONGSTATE,"Must call EPSSetUp() first");
979:   if (!ctx->sr) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONGSTATE,"Only available in interval computations, see EPSSetInterval()");
980:   if (k) *k = (ctx->npart==1)? 0: ctx->subc->color;
981:   if (n) *n = sr->numEigs;
982:   if (v) {
983:     BVCreateVec(sr->V,v);
984:   }
985:   return(0);
986: }

990: /*@C
991:    EPSKrylovSchurGetSubcommInfo - Gets information related to the case of
992:    doing spectrum slicing for a computational interval with multiple
993:    communicators.

995:    Collective on the subcommunicator (if v is given)

997:    Input Parameter:
998: .  eps - the eigenproblem solver context

1000:    Output Parameters:
1001: +  k - number of the subinterval for the calling process
1002: .  n - number of eigenvalues found in the k-th subinterval
1003: -  v - a vector owned by processes in the subcommunicator with dimensions
1004:        compatible for locally computed eigenvectors (or NULL)

1006:    Notes:
1007:    This function is only available for spectrum slicing runs.

1009:    The returned Vec should be destroyed by the user.

1011:    Level: advanced

1013: .seealso: EPSSetInterval(), EPSKrylovSchurSetPartitions(), EPSKrylovSchurGetSubcommPairs()
1014: @*/
1015: PetscErrorCode EPSKrylovSchurGetSubcommInfo(EPS eps,PetscInt *k,PetscInt *n,Vec *v)
1016: {

1021:   PetscTryMethod(eps,"EPSKrylovSchurGetSubcommInfo_C",(EPS,PetscInt*,PetscInt*,Vec*),(eps,k,n,v));
1022:   return(0);
1023: }

1027: static PetscErrorCode EPSKrylovSchurGetSubcommPairs_KrylovSchur(EPS eps,PetscInt i,PetscScalar *eig,Vec v)
1028: {
1029:   PetscErrorCode  ierr;
1030:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
1031:   EPS_SR          sr = ((EPS_KRYLOVSCHUR*)ctx->eps->data)->sr;

1034:   EPSCheckSolved(eps,1);
1035:   if (!ctx->sr) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONGSTATE,"Only available in interval computations, see EPSSetInterval()");
1036:   if (i<0 || i>=sr->numEigs) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Argument 2 out of range");
1037:   if (eig) *eig = sr->eigr[sr->perm[i]];
1038:   BVCopyVec(sr->V,sr->perm[i],v);
1039:   return(0);
1040: }

1044: /*@C
1045:    EPSKrylovSchurGetSubcommPairs - Gets the i-th eigenpair stored
1046:    internally in the multi-communicator to which the calling process belongs.

1048:    Collective on the subcommunicator (if v is given)

1050:    Input Parameter:
1051: +  eps - the eigenproblem solver context
1052: -  i   - index of the solution

1054:    Output Parameters:
1055: +  eig - the eigenvalue
1056: -  v   - the eigenvector

1058:    Notes:
1059:    It is allowed to pass NULL for v if the eigenvector is not required.
1060:    Otherwise, the caller must provide a valid Vec objects, i.e.,
1061:    it must be created by the calling program with EPSKrylovSchurGetSubcommInfo().

1063:    The index i should be a value between 0 and n-1, where n is the number of
1064:    vectors in the local subinterval, see EPSKrylovSchurGetSubcommInfo().

1066:    Level: advanced

1068: .seealso: EPSSetInterval(), EPSKrylovSchurSetPartitions(), EPSKrylovSchurGetSubcommInfo()
1069: @*/
1070: PetscErrorCode EPSKrylovSchurGetSubcommPairs(EPS eps,PetscInt i,PetscScalar *eig,Vec v)
1071: {

1077:   PetscTryMethod(eps,"EPSKrylovSchurGetSubcommPairs_C",(EPS,PetscInt,PetscScalar*,Vec),(eps,i,eig,v));
1078:   return(0);
1079: }

1083: PetscErrorCode EPSSetFromOptions_KrylovSchur(PetscOptions *PetscOptionsObject,EPS eps)
1084: {
1085:   PetscErrorCode  ierr;
1086:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
1087:   PetscBool       flg,lock,b;
1088:   PetscReal       keep;
1089:   PetscInt        i,j,k;

1092:   PetscOptionsHead(PetscOptionsObject,"EPS Krylov-Schur Options");
1093:   PetscOptionsReal("-eps_krylovschur_restart","Proportion of vectors kept after restart","EPSKrylovSchurSetRestart",0.5,&keep,&flg);
1094:   if (flg) {
1095:     EPSKrylovSchurSetRestart(eps,keep);
1096:   }
1097:   PetscOptionsBool("-eps_krylovschur_locking","Choose between locking and non-locking variants","EPSKrylovSchurSetLocking",PETSC_TRUE,&lock,&flg);
1098:   if (flg) {
1099:     EPSKrylovSchurSetLocking(eps,lock);
1100:   }
1101:   i = ctx->npart;
1102:   PetscOptionsInt("-eps_krylovschur_partitions","Number of partitions of the communicator for spectrum slicing","EPSKrylovSchurSetPartitions",ctx->npart,&i,&flg);
1103:   if (flg) {
1104:     EPSKrylovSchurSetPartitions(eps,i);
1105:   }
1106:   b = ctx->detect;
1107:   PetscOptionsBool("-eps_krylovschur_detect_zeros","Check zeros during factorizations at subinterval boundaries","EPSKrylovSchurSetDetectZeros",ctx->detect,&b,&flg);
1108:   if (flg) {
1109:     EPSKrylovSchurSetDetectZeros(eps,b);
1110:   }
1111:   i = 1;
1112:   j = k = PETSC_DECIDE;
1113:   PetscOptionsInt("-eps_krylovschur_nev","Number of eigenvalues to compute in each subsolve (only for spectrum slicing)","EPSKrylovSchurSetDimensions",40,&i,NULL);
1114:   PetscOptionsInt("-eps_krylovschur_ncv","Number of basis vectors in each subsolve (only for spectrum slicing)","EPSKrylovSchurSetDimensions",80,&j,NULL);
1115:   PetscOptionsInt("-eps_krylovschur_mpd","Maximum dimension of projected problem in each subsolve (only for spectrum slicing)","EPSKrylovSchurSetDimensions",80,&k,NULL);
1116:   EPSKrylovSchurSetDimensions(eps,i,j,k);
1117:   PetscOptionsTail();
1118:   return(0);
1119: }

1123: PetscErrorCode EPSView_KrylovSchur(EPS eps,PetscViewer viewer)
1124: {
1125:   PetscErrorCode  ierr;
1126:   EPS_KRYLOVSCHUR *ctx = (EPS_KRYLOVSCHUR*)eps->data;
1127:   PetscBool       isascii;

1130:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
1131:   if (isascii) {
1132:     PetscViewerASCIIPrintf(viewer,"  Krylov-Schur: %d%% of basis vectors kept after restart\n",(int)(100*ctx->keep));
1133:     PetscViewerASCIIPrintf(viewer,"  Krylov-Schur: using the %slocking variant\n",ctx->lock?"":"non-");
1134:     if (eps->which==EPS_ALL) {
1135:       PetscViewerASCIIPrintf(viewer,"  Krylov-Schur: doing spectrum slicing with nev=%D, ncv=%D, mpd=%D\n",ctx->nev,ctx->ncv,ctx->mpd);
1136:       if (ctx->npart>1) {
1137:         PetscViewerASCIIPrintf(viewer,"  Krylov-Schur: multi-communicator spectrum slicing with %d partitions\n",ctx->npart);
1138:         if (ctx->detect) { PetscViewerASCIIPrintf(viewer,"  Krylov-Schur: detecting zeros when factorizing at subinterval boundaries\n"); }
1139:       }
1140:     }
1141:   }
1142:   return(0);
1143: }

1147: PetscErrorCode EPSDestroy_KrylovSchur(EPS eps)
1148: {

1152:   PetscFree(eps->data);
1153:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetRestart_C",NULL);
1154:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetRestart_C",NULL);
1155:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetLocking_C",NULL);
1156:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetLocking_C",NULL);
1157:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetPartitions_C",NULL);
1158:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetPartitions_C",NULL);
1159:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetDetectZeros_C",NULL);
1160:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetDetectZeros_C",NULL);
1161:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetDimensions_C",NULL);
1162:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetDimensions_C",NULL);
1163:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetSubintervals_C",NULL);
1164:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetSubintervals_C",NULL);
1165:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetInertias_C",NULL);
1166:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetSubcommInfo_C",NULL);
1167:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetSubcommPairs_C",NULL);
1168:   return(0);
1169: }

1173: PetscErrorCode EPSReset_KrylovSchur(EPS eps)
1174: {
1177:   if (eps->which==EPS_ALL) {
1178:     EPSReset_KrylovSchur_Slice(eps);
1179:   }
1180:   return(0);
1181: }

1185: PETSC_EXTERN PetscErrorCode EPSCreate_KrylovSchur(EPS eps)
1186: {
1187:   EPS_KRYLOVSCHUR *ctx;
1188:   PetscErrorCode  ierr;

1191:   PetscNewLog(eps,&ctx);
1192:   eps->data   = (void*)ctx;
1193:   ctx->lock   = PETSC_TRUE;
1194:   ctx->nev    = 1;
1195:   ctx->npart  = 1;
1196:   ctx->detect = PETSC_FALSE;
1197:   ctx->global = PETSC_TRUE;

1199:   eps->ops->setup          = EPSSetUp_KrylovSchur;
1200:   eps->ops->setfromoptions = EPSSetFromOptions_KrylovSchur;
1201:   eps->ops->destroy        = EPSDestroy_KrylovSchur;
1202:   eps->ops->reset          = EPSReset_KrylovSchur;
1203:   eps->ops->view           = EPSView_KrylovSchur;
1204:   eps->ops->backtransform  = EPSBackTransform_Default;
1205:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetRestart_C",EPSKrylovSchurSetRestart_KrylovSchur);
1206:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetRestart_C",EPSKrylovSchurGetRestart_KrylovSchur);
1207:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetLocking_C",EPSKrylovSchurSetLocking_KrylovSchur);
1208:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetLocking_C",EPSKrylovSchurGetLocking_KrylovSchur);
1209:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetPartitions_C",EPSKrylovSchurSetPartitions_KrylovSchur);
1210:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetPartitions_C",EPSKrylovSchurGetPartitions_KrylovSchur);
1211:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetDetectZeros_C",EPSKrylovSchurSetDetectZeros_KrylovSchur);
1212:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetDetectZeros_C",EPSKrylovSchurGetDetectZeros_KrylovSchur);
1213:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetDimensions_C",EPSKrylovSchurSetDimensions_KrylovSchur);
1214:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetDimensions_C",EPSKrylovSchurGetDimensions_KrylovSchur);
1215:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurSetSubintervals_C",EPSKrylovSchurSetSubintervals_KrylovSchur);
1216:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetSubintervals_C",EPSKrylovSchurGetSubintervals_KrylovSchur);
1217:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetInertias_C",EPSKrylovSchurGetInertias_KrylovSchur);
1218:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetSubcommInfo_C",EPSKrylovSchurGetSubcommInfo_KrylovSchur);
1219:   PetscObjectComposeFunction((PetscObject)eps,"EPSKrylovSchurGetSubcommPairs_C",EPSKrylovSchurGetSubcommPairs_KrylovSchur);
1220:   return(0);
1221: }