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: }