Actual source code: stfunc.c

slepc-3.15.0 2021-03-31
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2021, 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:    The ST interface routines, callable by users
 12: */

 14: #include <slepc/private/stimpl.h>

 16: PetscClassId     ST_CLASSID = 0;
 17: PetscLogEvent    ST_SetUp = 0,ST_ComputeOperator = 0,ST_Apply = 0,ST_ApplyTranspose = 0,ST_MatSetUp = 0,ST_MatMult = 0,ST_MatMultTranspose = 0,ST_MatSolve = 0,ST_MatSolveTranspose = 0;
 18: static PetscBool STPackageInitialized = PETSC_FALSE;

 20: const char *STMatModes[] = {"COPY","INPLACE","SHELL","STMatMode","ST_MATMODE_",0};

 22: /*@C
 23:    STFinalizePackage - This function destroys everything in the Slepc interface
 24:    to the ST package. It is called from SlepcFinalize().

 26:    Level: developer

 28: .seealso: SlepcFinalize()
 29: @*/
 30: PetscErrorCode STFinalizePackage(void)
 31: {

 35:   PetscFunctionListDestroy(&STList);
 36:   STPackageInitialized = PETSC_FALSE;
 37:   STRegisterAllCalled  = PETSC_FALSE;
 38:   return(0);
 39: }

 41: /*@C
 42:    STInitializePackage - This function initializes everything in the ST package.
 43:    It is called from PetscDLLibraryRegister() when using dynamic libraries, and
 44:    on the first call to STCreate() when using static libraries.

 46:    Level: developer

 48: .seealso: SlepcInitialize()
 49: @*/
 50: PetscErrorCode STInitializePackage(void)
 51: {
 52:   char           logList[256];
 53:   PetscBool      opt,pkg;
 54:   PetscClassId   classids[1];

 58:   if (STPackageInitialized) return(0);
 59:   STPackageInitialized = PETSC_TRUE;
 60:   /* Register Classes */
 61:   PetscClassIdRegister("Spectral Transform",&ST_CLASSID);
 62:   /* Register Constructors */
 63:   STRegisterAll();
 64:   /* Register Events */
 65:   PetscLogEventRegister("STSetUp",ST_CLASSID,&ST_SetUp);
 66:   PetscLogEventRegister("STComputeOperatr",ST_CLASSID,&ST_ComputeOperator);
 67:   PetscLogEventRegister("STApply",ST_CLASSID,&ST_Apply);
 68:   PetscLogEventRegister("STApplyTranspose",ST_CLASSID,&ST_ApplyTranspose);
 69:   PetscLogEventRegister("STMatSetUp",ST_CLASSID,&ST_MatSetUp);
 70:   PetscLogEventRegister("STMatMult",ST_CLASSID,&ST_MatMult);
 71:   PetscLogEventRegister("STMatMultTranspose",ST_CLASSID,&ST_MatMultTranspose);
 72:   PetscLogEventRegister("STMatSolve",ST_CLASSID,&ST_MatSolve);
 73:   PetscLogEventRegister("STMatSolveTranspose",ST_CLASSID,&ST_MatSolveTranspose);
 74:   /* Process Info */
 75:   classids[0] = ST_CLASSID;
 76:   PetscInfoProcessClass("st",1,&classids[0]);
 77:   /* Process summary exclusions */
 78:   PetscOptionsGetString(NULL,NULL,"-log_exclude",logList,sizeof(logList),&opt);
 79:   if (opt) {
 80:     PetscStrInList("st",logList,',',&pkg);
 81:     if (pkg) { PetscLogEventDeactivateClass(ST_CLASSID); }
 82:   }
 83:   /* Register package finalizer */
 84:   PetscRegisterFinalize(STFinalizePackage);
 85:   return(0);
 86: }

 88: /*@
 89:    STReset - Resets the ST context to the initial state (prior to setup)
 90:    and destroys any allocated Vecs and Mats.

 92:    Collective on st

 94:    Input Parameter:
 95: .  st - the spectral transformation context

 97:    Level: advanced

 99: .seealso: STDestroy()
100: @*/
101: PetscErrorCode STReset(ST st)
102: {

107:   if (!st) return(0);
108:   STCheckNotSeized(st,1);
109:   if (st->ops->reset) { (*st->ops->reset)(st); }
110:   if (st->ksp) { KSPReset(st->ksp); }
111:   MatDestroyMatrices(PetscMax(2,st->nmat),&st->T);
112:   MatDestroyMatrices(PetscMax(2,st->nmat),&st->A);
113:   st->nmat = 0;
114:   PetscFree(st->Astate);
115:   MatDestroy(&st->Op);
116:   MatDestroy(&st->P);
117:   MatDestroy(&st->Pmat);
118:   VecDestroyVecs(st->nwork,&st->work);
119:   st->nwork = 0;
120:   VecDestroy(&st->wb);
121:   VecDestroy(&st->wht);
122:   VecDestroy(&st->D);
123:   st->state   = ST_STATE_INITIAL;
124:   st->opready = PETSC_FALSE;
125:   return(0);
126: }

128: /*@C
129:    STDestroy - Destroys ST context that was created with STCreate().

131:    Collective on st

133:    Input Parameter:
134: .  st - the spectral transformation context

136:    Level: beginner

138: .seealso: STCreate(), STSetUp()
139: @*/
140: PetscErrorCode STDestroy(ST *st)
141: {

145:   if (!*st) return(0);
147:   if (--((PetscObject)(*st))->refct > 0) { *st = 0; return(0); }
148:   STReset(*st);
149:   if ((*st)->ops->destroy) { (*(*st)->ops->destroy)(*st); }
150:   KSPDestroy(&(*st)->ksp);
151:   PetscHeaderDestroy(st);
152:   return(0);
153: }

155: /*@
156:    STCreate - Creates a spectral transformation context.

158:    Collective

160:    Input Parameter:
161: .  comm - MPI communicator

163:    Output Parameter:
164: .  st - location to put the spectral transformation context

166:    Level: beginner

168: .seealso: STSetUp(), STApply(), STDestroy(), ST
169: @*/
170: PetscErrorCode STCreate(MPI_Comm comm,ST *newst)
171: {
173:   ST             st;

177:   *newst = 0;
178:   STInitializePackage();
179:   SlepcHeaderCreate(st,ST_CLASSID,"ST","Spectral Transformation","ST",comm,STDestroy,STView);

181:   st->A            = NULL;
182:   st->nmat         = 0;
183:   st->sigma        = 0.0;
184:   st->defsigma     = 0.0;
185:   st->matmode      = ST_MATMODE_COPY;
186:   st->str          = UNKNOWN_NONZERO_PATTERN;
187:   st->transform    = PETSC_FALSE;
188:   st->D            = NULL;
189:   st->Pmat         = NULL;

191:   st->ksp          = NULL;
192:   st->usesksp      = PETSC_FALSE;
193:   st->nwork        = 0;
194:   st->work         = NULL;
195:   st->wb           = NULL;
196:   st->wht          = NULL;
197:   st->state        = ST_STATE_INITIAL;
198:   st->Astate       = NULL;
199:   st->T            = NULL;
200:   st->Op           = NULL;
201:   st->opseized     = PETSC_FALSE;
202:   st->opready      = PETSC_FALSE;
203:   st->P            = NULL;
204:   st->M            = NULL;
205:   st->sigma_set    = PETSC_FALSE;
206:   st->asymm        = PETSC_FALSE;
207:   st->aherm        = PETSC_FALSE;
208:   st->data         = NULL;

210:   *newst = st;
211:   return(0);
212: }

214: /*
215:    Checks whether the ST matrices are all symmetric or hermitian.
216: */
217: PETSC_STATIC_INLINE PetscErrorCode STMatIsSymmetricKnown(ST st,PetscBool *symm,PetscBool *herm)
218: {
220:   PetscInt       i;
221:   PetscBool      sbaij=PETSC_FALSE,set,flg=PETSC_FALSE;

224:   /* check if problem matrices are all sbaij */
225:   for (i=0;i<st->nmat;i++) {
226:     PetscObjectTypeCompareAny((PetscObject)st->A[i],&sbaij,MATSEQSBAIJ,MATMPISBAIJ,"");
227:     if (!sbaij) break;
228:   }
229:   /* check if user has set the symmetric flag */
230:   *symm = PETSC_TRUE;
231:   for (i=0;i<st->nmat;i++) {
232:     MatIsSymmetricKnown(st->A[i],&set,&flg);
233:     if (!set || !flg) { *symm = PETSC_FALSE; break; }
234:   }
235:   if (sbaij) *symm = PETSC_TRUE;
236: #if defined(PETSC_USE_COMPLEX)
237:   /* check if user has set the hermitian flag */
238:   *herm = PETSC_TRUE;
239:   for (i=0;i<st->nmat;i++) {
240:     MatIsHermitianKnown(st->A[i],&set,&flg);
241:     if (!set || !flg) { *herm = PETSC_FALSE; break; }
242:   }
243: #else
244:   *herm = *symm;
245: #endif
246:   return(0);
247: }

249: /*@
250:    STSetMatrices - Sets the matrices associated with the eigenvalue problem.

252:    Collective on st

254:    Input Parameters:
255: +  st - the spectral transformation context
256: .  n  - number of matrices in array A
257: -  A  - the array of matrices associated with the eigensystem

259:    Notes:
260:    It must be called before STSetUp(). If it is called again after STSetUp() then
261:    the ST object is reset.

263:    Level: intermediate

265: .seealso: STGetMatrix(), STGetNumMatrices(), STSetUp(), STReset()
266:  @*/
267: PetscErrorCode STSetMatrices(ST st,PetscInt n,Mat A[])
268: {
269:   PetscInt       i;
271:   PetscBool      same=PETSC_TRUE;

276:   if (n <= 0) SETERRQ1(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"Must have one or more matrices, you have %D",n);
279:   STCheckNotSeized(st,1);

281:   if (st->state) {
282:     if (n!=st->nmat) same = PETSC_FALSE;
283:     for (i=0;same&&i<n;i++) {
284:       if (A[i]!=st->A[i]) same = PETSC_FALSE;
285:     }
286:     if (!same) { STReset(st); }
287:   } else same = PETSC_FALSE;
288:   if (!same) {
289:     MatDestroyMatrices(PetscMax(2,st->nmat),&st->A);
290:     PetscCalloc1(PetscMax(2,n),&st->A);
291:     PetscLogObjectMemory((PetscObject)st,PetscMax(2,n)*sizeof(Mat));
292:     PetscFree(st->Astate);
293:     PetscMalloc1(PetscMax(2,n),&st->Astate);
294:     PetscLogObjectMemory((PetscObject)st,PetscMax(2,n)*sizeof(PetscInt));
295:   }
296:   for (i=0;i<n;i++) {
298:     PetscObjectReference((PetscObject)A[i]);
299:     MatDestroy(&st->A[i]);
300:     st->A[i] = A[i];
301:     st->Astate[i] = ((PetscObject)A[i])->state;
302:   }
303:   if (n==1) {
304:     st->A[1] = NULL;
305:     st->Astate[1] = 0;
306:   }
307:   st->nmat = n;
308:   if (same) st->state = ST_STATE_UPDATED;
309:   else st->state = ST_STATE_INITIAL;
310:   st->opready = PETSC_FALSE;
311:   if (!same) {
312:     STMatIsSymmetricKnown(st,&st->asymm,&st->aherm);
313:   }
314:   return(0);
315: }

317: /*@
318:    STGetMatrix - Gets the matrices associated with the original eigensystem.

320:    Not collective, though parallel Mats are returned if the ST is parallel

322:    Input Parameter:
323: +  st - the spectral transformation context
324: -  k  - the index of the requested matrix (starting in 0)

326:    Output Parameters:
327: .  A - the requested matrix

329:    Level: intermediate

331: .seealso: STSetMatrices(), STGetNumMatrices()
332: @*/
333: PetscErrorCode STGetMatrix(ST st,PetscInt k,Mat *A)
334: {
339:   STCheckMatrices(st,1);
340:   if (k<0 || k>=st->nmat) SETERRQ1(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %D",st->nmat-1);
341:   if (((PetscObject)st->A[k])->state!=st->Astate[k]) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"Cannot retrieve original matrices (have been modified)");
342:   *A = st->A[k];
343:   return(0);
344: }

346: /*@
347:    STGetMatrixTransformed - Gets the matrices associated with the transformed eigensystem.

349:    Not collective, though parallel Mats are returned if the ST is parallel

351:    Input Parameter:
352: +  st - the spectral transformation context
353: -  k  - the index of the requested matrix (starting in 0)

355:    Output Parameters:
356: .  T - the requested matrix

358:    Level: developer

360: .seealso: STGetMatrix(), STGetNumMatrices()
361: @*/
362: PetscErrorCode STGetMatrixTransformed(ST st,PetscInt k,Mat *T)
363: {
368:   STCheckMatrices(st,1);
369:   if (k<0 || k>=st->nmat) SETERRQ1(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %D",st->nmat-1);
370:   if (!st->T) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_POINTER,"There are no transformed matrices");
371:   *T = st->T[k];
372:   return(0);
373: }

375: /*@
376:    STGetNumMatrices - Returns the number of matrices stored in the ST.

378:    Not collective

380:    Input Parameter:
381: .  st - the spectral transformation context

383:    Output Parameters:
384: .  n - the number of matrices passed in STSetMatrices()

386:    Level: intermediate

388: .seealso: STSetMatrices()
389: @*/
390: PetscErrorCode STGetNumMatrices(ST st,PetscInt *n)
391: {
395:   *n = st->nmat;
396:   return(0);
397: }

399: /*@
400:    STResetMatrixState - Resets the stored state of the matrices in the ST.

402:    Logically Collective on st

404:    Input Parameter:
405: .  st - the spectral transformation context

407:    Note:
408:    This is useful in solvers where the user matrices are modified during
409:    the computation, as in nonlinear inverse iteration. The effect is that
410:    STGetMatrix() will retrieve the modified matrices as if they were
411:    the matrices originally provided by the user.

413:    Level: developer

415: .seealso: STGetMatrix(), EPSPowerSetNonlinear()
416: @*/
417: PetscErrorCode STResetMatrixState(ST st)
418: {
419:   PetscInt i;

423:   for (i=0;i<st->nmat;i++) st->Astate[i] = ((PetscObject)st->A[i])->state;
424:   return(0);
425: }

427: /*@
428:    STSetPreconditionerMat - Sets the matrix to be used to build the preconditioner.

430:    Collective on st

432:    Input Parameter:
433: +  st  - the spectral transformation context
434: -  mat - the matrix that will be used in constructing the preconditioner

436:    Notes:
437:    This matrix will be passed to the internal KSP object (via the last argument
438:    of KSPSetOperators()) as the matrix to be used when constructing the preconditioner.
439:    If no matrix is set or mat is set to NULL, A-sigma*B will be used
440:    to build the preconditioner, being sigma the value set by STSetShift().

442:    More precisely, this is relevant for spectral transformations that represent
443:    a rational matrix function, and use a KSP object for the denominator, called
444:    K in the description of STGetOperator(). It includes also the STPRECOND case.
445:    If the user has a good approximation to matrix K that can be used to build a
446:    cheap preconditioner, it can be passed with this function. Note that it affects
447:    only the Pmat argument of KSPSetOperators(), not the Amat argument.

449:    If a preconditioner matrix is set, the default is to use an iterative KSP
450:    rather than a direct method.

452:    Level: advanced

454: .seealso: STGetPreconditionerMat(), STSetShift(), STGetOperator()
455: @*/
456: PetscErrorCode STSetPreconditionerMat(ST st,Mat mat)
457: {

464:   STCheckNotSeized(st,1);
465:   PetscObjectReference((PetscObject)mat);
466:   MatDestroy(&st->Pmat);
467:   st->Pmat     = mat;
468:   st->state    = ST_STATE_INITIAL;
469:   st->opready  = PETSC_FALSE;
470:   return(0);
471: }

473: /*@
474:    STGetPreconditionerMat - Returns the matrix previously set by STSetPreconditionerMat().

476:    Collective on st

478:    Input Parameter:
479: .  st - the spectral transformation context

481:    Output Parameter:
482: .  mat - the matrix that will be used in constructing the preconditioner or
483:    NULL if no matrix was set by STSetPreconditionerMat().

485:    Level: advanced

487: .seealso: STSetPreconditionerMat()
488: @*/
489: PetscErrorCode STGetPreconditionerMat(ST st,Mat *mat)
490: {
494:   *mat = st->Pmat;
495:   return(0);
496: }

498: /*@
499:    STSetShift - Sets the shift associated with the spectral transformation.

501:    Logically Collective on st

503:    Input Parameters:
504: +  st - the spectral transformation context
505: -  shift - the value of the shift

507:    Notes:
508:    In some spectral transformations, changing the shift may have associated
509:    a lot of work, for example recomputing a factorization.

511:    This function is normally not directly called by users, since the shift is
512:    indirectly set by EPSSetTarget().

514:    Level: intermediate

516: .seealso: EPSSetTarget(), STGetShift(), STSetDefaultShift()
517: @*/
518: PetscErrorCode STSetShift(ST st,PetscScalar shift)
519: {

526:   if (st->sigma != shift) {
527:     STCheckNotSeized(st,1);
528:     if (st->state==ST_STATE_SETUP && st->ops->setshift) {
529:       (*st->ops->setshift)(st,shift);
530:     }
531:     st->sigma = shift;
532:   }
533:   st->sigma_set = PETSC_TRUE;
534:   return(0);
535: }

537: /*@
538:    STGetShift - Gets the shift associated with the spectral transformation.

540:    Not Collective

542:    Input Parameter:
543: .  st - the spectral transformation context

545:    Output Parameter:
546: .  shift - the value of the shift

548:    Level: intermediate

550: .seealso: STSetShift()
551: @*/
552: PetscErrorCode STGetShift(ST st,PetscScalar* shift)
553: {
557:   *shift = st->sigma;
558:   return(0);
559: }

561: /*@
562:    STSetDefaultShift - Sets the value of the shift that should be employed if
563:    the user did not specify one.

565:    Logically Collective on st

567:    Input Parameters:
568: +  st - the spectral transformation context
569: -  defaultshift - the default value of the shift

571:    Level: developer

573: .seealso: STSetShift()
574: @*/
575: PetscErrorCode STSetDefaultShift(ST st,PetscScalar defaultshift)
576: {
580:   if (st->defsigma != defaultshift) {
581:     st->defsigma = defaultshift;
582:     st->state    = ST_STATE_INITIAL;
583:     st->opready  = PETSC_FALSE;
584:   }
585:   return(0);
586: }

588: /*@
589:    STScaleShift - Multiply the shift with a given factor.

591:    Logically Collective on st

593:    Input Parameters:
594: +  st     - the spectral transformation context
595: -  factor - the scaling factor

597:    Note:
598:    This function does not update the transformation matrices, as opposed to
599:    STSetShift().

601:    Level: developer

603: .seealso: STSetShift()
604: @*/
605: PetscErrorCode STScaleShift(ST st,PetscScalar factor)
606: {
610:   st->sigma *= factor;
611:   return(0);
612: }

614: /*@
615:    STSetBalanceMatrix - Sets the diagonal matrix to be used for balancing.

617:    Collective on st

619:    Input Parameters:
620: +  st - the spectral transformation context
621: -  D  - the diagonal matrix (represented as a vector)

623:    Notes:
624:    If this matrix is set, STApply will effectively apply D*OP*D^{-1}. Use NULL
625:    to reset a previously passed D.

627:    Balancing is usually set via EPSSetBalance, but the advanced user may use
628:    this function to bypass the usual balancing methods.

630:    Level: developer

632: .seealso: EPSSetBalance(), STApply(), STGetBalanceMatrix()
633: @*/
634: PetscErrorCode STSetBalanceMatrix(ST st,Vec D)
635: {

640:   if (st->D == D) return(0);
641:   STCheckNotSeized(st,1);
642:   if (D) {
645:     PetscObjectReference((PetscObject)D);
646:   }
647:   VecDestroy(&st->D);
648:   st->D = D;
649:   st->state   = ST_STATE_INITIAL;
650:   st->opready = PETSC_FALSE;
651:   return(0);
652: }

654: /*@
655:    STGetBalanceMatrix - Gets the balance matrix used by the spectral transformation.

657:    Not collective, but vector is shared by all processors that share the ST

659:    Input Parameter:
660: .  st - the spectral transformation context

662:    Output Parameter:
663: .  D  - the diagonal matrix (represented as a vector)

665:    Note:
666:    If the matrix was not set, a null pointer will be returned.

668:    Level: developer

670: .seealso: STSetBalanceMatrix()
671: @*/
672: PetscErrorCode STGetBalanceMatrix(ST st,Vec *D)
673: {
677:   *D = st->D;
678:   return(0);
679: }

681: /*@C
682:    STMatCreateVecs - Get vector(s) compatible with the ST matrices.

684:    Collective on st

686:    Input Parameter:
687: .  st - the spectral transformation context

689:    Output Parameters:
690: +  right - (optional) vector that the matrix can be multiplied against
691: -  left  - (optional) vector that the matrix vector product can be stored in

693:    Level: developer
694: @*/
695: PetscErrorCode STMatCreateVecs(ST st,Vec *right,Vec *left)
696: {

700:   STCheckMatrices(st,1);
701:   MatCreateVecs(st->A[0],right,left);
702:   return(0);
703: }

705: /*@C
706:    STMatCreateVecsEmpty - Get vector(s) compatible with the ST matrices, i.e. with the same
707:    parallel layout, but without internal array.

709:    Collective on st

711:    Input Parameter:
712: .  st - the spectral transformation context

714:    Output Parameters:
715: +  right - (optional) vector that the matrix can be multiplied against
716: -  left  - (optional) vector that the matrix vector product can be stored in

718:    Level: developer

720: .seealso: MatCreateVecsEmpty()
721: @*/
722: PetscErrorCode STMatCreateVecsEmpty(ST st,Vec *right,Vec *left)
723: {

727:   STCheckMatrices(st,1);
728:   MatCreateVecsEmpty(st->A[0],right,left);
729:   return(0);
730: }

732: /*@
733:    STMatGetSize - Returns the number of rows and columns of the ST matrices.

735:    Not Collective

737:    Input Parameter:
738: .  st - the spectral transformation context

740:    Output Parameters:
741: +  m - the number of global rows
742: -  n - the number of global columns

744:    Level: developer
745: @*/
746: PetscErrorCode STMatGetSize(ST st,PetscInt *m,PetscInt *n)
747: {

751:   STCheckMatrices(st,1);
752:   MatGetSize(st->A[0],m,n);
753:   return(0);
754: }

756: /*@
757:    STMatGetLocalSize - Returns the number of local rows and columns of the ST matrices.

759:    Not Collective

761:    Input Parameter:
762: .  st - the spectral transformation context

764:    Output Parameters:
765: +  m - the number of local rows
766: -  n - the number of local columns

768:    Level: developer
769: @*/
770: PetscErrorCode STMatGetLocalSize(ST st,PetscInt *m,PetscInt *n)
771: {

775:   STCheckMatrices(st,1);
776:   MatGetLocalSize(st->A[0],m,n);
777:   return(0);
778: }

780: /*@C
781:    STSetOptionsPrefix - Sets the prefix used for searching for all
782:    ST options in the database.

784:    Logically Collective on st

786:    Input Parameters:
787: +  st     - the spectral transformation context
788: -  prefix - the prefix string to prepend to all ST option requests

790:    Notes:
791:    A hyphen (-) must NOT be given at the beginning of the prefix name.
792:    The first character of all runtime options is AUTOMATICALLY the
793:    hyphen.

795:    Level: advanced

797: .seealso: STAppendOptionsPrefix(), STGetOptionsPrefix()
798: @*/
799: PetscErrorCode STSetOptionsPrefix(ST st,const char *prefix)
800: {

805:   if (!st->ksp) { STGetKSP(st,&st->ksp); }
806:   KSPSetOptionsPrefix(st->ksp,prefix);
807:   KSPAppendOptionsPrefix(st->ksp,"st_");
808:   PetscObjectSetOptionsPrefix((PetscObject)st,prefix);
809:   return(0);
810: }

812: /*@C
813:    STAppendOptionsPrefix - Appends to the prefix used for searching for all
814:    ST options in the database.

816:    Logically Collective on st

818:    Input Parameters:
819: +  st     - the spectral transformation context
820: -  prefix - the prefix string to prepend to all ST option requests

822:    Notes:
823:    A hyphen (-) must NOT be given at the beginning of the prefix name.
824:    The first character of all runtime options is AUTOMATICALLY the
825:    hyphen.

827:    Level: advanced

829: .seealso: STSetOptionsPrefix(), STGetOptionsPrefix()
830: @*/
831: PetscErrorCode STAppendOptionsPrefix(ST st,const char *prefix)
832: {

837:   PetscObjectAppendOptionsPrefix((PetscObject)st,prefix);
838:   if (!st->ksp) { STGetKSP(st,&st->ksp); }
839:   KSPSetOptionsPrefix(st->ksp,((PetscObject)st)->prefix);
840:   KSPAppendOptionsPrefix(st->ksp,"st_");
841:   return(0);
842: }

844: /*@C
845:    STGetOptionsPrefix - Gets the prefix used for searching for all
846:    ST options in the database.

848:    Not Collective

850:    Input Parameters:
851: .  st - the spectral transformation context

853:    Output Parameters:
854: .  prefix - pointer to the prefix string used, is returned

856:    Note:
857:    On the Fortran side, the user should pass in a string 'prefix' of
858:    sufficient length to hold the prefix.

860:    Level: advanced

862: .seealso: STSetOptionsPrefix(), STAppendOptionsPrefix()
863: @*/
864: PetscErrorCode STGetOptionsPrefix(ST st,const char *prefix[])
865: {

871:   PetscObjectGetOptionsPrefix((PetscObject)st,prefix);
872:   return(0);
873: }

875: /*@C
876:    STView - Prints the ST data structure.

878:    Collective on st

880:    Input Parameters:
881: +  st - the ST context
882: -  viewer - optional visualization context

884:    Note:
885:    The available visualization contexts include
886: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
887: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
888:          output where only the first processor opens
889:          the file.  All other processors send their
890:          data to the first processor to print.

892:    The user can open an alternative visualization contexts with
893:    PetscViewerASCIIOpen() (output to a specified file).

895:    Level: beginner

897: .seealso: EPSView()
898: @*/
899: PetscErrorCode STView(ST st,PetscViewer viewer)
900: {
902:   STType         cstr;
903:   char           str[50];
904:   PetscBool      isascii,isstring;

908:   if (!viewer) {
909:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)st),&viewer);
910:   }

914:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
915:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
916:   if (isascii) {
917:     PetscObjectPrintClassNamePrefixType((PetscObject)st,viewer);
918:     if (st->ops->view) {
919:       PetscViewerASCIIPushTab(viewer);
920:       (*st->ops->view)(st,viewer);
921:       PetscViewerASCIIPopTab(viewer);
922:     }
923:     SlepcSNPrintfScalar(str,sizeof(str),st->sigma,PETSC_FALSE);
924:     PetscViewerASCIIPrintf(viewer,"  shift: %s\n",str);
925:     PetscViewerASCIIPrintf(viewer,"  number of matrices: %D\n",st->nmat);
926:     switch (st->matmode) {
927:     case ST_MATMODE_COPY:
928:       break;
929:     case ST_MATMODE_INPLACE:
930:       PetscViewerASCIIPrintf(viewer,"  shifting the matrix and unshifting at exit\n");
931:       break;
932:     case ST_MATMODE_SHELL:
933:       PetscViewerASCIIPrintf(viewer,"  using a shell matrix\n");
934:       break;
935:     }
936:     if (st->nmat>1 && st->matmode != ST_MATMODE_SHELL) {
937:       PetscViewerASCIIPrintf(viewer,"  all matrices have %s\n",MatStructures[st->str]);
938:     }
939:     if (st->transform && st->nmat>2) {
940:       PetscViewerASCIIPrintf(viewer,"  computing transformed matrices\n");
941:     }
942:   } else if (isstring) {
943:     STGetType(st,&cstr);
944:     PetscViewerStringSPrintf(viewer," %-7.7s",cstr);
945:     if (st->ops->view) { (*st->ops->view)(st,viewer); }
946:   }
947:   if (st->usesksp) {
948:     if (!st->ksp) { STGetKSP(st,&st->ksp); }
949:     PetscViewerASCIIPushTab(viewer);
950:     KSPView(st->ksp,viewer);
951:     PetscViewerASCIIPopTab(viewer);
952:   }
953:   return(0);
954: }

956: /*@C
957:    STViewFromOptions - View from options

959:    Collective on ST

961:    Input Parameters:
962: +  st   - the spectral transformation context
963: .  obj  - optional object
964: -  name - command line option

966:    Level: intermediate

968: .seealso: STView(), STCreate()
969: @*/
970: PetscErrorCode STViewFromOptions(ST st,PetscObject obj,const char name[])
971: {

976:   PetscObjectViewFromOptions((PetscObject)st,obj,name);
977:   return(0);
978: }

980: /*@C
981:    STRegister - Adds a method to the spectral transformation package.

983:    Not collective

985:    Input Parameters:
986: +  name - name of a new user-defined transformation
987: -  function - routine to create method context

989:    Notes:
990:    STRegister() may be called multiple times to add several user-defined
991:    spectral transformations.

993:    Sample usage:
994: .vb
995:     STRegister("my_transform",MyTransformCreate);
996: .ve

998:    Then, your spectral transform can be chosen with the procedural interface via
999: $     STSetType(st,"my_transform")
1000:    or at runtime via the option
1001: $     -st_type my_transform

1003:    Level: advanced

1005: .seealso: STRegisterAll()
1006: @*/
1007: PetscErrorCode STRegister(const char *name,PetscErrorCode (*function)(ST))
1008: {

1012:   STInitializePackage();
1013:   PetscFunctionListAdd(&STList,name,function);
1014:   return(0);
1015: }