Actual source code: nepbasic.c

slepc-3.11.2 2019-07-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2019, 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 NEP routines
 12: */

 14: #include <slepc/private/nepimpl.h>      /*I "slepcnep.h" I*/

 16: PetscFunctionList NEPList = 0;
 17: PetscBool         NEPRegisterAllCalled = PETSC_FALSE;
 18: PetscClassId      NEP_CLASSID = 0;
 19: PetscLogEvent     NEP_SetUp = 0,NEP_Solve = 0,NEP_Refine = 0,NEP_FunctionEval = 0,NEP_JacobianEval = 0,NEP_DerivativesEval = 0,NEP_Resolvent = 0;

 21: /*@
 22:    NEPCreate - Creates the default NEP context.

 24:    Collective on MPI_Comm

 26:    Input Parameter:
 27: .  comm - MPI communicator

 29:    Output Parameter:
 30: .  nep - location to put the NEP context

 32:    Level: beginner

 34: .seealso: NEPSetUp(), NEPSolve(), NEPDestroy(), NEP
 35: @*/
 36: PetscErrorCode NEPCreate(MPI_Comm comm,NEP *outnep)
 37: {
 39:   NEP            nep;

 43:   *outnep = 0;
 44:   NEPInitializePackage();
 45:   SlepcHeaderCreate(nep,NEP_CLASSID,"NEP","Nonlinear Eigenvalue Problem","NEP",comm,NEPDestroy,NEPView);

 47:   nep->max_it          = 0;
 48:   nep->nev             = 1;
 49:   nep->ncv             = 0;
 50:   nep->mpd             = 0;
 51:   nep->nini            = 0;
 52:   nep->target          = 0.0;
 53:   nep->tol             = PETSC_DEFAULT;
 54:   nep->conv            = NEP_CONV_REL;
 55:   nep->stop            = NEP_STOP_BASIC;
 56:   nep->which           = (NEPWhich)0;
 57:   nep->problem_type    = (NEPProblemType)0;
 58:   nep->refine          = NEP_REFINE_NONE;
 59:   nep->npart           = 1;
 60:   nep->rtol            = PETSC_DEFAULT;
 61:   nep->rits            = PETSC_DEFAULT;
 62:   nep->scheme          = (NEPRefineScheme)0;
 63:   nep->trackall        = PETSC_FALSE;
 64:   nep->twosided        = PETSC_FALSE;

 66:   nep->computefunction = NULL;
 67:   nep->computejacobian = NULL;
 68:   nep->functionctx     = NULL;
 69:   nep->jacobianctx     = NULL;
 70:   nep->computederivatives = NULL;
 71:   nep->derivativesctx  = NULL;
 72:   nep->converged       = NEPConvergedRelative;
 73:   nep->convergeduser   = NULL;
 74:   nep->convergeddestroy= NULL;
 75:   nep->stopping        = NEPStoppingBasic;
 76:   nep->stoppinguser    = NULL;
 77:   nep->stoppingdestroy = NULL;
 78:   nep->convergedctx    = NULL;
 79:   nep->stoppingctx     = NULL;
 80:   nep->numbermonitors  = 0;

 82:   nep->ds              = NULL;
 83:   nep->V               = NULL;
 84:   nep->W               = NULL;
 85:   nep->rg              = NULL;
 86:   nep->function        = NULL;
 87:   nep->function_pre    = NULL;
 88:   nep->jacobian        = NULL;
 89:   nep->derivatives     = NULL;
 90:   nep->A               = NULL;
 91:   nep->f               = NULL;
 92:   nep->nt              = 0;
 93:   nep->mstr            = DIFFERENT_NONZERO_PATTERN;
 94:   nep->IS              = NULL;
 95:   nep->eigr            = NULL;
 96:   nep->eigi            = NULL;
 97:   nep->errest          = NULL;
 98:   nep->perm            = NULL;
 99:   nep->nwork           = 0;
100:   nep->work            = NULL;
101:   nep->data            = NULL;

103:   nep->state           = NEP_STATE_INITIAL;
104:   nep->nconv           = 0;
105:   nep->its             = 0;
106:   nep->n               = 0;
107:   nep->nloc            = 0;
108:   nep->nrma            = NULL;
109:   nep->fui             = (NEPUserInterface)0;
110:   nep->useds           = PETSC_FALSE;
111:   nep->hasts           = PETSC_FALSE;
112:   nep->resolvent       = NULL;
113:   nep->reason          = NEP_CONVERGED_ITERATING;

115:   PetscNewLog(nep,&nep->sc);
116:   *outnep = nep;
117:   return(0);
118: }

120: /*@C
121:    NEPSetType - Selects the particular solver to be used in the NEP object.

123:    Logically Collective on NEP

125:    Input Parameters:
126: +  nep      - the nonlinear eigensolver context
127: -  type     - a known method

129:    Options Database Key:
130: .  -nep_type <method> - Sets the method; use -help for a list
131:     of available methods

133:    Notes:
134:    See "slepc/include/slepcnep.h" for available methods.

136:    Normally, it is best to use the NEPSetFromOptions() command and
137:    then set the NEP 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 NEPSetType() 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: NEPType
147: @*/
148: PetscErrorCode NEPSetType(NEP nep,NEPType type)
149: {
150:   PetscErrorCode ierr,(*r)(NEP);
151:   PetscBool      match;


157:   PetscObjectTypeCompare((PetscObject)nep,type,&match);
158:   if (match) return(0);

160:   PetscFunctionListFind(NEPList,type,&r);
161:   if (!r) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown NEP type given: %s",type);

163:   if (nep->ops->destroy) { (*nep->ops->destroy)(nep); }
164:   PetscMemzero(nep->ops,sizeof(struct _NEPOps));

166:   nep->state = NEP_STATE_INITIAL;
167:   PetscObjectChangeTypeName((PetscObject)nep,type);
168:   (*r)(nep);
169:   return(0);
170: }

172: /*@C
173:    NEPGetType - Gets the NEP type as a string from the NEP object.

175:    Not Collective

177:    Input Parameter:
178: .  nep - the eigensolver context

180:    Output Parameter:
181: .  name - name of NEP method

183:    Level: intermediate

185: .seealso: NEPSetType()
186: @*/
187: PetscErrorCode NEPGetType(NEP nep,NEPType *type)
188: {
192:   *type = ((PetscObject)nep)->type_name;
193:   return(0);
194: }

196: /*@C
197:    NEPRegister - Adds a method to the nonlinear 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:    NEPRegister() may be called multiple times to add several user-defined solvers.

208:    Sample usage:
209: .vb
210:     NEPRegister("my_solver",MySolverCreate);
211: .ve

213:    Then, your solver can be chosen with the procedural interface via
214: $     NEPSetType(nep,"my_solver")
215:    or at runtime via the option
216: $     -nep_type my_solver

218:    Level: advanced

220: .seealso: NEPRegisterAll()
221: @*/
222: PetscErrorCode NEPRegister(const char *name,PetscErrorCode (*function)(NEP))
223: {

227:   NEPInitializePackage();
228:   PetscFunctionListAdd(&NEPList,name,function);
229:   return(0);
230: }

232: /*
233:    NEPReset_Problem - Destroys the problem matrices.
234: @*/
235: PetscErrorCode NEPReset_Problem(NEP nep)
236: {
238:   PetscInt       i;

242:   MatDestroy(&nep->function);
243:   MatDestroy(&nep->function_pre);
244:   MatDestroy(&nep->jacobian);
245:   MatDestroy(&nep->derivatives);
246:   if (nep->fui==NEP_USER_INTERFACE_SPLIT) {
247:     MatDestroyMatrices(nep->nt,&nep->A);
248:     for (i=0;i<nep->nt;i++) {
249:       FNDestroy(&nep->f[i]);
250:     }
251:     PetscFree(nep->f);
252:     PetscFree(nep->nrma);
253:     nep->nt = 0;
254:   }
255:   return(0);
256: }
257: /*@
258:    NEPReset - Resets the NEP context to the initial state (prior to setup)
259:    and destroys any allocated Vecs and Mats.

261:    Collective on NEP

263:    Input Parameter:
264: .  nep - eigensolver context obtained from NEPCreate()

266:    Level: advanced

268: .seealso: NEPDestroy()
269: @*/
270: PetscErrorCode NEPReset(NEP nep)
271: {

276:   if (!nep) return(0);
277:   if (nep->ops->reset) { (nep->ops->reset)(nep); }
278:   if (nep->refineksp) { KSPReset(nep->refineksp); }
279:   NEPReset_Problem(nep);
280:   BVDestroy(&nep->V);
281:   BVDestroy(&nep->W);
282:   VecDestroyVecs(nep->nwork,&nep->work);
283:   MatDestroy(&nep->resolvent);
284:   nep->nwork = 0;
285:   nep->state = NEP_STATE_INITIAL;
286:   return(0);
287: }

289: /*@
290:    NEPDestroy - Destroys the NEP context.

292:    Collective on NEP

294:    Input Parameter:
295: .  nep - eigensolver context obtained from NEPCreate()

297:    Level: beginner

299: .seealso: NEPCreate(), NEPSetUp(), NEPSolve()
300: @*/
301: PetscErrorCode NEPDestroy(NEP *nep)
302: {

306:   if (!*nep) return(0);
308:   if (--((PetscObject)(*nep))->refct > 0) { *nep = 0; return(0); }
309:   NEPReset(*nep);
310:   if ((*nep)->ops->destroy) { (*(*nep)->ops->destroy)(*nep); }
311:   if ((*nep)->eigr) {
312:     PetscFree4((*nep)->eigr,(*nep)->eigi,(*nep)->errest,(*nep)->perm);
313:   }
314:   RGDestroy(&(*nep)->rg);
315:   DSDestroy(&(*nep)->ds);
316:   KSPDestroy(&(*nep)->refineksp);
317:   PetscSubcommDestroy(&(*nep)->refinesubc);
318:   PetscFree((*nep)->sc);
319:   /* just in case the initial vectors have not been used */
320:   SlepcBasisDestroy_Private(&(*nep)->nini,&(*nep)->IS);
321:   if ((*nep)->convergeddestroy) {
322:     (*(*nep)->convergeddestroy)((*nep)->convergedctx);
323:   }
324:   NEPMonitorCancel(*nep);
325:   PetscHeaderDestroy(nep);
326:   return(0);
327: }

329: /*@
330:    NEPSetBV - Associates a basis vectors object to the nonlinear eigensolver.

332:    Collective on NEP

334:    Input Parameters:
335: +  nep - eigensolver context obtained from NEPCreate()
336: -  bv  - the basis vectors object

338:    Note:
339:    Use NEPGetBV() to retrieve the basis vectors context (for example,
340:    to free it at the end of the computations).

342:    Level: advanced

344: .seealso: NEPGetBV()
345: @*/
346: PetscErrorCode NEPSetBV(NEP nep,BV bv)
347: {

354:   PetscObjectReference((PetscObject)bv);
355:   BVDestroy(&nep->V);
356:   nep->V = bv;
357:   PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->V);
358:   return(0);
359: }

361: /*@
362:    NEPGetBV - Obtain the basis vectors object associated to the nonlinear
363:    eigensolver object.

365:    Not Collective

367:    Input Parameters:
368: .  nep - eigensolver context obtained from NEPCreate()

370:    Output Parameter:
371: .  bv - basis vectors context

373:    Level: advanced

375: .seealso: NEPSetBV()
376: @*/
377: PetscErrorCode NEPGetBV(NEP nep,BV *bv)
378: {

384:   if (!nep->V) {
385:     BVCreate(PetscObjectComm((PetscObject)nep),&nep->V);
386:     PetscObjectIncrementTabLevel((PetscObject)nep->V,(PetscObject)nep,0);
387:     PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->V);
388:     PetscObjectSetOptions((PetscObject)nep->V,((PetscObject)nep)->options);
389:   }
390:   *bv = nep->V;
391:   return(0);
392: }

394: /*@
395:    NEPSetRG - Associates a region object to the nonlinear eigensolver.

397:    Collective on NEP

399:    Input Parameters:
400: +  nep - eigensolver context obtained from NEPCreate()
401: -  rg  - the region object

403:    Note:
404:    Use NEPGetRG() to retrieve the region context (for example,
405:    to free it at the end of the computations).

407:    Level: advanced

409: .seealso: NEPGetRG()
410: @*/
411: PetscErrorCode NEPSetRG(NEP nep,RG rg)
412: {

419:   PetscObjectReference((PetscObject)rg);
420:   RGDestroy(&nep->rg);
421:   nep->rg = rg;
422:   PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->rg);
423:   return(0);
424: }

426: /*@
427:    NEPGetRG - Obtain the region object associated to the
428:    nonlinear eigensolver object.

430:    Not Collective

432:    Input Parameters:
433: .  nep - eigensolver context obtained from NEPCreate()

435:    Output Parameter:
436: .  rg - region context

438:    Level: advanced

440: .seealso: NEPSetRG()
441: @*/
442: PetscErrorCode NEPGetRG(NEP nep,RG *rg)
443: {

449:   if (!nep->rg) {
450:     RGCreate(PetscObjectComm((PetscObject)nep),&nep->rg);
451:     PetscObjectIncrementTabLevel((PetscObject)nep->rg,(PetscObject)nep,0);
452:     PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->rg);
453:     PetscObjectSetOptions((PetscObject)nep->rg,((PetscObject)nep)->options);
454:   }
455:   *rg = nep->rg;
456:   return(0);
457: }

459: /*@
460:    NEPSetDS - Associates a direct solver object to the nonlinear eigensolver.

462:    Collective on NEP

464:    Input Parameters:
465: +  nep - eigensolver context obtained from NEPCreate()
466: -  ds  - the direct solver object

468:    Note:
469:    Use NEPGetDS() to retrieve the direct solver context (for example,
470:    to free it at the end of the computations).

472:    Level: advanced

474: .seealso: NEPGetDS()
475: @*/
476: PetscErrorCode NEPSetDS(NEP nep,DS ds)
477: {

484:   PetscObjectReference((PetscObject)ds);
485:   DSDestroy(&nep->ds);
486:   nep->ds = ds;
487:   PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->ds);
488:   return(0);
489: }

491: /*@
492:    NEPGetDS - Obtain the direct solver object associated to the
493:    nonlinear eigensolver object.

495:    Not Collective

497:    Input Parameters:
498: .  nep - eigensolver context obtained from NEPCreate()

500:    Output Parameter:
501: .  ds - direct solver context

503:    Level: advanced

505: .seealso: NEPSetDS()
506: @*/
507: PetscErrorCode NEPGetDS(NEP nep,DS *ds)
508: {

514:   if (!nep->ds) {
515:     DSCreate(PetscObjectComm((PetscObject)nep),&nep->ds);
516:     PetscObjectIncrementTabLevel((PetscObject)nep->ds,(PetscObject)nep,0);
517:     PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->ds);
518:     PetscObjectSetOptions((PetscObject)nep->ds,((PetscObject)nep)->options);
519:   }
520:   *ds = nep->ds;
521:   return(0);
522: }

524: /*@
525:    NEPRefineGetKSP - Obtain the ksp object used by the eigensolver
526:    object in the refinement phase.

528:    Not Collective

530:    Input Parameters:
531: .  nep - eigensolver context obtained from NEPCreate()

533:    Output Parameter:
534: .  ksp - ksp context

536:    Level: advanced

538: .seealso: NEPSetRefine()
539: @*/
540: PetscErrorCode NEPRefineGetKSP(NEP nep,KSP *ksp)
541: {

547:   if (!nep->refineksp) {
548:     if (nep->npart>1) {
549:       /* Split in subcomunicators */
550:       PetscSubcommCreate(PetscObjectComm((PetscObject)nep),&nep->refinesubc);
551:       PetscSubcommSetNumber(nep->refinesubc,nep->npart);
552:       PetscSubcommSetType(nep->refinesubc,PETSC_SUBCOMM_CONTIGUOUS);
553:       PetscLogObjectMemory((PetscObject)nep,sizeof(PetscSubcomm));
554:     }
555:     KSPCreate((nep->npart==1)?PetscObjectComm((PetscObject)nep):PetscSubcommChild(nep->refinesubc),&nep->refineksp);
556:     PetscObjectIncrementTabLevel((PetscObject)nep->refineksp,(PetscObject)nep,0);
557:     PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->refineksp);
558:     PetscObjectSetOptions((PetscObject)nep->refineksp,((PetscObject)nep)->options);
559:     KSPSetOptionsPrefix(*ksp,((PetscObject)nep)->prefix);
560:     KSPAppendOptionsPrefix(*ksp,"nep_refine_");
561:     KSPSetTolerances(nep->refineksp,SLEPC_DEFAULT_TOL,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);
562:   }
563:   *ksp = nep->refineksp;
564:   return(0);
565: }

567: /*@
568:    NEPSetTarget - Sets the value of the target.

570:    Logically Collective on NEP

572:    Input Parameters:
573: +  nep    - eigensolver context
574: -  target - the value of the target

576:    Options Database Key:
577: .  -nep_target <scalar> - the value of the target

579:    Notes:
580:    The target is a scalar value used to determine the portion of the spectrum
581:    of interest. It is used in combination with NEPSetWhichEigenpairs().

583:    In the case of complex scalars, a complex value can be provided in the
584:    command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
585:    -nep_target 1.0+2.0i

587:    Level: intermediate

589: .seealso: NEPGetTarget(), NEPSetWhichEigenpairs()
590: @*/
591: PetscErrorCode NEPSetTarget(NEP nep,PetscScalar target)
592: {
596:   nep->target = target;
597:   return(0);
598: }

600: /*@
601:    NEPGetTarget - Gets the value of the target.

603:    Not Collective

605:    Input Parameter:
606: .  nep - eigensolver context

608:    Output Parameter:
609: .  target - the value of the target

611:    Note:
612:    If the target was not set by the user, then zero is returned.

614:    Level: intermediate

616: .seealso: NEPSetTarget()
617: @*/
618: PetscErrorCode NEPGetTarget(NEP nep,PetscScalar* target)
619: {
623:   *target = nep->target;
624:   return(0);
625: }

627: /*@C
628:    NEPSetFunction - Sets the function to compute the nonlinear Function T(lambda)
629:    as well as the location to store the matrix.

631:    Logically Collective on NEP and Mat

633:    Input Parameters:
634: +  nep - the NEP context
635: .  A   - Function matrix
636: .  B   - preconditioner matrix (usually same as the Function)
637: .  fun - Function evaluation routine (if NULL then NEP retains any
638:          previously set value)
639: -  ctx - [optional] user-defined context for private data for the Function
640:          evaluation routine (may be NULL) (if NULL then NEP retains any
641:          previously set value)

643:    Calling Sequence of fun:
644: $   fun(NEP nep,PetscScalar lambda,Mat F,Mat P,void *ctx)

646: +  nep    - the NEP context
647: .  lambda - the scalar argument where T(.) must be evaluated
648: .  T      - matrix that will contain T(lambda)
649: .  P      - (optional) different matrix to build the preconditioner
650: -  ctx    - (optional) user-defined context, as set by NEPSetFunction()

652:    Level: beginner

654: .seealso: NEPGetFunction(), NEPSetJacobian()
655: @*/
656: PetscErrorCode NEPSetFunction(NEP nep,Mat A,Mat B,PetscErrorCode (*fun)(NEP,PetscScalar,Mat,Mat,void*),void *ctx)
657: {


667:   if (nep->state) { NEPReset(nep); }
668:   else if (nep->fui && nep->fui!=NEP_USER_INTERFACE_CALLBACK) { NEPReset_Problem(nep); }

670:   if (fun) nep->computefunction = fun;
671:   if (ctx) nep->functionctx     = ctx;
672:   if (A) {
673:     PetscObjectReference((PetscObject)A);
674:     MatDestroy(&nep->function);
675:     nep->function = A;
676:   }
677:   if (B) {
678:     PetscObjectReference((PetscObject)B);
679:     MatDestroy(&nep->function_pre);
680:     nep->function_pre = B;
681:   }
682:   nep->fui   = NEP_USER_INTERFACE_CALLBACK;
683:   nep->state = NEP_STATE_INITIAL;
684:   return(0);
685: }

687: /*@C
688:    NEPGetFunction - Returns the Function matrix and optionally the user
689:    provided context for evaluating the Function.

691:    Not Collective, but Mat object will be parallel if NEP object is

693:    Input Parameter:
694: .  nep - the nonlinear eigensolver context

696:    Output Parameters:
697: +  A   - location to stash Function matrix (or NULL)
698: .  B   - location to stash preconditioner matrix (or NULL)
699: .  fun - location to put Function function (or NULL)
700: -  ctx - location to stash Function context (or NULL)

702:    Level: advanced

704: .seealso: NEPSetFunction()
705: @*/
706: PetscErrorCode NEPGetFunction(NEP nep,Mat *A,Mat *B,PetscErrorCode (**fun)(NEP,PetscScalar,Mat,Mat,void*),void **ctx)
707: {
710:   NEPCheckCallback(nep,1);
711:   if (A)   *A   = nep->function;
712:   if (B)   *B   = nep->function_pre;
713:   if (fun) *fun = nep->computefunction;
714:   if (ctx) *ctx = nep->functionctx;
715:   return(0);
716: }

718: /*@C
719:    NEPSetJacobian - Sets the function to compute Jacobian T'(lambda) as well
720:    as the location to store the matrix.

722:    Logically Collective on NEP and Mat

724:    Input Parameters:
725: +  nep - the NEP context
726: .  A   - Jacobian matrix
727: .  jac - Jacobian evaluation routine (if NULL then NEP retains any
728:          previously set value)
729: -  ctx - [optional] user-defined context for private data for the Jacobian
730:          evaluation routine (may be NULL) (if NULL then NEP retains any
731:          previously set value)

733:    Calling Sequence of jac:
734: $   jac(NEP nep,PetscScalar lambda,Mat J,void *ctx)

736: +  nep    - the NEP context
737: .  lambda - the scalar argument where T'(.) must be evaluated
738: .  J      - matrix that will contain T'(lambda)
739: -  ctx    - (optional) user-defined context, as set by NEPSetJacobian()

741:    Level: beginner

743: .seealso: NEPSetFunction(), NEPGetJacobian()
744: @*/
745: PetscErrorCode NEPSetJacobian(NEP nep,Mat A,PetscErrorCode (*jac)(NEP,PetscScalar,Mat,void*),void *ctx)
746: {


754:   if (nep->state) { NEPReset(nep); }
755:   else if (nep->fui && nep->fui!=NEP_USER_INTERFACE_CALLBACK) { NEPReset_Problem(nep); }

757:   if (jac) nep->computejacobian = jac;
758:   if (ctx) nep->jacobianctx     = ctx;
759:   if (A) {
760:     PetscObjectReference((PetscObject)A);
761:     MatDestroy(&nep->jacobian);
762:     nep->jacobian = A;
763:   }
764:   nep->fui   = NEP_USER_INTERFACE_CALLBACK;
765:   nep->state = NEP_STATE_INITIAL;
766:   return(0);
767: }

769: /*@C
770:    NEPGetJacobian - Returns the Jacobian matrix and optionally the user
771:    provided routine and context for evaluating the Jacobian.

773:    Not Collective, but Mat object will be parallel if NEP object is

775:    Input Parameter:
776: .  nep - the nonlinear eigensolver context

778:    Output Parameters:
779: +  A   - location to stash Jacobian matrix (or NULL)
780: .  jac - location to put Jacobian function (or NULL)
781: -  ctx - location to stash Jacobian context (or NULL)

783:    Level: advanced

785: .seealso: NEPSetJacobian()
786: @*/
787: PetscErrorCode NEPGetJacobian(NEP nep,Mat *A,PetscErrorCode (**jac)(NEP,PetscScalar,Mat,void*),void **ctx)
788: {
791:   NEPCheckCallback(nep,1);
792:   if (A)   *A   = nep->jacobian;
793:   if (jac) *jac = nep->computejacobian;
794:   if (ctx) *ctx = nep->jacobianctx;
795:   return(0);
796: }

798: /*@
799:    NEPSetSplitOperator - Sets the operator of the nonlinear eigenvalue problem
800:    in split form.

802:    Collective on NEP, Mat and FN

804:    Input Parameters:
805: +  nep - the nonlinear eigensolver context
806: .  n   - number of terms in the split form
807: .  A   - array of matrices
808: .  f   - array of functions
809: -  str - structure flag for matrices

811:    Notes:
812:    The nonlinear operator is written as T(lambda) = sum_i A_i*f_i(lambda),
813:    for i=1,...,n. The derivative T'(lambda) can be obtained using the
814:    derivatives of f_i.

816:    The structure flag provides information about A_i's nonzero pattern
817:    (see MatStructure enum). If all matrices have the same pattern, then
818:    use SAME_NONZERO_PATTERN. If the patterns are different but contained
819:    in the pattern of the first one, then use SUBSET_NONZERO_PATTERN.
820:    Otherwise use DIFFERENT_NONZERO_PATTERN.

822:    This function must be called before NEPSetUp(). If it is called again
823:    after NEPSetUp() then the NEP object is reset.

825:    Level: beginner

827: .seealso: NEPGetSplitOperatorTerm(), NEPGetSplitOperatorInfo()
828:  @*/
829: PetscErrorCode NEPSetSplitOperator(NEP nep,PetscInt n,Mat A[],FN f[],MatStructure str)
830: {
831:   PetscInt       i;

837:   if (n <= 0) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Must have one or more terms, you have %D",n);

843:   for (i=0;i<n;i++) {
846:     PetscObjectReference((PetscObject)A[i]);
847:     PetscObjectReference((PetscObject)f[i]);
848:   }

850:   if (nep->state) { NEPReset(nep); }
851:   else { NEPReset_Problem(nep); }

853:   /* allocate space and copy matrices and functions */
854:   PetscMalloc1(n,&nep->A);
855:   PetscLogObjectMemory((PetscObject)nep,n*sizeof(Mat));
856:   for (i=0;i<n;i++) nep->A[i] = A[i];
857:   PetscMalloc1(n,&nep->f);
858:   PetscLogObjectMemory((PetscObject)nep,n*sizeof(FN));
859:   for (i=0;i<n;i++) nep->f[i] = f[i];
860:   PetscCalloc1(n,&nep->nrma);
861:   PetscLogObjectMemory((PetscObject)nep,n*sizeof(PetscReal));
862:   nep->nt    = n;
863:   nep->mstr  = str;
864:   nep->fui   = NEP_USER_INTERFACE_SPLIT;
865:   nep->state = NEP_STATE_INITIAL;
866:   return(0);
867: }

869: /*@
870:    NEPGetSplitOperatorTerm - Gets the matrices and functions associated with
871:    the nonlinear operator in split form.

873:    Not collective, though parallel Mats and FNs are returned if the NEP is parallel

875:    Input Parameter:
876: +  nep - the nonlinear eigensolver context
877: -  k   - the index of the requested term (starting in 0)

879:    Output Parameters:
880: +  A - the matrix of the requested term
881: -  f - the function of the requested term

883:    Level: intermediate

885: .seealso: NEPSetSplitOperator(), NEPGetSplitOperatorInfo()
886: @*/
887: PetscErrorCode NEPGetSplitOperatorTerm(NEP nep,PetscInt k,Mat *A,FN *f)
888: {
891:   NEPCheckSplit(nep,1);
892:   if (k<0 || k>=nep->nt) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %D",nep->nt-1);
893:   if (A) *A = nep->A[k];
894:   if (f) *f = nep->f[k];
895:   return(0);
896: }

898: /*@
899:    NEPGetSplitOperatorInfo - Returns the number of terms of the split form of
900:    the nonlinear operator, as well as the structure flag for matrices.

902:    Not collective

904:    Input Parameter:
905: .  nep - the nonlinear eigensolver context

907:    Output Parameters:
908: +  n   - the number of terms passed in NEPSetSplitOperator()
909: -  str - the matrix structure flag passed in NEPSetSplitOperator()

911:    Level: intermediate

913: .seealso: NEPSetSplitOperator(), NEPGetSplitOperatorTerm()
914: @*/
915: PetscErrorCode NEPGetSplitOperatorInfo(NEP nep,PetscInt *n,MatStructure *str)
916: {
919:   NEPCheckSplit(nep,1);
920:   if (n)   *n = nep->nt;
921:   if (str) *str = nep->mstr;
922:   return(0);
923: }

925: /*@C
926:    NEPSetDerivatives - Sets the function to compute the k-th derivative T^(k)(lambda)
927:    for any value of k (including 0), as well as the location to store the matrix.

929:    Logically Collective on NEP and Mat

931:    Input Parameters:
932: +  nep - the NEP context
933: .  A   - the matrix to store the computed derivative
934: .  der - routing to evaluate the k-th derivative (if NULL then NEP retains any
935:          previously set value)
936: -  ctx - [optional] user-defined context for private data for the derivatives
937:          evaluation routine (may be NULL) (if NULL then NEP retains any
938:          previously set value)

940:    Level: beginner

942: .seealso: NEPSetFunction(), NEPGetDerivatives()
943: @*/
944: PetscErrorCode NEPSetDerivatives(NEP nep,Mat A,PetscErrorCode (*der)(NEP,PetscScalar,PetscInt,Mat,void*),void *ctx)
945: {


953:   if (nep->state) { NEPReset(nep); }
954:   else { NEPReset_Problem(nep); }

956:   if (der) nep->computederivatives = der;
957:   if (ctx) nep->derivativesctx     = ctx;
958:   if (A) {
959:     PetscObjectReference((PetscObject)A);
960:     MatDestroy(&nep->derivatives);
961:     nep->derivatives = A;
962:   }
963:   nep->fui   = NEP_USER_INTERFACE_DERIVATIVES;
964:   nep->state = NEP_STATE_INITIAL;
965:   return(0);
966: }

968: /*@C
969:    NEPGetDerivatives - Returns the derivatives matrix and optionally the user
970:    provided routine and context for evaluating the derivatives.

972:    Not Collective, but Mat object will be parallel if NEP object is

974:    Input Parameter:
975: .  nep - the nonlinear eigensolver context

977:    Output Parameters:
978: +  A   - location to stash the derivatives matrix (or NULL)
979: .  der - location to put derivatives function (or NULL)
980: -  ctx - location to stash derivatives context (or NULL)

982:    Level: advanced

984: .seealso: NEPSetDerivatives()
985: @*/
986: PetscErrorCode NEPGetDerivatives(NEP nep,Mat *A,PetscErrorCode (**der)(NEP,PetscScalar,PetscInt,Mat,void*),void **ctx)
987: {
990:   NEPCheckDerivatives(nep,1);
991:   if (A)   *A   = nep->derivatives;
992:   if (der) *der = nep->computederivatives;
993:   if (ctx) *ctx = nep->derivativesctx;
994:   return(0);
995: }