1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, 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("STComputeOperator",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: PetscFree(st->Astate);
114: MatDestroy(&st->Op);
115: MatDestroy(&st->P);
116: VecDestroyVecs(st->nwork,&st->work);
117: st->nwork = 0;
118: VecDestroy(&st->wb);
119: VecDestroy(&st->wht);
120: VecDestroy(&st->D);
121: st->state = ST_STATE_INITIAL;
122: st->opready = PETSC_FALSE;
123: return(0);
124: }
126: /*@
127: STDestroy - Destroys ST context that was created with STCreate().
129: Collective on st
131: Input Parameter:
132: . st - the spectral transformation context
134: Level: beginner
136: .seealso: STCreate(), STSetUp()
137: @*/
138: PetscErrorCode STDestroy(ST *st)139: {
143: if (!*st) return(0);
145: if (--((PetscObject)(*st))->refct > 0) { *st = 0; return(0); }
146: STReset(*st);
147: if ((*st)->ops->destroy) { (*(*st)->ops->destroy)(*st); }
148: KSPDestroy(&(*st)->ksp);
149: PetscHeaderDestroy(st);
150: return(0);
151: }
153: /*@
154: STCreate - Creates a spectral transformation context.
156: Collective
158: Input Parameter:
159: . comm - MPI communicator
161: Output Parameter:
162: . st - location to put the spectral transformation context
164: Level: beginner
166: .seealso: STSetUp(), STApply(), STDestroy(), ST167: @*/
168: PetscErrorCode STCreate(MPI_Comm comm,ST *newst)169: {
171: ST st;
175: *newst = 0;
176: STInitializePackage();
177: SlepcHeaderCreate(st,ST_CLASSID,"ST","Spectral Transformation","ST",comm,STDestroy,STView);
179: st->A = NULL;
180: st->nmat = 0;
181: st->sigma = 0.0;
182: st->defsigma = 0.0;
183: st->matmode = ST_MATMODE_COPY;
184: st->str = DIFFERENT_NONZERO_PATTERN;
185: st->transform = PETSC_FALSE;
186: st->D = NULL;
188: st->ksp = NULL;
189: st->usesksp = PETSC_FALSE;
190: st->nwork = 0;
191: st->work = NULL;
192: st->wb = NULL;
193: st->wht = NULL;
194: st->state = ST_STATE_INITIAL;
195: st->Astate = NULL;
196: st->T = NULL;
197: st->Op = NULL;
198: st->opseized = PETSC_FALSE;
199: st->opready = PETSC_FALSE;
200: st->P = NULL;
201: st->M = NULL;
202: st->sigma_set = PETSC_FALSE;
203: st->asymm = PETSC_FALSE;
204: st->aherm = PETSC_FALSE;
205: st->data = NULL;
207: *newst = st;
208: return(0);
209: }
211: /*
212: Checks whether the ST matrices are all symmetric or hermitian.
213: */
214: PETSC_STATIC_INLINE PetscErrorCode STMatIsSymmetricKnown(ST st,PetscBool *symm,PetscBool *herm)215: {
217: PetscInt i;
218: PetscBool sbaij=PETSC_FALSE,set,flg=PETSC_FALSE;
221: /* check if problem matrices are all sbaij */
222: for (i=0;i<st->nmat;i++) {
223: PetscObjectTypeCompareAny((PetscObject)st->A[i],&sbaij,MATSEQSBAIJ,MATMPISBAIJ,"");
224: if (!sbaij) break;
225: }
226: /* check if user has set the symmetric flag */
227: *symm = PETSC_TRUE;
228: for (i=0;i<st->nmat;i++) {
229: MatIsSymmetricKnown(st->A[i],&set,&flg);
230: if (!set || !flg) { *symm = PETSC_FALSE; break; }
231: }
232: if (sbaij) *symm = PETSC_TRUE;
233: #if defined(PETSC_USE_COMPLEX)
234: /* check if user has set the hermitian flag */
235: *herm = PETSC_TRUE;
236: for (i=0;i<st->nmat;i++) {
237: MatIsHermitianKnown(st->A[i],&set,&flg);
238: if (!set || !flg) { *herm = PETSC_FALSE; break; }
239: }
240: #else
241: *herm = *symm;
242: #endif
243: return(0);
244: }
246: /*@
247: STSetMatrices - Sets the matrices associated with the eigenvalue problem.
249: Collective on st
251: Input Parameters:
252: + st - the spectral transformation context
253: . n - number of matrices in array A
254: - A - the array of matrices associated with the eigensystem
256: Notes:
257: It must be called before STSetUp(). If it is called again after STSetUp() then
258: the ST object is reset.
260: Level: intermediate
262: .seealso: STGetMatrix(), STGetNumMatrices(), STSetUp(), STReset()
263: @*/
264: PetscErrorCode STSetMatrices(ST st,PetscInt n,Mat A[])265: {
266: PetscInt i;
268: PetscBool same=PETSC_TRUE;
273: if (n <= 0) SETERRQ1(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"Must have one or more matrices, you have %D",n);
276: STCheckNotSeized(st,1);
278: if (st->state) {
279: if (n!=st->nmat) same = PETSC_FALSE;
280: for (i=0;same&&i<n;i++) {
281: if (A[i]!=st->A[i]) same = PETSC_FALSE;
282: }
283: if (!same) { STReset(st); }
284: } else same = PETSC_FALSE;
285: if (!same) {
286: MatDestroyMatrices(PetscMax(2,st->nmat),&st->A);
287: PetscCalloc1(PetscMax(2,n),&st->A);
288: PetscLogObjectMemory((PetscObject)st,PetscMax(2,n)*sizeof(Mat));
289: PetscFree(st->Astate);
290: PetscMalloc1(PetscMax(2,n),&st->Astate);
291: PetscLogObjectMemory((PetscObject)st,PetscMax(2,n)*sizeof(PetscInt));
292: }
293: for (i=0;i<n;i++) {
295: PetscObjectReference((PetscObject)A[i]);
296: MatDestroy(&st->A[i]);
297: st->A[i] = A[i];
298: st->Astate[i] = ((PetscObject)A[i])->state;
299: }
300: if (n==1) {
301: st->A[1] = NULL;
302: st->Astate[1] = 0;
303: }
304: st->nmat = n;
305: if (same) st->state = ST_STATE_UPDATED;
306: else st->state = ST_STATE_INITIAL;
307: st->opready = PETSC_FALSE;
308: if (!same) {
309: STMatIsSymmetricKnown(st,&st->asymm,&st->aherm);
310: }
311: return(0);
312: }
314: /*@
315: STGetMatrix - Gets the matrices associated with the original eigensystem.
317: Not collective, though parallel Mats are returned if the ST is parallel
319: Input Parameter:
320: + st - the spectral transformation context
321: - k - the index of the requested matrix (starting in 0)
323: Output Parameters:
324: . A - the requested matrix
326: Level: intermediate
328: .seealso: STSetMatrices(), STGetNumMatrices()
329: @*/
330: PetscErrorCode STGetMatrix(ST st,PetscInt k,Mat *A)331: {
336: STCheckMatrices(st,1);
337: if (k<0 || k>=st->nmat) SETERRQ1(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %D",st->nmat-1);
338: if (((PetscObject)st->A[k])->state!=st->Astate[k]) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_SUP,"Cannot retrieve original matrices (have been modified)");
339: *A = st->A[k];
340: return(0);
341: }
343: /*@
344: STGetMatrixTransformed - Gets the matrices associated with the transformed eigensystem.
346: Not collective, though parallel Mats are returned if the ST is parallel
348: Input Parameter:
349: + st - the spectral transformation context
350: - k - the index of the requested matrix (starting in 0)
352: Output Parameters:
353: . T - the requested matrix
355: Level: developer
357: .seealso: STGetMatrix(), STGetNumMatrices()
358: @*/
359: PetscErrorCode STGetMatrixTransformed(ST st,PetscInt k,Mat *T)360: {
365: STCheckMatrices(st,1);
366: if (k<0 || k>=st->nmat) SETERRQ1(PetscObjectComm((PetscObject)st),PETSC_ERR_ARG_OUTOFRANGE,"k must be between 0 and %D",st->nmat-1);
367: if (!st->T) SETERRQ(PetscObjectComm((PetscObject)st),PETSC_ERR_POINTER,"There are no transformed matrices");
368: *T = st->T[k];
369: return(0);
370: }
372: /*@
373: STGetNumMatrices - Returns the number of matrices stored in the ST.
375: Not collective
377: Input Parameter:
378: . st - the spectral transformation context
380: Output Parameters:
381: . n - the number of matrices passed in STSetMatrices()
383: Level: intermediate
385: .seealso: STSetMatrices()
386: @*/
387: PetscErrorCode STGetNumMatrices(ST st,PetscInt *n)388: {
392: *n = st->nmat;
393: return(0);
394: }
396: /*@
397: STResetMatrixState - Resets the stored state of the matrices in the ST.
399: Logically Collective on st
401: Input Parameter:
402: . st - the spectral transformation context
404: Note:
405: This is useful in solvers where the user matrices are modified during
406: the computation, as in nonlinear inverse iteration. The effect is that
407: STGetMatrix() will retrieve the modified matrices as if they were
408: the matrices originally provided by the user.
410: Level: developer
412: .seealso: STGetMatrix(), EPSPowerSetNonlinear()
413: @*/
414: PetscErrorCode STResetMatrixState(ST st)415: {
416: PetscInt i;
420: for (i=0;i<st->nmat;i++) st->Astate[i] = ((PetscObject)st->A[i])->state;
421: return(0);
422: }
424: /*@
425: STSetShift - Sets the shift associated with the spectral transformation.
427: Logically Collective on st
429: Input Parameters:
430: + st - the spectral transformation context
431: - shift - the value of the shift
433: Notes:
434: In some spectral transformations, changing the shift may have associated
435: a lot of work, for example recomputing a factorization.
437: This function is normally not directly called by users, since the shift is
438: indirectly set by EPSSetTarget().
440: Level: intermediate
442: .seealso: EPSSetTarget(), STGetShift(), STSetDefaultShift()
443: @*/
444: PetscErrorCode STSetShift(ST st,PetscScalar shift)445: {
452: if (st->sigma != shift) {
453: STCheckNotSeized(st,1);
454: if (st->state==ST_STATE_SETUP && st->ops->setshift) {
455: (*st->ops->setshift)(st,shift);
456: }
457: st->sigma = shift;
458: }
459: st->sigma_set = PETSC_TRUE;
460: return(0);
461: }
463: /*@
464: STGetShift - Gets the shift associated with the spectral transformation.
466: Not Collective
468: Input Parameter:
469: . st - the spectral transformation context
471: Output Parameter:
472: . shift - the value of the shift
474: Level: intermediate
476: .seealso: STSetShift()
477: @*/
478: PetscErrorCode STGetShift(ST st,PetscScalar* shift)479: {
483: *shift = st->sigma;
484: return(0);
485: }
487: /*@
488: STSetDefaultShift - Sets the value of the shift that should be employed if
489: the user did not specify one.
491: Logically Collective on st
493: Input Parameters:
494: + st - the spectral transformation context
495: - defaultshift - the default value of the shift
497: Level: developer
499: .seealso: STSetShift()
500: @*/
501: PetscErrorCode STSetDefaultShift(ST st,PetscScalar defaultshift)502: {
506: if (st->defsigma != defaultshift) {
507: st->defsigma = defaultshift;
508: st->state = ST_STATE_INITIAL;
509: st->opready = PETSC_FALSE;
510: }
511: return(0);
512: }
514: /*@
515: STScaleShift - Multiply the shift with a given factor.
517: Logically Collective on st
519: Input Parameters:
520: + st - the spectral transformation context
521: - factor - the scaling factor
523: Note:
524: This function does not update the transformation matrices, as opposed to
525: STSetShift().
527: Level: developer
529: .seealso: STSetShift()
530: @*/
531: PetscErrorCode STScaleShift(ST st,PetscScalar factor)532: {
536: st->sigma *= factor;
537: return(0);
538: }
540: /*@
541: STSetBalanceMatrix - Sets the diagonal matrix to be used for balancing.
543: Collective on st
545: Input Parameters:
546: + st - the spectral transformation context
547: - D - the diagonal matrix (represented as a vector)
549: Notes:
550: If this matrix is set, STApply will effectively apply D*OP*D^{-1}. Use NULL
551: to reset a previously passed D.
553: Balancing is usually set via EPSSetBalance, but the advanced user may use
554: this function to bypass the usual balancing methods.
556: Level: developer
558: .seealso: EPSSetBalance(), STApply(), STGetBalanceMatrix()
559: @*/
560: PetscErrorCode STSetBalanceMatrix(ST st,Vec D)561: {
566: if (st->D == D) return(0);
567: STCheckNotSeized(st,1);
568: if (D) {
571: PetscObjectReference((PetscObject)D);
572: }
573: VecDestroy(&st->D);
574: st->D = D;
575: st->state = ST_STATE_INITIAL;
576: st->opready = PETSC_FALSE;
577: return(0);
578: }
580: /*@
581: STGetBalanceMatrix - Gets the balance matrix used by the spectral transformation.
583: Not collective, but vector is shared by all processors that share the ST585: Input Parameter:
586: . st - the spectral transformation context
588: Output Parameter:
589: . D - the diagonal matrix (represented as a vector)
591: Note:
592: If the matrix was not set, a null pointer will be returned.
594: Level: developer
596: .seealso: STSetBalanceMatrix()
597: @*/
598: PetscErrorCode STGetBalanceMatrix(ST st,Vec *D)599: {
603: *D = st->D;
604: return(0);
605: }
607: /*@C
608: STMatCreateVecs - Get vector(s) compatible with the ST matrices.
610: Collective on st
612: Input Parameter:
613: . st - the spectral transformation context
615: Output Parameters:
616: + right - (optional) vector that the matrix can be multiplied against
617: - left - (optional) vector that the matrix vector product can be stored in
619: Level: developer
620: @*/
621: PetscErrorCode STMatCreateVecs(ST st,Vec *right,Vec *left)622: {
626: STCheckMatrices(st,1);
627: MatCreateVecs(st->A[0],right,left);
628: return(0);
629: }
631: /*@C
632: STMatCreateVecsEmpty - Get vector(s) compatible with the ST matrices, i.e. with the same
633: parallel layout, but without internal array.
635: Collective on st
637: Input Parameter:
638: . st - the spectral transformation context
640: Output Parameters:
641: + right - (optional) vector that the matrix can be multiplied against
642: - left - (optional) vector that the matrix vector product can be stored in
644: Level: developer
646: .seealso: MatCreateVecsEmpty()
647: @*/
648: PetscErrorCode STMatCreateVecsEmpty(ST st,Vec *right,Vec *left)649: {
653: STCheckMatrices(st,1);
654: MatCreateVecsEmpty(st->A[0],right,left);
655: return(0);
656: }
658: /*@
659: STMatGetSize - Returns the number of rows and columns of the ST matrices.
661: Not Collective
663: Input Parameter:
664: . st - the spectral transformation context
666: Output Parameters:
667: + m - the number of global rows
668: - n - the number of global columns
670: Level: developer
671: @*/
672: PetscErrorCode STMatGetSize(ST st,PetscInt *m,PetscInt *n)673: {
677: STCheckMatrices(st,1);
678: MatGetSize(st->A[0],m,n);
679: return(0);
680: }
682: /*@
683: STMatGetLocalSize - Returns the number of local rows and columns of the ST matrices.
685: Not Collective
687: Input Parameter:
688: . st - the spectral transformation context
690: Output Parameters:
691: + m - the number of local rows
692: - n - the number of local columns
694: Level: developer
695: @*/
696: PetscErrorCode STMatGetLocalSize(ST st,PetscInt *m,PetscInt *n)697: {
701: STCheckMatrices(st,1);
702: MatGetLocalSize(st->A[0],m,n);
703: return(0);
704: }
706: /*@C
707: STSetOptionsPrefix - Sets the prefix used for searching for all
708: ST options in the database.
710: Logically Collective on st
712: Input Parameters:
713: + st - the spectral transformation context
714: - prefix - the prefix string to prepend to all ST option requests
716: Notes:
717: A hyphen (-) must NOT be given at the beginning of the prefix name.
718: The first character of all runtime options is AUTOMATICALLY the
719: hyphen.
721: Level: advanced
723: .seealso: STAppendOptionsPrefix(), STGetOptionsPrefix()
724: @*/
725: PetscErrorCode STSetOptionsPrefix(ST st,const char *prefix)726: {
731: if (!st->ksp) { STGetKSP(st,&st->ksp); }
732: KSPSetOptionsPrefix(st->ksp,prefix);
733: KSPAppendOptionsPrefix(st->ksp,"st_");
734: PetscObjectSetOptionsPrefix((PetscObject)st,prefix);
735: return(0);
736: }
738: /*@C
739: STAppendOptionsPrefix - Appends to the prefix used for searching for all
740: ST options in the database.
742: Logically Collective on st
744: Input Parameters:
745: + st - the spectral transformation context
746: - prefix - the prefix string to prepend to all ST option requests
748: Notes:
749: A hyphen (-) must NOT be given at the beginning of the prefix name.
750: The first character of all runtime options is AUTOMATICALLY the
751: hyphen.
753: Level: advanced
755: .seealso: STSetOptionsPrefix(), STGetOptionsPrefix()
756: @*/
757: PetscErrorCode STAppendOptionsPrefix(ST st,const char *prefix)758: {
763: PetscObjectAppendOptionsPrefix((PetscObject)st,prefix);
764: if (!st->ksp) { STGetKSP(st,&st->ksp); }
765: KSPSetOptionsPrefix(st->ksp,((PetscObject)st)->prefix);
766: KSPAppendOptionsPrefix(st->ksp,"st_");
767: return(0);
768: }
770: /*@C
771: STGetOptionsPrefix - Gets the prefix used for searching for all
772: ST options in the database.
774: Not Collective
776: Input Parameters:
777: . st - the spectral transformation context
779: Output Parameters:
780: . prefix - pointer to the prefix string used, is returned
782: Note:
783: On the Fortran side, the user should pass in a string 'prefix' of
784: sufficient length to hold the prefix.
786: Level: advanced
788: .seealso: STSetOptionsPrefix(), STAppendOptionsPrefix()
789: @*/
790: PetscErrorCode STGetOptionsPrefix(ST st,const char *prefix[])791: {
797: PetscObjectGetOptionsPrefix((PetscObject)st,prefix);
798: return(0);
799: }
801: /*@C
802: STView - Prints the ST data structure.
804: Collective on st
806: Input Parameters:
807: + st - the ST context
808: - viewer - optional visualization context
810: Note:
811: The available visualization contexts include
812: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
813: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
814: output where only the first processor opens
815: the file. All other processors send their
816: data to the first processor to print.
818: The user can open an alternative visualization contexts with
819: PetscViewerASCIIOpen() (output to a specified file).
821: Level: beginner
823: .seealso: EPSView()
824: @*/
825: PetscErrorCode STView(ST st,PetscViewer viewer)826: {
828: STType cstr;
829: const char* pat=NULL;
830: char str[50];
831: PetscBool isascii,isstring;
835: if (!viewer) {
836: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)st),&viewer);
837: }
841: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
842: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERSTRING,&isstring);
843: if (isascii) {
844: PetscObjectPrintClassNamePrefixType((PetscObject)st,viewer);
845: if (st->ops->view) {
846: PetscViewerASCIIPushTab(viewer);
847: (*st->ops->view)(st,viewer);
848: PetscViewerASCIIPopTab(viewer);
849: }
850: SlepcSNPrintfScalar(str,50,st->sigma,PETSC_FALSE);
851: PetscViewerASCIIPrintf(viewer," shift: %s\n",str);
852: PetscViewerASCIIPrintf(viewer," number of matrices: %D\n",st->nmat);
853: switch (st->matmode) {
854: case ST_MATMODE_COPY:
855: break;
856: case ST_MATMODE_INPLACE:
857: PetscViewerASCIIPrintf(viewer," shifting the matrix and unshifting at exit\n");
858: break;
859: case ST_MATMODE_SHELL:
860: PetscViewerASCIIPrintf(viewer," using a shell matrix\n");
861: break;
862: }
863: if (st->nmat>1 && st->matmode != ST_MATMODE_SHELL) {
864: switch (st->str) {
865: case SAME_NONZERO_PATTERN: pat = "same nonzero pattern";break;
866: case DIFFERENT_NONZERO_PATTERN: pat = "different nonzero pattern";break;
867: case SUBSET_NONZERO_PATTERN: pat = "subset nonzero pattern";break;
868: }
869: PetscViewerASCIIPrintf(viewer," all matrices have %s\n",pat);
870: }
871: if (st->transform && st->nmat>2) {
872: PetscViewerASCIIPrintf(viewer," computing transformed matrices\n");
873: }
874: } else if (isstring) {
875: STGetType(st,&cstr);
876: PetscViewerStringSPrintf(viewer," %-7.7s",cstr);
877: if (st->ops->view) { (*st->ops->view)(st,viewer); }
878: }
879: if (st->usesksp) {
880: if (!st->ksp) { STGetKSP(st,&st->ksp); }
881: PetscViewerASCIIPushTab(viewer);
882: KSPView(st->ksp,viewer);
883: PetscViewerASCIIPopTab(viewer);
884: }
885: return(0);
886: }
888: /*@C
889: STViewFromOptions - View from options
891: Collective on ST893: Input Parameters:
894: + st - the spectral transformation context
895: . obj - optional object
896: - name - command line option
898: Level: intermediate
900: .seealso: STView(), STCreate()
901: @*/
902: PetscErrorCode STViewFromOptions(ST st,PetscObject obj,const char name[])903: {
908: PetscObjectViewFromOptions((PetscObject)st,obj,name);
909: return(0);
910: }
912: /*@C
913: STRegister - Adds a method to the spectral transformation package.
915: Not collective
917: Input Parameters:
918: + name - name of a new user-defined transformation
919: - function - routine to create method context
921: Notes:
922: STRegister() may be called multiple times to add several user-defined
923: spectral transformations.
925: Sample usage:
926: .vb
927: STRegister("my_transform",MyTransformCreate);
928: .ve
930: Then, your spectral transform can be chosen with the procedural interface via
931: $ STSetType(st,"my_transform")
932: or at runtime via the option
933: $ -st_type my_transform
935: Level: advanced
937: .seealso: STRegisterAll()
938: @*/
939: PetscErrorCode STRegister(const char *name,PetscErrorCode (*function)(ST))940: {
944: STInitializePackage();
945: PetscFunctionListAdd(&STList,name,function);
946: return(0);
947: }