Actual source code: primmes.c
1: /*
2: This file implements a wrapper to the PRIMME library
3: Contributed by Eloy Romero
5: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
6: SLEPc - Scalable Library for Eigenvalue Problem Computations
7: Copyright (c) 2002-2007, Universidad Politecnica de Valencia, Spain
9: This file is part of SLEPc. See the README file for conditions of use
10: and additional information.
11: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
12: */
14: #include src/eps/epsimpl.h
15: #include src/st/stimpl.h
18: #include "primme.h"
21: typedef struct {
22: primme_params primme; /* param struc */
23: primme_preset_method method; /* primme method */
24: Mat A; /* problem matrix */
25: Vec w; /* store reciprocal A diagonal */
26: Vec x,y; /* auxiliar vectors */
27: } EPS_PRIMME;
29: const char *methodList[] = {
30: "dynamic",
31: "default_min_time",
32: "default_min_matvecs",
33: "arnoldi",
34: "gd",
35: "gd_plusk",
36: "gd_olsen_plusk",
37: "jd_olsen_plusk",
38: "rqi",
39: "jdqr",
40: "jdqmr",
41: "jdqmr_etol",
42: "subspace_iteration",
43: "lobpcg_orthobasis",
44: "lobpcg_orthobasis_window"
45: };
46: const char *precondList[] = {"none", "diagonal"};
47: EPSPRIMMEMethod methodN[] = {
48: EPSPRIMME_DYNAMIC,
49: EPSPRIMME_DEFAULT_MIN_TIME,
50: EPSPRIMME_DEFAULT_MIN_MATVECS,
51: EPSPRIMME_ARNOLDI,
52: EPSPRIMME_GD,
53: EPSPRIMME_GD_PLUSK,
54: EPSPRIMME_GD_OLSEN_PLUSK,
55: EPSPRIMME_JD_OLSEN_PLUSK,
56: EPSPRIMME_RQI,
57: EPSPRIMME_JDQR,
58: EPSPRIMME_JDQMR,
59: EPSPRIMME_JDQMR_ETOL,
60: EPSPRIMME_SUBSPACE_ITERATION,
61: EPSPRIMME_LOBPCG_ORTHOBASIS,
62: EPSPRIMME_LOBPCG_ORTHOBASIS_WINDOW
63: };
64: EPSPRIMMEPrecond precondN[] = {EPSPRIMME_NONE, EPSPRIMME_DIAGONAL};
66: static void multMatvec_PRIMME(void *in, void *out, int *blockSize, primme_params *primme);
67: static void applyPreconditioner_PRIMME(void *in, void *out, int *blockSize, struct primme_params *primme);
69: static void par_GlobalSumDouble(void *sendBuf, void *recvBuf, int *count, primme_params *primme) {
70: MPI_Allreduce((double*)sendBuf, (double*)recvBuf, *count, MPI_DOUBLE, MPI_SUM, ((EPS)(primme->commInfo))->comm);
71: }
75: PetscErrorCode EPSSetUp_PRIMME(EPS eps)
76: {
78: PetscInt N, n;
79: int numProcs, procID;
80: EPS_PRIMME *ops = (EPS_PRIMME *)eps->data;
81: primme_params *primme = &(((EPS_PRIMME *)eps->data)->primme);
85: MPI_Comm_size(eps->comm,&numProcs);
86: MPI_Comm_rank(eps->comm,&procID);
87:
88: /* Check some constraints and set some default values */
89: VecGetSize(eps->vec_initial,&N);
90: VecGetLocalSize(eps->vec_initial,&n);
92: if (!eps->max_it) eps->max_it = PetscMax(1000,N);
93: STGetOperators(eps->OP, &ops->A, PETSC_NULL);
94: if (!ops->A) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"The problem matrix has to be specified first");
95: if (!eps->ishermitian)
96: SETERRQ(PETSC_ERR_SUP,"PRIMME is only available for Hermitian problems");
97: if (eps->isgeneralized)
98: SETERRQ(PETSC_ERR_SUP,"PRIMME is not available for generalized problems");
100: /* Transfer SLEPc options to PRIMME options */
101: primme->n = N;
102: primme->nLocal = n;
103: primme->numEvals = eps->nev;
104: primme->matrix = ops;
105: primme->commInfo = eps;
106: primme->maxMatvecs = eps->max_it;
107: primme->eps = eps->tol;
108: primme->numProcs = numProcs;
109: primme->procID = procID;
111: switch(eps->which) {
112: case EPS_LARGEST_REAL:
113: primme->target = primme_largest;
114: break;
116: case EPS_SMALLEST_REAL:
117: primme->target = primme_smallest;
118: break;
119:
120: default:
121: SETERRQ(PETSC_ERR_SUP,"PRIMME only allows EPS_LARGEST_REAL and EPS_SMALLEST_REAL for 'which' value");
122: break;
123: }
124:
125: if (primme_set_method(ops->method, primme) < 0)
126: SETERRQ(PETSC_ERR_SUP,"PRIMME method not valid");
127:
128: /* If user sets ncv, maxBasisSize is modified. If not, ncv is set as maxBasisSize */
129: if (eps->ncv) primme->maxBasisSize = eps->ncv;
130: else eps->ncv = primme->maxBasisSize;
131: if (eps->ncv < eps->nev+primme->maxBlockSize)
132: SETERRQ(PETSC_ERR_SUP,"PRIMME needs ncv >= nev+maxBlockSize");
134: /* Set workspace */
135: EPSAllocateSolutionContiguous(eps);
137: if (primme->correctionParams.precondition) {
138: /* Calc reciprocal A diagonal */
139: VecDuplicate(eps->vec_initial, &ops->w);
140: MatGetDiagonal(ops->A, ops->w);
141: VecReciprocal(ops->w);
142: primme->preconditioner = PETSC_NULL;
143: primme->applyPreconditioner = applyPreconditioner_PRIMME;
144: }
146: /* Prepare auxiliary vectors */
147: VecCreateMPIWithArray(PETSC_COMM_WORLD,n,N,PETSC_NULL,&ops->x);
148: VecCreateMPIWithArray(PETSC_COMM_WORLD,n,N,PETSC_NULL,&ops->y);
149:
150: return(0);
151: }
155: PetscErrorCode EPSSolve_PRIMME(EPS eps)
156: {
158: EPS_PRIMME *ops = (EPS_PRIMME *)eps->data;
159: PetscScalar *a;
160: #ifdef PETSC_USE_COMPLEX
161: int i;
162: PetscReal *evals;
163: #endif
167: /* Reset some parameters left from previous runs */
168: ops->primme.aNorm = 0.0;
169: ops->primme.initSize = 1;
170: ops->primme.iseed[0] = -1;
172: /* Copy vec_initial to V[0] vector */
173: VecCopy(eps->vec_initial, eps->V[0]);
174:
175: /* Call PRIMME solver */
176: VecGetArray(eps->V[0], &a);
177: #ifndef PETSC_USE_COMPLEX
178: dprimme(eps->eigr, a, eps->errest, &ops->primme);
179: #else
180: /* PRIMME returns real eigenvalues, but SLEPc works with complex ones */
181: PetscMalloc(eps->ncv*sizeof(PetscReal),&evals);
182: zprimme(evals, (Complex_Z*)a, eps->errest, &ops->primme);
183: for (i=0;i<eps->ncv;i++)
184: eps->eigr[i] = evals[i];
185: PetscFree(evals);
186: #endif
187: VecRestoreArray(eps->V[0], &a);
188:
189: switch(ierr) {
190: case 0: /* Successful */
191: break;
193: case -1:
194: SETERRQ(PETSC_ERR_SUP,"PRIMME: Failed to open output file");
195: break;
197: case -2:
198: SETERRQ(PETSC_ERR_SUP,"PRIMME: Insufficient integer or real workspace allocated");
199: break;
201: case -3:
202: SETERRQ(PETSC_ERR_SUP,"PRIMME: main_iter encountered a problem");
203: break;
205: default:
206: SETERRQ(PETSC_ERR_SUP,"PRIMME: some parameters wrong configured");
207: break;
208: }
210: eps->nconv = ops->primme.initSize >= 0 ? ops->primme.initSize : 0;
211: eps->reason = eps->ncv >= eps->nev ? EPS_CONVERGED_TOL : EPS_DIVERGED_ITS;
212: eps->its = ops->primme.stats.numOuterIterations;
213: eps->OP->applys = ops->primme.stats.numMatvecs;
215: return(0);
216: }
220: static void multMatvec_PRIMME(void *in, void *out, int *blockSize, primme_params *primme)
221: {
223: int i, N = primme->n;
224: EPS_PRIMME *ops = (EPS_PRIMME *)primme->matrix;
225: Vec x = ops->x;
226: Vec y = ops->y;
227: Mat A = ops->A;
230: for (i=0;i<*blockSize;i++) {
231: /* build vectors using 'in' an 'out' workspace */
232: VecPlaceArray(x, (PetscScalar*)in+N*i ); CHKERRABORT(PETSC_COMM_WORLD,ierr);
233: VecPlaceArray(y, (PetscScalar*)out+N*i ); CHKERRABORT(PETSC_COMM_WORLD,ierr);
235: MatMult(A, x, y); CHKERRABORT(PETSC_COMM_WORLD,ierr);
236:
237: VecResetArray(x); CHKERRABORT(PETSC_COMM_WORLD,ierr);
238: VecResetArray(y); CHKERRABORT(PETSC_COMM_WORLD,ierr);
239: }
240: PetscFunctionReturnVoid();
241: }
245: static void applyPreconditioner_PRIMME(void *in, void *out, int *blockSize, struct primme_params *primme)
246: {
248: int i, N = primme->n;
249: EPS_PRIMME *ops = (EPS_PRIMME *)primme->matrix;
250: Vec x = ops->x;
251: Vec y = ops->y;
252: Vec w = ops->w;
253:
255: for (i=0;i<*blockSize;i++) {
256: /* build vectors using 'in' an 'out' workspace */
257: VecPlaceArray(x, (PetscScalar*)in+N*i ); CHKERRABORT(PETSC_COMM_WORLD,ierr);
258: VecPlaceArray(y, (PetscScalar*)out+N*i ); CHKERRABORT(PETSC_COMM_WORLD,ierr);
260: VecPointwiseMult(y, w, x); CHKERRABORT(PETSC_COMM_WORLD,ierr);
261:
262: VecResetArray(x); CHKERRABORT(PETSC_COMM_WORLD,ierr);
263: VecResetArray(y); CHKERRABORT(PETSC_COMM_WORLD,ierr);
264: }
265: PetscFunctionReturnVoid();
266: }
271: PetscErrorCode EPSDestroy_PRIMME(EPS eps)
272: {
274: EPS_PRIMME *ops = (EPS_PRIMME *)eps->data;
278:
279: if (ops->primme.correctionParams.precondition) {
280: VecDestroy(ops->w);
281: }
282: primme_Free(&ops->primme);
283: VecDestroy(ops->x);
284: VecDestroy(ops->y);
285: PetscFree(eps->data);
286: EPSFreeSolutionContiguous(eps);
287:
288: return(0);
289: }
293: PetscErrorCode EPSView_PRIMME(EPS eps,PetscViewer viewer)
294: {
296: PetscTruth isascii;
297: primme_params *primme = &((EPS_PRIMME *)eps->data)->primme;
298: EPSPRIMMEMethod methodn;
299: EPSPRIMMEPrecond precondn;
302: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
303: if (!isascii) {
304: SETERRQ1(1,"Viewer type %s not supported for EPSPRIMME",((PetscObject)viewer)->type_name);
305: }
306:
307: PetscViewerASCIIPrintf(viewer,"PRIMME solver block size: %d\n",primme->maxBlockSize);
308: EPSPRIMMEGetMethod(eps, &methodn);
309: PetscViewerASCIIPrintf(viewer,"PRIMME solver method: %s\n", methodList[methodn]);
310: EPSPRIMMEGetPrecond(eps, &precondn);
311: PetscViewerASCIIPrintf(viewer,"PRIMME solver preconditioning: %s\n", precondList[precondn]);
313: return(0);
314: }
318: PetscErrorCode EPSSetFromOptions_PRIMME(EPS eps)
319: {
321: EPS_PRIMME *ops = (EPS_PRIMME *)eps->data;
322: PetscInt op;
323: PetscTruth flg;
326:
327: PetscOptionsHead("PRIMME options");
329: op = ops->primme.maxBlockSize;
330: PetscOptionsInt("-eps_primme_block_size"," maximum block size","EPSPRIMMESetBlockSize",op,&op,&flg);
331: if (flg) {EPSPRIMMESetBlockSize(eps,op);}
332: op = 0;
333: PetscOptionsEList("-eps_primme_method","set method for solving the eigenproblem",
334: "EPSPRIMMESetMethod",methodList,15,methodList[0],&op,&flg);
335: if (flg) {EPSPRIMMESetMethod(eps, methodN[op]);}
336: PetscOptionsEList("-eps_primme_precond","set preconditioner type",
337: "EPSPRIMMESetPrecond",precondList,2,precondList[0],&op,&flg);
338: if (flg) {EPSPRIMMESetPrecond(eps, precondN[op]);}
339:
340: PetscOptionsTail();
341:
342: return(0);
343: }
348: PetscErrorCode EPSPRIMMESetBlockSize_PRIMME(EPS eps,int bs)
349: {
350: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
354: if (bs == PETSC_DEFAULT) ops->primme.maxBlockSize = 1;
355: else if (bs <= 0) {
356: SETERRQ(1, "PRIMME: wrong block size");
357: } else ops->primme.maxBlockSize = bs;
358:
359: return(0);
360: }
365: /*@
366: EPSPRIMMESetBlockSize - The maximum block size the code will try to use.
367: The user should set
368: this based on the architecture specifics of the target computer,
369: as well as any a priori knowledge of multiplicities. The code does
370: NOT require BlockSize > 1 to find multiple eigenvalues. For some
371: methods, keeping BlockSize = 1 yields the best overall performance.
373: Collective on EPS
375: Input Parameters:
376: + eps - the eigenproblem solver context
377: - bs - block size
379: Options Database Key:
380: . -eps_primme_block_size - Sets the max allowed block size value
382: Notes:
383: If the block size is not set, the value established by primme_initialize
384: is used.
386: Level: advanced
387: .seealso: EPSPRIMMEGetBlockSize()
388: @*/
389: PetscErrorCode EPSPRIMMESetBlockSize(EPS eps,int bs)
390: {
391: PetscErrorCode ierr, (*f)(EPS,int);
394:
396: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMESetBlockSize_C",(void (**)())&f);
397: if (f) {
398: (*f)(eps,bs);
399: }
400:
401: return(0);
402: }
407: PetscErrorCode EPSPRIMMEGetBlockSize_PRIMME(EPS eps,int *bs)
408: {
409: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
413: if (bs) *bs = ops->primme.maxBlockSize;
414:
415: return(0);
416: }
421: /*@
422: EPSPRIMMEGetBlockSize - Get the maximum block size the code will try to use.
423:
424: Collective on EPS
426: Input Parameters:
427: . eps - the eigenproblem solver context
428:
429: Output Parameters:
430: . bs - returned block size
432: Level: advanced
433: .seealso: EPSPRIMMESetBlockSize()
434: @*/
435: PetscErrorCode EPSPRIMMEGetBlockSize(EPS eps,int *bs)
436: {
437: PetscErrorCode ierr, (*f)(EPS,int*);
440:
442: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMEGetBlockSize_C",(void (**)())&f);
443: if (f) {
444: (*f)(eps,bs);
445: }
446:
447: return(0);
448: }
453: PetscErrorCode EPSPRIMMESetMethod_PRIMME(EPS eps, EPSPRIMMEMethod method)
454: {
455: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
459: if (method == PETSC_DEFAULT) ops->method = DYNAMIC;
460: else ops->method = (primme_preset_method)method;
462: return(0);
463: }
468: /*@
469: EPSPRIMMESetMethod - Sets the method for the PRIMME library.
471: Collective on EPS
473: Input Parameters:
474: + eps - the eigenproblem solver context
475: - method - method that will be used by PRIMME. It must be one of:
476: EPSPRIMME_DYNAMIC, EPSPRIMME_DEFAULT_MIN_TIME(EPSPRIMME_JDQMR_ETOL),
477: EPSPRIMME_DEFAULT_MIN_MATVECS(EPSPRIMME_GD_OLSEN_PLUSK), EPSPRIMME_ARNOLDI,
478: EPSPRIMME_GD, EPSPRIMME_GD_PLUSK, EPSPRIMME_GD_OLSEN_PLUSK,
479: EPSPRIMME_JD_OLSEN_PLUSK, EPSPRIMME_RQI, EPSPRIMME_JDQR, EPSPRIMME_JDQMR,
480: EPSPRIMME_JDQMR_ETOL, EPSPRIMME_SUBSPACE_ITERATION,
481: EPSPRIMME_LOBPCG_ORTHOBASIS, EPSPRIMME_LOBPCG_ORTHOBASIS_WINDOW
483: Options Database Key:
484: . -eps_primme_set_method - Sets the method for the PRIMME library.
486: Note:
487: If not set, the method defaults to EPSPRIMME_DYNAMIC.
489: Level: advanced
491: .seealso: EPSPRIMMEGetMethod(), EPSPRIMMEMethod
492: @*/
493: PetscErrorCode EPSPRIMMESetMethod(EPS eps, EPSPRIMMEMethod method)
494: {
495: PetscErrorCode ierr, (*f)(EPS,EPSPRIMMEMethod);
498:
500: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMESetMethod_C",(void (**)())&f);
501: if (f) {
502: (*f)(eps,method);
503: }
504:
505: return(0);
506: }
511: PetscErrorCode EPSPRIMMEGetMethod_PRIMME(EPS eps, EPSPRIMMEMethod *method)
512: {
513: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
517: if (method)
518: *method = (EPSPRIMMEMethod)ops->method;
520: return(0);
521: }
526: /*@C
527: EPSPRIMMEGetMethod - Gets the method for the PRIMME library.
529: Mon Collective on EPS
531: Input Parameters:
532: . eps - the eigenproblem solver context
533:
534: Output Parameters:
535: . method - method that will be used by PRIMME. It must be one of:
536: EPSPRIMME_DYNAMIC, EPSPRIMME_DEFAULT_MIN_TIME(EPSPRIMME_JDQMR_ETOL),
537: EPSPRIMME_DEFAULT_MIN_MATVECS(EPSPRIMME_GD_OLSEN_PLUSK), EPSPRIMME_ARNOLDI,
538: EPSPRIMME_GD, EPSPRIMME_GD_PLUSK, EPSPRIMME_GD_OLSEN_PLUSK,
539: EPSPRIMME_JD_OLSEN_PLUSK, EPSPRIMME_RQI, EPSPRIMME_JDQR, EPSPRIMME_JDQMR,
540: EPSPRIMME_JDQMR_ETOL, EPSPRIMME_SUBSPACE_ITERATION,
541: EPSPRIMME_LOBPCG_ORTHOBASIS, EPSPRIMME_LOBPCG_ORTHOBASIS_WINDOW
543: Level: advanced
545: .seealso: EPSPRIMMESetMethod(), EPSPRIMMEMethod
546: @*/
547: PetscErrorCode EPSPRIMMEGetMethod(EPS eps, EPSPRIMMEMethod *method)
548: {
549: PetscErrorCode ierr, (*f)(EPS,EPSPRIMMEMethod*);
552:
554: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMEGetMethod_C",(void (**)())&f);
555: if (f) {
556: (*f)(eps,method);
557: }
558:
559: return(0);
560: }
566: PetscErrorCode EPSPRIMMESetPrecond_PRIMME(EPS eps, EPSPRIMMEPrecond precond)
567: {
568: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
572: if (precond == EPSPRIMME_NONE) ops->primme.correctionParams.precondition = 0;
573: else ops->primme.correctionParams.precondition = 1;
574:
575: return(0);
576: }
581: /*@
582: EPSPRIMMESetPrecond - Sets the preconditioner to be used in the PRIMME
583: library. Currently, only diagonal preconditioning is supported.
585: Collective on EPS
587: Input Parameters:
588: + eps - the eigenproblem solver context
589: - precond - either EPSPRIMME_NONE (no preconditioning) or EPSPRIMME_DIAGONAL
590: (diagonal preconditioning)
592: Options Database Key:
593: . -eps_primme_precond - Sets either 'none' or 'diagonal' preconditioner
595: Note:
596: The default is no preconditioning.
597:
598: Level: advanced
600: .seealso: EPSPRIMMEGetPrecond(), EPSPRIMMEPrecond
601: @*/
602: PetscErrorCode EPSPRIMMESetPrecond(EPS eps, EPSPRIMMEPrecond precond)
603: {
604: PetscErrorCode ierr, (*f)(EPS, EPSPRIMMEPrecond);
607:
609: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMESetPrecond_C",(void (**)())&f);
610: if (f) {
611: (*f)(eps, precond);
612: }
613:
614: return(0);
615: }
621: PetscErrorCode EPSPRIMMEGetPrecond_PRIMME(EPS eps, EPSPRIMMEPrecond *precond)
622: {
623: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
627: if (precond)
628: *precond = ops->primme.correctionParams.precondition ? EPSPRIMME_DIAGONAL : EPSPRIMME_NONE;
629:
630: return(0);
631: }
636: /*@C
637: EPSPRIMMEGetPrecond - Gets the preconditioner to be used in the PRIMME
638: library.
640: Collective on EPS
642: Input Parameters:
643: . eps - the eigenproblem solver context
644:
645: Output Parameters:
646: . precond - either EPSPRIMME_NONE (no preconditioning) or EPSPRIMME_DIAGONAL
647: (diagonal preconditioning)
649: Level: advanced
651: .seealso: EPSPRIMMESetPrecond(), EPSPRIMMEPrecond
652: @*/
653: PetscErrorCode EPSPRIMMEGetPrecond(EPS eps, EPSPRIMMEPrecond *precond)
654: {
655: PetscErrorCode ierr, (*f)(EPS, EPSPRIMMEPrecond*);
658:
660: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMEGetPrecond_C",(void (**)())&f);
661: if (f) {
662: (*f)(eps, precond);
663: }
664:
665: return(0);
666: }
672: PetscErrorCode EPSCreate_PRIMME(EPS eps)
673: {
675: EPS_PRIMME *primme;
678:
679: PetscNew(EPS_PRIMME,&primme);
680: PetscLogObjectMemory(eps,sizeof(EPS_PRIMME));
681: eps->data = (void *) primme;
682: eps->ops->solve = EPSSolve_PRIMME;
683: eps->ops->setup = EPSSetUp_PRIMME;
684: eps->ops->setfromoptions = EPSSetFromOptions_PRIMME;
685: eps->ops->destroy = EPSDestroy_PRIMME;
686: eps->ops->view = EPSView_PRIMME;
687: eps->ops->backtransform = EPSBackTransform_Default;
688: eps->ops->computevectors = EPSComputeVectors_Default;
690: primme_initialize(&primme->primme);
691: primme->primme.matrixMatvec = multMatvec_PRIMME;
692: primme->primme.globalSumDouble = par_GlobalSumDouble;
693: primme->method = DYNAMIC;
694: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMESetBlockSize_C","EPSPRIMMESetBlockSize_PRIMME",EPSPRIMMESetBlockSize_PRIMME);
695: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMESetMethod_C","EPSPRIMMESetMethod_PRIMME",EPSPRIMMESetMethod_PRIMME);
696: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMESetPrecond_C","EPSPRIMMESetPrecond_PRIMME",EPSPRIMMESetPrecond_PRIMME);
697:
698: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMEGetBlockSize_C","EPSPRIMMEGetBlockSize_PRIMME",EPSPRIMMEGetBlockSize_PRIMME);
699: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMEGetMethod_C","EPSPRIMMEGetMethod_PRIMME",EPSPRIMMEGetMethod_PRIMME);
700: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMEGetPrecond_C","EPSPRIMMEGetPrecond_PRIMME",EPSPRIMMEGetPrecond_PRIMME);
702: eps->which = EPS_LARGEST_REAL;
704: return(0);
705: }