Actual source code: epsbasic.c
slepc-3.10.2 2019-02-11
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2018, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: Basic EPS routines
12: */
14: #include <slepc/private/epsimpl.h> /*I "slepceps.h" I*/
16: PetscFunctionList EPSList = 0;
17: PetscBool EPSRegisterAllCalled = PETSC_FALSE;
18: PetscClassId EPS_CLASSID = 0;
19: PetscLogEvent EPS_SetUp = 0,EPS_Solve = 0;
21: /*@
22: EPSCreate - Creates the default EPS context.
24: Collective on MPI_Comm
26: Input Parameter:
27: . comm - MPI communicator
29: Output Parameter:
30: . eps - location to put the EPS context
32: Note:
33: The default EPS type is EPSKRYLOVSCHUR
35: Level: beginner
37: .seealso: EPSSetUp(), EPSSolve(), EPSDestroy(), EPS
38: @*/
39: PetscErrorCode EPSCreate(MPI_Comm comm,EPS *outeps)
40: {
42: EPS eps;
46: *outeps = 0;
47: EPSInitializePackage();
48: SlepcHeaderCreate(eps,EPS_CLASSID,"EPS","Eigenvalue Problem Solver","EPS",comm,EPSDestroy,EPSView);
50: eps->max_it = 0;
51: eps->nev = 1;
52: eps->ncv = 0;
53: eps->mpd = 0;
54: eps->nini = 0;
55: eps->nds = 0;
56: eps->target = 0.0;
57: eps->tol = PETSC_DEFAULT;
58: eps->conv = EPS_CONV_REL;
59: eps->stop = EPS_STOP_BASIC;
60: eps->which = (EPSWhich)0;
61: eps->inta = 0.0;
62: eps->intb = 0.0;
63: eps->problem_type = (EPSProblemType)0;
64: eps->extraction = EPS_RITZ;
65: eps->balance = EPS_BALANCE_NONE;
66: eps->balance_its = 5;
67: eps->balance_cutoff = 1e-8;
68: eps->trueres = PETSC_FALSE;
69: eps->trackall = PETSC_FALSE;
70: eps->purify = PETSC_TRUE;
72: eps->converged = EPSConvergedRelative;
73: eps->convergeduser = NULL;
74: eps->convergeddestroy= NULL;
75: eps->stopping = EPSStoppingBasic;
76: eps->stoppinguser = NULL;
77: eps->stoppingdestroy = NULL;
78: eps->arbitrary = NULL;
79: eps->convergedctx = NULL;
80: eps->stoppingctx = NULL;
81: eps->arbitraryctx = NULL;
82: eps->numbermonitors = 0;
84: eps->st = NULL;
85: eps->ds = NULL;
86: eps->V = NULL;
87: eps->rg = NULL;
88: eps->D = NULL;
89: eps->IS = NULL;
90: eps->defl = NULL;
91: eps->eigr = NULL;
92: eps->eigi = NULL;
93: eps->errest = NULL;
94: eps->rr = NULL;
95: eps->ri = NULL;
96: eps->perm = NULL;
97: eps->nwork = 0;
98: eps->work = NULL;
99: eps->data = NULL;
101: eps->state = EPS_STATE_INITIAL;
102: eps->categ = EPS_CATEGORY_KRYLOV;
103: eps->nconv = 0;
104: eps->its = 0;
105: eps->nloc = 0;
106: eps->nrma = 0.0;
107: eps->nrmb = 0.0;
108: eps->useds = PETSC_FALSE;
109: eps->isgeneralized = PETSC_FALSE;
110: eps->ispositive = PETSC_FALSE;
111: eps->ishermitian = PETSC_FALSE;
112: eps->reason = EPS_CONVERGED_ITERATING;
114: PetscNewLog(eps,&eps->sc);
115: *outeps = eps;
116: return(0);
117: }
119: /*@C
120: EPSSetType - Selects the particular solver to be used in the EPS object.
122: Logically Collective on EPS
124: Input Parameters:
125: + eps - the eigensolver context
126: - type - a known method
128: Options Database Key:
129: . -eps_type <method> - Sets the method; use -help for a list
130: of available methods
132: Notes:
133: See "slepc/include/slepceps.h" for available methods. The default
134: is EPSKRYLOVSCHUR.
136: Normally, it is best to use the EPSSetFromOptions() command and
137: then set the EPS type from the options database rather than by using
138: this routine. Using the options database provides the user with
139: maximum flexibility in evaluating the different available methods.
140: The EPSSetType() routine is provided for those situations where it
141: is necessary to set the iterative solver independently of the command
142: line or options database.
144: Level: intermediate
146: .seealso: STSetType(), EPSType
147: @*/
148: PetscErrorCode EPSSetType(EPS eps,EPSType type)
149: {
150: PetscErrorCode ierr,(*r)(EPS);
151: PetscBool match;
157: PetscObjectTypeCompare((PetscObject)eps,type,&match);
158: if (match) return(0);
160: PetscFunctionListFind(EPSList,type,&r);
161: if (!r) SETERRQ1(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown EPS type given: %s",type);
163: if (eps->ops->destroy) { (*eps->ops->destroy)(eps); }
164: PetscMemzero(eps->ops,sizeof(struct _EPSOps));
166: eps->state = EPS_STATE_INITIAL;
167: PetscObjectChangeTypeName((PetscObject)eps,type);
168: (*r)(eps);
169: return(0);
170: }
172: /*@C
173: EPSGetType - Gets the EPS type as a string from the EPS object.
175: Not Collective
177: Input Parameter:
178: . eps - the eigensolver context
180: Output Parameter:
181: . name - name of EPS method
183: Level: intermediate
185: .seealso: EPSSetType()
186: @*/
187: PetscErrorCode EPSGetType(EPS eps,EPSType *type)
188: {
192: *type = ((PetscObject)eps)->type_name;
193: return(0);
194: }
196: /*@C
197: EPSRegister - Adds a method to the eigenproblem solver package.
199: Not Collective
201: Input Parameters:
202: + name - name of a new user-defined solver
203: - function - routine to create the solver context
205: Notes:
206: EPSRegister() may be called multiple times to add several user-defined solvers.
208: Sample usage:
209: .vb
210: EPSRegister("my_solver",MySolverCreate);
211: .ve
213: Then, your solver can be chosen with the procedural interface via
214: $ EPSSetType(eps,"my_solver")
215: or at runtime via the option
216: $ -eps_type my_solver
218: Level: advanced
220: .seealso: EPSRegisterAll()
221: @*/
222: PetscErrorCode EPSRegister(const char *name,PetscErrorCode (*function)(EPS))
223: {
227: EPSInitializePackage();
228: PetscFunctionListAdd(&EPSList,name,function);
229: return(0);
230: }
232: /*@
233: EPSReset - Resets the EPS context to the initial state (prior to setup)
234: and destroys any allocated Vecs and Mats.
236: Collective on EPS
238: Input Parameter:
239: . eps - eigensolver context obtained from EPSCreate()
241: Note:
242: This can be used when a problem of different matrix size wants to be solved.
243: All options that have previously been set are preserved, so in a next use
244: the solver configuration is the same, but new sizes for matrices and vectors
245: are allowed.
247: Level: advanced
249: .seealso: EPSDestroy()
250: @*/
251: PetscErrorCode EPSReset(EPS eps)
252: {
257: if (!eps) return(0);
258: if (eps->ops->reset) { (eps->ops->reset)(eps); }
259: if (eps->st) { STReset(eps->st); }
260: VecDestroy(&eps->D);
261: BVDestroy(&eps->V);
262: VecDestroyVecs(eps->nwork,&eps->work);
263: eps->nwork = 0;
264: eps->state = EPS_STATE_INITIAL;
265: return(0);
266: }
268: /*@
269: EPSDestroy - Destroys the EPS context.
271: Collective on EPS
273: Input Parameter:
274: . eps - eigensolver context obtained from EPSCreate()
276: Level: beginner
278: .seealso: EPSCreate(), EPSSetUp(), EPSSolve()
279: @*/
280: PetscErrorCode EPSDestroy(EPS *eps)
281: {
285: if (!*eps) return(0);
287: if (--((PetscObject)(*eps))->refct > 0) { *eps = 0; return(0); }
288: EPSReset(*eps);
289: if ((*eps)->ops->destroy) { (*(*eps)->ops->destroy)(*eps); }
290: if ((*eps)->eigr) {
291: PetscFree4((*eps)->eigr,(*eps)->eigi,(*eps)->errest,(*eps)->perm);
292: }
293: if ((*eps)->rr) {
294: PetscFree2((*eps)->rr,(*eps)->ri);
295: }
296: STDestroy(&(*eps)->st);
297: RGDestroy(&(*eps)->rg);
298: DSDestroy(&(*eps)->ds);
299: PetscFree((*eps)->sc);
300: /* just in case the initial vectors have not been used */
301: SlepcBasisDestroy_Private(&(*eps)->nds,&(*eps)->defl);
302: SlepcBasisDestroy_Private(&(*eps)->nini,&(*eps)->IS);
303: if ((*eps)->convergeddestroy) {
304: (*(*eps)->convergeddestroy)((*eps)->convergedctx);
305: }
306: EPSMonitorCancel(*eps);
307: PetscHeaderDestroy(eps);
308: return(0);
309: }
311: /*@
312: EPSSetTarget - Sets the value of the target.
314: Logically Collective on EPS
316: Input Parameters:
317: + eps - eigensolver context
318: - target - the value of the target
320: Options Database Key:
321: . -eps_target <scalar> - the value of the target
323: Notes:
324: The target is a scalar value used to determine the portion of the spectrum
325: of interest. It is used in combination with EPSSetWhichEigenpairs().
327: In the case of complex scalars, a complex value can be provided in the
328: command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
329: -eps_target 1.0+2.0i
331: Level: intermediate
333: .seealso: EPSGetTarget(), EPSSetWhichEigenpairs()
334: @*/
335: PetscErrorCode EPSSetTarget(EPS eps,PetscScalar target)
336: {
342: eps->target = target;
343: if (!eps->st) { EPSGetST(eps,&eps->st); }
344: STSetDefaultShift(eps->st,target);
345: return(0);
346: }
348: /*@
349: EPSGetTarget - Gets the value of the target.
351: Not Collective
353: Input Parameter:
354: . eps - eigensolver context
356: Output Parameter:
357: . target - the value of the target
359: Note:
360: If the target was not set by the user, then zero is returned.
362: Level: intermediate
364: .seealso: EPSSetTarget()
365: @*/
366: PetscErrorCode EPSGetTarget(EPS eps,PetscScalar* target)
367: {
371: *target = eps->target;
372: return(0);
373: }
375: /*@
376: EPSSetInterval - Defines the computational interval for spectrum slicing.
378: Logically Collective on EPS
380: Input Parameters:
381: + eps - eigensolver context
382: . inta - left end of the interval
383: - intb - right end of the interval
385: Options Database Key:
386: . -eps_interval <a,b> - set [a,b] as the interval of interest
388: Notes:
389: Spectrum slicing is a technique employed for computing all eigenvalues of
390: symmetric eigenproblems in a given interval. This function provides the
391: interval to be considered. It must be used in combination with EPS_ALL, see
392: EPSSetWhichEigenpairs().
394: In the command-line option, two values must be provided. For an open interval,
395: one can give an infinite, e.g., -eps_interval 1.0,inf or -eps_interval -inf,1.0.
396: An open interval in the programmatic interface can be specified with
397: PETSC_MAX_REAL and -PETSC_MAX_REAL.
399: Level: intermediate
401: .seealso: EPSGetInterval(), EPSSetWhichEigenpairs()
402: @*/
403: PetscErrorCode EPSSetInterval(EPS eps,PetscReal inta,PetscReal intb)
404: {
409: if (inta>intb) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be inta<intb");
410: if (eps->inta != inta || eps->intb != intb) {
411: eps->inta = inta;
412: eps->intb = intb;
413: eps->state = EPS_STATE_INITIAL;
414: }
415: return(0);
416: }
418: /*@
419: EPSGetInterval - Gets the computational interval for spectrum slicing.
421: Not Collective
423: Input Parameter:
424: . eps - eigensolver context
426: Output Parameters:
427: + inta - left end of the interval
428: - intb - right end of the interval
430: Level: intermediate
432: Note:
433: If the interval was not set by the user, then zeros are returned.
435: .seealso: EPSSetInterval()
436: @*/
437: PetscErrorCode EPSGetInterval(EPS eps,PetscReal* inta,PetscReal* intb)
438: {
441: if (inta) *inta = eps->inta;
442: if (intb) *intb = eps->intb;
443: return(0);
444: }
446: /*@
447: EPSSetST - Associates a spectral transformation object to the eigensolver.
449: Collective on EPS
451: Input Parameters:
452: + eps - eigensolver context obtained from EPSCreate()
453: - st - the spectral transformation object
455: Note:
456: Use EPSGetST() to retrieve the spectral transformation context (for example,
457: to free it at the end of the computations).
459: Level: advanced
461: .seealso: EPSGetST()
462: @*/
463: PetscErrorCode EPSSetST(EPS eps,ST st)
464: {
471: PetscObjectReference((PetscObject)st);
472: STDestroy(&eps->st);
473: eps->st = st;
474: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
475: return(0);
476: }
478: /*@
479: EPSGetST - Obtain the spectral transformation (ST) object associated
480: to the eigensolver object.
482: Not Collective
484: Input Parameters:
485: . eps - eigensolver context obtained from EPSCreate()
487: Output Parameter:
488: . st - spectral transformation context
490: Level: intermediate
492: .seealso: EPSSetST()
493: @*/
494: PetscErrorCode EPSGetST(EPS eps,ST *st)
495: {
501: if (!eps->st) {
502: STCreate(PetscObjectComm((PetscObject)eps),&eps->st);
503: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->st);
504: }
505: *st = eps->st;
506: return(0);
507: }
509: /*@
510: EPSSetBV - Associates a basis vectors object to the eigensolver.
512: Collective on EPS
514: Input Parameters:
515: + eps - eigensolver context obtained from EPSCreate()
516: - V - the basis vectors object
518: Note:
519: Use EPSGetBV() to retrieve the basis vectors context (for example,
520: to free them at the end of the computations).
522: Level: advanced
524: .seealso: EPSGetBV()
525: @*/
526: PetscErrorCode EPSSetBV(EPS eps,BV V)
527: {
534: PetscObjectReference((PetscObject)V);
535: BVDestroy(&eps->V);
536: eps->V = V;
537: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
538: return(0);
539: }
541: /*@
542: EPSGetBV - Obtain the basis vectors object associated to the eigensolver object.
544: Not Collective
546: Input Parameters:
547: . eps - eigensolver context obtained from EPSCreate()
549: Output Parameter:
550: . V - basis vectors context
552: Level: advanced
554: .seealso: EPSSetBV()
555: @*/
556: PetscErrorCode EPSGetBV(EPS eps,BV *V)
557: {
563: if (!eps->V) {
564: BVCreate(PetscObjectComm((PetscObject)eps),&eps->V);
565: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->V);
566: }
567: *V = eps->V;
568: return(0);
569: }
571: /*@
572: EPSSetRG - Associates a region object to the eigensolver.
574: Collective on EPS
576: Input Parameters:
577: + eps - eigensolver context obtained from EPSCreate()
578: - rg - the region object
580: Note:
581: Use EPSGetRG() to retrieve the region context (for example,
582: to free it at the end of the computations).
584: Level: advanced
586: .seealso: EPSGetRG()
587: @*/
588: PetscErrorCode EPSSetRG(EPS eps,RG rg)
589: {
596: PetscObjectReference((PetscObject)rg);
597: RGDestroy(&eps->rg);
598: eps->rg = rg;
599: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
600: return(0);
601: }
603: /*@
604: EPSGetRG - Obtain the region object associated to the eigensolver.
606: Not Collective
608: Input Parameters:
609: . eps - eigensolver context obtained from EPSCreate()
611: Output Parameter:
612: . rg - region context
614: Level: advanced
616: .seealso: EPSSetRG()
617: @*/
618: PetscErrorCode EPSGetRG(EPS eps,RG *rg)
619: {
625: if (!eps->rg) {
626: RGCreate(PetscObjectComm((PetscObject)eps),&eps->rg);
627: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->rg);
628: }
629: *rg = eps->rg;
630: return(0);
631: }
633: /*@
634: EPSSetDS - Associates a direct solver object to the eigensolver.
636: Collective on EPS
638: Input Parameters:
639: + eps - eigensolver context obtained from EPSCreate()
640: - ds - the direct solver object
642: Note:
643: Use EPSGetDS() to retrieve the direct solver context (for example,
644: to free it at the end of the computations).
646: Level: advanced
648: .seealso: EPSGetDS()
649: @*/
650: PetscErrorCode EPSSetDS(EPS eps,DS ds)
651: {
658: PetscObjectReference((PetscObject)ds);
659: DSDestroy(&eps->ds);
660: eps->ds = ds;
661: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
662: return(0);
663: }
665: /*@
666: EPSGetDS - Obtain the direct solver object associated to the eigensolver object.
668: Not Collective
670: Input Parameters:
671: . eps - eigensolver context obtained from EPSCreate()
673: Output Parameter:
674: . ds - direct solver context
676: Level: advanced
678: .seealso: EPSSetDS()
679: @*/
680: PetscErrorCode EPSGetDS(EPS eps,DS *ds)
681: {
687: if (!eps->ds) {
688: DSCreate(PetscObjectComm((PetscObject)eps),&eps->ds);
689: PetscLogObjectParent((PetscObject)eps,(PetscObject)eps->ds);
690: }
691: *ds = eps->ds;
692: return(0);
693: }
695: /*@
696: EPSIsGeneralized - Ask if the EPS object corresponds to a generalized
697: eigenvalue problem.
699: Not collective
701: Input Parameter:
702: . eps - the eigenproblem solver context
704: Output Parameter:
705: . is - the answer
707: Level: intermediate
709: .seealso: EPSIsHermitian(), EPSIsPositive()
710: @*/
711: PetscErrorCode EPSIsGeneralized(EPS eps,PetscBool* is)
712: {
716: *is = eps->isgeneralized;
717: return(0);
718: }
720: /*@
721: EPSIsHermitian - Ask if the EPS object corresponds to a Hermitian
722: eigenvalue problem.
724: Not collective
726: Input Parameter:
727: . eps - the eigenproblem solver context
729: Output Parameter:
730: . is - the answer
732: Level: intermediate
734: .seealso: EPSIsGeneralized(), EPSIsPositive()
735: @*/
736: PetscErrorCode EPSIsHermitian(EPS eps,PetscBool* is)
737: {
741: *is = eps->ishermitian;
742: return(0);
743: }
745: /*@
746: EPSIsPositive - Ask if the EPS object corresponds to an eigenvalue
747: problem type that requires a positive (semi-) definite matrix B.
749: Not collective
751: Input Parameter:
752: . eps - the eigenproblem solver context
754: Output Parameter:
755: . is - the answer
757: Level: intermediate
759: .seealso: EPSIsGeneralized(), EPSIsHermitian()
760: @*/
761: PetscErrorCode EPSIsPositive(EPS eps,PetscBool* is)
762: {
766: *is = eps->ispositive;
767: return(0);
768: }