Line data Source code
1 : /*
2 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3 : SLEPc - Scalable Library for Eigenvalue Problem Computations
4 : Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5 :
6 : This file is part of SLEPc.
7 : SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9 : */
10 : /*
11 : Basic NEP routines
12 : */
13 :
14 : #include <slepc/private/nepimpl.h> /*I "slepcnep.h" I*/
15 :
16 : /* Logging support */
17 : PetscClassId NEP_CLASSID = 0;
18 : PetscLogEvent NEP_SetUp = 0,NEP_Solve = 0,NEP_Refine = 0,NEP_FunctionEval = 0,NEP_JacobianEval = 0,NEP_Resolvent = 0,NEP_CISS_SVD = 0;
19 :
20 : /* List of registered NEP routines */
21 : PetscFunctionList NEPList = NULL;
22 : PetscBool NEPRegisterAllCalled = PETSC_FALSE;
23 :
24 : /* List of registered NEP monitors */
25 : PetscFunctionList NEPMonitorList = NULL;
26 : PetscFunctionList NEPMonitorCreateList = NULL;
27 : PetscFunctionList NEPMonitorDestroyList = NULL;
28 : PetscBool NEPMonitorRegisterAllCalled = PETSC_FALSE;
29 :
30 : /*@
31 : NEPCreate - Creates the default NEP context.
32 :
33 : Collective
34 :
35 : Input Parameter:
36 : . comm - MPI communicator
37 :
38 : Output Parameter:
39 : . outnep - location to put the NEP context
40 :
41 : Level: beginner
42 :
43 : .seealso: NEPSetUp(), NEPSolve(), NEPDestroy(), NEP
44 : @*/
45 102 : PetscErrorCode NEPCreate(MPI_Comm comm,NEP *outnep)
46 : {
47 102 : NEP nep;
48 :
49 102 : PetscFunctionBegin;
50 102 : PetscAssertPointer(outnep,2);
51 102 : *outnep = NULL;
52 102 : PetscCall(NEPInitializePackage());
53 102 : PetscCall(SlepcHeaderCreate(nep,NEP_CLASSID,"NEP","Nonlinear Eigenvalue Problem","NEP",comm,NEPDestroy,NEPView));
54 :
55 102 : nep->max_it = PETSC_DEFAULT;
56 102 : nep->nev = 1;
57 102 : nep->ncv = PETSC_DEFAULT;
58 102 : nep->mpd = PETSC_DEFAULT;
59 102 : nep->nini = 0;
60 102 : nep->target = 0.0;
61 102 : nep->tol = PETSC_DEFAULT;
62 102 : nep->conv = NEP_CONV_REL;
63 102 : nep->stop = NEP_STOP_BASIC;
64 102 : nep->which = (NEPWhich)0;
65 102 : nep->problem_type = (NEPProblemType)0;
66 102 : nep->refine = NEP_REFINE_NONE;
67 102 : nep->npart = 1;
68 102 : nep->rtol = PETSC_DEFAULT;
69 102 : nep->rits = PETSC_DEFAULT;
70 102 : nep->scheme = (NEPRefineScheme)0;
71 102 : nep->trackall = PETSC_FALSE;
72 102 : nep->twosided = PETSC_FALSE;
73 :
74 102 : nep->computefunction = NULL;
75 102 : nep->computejacobian = NULL;
76 102 : nep->functionctx = NULL;
77 102 : nep->jacobianctx = NULL;
78 102 : nep->converged = NEPConvergedRelative;
79 102 : nep->convergeduser = NULL;
80 102 : nep->convergeddestroy= NULL;
81 102 : nep->stopping = NEPStoppingBasic;
82 102 : nep->stoppinguser = NULL;
83 102 : nep->stoppingdestroy = NULL;
84 102 : nep->convergedctx = NULL;
85 102 : nep->stoppingctx = NULL;
86 102 : nep->numbermonitors = 0;
87 :
88 102 : nep->ds = NULL;
89 102 : nep->V = NULL;
90 102 : nep->W = NULL;
91 102 : nep->rg = NULL;
92 102 : nep->function = NULL;
93 102 : nep->function_pre = NULL;
94 102 : nep->jacobian = NULL;
95 102 : nep->A = NULL;
96 102 : nep->f = NULL;
97 102 : nep->nt = 0;
98 102 : nep->mstr = UNKNOWN_NONZERO_PATTERN;
99 102 : nep->P = NULL;
100 102 : nep->mstrp = UNKNOWN_NONZERO_PATTERN;
101 102 : nep->IS = NULL;
102 102 : nep->eigr = NULL;
103 102 : nep->eigi = NULL;
104 102 : nep->errest = NULL;
105 102 : nep->perm = NULL;
106 102 : nep->nwork = 0;
107 102 : nep->work = NULL;
108 102 : nep->data = NULL;
109 :
110 102 : nep->state = NEP_STATE_INITIAL;
111 102 : nep->nconv = 0;
112 102 : nep->its = 0;
113 102 : nep->n = 0;
114 102 : nep->nloc = 0;
115 102 : nep->nrma = NULL;
116 102 : nep->fui = (NEPUserInterface)0;
117 102 : nep->useds = PETSC_FALSE;
118 102 : nep->resolvent = NULL;
119 102 : nep->reason = NEP_CONVERGED_ITERATING;
120 :
121 102 : PetscCall(PetscNew(&nep->sc));
122 102 : *outnep = nep;
123 102 : PetscFunctionReturn(PETSC_SUCCESS);
124 : }
125 :
126 : /*@C
127 : NEPSetType - Selects the particular solver to be used in the NEP object.
128 :
129 : Logically Collective
130 :
131 : Input Parameters:
132 : + nep - the nonlinear eigensolver context
133 : - type - a known method
134 :
135 : Options Database Key:
136 : . -nep_type <method> - Sets the method; use -help for a list
137 : of available methods
138 :
139 : Notes:
140 : See "slepc/include/slepcnep.h" for available methods.
141 :
142 : Normally, it is best to use the NEPSetFromOptions() command and
143 : then set the NEP type from the options database rather than by using
144 : this routine. Using the options database provides the user with
145 : maximum flexibility in evaluating the different available methods.
146 : The NEPSetType() routine is provided for those situations where it
147 : is necessary to set the iterative solver independently of the command
148 : line or options database.
149 :
150 : Level: intermediate
151 :
152 : .seealso: NEPType
153 : @*/
154 108 : PetscErrorCode NEPSetType(NEP nep,NEPType type)
155 : {
156 108 : PetscErrorCode (*r)(NEP);
157 108 : PetscBool match;
158 :
159 108 : PetscFunctionBegin;
160 108 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
161 108 : PetscAssertPointer(type,2);
162 :
163 108 : PetscCall(PetscObjectTypeCompare((PetscObject)nep,type,&match));
164 108 : if (match) PetscFunctionReturn(PETSC_SUCCESS);
165 :
166 105 : PetscCall(PetscFunctionListFind(NEPList,type,&r));
167 105 : PetscCheck(r,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown NEP type given: %s",type);
168 :
169 105 : PetscTryTypeMethod(nep,destroy);
170 105 : PetscCall(PetscMemzero(nep->ops,sizeof(struct _NEPOps)));
171 :
172 105 : nep->state = NEP_STATE_INITIAL;
173 105 : PetscCall(PetscObjectChangeTypeName((PetscObject)nep,type));
174 105 : PetscCall((*r)(nep));
175 105 : PetscFunctionReturn(PETSC_SUCCESS);
176 : }
177 :
178 : /*@C
179 : NEPGetType - Gets the NEP type as a string from the NEP object.
180 :
181 : Not Collective
182 :
183 : Input Parameter:
184 : . nep - the eigensolver context
185 :
186 : Output Parameter:
187 : . type - name of NEP method
188 :
189 : Level: intermediate
190 :
191 : .seealso: NEPSetType()
192 : @*/
193 17 : PetscErrorCode NEPGetType(NEP nep,NEPType *type)
194 : {
195 17 : PetscFunctionBegin;
196 17 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
197 17 : PetscAssertPointer(type,2);
198 17 : *type = ((PetscObject)nep)->type_name;
199 17 : PetscFunctionReturn(PETSC_SUCCESS);
200 : }
201 :
202 : /*@C
203 : NEPRegister - Adds a method to the nonlinear eigenproblem solver package.
204 :
205 : Not Collective
206 :
207 : Input Parameters:
208 : + name - name of a new user-defined solver
209 : - function - routine to create the solver context
210 :
211 : Notes:
212 : NEPRegister() may be called multiple times to add several user-defined solvers.
213 :
214 : Example Usage:
215 : .vb
216 : NEPRegister("my_solver",MySolverCreate);
217 : .ve
218 :
219 : Then, your solver can be chosen with the procedural interface via
220 : $ NEPSetType(nep,"my_solver")
221 : or at runtime via the option
222 : $ -nep_type my_solver
223 :
224 : Level: advanced
225 :
226 : .seealso: NEPRegisterAll()
227 : @*/
228 515 : PetscErrorCode NEPRegister(const char *name,PetscErrorCode (*function)(NEP))
229 : {
230 515 : PetscFunctionBegin;
231 515 : PetscCall(NEPInitializePackage());
232 515 : PetscCall(PetscFunctionListAdd(&NEPList,name,function));
233 515 : PetscFunctionReturn(PETSC_SUCCESS);
234 : }
235 :
236 : /*@C
237 : NEPMonitorRegister - Adds NEP monitor routine.
238 :
239 : Not Collective
240 :
241 : Input Parameters:
242 : + name - name of a new monitor routine
243 : . vtype - a PetscViewerType for the output
244 : . format - a PetscViewerFormat for the output
245 : . monitor - monitor routine
246 : . create - creation routine, or NULL
247 : - destroy - destruction routine, or NULL
248 :
249 : Notes:
250 : NEPMonitorRegister() may be called multiple times to add several user-defined monitors.
251 :
252 : Example Usage:
253 : .vb
254 : NEPMonitorRegister("my_monitor",PETSCVIEWERASCII,PETSC_VIEWER_ASCII_INFO_DETAIL,MyMonitor,NULL,NULL);
255 : .ve
256 :
257 : Then, your monitor can be chosen with the procedural interface via
258 : $ NEPMonitorSetFromOptions(nep,"-nep_monitor_my_monitor","my_monitor",NULL)
259 : or at runtime via the option
260 : $ -nep_monitor_my_monitor
261 :
262 : Level: advanced
263 :
264 : .seealso: NEPMonitorRegisterAll()
265 : @*/
266 618 : PetscErrorCode NEPMonitorRegister(const char name[],PetscViewerType vtype,PetscViewerFormat format,PetscErrorCode (*monitor)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,PetscViewerAndFormat*),PetscErrorCode (*create)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**),PetscErrorCode (*destroy)(PetscViewerAndFormat**))
267 : {
268 618 : char key[PETSC_MAX_PATH_LEN];
269 :
270 618 : PetscFunctionBegin;
271 618 : PetscCall(NEPInitializePackage());
272 618 : PetscCall(SlepcMonitorMakeKey_Internal(name,vtype,format,key));
273 618 : PetscCall(PetscFunctionListAdd(&NEPMonitorList,key,monitor));
274 618 : if (create) PetscCall(PetscFunctionListAdd(&NEPMonitorCreateList,key,create));
275 618 : if (destroy) PetscCall(PetscFunctionListAdd(&NEPMonitorDestroyList,key,destroy));
276 618 : PetscFunctionReturn(PETSC_SUCCESS);
277 : }
278 :
279 : /*
280 : NEPReset_Problem - Destroys the problem matrices.
281 : */
282 192 : PetscErrorCode NEPReset_Problem(NEP nep)
283 : {
284 192 : PetscInt i;
285 :
286 192 : PetscFunctionBegin;
287 192 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
288 192 : PetscCall(MatDestroy(&nep->function));
289 192 : PetscCall(MatDestroy(&nep->function_pre));
290 192 : PetscCall(MatDestroy(&nep->jacobian));
291 192 : if (nep->fui==NEP_USER_INTERFACE_SPLIT) {
292 85 : PetscCall(MatDestroyMatrices(nep->nt,&nep->A));
293 332 : for (i=0;i<nep->nt;i++) PetscCall(FNDestroy(&nep->f[i]));
294 85 : PetscCall(PetscFree(nep->f));
295 85 : PetscCall(PetscFree(nep->nrma));
296 85 : if (nep->P) PetscCall(MatDestroyMatrices(nep->nt,&nep->P));
297 85 : nep->nt = 0;
298 : }
299 192 : PetscFunctionReturn(PETSC_SUCCESS);
300 : }
301 : /*@
302 : NEPReset - Resets the NEP context to the initial state (prior to setup)
303 : and destroys any allocated Vecs and Mats.
304 :
305 : Collective
306 :
307 : Input Parameter:
308 : . nep - eigensolver context obtained from NEPCreate()
309 :
310 : Level: advanced
311 :
312 : .seealso: NEPDestroy()
313 : @*/
314 121 : PetscErrorCode NEPReset(NEP nep)
315 : {
316 121 : PetscFunctionBegin;
317 121 : if (nep) PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
318 121 : if (!nep) PetscFunctionReturn(PETSC_SUCCESS);
319 121 : PetscTryTypeMethod(nep,reset);
320 121 : if (nep->refineksp) PetscCall(KSPReset(nep->refineksp));
321 121 : PetscCall(NEPReset_Problem(nep));
322 121 : PetscCall(BVDestroy(&nep->V));
323 121 : PetscCall(BVDestroy(&nep->W));
324 121 : PetscCall(VecDestroyVecs(nep->nwork,&nep->work));
325 121 : PetscCall(MatDestroy(&nep->resolvent));
326 121 : nep->nwork = 0;
327 121 : nep->state = NEP_STATE_INITIAL;
328 121 : PetscFunctionReturn(PETSC_SUCCESS);
329 : }
330 :
331 : /*@C
332 : NEPDestroy - Destroys the NEP context.
333 :
334 : Collective
335 :
336 : Input Parameter:
337 : . nep - eigensolver context obtained from NEPCreate()
338 :
339 : Level: beginner
340 :
341 : .seealso: NEPCreate(), NEPSetUp(), NEPSolve()
342 : @*/
343 102 : PetscErrorCode NEPDestroy(NEP *nep)
344 : {
345 102 : PetscFunctionBegin;
346 102 : if (!*nep) PetscFunctionReturn(PETSC_SUCCESS);
347 102 : PetscValidHeaderSpecific(*nep,NEP_CLASSID,1);
348 102 : if (--((PetscObject)*nep)->refct > 0) { *nep = NULL; PetscFunctionReturn(PETSC_SUCCESS); }
349 102 : PetscCall(NEPReset(*nep));
350 102 : PetscTryTypeMethod(*nep,destroy);
351 102 : if ((*nep)->eigr) PetscCall(PetscFree4((*nep)->eigr,(*nep)->eigi,(*nep)->errest,(*nep)->perm));
352 102 : PetscCall(RGDestroy(&(*nep)->rg));
353 102 : PetscCall(DSDestroy(&(*nep)->ds));
354 102 : PetscCall(KSPDestroy(&(*nep)->refineksp));
355 102 : PetscCall(PetscSubcommDestroy(&(*nep)->refinesubc));
356 102 : PetscCall(PetscFree((*nep)->sc));
357 : /* just in case the initial vectors have not been used */
358 102 : PetscCall(SlepcBasisDestroy_Private(&(*nep)->nini,&(*nep)->IS));
359 102 : if ((*nep)->convergeddestroy) PetscCall((*(*nep)->convergeddestroy)((*nep)->convergedctx));
360 102 : PetscCall(NEPMonitorCancel(*nep));
361 102 : PetscCall(PetscHeaderDestroy(nep));
362 102 : PetscFunctionReturn(PETSC_SUCCESS);
363 : }
364 :
365 : /*@
366 : NEPSetBV - Associates a basis vectors object to the nonlinear eigensolver.
367 :
368 : Collective
369 :
370 : Input Parameters:
371 : + nep - eigensolver context obtained from NEPCreate()
372 : - bv - the basis vectors object
373 :
374 : Note:
375 : Use NEPGetBV() to retrieve the basis vectors context (for example,
376 : to free it at the end of the computations).
377 :
378 : Level: advanced
379 :
380 : .seealso: NEPGetBV()
381 : @*/
382 0 : PetscErrorCode NEPSetBV(NEP nep,BV bv)
383 : {
384 0 : PetscFunctionBegin;
385 0 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
386 0 : PetscValidHeaderSpecific(bv,BV_CLASSID,2);
387 0 : PetscCheckSameComm(nep,1,bv,2);
388 0 : PetscCall(PetscObjectReference((PetscObject)bv));
389 0 : PetscCall(BVDestroy(&nep->V));
390 0 : nep->V = bv;
391 0 : PetscFunctionReturn(PETSC_SUCCESS);
392 : }
393 :
394 : /*@
395 : NEPGetBV - Obtain the basis vectors object associated to the nonlinear
396 : eigensolver object.
397 :
398 : Not Collective
399 :
400 : Input Parameters:
401 : . nep - eigensolver context obtained from NEPCreate()
402 :
403 : Output Parameter:
404 : . bv - basis vectors context
405 :
406 : Level: advanced
407 :
408 : .seealso: NEPSetBV()
409 : @*/
410 122 : PetscErrorCode NEPGetBV(NEP nep,BV *bv)
411 : {
412 122 : PetscFunctionBegin;
413 122 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
414 122 : PetscAssertPointer(bv,2);
415 122 : if (!nep->V) {
416 121 : PetscCall(BVCreate(PetscObjectComm((PetscObject)nep),&nep->V));
417 121 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)nep->V,(PetscObject)nep,0));
418 121 : PetscCall(PetscObjectSetOptions((PetscObject)nep->V,((PetscObject)nep)->options));
419 : }
420 122 : *bv = nep->V;
421 122 : PetscFunctionReturn(PETSC_SUCCESS);
422 : }
423 :
424 : /*@
425 : NEPSetRG - Associates a region object to the nonlinear eigensolver.
426 :
427 : Collective
428 :
429 : Input Parameters:
430 : + nep - eigensolver context obtained from NEPCreate()
431 : - rg - the region object
432 :
433 : Note:
434 : Use NEPGetRG() to retrieve the region context (for example,
435 : to free it at the end of the computations).
436 :
437 : Level: advanced
438 :
439 : .seealso: NEPGetRG()
440 : @*/
441 1 : PetscErrorCode NEPSetRG(NEP nep,RG rg)
442 : {
443 1 : PetscFunctionBegin;
444 1 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
445 1 : if (rg) {
446 1 : PetscValidHeaderSpecific(rg,RG_CLASSID,2);
447 1 : PetscCheckSameComm(nep,1,rg,2);
448 : }
449 1 : PetscCall(PetscObjectReference((PetscObject)rg));
450 1 : PetscCall(RGDestroy(&nep->rg));
451 1 : nep->rg = rg;
452 1 : PetscFunctionReturn(PETSC_SUCCESS);
453 : }
454 :
455 : /*@
456 : NEPGetRG - Obtain the region object associated to the
457 : nonlinear eigensolver object.
458 :
459 : Not Collective
460 :
461 : Input Parameters:
462 : . nep - eigensolver context obtained from NEPCreate()
463 :
464 : Output Parameter:
465 : . rg - region context
466 :
467 : Level: advanced
468 :
469 : .seealso: NEPSetRG()
470 : @*/
471 102 : PetscErrorCode NEPGetRG(NEP nep,RG *rg)
472 : {
473 102 : PetscFunctionBegin;
474 102 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
475 102 : PetscAssertPointer(rg,2);
476 102 : if (!nep->rg) {
477 101 : PetscCall(RGCreate(PetscObjectComm((PetscObject)nep),&nep->rg));
478 101 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)nep->rg,(PetscObject)nep,0));
479 101 : PetscCall(PetscObjectSetOptions((PetscObject)nep->rg,((PetscObject)nep)->options));
480 : }
481 102 : *rg = nep->rg;
482 102 : PetscFunctionReturn(PETSC_SUCCESS);
483 : }
484 :
485 : /*@
486 : NEPSetDS - Associates a direct solver object to the nonlinear eigensolver.
487 :
488 : Collective
489 :
490 : Input Parameters:
491 : + nep - eigensolver context obtained from NEPCreate()
492 : - ds - the direct solver object
493 :
494 : Note:
495 : Use NEPGetDS() to retrieve the direct solver context (for example,
496 : to free it at the end of the computations).
497 :
498 : Level: advanced
499 :
500 : .seealso: NEPGetDS()
501 : @*/
502 0 : PetscErrorCode NEPSetDS(NEP nep,DS ds)
503 : {
504 0 : PetscFunctionBegin;
505 0 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
506 0 : PetscValidHeaderSpecific(ds,DS_CLASSID,2);
507 0 : PetscCheckSameComm(nep,1,ds,2);
508 0 : PetscCall(PetscObjectReference((PetscObject)ds));
509 0 : PetscCall(DSDestroy(&nep->ds));
510 0 : nep->ds = ds;
511 0 : PetscFunctionReturn(PETSC_SUCCESS);
512 : }
513 :
514 : /*@
515 : NEPGetDS - Obtain the direct solver object associated to the
516 : nonlinear eigensolver object.
517 :
518 : Not Collective
519 :
520 : Input Parameters:
521 : . nep - eigensolver context obtained from NEPCreate()
522 :
523 : Output Parameter:
524 : . ds - direct solver context
525 :
526 : Level: advanced
527 :
528 : .seealso: NEPSetDS()
529 : @*/
530 80 : PetscErrorCode NEPGetDS(NEP nep,DS *ds)
531 : {
532 80 : PetscFunctionBegin;
533 80 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
534 80 : PetscAssertPointer(ds,2);
535 80 : if (!nep->ds) {
536 79 : PetscCall(DSCreate(PetscObjectComm((PetscObject)nep),&nep->ds));
537 79 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)nep->ds,(PetscObject)nep,0));
538 79 : PetscCall(PetscObjectSetOptions((PetscObject)nep->ds,((PetscObject)nep)->options));
539 : }
540 80 : *ds = nep->ds;
541 80 : PetscFunctionReturn(PETSC_SUCCESS);
542 : }
543 :
544 : /*@
545 : NEPRefineGetKSP - Obtain the ksp object used by the eigensolver
546 : object in the refinement phase.
547 :
548 : Collective
549 :
550 : Input Parameters:
551 : . nep - eigensolver context obtained from NEPCreate()
552 :
553 : Output Parameter:
554 : . ksp - ksp context
555 :
556 : Level: advanced
557 :
558 : .seealso: NEPSetRefine()
559 : @*/
560 106 : PetscErrorCode NEPRefineGetKSP(NEP nep,KSP *ksp)
561 : {
562 106 : MPI_Comm comm;
563 :
564 106 : PetscFunctionBegin;
565 106 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
566 106 : PetscAssertPointer(ksp,2);
567 106 : if (!nep->refineksp) {
568 102 : if (nep->npart>1) {
569 : /* Split in subcomunicators */
570 8 : PetscCall(PetscSubcommCreate(PetscObjectComm((PetscObject)nep),&nep->refinesubc));
571 8 : PetscCall(PetscSubcommSetNumber(nep->refinesubc,nep->npart));
572 8 : PetscCall(PetscSubcommSetType(nep->refinesubc,PETSC_SUBCOMM_CONTIGUOUS));
573 8 : PetscCall(PetscSubcommGetChild(nep->refinesubc,&comm));
574 94 : } else PetscCall(PetscObjectGetComm((PetscObject)nep,&comm));
575 102 : PetscCall(KSPCreate(comm,&nep->refineksp));
576 102 : PetscCall(PetscObjectIncrementTabLevel((PetscObject)nep->refineksp,(PetscObject)nep,0));
577 102 : PetscCall(PetscObjectSetOptions((PetscObject)nep->refineksp,((PetscObject)nep)->options));
578 102 : PetscCall(KSPSetOptionsPrefix(*ksp,((PetscObject)nep)->prefix));
579 102 : PetscCall(KSPAppendOptionsPrefix(*ksp,"nep_refine_"));
580 203 : PetscCall(KSPSetTolerances(nep->refineksp,SlepcDefaultTol(nep->rtol),PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT));
581 : }
582 106 : *ksp = nep->refineksp;
583 106 : PetscFunctionReturn(PETSC_SUCCESS);
584 : }
585 :
586 : /*@
587 : NEPSetTarget - Sets the value of the target.
588 :
589 : Logically Collective
590 :
591 : Input Parameters:
592 : + nep - eigensolver context
593 : - target - the value of the target
594 :
595 : Options Database Key:
596 : . -nep_target <scalar> - the value of the target
597 :
598 : Notes:
599 : The target is a scalar value used to determine the portion of the spectrum
600 : of interest. It is used in combination with NEPSetWhichEigenpairs().
601 :
602 : In the case of complex scalars, a complex value can be provided in the
603 : command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
604 : -nep_target 1.0+2.0i
605 :
606 : Level: intermediate
607 :
608 : .seealso: NEPGetTarget(), NEPSetWhichEigenpairs()
609 : @*/
610 89 : PetscErrorCode NEPSetTarget(NEP nep,PetscScalar target)
611 : {
612 89 : PetscFunctionBegin;
613 89 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
614 356 : PetscValidLogicalCollectiveScalar(nep,target,2);
615 89 : nep->target = target;
616 89 : PetscFunctionReturn(PETSC_SUCCESS);
617 : }
618 :
619 : /*@
620 : NEPGetTarget - Gets the value of the target.
621 :
622 : Not Collective
623 :
624 : Input Parameter:
625 : . nep - eigensolver context
626 :
627 : Output Parameter:
628 : . target - the value of the target
629 :
630 : Note:
631 : If the target was not set by the user, then zero is returned.
632 :
633 : Level: intermediate
634 :
635 : .seealso: NEPSetTarget()
636 : @*/
637 2 : PetscErrorCode NEPGetTarget(NEP nep,PetscScalar* target)
638 : {
639 2 : PetscFunctionBegin;
640 2 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
641 2 : PetscAssertPointer(target,2);
642 2 : *target = nep->target;
643 2 : PetscFunctionReturn(PETSC_SUCCESS);
644 : }
645 :
646 : /*@C
647 : NEPSetFunction - Sets the function to compute the nonlinear Function T(lambda)
648 : as well as the location to store the matrix.
649 :
650 : Collective
651 :
652 : Input Parameters:
653 : + nep - the NEP context
654 : . A - Function matrix
655 : . B - preconditioner matrix (usually same as A)
656 : . fun - Function evaluation routine (if NULL then NEP retains any
657 : previously set value)
658 : - ctx - [optional] user-defined context for private data for the Function
659 : evaluation routine (may be NULL) (if NULL then NEP retains any
660 : previously set value)
661 :
662 : Calling sequence of fun:
663 : $ PetscErrorCode fun(NEP nep,PetscScalar lambda,Mat T,Mat P,void *ctx)
664 : + nep - the NEP context
665 : . lambda - the scalar argument where T(.) must be evaluated
666 : . T - matrix that will contain T(lambda)
667 : . P - (optional) different matrix to build the preconditioner
668 : - ctx - (optional) user-defined context, as set by NEPSetFunction()
669 :
670 : Level: beginner
671 :
672 : .seealso: NEPGetFunction(), NEPSetJacobian()
673 : @*/
674 36 : PetscErrorCode NEPSetFunction(NEP nep,Mat A,Mat B,PetscErrorCode (*fun)(NEP nep,PetscScalar lambda,Mat T,Mat P,void *ctx),void *ctx)
675 : {
676 36 : PetscFunctionBegin;
677 36 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
678 36 : if (A) PetscValidHeaderSpecific(A,MAT_CLASSID,2);
679 36 : if (B) PetscValidHeaderSpecific(B,MAT_CLASSID,3);
680 36 : if (A) PetscCheckSameComm(nep,1,A,2);
681 36 : if (B) PetscCheckSameComm(nep,1,B,3);
682 :
683 36 : if (nep->state) PetscCall(NEPReset(nep));
684 31 : else if (nep->fui && nep->fui!=NEP_USER_INTERFACE_CALLBACK) PetscCall(NEPReset_Problem(nep));
685 :
686 36 : if (fun) nep->computefunction = fun;
687 36 : if (ctx) nep->functionctx = ctx;
688 36 : if (A) {
689 36 : PetscCall(PetscObjectReference((PetscObject)A));
690 36 : PetscCall(MatDestroy(&nep->function));
691 36 : nep->function = A;
692 : }
693 36 : if (B) {
694 36 : PetscCall(PetscObjectReference((PetscObject)B));
695 36 : PetscCall(MatDestroy(&nep->function_pre));
696 36 : nep->function_pre = B;
697 : }
698 36 : nep->fui = NEP_USER_INTERFACE_CALLBACK;
699 36 : nep->state = NEP_STATE_INITIAL;
700 36 : PetscFunctionReturn(PETSC_SUCCESS);
701 : }
702 :
703 : /*@C
704 : NEPGetFunction - Returns the Function matrix and optionally the user
705 : provided context for evaluating the Function.
706 :
707 : Not Collective
708 :
709 : Input Parameter:
710 : . nep - the nonlinear eigensolver context
711 :
712 : Output Parameters:
713 : + A - location to stash Function matrix (or NULL)
714 : . B - location to stash preconditioner matrix (or NULL)
715 : . fun - location to put Function function (or NULL)
716 : - ctx - location to stash Function context (or NULL)
717 :
718 : Level: advanced
719 :
720 : .seealso: NEPSetFunction()
721 : @*/
722 72 : PetscErrorCode NEPGetFunction(NEP nep,Mat *A,Mat *B,PetscErrorCode (**fun)(NEP,PetscScalar,Mat,Mat,void*),void **ctx)
723 : {
724 72 : PetscFunctionBegin;
725 72 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
726 72 : NEPCheckCallback(nep,1);
727 72 : if (A) *A = nep->function;
728 72 : if (B) *B = nep->function_pre;
729 72 : if (fun) *fun = nep->computefunction;
730 72 : if (ctx) *ctx = nep->functionctx;
731 72 : PetscFunctionReturn(PETSC_SUCCESS);
732 : }
733 :
734 : /*@C
735 : NEPSetJacobian - Sets the function to compute the Jacobian T'(lambda) as well
736 : as the location to store the matrix.
737 :
738 : Collective
739 :
740 : Input Parameters:
741 : + nep - the NEP context
742 : . A - Jacobian matrix
743 : . jac - Jacobian evaluation routine (if NULL then NEP retains any
744 : previously set value)
745 : - ctx - [optional] user-defined context for private data for the Jacobian
746 : evaluation routine (may be NULL) (if NULL then NEP retains any
747 : previously set value)
748 :
749 : Calling sequence of jac:
750 : $ PetscErrorCode jac(NEP nep,PetscScalar lambda,Mat J,void *ctx)
751 : + nep - the NEP context
752 : . lambda - the scalar argument where T'(.) must be evaluated
753 : . J - matrix that will contain T'(lambda)
754 : - ctx - (optional) user-defined context, as set by NEPSetJacobian()
755 :
756 : Level: beginner
757 :
758 : .seealso: NEPSetFunction(), NEPGetJacobian()
759 : @*/
760 33 : PetscErrorCode NEPSetJacobian(NEP nep,Mat A,PetscErrorCode (*jac)(NEP nep,PetscScalar lambda,Mat J,void *ctx),void *ctx)
761 : {
762 33 : PetscFunctionBegin;
763 33 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
764 33 : if (A) PetscValidHeaderSpecific(A,MAT_CLASSID,2);
765 33 : if (A) PetscCheckSameComm(nep,1,A,2);
766 :
767 33 : if (nep->state) PetscCall(NEPReset(nep));
768 33 : else if (nep->fui && nep->fui!=NEP_USER_INTERFACE_CALLBACK) PetscCall(NEPReset_Problem(nep));
769 :
770 33 : if (jac) nep->computejacobian = jac;
771 33 : if (ctx) nep->jacobianctx = ctx;
772 33 : if (A) {
773 33 : PetscCall(PetscObjectReference((PetscObject)A));
774 33 : PetscCall(MatDestroy(&nep->jacobian));
775 33 : nep->jacobian = A;
776 : }
777 33 : nep->fui = NEP_USER_INTERFACE_CALLBACK;
778 33 : nep->state = NEP_STATE_INITIAL;
779 33 : PetscFunctionReturn(PETSC_SUCCESS);
780 : }
781 :
782 : /*@C
783 : NEPGetJacobian - Returns the Jacobian matrix and optionally the user
784 : provided routine and context for evaluating the Jacobian.
785 :
786 : Not Collective
787 :
788 : Input Parameter:
789 : . nep - the nonlinear eigensolver context
790 :
791 : Output Parameters:
792 : + A - location to stash Jacobian matrix (or NULL)
793 : . jac - location to put Jacobian function (or NULL)
794 : - ctx - location to stash Jacobian context (or NULL)
795 :
796 : Level: advanced
797 :
798 : .seealso: NEPSetJacobian()
799 : @*/
800 0 : PetscErrorCode NEPGetJacobian(NEP nep,Mat *A,PetscErrorCode (**jac)(NEP,PetscScalar,Mat,void*),void **ctx)
801 : {
802 0 : PetscFunctionBegin;
803 0 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
804 0 : NEPCheckCallback(nep,1);
805 0 : if (A) *A = nep->jacobian;
806 0 : if (jac) *jac = nep->computejacobian;
807 0 : if (ctx) *ctx = nep->jacobianctx;
808 0 : PetscFunctionReturn(PETSC_SUCCESS);
809 : }
810 :
811 : /*@
812 : NEPSetSplitOperator - Sets the operator of the nonlinear eigenvalue problem
813 : in split form.
814 :
815 : Collective
816 :
817 : Input Parameters:
818 : + nep - the nonlinear eigensolver context
819 : . nt - number of terms in the split form
820 : . A - array of matrices
821 : . f - array of functions
822 : - str - structure flag for matrices
823 :
824 : Notes:
825 : The nonlinear operator is written as T(lambda) = sum_i A_i*f_i(lambda),
826 : for i=1,...,n. The derivative T'(lambda) can be obtained using the
827 : derivatives of f_i.
828 :
829 : The structure flag provides information about A_i's nonzero pattern
830 : (see MatStructure enum). If all matrices have the same pattern, then
831 : use SAME_NONZERO_PATTERN. If the patterns are different but contained
832 : in the pattern of the first one, then use SUBSET_NONZERO_PATTERN. If
833 : patterns are known to be different, use DIFFERENT_NONZERO_PATTERN.
834 : If set to UNKNOWN_NONZERO_PATTERN, the patterns will be compared to
835 : determine if they are equal.
836 :
837 : This function must be called before NEPSetUp(). If it is called again
838 : after NEPSetUp() then the NEP object is reset.
839 :
840 : Level: beginner
841 :
842 : .seealso: NEPGetSplitOperatorTerm(), NEPGetSplitOperatorInfo(), NEPSetSplitPreconditioner()
843 : @*/
844 85 : PetscErrorCode NEPSetSplitOperator(NEP nep,PetscInt nt,Mat A[],FN f[],MatStructure str)
845 : {
846 85 : PetscInt i,n=0,m,m0=0,mloc,nloc,mloc0=0;
847 :
848 85 : PetscFunctionBegin;
849 85 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
850 340 : PetscValidLogicalCollectiveInt(nep,nt,2);
851 85 : PetscCheck(nt>0,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Must have one or more terms, you have %" PetscInt_FMT,nt);
852 85 : PetscAssertPointer(A,3);
853 85 : PetscAssertPointer(f,4);
854 340 : PetscValidLogicalCollectiveEnum(nep,str,5);
855 :
856 332 : for (i=0;i<nt;i++) {
857 247 : PetscValidHeaderSpecific(A[i],MAT_CLASSID,3);
858 247 : PetscCheckSameComm(nep,1,A[i],3);
859 247 : PetscValidHeaderSpecific(f[i],FN_CLASSID,4);
860 247 : PetscCheckSameComm(nep,1,f[i],4);
861 247 : PetscCall(MatGetSize(A[i],&m,&n));
862 247 : PetscCall(MatGetLocalSize(A[i],&mloc,&nloc));
863 247 : PetscCheck(m==n,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONG,"A[%" PetscInt_FMT "] is a non-square matrix (%" PetscInt_FMT " rows, %" PetscInt_FMT " cols)",i,m,n);
864 247 : PetscCheck(mloc==nloc,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONG,"A[%" PetscInt_FMT "] does not have equal row and column local sizes (%" PetscInt_FMT ", %" PetscInt_FMT ")",i,mloc,nloc);
865 247 : if (!i) { m0 = m; mloc0 = mloc; }
866 247 : PetscCheck(m==m0,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_INCOMP,"Dimensions of A[%" PetscInt_FMT "] do not match with previous matrices (%" PetscInt_FMT ", %" PetscInt_FMT ")",i,m,m0);
867 247 : PetscCheck(mloc==mloc0,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_INCOMP,"Local dimensions of A[%" PetscInt_FMT "] do not match with previous matrices (%" PetscInt_FMT ", %" PetscInt_FMT ")",i,mloc,mloc0);
868 247 : PetscCall(PetscObjectReference((PetscObject)A[i]));
869 247 : PetscCall(PetscObjectReference((PetscObject)f[i]));
870 : }
871 :
872 85 : if (nep->state && (n!=nep->n || nloc!=nep->nloc)) PetscCall(NEPReset(nep));
873 71 : else PetscCall(NEPReset_Problem(nep));
874 :
875 : /* allocate space and copy matrices and functions */
876 85 : PetscCall(PetscMalloc1(nt,&nep->A));
877 332 : for (i=0;i<nt;i++) nep->A[i] = A[i];
878 85 : PetscCall(PetscMalloc1(nt,&nep->f));
879 332 : for (i=0;i<nt;i++) nep->f[i] = f[i];
880 85 : PetscCall(PetscCalloc1(nt,&nep->nrma));
881 85 : nep->nt = nt;
882 85 : nep->mstr = str;
883 85 : nep->fui = NEP_USER_INTERFACE_SPLIT;
884 85 : nep->state = NEP_STATE_INITIAL;
885 85 : PetscFunctionReturn(PETSC_SUCCESS);
886 : }
887 :
888 : /*@
889 : NEPGetSplitOperatorTerm - Gets the matrices and functions associated with
890 : the nonlinear operator in split form.
891 :
892 : Not Collective
893 :
894 : Input Parameters:
895 : + nep - the nonlinear eigensolver context
896 : - k - the index of the requested term (starting in 0)
897 :
898 : Output Parameters:
899 : + A - the matrix of the requested term
900 : - f - the function of the requested term
901 :
902 : Level: intermediate
903 :
904 : .seealso: NEPSetSplitOperator(), NEPGetSplitOperatorInfo()
905 : @*/
906 25 : PetscErrorCode NEPGetSplitOperatorTerm(NEP nep,PetscInt k,Mat *A,FN *f)
907 : {
908 25 : PetscFunctionBegin;
909 25 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
910 100 : PetscValidLogicalCollectiveInt(nep,k,2);
911 25 : NEPCheckSplit(nep,1);
912 25 : PetscCheck(k>=0 && k<nep->nt,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %" PetscInt_FMT,nep->nt-1);
913 25 : if (A) *A = nep->A[k];
914 25 : if (f) *f = nep->f[k];
915 25 : PetscFunctionReturn(PETSC_SUCCESS);
916 : }
917 :
918 : /*@
919 : NEPGetSplitOperatorInfo - Returns the number of terms of the split form of
920 : the nonlinear operator, as well as the structure flag for matrices.
921 :
922 : Not Collective
923 :
924 : Input Parameter:
925 : . nep - the nonlinear eigensolver context
926 :
927 : Output Parameters:
928 : + n - the number of terms passed in NEPSetSplitOperator()
929 : - str - the matrix structure flag passed in NEPSetSplitOperator()
930 :
931 : Level: intermediate
932 :
933 : .seealso: NEPSetSplitOperator(), NEPGetSplitOperatorTerm()
934 : @*/
935 9 : PetscErrorCode NEPGetSplitOperatorInfo(NEP nep,PetscInt *n,MatStructure *str)
936 : {
937 9 : PetscFunctionBegin;
938 9 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
939 9 : NEPCheckSplit(nep,1);
940 9 : if (n) *n = nep->nt;
941 9 : if (str) *str = nep->mstr;
942 9 : PetscFunctionReturn(PETSC_SUCCESS);
943 : }
944 :
945 : /*@
946 : NEPSetSplitPreconditioner - Sets an operator in split form from which
947 : to build the preconditioner to be used when solving the nonlinear
948 : eigenvalue problem in split form.
949 :
950 : Collective
951 :
952 : Input Parameters:
953 : + nep - the nonlinear eigensolver context
954 : . ntp - number of terms in the split preconditioner
955 : . P - array of matrices
956 : - strp - structure flag for matrices
957 :
958 : Notes:
959 : The matrix for the preconditioner is expressed as P(lambda) =
960 : sum_i P_i*f_i(lambda), for i=1,...,n, where the f_i functions
961 : are the same as in NEPSetSplitOperator(). It is not necessary to call
962 : this function. If it is not invoked, then the preconditioner is
963 : built from T(lambda), i.e., both matrices and functions passed in
964 : NEPSetSplitOperator().
965 :
966 : The structure flag provides information about P_i's nonzero pattern
967 : in the same way as in NEPSetSplitOperator().
968 :
969 : If the functions defining the preconditioner operator were different
970 : from the ones given in NEPSetSplitOperator(), then the split form
971 : cannot be used. Use the callback interface instead.
972 :
973 : Use ntp=0 to reset a previously set split preconditioner.
974 :
975 : Level: advanced
976 :
977 : .seealso: NEPGetSplitPreconditionerTerm(), NEPGetSplitPreconditionerInfo(), NEPSetSplitOperator()
978 : @*/
979 4 : PetscErrorCode NEPSetSplitPreconditioner(NEP nep,PetscInt ntp,Mat P[],MatStructure strp)
980 : {
981 4 : PetscInt i,n=0,m,m0=0,mloc,nloc,mloc0=0;
982 :
983 4 : PetscFunctionBegin;
984 4 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
985 16 : PetscValidLogicalCollectiveInt(nep,ntp,2);
986 4 : PetscCheck(ntp>=0,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Negative value of ntp = %" PetscInt_FMT,ntp);
987 4 : PetscCheck(nep->fui==NEP_USER_INTERFACE_SPLIT,PetscObjectComm((PetscObject)nep),PETSC_ERR_ORDER,"Must call NEPSetSplitOperator first");
988 4 : PetscCheck(ntp==0 || nep->nt==ntp,PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"The number of terms must be the same as in NEPSetSplitOperator()");
989 4 : if (ntp) PetscAssertPointer(P,3);
990 16 : PetscValidLogicalCollectiveEnum(nep,strp,4);
991 :
992 16 : for (i=0;i<ntp;i++) {
993 12 : PetscValidHeaderSpecific(P[i],MAT_CLASSID,3);
994 12 : PetscCheckSameComm(nep,1,P[i],3);
995 12 : PetscCall(MatGetSize(P[i],&m,&n));
996 12 : PetscCall(MatGetLocalSize(P[i],&mloc,&nloc));
997 12 : PetscCheck(m==n,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONG,"P[%" PetscInt_FMT "] is a non-square matrix (%" PetscInt_FMT " rows, %" PetscInt_FMT " cols)",i,m,n);
998 12 : PetscCheck(mloc==nloc,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONG,"P[%" PetscInt_FMT "] does not have equal row and column local sizes (%" PetscInt_FMT ", %" PetscInt_FMT ")",i,mloc,nloc);
999 12 : if (!i) { m0 = m; mloc0 = mloc; }
1000 12 : PetscCheck(m==m0,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_INCOMP,"Dimensions of P[%" PetscInt_FMT "] do not match with previous matrices (%" PetscInt_FMT ", %" PetscInt_FMT ")",i,m,m0);
1001 12 : PetscCheck(mloc==mloc0,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_INCOMP,"Local dimensions of P[%" PetscInt_FMT "] do not match with previous matrices (%" PetscInt_FMT ", %" PetscInt_FMT ")",i,mloc,mloc0);
1002 12 : PetscCall(PetscObjectReference((PetscObject)P[i]));
1003 : }
1004 :
1005 4 : PetscCheck(!nep->state,PetscObjectComm((PetscObject)nep),PETSC_ERR_ORDER,"To call this function after NEPSetUp(), you must call NEPSetSplitOperator() again");
1006 4 : if (nep->P) PetscCall(MatDestroyMatrices(nep->nt,&nep->P));
1007 :
1008 : /* allocate space and copy matrices */
1009 4 : if (ntp) {
1010 4 : PetscCall(PetscMalloc1(ntp,&nep->P));
1011 16 : for (i=0;i<ntp;i++) nep->P[i] = P[i];
1012 : }
1013 4 : nep->mstrp = strp;
1014 4 : nep->state = NEP_STATE_INITIAL;
1015 4 : PetscFunctionReturn(PETSC_SUCCESS);
1016 : }
1017 :
1018 : /*@
1019 : NEPGetSplitPreconditionerTerm - Gets the matrices associated with
1020 : the split preconditioner.
1021 :
1022 : Not Collective
1023 :
1024 : Input Parameters:
1025 : + nep - the nonlinear eigensolver context
1026 : - k - the index of the requested term (starting in 0)
1027 :
1028 : Output Parameter:
1029 : . P - the matrix of the requested term
1030 :
1031 : Level: advanced
1032 :
1033 : .seealso: NEPSetSplitPreconditioner(), NEPGetSplitPreconditionerInfo()
1034 : @*/
1035 4 : PetscErrorCode NEPGetSplitPreconditionerTerm(NEP nep,PetscInt k,Mat *P)
1036 : {
1037 4 : PetscFunctionBegin;
1038 4 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
1039 16 : PetscValidLogicalCollectiveInt(nep,k,2);
1040 4 : PetscAssertPointer(P,3);
1041 4 : NEPCheckSplit(nep,1);
1042 4 : PetscCheck(k>=0 && k<nep->nt,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %" PetscInt_FMT,nep->nt-1);
1043 4 : PetscCheck(nep->P,PetscObjectComm((PetscObject)nep),PETSC_ERR_ORDER,"You have not called NEPSetSplitPreconditioner()");
1044 4 : *P = nep->P[k];
1045 4 : PetscFunctionReturn(PETSC_SUCCESS);
1046 : }
1047 :
1048 : /*@
1049 : NEPGetSplitPreconditionerInfo - Returns the number of terms of the split
1050 : preconditioner, as well as the structure flag for matrices.
1051 :
1052 : Not Collective
1053 :
1054 : Input Parameter:
1055 : . nep - the nonlinear eigensolver context
1056 :
1057 : Output Parameters:
1058 : + n - the number of terms passed in NEPSetSplitPreconditioner()
1059 : - strp - the matrix structure flag passed in NEPSetSplitPreconditioner()
1060 :
1061 : Level: advanced
1062 :
1063 : .seealso: NEPSetSplitPreconditioner(), NEPGetSplitPreconditionerTerm()
1064 : @*/
1065 4 : PetscErrorCode NEPGetSplitPreconditionerInfo(NEP nep,PetscInt *n,MatStructure *strp)
1066 : {
1067 4 : PetscFunctionBegin;
1068 4 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
1069 4 : NEPCheckSplit(nep,1);
1070 4 : if (n) *n = nep->P? nep->nt: 0;
1071 4 : if (strp) *strp = nep->mstrp;
1072 4 : PetscFunctionReturn(PETSC_SUCCESS);
1073 : }
|