Actual source code: primmes.c
1: /*
2: This file implements a wrapper to the PRIMME library
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2009, Universidad Politecnica de Valencia, Spain
8: This file is part of SLEPc.
9:
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include private/epsimpl.h
25: #include private/stimpl.h
28: #include "primme.h"
31: typedef struct {
32: primme_params primme; /* param struc */
33: primme_preset_method method; /* primme method */
34: Mat A; /* problem matrix */
35: Vec w; /* store reciprocal A diagonal */
36: Vec x,y; /* auxiliar vectors */
37: } EPS_PRIMME;
39: const char *methodList[] = {
40: "dynamic",
41: "default_min_time",
42: "default_min_matvecs",
43: "arnoldi",
44: "gd",
45: "gd_plusk",
46: "gd_olsen_plusk",
47: "jd_olsen_plusk",
48: "rqi",
49: "jdqr",
50: "jdqmr",
51: "jdqmr_etol",
52: "subspace_iteration",
53: "lobpcg_orthobasis",
54: "lobpcg_orthobasis_window"
55: };
56: const char *precondList[] = {"none", "diagonal"};
57: EPSPRIMMEMethod methodN[] = {
58: EPSPRIMME_DYNAMIC,
59: EPSPRIMME_DEFAULT_MIN_TIME,
60: EPSPRIMME_DEFAULT_MIN_MATVECS,
61: EPSPRIMME_ARNOLDI,
62: EPSPRIMME_GD,
63: EPSPRIMME_GD_PLUSK,
64: EPSPRIMME_GD_OLSEN_PLUSK,
65: EPSPRIMME_JD_OLSEN_PLUSK,
66: EPSPRIMME_RQI,
67: EPSPRIMME_JDQR,
68: EPSPRIMME_JDQMR,
69: EPSPRIMME_JDQMR_ETOL,
70: EPSPRIMME_SUBSPACE_ITERATION,
71: EPSPRIMME_LOBPCG_ORTHOBASIS,
72: EPSPRIMME_LOBPCG_ORTHOBASISW
73: };
74: EPSPRIMMEPrecond precondN[] = {EPSPRIMME_NONE, EPSPRIMME_DIAGONAL};
76: static void multMatvec_PRIMME(void *in, void *out, int *blockSize, primme_params *primme);
77: static void applyPreconditioner_PRIMME(void *in, void *out, int *blockSize, struct primme_params *primme);
79: void par_GlobalSumDouble(void *sendBuf, void *recvBuf, int *count, primme_params *primme) {
81: MPI_Allreduce((double*)sendBuf, (double*)recvBuf, *count, MPI_DOUBLE, MPI_SUM, ((PetscObject)(primme->commInfo))->comm);CHKERRABORT(((PetscObject)(primme->commInfo))->comm,ierr);
82: }
86: PetscErrorCode EPSSetUp_PRIMME(EPS eps)
87: {
89: PetscInt N, n;
90: PetscMPIInt numProcs, procID;
91: EPS_PRIMME *ops = (EPS_PRIMME *)eps->data;
92: primme_params *primme = &(((EPS_PRIMME *)eps->data)->primme);
96: MPI_Comm_size(((PetscObject)eps)->comm,&numProcs);
97: MPI_Comm_rank(((PetscObject)eps)->comm,&procID);
98:
99: /* Check some constraints and set some default values */
100: VecGetSize(eps->vec_initial,&N);
101: VecGetLocalSize(eps->vec_initial,&n);
103: if (!eps->max_it) eps->max_it = PetscMax(1000,N);
104: STGetOperators(eps->OP, &ops->A, PETSC_NULL);
105: if (!ops->A) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"The problem matrix has to be specified first");
106: if (!eps->ishermitian)
107: SETERRQ(PETSC_ERR_SUP,"PRIMME is only available for Hermitian problems");
108: if (eps->isgeneralized)
109: SETERRQ(PETSC_ERR_SUP,"PRIMME is not available for generalized problems");
111: /* Transfer SLEPc options to PRIMME options */
112: primme->n = N;
113: primme->nLocal = n;
114: primme->numEvals = eps->nev;
115: primme->matrix = ops;
116: primme->commInfo = eps;
117: primme->maxMatvecs = eps->max_it;
118: primme->eps = eps->tol;
119: primme->numProcs = numProcs;
120: primme->procID = procID;
121: primme->printLevel = 0;
123: switch(eps->which) {
124: case EPS_LARGEST_REAL:
125: primme->target = primme_largest;
126: break;
128: case EPS_SMALLEST_REAL:
129: primme->target = primme_smallest;
130: break;
131:
132: default:
133: SETERRQ(PETSC_ERR_SUP,"PRIMME only allows EPS_LARGEST_REAL and EPS_SMALLEST_REAL for 'which' value");
134: break;
135: }
136:
137: if (primme_set_method(ops->method, primme) < 0)
138: SETERRQ(PETSC_ERR_SUP,"PRIMME method not valid");
139:
140: /* If user sets ncv, maxBasisSize is modified. If not, ncv is set as maxBasisSize */
141: if (eps->ncv) primme->maxBasisSize = eps->ncv;
142: else eps->ncv = primme->maxBasisSize;
143: if (eps->ncv < eps->nev+primme->maxBlockSize)
144: SETERRQ(PETSC_ERR_SUP,"PRIMME needs ncv >= nev+maxBlockSize");
145: if (eps->mpd) PetscInfo(eps,"Warning: parameter mpd ignored\n");
147: if (eps->extraction) {
148: PetscInfo(eps,"Warning: extraction type ignored\n");
149: }
151: /* Set workspace */
152: EPSAllocateSolution(eps);
154: if (primme->correctionParams.precondition) {
155: /* Calc reciprocal A diagonal */
156: VecDuplicate(eps->vec_initial, &ops->w);
157: MatGetDiagonal(ops->A, ops->w);
158: VecReciprocal(ops->w);
159: primme->preconditioner = PETSC_NULL;
160: primme->applyPreconditioner = applyPreconditioner_PRIMME;
161: }
163: /* Prepare auxiliary vectors */
164: VecCreateMPIWithArray(PETSC_COMM_WORLD,n,N,PETSC_NULL,&ops->x);
165: VecCreateMPIWithArray(PETSC_COMM_WORLD,n,N,PETSC_NULL,&ops->y);
166:
167: return(0);
168: }
172: PetscErrorCode EPSSolve_PRIMME(EPS eps)
173: {
175: EPS_PRIMME *ops = (EPS_PRIMME *)eps->data;
176: PetscScalar *a;
177: #ifdef PETSC_USE_COMPLEX
178: PetscInt i;
179: PetscReal *evals;
180: #endif
184: /* Reset some parameters left from previous runs */
185: ops->primme.aNorm = 0.0;
186: ops->primme.initSize = 1;
187: ops->primme.iseed[0] = -1;
189: /* Copy vec_initial to V[0] vector */
190: VecCopy(eps->vec_initial, eps->V[0]);
191:
192: /* Call PRIMME solver */
193: VecGetArray(eps->V[0], &a);
194: #ifndef PETSC_USE_COMPLEX
195: dprimme(eps->eigr, a, eps->errest, &ops->primme);
196: #else
197: /* PRIMME returns real eigenvalues, but SLEPc works with complex ones */
198: PetscMalloc(eps->ncv*sizeof(PetscReal),&evals);
199: zprimme(evals, (Complex_Z*)a, eps->errest, &ops->primme);
200: for (i=0;i<eps->ncv;i++)
201: eps->eigr[i] = evals[i];
202: PetscFree(evals);
203: #endif
204: VecRestoreArray(eps->V[0], &a);
205:
206: switch(ierr) {
207: case 0: /* Successful */
208: break;
210: case -1:
211: SETERRQ(PETSC_ERR_SUP,"PRIMME: Failed to open output file");
212: break;
214: case -2:
215: SETERRQ(PETSC_ERR_SUP,"PRIMME: Insufficient integer or real workspace allocated");
216: break;
218: case -3:
219: SETERRQ(PETSC_ERR_SUP,"PRIMME: main_iter encountered a problem");
220: break;
222: default:
223: SETERRQ(PETSC_ERR_SUP,"PRIMME: some parameters wrong configured");
224: break;
225: }
227: eps->nconv = ops->primme.initSize >= 0 ? ops->primme.initSize : 0;
228: eps->reason = eps->ncv >= eps->nev ? EPS_CONVERGED_TOL : EPS_DIVERGED_ITS;
229: eps->its = ops->primme.stats.numOuterIterations;
230: eps->OP->applys = ops->primme.stats.numMatvecs;
232: return(0);
233: }
237: static void multMatvec_PRIMME(void *in, void *out, int *blockSize, primme_params *primme)
238: {
240: PetscInt i, N = primme->n;
241: EPS_PRIMME *ops = (EPS_PRIMME *)primme->matrix;
242: Vec x = ops->x;
243: Vec y = ops->y;
244: Mat A = ops->A;
247: for (i=0;i<*blockSize;i++) {
248: /* build vectors using 'in' an 'out' workspace */
249: VecPlaceArray(x, (PetscScalar*)in+N*i ); CHKERRABORT(PETSC_COMM_WORLD,ierr);
250: VecPlaceArray(y, (PetscScalar*)out+N*i ); CHKERRABORT(PETSC_COMM_WORLD,ierr);
252: MatMult(A, x, y); CHKERRABORT(PETSC_COMM_WORLD,ierr);
253:
254: VecResetArray(x); CHKERRABORT(PETSC_COMM_WORLD,ierr);
255: VecResetArray(y); CHKERRABORT(PETSC_COMM_WORLD,ierr);
256: }
257: PetscFunctionReturnVoid();
258: }
262: static void applyPreconditioner_PRIMME(void *in, void *out, int *blockSize, struct primme_params *primme)
263: {
265: PetscInt i, N = primme->n;
266: EPS_PRIMME *ops = (EPS_PRIMME *)primme->matrix;
267: Vec x = ops->x;
268: Vec y = ops->y;
269: Vec w = ops->w;
270:
272: for (i=0;i<*blockSize;i++) {
273: /* build vectors using 'in' an 'out' workspace */
274: VecPlaceArray(x, (PetscScalar*)in+N*i ); CHKERRABORT(PETSC_COMM_WORLD,ierr);
275: VecPlaceArray(y, (PetscScalar*)out+N*i ); CHKERRABORT(PETSC_COMM_WORLD,ierr);
277: VecPointwiseMult(y, w, x); CHKERRABORT(PETSC_COMM_WORLD,ierr);
278:
279: VecResetArray(x); CHKERRABORT(PETSC_COMM_WORLD,ierr);
280: VecResetArray(y); CHKERRABORT(PETSC_COMM_WORLD,ierr);
281: }
282: PetscFunctionReturnVoid();
283: }
288: PetscErrorCode EPSDestroy_PRIMME(EPS eps)
289: {
291: EPS_PRIMME *ops = (EPS_PRIMME *)eps->data;
295:
296: if (ops->primme.correctionParams.precondition) {
297: VecDestroy(ops->w);
298: }
299: primme_Free(&ops->primme);
300: VecDestroy(ops->x);
301: VecDestroy(ops->y);
302: PetscFree(eps->data);
303: EPSFreeSolution(eps);
304:
305: return(0);
306: }
310: PetscErrorCode EPSView_PRIMME(EPS eps,PetscViewer viewer)
311: {
313: PetscTruth isascii;
314: primme_params *primme = &((EPS_PRIMME *)eps->data)->primme;
315: EPSPRIMMEMethod methodn;
316: EPSPRIMMEPrecond precondn;
319: PetscTypeCompare((PetscObject)viewer,PETSC_VIEWER_ASCII,&isascii);
320: if (!isascii) {
321: SETERRQ1(1,"Viewer type %s not supported for EPSPRIMME",((PetscObject)viewer)->type_name);
322: }
323:
324: PetscViewerASCIIPrintf(viewer,"PRIMME solver block size: %d\n",primme->maxBlockSize);
325: EPSPRIMMEGetMethod(eps, &methodn);
326: PetscViewerASCIIPrintf(viewer,"PRIMME solver method: %s\n", methodList[methodn]);
327: EPSPRIMMEGetPrecond(eps, &precondn);
328: PetscViewerASCIIPrintf(viewer,"PRIMME solver preconditioning: %s\n", precondList[precondn]);
330: /* Display PRIMME params */
331: primme_display_params(*primme);
333: return(0);
334: }
338: PetscErrorCode EPSSetFromOptions_PRIMME(EPS eps)
339: {
341: EPS_PRIMME *ops = (EPS_PRIMME *)eps->data;
342: PetscInt op;
343: PetscTruth flg;
346:
347: PetscOptionsHead("PRIMME options");
349: op = ops->primme.maxBlockSize;
350: PetscOptionsInt("-eps_primme_block_size"," maximum block size","EPSPRIMMESetBlockSize",op,&op,&flg);
351: if (flg) {EPSPRIMMESetBlockSize(eps,op);}
352: op = 0;
353: PetscOptionsEList("-eps_primme_method","set method for solving the eigenproblem",
354: "EPSPRIMMESetMethod",methodList,15,methodList[1],&op,&flg);
355: if (flg) {EPSPRIMMESetMethod(eps, methodN[op]);}
356: PetscOptionsEList("-eps_primme_precond","set preconditioner type",
357: "EPSPRIMMESetPrecond",precondList,2,precondList[0],&op,&flg);
358: if (flg) {EPSPRIMMESetPrecond(eps, precondN[op]);}
359:
360: PetscOptionsTail();
361:
362: return(0);
363: }
368: PetscErrorCode EPSPRIMMESetBlockSize_PRIMME(EPS eps,PetscInt bs)
369: {
370: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
374: if (bs == PETSC_DEFAULT) ops->primme.maxBlockSize = 1;
375: else if (bs <= 0) {
376: SETERRQ(1, "PRIMME: wrong block size");
377: } else ops->primme.maxBlockSize = bs;
378:
379: return(0);
380: }
385: /*@
386: EPSPRIMMESetBlockSize - The maximum block size the code will try to use.
387: The user should set
388: this based on the architecture specifics of the target computer,
389: as well as any a priori knowledge of multiplicities. The code does
390: NOT require BlockSize > 1 to find multiple eigenvalues. For some
391: methods, keeping BlockSize = 1 yields the best overall performance.
393: Collective on EPS
395: Input Parameters:
396: + eps - the eigenproblem solver context
397: - bs - block size
399: Options Database Key:
400: . -eps_primme_block_size - Sets the max allowed block size value
402: Notes:
403: If the block size is not set, the value established by primme_initialize
404: is used.
406: Level: advanced
407: .seealso: EPSPRIMMEGetBlockSize()
408: @*/
409: PetscErrorCode EPSPRIMMESetBlockSize(EPS eps,PetscInt bs)
410: {
411: PetscErrorCode ierr, (*f)(EPS,PetscInt);
414:
416: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMESetBlockSize_C",(void (**)())&f);
417: if (f) {
418: (*f)(eps,bs);
419: }
420:
421: return(0);
422: }
427: PetscErrorCode EPSPRIMMEGetBlockSize_PRIMME(EPS eps,PetscInt *bs)
428: {
429: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
433: if (bs) *bs = ops->primme.maxBlockSize;
434:
435: return(0);
436: }
441: /*@
442: EPSPRIMMEGetBlockSize - Get the maximum block size the code will try to use.
443:
444: Collective on EPS
446: Input Parameters:
447: . eps - the eigenproblem solver context
448:
449: Output Parameters:
450: . bs - returned block size
452: Level: advanced
453: .seealso: EPSPRIMMESetBlockSize()
454: @*/
455: PetscErrorCode EPSPRIMMEGetBlockSize(EPS eps,PetscInt *bs)
456: {
457: PetscErrorCode ierr, (*f)(EPS,PetscInt*);
460:
462: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMEGetBlockSize_C",(void (**)())&f);
463: if (f) {
464: (*f)(eps,bs);
465: }
466:
467: return(0);
468: }
473: PetscErrorCode EPSPRIMMESetMethod_PRIMME(EPS eps, EPSPRIMMEMethod method)
474: {
475: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
479: if (method == PETSC_DEFAULT) ops->method = DEFAULT_MIN_TIME;
480: else ops->method = (primme_preset_method)method;
482: return(0);
483: }
488: /*@
489: EPSPRIMMESetMethod - Sets the method for the PRIMME library.
491: Collective on EPS
493: Input Parameters:
494: + eps - the eigenproblem solver context
495: - method - method that will be used by PRIMME. It must be one of:
496: EPSPRIMME_DYNAMIC, EPSPRIMME_DEFAULT_MIN_TIME(EPSPRIMME_JDQMR_ETOL),
497: EPSPRIMME_DEFAULT_MIN_MATVECS(EPSPRIMME_GD_OLSEN_PLUSK), EPSPRIMME_ARNOLDI,
498: EPSPRIMME_GD, EPSPRIMME_GD_PLUSK, EPSPRIMME_GD_OLSEN_PLUSK,
499: EPSPRIMME_JD_OLSEN_PLUSK, EPSPRIMME_RQI, EPSPRIMME_JDQR, EPSPRIMME_JDQMR,
500: EPSPRIMME_JDQMR_ETOL, EPSPRIMME_SUBSPACE_ITERATION,
501: EPSPRIMME_LOBPCG_ORTHOBASIS, EPSPRIMME_LOBPCG_ORTHOBASISW
503: Options Database Key:
504: . -eps_primme_set_method - Sets the method for the PRIMME library.
506: Note:
507: If not set, the method defaults to EPSPRIMME_DEFAULT_MIN_TIME.
509: Level: advanced
511: .seealso: EPSPRIMMEGetMethod(), EPSPRIMMEMethod
512: @*/
513: PetscErrorCode EPSPRIMMESetMethod(EPS eps, EPSPRIMMEMethod method)
514: {
515: PetscErrorCode ierr, (*f)(EPS,EPSPRIMMEMethod);
518:
520: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMESetMethod_C",(void (**)())&f);
521: if (f) {
522: (*f)(eps,method);
523: }
524:
525: return(0);
526: }
531: PetscErrorCode EPSPRIMMEGetMethod_PRIMME(EPS eps, EPSPRIMMEMethod *method)
532: {
533: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
537: if (method)
538: *method = (EPSPRIMMEMethod)ops->method;
540: return(0);
541: }
546: /*@C
547: EPSPRIMMEGetMethod - Gets the method for the PRIMME library.
549: Mon Collective on EPS
551: Input Parameters:
552: . eps - the eigenproblem solver context
553:
554: Output Parameters:
555: . method - method that will be used by PRIMME. It must be one of:
556: EPSPRIMME_DYNAMIC, EPSPRIMME_DEFAULT_MIN_TIME(EPSPRIMME_JDQMR_ETOL),
557: EPSPRIMME_DEFAULT_MIN_MATVECS(EPSPRIMME_GD_OLSEN_PLUSK), EPSPRIMME_ARNOLDI,
558: EPSPRIMME_GD, EPSPRIMME_GD_PLUSK, EPSPRIMME_GD_OLSEN_PLUSK,
559: EPSPRIMME_JD_OLSEN_PLUSK, EPSPRIMME_RQI, EPSPRIMME_JDQR, EPSPRIMME_JDQMR,
560: EPSPRIMME_JDQMR_ETOL, EPSPRIMME_SUBSPACE_ITERATION,
561: EPSPRIMME_LOBPCG_ORTHOBASIS, EPSPRIMME_LOBPCG_ORTHOBASISW
563: Level: advanced
565: .seealso: EPSPRIMMESetMethod(), EPSPRIMMEMethod
566: @*/
567: PetscErrorCode EPSPRIMMEGetMethod(EPS eps, EPSPRIMMEMethod *method)
568: {
569: PetscErrorCode ierr, (*f)(EPS,EPSPRIMMEMethod*);
572:
574: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMEGetMethod_C",(void (**)())&f);
575: if (f) {
576: (*f)(eps,method);
577: }
578:
579: return(0);
580: }
586: PetscErrorCode EPSPRIMMESetPrecond_PRIMME(EPS eps, EPSPRIMMEPrecond precond)
587: {
588: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
592: if (precond == EPSPRIMME_NONE) ops->primme.correctionParams.precondition = 0;
593: else ops->primme.correctionParams.precondition = 1;
594:
595: return(0);
596: }
601: /*@
602: EPSPRIMMESetPrecond - Sets the preconditioner to be used in the PRIMME
603: library. Currently, only diagonal preconditioning is supported.
605: Collective on EPS
607: Input Parameters:
608: + eps - the eigenproblem solver context
609: - precond - either EPSPRIMME_NONE (no preconditioning) or EPSPRIMME_DIAGONAL
610: (diagonal preconditioning)
612: Options Database Key:
613: . -eps_primme_precond - Sets either 'none' or 'diagonal' preconditioner
615: Note:
616: The default is no preconditioning.
617:
618: Level: advanced
620: .seealso: EPSPRIMMEGetPrecond(), EPSPRIMMEPrecond
621: @*/
622: PetscErrorCode EPSPRIMMESetPrecond(EPS eps, EPSPRIMMEPrecond precond)
623: {
624: PetscErrorCode ierr, (*f)(EPS, EPSPRIMMEPrecond);
627:
629: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMESetPrecond_C",(void (**)())&f);
630: if (f) {
631: (*f)(eps, precond);
632: }
633:
634: return(0);
635: }
641: PetscErrorCode EPSPRIMMEGetPrecond_PRIMME(EPS eps, EPSPRIMMEPrecond *precond)
642: {
643: EPS_PRIMME *ops = (EPS_PRIMME *) eps->data;
647: if (precond)
648: *precond = ops->primme.correctionParams.precondition ? EPSPRIMME_DIAGONAL : EPSPRIMME_NONE;
649:
650: return(0);
651: }
656: /*@C
657: EPSPRIMMEGetPrecond - Gets the preconditioner to be used in the PRIMME
658: library.
660: Collective on EPS
662: Input Parameters:
663: . eps - the eigenproblem solver context
664:
665: Output Parameters:
666: . precond - either EPSPRIMME_NONE (no preconditioning) or EPSPRIMME_DIAGONAL
667: (diagonal preconditioning)
669: Level: advanced
671: .seealso: EPSPRIMMESetPrecond(), EPSPRIMMEPrecond
672: @*/
673: PetscErrorCode EPSPRIMMEGetPrecond(EPS eps, EPSPRIMMEPrecond *precond)
674: {
675: PetscErrorCode ierr, (*f)(EPS, EPSPRIMMEPrecond*);
678:
680: PetscObjectQueryFunction((PetscObject)eps,"EPSPRIMMEGetPrecond_C",(void (**)())&f);
681: if (f) {
682: (*f)(eps, precond);
683: }
684:
685: return(0);
686: }
692: PetscErrorCode EPSCreate_PRIMME(EPS eps)
693: {
695: EPS_PRIMME *primme;
698:
699: PetscNew(EPS_PRIMME,&primme);
700: PetscLogObjectMemory(eps,sizeof(EPS_PRIMME));
701: eps->data = (void *) primme;
702: eps->ops->solve = EPSSolve_PRIMME;
703: eps->ops->setup = EPSSetUp_PRIMME;
704: eps->ops->setfromoptions = EPSSetFromOptions_PRIMME;
705: eps->ops->destroy = EPSDestroy_PRIMME;
706: eps->ops->view = EPSView_PRIMME;
707: eps->ops->backtransform = EPSBackTransform_Default;
708: eps->ops->computevectors = EPSComputeVectors_Default;
710: primme_initialize(&primme->primme);
711: primme->primme.matrixMatvec = multMatvec_PRIMME;
712: primme->primme.globalSumDouble = par_GlobalSumDouble;
713: primme->method = (primme_preset_method)EPSPRIMME_DEFAULT_MIN_TIME;
714: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMESetBlockSize_C","EPSPRIMMESetBlockSize_PRIMME",EPSPRIMMESetBlockSize_PRIMME);
715: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMESetMethod_C","EPSPRIMMESetMethod_PRIMME",EPSPRIMMESetMethod_PRIMME);
716: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMESetPrecond_C","EPSPRIMMESetPrecond_PRIMME",EPSPRIMMESetPrecond_PRIMME);
717:
718: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMEGetBlockSize_C","EPSPRIMMEGetBlockSize_PRIMME",EPSPRIMMEGetBlockSize_PRIMME);
719: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMEGetMethod_C","EPSPRIMMEGetMethod_PRIMME",EPSPRIMMEGetMethod_PRIMME);
720: PetscObjectComposeFunctionDynamic((PetscObject)eps,"EPSPRIMMEGetPrecond_C","EPSPRIMMEGetPrecond_PRIMME",EPSPRIMMEGetPrecond_PRIMME);
722: eps->which = EPS_LARGEST_REAL;
724: return(0);
725: }