Actual source code: lobpcg.c
slepc-3.6.1 2015-09-03
1: /*
3: SLEPc eigensolver: "lobpcg"
5: Method: Locally Optimal Block Preconditioned Conjugate Gradient
7: Algorithm:
9: LOBPCG with soft and hard locking. Follows the implementation
10: in BLOPEX [2].
12: References:
14: [1] A. V. Knyazev, "Toward the optimal preconditioned eigensolver:
15: locally optimal block preconditioned conjugate gradient method",
16: SIAM J. Sci. Comput. 23(2):517-541, 2001.
18: [2] A. V. Knyazev et al., "Block Locally Optimal Preconditioned
19: Eigenvalue Xolvers (BLOPEX) in Hypre and PETSc", SIAM J. Sci.
20: Comput. 29(5):2224-2239, 2007.
22: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
23: SLEPc - Scalable Library for Eigenvalue Problem Computations
24: Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain
26: This file is part of SLEPc.
28: SLEPc is free software: you can redistribute it and/or modify it under the
29: terms of version 3 of the GNU Lesser General Public License as published by
30: the Free Software Foundation.
32: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
33: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
34: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
35: more details.
37: You should have received a copy of the GNU Lesser General Public License
38: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
39: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
40: */
42: #include <slepc/private/epsimpl.h> /*I "slepceps.h" I*/
43: #include <slepc/private/dsimpl.h> /*I "slepcds.h" I*/
45: typedef struct {
46: PetscInt bs; /* block size */
47: PetscBool lock; /* soft locking active/inactive */
48: } EPS_LOBPCG;
52: PetscErrorCode EPSSetDimensions_LOBPCG(EPS eps,PetscInt nev,PetscInt *ncv,PetscInt *mpd)
53: {
54: EPS_LOBPCG *ctx = (EPS_LOBPCG*)eps->data;
55: PetscInt k;
58: k = PetscMax(3*ctx->bs,((eps->nev-1)/ctx->bs+3)*ctx->bs);
59: if (*ncv) { /* ncv set */
60: if (*ncv<k) SETERRQ(PetscObjectComm((PetscObject)eps),1,"The value of ncv is not sufficiently large");
61: } else *ncv = k;
62: if (!*mpd) *mpd = 3*ctx->bs;
63: else if (*mpd!=3*ctx->bs) SETERRQ(PetscObjectComm((PetscObject)eps),1,"This solver does not allow a value of mpd different from 3*blocksize");
64: return(0);
65: }
69: PetscErrorCode EPSSetUp_LOBPCG(EPS eps)
70: {
72: EPS_LOBPCG *ctx = (EPS_LOBPCG*)eps->data;
73: PetscBool precond,istrivial;
76: if (!eps->ishermitian || (eps->isgeneralized && !eps->ispositive)) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"LOBPCG only works for Hermitian problems");
77: if (!ctx->bs) ctx->bs = PetscMin(16,eps->nev);
78: if (eps->n-eps->nds<5*ctx->bs) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"The problem size is too small relative to the block size");
79: EPSSetDimensions_LOBPCG(eps,eps->nev,&eps->ncv,&eps->mpd);
80: if (!eps->max_it) eps->max_it = PetscMax(100,2*eps->n/eps->ncv);
81: if (!eps->which) eps->which = EPS_SMALLEST_REAL;
82: if (eps->which!=EPS_SMALLEST_REAL) SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->which");
83: if (eps->arbitrary) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Arbitrary selection of eigenpairs not supported in this solver");
84: if (eps->extraction) { PetscInfo(eps,"Warning: extraction type ignored\n"); }
85: RGIsTrivial(eps->rg,&istrivial);
86: if (!istrivial) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"This solver does not support region filtering");
88: STSetUp(eps->st);
89: PetscObjectTypeCompare((PetscObject)eps->st,STPRECOND,&precond);
90: if (!precond) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"LOBPCG only works with precond ST");
92: EPSAllocateSolution(eps,0);
93: EPS_SetInnerProduct(eps);
94: DSSetType(eps->ds,DSGHEP);
95: DSAllocate(eps->ds,eps->mpd);
96: EPSSetWorkVecs(eps,1);
97: return(0);
98: }
102: PetscErrorCode EPSSolve_LOBPCG(EPS eps)
103: {
105: EPS_LOBPCG *ctx = (EPS_LOBPCG*)eps->data;
106: PetscInt i,j,k,ld,nv,ini,kini,nmat,nc,nconv,bdone,its;
107: PetscReal norm;
108: PetscBool breakdown,countc;
109: Mat A,B,M;
110: Vec v,w=eps->work[0];
111: BV X,Y,Z,R,P,AX,BX;
114: DSGetLeadingDimension(eps->ds,&ld);
115: STGetNumMatrices(eps->st,&nmat);
116: STGetOperators(eps->st,0,&A);
117: if (nmat>1) { STGetOperators(eps->st,1,&B); }
118: else B = NULL;
120: /* 1. Allocate memory */
121: BVDuplicateResize(eps->V,3*ctx->bs,&Z);
122: BVDuplicateResize(eps->V,ctx->bs,&X);
123: BVDuplicateResize(eps->V,ctx->bs,&R);
124: BVDuplicateResize(eps->V,ctx->bs,&P);
125: BVDuplicateResize(eps->V,ctx->bs,&AX);
126: if (B) {
127: BVDuplicateResize(eps->V,ctx->bs,&BX);
128: }
129: nc = eps->nds;
130: if (nc>0 || eps->nev>ctx->bs) {
131: BVDuplicateResize(eps->V,nc+eps->nev,&Y);
132: }
133: if (nc>0) {
134: for (j=0;j<nc;j++) {
135: BVGetColumn(eps->V,-nc+j,&v);
136: BVInsertVec(Y,j,v);
137: BVRestoreColumn(eps->V,-nc+j,&v);
138: }
139: BVSetActiveColumns(Y,0,nc);
140: }
142: /* 2. Apply the constraints to the initial vectors */
143: kini = eps->nini;
144: while (kini<eps->ncv-2*ctx->bs) { /* Generate more initial vectors if necessary */
145: BVSetRandomColumn(eps->V,kini,eps->rand);
146: BVOrthogonalizeColumn(eps->V,kini,NULL,&norm,&breakdown);
147: if (norm>0.0 && !breakdown) {
148: BVScaleColumn(eps->V,kini,1.0/norm);
149: kini++;
150: }
151: }
152: nv = ctx->bs;
153: BVSetActiveColumns(eps->V,0,nv);
154: BVSetActiveColumns(Z,0,nv);
155: BVCopy(eps->V,Z);
156: BVCopy(Z,X);
158: /* 3. B-orthogonalize initial vectors */
159: if (B) {
160: BVMatMult(X,B,BX);
161: }
163: /* 4. Compute initial Ritz vectors */
164: BVMatMult(X,A,AX);
165: DSSetDimensions(eps->ds,nv,0,0,0);
166: DSGetMat(eps->ds,DS_MAT_A,&M);
167: BVMatProject(AX,NULL,X,M);
168: DSRestoreMat(eps->ds,DS_MAT_A,&M);
169: DSSetIdentity(eps->ds,DS_MAT_B);
170: DSSetState(eps->ds,DS_STATE_RAW);
171: DSSolve(eps->ds,eps->eigr,eps->eigi);
172: DSSort(eps->ds,eps->eigr,eps->eigi,NULL,NULL,NULL);
173: DSVectors(eps->ds,DS_MAT_X,NULL,NULL);
174: DSGetMat(eps->ds,DS_MAT_X,&M);
175: BVMultInPlace(X,M,0,nv);
176: BVMultInPlace(AX,M,0,nv);
177: DSRestoreMat(eps->ds,DS_MAT_X,&M);
179: /* 5. Initialize range of active iterates */
180: bdone = 0; /* finished blocks, the leading bdone*bs columns of V are eigenvectors */
181: nconv = 0; /* number of converged eigenvalues in the current block */
182: its = 0; /* iterations for the current block */
184: /* 6. Main loop */
185: while (eps->reason == EPS_CONVERGED_ITERATING) {
187: if (ctx->lock) {
188: BVSetActiveColumns(R,nconv,ctx->bs);
189: BVSetActiveColumns(AX,nconv,ctx->bs);
190: if (B) {
191: BVSetActiveColumns(BX,nconv,ctx->bs);
192: }
193: }
195: /* 7. Compute residuals */
196: DSGetMat(eps->ds,DS_MAT_A,&M);
197: BVCopy(AX,R);
198: if (B) {
199: BVMult(R,-1.0,1.0,BX,M);
200: } else {
201: BVMult(R,-1.0,1.0,X,M);
202: }
203: DSRestoreMat(eps->ds,DS_MAT_A,&M);
205: /* 8. Compute residual norms and update index set of active iterates */
206: ini = (ctx->lock)? nconv: 0;
207: k = ini;
208: countc = PETSC_TRUE;
209: for (j=ini;j<ctx->bs;j++) {
210: i = bdone*ctx->bs+j;
211: BVGetColumn(R,j,&v);
212: VecNorm(v,NORM_2,&norm);
213: BVRestoreColumn(R,j,&v);
214: (*eps->converged)(eps,eps->eigr[i],eps->eigi[i],norm,&eps->errest[i],eps->convergedctx);
215: if (countc) {
216: if (eps->errest[i] < eps->tol) k++;
217: else countc = PETSC_FALSE;
218: }
219: if (!countc && !eps->trackall) break;
220: }
221: nconv = k;
222: eps->nconv = bdone*ctx->bs + nconv;
223: if (eps->its+its) {
224: EPSMonitor(eps,eps->its+its,eps->nconv,eps->eigr,eps->eigi,eps->errest,(bdone+1)*ctx->bs);
225: }
226: if (eps->nconv >= eps->nev || nconv == ctx->bs) {
227: BVSetActiveColumns(eps->V,bdone*ctx->bs,bdone*ctx->bs+nconv);
228: BVSetActiveColumns(Z,0,nconv);
229: BVSetActiveColumns(X,0,nconv);
230: BVCopy(X,eps->V);
231: BVCopy(X,Z);
232: }
233: if (eps->nconv >= eps->nev) eps->reason = EPS_CONVERGED_TOL;
234: else if (eps->its+its >= eps->max_it) eps->reason = EPS_DIVERGED_ITS;
235: if (eps->nconv >= eps->nev || nconv == ctx->bs) {
236: eps->its += its;
237: its = 0;
238: }
239: its++;
240: if (eps->reason != EPS_CONVERGED_ITERATING) break;
242: if (nconv == ctx->bs) { /* block finished: lock eigenvectors and compute new R */
244: /* extend constraints */
245: BVSetActiveColumns(Y,nc+bdone*ctx->bs,nc+(bdone+1)*ctx->bs);
246: BVCopy(Z,Y);
247: for (j=0;j<ctx->bs;j++) {
248: BVSetActiveColumns(Y,0,nc+bdone*ctx->bs+j);
249: BVGetColumn(Y,nc+bdone*ctx->bs+j,&v);
250: BVOrthogonalizeVec(Y,v,NULL,NULL,NULL);
251: VecNormalize(v,NULL);
252: BVRestoreColumn(Y,nc+bdone*ctx->bs+j,&v);
253: }
254: BVSetActiveColumns(Y,0,nc+(bdone+1)*ctx->bs);
256: bdone++;
257: nconv = 0;
259: /* set new initial vectors */
260: BVSetActiveColumns(eps->V,bdone*ctx->bs,(bdone+1)*ctx->bs);
261: BVCopy(eps->V,Z);
262: for (j=0;j<ctx->bs;j++) {
263: BVGetColumn(Z,j,&v);
264: BVOrthogonalizeVec(Y,v,NULL,NULL,NULL);
265: VecNormalize(v,NULL);
266: BVRestoreColumn(Z,j,&v);
267: }
268: BVSetActiveColumns(X,nconv,ctx->bs);
269: BVSetActiveColumns(AX,nconv,ctx->bs);
270: BVCopy(Z,X);
272: /* B-orthogonalize initial vectors */
273: if (B) {
274: BVMatMult(X,B,BX);
275: }
277: /* Compute initial Ritz vectors */
278: nv = ctx->bs;
279: BVMatMult(X,A,AX);
280: DSSetDimensions(eps->ds,nv,0,0,0);
281: DSGetMat(eps->ds,DS_MAT_A,&M);
282: BVMatProject(AX,NULL,X,M);
283: DSRestoreMat(eps->ds,DS_MAT_A,&M);
284: DSGetMat(eps->ds,DS_MAT_B,&M);
285: if (B) {
286: BVMatProject(Z,B,Z,M);
287: } else {
288: BVDot(Z,Z,M);
289: }
290: DSRestoreMat(eps->ds,DS_MAT_B,&M);
291: DSSetState(eps->ds,DS_STATE_RAW);
292: DSSolve(eps->ds,eps->eigr+bdone*ctx->bs,eps->eigi);
293: DSSort(eps->ds,eps->eigr+bdone*ctx->bs,eps->eigi,NULL,NULL,NULL);
294: DSVectors(eps->ds,DS_MAT_X,NULL,NULL);
295: DSGetMat(eps->ds,DS_MAT_X,&M);
296: BVMultInPlace(Z,M,0,nv);
297: BVMultInPlace(X,M,0,nv);
298: BVMultInPlace(AX,M,0,nv);
299: DSRestoreMat(eps->ds,DS_MAT_X,&M);
301: EPSMonitor(eps,eps->its+its-1,eps->nconv,eps->eigr,eps->eigi,eps->errest,(bdone+1)*ctx->bs);
303: if (ctx->lock) {
304: BVSetActiveColumns(R,nconv,ctx->bs);
305: BVSetActiveColumns(AX,nconv,ctx->bs);
306: if (B) {
307: BVSetActiveColumns(BX,nconv,ctx->bs);
308: }
309: }
311: /* compute residuals */
312: DSGetMat(eps->ds,DS_MAT_A,&M);
313: BVCopy(AX,R);
314: if (B) {
315: BVMult(R,-1.0,1.0,BX,M);
316: } else {
317: BVMult(R,-1.0,1.0,X,M);
318: }
319: DSRestoreMat(eps->ds,DS_MAT_A,&M);
320: }
322: ini = (ctx->lock)? nconv: 0;
323: if (ctx->lock) {
324: BVSetActiveColumns(R,nconv,ctx->bs);
325: BVSetActiveColumns(P,nconv,ctx->bs);
326: BVSetActiveColumns(AX,nconv,ctx->bs);
327: if (B) {
328: BVSetActiveColumns(BX,nconv,ctx->bs);
329: }
330: }
332: /* 9. Apply preconditioner to the residuals */
333: for (j=ini;j<ctx->bs;j++) {
334: BVGetColumn(R,j,&v);
335: STMatSolve(eps->st,v,w);
336: if (nc+bdone*ctx->bs>0) {
337: BVOrthogonalizeVec(Y,w,NULL,NULL,NULL);
338: VecNormalize(w,NULL);
339: }
340: VecCopy(w,v);
341: BVRestoreColumn(R,j,&v);
342: }
344: /* 11. B-orthonormalize preconditioned residuals */
345: BVOrthogonalize(R,NULL);
347: /* 13-16. B-orthonormalize conjugate directions */
348: if (its>1) {
349: BVOrthogonalize(P,NULL);
350: }
352: /* 17-23. Compute symmetric Gram matrices */
353: BVSetActiveColumns(Z,0,ctx->bs);
354: BVSetActiveColumns(X,0,ctx->bs);
355: BVCopy(X,Z);
356: BVSetActiveColumns(Z,ctx->bs,2*ctx->bs-ini);
357: BVCopy(R,Z);
358: if (its>1) {
359: BVSetActiveColumns(Z,2*ctx->bs-ini,3*ctx->bs-2*ini);
360: BVCopy(P,Z);
361: }
363: if (its>1) nv = 3*ctx->bs-2*ini;
364: else nv = 2*ctx->bs-ini;
366: BVSetActiveColumns(Z,0,nv);
367: DSSetDimensions(eps->ds,nv,0,0,0);
368: DSGetMat(eps->ds,DS_MAT_A,&M);
369: BVMatProject(Z,A,Z,M);
370: DSRestoreMat(eps->ds,DS_MAT_A,&M);
371: DSGetMat(eps->ds,DS_MAT_B,&M);
372: if (B) {
373: BVMatProject(Z,B,Z,M);
374: } else {
375: BVDot(Z,Z,M);
376: }
377: DSRestoreMat(eps->ds,DS_MAT_B,&M);
378:
379: /* 24. Solve the generalized eigenvalue problem */
380: DSSetState(eps->ds,DS_STATE_RAW);
381: DSSolve(eps->ds,eps->eigr+bdone*ctx->bs,eps->eigi);
382: DSSort(eps->ds,eps->eigr+bdone*ctx->bs,eps->eigi,NULL,NULL,NULL);
383: DSVectors(eps->ds,DS_MAT_X,NULL,NULL);
384:
385: /* 25-33. Compute Ritz vectors */
386: DSGetMat(eps->ds,DS_MAT_X,&M);
387: BVSetActiveColumns(Z,ctx->bs,nv);
388: if (ctx->lock) {
389: BVSetActiveColumns(P,0,ctx->bs);
390: }
391: BVMult(P,1.0,0.0,Z,M);
392: BVCopy(P,X);
393: if (ctx->lock) {
394: BVSetActiveColumns(P,nconv,ctx->bs);
395: }
396: BVSetActiveColumns(Z,0,ctx->bs);
397: BVMult(X,1.0,1.0,Z,M);
398: if (ctx->lock) {
399: BVSetActiveColumns(X,nconv,ctx->bs);
400: }
401: BVMatMult(X,A,AX);
402: if (B) {
403: BVMatMult(X,B,BX);
404: }
405: DSRestoreMat(eps->ds,DS_MAT_X,&M);
406: }
408: BVDestroy(&Z);
409: BVDestroy(&X);
410: BVDestroy(&R);
411: BVDestroy(&P);
412: BVDestroy(&AX);
413: if (B) {
414: BVDestroy(&BX);
415: }
416: if (nc>0 || eps->nev>ctx->bs) {
417: BVDestroy(&Y);
418: }
419: return(0);
420: }
424: static PetscErrorCode EPSLOBPCGSetBlockSize_LOBPCG(EPS eps,PetscInt bs)
425: {
426: EPS_LOBPCG *ctx = (EPS_LOBPCG*)eps->data;
429: ctx->bs = bs;
430: return(0);
431: }
435: /*@
436: EPSLOBPCGSetBlockSize - Sets the block size of the LOBPCG method.
438: Logically Collective on EPS
440: Input Parameters:
441: + eps - the eigenproblem solver context
442: - bs - the block size
444: Options Database Key:
445: . -eps_lobpcg_blocksize - Sets the block size
447: Level: advanced
449: .seealso: EPSLOBPCGGetBlockSize()
450: @*/
451: PetscErrorCode EPSLOBPCGSetBlockSize(EPS eps,PetscInt bs)
452: {
458: PetscTryMethod(eps,"EPSLOBPCGSetBlockSize_C",(EPS,PetscInt),(eps,bs));
459: return(0);
460: }
464: static PetscErrorCode EPSLOBPCGGetBlockSize_LOBPCG(EPS eps,PetscInt *bs)
465: {
466: EPS_LOBPCG *ctx = (EPS_LOBPCG*)eps->data;
469: *bs = ctx->bs;
470: return(0);
471: }
475: /*@
476: EPSLOBPCGGetBlockSize - Gets the block size used in the LOBPCG method.
478: Not Collective
480: Input Parameter:
481: . eps - the eigenproblem solver context
483: Output Parameter:
484: . bs - the block size
486: Level: advanced
488: .seealso: EPSLOBPCGSetBlockSize()
489: @*/
490: PetscErrorCode EPSLOBPCGGetBlockSize(EPS eps,PetscInt *bs)
491: {
497: PetscTryMethod(eps,"EPSLOBPCGGetBlockSize_C",(EPS,PetscInt*),(eps,bs));
498: return(0);
499: }
503: static PetscErrorCode EPSLOBPCGSetLocking_LOBPCG(EPS eps,PetscBool lock)
504: {
505: EPS_LOBPCG *ctx = (EPS_LOBPCG*)eps->data;
508: ctx->lock = lock;
509: return(0);
510: }
514: /*@
515: EPSLOBPCGSetLocking - Choose between locking and non-locking variants of
516: the LOBPCG method.
518: Logically Collective on EPS
520: Input Parameters:
521: + eps - the eigenproblem solver context
522: - lock - true if the locking variant must be selected
524: Options Database Key:
525: . -eps_lobpcg_locking - Sets the locking flag
527: Notes:
528: This flag refers to soft locking (converged vectors within the current
529: block iterate), since hard locking is always used (when nev is larger
530: than the block size).
532: Level: advanced
534: .seealso: EPSLOBPCGGetLocking()
535: @*/
536: PetscErrorCode EPSLOBPCGSetLocking(EPS eps,PetscBool lock)
537: {
543: PetscTryMethod(eps,"EPSLOBPCGSetLocking_C",(EPS,PetscBool),(eps,lock));
544: return(0);
545: }
549: static PetscErrorCode EPSLOBPCGGetLocking_LOBPCG(EPS eps,PetscBool *lock)
550: {
551: EPS_LOBPCG *ctx = (EPS_LOBPCG*)eps->data;
554: *lock = ctx->lock;
555: return(0);
556: }
560: /*@
561: EPSLOBPCGGetLocking - Gets the locking flag used in the LOBPCG method.
563: Not Collective
565: Input Parameter:
566: . eps - the eigenproblem solver context
568: Output Parameter:
569: . lock - the locking flag
571: Level: advanced
573: .seealso: EPSLOBPCGSetLocking()
574: @*/
575: PetscErrorCode EPSLOBPCGGetLocking(EPS eps,PetscBool *lock)
576: {
582: PetscTryMethod(eps,"EPSLOBPCGGetLocking_C",(EPS,PetscBool*),(eps,lock));
583: return(0);
584: }
588: PetscErrorCode EPSView_LOBPCG(EPS eps,PetscViewer viewer)
589: {
591: EPS_LOBPCG *ctx = (EPS_LOBPCG*)eps->data;
592: PetscBool isascii;
595: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
596: if (isascii) {
597: PetscViewerASCIIPrintf(viewer," LOBPCG: block size %D\n",ctx->bs);
598: PetscViewerASCIIPrintf(viewer," LOBPCG: soft locking %sactivated\n",ctx->lock?"":"de");
599: }
600: return(0);
601: }
605: PetscErrorCode EPSSetFromOptions_LOBPCG(PetscOptions *PetscOptionsObject,EPS eps)
606: {
608: PetscBool lock,flg;
609: PetscInt bs;
610: KSP ksp;
613: PetscOptionsHead(PetscOptionsObject,"EPS LOBPCG Options");
614: PetscOptionsInt("-eps_lobpcg_blocksize","LOBPCG block size","EPSLOBPCGSetBlockSize",20,&bs,&flg);
615: if (flg) {
616: EPSLOBPCGSetBlockSize(eps,bs);
617: }
618: PetscOptionsBool("-eps_lobpcg_locking","Choose between locking and non-locking variants","EPSLOBPCGSetLocking",PETSC_TRUE,&lock,&flg);
619: if (flg) {
620: EPSLOBPCGSetLocking(eps,lock);
621: }
623: /* Set STPrecond as the default ST */
624: if (!((PetscObject)eps->st)->type_name) {
625: STSetType(eps->st,STPRECOND);
626: }
628: /* Set the default options of the KSP */
629: STGetKSP(eps->st,&ksp);
630: if (!((PetscObject)ksp)->type_name) {
631: KSPSetType(ksp,KSPPREONLY);
632: }
633: PetscOptionsTail();
634: return(0);
635: }
639: PetscErrorCode EPSDestroy_LOBPCG(EPS eps)
640: {
644: PetscFree(eps->data);
645: PetscObjectComposeFunction((PetscObject)eps,"EPSLOBPCGSetBlockSize_C",NULL);
646: PetscObjectComposeFunction((PetscObject)eps,"EPSLOBPCGGetBlockSize_C",NULL);
647: PetscObjectComposeFunction((PetscObject)eps,"EPSLOBPCGSetLocking_C",NULL);
648: PetscObjectComposeFunction((PetscObject)eps,"EPSLOBPCGGetLocking_C",NULL);
649: return(0);
650: }
654: PETSC_EXTERN PetscErrorCode EPSCreate_LOBPCG(EPS eps)
655: {
656: EPS_LOBPCG *lobpcg;
660: PetscNewLog(eps,&lobpcg);
661: eps->data = (void*)lobpcg;
662: lobpcg->lock = PETSC_TRUE;
664: eps->ops->setup = EPSSetUp_LOBPCG;
665: eps->ops->solve = EPSSolve_LOBPCG;
666: eps->ops->setfromoptions = EPSSetFromOptions_LOBPCG;
667: eps->ops->destroy = EPSDestroy_LOBPCG;
668: eps->ops->view = EPSView_LOBPCG;
669: eps->ops->backtransform = EPSBackTransform_Default;
670: STSetType(eps->st,STPRECOND);
671: STPrecondSetKSPHasMat(eps->st,PETSC_TRUE);
672: PetscObjectComposeFunction((PetscObject)eps,"EPSLOBPCGSetBlockSize_C",EPSLOBPCGSetBlockSize_LOBPCG);
673: PetscObjectComposeFunction((PetscObject)eps,"EPSLOBPCGGetBlockSize_C",EPSLOBPCGGetBlockSize_LOBPCG);
674: PetscObjectComposeFunction((PetscObject)eps,"EPSLOBPCGSetLocking_C",EPSLOBPCGSetLocking_LOBPCG);
675: PetscObjectComposeFunction((PetscObject)eps,"EPSLOBPCGGetLocking_C",EPSLOBPCGGetLocking_LOBPCG);
676: return(0);
677: }