1: /*
2: Basic NEP routines, Create, View, etc.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2016, 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/nepimpl.h> /*I "slepcnep.h" I*/
26: PetscFunctionList NEPList = 0;
27: PetscBool NEPRegisterAllCalled = PETSC_FALSE;
28: PetscClassId NEP_CLASSID = 0;
29: PetscLogEvent NEP_SetUp = 0,NEP_Solve = 0,NEP_Refine = 0,NEP_FunctionEval = 0,NEP_JacobianEval = 0,NEP_DerivativesEval = 0;
33: /*@
34: NEPCreate - Creates the default NEP context.
36: Collective on MPI_Comm
38: Input Parameter:
39: . comm - MPI communicator
41: Output Parameter:
42: . nep - location to put the NEP context
44: Level: beginner
46: .seealso: NEPSetUp(), NEPSolve(), NEPDestroy(), NEP 47: @*/
48: PetscErrorCode NEPCreate(MPI_Comm comm,NEP *outnep) 49: {
51: NEP nep;
55: *outnep = 0;
56: NEPInitializePackage();
57: SlepcHeaderCreate(nep,NEP_CLASSID,"NEP","Nonlinear Eigenvalue Problem","NEP",comm,NEPDestroy,NEPView);
59: nep->max_it = 0;
60: nep->nev = 1;
61: nep->ncv = 0;
62: nep->mpd = 0;
63: nep->nini = 0;
64: nep->target = 0.0;
65: nep->tol = PETSC_DEFAULT;
66: nep->conv = NEP_CONV_REL;
67: nep->stop = NEP_STOP_BASIC;
68: nep->which = (NEPWhich)0;
69: nep->refine = NEP_REFINE_NONE;
70: nep->npart = 1;
71: nep->rtol = PETSC_DEFAULT;
72: nep->rits = PETSC_DEFAULT;
73: nep->scheme = (NEPRefineScheme)0;
74: nep->trackall = PETSC_FALSE;
76: nep->computefunction = NULL;
77: nep->computejacobian = NULL;
78: nep->functionctx = NULL;
79: nep->jacobianctx = NULL;
80: nep->computederivatives = NULL;
81: nep->derivativesctx = NULL;
82: nep->converged = NEPConvergedRelative;
83: nep->convergeddestroy= NULL;
84: nep->stopping = NEPStoppingBasic;
85: nep->stoppingdestroy = NULL;
86: nep->convergedctx = NULL;
87: nep->stoppingctx = NULL;
88: nep->numbermonitors = 0;
90: nep->ds = NULL;
91: nep->V = NULL;
92: nep->rg = NULL;
93: nep->function = NULL;
94: nep->function_pre = NULL;
95: nep->jacobian = NULL;
96: nep->derivatives = NULL;
97: nep->A = NULL;
98: nep->f = NULL;
99: nep->nt = 0;
100: nep->mstr = DIFFERENT_NONZERO_PATTERN;
101: nep->IS = NULL;
102: nep->eigr = NULL;
103: nep->eigi = NULL;
104: nep->errest = NULL;
105: nep->perm = NULL;
106: nep->nwork = 0;
107: nep->work = NULL;
108: nep->data = NULL;
110: nep->state = NEP_STATE_INITIAL;
111: nep->nconv = 0;
112: nep->its = 0;
113: nep->n = 0;
114: nep->nloc = 0;
115: nep->nrma = NULL;
116: nep->fui = (NEPUserInterface)0;
117: nep->reason = NEP_CONVERGED_ITERATING;
119: PetscNewLog(nep,&nep->sc);
120: *outnep = nep;
121: return(0);
122: }
126: /*@C
127: NEPSetType - Selects the particular solver to be used in the NEP object.
129: Logically Collective on NEP131: Input Parameters:
132: + nep - the nonlinear eigensolver context
133: - type - a known method
135: Options Database Key:
136: . -nep_type <method> - Sets the method; use -help for a list
137: of available methods
139: Notes:
140: See "slepc/include/slepcnep.h" for available methods.
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.
150: Level: intermediate
152: .seealso: NEPType153: @*/
154: PetscErrorCode NEPSetType(NEP nep,NEPType type)155: {
156: PetscErrorCode ierr,(*r)(NEP);
157: PetscBool match;
163: PetscObjectTypeCompare((PetscObject)nep,type,&match);
164: if (match) return(0);
166: PetscFunctionListFind(NEPList,type,&r);
167: if (!r) SETERRQ1(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_UNKNOWN_TYPE,"Unknown NEP type given: %s",type);
169: if (nep->ops->destroy) { (*nep->ops->destroy)(nep); }
170: PetscMemzero(nep->ops,sizeof(struct _NEPOps));
172: nep->state = NEP_STATE_INITIAL;
173: PetscObjectChangeTypeName((PetscObject)nep,type);
174: (*r)(nep);
175: return(0);
176: }
180: /*@C
181: NEPGetType - Gets the NEP type as a string from the NEP object.
183: Not Collective
185: Input Parameter:
186: . nep - the eigensolver context
188: Output Parameter:
189: . name - name of NEP method
191: Level: intermediate
193: .seealso: NEPSetType()
194: @*/
195: PetscErrorCode NEPGetType(NEP nep,NEPType *type)196: {
200: *type = ((PetscObject)nep)->type_name;
201: return(0);
202: }
206: /*@C
207: NEPRegister - Adds a method to the nonlinear eigenproblem solver package.
209: Not Collective
211: Input Parameters:
212: + name - name of a new user-defined solver
213: - function - routine to create the solver context
215: Notes:
216: NEPRegister() may be called multiple times to add several user-defined solvers.
218: Sample usage:
219: .vb
220: NEPRegister("my_solver",MySolverCreate);
221: .ve
223: Then, your solver can be chosen with the procedural interface via
224: $ NEPSetType(nep,"my_solver")
225: or at runtime via the option
226: $ -nep_type my_solver
228: Level: advanced
230: .seealso: NEPRegisterAll()
231: @*/
232: PetscErrorCode NEPRegister(const char *name,PetscErrorCode (*function)(NEP))233: {
237: PetscFunctionListAdd(&NEPList,name,function);
238: return(0);
239: }
243: /*
244: NEPReset_Problem - Destroys the problem matrices.
245: @*/
246: PetscErrorCode NEPReset_Problem(NEP nep)247: {
249: PetscInt i;
253: MatDestroy(&nep->function);
254: MatDestroy(&nep->function_pre);
255: MatDestroy(&nep->jacobian);
256: MatDestroy(&nep->derivatives);
257: if (nep->fui==NEP_USER_INTERFACE_SPLIT) {
258: MatDestroyMatrices(nep->nt,&nep->A);
259: for (i=0;i<nep->nt;i++) {
260: FNDestroy(&nep->f[i]);
261: }
262: PetscFree(nep->f);
263: PetscFree(nep->nrma);
264: }
265: return(0);
266: }
269: /*@
270: NEPReset - Resets the NEP context to the initial state and removes any
271: allocated objects.
273: Collective on NEP275: Input Parameter:
276: . nep - eigensolver context obtained from NEPCreate()
278: Level: advanced
280: .seealso: NEPDestroy()
281: @*/
282: PetscErrorCode NEPReset(NEP nep)283: {
285: PetscInt ncols;
289: if (nep->ops->reset) { (nep->ops->reset)(nep); }
290: if (nep->ds) { DSReset(nep->ds); }
291: NEPReset_Problem(nep);
292: BVGetSizes(nep->V,NULL,NULL,&ncols);
293: if (ncols) {
294: PetscFree4(nep->eigr,nep->eigi,nep->errest,nep->perm);
295: }
296: BVDestroy(&nep->V);
297: VecDestroyVecs(nep->nwork,&nep->work);
298: KSPDestroy(&nep->refineksp);
299: PetscSubcommDestroy(&nep->refinesubc);
300: nep->nwork = 0;
301: nep->state = NEP_STATE_INITIAL;
302: return(0);
303: }
307: /*@
308: NEPDestroy - Destroys the NEP context.
310: Collective on NEP312: Input Parameter:
313: . nep - eigensolver context obtained from NEPCreate()
315: Level: beginner
317: .seealso: NEPCreate(), NEPSetUp(), NEPSolve()
318: @*/
319: PetscErrorCode NEPDestroy(NEP *nep)320: {
324: if (!*nep) return(0);
326: if (--((PetscObject)(*nep))->refct > 0) { *nep = 0; return(0); }
327: NEPReset(*nep);
328: if ((*nep)->ops->destroy) { (*(*nep)->ops->destroy)(*nep); }
329: RGDestroy(&(*nep)->rg);
330: DSDestroy(&(*nep)->ds);
331: PetscFree((*nep)->sc);
332: /* just in case the initial vectors have not been used */
333: SlepcBasisDestroy_Private(&(*nep)->nini,&(*nep)->IS);
334: if ((*nep)->convergeddestroy) {
335: (*(*nep)->convergeddestroy)((*nep)->convergedctx);
336: }
337: NEPMonitorCancel(*nep);
338: PetscHeaderDestroy(nep);
339: return(0);
340: }
344: /*@
345: NEPSetBV - Associates a basis vectors object to the nonlinear eigensolver.
347: Collective on NEP349: Input Parameters:
350: + nep - eigensolver context obtained from NEPCreate()
351: - bv - the basis vectors object
353: Note:
354: Use NEPGetBV() to retrieve the basis vectors context (for example,
355: to free it at the end of the computations).
357: Level: advanced
359: .seealso: NEPGetBV()
360: @*/
361: PetscErrorCode NEPSetBV(NEP nep,BV bv)362: {
369: PetscObjectReference((PetscObject)bv);
370: BVDestroy(&nep->V);
371: nep->V = bv;
372: PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->V);
373: return(0);
374: }
378: /*@
379: NEPGetBV - Obtain the basis vectors object associated to the nonlinear
380: eigensolver object.
382: Not Collective
384: Input Parameters:
385: . nep - eigensolver context obtained from NEPCreate()
387: Output Parameter:
388: . bv - basis vectors context
390: Level: advanced
392: .seealso: NEPSetBV()
393: @*/
394: PetscErrorCode NEPGetBV(NEP nep,BV *bv)395: {
401: if (!nep->V) {
402: BVCreate(PetscObjectComm((PetscObject)nep),&nep->V);
403: PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->V);
404: }
405: *bv = nep->V;
406: return(0);
407: }
411: /*@
412: NEPSetRG - Associates a region object to the nonlinear eigensolver.
414: Collective on NEP416: Input Parameters:
417: + nep - eigensolver context obtained from NEPCreate()
418: - rg - the region object
420: Note:
421: Use NEPGetRG() to retrieve the region context (for example,
422: to free it at the end of the computations).
424: Level: advanced
426: .seealso: NEPGetRG()
427: @*/
428: PetscErrorCode NEPSetRG(NEP nep,RG rg)429: {
436: PetscObjectReference((PetscObject)rg);
437: RGDestroy(&nep->rg);
438: nep->rg = rg;
439: PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->rg);
440: return(0);
441: }
445: /*@
446: NEPGetRG - Obtain the region object associated to the
447: nonlinear eigensolver object.
449: Not Collective
451: Input Parameters:
452: . nep - eigensolver context obtained from NEPCreate()
454: Output Parameter:
455: . rg - region context
457: Level: advanced
459: .seealso: NEPSetRG()
460: @*/
461: PetscErrorCode NEPGetRG(NEP nep,RG *rg)462: {
468: if (!nep->rg) {
469: RGCreate(PetscObjectComm((PetscObject)nep),&nep->rg);
470: PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->rg);
471: }
472: *rg = nep->rg;
473: return(0);
474: }
478: /*@
479: NEPSetDS - Associates a direct solver object to the nonlinear eigensolver.
481: Collective on NEP483: Input Parameters:
484: + nep - eigensolver context obtained from NEPCreate()
485: - ds - the direct solver object
487: Note:
488: Use NEPGetDS() to retrieve the direct solver context (for example,
489: to free it at the end of the computations).
491: Level: advanced
493: .seealso: NEPGetDS()
494: @*/
495: PetscErrorCode NEPSetDS(NEP nep,DS ds)496: {
503: PetscObjectReference((PetscObject)ds);
504: DSDestroy(&nep->ds);
505: nep->ds = ds;
506: PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->ds);
507: return(0);
508: }
512: /*@
513: NEPGetDS - Obtain the direct solver object associated to the
514: nonlinear eigensolver object.
516: Not Collective
518: Input Parameters:
519: . nep - eigensolver context obtained from NEPCreate()
521: Output Parameter:
522: . ds - direct solver context
524: Level: advanced
526: .seealso: NEPSetDS()
527: @*/
528: PetscErrorCode NEPGetDS(NEP nep,DS *ds)529: {
535: if (!nep->ds) {
536: DSCreate(PetscObjectComm((PetscObject)nep),&nep->ds);
537: PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->ds);
538: }
539: *ds = nep->ds;
540: return(0);
541: }
545: /*@
546: NEPRefineGetKSP - Obtain the ksp object used by the eigensolver
547: object in the refinement phase.
549: Not Collective
551: Input Parameters:
552: . nep - eigensolver context obtained from NEPCreate()
554: Output Parameter:
555: . ksp - ksp context
557: Level: advanced
559: .seealso: NEPSetRefine()
560: @*/
561: PetscErrorCode NEPRefineGetKSP(NEP nep,KSP *ksp)562: {
568: if (!nep->refineksp) {
569: if (nep->npart>1) {
570: /* Split in subcomunicators */
571: PetscSubcommCreate(PetscObjectComm((PetscObject)nep),&nep->refinesubc);
572: PetscSubcommSetNumber(nep->refinesubc,nep->npart);
573: PetscSubcommSetType(nep->refinesubc,PETSC_SUBCOMM_CONTIGUOUS);
574: PetscLogObjectMemory((PetscObject)nep,sizeof(PetscSubcomm));
575: }
576: KSPCreate((nep->npart==1)?PetscObjectComm((PetscObject)nep):PetscSubcommChild(nep->refinesubc),&nep->refineksp);
577: PetscLogObjectParent((PetscObject)nep,(PetscObject)nep->refineksp);
578: KSPSetOptionsPrefix(*ksp,((PetscObject)nep)->prefix);
579: KSPAppendOptionsPrefix(*ksp,"nep_refine_");
580: KSPSetErrorIfNotConverged(*ksp,PETSC_TRUE);
581: }
582: *ksp = nep->refineksp;
583: return(0);
584: }
588: /*@
589: NEPSetTarget - Sets the value of the target.
591: Logically Collective on NEP593: Input Parameters:
594: + nep - eigensolver context
595: - target - the value of the target
597: Options Database Key:
598: . -nep_target <scalar> - the value of the target
600: Notes:
601: The target is a scalar value used to determine the portion of the spectrum
602: of interest. It is used in combination with NEPSetWhichEigenpairs().
604: In the case of complex scalars, a complex value can be provided in the
605: command line with [+/-][realnumber][+/-]realnumberi with no spaces, e.g.
606: -nep_target 1.0+2.0i
608: Level: intermediate
610: .seealso: NEPGetTarget(), NEPSetWhichEigenpairs()
611: @*/
612: PetscErrorCode NEPSetTarget(NEP nep,PetscScalar target)613: {
617: nep->target = target;
618: return(0);
619: }
623: /*@
624: NEPGetTarget - Gets the value of the target.
626: Not Collective
628: Input Parameter:
629: . nep - eigensolver context
631: Output Parameter:
632: . target - the value of the target
634: Note:
635: If the target was not set by the user, then zero is returned.
637: Level: intermediate
639: .seealso: NEPSetTarget()
640: @*/
641: PetscErrorCode NEPGetTarget(NEP nep,PetscScalar* target)642: {
646: *target = nep->target;
647: return(0);
648: }
652: /*@C
653: NEPSetFunction - Sets the function to compute the nonlinear Function T(lambda)
654: as well as the location to store the matrix.
656: Logically Collective on NEP and Mat
658: Input Parameters:
659: + nep - the NEP context
660: . A - Function matrix
661: . B - preconditioner matrix (usually same as the Function)
662: . fun - Function evaluation routine (if NULL then NEP retains any
663: previously set value)
664: - ctx - [optional] user-defined context for private data for the Function
665: evaluation routine (may be NULL) (if NULL then NEP retains any
666: previously set value)
668: Calling Sequence of fun:
669: $ fun(NEP nep,PetscScalar lambda,Mat F,Mat P,void *ctx)
671: + nep - the NEP context
672: . lambda - the scalar argument where T(.) must be evaluated
673: . T - matrix that will contain T(lambda)
674: . P - (optional) different matrix to build the preconditioner
675: - ctx - (optional) user-defined context, as set by NEPSetFunction()
677: Level: beginner
679: .seealso: NEPGetFunction(), NEPSetJacobian()
680: @*/
681: PetscErrorCode NEPSetFunction(NEP nep,Mat A,Mat B,PetscErrorCode (*fun)(NEP,PetscScalar,Mat,Mat,void*),void *ctx)682: {
692: if (nep->fui && nep->fui!=NEP_USER_INTERFACE_CALLBACK) { /* clean previous user info */
693: NEPReset_Problem(nep);
694: }
696: if (fun) nep->computefunction = fun;
697: if (ctx) nep->functionctx = ctx;
698: if (A) {
699: PetscObjectReference((PetscObject)A);
700: MatDestroy(&nep->function);
701: nep->function = A;
702: }
703: if (B) {
704: PetscObjectReference((PetscObject)B);
705: MatDestroy(&nep->function_pre);
706: nep->function_pre = B;
707: }
708: nep->fui = NEP_USER_INTERFACE_CALLBACK;
709: return(0);
710: }
714: /*@C
715: NEPGetFunction - Returns the Function matrix and optionally the user
716: provided context for evaluating the Function.
718: Not Collective, but Mat object will be parallel if NEP object is
720: Input Parameter:
721: . nep - the nonlinear eigensolver context
723: Output Parameters:
724: + A - location to stash Function matrix (or NULL)
725: . B - location to stash preconditioner matrix (or NULL)
726: . fun - location to put Function function (or NULL)
727: - ctx - location to stash Function context (or NULL)
729: Level: advanced
731: .seealso: NEPSetFunction()
732: @*/
733: PetscErrorCode NEPGetFunction(NEP nep,Mat *A,Mat *B,PetscErrorCode (**fun)(NEP,PetscScalar,Mat,Mat,void*),void **ctx)734: {
737: NEPCheckCallback(nep,1);
738: if (A) *A = nep->function;
739: if (B) *B = nep->function_pre;
740: if (fun) *fun = nep->computefunction;
741: if (ctx) *ctx = nep->functionctx;
742: return(0);
743: }
747: /*@C
748: NEPSetJacobian - Sets the function to compute Jacobian T'(lambda) as well
749: as the location to store the matrix.
751: Logically Collective on NEP and Mat
753: Input Parameters:
754: + nep - the NEP context
755: . A - Jacobian matrix
756: . jac - Jacobian evaluation routine (if NULL then NEP retains any
757: previously set value)
758: - ctx - [optional] user-defined context for private data for the Jacobian
759: evaluation routine (may be NULL) (if NULL then NEP retains any
760: previously set value)
762: Calling Sequence of jac:
763: $ jac(NEP nep,PetscScalar lambda,Mat J,void *ctx)
765: + nep - the NEP context
766: . lambda - the scalar argument where T'(.) must be evaluated
767: . J - matrix that will contain T'(lambda)
768: - ctx - (optional) user-defined context, as set by NEPSetJacobian()
770: Level: beginner
772: .seealso: NEPSetFunction(), NEPGetJacobian()
773: @*/
774: PetscErrorCode NEPSetJacobian(NEP nep,Mat A,PetscErrorCode (*jac)(NEP,PetscScalar,Mat,void*),void *ctx)775: {
783: if (nep->fui && nep->fui!=NEP_USER_INTERFACE_CALLBACK) { /* clean previous user info */
784: NEPReset_Problem(nep);
785: }
787: if (jac) nep->computejacobian = jac;
788: if (ctx) nep->jacobianctx = ctx;
789: if (A) {
790: PetscObjectReference((PetscObject)A);
791: MatDestroy(&nep->jacobian);
792: nep->jacobian = A;
793: }
794: nep->fui = NEP_USER_INTERFACE_CALLBACK;
795: return(0);
796: }
800: /*@C
801: NEPGetJacobian - Returns the Jacobian matrix and optionally the user
802: provided routine and context for evaluating the Jacobian.
804: Not Collective, but Mat object will be parallel if NEP object is
806: Input Parameter:
807: . nep - the nonlinear eigensolver context
809: Output Parameters:
810: + A - location to stash Jacobian matrix (or NULL)
811: . jac - location to put Jacobian function (or NULL)
812: - ctx - location to stash Jacobian context (or NULL)
814: Level: advanced
816: .seealso: NEPSetJacobian()
817: @*/
818: PetscErrorCode NEPGetJacobian(NEP nep,Mat *A,PetscErrorCode (**jac)(NEP,PetscScalar,Mat,void*),void **ctx)819: {
822: NEPCheckCallback(nep,1);
823: if (A) *A = nep->jacobian;
824: if (jac) *jac = nep->computejacobian;
825: if (ctx) *ctx = nep->jacobianctx;
826: return(0);
827: }
831: /*@
832: NEPSetSplitOperator - Sets the operator of the nonlinear eigenvalue problem
833: in split form.
835: Collective on NEP, Mat and FN837: Input Parameters:
838: + nep - the nonlinear eigensolver context
839: . n - number of terms in the split form
840: . A - array of matrices
841: . f - array of functions
842: - str - structure flag for matrices
844: Notes:
845: The nonlinear operator is written as T(lambda) = sum_i A_i*f_i(lambda),
846: for i=1,...,n. The derivative T'(lambda) can be obtained using the
847: derivatives of f_i.
849: The structure flag provides information about A_i's nonzero pattern
850: (see MatStructure enum). If all matrices have the same pattern, then
851: use SAME_NONZERO_PATTERN. If the patterns are different but contained
852: in the pattern of the first one, then use SUBSET_NONZERO_PATTERN.
853: Otherwise use DIFFERENT_NONZERO_PATTERN.
855: This function must be called before NEPSetUp(). If it is called again
856: after NEPSetUp() then the NEP object is reset.
858: Level: beginner
860: .seealso: NEPGetSplitOperatorTerm(), NEPGetSplitOperatorInfo()
861: @*/
862: PetscErrorCode NEPSetSplitOperator(NEP nep,PetscInt n,Mat A[],FN f[],MatStructure str)863: {
864: PetscInt i;
870: if (n <= 0) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must have one or more terms, you have %D",n);
875: if (nep->state) { NEPReset(nep); }
876: /* clean previously stored information */
877: NEPReset_Problem(nep);
878: /* allocate space and copy matrices and functions */
879: PetscMalloc1(n,&nep->A);
880: PetscLogObjectMemory((PetscObject)nep,n*sizeof(Mat));
881: for (i=0;i<n;i++) {
883: PetscObjectReference((PetscObject)A[i]);
884: nep->A[i] = A[i];
885: }
886: PetscMalloc1(n,&nep->f);
887: PetscLogObjectMemory((PetscObject)nep,n*sizeof(FN));
888: for (i=0;i<n;i++) {
890: PetscObjectReference((PetscObject)f[i]);
891: nep->f[i] = f[i];
892: }
893: PetscCalloc1(n,&nep->nrma);
894: PetscLogObjectMemory((PetscObject)nep,n*sizeof(PetscReal));
895: nep->nt = n;
896: nep->mstr = str;
897: nep->fui = NEP_USER_INTERFACE_SPLIT;
898: return(0);
899: }
903: /*@
904: NEPGetSplitOperatorTerm - Gets the matrices and functions associated with
905: the nonlinear operator in split form.
907: Not collective, though parallel Mats and FNs are returned if the NEP is parallel
909: Input Parameter:
910: + nep - the nonlinear eigensolver context
911: - k - the index of the requested term (starting in 0)
913: Output Parameters:
914: + A - the matrix of the requested term
915: - f - the function of the requested term
917: Level: intermediate
919: .seealso: NEPSetSplitOperator(), NEPGetSplitOperatorInfo()
920: @*/
921: PetscErrorCode NEPGetSplitOperatorTerm(NEP nep,PetscInt k,Mat *A,FN *f)922: {
925: NEPCheckSplit(nep,1);
926: if (k<0 || k>=nep->nt) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %D",nep->nt-1);
927: if (A) *A = nep->A[k];
928: if (f) *f = nep->f[k];
929: return(0);
930: }
934: /*@
935: NEPGetSplitOperatorInfo - Returns the number of terms of the split form of
936: the nonlinear operator, as well as the structure flag for matrices.
938: Not collective
940: Input Parameter:
941: . nep - the nonlinear eigensolver context
943: Output Parameters:
944: + n - the number of terms passed in NEPSetSplitOperator()
945: - str - the matrix structure flag passed in NEPSetSplitOperator()
947: Level: intermediate
949: .seealso: NEPSetSplitOperator(), NEPGetSplitOperatorTerm()
950: @*/
951: PetscErrorCode NEPGetSplitOperatorInfo(NEP nep,PetscInt *n,MatStructure *str)952: {
955: NEPCheckSplit(nep,1);
956: if (n) *n = nep->nt;
957: if (str) *str = nep->mstr;
958: return(0);
959: }
963: /*@C
964: NEPSetDerivatives - Sets the function to compute the k-th derivative T^(k)(lambda)
965: for any value of k (including 0), as well as the location to store the matrix.
967: Logically Collective on NEP and Mat
969: Input Parameters:
970: + nep - the NEP context
971: . A - the matrix to store the computed derivative
972: . der - routing to evaluate the k-th derivative (if NULL then NEP retains any
973: previously set value)
974: - ctx - [optional] user-defined context for private data for the derivatives
975: evaluation routine (may be NULL) (if NULL then NEP retains any
976: previously set value)
978: Level: beginner
980: .seealso: NEPSetFunction(), NEPGetDerivatives()
981: @*/
982: PetscErrorCode NEPSetDerivatives(NEP nep,Mat A,PetscErrorCode (*der)(NEP,PetscScalar,PetscInt,Mat,void*),void *ctx)983: {
991: if (nep->fui && nep->fui!=NEP_USER_INTERFACE_DERIVATIVES) { /* clean previous user info */
992: NEPReset_Problem(nep);
993: }
995: if (der) nep->computederivatives = der;
996: if (ctx) nep->derivativesctx = ctx;
997: if (A) {
998: PetscObjectReference((PetscObject)A);
999: MatDestroy(&nep->derivatives);
1000: nep->derivatives = A;
1001: }
1002: nep->fui = NEP_USER_INTERFACE_DERIVATIVES;
1003: return(0);
1004: }
1008: /*@C
1009: NEPGetDerivatives - Returns the derivatives matrix and optionally the user
1010: provided routine and context for evaluating the derivatives.
1012: Not Collective, but Mat object will be parallel if NEP object is
1014: Input Parameter:
1015: . nep - the nonlinear eigensolver context
1017: Output Parameters:
1018: + A - location to stash the derivatives matrix (or NULL)
1019: . der - location to put derivatives function (or NULL)
1020: - ctx - location to stash derivatives context (or NULL)
1022: Level: advanced
1024: .seealso: NEPSetDerivatives()
1025: @*/
1026: PetscErrorCode NEPGetDerivatives(NEP nep,Mat *A,PetscErrorCode (**der)(NEP,PetscScalar,PetscInt,Mat,void*),void **ctx)1027: {
1030: NEPCheckDerivatives(nep,1);
1031: if (A) *A = nep->derivatives;
1032: if (der) *der = nep->computederivatives;
1033: if (ctx) *ctx = nep->derivativesctx;
1034: return(0);
1035: }