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(), ST169: @*/
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 Parameters:
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 Parameters:
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 Parameters:
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 ST659: 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 ST961: 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: }