Actual source code: primme.c

slepc-3.6.1 2015-09-03
Report Typos and Errors
  1: /*
  2:    This file implements a wrapper to the PRIMME package

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain

  8:    This file is part of SLEPc.

 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 <slepc/private/epsimpl.h>    /*I "slepceps.h" I*/

 26: PetscErrorCode EPSSolve_PRIMME(EPS);

 28: EXTERN_C_BEGIN
 29: #include <primme.h>
 30: EXTERN_C_END

 32: typedef struct {
 33:   primme_params primme;           /* param struc */
 34:   primme_preset_method method;    /* primme method */
 35:   Mat       A;                    /* problem matrix */
 36:   EPS       eps;                  /* EPS current context */
 37:   KSP       ksp;                  /* linear solver and preconditioner */
 38:   Vec       x,y;                  /* auxiliary vectors */
 39:   PetscReal target;               /* a copy of eps's target */
 40: } EPS_PRIMME;

 42: EPSPRIMMEMethod methodN[] = {
 43:   EPS_PRIMME_DYNAMIC,
 44:   EPS_PRIMME_DEFAULT_MIN_TIME,
 45:   EPS_PRIMME_DEFAULT_MIN_MATVECS,
 46:   EPS_PRIMME_ARNOLDI,
 47:   EPS_PRIMME_GD,
 48:   EPS_PRIMME_GD_PLUSK,
 49:   EPS_PRIMME_GD_OLSEN_PLUSK,
 50:   EPS_PRIMME_JD_OLSEN_PLUSK,
 51:   EPS_PRIMME_RQI,
 52:   EPS_PRIMME_JDQR,
 53:   EPS_PRIMME_JDQMR,
 54:   EPS_PRIMME_JDQMR_ETOL,
 55:   EPS_PRIMME_SUBSPACE_ITERATION,
 56:   EPS_PRIMME_LOBPCG_ORTHOBASIS,
 57:   EPS_PRIMME_LOBPCG_ORTHOBASISW
 58: };

 60: static void multMatvec_PRIMME(void *in,void *out,int *blockSize,primme_params *primme);
 61: static void applyPreconditioner_PRIMME(void *in,void *out,int *blockSize,struct primme_params *primme);

 63: static void par_GlobalSumDouble(void *sendBuf,void *recvBuf,int *count,primme_params *primme)
 64: {
 66:   MPI_Allreduce((double*)sendBuf,(double*)recvBuf,*count,MPI_DOUBLE,MPI_SUM,PetscObjectComm((PetscObject)primme->commInfo));CHKERRABORT(PetscObjectComm((PetscObject)primme->commInfo),ierr);
 67: }

 71: PetscErrorCode EPSSetUp_PRIMME(EPS eps)
 72: {
 74:   PetscMPIInt    numProcs,procID;
 75:   EPS_PRIMME     *ops = (EPS_PRIMME*)eps->data;
 76:   primme_params  *primme = &ops->primme;
 77:   PetscBool      istrivial,flg;

 80:   MPI_Comm_size(PetscObjectComm((PetscObject)eps),&numProcs);
 81:   MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&procID);

 83:   /* Check some constraints and set some default values */
 84:   if (!eps->max_it) eps->max_it = PetscMax(1000,eps->n);
 85:   STGetOperators(eps->st,0,&ops->A);
 86:   if (!eps->ishermitian) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"PRIMME is only available for Hermitian problems");
 87:   if (eps->isgeneralized) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"PRIMME is not available for generalized problems");
 88:   if (eps->arbitrary) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Arbitrary selection of eigenpairs not supported in this solver");
 89:   if (!eps->which) eps->which = EPS_LARGEST_REAL;
 90:   if (eps->converged != EPSConvergedAbsolute) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"PRIMME only supports absolute convergence test");
 91:   RGIsTrivial(eps->rg,&istrivial);
 92:   if (!istrivial) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"This solver does not support region filtering");

 94:   STSetUp(eps->st);
 95:   PetscObjectTypeCompare((PetscObject)eps->st,STPRECOND,&flg);
 96:   if (!flg) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"PRIMME only works with STPRECOND");

 98:   /* Transfer SLEPc options to PRIMME options */
 99:   primme->n          = eps->n;
100:   primme->nLocal     = eps->nloc;
101:   primme->numEvals   = eps->nev;
102:   primme->matrix     = ops;
103:   primme->commInfo   = eps;
104:   primme->maxMatvecs = eps->max_it;
105:   primme->eps        = eps->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:eps->tol;
106:   primme->numProcs   = numProcs;
107:   primme->procID     = procID;
108:   primme->printLevel = 0;
109:   primme->correctionParams.precondition = 1;

111:   switch (eps->which) {
112:     case EPS_LARGEST_REAL:
113:       primme->target = primme_largest;
114:       break;
115:     case EPS_SMALLEST_REAL:
116:       primme->target = primme_smallest;
117:       break;
118:     case EPS_TARGET_MAGNITUDE:
119:     case EPS_TARGET_REAL:
120:       primme->target = primme_closest_abs;
121:       primme->numTargetShifts = 1;
122:       ops->target = PetscRealPart(eps->target);
123:       primme->targetShifts = &ops->target;
124:       break;
125:     default:
126:       SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"'which' value not supported by PRIMME");
127:       break;
128:   }

130:   if (primme_set_method(ops->method,primme) < 0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"PRIMME method not valid");

132:   /* If user sets ncv, maxBasisSize is modified. If not, ncv is set as maxBasisSize */
133:   if (eps->ncv) primme->maxBasisSize = eps->ncv;
134:   else eps->ncv = primme->maxBasisSize;
135:   if (eps->ncv < eps->nev+primme->maxBlockSize)  SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"PRIMME needs ncv >= nev+maxBlockSize");
136:   if (eps->mpd) { PetscInfo(eps,"Warning: parameter mpd ignored\n"); }

138:   if (eps->extraction) { PetscInfo(eps,"Warning: extraction type ignored\n"); }

140:   /* Set workspace */
141:   EPSAllocateSolution(eps,0);
142:   PetscObjectTypeCompare((PetscObject)eps->V,BVVECS,&flg);
143:   if (flg) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"This solver requires a BV with contiguous storage");

145:   /* Setup the preconditioner */
146:   ops->eps = eps;
147:   if (primme->correctionParams.precondition) {
148:     STGetKSP(eps->st,&ops->ksp);
149:     PetscObjectTypeCompare((PetscObject)ops->ksp,KSPPREONLY,&flg);
150:     if (!flg) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"PRIMME only works with KSPPREONLY");
151:     primme->preconditioner = NULL;
152:     primme->applyPreconditioner = applyPreconditioner_PRIMME;
153:   }

155:   /* Prepare auxiliary vectors */
156:   VecCreateMPIWithArray(PetscObjectComm((PetscObject)eps),1,eps->nloc,eps->n,NULL,&ops->x);
157:   VecCreateMPIWithArray(PetscObjectComm((PetscObject)eps),1,eps->nloc,eps->n,NULL,&ops->y);
158:   PetscLogObjectParent((PetscObject)eps,(PetscObject)ops->x);
159:   PetscLogObjectParent((PetscObject)eps,(PetscObject)ops->y);

161:   /* dispatch solve method */
162:   eps->ops->solve = EPSSolve_PRIMME;
163:   return(0);
164: }

168: PetscErrorCode EPSSolve_PRIMME(EPS eps)
169: {
171:   EPS_PRIMME     *ops = (EPS_PRIMME*)eps->data;
172:   PetscScalar    *a;
173:   Vec            v0;
174: #if defined(PETSC_USE_COMPLEX)
175:   PetscInt       i;
176:   PetscReal      *evals;
177: #endif

180:   /* Reset some parameters left from previous runs */
181:   ops->primme.aNorm    = 1.0;
182:   ops->primme.initSize = eps->nini;
183:   ops->primme.iseed[0] = -1;

185:   /* Call PRIMME solver */
186:   BVGetColumn(eps->V,0,&v0);
187:   VecGetArray(v0,&a);
188: #if !defined(PETSC_USE_COMPLEX)
189:   dprimme(eps->eigr,a,eps->errest,&ops->primme);
190:   if (ierr) SETERRQ1(PetscObjectComm((PetscObject)eps),PETSC_ERR_LIB,"PRIMME library failed with error code=%d",ierr);
191: #else
192:   /* PRIMME returns real eigenvalues, but SLEPc works with complex ones */
193:   PetscMalloc1(eps->ncv,&evals);
194:   zprimme(evals,(Complex_Z*)a,eps->errest,&ops->primme);
195:   if (ierr) SETERRQ1(PetscObjectComm((PetscObject)eps),PETSC_ERR_LIB,"PRIMME library failed with error code=%d",ierr);
196:   for (i=0;i<eps->ncv;i++) eps->eigr[i] = evals[i];
197:   PetscFree(evals);
198: #endif
199:   VecRestoreArray(v0,&a);
200:   BVRestoreColumn(eps->V,0,&v0);

202:   eps->nconv  = ops->primme.initSize >= 0 ? ops->primme.initSize : 0;
203:   eps->reason = eps->ncv >= eps->nev ? EPS_CONVERGED_TOL: EPS_DIVERGED_ITS;
204:   eps->its    = ops->primme.stats.numOuterIterations;
205:   return(0);
206: }

210: static void multMatvec_PRIMME(void *in,void *out,int *blockSize,primme_params *primme)
211: {
213:   PetscInt       i,N = primme->n;
214:   EPS_PRIMME     *ops = (EPS_PRIMME*)primme->matrix;
215:   Vec            x = ops->x,y = ops->y;
216:   Mat            A = ops->A;

219:   for (i=0;i<*blockSize;i++) {
220:     /* build vectors using 'in' an 'out' workspace */
221:     VecPlaceArray(x,(PetscScalar*)in+N*i);CHKERRABORT(PetscObjectComm((PetscObject)A),ierr);
222:     VecPlaceArray(y,(PetscScalar*)out+N*i);CHKERRABORT(PetscObjectComm((PetscObject)A),ierr);

224:     MatMult(A,x,y);CHKERRABORT(PetscObjectComm((PetscObject)A),ierr);

226:     VecResetArray(x);CHKERRABORT(PetscObjectComm((PetscObject)A),ierr);
227:     VecResetArray(y);CHKERRABORT(PetscObjectComm((PetscObject)A),ierr);
228:   }
229:   PetscFunctionReturnVoid();
230: }

234: static void applyPreconditioner_PRIMME(void *in,void *out,int *blockSize,struct primme_params *primme)
235: {
237:   PetscInt       i,N = primme->n;
238:   EPS_PRIMME     *ops = (EPS_PRIMME*)primme->matrix;
239:   Vec            x = ops->x,y = ops->y;

242:   for (i=0;i<*blockSize;i++) {
243:     /* build vectors using 'in' an 'out' workspace */
244:     VecPlaceArray(x,(PetscScalar*)in+N*i);CHKERRABORT(PetscObjectComm((PetscObject)ops->ksp),ierr);
245:     VecPlaceArray(y,(PetscScalar*)out+N*i);CHKERRABORT(PetscObjectComm((PetscObject)ops->ksp),ierr);
246:     KSPSolve(ops->ksp,x,y);CHKERRABORT(PetscObjectComm((PetscObject)ops->ksp),ierr);
247:     VecResetArray(x);CHKERRABORT(PetscObjectComm((PetscObject)ops->ksp),ierr);
248:     VecResetArray(y);CHKERRABORT(PetscObjectComm((PetscObject)ops->ksp),ierr);
249:   }
250:   PetscFunctionReturnVoid();
251: }

255: PetscErrorCode EPSReset_PRIMME(EPS eps)
256: {
258:   EPS_PRIMME     *ops = (EPS_PRIMME*)eps->data;

261:   primme_Free(&ops->primme);
262:   VecDestroy(&ops->x);
263:   VecDestroy(&ops->y);
264:   return(0);
265: }

269: PetscErrorCode EPSDestroy_PRIMME(EPS eps)
270: {

274:   PetscFree(eps->data);
275:   PetscObjectComposeFunction((PetscObject)eps,"EPSPRIMMESetBlockSize_C",NULL);
276:   PetscObjectComposeFunction((PetscObject)eps,"EPSPRIMMESetMethod_C",NULL);
277:   PetscObjectComposeFunction((PetscObject)eps,"EPSPRIMMEGetBlockSize_C",NULL);
278:   PetscObjectComposeFunction((PetscObject)eps,"EPSPRIMMEGetMethod_C",NULL);
279:   return(0);
280: }

284: PetscErrorCode EPSView_PRIMME(EPS eps,PetscViewer viewer)
285: {
286:   PetscErrorCode  ierr;
287:   PetscBool       isascii;
288:   primme_params   *primme = &((EPS_PRIMME*)eps->data)->primme;
289:   EPSPRIMMEMethod methodn;
290:   PetscMPIInt     rank;

293:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
294:   if (isascii) {
295:     PetscViewerASCIIPrintf(viewer,"  PRIMME: block size=%d\n",primme->maxBlockSize);
296:     EPSPRIMMEGetMethod(eps,&methodn);
297:     PetscViewerASCIIPrintf(viewer,"  PRIMME: solver method: %s\n",EPSPRIMMEMethods[methodn]);

299:     /* Display PRIMME params */
300:     MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank);
301:     if (!rank) primme_display_params(*primme);
302:   }
303:   return(0);
304: }

308: PetscErrorCode EPSSetFromOptions_PRIMME(PetscOptions *PetscOptionsObject,EPS eps)
309: {
310:   PetscErrorCode  ierr;
311:   EPS_PRIMME      *ctx = (EPS_PRIMME*)eps->data;
312:   PetscInt        bs;
313:   EPSPRIMMEMethod meth;
314:   PetscBool       flg;
315:   KSP             ksp;

318:   PetscOptionsHead(PetscOptionsObject,"EPS PRIMME Options");
319:   PetscOptionsInt("-eps_primme_block_size","Maximum block size","EPSPRIMMESetBlockSize",ctx->primme.maxBlockSize,&bs,&flg);
320:   if (flg) {
321:     EPSPRIMMESetBlockSize(eps,bs);
322:   }
323:   PetscOptionsEnum("-eps_primme_method","Method for solving the eigenproblem","EPSPRIMMESetMethod",EPSPRIMMEMethods,(PetscEnum)ctx->method,(PetscEnum*)&meth,&flg);
324:   if (flg) {
325:     EPSPRIMMESetMethod(eps,meth);
326:   }

328:   /* Set STPrecond as the default ST */
329:   if (!((PetscObject)eps->st)->type_name) {
330:     STSetType(eps->st,STPRECOND);
331:   }
332:   STPrecondSetKSPHasMat(eps->st,PETSC_TRUE);

334:   /* Set the default options of the KSP */
335:   STGetKSP(eps->st,&ksp);
336:   if (!((PetscObject)ksp)->type_name) {
337:     KSPSetType(ksp,KSPPREONLY);
338:   }
339:   PetscOptionsTail();
340:   return(0);
341: }

345: static PetscErrorCode EPSPRIMMESetBlockSize_PRIMME(EPS eps,PetscInt bs)
346: {
347:   EPS_PRIMME *ops = (EPS_PRIMME*)eps->data;

350:   if (bs == PETSC_DEFAULT) ops->primme.maxBlockSize = 1;
351:   else if (bs <= 0) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"PRIMME: block size must be positive");
352:   else ops->primme.maxBlockSize = bs;
353:   return(0);
354: }

358: /*@
359:    EPSPRIMMESetBlockSize - The maximum block size the code will try to use.
360:    The user should set
361:    this based on the architecture specifics of the target computer,
362:    as well as any a priori knowledge of multiplicities. The code does
363:    NOT require BlockSize > 1 to find multiple eigenvalues.  For some
364:    methods, keeping BlockSize = 1 yields the best overall performance.

366:    Logically Collective on EPS

368:    Input Parameters:
369: +  eps - the eigenproblem solver context
370: -  bs - block size

372:    Options Database Key:
373: .  -eps_primme_block_size - Sets the max allowed block size value

375:    Notes:
376:    If the block size is not set, the value established by primme_initialize
377:    is used.

379:    Level: advanced

381: .seealso: EPSPRIMMEGetBlockSize()
382: @*/
383: PetscErrorCode EPSPRIMMESetBlockSize(EPS eps,PetscInt bs)
384: {

390:   PetscTryMethod(eps,"EPSPRIMMESetBlockSize_C",(EPS,PetscInt),(eps,bs));
391:   return(0);
392: }

396: static PetscErrorCode EPSPRIMMEGetBlockSize_PRIMME(EPS eps,PetscInt *bs)
397: {
398:   EPS_PRIMME *ops = (EPS_PRIMME*)eps->data;

401:   if (bs) *bs = ops->primme.maxBlockSize;
402:   return(0);
403: }

407: /*@
408:    EPSPRIMMEGetBlockSize - Get the maximum block size the code will try to use.

410:    Not Collective

412:    Input Parameter:
413: .  eps - the eigenproblem solver context

415:    Output Parameter:
416: .  bs - returned block size

418:    Level: advanced

420: .seealso: EPSPRIMMESetBlockSize()
421: @*/
422: PetscErrorCode EPSPRIMMEGetBlockSize(EPS eps,PetscInt *bs)
423: {

428:   PetscTryMethod(eps,"EPSPRIMMEGetBlockSize_C",(EPS,PetscInt*),(eps,bs));
429:   return(0);
430: }

434: static PetscErrorCode EPSPRIMMESetMethod_PRIMME(EPS eps,EPSPRIMMEMethod method)
435: {
436:   EPS_PRIMME *ops = (EPS_PRIMME*)eps->data;

439:   if (method == PETSC_DEFAULT) ops->method = DEFAULT_MIN_TIME;
440:   else ops->method = (primme_preset_method)method;
441:   return(0);
442: }

446: /*@
447:    EPSPRIMMESetMethod - Sets the method for the PRIMME library.

449:    Logically Collective on EPS

451:    Input Parameters:
452: +  eps - the eigenproblem solver context
453: -  method - method that will be used by PRIMME. It must be one of:
454:     EPS_PRIMME_DYNAMIC, EPS_PRIMME_DEFAULT_MIN_TIME(EPS_PRIMME_JDQMR_ETOL),
455:     EPS_PRIMME_DEFAULT_MIN_MATVECS(EPS_PRIMME_GD_OLSEN_PLUSK), EPS_PRIMME_ARNOLDI,
456:     EPS_PRIMME_GD, EPS_PRIMME_GD_PLUSK, EPS_PRIMME_GD_OLSEN_PLUSK,
457:     EPS_PRIMME_JD_OLSEN_PLUSK, EPS_PRIMME_RQI, EPS_PRIMME_JDQR, EPS_PRIMME_JDQMR,
458:     EPS_PRIMME_JDQMR_ETOL, EPS_PRIMME_SUBSPACE_ITERATION,
459:     EPS_PRIMME_LOBPCG_ORTHOBASIS, EPS_PRIMME_LOBPCG_ORTHOBASISW

461:    Options Database Key:
462: .  -eps_primme_method - Sets the method for the PRIMME library (one of 
463:     'dynamic', 'default_min_time', 'default_min_matvecs', 'arnoldi',
464:     'gd', 'gd_plusk', 'gd_olsen_plusk', 'jd_olsen_plusk', 'rqi', 'jdqr', 'jdqmr',
465:     'jdqmr_etol', 'subspace_iteration', 'lobpcg_orthobasis', 'lobpcg_orthobasisw').

467:    Note:
468:    If not set, the method defaults to EPS_PRIMME_DEFAULT_MIN_TIME.

470:    Level: advanced

472: .seealso: EPSPRIMMEGetMethod(), EPSPRIMMEMethod
473: @*/
474: PetscErrorCode EPSPRIMMESetMethod(EPS eps,EPSPRIMMEMethod method)
475: {

481:   PetscTryMethod(eps,"EPSPRIMMESetMethod_C",(EPS,EPSPRIMMEMethod),(eps,method));
482:   return(0);
483: }

487: static PetscErrorCode EPSPRIMMEGetMethod_PRIMME(EPS eps,EPSPRIMMEMethod *method)
488: {
489:   EPS_PRIMME *ops = (EPS_PRIMME*)eps->data;

492:   if (method) *method = (EPSPRIMMEMethod)ops->method;
493:   return(0);
494: }

498: /*@
499:    EPSPRIMMEGetMethod - Gets the method for the PRIMME library.

501:    Not Collective

503:    Input Parameter:
504: .  eps - the eigenproblem solver context

506:    Output Parameter:
507: .  method - method that will be used by PRIMME, one of
508:     EPS_PRIMME_DYNAMIC, EPS_PRIMME_DEFAULT_MIN_TIME(EPS_PRIMME_JDQMR_ETOL),
509:     EPS_PRIMME_DEFAULT_MIN_MATVECS(EPS_PRIMME_GD_OLSEN_PLUSK), EPS_PRIMME_ARNOLDI,
510:     EPS_PRIMME_GD, EPS_PRIMME_GD_PLUSK, EPS_PRIMME_GD_OLSEN_PLUSK,
511:     EPS_PRIMME_JD_OLSEN_PLUSK, EPS_PRIMME_RQI, EPS_PRIMME_JDQR, EPS_PRIMME_JDQMR,
512:     EPS_PRIMME_JDQMR_ETOL, EPS_PRIMME_SUBSPACE_ITERATION,
513:     EPS_PRIMME_LOBPCG_ORTHOBASIS, EPS_PRIMME_LOBPCG_ORTHOBASISW

515:    Level: advanced

517: .seealso: EPSPRIMMESetMethod(), EPSPRIMMEMethod
518: @*/
519: PetscErrorCode EPSPRIMMEGetMethod(EPS eps,EPSPRIMMEMethod *method)
520: {

525:   PetscTryMethod(eps,"EPSPRIMMEGetMethod_C",(EPS,EPSPRIMMEMethod*),(eps,method));
526:   return(0);
527: }

531: PETSC_EXTERN PetscErrorCode EPSCreate_PRIMME(EPS eps)
532: {
534:   EPS_PRIMME     *primme;

537:   PetscNewLog(eps,&primme);
538:   eps->data = (void*)primme;

540:   eps->ops->setup          = EPSSetUp_PRIMME;
541:   eps->ops->setfromoptions = EPSSetFromOptions_PRIMME;
542:   eps->ops->destroy        = EPSDestroy_PRIMME;
543:   eps->ops->reset          = EPSReset_PRIMME;
544:   eps->ops->view           = EPSView_PRIMME;
545:   eps->ops->backtransform  = EPSBackTransform_Default;

547:   primme_initialize(&primme->primme);
548:   primme->primme.matrixMatvec = multMatvec_PRIMME;
549:   primme->primme.globalSumDouble = par_GlobalSumDouble;
550:   primme->method = (primme_preset_method)EPS_PRIMME_DEFAULT_MIN_TIME;

552:   PetscObjectComposeFunction((PetscObject)eps,"EPSPRIMMESetBlockSize_C",EPSPRIMMESetBlockSize_PRIMME);
553:   PetscObjectComposeFunction((PetscObject)eps,"EPSPRIMMESetMethod_C",EPSPRIMMESetMethod_PRIMME);
554:   PetscObjectComposeFunction((PetscObject)eps,"EPSPRIMMEGetBlockSize_C",EPSPRIMMEGetBlockSize_PRIMME);
555:   PetscObjectComposeFunction((PetscObject)eps,"EPSPRIMMEGetMethod_C",EPSPRIMMEGetMethod_PRIMME);
556:   return(0);
557: }