1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, 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: ST interface routines, callable by users
12: */
14: #include <slepc/private/stimpl.h> 16: PetscErrorCode STApply_Generic(ST st,Vec x,Vec y) 17: {
18: if (st->M && st->P) {
19: MatMult(st->M,x,st->work[0]);
20: STMatSolve(st,st->work[0],y);
21: } else if (st->M) MatMult(st->M,x,y);
22: else STMatSolve(st,x,y);
23: PetscFunctionReturn(0);
24: }
26: /*@
27: STApply - Applies the spectral transformation operator to a vector, for
28: instance (A - sB)^-1 B in the case of the shift-and-invert transformation
29: and generalized eigenproblem.
31: Collective on st
33: Input Parameters:
34: + st - the spectral transformation context
35: - x - input vector
37: Output Parameter:
38: . y - output vector
40: Level: developer
42: .seealso: STApplyTranspose(), STApplyHermitianTranspose()
43: @*/
44: PetscErrorCode STApply(ST st,Vec x,Vec y) 45: {
46: Mat Op;
52: STCheckMatrices(st,1);
54: VecSetErrorIfLocked(y,3);
56: STGetOperator_Private(st,&Op);
57: MatMult(Op,x,y);
58: PetscFunctionReturn(0);
59: }
61: PetscErrorCode STApplyMat_Generic(ST st,Mat B,Mat C) 62: {
63: Mat work;
65: if (st->M && st->P) {
66: MatMatMult(st->M,B,MAT_INITIAL_MATRIX,PETSC_DEFAULT,&work);
67: STMatMatSolve(st,work,C);
68: MatDestroy(&work);
69: } else if (st->M) MatMatMult(st->M,B,MAT_REUSE_MATRIX,PETSC_DEFAULT,&C);
70: else STMatMatSolve(st,B,C);
71: PetscFunctionReturn(0);
72: }
74: /*@
75: STApplyMat - Applies the spectral transformation operator to a matrix, for
76: instance (A - sB)^-1 B in the case of the shift-and-invert transformation
77: and generalized eigenproblem.
79: Collective on st
81: Input Parameters:
82: + st - the spectral transformation context
83: - X - input matrix
85: Output Parameter:
86: . Y - output matrix
88: Level: developer
90: .seealso: STApply()
91: @*/
92: PetscErrorCode STApplyMat(ST st,Mat X,Mat Y) 93: {
98: STCheckMatrices(st,1);
101: (*st->ops->applymat)(st,X,Y);
102: PetscFunctionReturn(0);
103: }
105: PetscErrorCode STApplyTranspose_Generic(ST st,Vec x,Vec y)106: {
107: if (st->M && st->P) {
108: STMatSolveTranspose(st,x,st->work[0]);
109: MatMultTranspose(st->M,st->work[0],y);
110: } else if (st->M) MatMultTranspose(st->M,x,y);
111: else STMatSolveTranspose(st,x,y);
112: PetscFunctionReturn(0);
113: }
115: /*@
116: STApplyTranspose - Applies the transpose of the operator to a vector, for
117: instance B^T(A - sB)^-T in the case of the shift-and-invert transformation
118: and generalized eigenproblem.
120: Collective on st
122: Input Parameters:
123: + st - the spectral transformation context
124: - x - input vector
126: Output Parameter:
127: . y - output vector
129: Level: developer
131: .seealso: STApply(), STApplyHermitianTranspose()
132: @*/
133: PetscErrorCode STApplyTranspose(ST st,Vec x,Vec y)134: {
135: Mat Op;
141: STCheckMatrices(st,1);
143: VecSetErrorIfLocked(y,3);
145: STGetOperator_Private(st,&Op);
146: MatMultTranspose(Op,x,y);
147: PetscFunctionReturn(0);
148: }
150: /*@
151: STApplyHermitianTranspose - Applies the hermitian-transpose of the operator
152: to a vector, for instance B^H(A - sB)^-H in the case of the shift-and-invert
153: transformation and generalized eigenproblem.
155: Collective on st
157: Input Parameters:
158: + st - the spectral transformation context
159: - x - input vector
161: Output Parameter:
162: . y - output vector
164: Note:
165: Currently implemented via STApplyTranspose() with appropriate conjugation.
167: Level: developer
169: .seealso: STApply(), STApplyTranspose()
170: @*/
171: PetscErrorCode STApplyHermitianTranspose(ST st,Vec x,Vec y)172: {
173: Mat Op;
179: STCheckMatrices(st,1);
181: VecSetErrorIfLocked(y,3);
183: STGetOperator_Private(st,&Op);
184: MatMultHermitianTranspose(Op,x,y);
185: PetscFunctionReturn(0);
186: }
188: /*@
189: STGetBilinearForm - Returns the matrix used in the bilinear form with a
190: generalized problem with semi-definite B.
192: Not collective, though a parallel Mat may be returned
194: Input Parameters:
195: . st - the spectral transformation context
197: Output Parameter:
198: . B - output matrix
200: Notes:
201: The output matrix B must be destroyed after use. It will be NULL in
202: case of standard eigenproblems.
204: Level: developer
206: .seealso: BVSetMatrix()
207: @*/
208: PetscErrorCode STGetBilinearForm(ST st,Mat *B)209: {
213: STCheckMatrices(st,1);
214: (*st->ops->getbilinearform)(st,B);
215: PetscFunctionReturn(0);
216: }
218: PetscErrorCode STGetBilinearForm_Default(ST st,Mat *B)219: {
220: if (st->nmat==1) *B = NULL;
221: else {
222: *B = st->A[1];
223: PetscObjectReference((PetscObject)*B);
224: }
225: PetscFunctionReturn(0);
226: }
228: static PetscErrorCode MatMult_STOperator(Mat Op,Vec x,Vec y)229: {
230: ST st;
232: MatShellGetContext(Op,&st);
233: STSetUp(st);
234: PetscLogEventBegin(ST_Apply,st,x,y,0);
235: if (st->D) { /* with balancing */
236: VecPointwiseDivide(st->wb,x,st->D);
237: (*st->ops->apply)(st,st->wb,y);
238: VecPointwiseMult(y,y,st->D);
239: } else (*st->ops->apply)(st,x,y);
240: PetscLogEventEnd(ST_Apply,st,x,y,0);
241: PetscFunctionReturn(0);
242: }
244: static PetscErrorCode MatMultTranspose_STOperator(Mat Op,Vec x,Vec y)245: {
246: ST st;
248: MatShellGetContext(Op,&st);
249: STSetUp(st);
250: PetscLogEventBegin(ST_ApplyTranspose,st,x,y,0);
251: if (st->D) { /* with balancing */
252: VecPointwiseMult(st->wb,x,st->D);
253: (*st->ops->applytrans)(st,st->wb,y);
254: VecPointwiseDivide(y,y,st->D);
255: } else (*st->ops->applytrans)(st,x,y);
256: PetscLogEventEnd(ST_ApplyTranspose,st,x,y,0);
257: PetscFunctionReturn(0);
258: }
260: #if defined(PETSC_USE_COMPLEX)
261: static PetscErrorCode MatMultHermitianTranspose_STOperator(Mat Op,Vec x,Vec y)262: {
263: ST st;
265: MatShellGetContext(Op,&st);
266: STSetUp(st);
267: PetscLogEventBegin(ST_ApplyTranspose,st,x,y,0);
268: if (!st->wht) {
269: MatCreateVecs(st->A[0],&st->wht,NULL);
270: PetscLogObjectParent((PetscObject)st,(PetscObject)st->wht);
271: }
272: VecCopy(x,st->wht);
273: VecConjugate(st->wht);
274: if (st->D) { /* with balancing */
275: VecPointwiseMult(st->wb,st->wht,st->D);
276: (*st->ops->applytrans)(st,st->wb,y);
277: VecPointwiseDivide(y,y,st->D);
278: } else (*st->ops->applytrans)(st,st->wht,y);
279: VecConjugate(y);
280: PetscLogEventEnd(ST_ApplyTranspose,st,x,y,0);
281: PetscFunctionReturn(0);
282: }
283: #endif
285: static PetscErrorCode MatMatMult_STOperator(Mat Op,Mat B,Mat C,void *ctx)286: {
287: ST st;
289: MatShellGetContext(Op,&st);
290: STSetUp(st);
291: PetscLogEventBegin(ST_Apply,st,B,C,0);
292: STApplyMat_Generic(st,B,C);
293: PetscLogEventEnd(ST_Apply,st,B,C,0);
294: PetscFunctionReturn(0);
295: }
297: PetscErrorCode STGetOperator_Private(ST st,Mat *Op)298: {
299: PetscInt m,n,M,N;
300: Vec v;
301: VecType vtype;
303: if (!st->Op) {
304: if (Op) *Op = NULL;
305: /* create the shell matrix */
306: MatGetLocalSize(st->A[0],&m,&n);
307: MatGetSize(st->A[0],&M,&N);
308: MatCreateShell(PetscObjectComm((PetscObject)st),m,n,M,N,st,&st->Op);
309: MatShellSetOperation(st->Op,MATOP_MULT,(void(*)(void))MatMult_STOperator);
310: MatShellSetOperation(st->Op,MATOP_MULT_TRANSPOSE,(void(*)(void))MatMultTranspose_STOperator);
311: #if defined(PETSC_USE_COMPLEX)
312: MatShellSetOperation(st->Op,MATOP_MULT_HERMITIAN_TRANSPOSE,(void(*)(void))MatMultHermitianTranspose_STOperator);
313: #else
314: MatShellSetOperation(st->Op,MATOP_MULT_HERMITIAN_TRANSPOSE,(void(*)(void))MatMultTranspose_STOperator);
315: #endif
316: if (!st->D && st->ops->apply==STApply_Generic) {
317: MatShellSetMatProductOperation(st->Op,MATPRODUCT_AB,NULL,MatMatMult_STOperator,NULL,MATDENSE,MATDENSE);
318: MatShellSetMatProductOperation(st->Op,MATPRODUCT_AB,NULL,MatMatMult_STOperator,NULL,MATDENSECUDA,MATDENSECUDA);
319: }
320: /* make sure the shell matrix generates a vector of the same type as the problem matrices */
321: MatCreateVecs(st->A[0],&v,NULL);
322: VecGetType(v,&vtype);
323: MatShellSetVecType(st->Op,vtype);
324: VecDestroy(&v);
325: /* build the operator matrices */
326: STComputeOperator(st);
327: }
328: if (Op) *Op = st->Op;
329: PetscFunctionReturn(0);
330: }
332: /*@
333: STGetOperator - Returns a shell matrix that represents the operator of the
334: spectral transformation.
336: Collective on st
338: Input Parameter:
339: . st - the spectral transformation context
341: Output Parameter:
342: . Op - operator matrix
344: Notes:
345: The operator is defined in linear eigenproblems only, not in polynomial ones,
346: so the call will fail if more than 2 matrices were passed in STSetMatrices().
348: The returned shell matrix is essentially a wrapper to the STApply() and
349: STApplyTranspose() operations. The operator can often be expressed as
351: $ Op = D*inv(K)*M*inv(D)
353: where D is the balancing matrix, and M and K are two matrices corresponding
354: to the numerator and denominator for spectral transformations that represent
355: a rational matrix function. In the case of STSHELL, the inner part inv(K)*M
356: is replaced by the user-provided operation from STShellSetApply().
358: The preconditioner matrix K typically depends on the value of the shift, and
359: its inverse is handled via an internal KSP object. Normal usage does not
360: require explicitly calling STGetOperator(), but it can be used to force the
361: creation of K and M, and then K is passed to the KSP. This is useful for
362: setting options associated with the PCFactor (to set MUMPS options, for instance).
364: The returned matrix must NOT be destroyed by the user. Instead, when no
365: longer needed it must be returned with STRestoreOperator(). In particular,
366: this is required before modifying the ST matrices or the shift.
368: A NULL pointer can be passed in Op in case the matrix is not required but we
369: want to force its creation. In this case, STRestoreOperator() should not be
370: called.
372: Level: advanced
374: .seealso: STApply(), STApplyTranspose(), STSetBalanceMatrix(), STShellSetApply(),
375: STGetKSP(), STSetShift(), STRestoreOperator(), STSetMatrices()
376: @*/
377: PetscErrorCode STGetOperator(ST st,Mat *Op)378: {
381: STCheckMatrices(st,1);
382: STCheckNotSeized(st,1);
384: STGetOperator_Private(st,Op);
385: if (Op) st->opseized = PETSC_TRUE;
386: PetscFunctionReturn(0);
387: }
389: /*@
390: STRestoreOperator - Restore the previously seized operator matrix.
392: Collective on st
394: Input Parameters:
395: + st - the spectral transformation context
396: - Op - operator matrix
398: Notes:
399: The arguments must match the corresponding call to STGetOperator().
401: Level: advanced
403: .seealso: STGetOperator()
404: @*/
405: PetscErrorCode STRestoreOperator(ST st,Mat *Op)406: {
411: *Op = NULL;
412: st->opseized = PETSC_FALSE;
413: PetscFunctionReturn(0);
414: }
416: /*
417: STComputeOperator - Computes the matrices that constitute the operator
419: Op = D*inv(K)*M*inv(D).
421: K and M are computed here (D is user-provided) from the system matrices
422: and the shift sigma (whenever these are changed, this function recomputes
423: K and M). This is used only in linear eigenproblems (nmat<3).
425: K is the "preconditioner matrix": it is the denominator in rational operators,
426: e.g. (A-sigma*B) in shift-and-invert. In non-rational transformations such
427: as STFILTER, K=NULL which means identity. After computing K, it is passed to
428: the internal KSP object via KSPSetOperators.
430: M is the numerator in rational operators. If unused it is set to NULL (e.g.
431: in STPRECOND).
433: STSHELL does not compute anything here, but sets the flag as if it was ready.
434: */
435: PetscErrorCode STComputeOperator(ST st)436: {
437: PC pc;
441: if (!st->opready && st->ops->computeoperator) {
442: PetscInfo(st,"Building the operator matrices\n");
443: STCheckMatrices(st,1);
444: if (!st->T) {
445: PetscCalloc1(PetscMax(2,st->nmat),&st->T);
446: PetscLogObjectMemory((PetscObject)st,PetscMax(2,st->nmat)*sizeof(Mat));
447: }
448: PetscLogEventBegin(ST_ComputeOperator,st,0,0,0);
449: (*st->ops->computeoperator)(st);
450: PetscLogEventEnd(ST_ComputeOperator,st,0,0,0);
451: if (st->usesksp) {
452: if (!st->ksp) STGetKSP(st,&st->ksp);
453: if (st->P) {
454: STSetDefaultKSP(st);
455: ST_KSPSetOperators(st,st->P,st->Pmat?st->Pmat:st->P);
456: } else {
457: /* STPRECOND defaults to PCNONE if st->P is empty */
458: KSPGetPC(st->ksp,&pc);
459: PCSetType(pc,PCNONE);
460: }
461: }
462: }
463: st->opready = PETSC_TRUE;
464: PetscFunctionReturn(0);
465: }
467: /*@
468: STSetUp - Prepares for the use of a spectral transformation.
470: Collective on st
472: Input Parameter:
473: . st - the spectral transformation context
475: Level: advanced
477: .seealso: STCreate(), STApply(), STDestroy()
478: @*/
479: PetscErrorCode STSetUp(ST st)480: {
481: PetscInt i,n,k;
485: STCheckMatrices(st,1);
486: switch (st->state) {
487: case ST_STATE_INITIAL:
488: PetscInfo(st,"Setting up new ST\n");
489: if (!((PetscObject)st)->type_name) STSetType(st,STSHIFT);
490: break;
491: case ST_STATE_SETUP:
492: PetscFunctionReturn(0);
493: case ST_STATE_UPDATED:
494: PetscInfo(st,"Setting up updated ST\n");
495: break;
496: }
497: PetscLogEventBegin(ST_SetUp,st,0,0,0);
498: if (st->state!=ST_STATE_UPDATED) {
499: if (!(st->nmat<3 && st->opready)) {
500: if (st->T) {
501: for (i=0;i<PetscMax(2,st->nmat);i++) MatDestroy(&st->T[i]);
502: }
503: MatDestroy(&st->P);
504: }
505: }
506: if (st->D) {
507: MatGetLocalSize(st->A[0],NULL,&n);
508: VecGetLocalSize(st->D,&k);
510: if (!st->wb) {
511: VecDuplicate(st->D,&st->wb);
512: PetscLogObjectParent((PetscObject)st,(PetscObject)st->wb);
513: }
514: }
515: if (st->nmat<3 && st->transform) STComputeOperator(st);
516: else {
517: if (!st->T) {
518: PetscCalloc1(PetscMax(2,st->nmat),&st->T);
519: PetscLogObjectMemory((PetscObject)st,PetscMax(2,st->nmat)*sizeof(Mat));
520: }
521: }
522: if (st->ops->setup) (*st->ops->setup)(st);
523: st->state = ST_STATE_SETUP;
524: PetscLogEventEnd(ST_SetUp,st,0,0,0);
525: PetscFunctionReturn(0);
526: }
528: /*
529: Computes coefficients for the transformed polynomial,
530: and stores the result in argument S.
532: alpha - value of the parameter of the transformed polynomial
533: beta - value of the previous shift (only used in inplace mode)
534: k - index of first matrix included in the computation
535: coeffs - coefficients of the expansion
536: initial - true if this is the first time
537: precond - whether the preconditioner matrix must be computed
538: */
539: PetscErrorCode STMatMAXPY_Private(ST st,PetscScalar alpha,PetscScalar beta,PetscInt k,PetscScalar *coeffs,PetscBool initial,PetscBool precond,Mat *S)540: {
541: PetscInt *matIdx=NULL,nmat,i,ini=-1;
542: PetscScalar t=1.0,ta,gamma;
543: PetscBool nz=PETSC_FALSE;
544: Mat *A=precond?st->Psplit:st->A;
545: MatStructure str=precond?st->strp:st->str;
547: nmat = st->nmat-k;
548: switch (st->matmode) {
549: case ST_MATMODE_INPLACE:
552: if (initial) {
553: PetscObjectReference((PetscObject)A[0]);
554: *S = A[0];
555: gamma = alpha;
556: } else gamma = alpha-beta;
557: if (gamma != 0.0) {
558: if (st->nmat>1) MatAXPY(*S,gamma,A[1],str);
559: else MatShift(*S,gamma);
560: }
561: break;
562: case ST_MATMODE_SHELL:
564: if (initial) {
565: if (st->nmat>2) {
566: PetscMalloc1(nmat,&matIdx);
567: for (i=0;i<nmat;i++) matIdx[i] = k+i;
568: }
569: STMatShellCreate(st,alpha,nmat,matIdx,coeffs,S);
570: PetscLogObjectParent((PetscObject)st,(PetscObject)*S);
571: if (st->nmat>2) PetscFree(matIdx);
572: } else STMatShellShift(*S,alpha);
573: break;
574: case ST_MATMODE_COPY:
575: if (coeffs) {
576: for (i=0;i<nmat && ini==-1;i++) {
577: if (coeffs[i]!=0.0) ini = i;
578: else t *= alpha;
579: }
580: if (coeffs[ini] != 1.0) nz = PETSC_TRUE;
581: for (i=ini+1;i<nmat&&!nz;i++) if (coeffs[i]!=0.0) nz = PETSC_TRUE;
582: } else { nz = PETSC_TRUE; ini = 0; }
583: if ((alpha == 0.0 || !nz) && t==1.0) {
584: PetscObjectReference((PetscObject)A[k+ini]);
585: MatDestroy(S);
586: *S = A[k+ini];
587: } else {
588: if (*S && *S!=A[k+ini]) {
589: MatSetOption(*S,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
590: MatCopy(A[k+ini],*S,DIFFERENT_NONZERO_PATTERN);
591: } else {
592: MatDestroy(S);
593: MatDuplicate(A[k+ini],MAT_COPY_VALUES,S);
594: MatSetOption(*S,MAT_NEW_NONZERO_ALLOCATION_ERR,PETSC_FALSE);
595: PetscLogObjectParent((PetscObject)st,(PetscObject)*S);
596: }
597: if (coeffs && coeffs[ini]!=1.0) MatScale(*S,coeffs[ini]);
598: for (i=ini+k+1;i<PetscMax(2,st->nmat);i++) {
599: t *= alpha;
600: ta = t;
601: if (coeffs) ta *= coeffs[i-k];
602: if (ta!=0.0) {
603: if (st->nmat>1) MatAXPY(*S,ta,A[i],str);
604: else MatShift(*S,ta);
605: }
606: }
607: }
608: }
609: MatSetOption(*S,MAT_SYMMETRIC,st->asymm);
610: MatSetOption(*S,MAT_HERMITIAN,(PetscImaginaryPart(st->sigma)==0.0)?st->aherm:PETSC_FALSE);
611: PetscFunctionReturn(0);
612: }
614: /*
615: Computes the values of the coefficients required by STMatMAXPY_Private
616: for the case of monomial basis.
617: */
618: PetscErrorCode STCoeffs_Monomial(ST st, PetscScalar *coeffs)619: {
620: PetscInt k,i,ini,inip;
622: /* Compute binomial coefficients */
623: ini = (st->nmat*(st->nmat-1))/2;
624: for (i=0;i<st->nmat;i++) coeffs[ini+i]=1.0;
625: for (k=st->nmat-1;k>=1;k--) {
626: inip = ini+1;
627: ini = (k*(k-1))/2;
628: coeffs[ini] = 1.0;
629: for (i=1;i<k;i++) coeffs[ini+i] = coeffs[ini+i-1]+coeffs[inip+i-1];
630: }
631: PetscFunctionReturn(0);
632: }
634: /*@
635: STPostSolve - Optional post-solve phase, intended for any actions that must
636: be performed on the ST object after the eigensolver has finished.
638: Collective on st
640: Input Parameters:
641: . st - the spectral transformation context
643: Level: developer
645: .seealso: EPSSolve()
646: @*/
647: PetscErrorCode STPostSolve(ST st)648: {
651: if (st->ops->postsolve) (*st->ops->postsolve)(st);
652: PetscFunctionReturn(0);
653: }
655: /*@
656: STBackTransform - Back-transformation phase, intended for
657: spectral transformations which require to transform the computed
658: eigenvalues back to the original eigenvalue problem.
660: Not Collective
662: Input Parameters:
663: + st - the spectral transformation context
664: . n - number of eigenvalues
665: . eigr - real part of a computed eigenvalues
666: - eigi - imaginary part of a computed eigenvalues
668: Level: developer
670: .seealso: STIsInjective()
671: @*/
672: PetscErrorCode STBackTransform(ST st,PetscInt n,PetscScalar* eigr,PetscScalar* eigi)673: {
676: if (st->ops->backtransform) (*st->ops->backtransform)(st,n,eigr,eigi);
677: PetscFunctionReturn(0);
678: }
680: /*@
681: STIsInjective - Ask if this spectral transformation is injective or not
682: (that is, if it corresponds to a one-to-one mapping). If not, then it
683: does not make sense to call STBackTransform().
685: Not collective
687: Input Parameter:
688: . st - the spectral transformation context
690: Output Parameter:
691: . is - the answer
693: Level: developer
695: .seealso: STBackTransform()
696: @*/
697: PetscErrorCode STIsInjective(ST st,PetscBool* is)698: {
699: PetscBool shell;
705: PetscObjectTypeCompare((PetscObject)st,STSHELL,&shell);
706: if (shell) STIsInjective_Shell(st,is);
707: else *is = st->ops->backtransform? PETSC_TRUE: PETSC_FALSE;
708: PetscFunctionReturn(0);
709: }
711: /*@
712: STMatSetUp - Build the preconditioner matrix used in STMatSolve().
714: Collective on st
716: Input Parameters:
717: + st - the spectral transformation context
718: . sigma - the shift
719: - coeffs - the coefficients (may be NULL)
721: Note:
722: This function is not intended to be called by end users, but by SLEPc
723: solvers that use ST. It builds matrix st->P as follows, then calls KSPSetUp().
724: .vb
725: If (coeffs) st->P = Sum_{i=0..nmat-1} coeffs[i]*sigma^i*A_i
726: else st->P = Sum_{i=0..nmat-1} sigma^i*A_i
727: .ve
729: Level: developer
731: .seealso: STMatSolve()
732: @*/
733: PetscErrorCode STMatSetUp(ST st,PetscScalar sigma,PetscScalar *coeffs)734: {
737: STCheckMatrices(st,1);
739: PetscLogEventBegin(ST_MatSetUp,st,0,0,0);
740: STMatMAXPY_Private(st,sigma,0.0,0,coeffs,PETSC_TRUE,PETSC_FALSE,&st->P);
741: if (st->Psplit) STMatMAXPY_Private(st,sigma,0.0,0,coeffs,PETSC_TRUE,PETSC_TRUE,&st->Pmat);
742: ST_KSPSetOperators(st,st->P,st->Pmat?st->Pmat:st->P);
743: KSPSetUp(st->ksp);
744: PetscLogEventEnd(ST_MatSetUp,st,0,0,0);
745: PetscFunctionReturn(0);
746: }
748: /*@
749: STSetWorkVecs - Sets a number of work vectors into the ST object.
751: Collective on st
753: Input Parameters:
754: + st - the spectral transformation context
755: - nw - number of work vectors to allocate
757: Developer Notes:
758: This is SLEPC_EXTERN because it may be required by shell STs.
760: Level: developer
762: .seealso: STMatCreateVecs()
763: @*/
764: PetscErrorCode STSetWorkVecs(ST st,PetscInt nw)765: {
766: PetscInt i;
771: if (st->nwork < nw) {
772: VecDestroyVecs(st->nwork,&st->work);
773: st->nwork = nw;
774: PetscMalloc1(nw,&st->work);
775: for (i=0;i<nw;i++) STMatCreateVecs(st,&st->work[i],NULL);
776: PetscLogObjectParents(st,nw,st->work);
777: }
778: PetscFunctionReturn(0);
779: }