Actual source code: stoar.c
slepc-3.13.3 2020-06-14
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, 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: SLEPc polynomial eigensolver: "stoar"
13: Method: S-TOAR
15: Algorithm:
17: Symmetric Two-Level Orthogonal Arnoldi.
19: References:
21: [1] C. Campos and J.E. Roman, "Restarted Q-Arnoldi-type methods
22: exploiting symmetry in quadratic eigenvalue problems", BIT
23: Numer. Math. 56(4):1213-1236, 2016.
24: */
26: #include <slepc/private/pepimpl.h>
27: #include "../src/pep/impls/krylov/pepkrylov.h"
28: #include <slepcblaslapack.h>
30: static PetscBool cited = PETSC_FALSE;
31: static const char citation[] =
32: "@Article{slepc-stoar,\n"
33: " author = \"C. Campos and J. E. Roman\",\n"
34: " title = \"Restarted {Q-Arnoldi-type} methods exploiting symmetry in quadratic eigenvalue problems\",\n"
35: " journal = \"{BIT} Numer. Math.\",\n"
36: " volume = \"56\",\n"
37: " number = \"4\",\n"
38: " pages = \"1213--1236\",\n"
39: " year = \"2016,\"\n"
40: " doi = \"https://doi.org/10.1007/s10543-016-0601-5\"\n"
41: "}\n";
43: typedef struct {
44: PetscReal scal[2];
45: Mat A[2];
46: Vec t;
47: } PEP_STOAR_MATSHELL;
49: static PetscErrorCode MatMult_STOAR(Mat A,Vec x,Vec y)
50: {
51: PetscErrorCode ierr;
52: PEP_STOAR_MATSHELL *ctx;
55: MatShellGetContext(A,&ctx);
56: MatMult(ctx->A[0],x,y);
57: VecScale(y,ctx->scal[0]);
58: if (ctx->scal[1]) {
59: MatMult(ctx->A[1],x,ctx->t);
60: VecAXPY(y,ctx->scal[1],ctx->t);
61: }
62: return(0);
63: }
65: static PetscErrorCode MatDestroy_STOAR(Mat A)
66: {
67: PEP_STOAR_MATSHELL *ctx;
68: PetscErrorCode ierr;
71: MatShellGetContext(A,(void**)&ctx);
72: VecDestroy(&ctx->t);
73: PetscFree(ctx);
74: return(0);
75: }
77: PetscErrorCode PEPSTOARSetUpInnerMatrix(PEP pep,Mat *B)
78: {
79: Mat pB[4],Bs[3],D[3];
80: PetscInt i,j,n,m;
81: PEP_STOAR_MATSHELL *ctxMat[3];
82: PEP_STOAR *ctx=(PEP_STOAR*)pep->data;
83: PetscErrorCode ierr;
86: for (i=0;i<3;i++) {
87: STGetMatrixTransformed(pep->st,i,&D[i]); /* D[2] = M */
88: }
89: MatGetLocalSize(D[2],&m,&n);
91: for (j=0;j<3;j++) {
92: PetscNew(ctxMat+j);
93: MatCreateShell(PetscObjectComm((PetscObject)pep),m,n,PETSC_DETERMINE,PETSC_DETERMINE,ctxMat[j],&Bs[j]);
94: MatShellSetOperation(Bs[j],MATOP_MULT,(void(*)(void))MatMult_STOAR);
95: MatShellSetOperation(Bs[j],MATOP_DESTROY,(void(*)(void))MatDestroy_STOAR);
96: }
97: for (i=0;i<4;i++) pB[i] = NULL;
98: if (ctx->alpha) {
99: ctxMat[0]->A[0] = D[0]; ctxMat[0]->scal[0] = ctx->alpha; ctxMat[0]->scal[1] = 0.0;
100: ctxMat[2]->A[0] = D[2]; ctxMat[2]->scal[0] = -ctx->alpha*pep->sfactor*pep->sfactor; ctxMat[2]->scal[1] = 0.0;
101: pB[0] = Bs[0]; pB[3] = Bs[2];
102: }
103: if (ctx->beta) {
104: i = (ctx->alpha)?1:0;
105: ctxMat[0]->scal[1] = 0.0;
106: ctxMat[0]->A[i] = D[1]; ctxMat[0]->scal[i] = -ctx->beta*pep->sfactor;
107: ctxMat[1]->A[0] = D[2]; ctxMat[1]->scal[0] = -ctx->beta*pep->sfactor*pep->sfactor; ctxMat[1]->scal[1] = 0.0;
108: pB[0] = Bs[0]; pB[1] = pB[2] = Bs[1];
109: }
110: BVCreateVec(pep->V,&ctxMat[0]->t);
111: MatCreateNest(PetscObjectComm((PetscObject)pep),2,NULL,2,NULL,pB,B);
112: for (j=0;j<3;j++) { MatDestroy(&Bs[j]); }
113: return(0);
114: }
116: PetscErrorCode PEPSetUp_STOAR(PEP pep)
117: {
118: PetscErrorCode ierr;
119: PetscBool shift,sinv,flg;
120: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
121: PetscInt ld,i;
122: PetscReal eta;
123: BVOrthogType otype;
124: BVOrthogBlockType obtype;
127: if (pep->problem_type!=PEP_HERMITIAN && pep->problem_type!=PEP_HYPERBOLIC) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Requested method is only available for Hermitian (or hyperbolic) problems");
128: if (pep->nmat!=3) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Solver only available for quadratic problems");
129: if (pep->basis!=PEP_BASIS_MONOMIAL) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Solver not implemented for non-monomial bases");
130: /* spectrum slicing requires special treatment of default values */
131: PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv);
132: if (pep->which==PEP_ALL) {
133: if (pep->stopping!=PEPStoppingBasic) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Spectrum slicing does not support user-defined stopping test");
134: pep->ops->solve = PEPSolve_STOAR_QSlice;
135: pep->ops->extractvectors = NULL;
136: pep->ops->setdefaultst = NULL;
137: PEPSetUp_STOAR_QSlice(pep);
138: } else {
139: PEPSetDimensions_Default(pep,pep->nev,&pep->ncv,&pep->mpd);
140: if (!ctx->lock && pep->mpd<pep->ncv) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Should not use mpd parameter in non-locking variant");
141: if (pep->max_it==PETSC_DEFAULT) pep->max_it = PetscMax(100,2*(pep->nmat-1)*pep->n/pep->ncv);
142: pep->ops->solve = PEPSolve_STOAR;
143: ld = pep->ncv+2;
144: DSSetType(pep->ds,DSGHIEP);
145: DSSetCompact(pep->ds,PETSC_TRUE);
146: DSAllocate(pep->ds,ld);
147: PEPBasisCoefficients(pep,pep->pbc);
148: STGetTransform(pep->st,&flg);
149: if (!flg) {
150: PetscFree(pep->solvematcoeffs);
151: PetscMalloc1(pep->nmat,&pep->solvematcoeffs);
152: PetscLogObjectMemory((PetscObject)pep,pep->nmat*sizeof(PetscScalar));
153: if (sinv) {
154: PEPEvaluateBasis(pep,pep->target,0,pep->solvematcoeffs,NULL);
155: } else {
156: for (i=0;i<pep->nmat-1;i++) pep->solvematcoeffs[i] = 0.0;
157: pep->solvematcoeffs[pep->nmat-1] = 1.0;
158: }
159: }
160: }
162: pep->lineariz = PETSC_TRUE;
163: PetscObjectTypeCompare((PetscObject)pep->st,STSHIFT,&shift);
164: PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv);
165: if (!shift && !sinv) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"Only STSHIFT and STSINVERT spectral transformations can be used");
166: if (!pep->which) {
167: if (sinv) pep->which = PEP_TARGET_MAGNITUDE;
168: else pep->which = PEP_LARGEST_MAGNITUDE;
169: }
171: PEPAllocateSolution(pep,2);
172: PEPSetWorkVecs(pep,4);
173: BVDestroy(&ctx->V);
174: BVCreateTensor(pep->V,pep->nmat-1,&ctx->V);
175: BVGetOrthogonalization(pep->V,&otype,NULL,&eta,&obtype);
176: BVSetOrthogonalization(ctx->V,otype,BV_ORTHOG_REFINE_ALWAYS,eta,obtype);
177: return(0);
178: }
180: /*
181: Compute a run of Lanczos iterations. dim(work)=(ctx->ld)*4
182: */
183: static PetscErrorCode PEPSTOARrun(PEP pep,PetscReal *a,PetscReal *b,PetscReal *omega,PetscInt k,PetscInt *M,PetscBool *breakdown,PetscBool *symmlost,Vec *t_)
184: {
186: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
187: PetscInt i,j,m=*M,l,lock;
188: PetscInt lds,d,ld,offq,nqt,ldds;
189: Vec v=t_[0],t=t_[1],q=t_[2];
190: PetscReal norm,sym=0.0,fro=0.0,*f;
191: PetscScalar *y,*S,*x,sigma;
192: PetscBLASInt j_,one=1;
193: PetscBool lindep,flg,sinvert=PETSC_FALSE;
194: Mat MS;
197: PetscMalloc1(*M,&y);
198: BVGetSizes(pep->V,NULL,NULL,&ld);
199: BVTensorGetDegree(ctx->V,&d);
200: BVGetActiveColumns(pep->V,&lock,&nqt);
201: lds = d*ld;
202: offq = ld;
203: DSGetLeadingDimension(pep->ds,&ldds);
204: *breakdown = PETSC_FALSE; /* ----- */
205: DSGetDimensions(pep->ds,NULL,NULL,&l,NULL,NULL);
206: BVSetActiveColumns(ctx->V,0,m);
207: BVSetActiveColumns(pep->V,0,nqt);
208: STGetTransform(pep->st,&flg);
209: if (!flg) {
210: /* spectral transformation handled by the solver */
211: PetscObjectTypeCompareAny((PetscObject)pep->st,&flg,STSINVERT,STSHIFT,"");
212: if (!flg) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_SUP,"ST type not supported for TOAR without transforming matrices");
213: PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinvert);
214: STGetShift(pep->st,&sigma);
215: }
216: for (j=k;j<m;j++) {
217: /* apply operator */
218: BVTensorGetFactors(ctx->V,NULL,&MS);
219: MatDenseGetArray(MS,&S);
220: BVGetColumn(pep->V,nqt,&t);
221: BVMultVec(pep->V,1.0,0.0,v,S+j*lds);
222: if (!sinvert) {
223: STMatMult(pep->st,0,v,q);
224: BVMultVec(pep->V,1.0,0.0,v,S+offq+j*lds);
225: STMatMult(pep->st,1,v,t);
226: VecAXPY(q,pep->sfactor,t);
227: if (ctx->beta && ctx->alpha) {
228: STMatMult(pep->st,2,v,t);
229: VecAXPY(q,-pep->sfactor*pep->sfactor*ctx->beta/ctx->alpha,t);
230: }
231: STMatSolve(pep->st,q,t);
232: VecScale(t,-1.0/(pep->sfactor*pep->sfactor));
233: } else {
234: STMatMult(pep->st,1,v,q);
235: STMatMult(pep->st,2,v,t);
236: VecAXPY(q,sigma*pep->sfactor,t);
237: VecScale(q,pep->sfactor);
238: BVMultVec(pep->V,1.0,0.0,v,S+offq+j*lds);
239: STMatMult(pep->st,2,v,t);
240: VecAXPY(q,pep->sfactor*pep->sfactor,t);
241: STMatSolve(pep->st,q,t);
242: VecScale(t,-1.0);
243: }
244: BVRestoreColumn(pep->V,nqt,&t);
246: /* orthogonalize */
247: if (!sinvert) x = S+offq+(j+1)*lds;
248: else x = S+(j+1)*lds;
249: BVOrthogonalizeColumn(pep->V,nqt,x,&norm,&lindep);
251: if (!lindep) {
252: if (!sinvert) *(S+offq+(j+1)*lds+nqt) = norm;
253: else *(S+(j+1)*lds+nqt) = norm;
254: BVScaleColumn(pep->V,nqt,1.0/norm);
255: nqt++;
256: }
257: if (!sinvert) {
258: for (i=0;i<=nqt-1;i++) *(S+(j+1)*lds+i) = *(S+offq+j*lds+i);
259: if (ctx->beta && ctx->alpha) {
260: for (i=0;i<=nqt-1;i++) *(S+(j+1)*lds+offq+i) -= *(S+(j+1)*lds+i)*ctx->beta/ctx->alpha;
261: }
262: } else for (i=0;i<nqt;i++) *(S+(j+1)*lds+offq+i) = *(S+j*lds+i)+sigma*(*(S+(j+1)*lds+i));
263: BVSetActiveColumns(pep->V,0,nqt);
264: MatDenseRestoreArray(MS,&S);
265: BVTensorRestoreFactors(ctx->V,NULL,&MS);
267: /* level-2 orthogonalization */
268: BVOrthogonalizeColumn(ctx->V,j+1,y,&norm,&lindep);
269: a[j] = PetscRealPart(y[j]);
270: omega[j+1] = (norm > 0)?1.0:-1.0;
271: BVScaleColumn(ctx->V,j+1,1.0/norm);
272: b[j] = PetscAbsReal(norm);
274: /* check symmetry */
275: DSGetArrayReal(pep->ds,DS_MAT_T,&f);
276: if (j==k) {
277: for (i=l;i<j-1;i++) y[i] = PetscAbsScalar(y[i])-PetscAbsReal(f[2*ldds+i]);
278: for (i=0;i<l;i++) y[i] = 0.0;
279: }
280: DSRestoreArrayReal(pep->ds,DS_MAT_T,&f);
281: if (j>0) y[j-1] = PetscAbsScalar(y[j-1])-PetscAbsReal(b[j-1]);
282: PetscBLASIntCast(j,&j_);
283: sym = SlepcAbs(BLASnrm2_(&j_,y,&one),sym);
284: fro = SlepcAbs(fro,SlepcAbs(a[j],b[j]));
285: if (j>0) fro = SlepcAbs(fro,b[j-1]);
286: if (sym/fro>PetscMax(PETSC_SQRT_MACHINE_EPSILON,10*pep->tol)) {
287: *symmlost = PETSC_TRUE;
288: *M=j;
289: break;
290: }
291: }
292: BVSetActiveColumns(pep->V,lock,nqt);
293: BVSetActiveColumns(ctx->V,0,*M);
294: PetscFree(y);
295: return(0);
296: }
298: #if 0
299: static PetscErrorCode PEPSTOARpreKConvergence(PEP pep,PetscInt nv,PetscReal *norm,Vec *w)
300: {
302: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
303: PetscBLASInt n_,one=1;
304: PetscInt lds=2*ctx->ld;
305: PetscReal t1,t2;
306: PetscScalar *S=ctx->S;
309: PetscBLASIntCast(nv+2,&n_);
310: t1 = BLASnrm2_(&n_,S+nv*2*ctx->ld,&one);
311: t2 = BLASnrm2_(&n_,S+(nv*2+1)*ctx->ld,&one);
312: *norm = SlepcAbs(t1,t2);
313: BVSetActiveColumns(pep->V,0,nv+2);
314: BVMultVec(pep->V,1.0,0.0,w[1],S+nv*lds);
315: STMatMult(pep->st,0,w[1],w[2]);
316: VecNorm(w[2],NORM_2,&t1);
317: BVMultVec(pep->V,1.0,0.0,w[1],S+ctx->ld+nv*lds);
318: STMatMult(pep->st,2,w[1],w[2]);
319: VecNorm(w[2],NORM_2,&t2);
320: t2 *= pep->sfactor*pep->sfactor;
321: *norm = PetscMax(*norm,SlepcAbs(t1,t2));
322: return(0);
323: }
324: #endif
326: PetscErrorCode PEPSolve_STOAR(PEP pep)
327: {
329: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
330: PetscInt j,k,l,nv=0,ld,ldds,t,nq=0;
331: PetscInt nconv=0,deg=pep->nmat-1;
332: PetscScalar *Q,*om,sigma;
333: PetscReal beta,norm=1.0,*omega,*a,*b,*r;
334: PetscBool breakdown,symmlost=PETSC_FALSE,sinv=PETSC_FALSE,falselock=PETSC_TRUE,flg;
335: Mat MQ,A;
336: Vec vomega;
339: PetscCitationsRegister(citation,&cited);
340: PEPSTOARSetUpInnerMatrix(pep,&A);
341: BVSetMatrix(ctx->V,A,PETSC_TRUE);
342: MatDestroy(&A);
343: if (ctx->lock) {
344: /* undocumented option to use a cheaper locking instead of the true locking */
345: PetscOptionsGetBool(NULL,NULL,"-pep_stoar_falselocking",&falselock,NULL);
346: }
347: BVGetSizes(pep->V,NULL,NULL,&ld);
348: STGetShift(pep->st,&sigma);
349: STGetTransform(pep->st,&flg);
350: if (pep->sfactor!=1.0) {
351: if (!flg) {
352: pep->target /= pep->sfactor;
353: RGPushScale(pep->rg,1.0/pep->sfactor);
354: STScaleShift(pep->st,1.0/pep->sfactor);
355: sigma /= pep->sfactor;
356: } else {
357: PetscObjectTypeCompare((PetscObject)pep->st,STSINVERT,&sinv);
358: pep->target = sinv?pep->target*pep->sfactor:pep->target/pep->sfactor;
359: RGPushScale(pep->rg,sinv?pep->sfactor:1.0/pep->sfactor);
360: STScaleShift(pep->st,sinv?pep->sfactor:1.0/pep->sfactor);
361: }
362: }
363: if (flg) sigma = 0.0;
365: /* Get the starting Arnoldi vector */
366: BVTensorBuildFirstColumn(ctx->V,pep->nini);
367: DSGetArrayReal(pep->ds,DS_MAT_D,&omega);
368: VecCreateSeq(PETSC_COMM_SELF,1,&vomega);
369: BVSetActiveColumns(ctx->V,0,1);
370: BVGetSignature(ctx->V,vomega);
371: VecGetArray(vomega,&om);
372: omega[0] = PetscRealPart(om[0]);
373: VecRestoreArray(vomega,&om);
374: DSRestoreArrayReal(pep->ds,DS_MAT_D,&omega);
375: VecDestroy(&vomega);
377: /* Restart loop */
378: l = 0;
379: DSGetLeadingDimension(pep->ds,&ldds);
380: while (pep->reason == PEP_CONVERGED_ITERATING) {
381: pep->its++;
382: DSGetArrayReal(pep->ds,DS_MAT_T,&a);
383: b = a+ldds;
384: DSGetArrayReal(pep->ds,DS_MAT_D,&omega);
386: /* Compute an nv-step Lanczos factorization */
387: nv = PetscMin(pep->nconv+pep->mpd,pep->ncv);
388: PEPSTOARrun(pep,a,b,omega,pep->nconv+l,&nv,&breakdown,&symmlost,pep->work);
389: beta = b[nv-1];
390: if (symmlost && nv==pep->nconv+l) {
391: pep->reason = PEP_DIVERGED_SYMMETRY_LOST;
392: pep->nconv = nconv;
393: if (falselock || !ctx->lock) {
394: BVSetActiveColumns(ctx->V,0,pep->nconv);
395: BVTensorCompress(ctx->V,0);
396: }
397: break;
398: }
399: DSRestoreArrayReal(pep->ds,DS_MAT_T,&a);
400: DSRestoreArrayReal(pep->ds,DS_MAT_D,&omega);
401: DSSetDimensions(pep->ds,nv,0,pep->nconv,pep->nconv+l);
402: if (l==0) {
403: DSSetState(pep->ds,DS_STATE_INTERMEDIATE);
404: } else {
405: DSSetState(pep->ds,DS_STATE_RAW);
406: }
408: /* Solve projected problem */
409: DSSolve(pep->ds,pep->eigr,pep->eigi);
410: DSSort(pep->ds,pep->eigr,pep->eigi,NULL,NULL,NULL);
411: DSSynchronize(pep->ds,pep->eigr,pep->eigi);
413: /* Check convergence */
414: /* PEPSTOARpreKConvergence(pep,nv,&norm,pep->work);*/
415: norm = 1.0;
416: DSGetDimensions(pep->ds,NULL,NULL,NULL,NULL,&t);
417: PEPKrylovConvergence(pep,PETSC_FALSE,pep->nconv,t-pep->nconv,PetscAbsReal(beta)*norm,&k);
418: (*pep->stopping)(pep,pep->its,pep->max_it,k,pep->nev,&pep->reason,pep->stoppingctx);
420: /* Update l */
421: if (pep->reason != PEP_CONVERGED_ITERATING || breakdown) l = 0;
422: else {
423: l = PetscMax(1,(PetscInt)((nv-k)/2));
424: l = PetscMin(l,t);
425: if (!breakdown) {
426: DSGetArrayReal(pep->ds,DS_MAT_T,&a);
427: if (*(a+ldds+k+l-1)!=0) {
428: if (k+l<nv-1) l = l+1;
429: else l = l-1;
430: }
431: /* Prepare the Rayleigh quotient for restart */
432: DSGetArray(pep->ds,DS_MAT_Q,&Q);
433: DSGetArrayReal(pep->ds,DS_MAT_D,&omega);
434: r = a + 2*ldds;
435: for (j=k;j<k+l;j++) {
436: r[j] = PetscRealPart(Q[nv-1+j*ldds]*beta);
437: }
438: b = a+ldds;
439: b[k+l-1] = r[k+l-1];
440: omega[k+l] = omega[nv];
441: DSRestoreArray(pep->ds,DS_MAT_Q,&Q);
442: DSRestoreArrayReal(pep->ds,DS_MAT_T,&a);
443: DSRestoreArrayReal(pep->ds,DS_MAT_D,&omega);
444: }
445: }
446: nconv = k;
447: if (!ctx->lock && pep->reason == PEP_CONVERGED_ITERATING && !breakdown) { l += k; k = 0; } /* non-locking variant: reset no. of converged pairs */
449: /* Update S */
450: DSGetMat(pep->ds,DS_MAT_Q,&MQ);
451: BVMultInPlace(ctx->V,MQ,pep->nconv,k+l);
452: MatDestroy(&MQ);
454: /* Copy last column of S */
455: BVCopyColumn(ctx->V,nv,k+l);
456: DSGetArrayReal(pep->ds,DS_MAT_D,&omega);
457: VecCreateSeq(PETSC_COMM_SELF,k+l,&vomega);
458: VecGetArray(vomega,&om);
459: for (j=0;j<k+l;j++) om[j] = omega[j];
460: VecRestoreArray(vomega,&om);
461: BVSetActiveColumns(ctx->V,0,k+l);
462: BVSetSignature(ctx->V,vomega);
463: VecDestroy(&vomega);
464: DSRestoreArrayReal(pep->ds,DS_MAT_D,&omega);
466: if (breakdown && pep->reason == PEP_CONVERGED_ITERATING) {
467: /* stop if breakdown */
468: PetscInfo2(pep,"Breakdown TOAR method (it=%D norm=%g)\n",pep->its,(double)beta);
469: pep->reason = PEP_DIVERGED_BREAKDOWN;
470: }
471: if (pep->reason != PEP_CONVERGED_ITERATING) l--;
472: BVGetActiveColumns(pep->V,NULL,&nq);
473: if (k+l+deg<=nq) {
474: BVSetActiveColumns(ctx->V,pep->nconv,k+l+1);
475: if (!falselock && ctx->lock) {
476: BVTensorCompress(ctx->V,k-pep->nconv);
477: } else {
478: BVTensorCompress(ctx->V,0);
479: }
480: }
481: pep->nconv = k;
482: PEPMonitor(pep,pep->its,nconv,pep->eigr,pep->eigi,pep->errest,nv);
483: }
485: if (pep->nconv>0) {
486: BVSetActiveColumns(ctx->V,0,pep->nconv);
487: BVGetActiveColumns(pep->V,NULL,&nq);
488: BVSetActiveColumns(pep->V,0,nq);
489: if (nq>pep->nconv) {
490: BVTensorCompress(ctx->V,pep->nconv);
491: BVSetActiveColumns(pep->V,0,pep->nconv);
492: }
493: }
494: STGetTransform(pep->st,&flg);
495: if (!flg && pep->ops->backtransform) {
496: (*pep->ops->backtransform)(pep);
497: }
498: if (pep->sfactor!=1.0) {
499: for (j=0;j<pep->nconv;j++) {
500: pep->eigr[j] *= pep->sfactor;
501: pep->eigi[j] *= pep->sfactor;
502: }
503: }
504: /* restore original values */
505: if (!flg) {
506: pep->target *= pep->sfactor;
507: STScaleShift(pep->st,pep->sfactor);
508: } else {
509: STScaleShift(pep->st,sinv?1.0/pep->sfactor:pep->sfactor);
510: pep->target = (sinv)?pep->target/pep->sfactor:pep->target*pep->sfactor;
511: }
512: if (pep->sfactor!=1.0) { RGPopScale(pep->rg); }
514: /* truncate Schur decomposition and change the state to raw so that
515: DSVectors() computes eigenvectors from scratch */
516: DSSetDimensions(pep->ds,pep->nconv,0,0,0);
517: DSSetState(pep->ds,DS_STATE_RAW);
518: return(0);
519: }
521: PetscErrorCode PEPSetFromOptions_STOAR(PetscOptionItems *PetscOptionsObject,PEP pep)
522: {
524: PetscBool flg,lock,b,f1,f2,f3;
525: PetscInt i,j,k;
526: PetscReal array[2]={0,0};
527: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
530: PetscOptionsHead(PetscOptionsObject,"PEP STOAR Options");
532: PetscOptionsBool("-pep_stoar_locking","Choose between locking and non-locking variants","PEPSTOARSetLocking",PETSC_FALSE,&lock,&flg);
533: if (flg) { PEPSTOARSetLocking(pep,lock); }
535: b = ctx->detect;
536: PetscOptionsBool("-pep_stoar_detect_zeros","Check zeros during factorizations at interval boundaries","PEPSTOARSetDetectZeros",ctx->detect,&b,&flg);
537: if (flg) { PEPSTOARSetDetectZeros(pep,b); }
539: i = 1;
540: j = k = PETSC_DECIDE;
541: PetscOptionsInt("-pep_stoar_nev","Number of eigenvalues to compute in each subsolve (only for spectrum slicing)","PEPSTOARSetDimensions",20,&i,&f1);
542: PetscOptionsInt("-pep_stoar_ncv","Number of basis vectors in each subsolve (only for spectrum slicing)","PEPSTOARSetDimensions",40,&j,&f2);
543: PetscOptionsInt("-pep_stoar_mpd","Maximum dimension of projected problem in each subsolve (only for spectrum slicing)","PEPSTOARSetDimensions",40,&k,&f3);
544: if (f1 || f2 || f3) { PEPSTOARSetDimensions(pep,i,j,k); }
546: k = 2;
547: PetscOptionsRealArray("-pep_stoar_linearization","Parameters of the linearization","PEPSTOARSetLinearization",array,&k,&flg);
548: if (flg) {
549: PEPSTOARSetLinearization(pep,array[0],array[1]);
550: }
552: b = ctx->checket;
553: PetscOptionsBool("-pep_stoar_check_eigenvalue_type","Check eigenvalue type during spectrum slicing","PEPSTOARSetCheckEigenvalueType",ctx->checket,&b,&flg);
554: if (flg) { PEPSTOARSetCheckEigenvalueType(pep,b); }
556: PetscOptionsTail();
557: return(0);
558: }
560: static PetscErrorCode PEPSTOARSetLocking_STOAR(PEP pep,PetscBool lock)
561: {
562: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
565: ctx->lock = lock;
566: return(0);
567: }
569: /*@
570: PEPSTOARSetLocking - Choose between locking and non-locking variants of
571: the STOAR method.
573: Logically Collective on pep
575: Input Parameters:
576: + pep - the eigenproblem solver context
577: - lock - true if the locking variant must be selected
579: Options Database Key:
580: . -pep_stoar_locking - Sets the locking flag
582: Notes:
583: The default is to lock converged eigenpairs when the method restarts.
584: This behaviour can be changed so that all directions are kept in the
585: working subspace even if already converged to working accuracy (the
586: non-locking variant).
588: Level: advanced
590: .seealso: PEPSTOARGetLocking()
591: @*/
592: PetscErrorCode PEPSTOARSetLocking(PEP pep,PetscBool lock)
593: {
599: PetscTryMethod(pep,"PEPSTOARSetLocking_C",(PEP,PetscBool),(pep,lock));
600: return(0);
601: }
603: static PetscErrorCode PEPSTOARGetLocking_STOAR(PEP pep,PetscBool *lock)
604: {
605: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
608: *lock = ctx->lock;
609: return(0);
610: }
612: /*@
613: PEPSTOARGetLocking - Gets the locking flag used in the STOAR method.
615: Not Collective
617: Input Parameter:
618: . pep - the eigenproblem solver context
620: Output Parameter:
621: . lock - the locking flag
623: Level: advanced
625: .seealso: PEPSTOARSetLocking()
626: @*/
627: PetscErrorCode PEPSTOARGetLocking(PEP pep,PetscBool *lock)
628: {
634: PetscUseMethod(pep,"PEPSTOARGetLocking_C",(PEP,PetscBool*),(pep,lock));
635: return(0);
636: }
638: static PetscErrorCode PEPSTOARGetInertias_STOAR(PEP pep,PetscInt *n,PetscReal **shifts,PetscInt **inertias)
639: {
641: PetscInt i,numsh;
642: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
643: PEP_SR sr = ctx->sr;
646: if (!pep->state) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_WRONGSTATE,"Must call PEPSetUp() first");
647: if (!ctx->sr) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_WRONGSTATE,"Only available in interval computations, see PEPSetInterval()");
648: switch (pep->state) {
649: case PEP_STATE_INITIAL:
650: break;
651: case PEP_STATE_SETUP:
652: if (n) *n = 2;
653: if (shifts) {
654: PetscMalloc1(2,shifts);
655: (*shifts)[0] = pep->inta;
656: (*shifts)[1] = pep->intb;
657: }
658: if (inertias) {
659: PetscMalloc1(2,inertias);
660: (*inertias)[0] = (sr->dir==1)?sr->inertia0:sr->inertia1;
661: (*inertias)[1] = (sr->dir==1)?sr->inertia1:sr->inertia0;
662: }
663: break;
664: case PEP_STATE_SOLVED:
665: case PEP_STATE_EIGENVECTORS:
666: numsh = ctx->nshifts;
667: if (n) *n = numsh;
668: if (shifts) {
669: PetscMalloc1(numsh,shifts);
670: for (i=0;i<numsh;i++) (*shifts)[i] = ctx->shifts[i];
671: }
672: if (inertias) {
673: PetscMalloc1(numsh,inertias);
674: for (i=0;i<numsh;i++) (*inertias)[i] = ctx->inertias[i];
675: }
676: break;
677: }
678: return(0);
679: }
681: /*@C
682: PEPSTOARGetInertias - Gets the values of the shifts and their
683: corresponding inertias in case of doing spectrum slicing for a
684: computational interval.
686: Not Collective
688: Input Parameter:
689: . pep - the eigenproblem solver context
691: Output Parameters:
692: + n - number of shifts, including the endpoints of the interval
693: . shifts - the values of the shifts used internally in the solver
694: - inertias - the values of the inertia in each shift
696: Notes:
697: If called after PEPSolve(), all shifts used internally by the solver are
698: returned (including both endpoints and any intermediate ones). If called
699: before PEPSolve() and after PEPSetUp() then only the information of the
700: endpoints of subintervals is available.
702: This function is only available for spectrum slicing runs.
704: The returned arrays should be freed by the user. Can pass NULL in any of
705: the two arrays if not required.
707: Fortran Notes:
708: The calling sequence from Fortran is
709: .vb
710: PEPSTOARGetInertias(pep,n,shifts,inertias,ierr)
711: integer n
712: double precision shifts(*)
713: integer inertias(*)
714: .ve
715: The arrays should be at least of length n. The value of n can be determined
716: by an initial call
717: .vb
718: PEPSTOARGetInertias(pep,n,PETSC_NULL_REAL,PETSC_NULL_INTEGER,ierr)
719: .ve
721: Level: advanced
723: .seealso: PEPSetInterval()
724: @*/
725: PetscErrorCode PEPSTOARGetInertias(PEP pep,PetscInt *n,PetscReal **shifts,PetscInt **inertias)
726: {
732: PetscUseMethod(pep,"PEPSTOARGetInertias_C",(PEP,PetscInt*,PetscReal**,PetscInt**),(pep,n,shifts,inertias));
733: return(0);
734: }
736: static PetscErrorCode PEPSTOARSetDetectZeros_STOAR(PEP pep,PetscBool detect)
737: {
738: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
741: ctx->detect = detect;
742: pep->state = PEP_STATE_INITIAL;
743: return(0);
744: }
746: /*@
747: PEPSTOARSetDetectZeros - Sets a flag to enforce detection of
748: zeros during the factorizations throughout the spectrum slicing computation.
750: Logically Collective on pep
752: Input Parameters:
753: + pep - the eigenproblem solver context
754: - detect - check for zeros
756: Options Database Key:
757: . -pep_stoar_detect_zeros - Check for zeros; this takes an optional
758: bool value (0/1/no/yes/true/false)
760: Notes:
761: A zero in the factorization indicates that a shift coincides with an eigenvalue.
763: This flag is turned off by default, and may be necessary in some cases.
764: This feature currently requires an external package for factorizations
765: with support for zero detection, e.g. MUMPS.
767: Level: advanced
769: .seealso: PEPSetInterval()
770: @*/
771: PetscErrorCode PEPSTOARSetDetectZeros(PEP pep,PetscBool detect)
772: {
778: PetscTryMethod(pep,"PEPSTOARSetDetectZeros_C",(PEP,PetscBool),(pep,detect));
779: return(0);
780: }
782: static PetscErrorCode PEPSTOARGetDetectZeros_STOAR(PEP pep,PetscBool *detect)
783: {
784: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
787: *detect = ctx->detect;
788: return(0);
789: }
791: /*@
792: PEPSTOARGetDetectZeros - Gets the flag that enforces zero detection
793: in spectrum slicing.
795: Not Collective
797: Input Parameter:
798: . pep - the eigenproblem solver context
800: Output Parameter:
801: . detect - whether zeros detection is enforced during factorizations
803: Level: advanced
805: .seealso: PEPSTOARSetDetectZeros()
806: @*/
807: PetscErrorCode PEPSTOARGetDetectZeros(PEP pep,PetscBool *detect)
808: {
814: PetscUseMethod(pep,"PEPSTOARGetDetectZeros_C",(PEP,PetscBool*),(pep,detect));
815: return(0);
816: }
818: static PetscErrorCode PEPSTOARSetLinearization_STOAR(PEP pep,PetscReal alpha,PetscReal beta)
819: {
820: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
823: if (beta==0.0 && alpha==0.0) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_WRONG,"Parameters alpha and beta cannot be zero simultaneously");
824: ctx->alpha = alpha;
825: ctx->beta = beta;
826: return(0);
827: }
829: /*@
830: PEPSTOARSetLinearization - Set the coefficients that define
831: the linearization of a quadratic eigenproblem.
833: Logically Collective on pep
835: Input Parameters:
836: + pep - polynomial eigenvalue solver
837: . alpha - first parameter of the linearization
838: - beta - second parameter of the linearization
840: Options Database Key:
841: . -pep_stoar_linearization <alpha,beta> - Sets the coefficients
843: Notes:
844: Cannot pass zero for both alpha and beta. The default values are
845: alpha=1 and beta=0.
847: Level: advanced
849: .seealso: PEPSTOARGetLinearization()
850: @*/
851: PetscErrorCode PEPSTOARSetLinearization(PEP pep,PetscReal alpha,PetscReal beta)
852: {
859: PetscTryMethod(pep,"PEPSTOARSetLinearization_C",(PEP,PetscReal,PetscReal),(pep,alpha,beta));
860: return(0);
861: }
863: static PetscErrorCode PEPSTOARGetLinearization_STOAR(PEP pep,PetscReal *alpha,PetscReal *beta)
864: {
865: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
868: if (alpha) *alpha = ctx->alpha;
869: if (beta) *beta = ctx->beta;
870: return(0);
871: }
873: /*@
874: PEPSTOARGetLinearization - Returns the coefficients that define
875: the linearization of a quadratic eigenproblem.
877: Not Collective
879: Input Parameter:
880: . pep - polynomial eigenvalue solver
882: Output Parameters:
883: + alpha - the first parameter of the linearization
884: - beta - the second parameter of the linearization
886: Level: advanced
888: .seealso: PEPSTOARSetLinearization()
889: @*/
890: PetscErrorCode PEPSTOARGetLinearization(PEP pep,PetscReal *alpha,PetscReal *beta)
891: {
896: PetscUseMethod(pep,"PEPSTOARGetLinearization_C",(PEP,PetscReal*,PetscReal*),(pep,alpha,beta));
897: return(0);
898: }
900: static PetscErrorCode PEPSTOARSetDimensions_STOAR(PEP pep,PetscInt nev,PetscInt ncv,PetscInt mpd)
901: {
902: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
905: if (nev<1) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
906: ctx->nev = nev;
907: if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
908: ctx->ncv = PETSC_DEFAULT;
909: } else {
910: if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
911: ctx->ncv = ncv;
912: }
913: if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
914: ctx->mpd = PETSC_DEFAULT;
915: } else {
916: if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
917: ctx->mpd = mpd;
918: }
919: pep->state = PEP_STATE_INITIAL;
920: return(0);
921: }
923: /*@
924: PEPSTOARSetDimensions - Sets the dimensions used for each subsolve
925: step in case of doing spectrum slicing for a computational interval.
926: The meaning of the parameters is the same as in PEPSetDimensions().
928: Logically Collective on pep
930: Input Parameters:
931: + pep - the eigenproblem solver context
932: . nev - number of eigenvalues to compute
933: . ncv - the maximum dimension of the subspace to be used by the subsolve
934: - mpd - the maximum dimension allowed for the projected problem
936: Options Database Key:
937: + -eps_stoar_nev <nev> - Sets the number of eigenvalues
938: . -eps_stoar_ncv <ncv> - Sets the dimension of the subspace
939: - -eps_stoar_mpd <mpd> - Sets the maximum projected dimension
941: Level: advanced
943: .seealso: PEPSTOARGetDimensions(), PEPSetDimensions(), PEPSetInterval()
944: @*/
945: PetscErrorCode PEPSTOARSetDimensions(PEP pep,PetscInt nev,PetscInt ncv,PetscInt mpd)
946: {
954: PetscTryMethod(pep,"PEPSTOARSetDimensions_C",(PEP,PetscInt,PetscInt,PetscInt),(pep,nev,ncv,mpd));
955: return(0);
956: }
958: static PetscErrorCode PEPSTOARGetDimensions_STOAR(PEP pep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
959: {
960: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
963: if (nev) *nev = ctx->nev;
964: if (ncv) *ncv = ctx->ncv;
965: if (mpd) *mpd = ctx->mpd;
966: return(0);
967: }
969: /*@
970: PEPSTOARGetDimensions - Gets the dimensions used for each subsolve
971: step in case of doing spectrum slicing for a computational interval.
973: Not Collective
975: Input Parameter:
976: . pep - the eigenproblem solver context
978: Output Parameters:
979: + nev - number of eigenvalues to compute
980: . ncv - the maximum dimension of the subspace to be used by the subsolve
981: - mpd - the maximum dimension allowed for the projected problem
983: Level: advanced
985: .seealso: PEPSTOARSetDimensions()
986: @*/
987: PetscErrorCode PEPSTOARGetDimensions(PEP pep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
988: {
993: PetscUseMethod(pep,"PEPSTOARGetDimensions_C",(PEP,PetscInt*,PetscInt*,PetscInt*),(pep,nev,ncv,mpd));
994: return(0);
995: }
997: static PetscErrorCode PEPSTOARSetCheckEigenvalueType_STOAR(PEP pep,PetscBool checket)
998: {
999: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
1002: ctx->checket = checket;
1003: pep->state = PEP_STATE_INITIAL;
1004: return(0);
1005: }
1007: /*@
1008: PEPSTOARSetCheckEigenvalueType - Sets a flag to check that all the eigenvalues
1009: obtained throughout the spectrum slicing computation have the same definite type.
1011: Logically Collective on pep
1013: Input Parameters:
1014: + pep - the eigenproblem solver context
1015: - checket - check eigenvalue type
1017: Options Database Key:
1018: . -pep_stoar_check_eigenvalue_type - Check eigenvalue type; this takes an optional
1019: bool value (0/1/no/yes/true/false)
1021: Notes:
1022: This option is relevant only for spectrum slicing computations, but it is
1023: ignored if the problem type is PEP_HYPERBOLIC.
1025: This flag is turned on by default, to guarantee that the computed eigenvalues
1026: have the same type (otherwise the computed solution might be wrong). But since
1027: the check is computationally quite expensive, the check may be turned off if
1028: the user knows for sure that all eigenvalues in the requested interval have
1029: the same type.
1031: Level: advanced
1033: .seealso: PEPSetProblemType(), PEPSetInterval()
1034: @*/
1035: PetscErrorCode PEPSTOARSetCheckEigenvalueType(PEP pep,PetscBool checket)
1036: {
1042: PetscTryMethod(pep,"PEPSTOARSetCheckEigenvalueType_C",(PEP,PetscBool),(pep,checket));
1043: return(0);
1044: }
1046: static PetscErrorCode PEPSTOARGetCheckEigenvalueType_STOAR(PEP pep,PetscBool *checket)
1047: {
1048: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
1051: *checket = ctx->checket;
1052: return(0);
1053: }
1055: /*@
1056: PEPSTOARGetCheckEigenvalueType - Gets the flag for the eigenvalue type
1057: check in spectrum slicing.
1059: Not Collective
1061: Input Parameter:
1062: . pep - the eigenproblem solver context
1064: Output Parameter:
1065: . checket - whether eigenvalue type must be checked during spectrum slcing
1067: Level: advanced
1069: .seealso: PEPSTOARSetCheckEigenvalueType()
1070: @*/
1071: PetscErrorCode PEPSTOARGetCheckEigenvalueType(PEP pep,PetscBool *checket)
1072: {
1078: PetscUseMethod(pep,"PEPSTOARGetCheckEigenvalueType_C",(PEP,PetscBool*),(pep,checket));
1079: return(0);
1080: }
1082: PetscErrorCode PEPView_STOAR(PEP pep,PetscViewer viewer)
1083: {
1085: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
1086: PetscBool isascii;
1089: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
1090: if (isascii) {
1091: PetscViewerASCIIPrintf(viewer," using the %slocking variant\n",ctx->lock?"":"non-");
1092: PetscViewerASCIIPrintf(viewer," linearization parameters: alpha=%g beta=%g\n",(double)ctx->alpha,(double)ctx->beta);
1093: if (pep->which==PEP_ALL && !ctx->hyperbolic) {
1094: PetscViewerASCIIPrintf(viewer," checking eigenvalue type: %s\n",ctx->checket?"enabled":"disabled");
1095: }
1096: }
1097: return(0);
1098: }
1100: PetscErrorCode PEPReset_STOAR(PEP pep)
1101: {
1105: if (pep->which==PEP_ALL) {
1106: PEPReset_STOAR_QSlice(pep);
1107: }
1108: return(0);
1109: }
1111: PetscErrorCode PEPDestroy_STOAR(PEP pep)
1112: {
1114: PEP_STOAR *ctx = (PEP_STOAR*)pep->data;
1117: BVDestroy(&ctx->V);
1118: PetscFree(pep->data);
1119: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetLocking_C",NULL);
1120: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetLocking_C",NULL);
1121: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetDetectZeros_C",NULL);
1122: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetDetectZeros_C",NULL);
1123: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetInertias_C",NULL);
1124: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetDimensions_C",NULL);
1125: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetDimensions_C",NULL);
1126: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetLinearization_C",NULL);
1127: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetLinearization_C",NULL);
1128: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetCheckEigenvalueType_C",NULL);
1129: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetCheckEigenvalueType_C",NULL);
1130: return(0);
1131: }
1133: SLEPC_EXTERN PetscErrorCode PEPCreate_STOAR(PEP pep)
1134: {
1136: PEP_STOAR *ctx;
1139: PetscNewLog(pep,&ctx);
1140: pep->data = (void*)ctx;
1141: ctx->lock = PETSC_TRUE;
1142: ctx->nev = 1;
1143: ctx->ncv = PETSC_DEFAULT;
1144: ctx->mpd = PETSC_DEFAULT;
1145: ctx->alpha = 1.0;
1146: ctx->beta = 0.0;
1147: ctx->checket = PETSC_TRUE;
1149: pep->ops->setup = PEPSetUp_STOAR;
1150: pep->ops->setfromoptions = PEPSetFromOptions_STOAR;
1151: pep->ops->destroy = PEPDestroy_STOAR;
1152: pep->ops->view = PEPView_STOAR;
1153: pep->ops->backtransform = PEPBackTransform_Default;
1154: pep->ops->computevectors = PEPComputeVectors_Default;
1155: pep->ops->extractvectors = PEPExtractVectors_TOAR;
1156: pep->ops->reset = PEPReset_STOAR;
1158: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetLocking_C",PEPSTOARSetLocking_STOAR);
1159: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetLocking_C",PEPSTOARGetLocking_STOAR);
1160: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetDetectZeros_C",PEPSTOARSetDetectZeros_STOAR);
1161: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetDetectZeros_C",PEPSTOARGetDetectZeros_STOAR);
1162: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetInertias_C",PEPSTOARGetInertias_STOAR);
1163: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetDimensions_C",PEPSTOARGetDimensions_STOAR);
1164: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetDimensions_C",PEPSTOARSetDimensions_STOAR);
1165: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetLinearization_C",PEPSTOARSetLinearization_STOAR);
1166: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetLinearization_C",PEPSTOARGetLinearization_STOAR);
1167: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARSetCheckEigenvalueType_C",PEPSTOARSetCheckEigenvalueType_STOAR);
1168: PetscObjectComposeFunction((PetscObject)pep,"PEPSTOARGetCheckEigenvalueType_C",PEPSTOARGetCheckEigenvalueType_STOAR);
1169: return(0);
1170: }