Actual source code: stfunc.c

slepc-3.12.2 2020-01-13
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2019, 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>            /*I "slepcst.h" I*/

 16: PetscClassId     ST_CLASSID = 0;
 17: PetscLogEvent    ST_SetUp = 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;

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

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

 93:    Collective on st

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

 98:    Level: advanced

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

108:   if (!st) return(0);
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:   PetscFree(st->Astate);
114:   MatDestroy(&st->P);
115:   VecDestroyVecs(st->nwork,&st->work);
116:   st->nwork = 0;
117:   VecDestroy(&st->wb);
118:   VecDestroy(&st->D);
119:   st->state = ST_STATE_INITIAL;
120:   return(0);
121: }

123: /*@
124:    STDestroy - Destroys ST context that was created with STCreate().

126:    Collective on st

128:    Input Parameter:
129: .  st - the spectral transformation context

131:    Level: beginner

133: .seealso: STCreate(), STSetUp()
134: @*/
135: PetscErrorCode STDestroy(ST *st)
136: {

140:   if (!*st) return(0);
142:   if (--((PetscObject)(*st))->refct > 0) { *st = 0; return(0); }
143:   STReset(*st);
144:   if ((*st)->ops->destroy) { (*(*st)->ops->destroy)(*st); }
145:   KSPDestroy(&(*st)->ksp);
146:   PetscHeaderDestroy(st);
147:   return(0);
148: }

150: /*@
151:    STCreate - Creates a spectral transformation context.

153:    Collective

155:    Input Parameter:
156: .  comm - MPI communicator

158:    Output Parameter:
159: .  st - location to put the spectral transformation context

161:    Level: beginner

163: .seealso: STSetUp(), STApply(), STDestroy(), ST
164: @*/
165: PetscErrorCode STCreate(MPI_Comm comm,ST *newst)
166: {
168:   ST             st;

172:   *newst = 0;
173:   STInitializePackage();
174:   SlepcHeaderCreate(st,ST_CLASSID,"ST","Spectral Transformation","ST",comm,STDestroy,STView);

176:   st->A            = NULL;
177:   st->nmat         = 0;
178:   st->sigma        = 0.0;
179:   st->defsigma     = 0.0;
180:   st->matmode      = ST_MATMODE_COPY;
181:   st->str          = DIFFERENT_NONZERO_PATTERN;
182:   st->transform    = PETSC_FALSE;
183:   st->D            = NULL;

185:   st->ksp          = NULL;
186:   st->usesksp      = PETSC_FALSE;
187:   st->nwork        = 0;
188:   st->work         = NULL;
189:   st->wb           = NULL;
190:   st->state        = ST_STATE_INITIAL;
191:   st->Astate       = NULL;
192:   st->T            = NULL;
193:   st->P            = NULL;
194:   st->sigma_set    = PETSC_FALSE;
195:   st->data         = NULL;

197:   *newst = st;
198:   return(0);
199: }

201: /*@
202:    STSetMatrices - Sets the matrices associated with the eigenvalue problem.

204:    Collective on st

206:    Input Parameters:
207: +  st - the spectral transformation context
208: .  n  - number of matrices in array A
209: -  A  - the array of matrices associated with the eigensystem

211:    Notes:
212:    It must be called before STSetUp(). If it is called again after STSetUp() then
213:    the ST object is reset.

215:    Level: intermediate

217: .seealso: STGetMatrix(), STGetNumMatrices(), STSetUp(), STReset()
218:  @*/
219: PetscErrorCode STSetMatrices(ST st,PetscInt n,Mat A[])
220: {
221:   PetscInt       i;
223:   PetscBool      same=PETSC_TRUE;

228:   if (n <= 0) SETERRQ1(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"Must have one or more matrices, you have %D",n);
231:   if (st->state) {
232:     if (n!=st->nmat) same = PETSC_FALSE;
233:     for (i=0;same&&i<n;i++) {
234:       if (A[i]!=st->A[i]) same = PETSC_FALSE;
235:     }
236:     if (!same) { STReset(st); }
237:   } else same = PETSC_FALSE;
238:   if (!same) {
239:     MatDestroyMatrices(PetscMax(2,st->nmat),&st->A);
240:     PetscCalloc1(PetscMax(2,n),&st->A);
241:     PetscLogObjectMemory((PetscObject)st,PetscMax(2,n)*sizeof(Mat));
242:     PetscFree(st->Astate);
243:     PetscMalloc1(PetscMax(2,n),&st->Astate);
244:     PetscLogObjectMemory((PetscObject)st,PetscMax(2,n)*sizeof(PetscInt));
245:   }
246:   for (i=0;i<n;i++) {
248:     PetscObjectReference((PetscObject)A[i]);
249:     MatDestroy(&st->A[i]);
250:     st->A[i] = A[i];
251:     st->Astate[i] = ((PetscObject)A[i])->state;
252:   }
253:   if (n==1) {
254:     st->A[1] = NULL;
255:     st->Astate[1] = 0;
256:   }
257:   st->nmat = n;
258:   if (same) st->state = ST_STATE_UPDATED;
259:   else st->state = ST_STATE_INITIAL;
260:   return(0);
261: }

263: /*@
264:    STGetMatrix - Gets the matrices associated with the original eigensystem.

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

268:    Input Parameter:
269: +  st - the spectral transformation context
270: -  k  - the index of the requested matrix (starting in 0)

272:    Output Parameters:
273: .  A - the requested matrix

275:    Level: intermediate

277: .seealso: STSetMatrices(), STGetNumMatrices()
278: @*/
279: PetscErrorCode STGetMatrix(ST st,PetscInt k,Mat *A)
280: {
285:   STCheckMatrices(st,1);
286:   if (k<0 || k>=st->nmat) SETERRQ1(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %D",st->nmat-1);
287:   if (((PetscObject)st->A[k])->state!=st->Astate[k]) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"Cannot retrieve original matrices (have been modified)");
288:   *A = st->A[k];
289:   return(0);
290: }

292: /*@
293:    STGetMatrixTransformed - Gets the matrices associated with the transformed eigensystem.

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

297:    Input Parameter:
298: +  st - the spectral transformation context
299: -  k  - the index of the requested matrix (starting in 0)

301:    Output Parameters:
302: .  T - the requested matrix

304:    Level: developer

306: .seealso: STGetMatrix(), STGetNumMatrices()
307: @*/
308: PetscErrorCode STGetMatrixTransformed(ST st,PetscInt k,Mat *T)
309: {
314:   STCheckMatrices(st,1);
315:   if (k<0 || k>=st->nmat) SETERRQ1(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %D",st->nmat-1);
316:   if (!st->T) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_POINTER,"There are no transformed matrices");
317:   *T = st->T[k];
318:   return(0);
319: }

321: /*@
322:    STGetNumMatrices - Returns the number of matrices stored in the ST.

324:    Not collective

326:    Input Parameter:
327: .  st - the spectral transformation context

329:    Output Parameters:
330: .  n - the number of matrices passed in STSetMatrices()

332:    Level: intermediate

334: .seealso: STSetMatrices()
335: @*/
336: PetscErrorCode STGetNumMatrices(ST st,PetscInt *n)
337: {
341:   *n = st->nmat;
342:   return(0);
343: }

345: /*@
346:    STResetMatrixState - Resets the stored state of the matrices in the ST.

348:    Logically Collective on st

350:    Input Parameter:
351: .  st - the spectral transformation context

353:    Note:
354:    This is useful in solvers where the user matrices are modified during
355:    the computation, as in nonlinear inverse iteration. The effect is that
356:    STGetMatrix() will retrieve the modified matrices as if they were
357:    the matrices originally provided by the user.

359:    Level: developer

361: .seealso: STGetMatrix(), EPSPowerSetNonlinear()
362: @*/
363: PetscErrorCode STResetMatrixState(ST st)
364: {
365:   PetscInt i;

369:   for (i=0;i<st->nmat;i++) st->Astate[i] = ((PetscObject)st->A[i])->state;
370:   return(0);
371: }

373: /*@
374:    STSetShift - Sets the shift associated with the spectral transformation.

376:    Logically Collective on st

378:    Input Parameters:
379: +  st - the spectral transformation context
380: -  shift - the value of the shift

382:    Notes:
383:    In some spectral transformations, changing the shift may have associated
384:    a lot of work, for example recomputing a factorization.

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

389:    Level: intermediate

391: .seealso: EPSSetTarget(), STGetShift(), STSetDefaultShift()
392: @*/
393: PetscErrorCode STSetShift(ST st,PetscScalar shift)
394: {

401:   if (st->state==ST_STATE_SETUP && st->sigma != shift) {
402:     if (st->ops->setshift) {
403:       (*st->ops->setshift)(st,shift);
404:     }
405:   }
406:   st->sigma = shift;
407:   st->sigma_set = PETSC_TRUE;
408:   return(0);
409: }

411: /*@
412:    STGetShift - Gets the shift associated with the spectral transformation.

414:    Not Collective

416:    Input Parameter:
417: .  st - the spectral transformation context

419:    Output Parameter:
420: .  shift - the value of the shift

422:    Level: intermediate

424: .seealso: STSetShift()
425: @*/
426: PetscErrorCode STGetShift(ST st,PetscScalar* shift)
427: {
431:   *shift = st->sigma;
432:   return(0);
433: }

435: /*@
436:    STSetDefaultShift - Sets the value of the shift that should be employed if
437:    the user did not specify one.

439:    Logically Collective on st

441:    Input Parameters:
442: +  st - the spectral transformation context
443: -  defaultshift - the default value of the shift

445:    Level: developer

447: .seealso: STSetShift()
448: @*/
449: PetscErrorCode STSetDefaultShift(ST st,PetscScalar defaultshift)
450: {
454:   st->defsigma = defaultshift;
455:   return(0);
456: }

458: /*@
459:    STScaleShift - Multiply the shift with a given factor.

461:    Logically Collective on st

463:    Input Parameters:
464: +  st     - the spectral transformation context
465: -  factor - the scaling factor

467:    Note:
468:    This function does not update the transformation matrices, as opposed to
469:    STSetShift().

471:    Level: developer

473: .seealso: STSetShift()
474: @*/
475: PetscErrorCode STScaleShift(ST st,PetscScalar factor)
476: {
480:   st->sigma *= factor;
481:   return(0);
482: }

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

487:    Collective on st

489:    Input Parameters:
490: +  st - the spectral transformation context
491: -  D  - the diagonal matrix (represented as a vector)

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

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

500:    Level: developer

502: .seealso: EPSSetBalance(), STApply(), STGetBalanceMatrix()
503: @*/
504: PetscErrorCode STSetBalanceMatrix(ST st,Vec D)
505: {

510:   if (D) {
513:     PetscObjectReference((PetscObject)D);
514:   }
515:   VecDestroy(&st->D);
516:   st->D = D;
517:   st->state = ST_STATE_INITIAL;
518:   return(0);
519: }

521: /*@
522:    STGetBalanceMatrix - Gets the balance matrix used by the spectral transformation.

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

526:    Input Parameter:
527: .  st - the spectral transformation context

529:    Output Parameter:
530: .  D  - the diagonal matrix (represented as a vector)

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

535:    Level: developer

537: .seealso: STSetBalanceMatrix()
538: @*/
539: PetscErrorCode STGetBalanceMatrix(ST st,Vec *D)
540: {
544:   *D = st->D;
545:   return(0);
546: }

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

551:    Collective on st

553:    Input Parameter:
554: .  st - the spectral transformation context

556:    Output Parameters:
557: +  right - (optional) vector that the matrix can be multiplied against
558: -  left  - (optional) vector that the matrix vector product can be stored in

560:    Level: developer
561: @*/
562: PetscErrorCode STMatCreateVecs(ST st,Vec *right,Vec *left)
563: {

567:   STCheckMatrices(st,1);
568:   MatCreateVecs(st->A[0],right,left);
569:   return(0);
570: }

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

576:    Collective on st

578:    Input Parameter:
579: .  st - the spectral transformation context

581:    Output Parameters:
582: +  right - (optional) vector that the matrix can be multiplied against
583: -  left  - (optional) vector that the matrix vector product can be stored in

585:    Level: developer

587: .seealso: MatCreateVecsEmpty()
588: @*/
589: PetscErrorCode STMatCreateVecsEmpty(ST st,Vec *right,Vec *left)
590: {

594:   STCheckMatrices(st,1);
595:   MatCreateVecsEmpty(st->A[0],right,left);
596:   return(0);
597: }

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

602:    Not Collective

604:    Input Parameter:
605: .  st - the spectral transformation context

607:    Output Parameters:
608: +  m - the number of global rows
609: -  n - the number of global columns

611:    Level: developer
612: @*/
613: PetscErrorCode STMatGetSize(ST st,PetscInt *m,PetscInt *n)
614: {

618:   STCheckMatrices(st,1);
619:   MatGetSize(st->A[0],m,n);
620:   return(0);
621: }

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

626:    Not Collective

628:    Input Parameter:
629: .  st - the spectral transformation context

631:    Output Parameters:
632: +  m - the number of local rows
633: -  n - the number of local columns

635:    Level: developer
636: @*/
637: PetscErrorCode STMatGetLocalSize(ST st,PetscInt *m,PetscInt *n)
638: {

642:   STCheckMatrices(st,1);
643:   MatGetLocalSize(st->A[0],m,n);
644:   return(0);
645: }

647: /*@C
648:    STSetOptionsPrefix - Sets the prefix used for searching for all
649:    ST options in the database.

651:    Logically Collective on st

653:    Input Parameters:
654: +  st     - the spectral transformation context
655: -  prefix - the prefix string to prepend to all ST option requests

657:    Notes:
658:    A hyphen (-) must NOT be given at the beginning of the prefix name.
659:    The first character of all runtime options is AUTOMATICALLY the
660:    hyphen.

662:    Level: advanced

664: .seealso: STAppendOptionsPrefix(), STGetOptionsPrefix()
665: @*/
666: PetscErrorCode STSetOptionsPrefix(ST st,const char *prefix)
667: {

672:   if (!st->ksp) { STGetKSP(st,&st->ksp); }
673:   KSPSetOptionsPrefix(st->ksp,prefix);
674:   KSPAppendOptionsPrefix(st->ksp,"st_");
675:   PetscObjectSetOptionsPrefix((PetscObject)st,prefix);
676:   return(0);
677: }

679: /*@C
680:    STAppendOptionsPrefix - Appends to the prefix used for searching for all
681:    ST options in the database.

683:    Logically Collective on st

685:    Input Parameters:
686: +  st     - the spectral transformation context
687: -  prefix - the prefix string to prepend to all ST option requests

689:    Notes:
690:    A hyphen (-) must NOT be given at the beginning of the prefix name.
691:    The first character of all runtime options is AUTOMATICALLY the
692:    hyphen.

694:    Level: advanced

696: .seealso: STSetOptionsPrefix(), STGetOptionsPrefix()
697: @*/
698: PetscErrorCode STAppendOptionsPrefix(ST st,const char *prefix)
699: {

704:   PetscObjectAppendOptionsPrefix((PetscObject)st,prefix);
705:   if (!st->ksp) { STGetKSP(st,&st->ksp); }
706:   KSPSetOptionsPrefix(st->ksp,((PetscObject)st)->prefix);
707:   KSPAppendOptionsPrefix(st->ksp,"st_");
708:   return(0);
709: }

711: /*@C
712:    STGetOptionsPrefix - Gets the prefix used for searching for all
713:    ST options in the database.

715:    Not Collective

717:    Input Parameters:
718: .  st - the spectral transformation context

720:    Output Parameters:
721: .  prefix - pointer to the prefix string used, is returned

723:    Note:
724:    On the Fortran side, the user should pass in a string 'prefix' of
725:    sufficient length to hold the prefix.

727:    Level: advanced

729: .seealso: STSetOptionsPrefix(), STAppendOptionsPrefix()
730: @*/
731: PetscErrorCode STGetOptionsPrefix(ST st,const char *prefix[])
732: {

738:   PetscObjectGetOptionsPrefix((PetscObject)st,prefix);
739:   return(0);
740: }

742: /*@C
743:    STView - Prints the ST data structure.

745:    Collective on st

747:    Input Parameters:
748: +  st - the ST context
749: -  viewer - optional visualization context

751:    Note:
752:    The available visualization contexts include
753: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
754: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
755:          output where only the first processor opens
756:          the file.  All other processors send their
757:          data to the first processor to print.

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

762:    Level: beginner

764: .seealso: EPSView(), PetscViewerASCIIOpen()
765: @*/
766: PetscErrorCode STView(ST st,PetscViewer viewer)
767: {
769:   STType         cstr;
770:   const char*    pat=NULL;
771:   char           str[50];
772:   PetscBool      isascii,isstring;

776:   if (!viewer) {
777:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)st),&viewer);
778:   }

782:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
783:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
784:   if (isascii) {
785:     PetscObjectPrintClassNamePrefixType((PetscObject)st,viewer);
786:     if (st->ops->view) {
787:       PetscViewerASCIIPushTab(viewer);
788:       (*st->ops->view)(st,viewer);
789:       PetscViewerASCIIPopTab(viewer);
790:     }
791:     SlepcSNPrintfScalar(str,50,st->sigma,PETSC_FALSE);
792:     PetscViewerASCIIPrintf(viewer,"  shift: %s\n",str);
793:     PetscViewerASCIIPrintf(viewer,"  number of matrices: %D\n",st->nmat);
794:     switch (st->matmode) {
795:     case ST_MATMODE_COPY:
796:       break;
797:     case ST_MATMODE_INPLACE:
798:       PetscViewerASCIIPrintf(viewer,"  shifting the matrix and unshifting at exit\n");
799:       break;
800:     case ST_MATMODE_SHELL:
801:       PetscViewerASCIIPrintf(viewer,"  using a shell matrix\n");
802:       break;
803:     }
804:     if (st->nmat>1 && st->matmode != ST_MATMODE_SHELL) {
805:       switch (st->str) {
806:         case SAME_NONZERO_PATTERN:      pat = "same nonzero pattern";break;
807:         case DIFFERENT_NONZERO_PATTERN: pat = "different nonzero pattern";break;
808:         case SUBSET_NONZERO_PATTERN:    pat = "subset nonzero pattern";break;
809:       }
810:       PetscViewerASCIIPrintf(viewer,"  all matrices have %s\n",pat);
811:     }
812:     if (st->transform && st->nmat>2) {
813:       PetscViewerASCIIPrintf(viewer,"  computing transformed matrices\n");
814:     }
815:   } else if (isstring) {
816:     STGetType(st,&cstr);
817:     PetscViewerStringSPrintf(viewer," %-7.7s",cstr);
818:     if (st->ops->view) { (*st->ops->view)(st,viewer); }
819:   }
820:   if (st->usesksp) {
821:     if (!st->ksp) { STGetKSP(st,&st->ksp); }
822:     PetscViewerASCIIPushTab(viewer);
823:     KSPView(st->ksp,viewer);
824:     PetscViewerASCIIPopTab(viewer);
825:   }
826:   return(0);
827: }

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

832:    Not collective

834:    Input Parameters:
835: +  name - name of a new user-defined transformation
836: -  function - routine to create method context

838:    Notes:
839:    STRegister() may be called multiple times to add several user-defined
840:    spectral transformations.

842:    Sample usage:
843: .vb
844:     STRegister("my_transform",MyTransformCreate);
845: .ve

847:    Then, your spectral transform can be chosen with the procedural interface via
848: $     STSetType(st,"my_transform")
849:    or at runtime via the option
850: $     -st_type my_transform

852:    Level: advanced

854: .seealso: STRegisterAll()
855: @*/
856: PetscErrorCode STRegister(const char *name,PetscErrorCode (*function)(ST))
857: {

861:   STInitializePackage();
862:   PetscFunctionListAdd(&STList,name,function);
863:   return(0);
864: }