Actual source code: arnoldi.c

  1: /*                       

  3:    SLEPc eigensolver: "arnoldi"

  5:    Method: Explicitly Restarted Arnoldi

  7:    Algorithm:

  9:        Arnoldi method with explicit restart and deflation.

 11:    References:

 13:        [1] "Arnoldi Methods in SLEPc", SLEPc Technical Report STR-4, 
 14:            available at http://www.grycap.upv.es/slepc.

 16:    Last update: Feb 2009

 18:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 19:    SLEPc - Scalable Library for Eigenvalue Problem Computations
 20:    Copyright (c) 2002-2010, Universidad Politecnica de Valencia, Spain

 22:    This file is part of SLEPc.
 23:       
 24:    SLEPc is free software: you can redistribute it and/or modify it under  the
 25:    terms of version 3 of the GNU Lesser General Public License as published by
 26:    the Free Software Foundation.

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

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

 38:  #include private/epsimpl.h
 39:  #include slepcblaslapack.h

 41: PetscErrorCode EPSSolve_ARNOLDI(EPS);

 43: typedef struct {
 44:   PetscTruth delayed;
 45: } EPS_ARNOLDI;

 49: PetscErrorCode EPSSetUp_ARNOLDI(EPS eps)
 50: {

 54:   if (eps->ncv) { /* ncv set */
 55:     if (eps->ncv<eps->nev) SETERRQ(1,"The value of ncv must be at least nev");
 56:   }
 57:   else if (eps->mpd) { /* mpd set */
 58:     eps->ncv = PetscMin(eps->n,eps->nev+eps->mpd);
 59:   }
 60:   else { /* neither set: defaults depend on nev being small or large */
 61:     if (eps->nev<500) eps->ncv = PetscMin(eps->n,PetscMax(2*eps->nev,eps->nev+15));
 62:     else { eps->mpd = 500; eps->ncv = PetscMin(eps->n,eps->nev+eps->mpd); }
 63:   }
 64:   if (!eps->mpd) eps->mpd = eps->ncv;
 65:   if (eps->ncv>eps->nev+eps->mpd) SETERRQ(1,"The value of ncv must not be larger than nev+mpd");
 66:   if (!eps->max_it) eps->max_it = PetscMax(100,2*eps->n/eps->ncv);
 67:   if (!eps->which) eps->which = EPS_LARGEST_MAGNITUDE;
 68:   if (eps->ishermitian && (eps->which==EPS_LARGEST_IMAGINARY || eps->which==EPS_SMALLEST_IMAGINARY))
 69:     SETERRQ(1,"Wrong value of eps->which");

 71:   if (!eps->extraction) {
 72:     EPSSetExtraction(eps,EPS_RITZ);
 73:   }

 75:   EPSAllocateSolution(eps);
 76:   PetscFree(eps->T);
 77:   PetscMalloc(eps->ncv*eps->ncv*sizeof(PetscScalar),&eps->T);
 78:   if (eps->leftvecs) {
 79:     PetscFree(eps->Tl);
 80:     PetscMalloc(eps->ncv*eps->ncv*sizeof(PetscScalar),&eps->Tl);
 81:     PetscInfo(eps,"Warning: parameter mpd ignored\n");
 82:     EPSDefaultGetWork(eps,2);
 83:   } else {
 84:     EPSDefaultGetWork(eps,1);
 85:   }

 87:   /* dispatch solve method */
 88:   if (eps->leftvecs) SETERRQ(PETSC_ERR_SUP,"Left vectors not supported in this solver");
 89:   eps->ops->solve = EPSSolve_ARNOLDI;
 90:   return(0);
 91: }

 95: /*
 96:    EPSDelayedArnoldi - This function is equivalent to EPSBasicArnoldi but
 97:    performs the computation in a different way. The main idea is that
 98:    reorthogonalization is delayed to the next Arnoldi step. This version is
 99:    more scalable but in some cases convergence may stagnate.
100: */
101: PetscErrorCode EPSDelayedArnoldi(EPS eps,PetscScalar *H,PetscInt ldh,Vec *V,PetscInt k,PetscInt *M,Vec f,PetscReal *beta,PetscTruth *breakdown)
102: {
104:   PetscInt       i,j,m=*M;
105:   Vec            u,t;
106:   PetscScalar    shh[100],*lhh,dot,dot2;
107:   PetscReal      norm1=0.0,norm2;

110:   if (m<=100) lhh = shh;
111:   else { PetscMalloc(m*sizeof(PetscScalar),&lhh); }
112:   VecDuplicate(f,&u);
113:   VecDuplicate(f,&t);

115:   for (j=k;j<m;j++) {
116:     STApply(eps->OP,V[j],f);
117:     IPOrthogonalize(eps->ip,0,PETSC_NULL,eps->nds,PETSC_NULL,eps->DS,f,PETSC_NULL,PETSC_NULL,PETSC_NULL);

119:     IPMInnerProductBegin(eps->ip,f,j+1,V,H+ldh*j);
120:     if (j>k) {
121:       IPMInnerProductBegin(eps->ip,V[j],j,V,lhh);
122:       IPInnerProductBegin(eps->ip,V[j],V[j],&dot);
123:     }
124:     if (j>k+1) {
125:       IPNormBegin(eps->ip,u,&norm2);
126:       VecDotBegin(u,V[j-2],&dot2);
127:     }
128: 
129:     IPMInnerProductEnd(eps->ip,f,j+1,V,H+ldh*j);
130:     if (j>k) {
131:       IPMInnerProductEnd(eps->ip,V[j],j,V,lhh);
132:       IPInnerProductEnd(eps->ip,V[j],V[j],&dot);
133:     }
134:     if (j>k+1) {
135:       IPNormEnd(eps->ip,u,&norm2);
136:       VecDotEnd(u,V[j-2],&dot2);
137:       if (PetscAbsScalar(dot2/norm2) > PETSC_MACHINE_EPSILON) {
138:         *breakdown = PETSC_TRUE;
139:         *M = j-1;
140:         *beta = norm2;

142:         if (m>100) { PetscFree(lhh); }
143:         VecDestroy(u);
144:         VecDestroy(t);
145:         return(0);
146:       }
147:     }
148: 
149:     if (j>k) {
150:       norm1 = sqrt(PetscRealPart(dot));
151:       for (i=0;i<j;i++)
152:         H[ldh*j+i] = H[ldh*j+i]/norm1;
153:       H[ldh*j+j] = H[ldh*j+j]/dot;
154: 
155:       VecCopy(V[j],t);
156:       VecScale(V[j],1.0/norm1);
157:       VecScale(f,1.0/norm1);
158:     }

160:     SlepcVecMAXPBY(f,1.0,-1.0,j+1,H+ldh*j,V);

162:     if (j>k) {
163:       SlepcVecMAXPBY(t,1.0,-1.0,j,lhh,V);
164:       for (i=0;i<j;i++)
165:         H[ldh*(j-1)+i] += lhh[i];
166:     }

168:     if (j>k+1) {
169:       VecCopy(u,V[j-1]);
170:       VecScale(V[j-1],1.0/norm2);
171:       H[ldh*(j-2)+j-1] = norm2;
172:     }

174:     if (j<m-1) {
175:       VecCopy(f,V[j+1]);
176:       VecCopy(t,u);
177:     }
178:   }

180:   IPNorm(eps->ip,t,&norm2);
181:   VecScale(t,1.0/norm2);
182:   VecCopy(t,V[m-1]);
183:   H[ldh*(m-2)+m-1] = norm2;

185:   IPMInnerProduct(eps->ip,f,m,V,lhh);
186: 
187:   SlepcVecMAXPBY(f,1.0,-1.0,m,lhh,V);
188:   for (i=0;i<m;i++)
189:     H[ldh*(m-1)+i] += lhh[i];

191:   IPNorm(eps->ip,f,beta);
192:   VecScale(f,1.0 / *beta);
193:   *breakdown = PETSC_FALSE;
194: 
195:   if (m>100) { PetscFree(lhh); }
196:   VecDestroy(u);
197:   VecDestroy(t);

199:   return(0);
200: }

204: /*
205:    EPSDelayedArnoldi1 - This function is similar to EPSDelayedArnoldi1,
206:    but without reorthogonalization (only delayed normalization).
207: */
208: PetscErrorCode EPSDelayedArnoldi1(EPS eps,PetscScalar *H,PetscInt ldh,Vec *V,PetscInt k,PetscInt *M,Vec f,PetscReal *beta,PetscTruth *breakdown)
209: {
211:   PetscInt       i,j,m=*M;
212:   PetscScalar    dot;
213:   PetscReal      norm=0.0;


217:   for (j=k;j<m;j++) {
218:     STApply(eps->OP,V[j],f);
219:     IPOrthogonalize(eps->ip,0,PETSC_NULL,eps->nds,PETSC_NULL,eps->DS,f,PETSC_NULL,PETSC_NULL,PETSC_NULL);

221:     IPMInnerProductBegin(eps->ip,f,j+1,V,H+ldh*j);
222:     if (j>k) {
223:       IPInnerProductBegin(eps->ip,V[j],V[j],&dot);
224:     }
225: 
226:     IPMInnerProductEnd(eps->ip,f,j+1,V,H+ldh*j);
227:     if (j>k) {
228:       IPInnerProductEnd(eps->ip,V[j],V[j],&dot);
229:     }
230: 
231:     if (j>k) {
232:       norm = sqrt(PetscRealPart(dot));
233:       VecScale(V[j],1.0/norm);
234:       H[ldh*(j-1)+j] = norm;

236:       for (i=0;i<j;i++)
237:         H[ldh*j+i] = H[ldh*j+i]/norm;
238:       H[ldh*j+j] = H[ldh*j+j]/dot;
239:       VecScale(f,1.0/norm);
240:     }

242:     SlepcVecMAXPBY(f,1.0,-1.0,j+1,H+ldh*j,V);

244:     if (j<m-1) {
245:       VecCopy(f,V[j+1]);
246:     }
247:   }

249:   IPNorm(eps->ip,f,beta);
250:   VecScale(f,1.0 / *beta);
251:   *breakdown = PETSC_FALSE;
252: 
253:   return(0);
254: }

258: /*
259:    EPSProjectedArnoldi - Solves the projected eigenproblem.

261:    On input:
262:      S is the projected matrix (leading dimension is lds)

264:    On output:
265:      S has (real) Schur form with diagonal blocks sorted appropriately
266:      Q contains the corresponding Schur vectors (order n, leading dimension n)
267: */
268: PetscErrorCode EPSProjectedArnoldi(EPS eps,PetscScalar *S,PetscInt lds,PetscScalar *Q,PetscInt n)
269: {
271:   PetscInt       i;

274:   /* Initialize orthogonal matrix */
275:   PetscMemzero(Q,n*n*sizeof(PetscScalar));
276:   for (i=0;i<n;i++)
277:     Q[i*(n+1)] = 1.0;
278:   /* Reduce S to (quasi-)triangular form, S <- Q S Q' */
279:   EPSDenseSchur(n,eps->nconv,S,lds,Q,eps->eigr,eps->eigi);
280:   /* Sort the remaining columns of the Schur form */
281:   EPSSortDenseSchur(eps,n,eps->nconv,S,lds,Q,eps->eigr,eps->eigi);
282:   return(0);
283: }

287: /*
288:    EPSUpdateVectors - Computes approximate Schur vectors (or eigenvectors) by
289:    either Ritz extraction (U=U*Q) or refined Ritz extraction 

291:    On input:
292:      n is the size of U
293:      U is the orthogonal basis of the subspace used for projecting
294:      s is the index of the first vector computed
295:      e+1 is the index of the last vector computed
296:      Q contains the corresponding Schur vectors of the projected matrix (size n x n, leading dimension ldq)
297:      H is the (extended) projected matrix (size n+1 x n, leading dimension ldh)

299:    On output:
300:      v is the resulting vector
301: */
302: PetscErrorCode EPSUpdateVectors(EPS eps,PetscInt n_,Vec *U,PetscInt s,PetscInt e,PetscScalar *Q,PetscInt ldq,PetscScalar *H,PetscInt ldh_)
303: {
304: #if defined(PETSC_MISSING_LAPACK_GESVD) 
305:   SETERRQ(PETSC_ERR_SUP,"GESVD - Lapack routine is unavailable.");
306: #else
308:   PetscTruth     isrefined;
309:   PetscInt       i,j,k;
310:   PetscBLASInt   n1,lwork,idummy=1,info,n=n_,ldh=ldh_;
311:   PetscScalar    *B,sdummy,*work;
312:   PetscReal      *sigma;

315:   isrefined = (eps->extraction==EPS_REFINED || eps->extraction==EPS_REFINED_HARMONIC)?PETSC_TRUE:PETSC_FALSE;
316:   if (isrefined) {
317:     /* Refined Ritz extraction */
318:     n1 = n+1;
319:     PetscMalloc(n1*n*sizeof(PetscScalar),&B);
320:     PetscMalloc(6*n*sizeof(PetscReal),&sigma);
321:     lwork = 10*n;
322:     PetscMalloc(lwork*sizeof(PetscScalar),&work);
323: 
324:     for (k=s;k<e;k++) {
325:       /* copy H to B */
326:       for (i=0;i<=n;i++) {
327:         for (j=0;j<n;j++) {
328:           B[i+j*n1] = H[i+j*ldh];
329:         }
330:       }
331:       /* subtract ritz value from diagonal of B^ */
332:       for (i=0;i<n;i++) {
333:         B[i+i*n1] -= eps->eigr[k];  /* MISSING: complex case */
334:       }
335:       /* compute SVD of [H-mu*I] */
336:   #if !defined(PETSC_USE_COMPLEX)
337:       LAPACKgesvd_("N","O",&n1,&n,B,&n1,sigma,&sdummy,&idummy,&sdummy,&idummy,work,&lwork,&info);
338:   #else
339:       LAPACKgesvd_("N","O",&n1,&n,B,&n1,sigma,&sdummy,&idummy,&sdummy,&idummy,work,&lwork,sigma+n,&info);
340:   #endif
341:       if (info) SETERRQ1(PETSC_ERR_LIB,"Error in Lapack xGESVD %d",info);
342:       /* the smallest singular value is the new error estimate */
343:       eps->errest[k] = sigma[n-1];
344:       /* update vector with right singular vector associated to smallest singular value */
345:       for (i=0;i<n;i++)
346:         Q[k*ldq+i] = B[n-1+i*n1];
347:     }
348:     /* free workspace */
349:     PetscFree(B);
350:     PetscFree(sigma);
351:     PetscFree(work);
352:   }
353:   /* Ritz extraction: v = U*q */
354:   SlepcUpdateVectors(n_,U,s,e,Q,ldq,PETSC_FALSE);
355:   return(0);
356: #endif
357: }

361: PetscErrorCode EPSSolve_ARNOLDI(EPS eps)
362: {
364:   PetscInt       i,k,lwork,nv;
365:   Vec            f=eps->work[0];
366:   PetscScalar    *H=eps->T,*U,*g,*work,*Hcopy;
367:   PetscReal      beta,gnorm,corrf=1.0;
368:   PetscTruth     breakdown;
369:   IPOrthogonalizationRefinementType orthog_ref;
370:   EPS_ARNOLDI    *arnoldi = (EPS_ARNOLDI *)eps->data;

373:   PetscMemzero(eps->T,eps->ncv*eps->ncv*sizeof(PetscScalar));
374:   PetscMalloc(eps->ncv*eps->ncv*sizeof(PetscScalar),&U);
375:   lwork = PetscMax((eps->ncv+1)*eps->ncv,7*eps->ncv);
376:   PetscMalloc(lwork*sizeof(PetscScalar),&work);
377:   if (eps->extraction==EPS_HARMONIC || eps->extraction==EPS_REFINED_HARMONIC) {
378:     PetscMalloc(eps->ncv*sizeof(PetscScalar),&g);
379:   }
380:   if (eps->extraction==EPS_REFINED || eps->extraction==EPS_REFINED_HARMONIC) {
381:     PetscMalloc((eps->ncv+1)*eps->ncv*sizeof(PetscScalar),&Hcopy);
382:   }
383: 
384:   IPGetOrthogonalization(eps->ip,PETSC_NULL,&orthog_ref,PETSC_NULL);

386:   /* Get the starting Arnoldi vector */
387:   EPSGetStartVector(eps,0,eps->V[0],PETSC_NULL);
388: 
389:   /* Restart loop */
390:   while (eps->reason == EPS_CONVERGED_ITERATING) {
391:     eps->its++;

393:     /* Compute an nv-step Arnoldi factorization */
394:     nv = PetscMin(eps->nconv+eps->mpd,eps->ncv);
395:     if (!arnoldi->delayed) {
396:       EPSBasicArnoldi(eps,PETSC_FALSE,H,eps->ncv,eps->V,eps->nconv,&nv,f,&beta,&breakdown);
397:     } else if (orthog_ref == IP_ORTH_REFINE_NEVER) {
398:       EPSDelayedArnoldi1(eps,H,eps->ncv,eps->V,eps->nconv,&nv,f,&beta,&breakdown);
399:     } else {
400:       EPSDelayedArnoldi(eps,H,eps->ncv,eps->V,eps->nconv,&nv,f,&beta,&breakdown);
401:     }

403:     if (eps->extraction==EPS_REFINED || eps->extraction==EPS_REFINED_HARMONIC) {
404:       PetscMemcpy(Hcopy,H,eps->ncv*eps->ncv*sizeof(PetscScalar));
405:       for (i=0;i<nv-1;i++) Hcopy[nv+i*eps->ncv] = 0.0;
406:       Hcopy[nv+(nv-1)*eps->ncv] = beta;
407:     }

409:     /* Compute translation of Krylov decomposition if harmonic extraction used */
410:     if (eps->extraction==EPS_HARMONIC || eps->extraction==EPS_REFINED_HARMONIC) {
411:       EPSTranslateHarmonic(nv,H,eps->ncv,eps->target,(PetscScalar)beta,g,work);
412:       gnorm = 0.0;
413:       for (i=0;i<nv;i++)
414:         gnorm = gnorm + PetscRealPart(g[i]*PetscConj(g[i]));
415:       corrf = sqrt(1.0+gnorm);
416:     }

418:     /* Solve projected problem */
419:     EPSProjectedArnoldi(eps,H,eps->ncv,U,nv);

421:     /* Check convergence */
422:     EPSKrylovConvergence(eps,PETSC_FALSE,eps->nconv,nv-eps->nconv,H,eps->ncv,U,eps->V,nv,beta,corrf,&k,work);

424:     EPSUpdateVectors(eps,nv,eps->V,eps->nconv,PetscMin(k+1,nv),U,nv,Hcopy,eps->ncv);
425:     eps->nconv = k;

427:     EPSMonitor(eps,eps->its,eps->nconv,eps->eigr,eps->eigi,eps->errest,nv);
428:     if (breakdown) {
429:       PetscInfo2(eps,"Breakdown in Arnoldi method (it=%i norm=%g)\n",eps->its,beta);
430:       EPSGetStartVector(eps,k,eps->V[k],&breakdown);
431:       if (breakdown) {
432:         eps->reason = EPS_DIVERGED_BREAKDOWN;
433:         PetscInfo(eps,"Unable to generate more start vectors\n");
434:       }
435:     }
436:     if (eps->its >= eps->max_it) eps->reason = EPS_DIVERGED_ITS;
437:     if (eps->nconv >= eps->nev) eps->reason = EPS_CONVERGED_TOL;
438:   }
439: 
440:   PetscFree(U);
441:   PetscFree(work);
442:   if (eps->extraction==EPS_HARMONIC || eps->extraction==EPS_REFINED_HARMONIC) {
443:     PetscFree(g);
444:   }
445:   if (eps->extraction==EPS_REFINED || eps->extraction==EPS_REFINED_HARMONIC) {
446:     PetscFree(Hcopy);
447:   }
448:   return(0);
449: }

453: PetscErrorCode EPSSetFromOptions_ARNOLDI(EPS eps)
454: {
456:   PetscTruth     set,val;
457:   EPS_ARNOLDI    *arnoldi = (EPS_ARNOLDI *)eps->data;

460:   PetscOptionsBegin(((PetscObject)eps)->comm,((PetscObject)eps)->prefix,"ARNOLDI Options","EPS");
461:   PetscOptionsTruth("-eps_arnoldi_delayed","Arnoldi with delayed reorthogonalization","EPSArnoldiSetDelayed",arnoldi->delayed,&val,&set);
462:   if (set) {
463:     EPSArnoldiSetDelayed(eps,val);
464:   }
465:   PetscOptionsEnd();
466:   return(0);
467: }

472: PetscErrorCode EPSArnoldiSetDelayed_ARNOLDI(EPS eps,PetscTruth delayed)
473: {
474:   EPS_ARNOLDI    *arnoldi = (EPS_ARNOLDI *)eps->data;

477:   arnoldi->delayed = delayed;
478:   return(0);
479: }

484: /*@
485:    EPSArnoldiSetDelayed - Activates or deactivates delayed reorthogonalization 
486:    in the Arnoldi iteration. 

488:    Collective on EPS

490:    Input Parameters:
491: +  eps - the eigenproblem solver context
492: -  delayed - boolean flag

494:    Options Database Key:
495: .  -eps_arnoldi_delayed - Activates delayed reorthogonalization in Arnoldi
496:    
497:    Note:
498:    Delayed reorthogonalization is an aggressive optimization for the Arnoldi
499:    eigensolver than may provide better scalability, but sometimes makes the
500:    solver converge less than the default algorithm.

502:    Level: advanced

504: .seealso: EPSArnoldiGetDelayed()
505: @*/
506: PetscErrorCode EPSArnoldiSetDelayed(EPS eps,PetscTruth delayed)
507: {
508:   PetscErrorCode ierr, (*f)(EPS,PetscTruth);

512:   PetscObjectQueryFunction((PetscObject)eps,"EPSArnoldiSetDelayed_C",(void (**)())&f);
513:   if (f) {
514:     (*f)(eps,delayed);
515:   }
516:   return(0);
517: }

522: PetscErrorCode EPSArnoldiGetDelayed_ARNOLDI(EPS eps,PetscTruth *delayed)
523: {
524:   EPS_ARNOLDI    *arnoldi = (EPS_ARNOLDI *)eps->data;

527:   *delayed = arnoldi->delayed;
528:   return(0);
529: }

534: /*@C
535:    EPSArnoldiGetDelayed - Gets the type of reorthogonalization used during the Arnoldi
536:    iteration. 

538:    Collective on EPS

540:    Input Parameter:
541: .  eps - the eigenproblem solver context

543:    Input Parameter:
544: .  delayed - boolean flag indicating if delayed reorthogonalization has been enabled

546:    Level: advanced

548: .seealso: EPSArnoldiSetDelayed()
549: @*/
550: PetscErrorCode EPSArnoldiGetDelayed(EPS eps,PetscTruth *delayed)
551: {
552:   PetscErrorCode ierr, (*f)(EPS,PetscTruth*);

556:   PetscObjectQueryFunction((PetscObject)eps,"EPSArnoldiGetDelayed_C",(void (**)())&f);
557:   if (f) {
558:     (*f)(eps,delayed);
559:   }
560:   return(0);
561: }

565: PetscErrorCode EPSDestroy_ARNOLDI(EPS eps)
566: {

571:   EPSDestroy_Default(eps);
572:   PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSArnoldiSetDelayed_C","",PETSC_NULL);
573:   PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSArnoldiGetDelayed_C","",PETSC_NULL);
574:   return(0);
575: }

579: PetscErrorCode EPSView_ARNOLDI(EPS eps,PetscViewer viewer)
580: {
582:   PetscTruth     isascii;
583:   EPS_ARNOLDI    *arnoldi = (EPS_ARNOLDI *)eps->data;

586:   PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
587:   if (!isascii) {
588:     SETERRQ1(1,"Viewer type %s not supported for EPSARNOLDI",((PetscObject)viewer)->type_name);
589:   }
590:   if (arnoldi->delayed) {
591:     PetscViewerASCIIPrintf(viewer,"using delayed reorthogonalization\n");
592:   }
593:   return(0);
594: }

596: EXTERN PetscErrorCode EPSSolve_TS_ARNOLDI(EPS);

601: PetscErrorCode EPSCreate_ARNOLDI(EPS eps)
602: {
604:   EPS_ARNOLDI    *arnoldi;
605: 
607:   PetscNew(EPS_ARNOLDI,&arnoldi);
608:   PetscLogObjectMemory(eps,sizeof(EPS_ARNOLDI));
609:   eps->data                      = (void *)arnoldi;
610:   eps->ops->setup                = EPSSetUp_ARNOLDI;
611:   eps->ops->setfromoptions       = EPSSetFromOptions_ARNOLDI;
612:   eps->ops->destroy              = EPSDestroy_ARNOLDI;
613:   eps->ops->view                 = EPSView_ARNOLDI;
614:   eps->ops->backtransform        = EPSBackTransform_Default;
615:   eps->ops->computevectors       = EPSComputeVectors_Schur;
616:   arnoldi->delayed               = PETSC_FALSE;
617:   PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSArnoldiSetDelayed_C","EPSArnoldiSetDelayed_ARNOLDI",EPSArnoldiSetDelayed_ARNOLDI);
618:   PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSArnoldiGetDelayed_C","EPSArnoldiGetDelayed_ARNOLDI",EPSArnoldiGetDelayed_ARNOLDI);
619:   return(0);
620: }