Actual source code: slepcst.h
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: Spectral transformation module for eigenvalue problems
12: */
14: #pragma once
16: #include <slepcsys.h>
17: #include <slepcbv.h>
18: #include <petscksp.h>
20: /* SUBMANSEC = ST */
22: SLEPC_EXTERN PetscErrorCode STInitializePackage(void);
23: SLEPC_EXTERN PetscErrorCode STFinalizePackage(void);
25: /*S
26: ST - Spectral transformation, encapsulates the functionality required
27: for acceleration techniques based on the transformation of the spectrum,
28: e.g., shift-and-invert.
30: Level: beginner
32: .seealso: [](ch:st), `STCreate()`
33: S*/
34: typedef struct _p_ST* ST;
36: /*J
37: STType - String with the name of the spectral transformation type.
39: Level: beginner
41: .seealso: [](ch:st), `ST`, `STSetType()`
42: J*/
43: typedef const char *STType;
44: #define STSHIFT "shift"
45: #define STSINVERT "sinvert"
46: #define STCAYLEY "cayley"
47: #define STPRECOND "precond"
48: #define STFILTER "filter"
49: #define STSHELL "shell"
51: /* Logging support */
52: SLEPC_EXTERN PetscClassId ST_CLASSID;
54: SLEPC_EXTERN PetscErrorCode STCreate(MPI_Comm,ST*);
55: SLEPC_EXTERN PetscErrorCode STDestroy(ST*);
56: SLEPC_EXTERN PetscErrorCode STReset(ST);
57: SLEPC_EXTERN PetscErrorCode STSetType(ST,STType);
58: SLEPC_EXTERN PetscErrorCode STGetType(ST,STType*);
59: SLEPC_EXTERN PetscErrorCode STSetMatrices(ST,PetscInt,Mat[]);
60: SLEPC_EXTERN PetscErrorCode STGetMatrix(ST,PetscInt,Mat*);
61: SLEPC_EXTERN PetscErrorCode STGetMatrixTransformed(ST,PetscInt,Mat*);
62: SLEPC_EXTERN PetscErrorCode STGetNumMatrices(ST,PetscInt*);
63: SLEPC_EXTERN PetscErrorCode STGetOperator(ST,Mat*);
64: SLEPC_EXTERN PetscErrorCode STRestoreOperator(ST,Mat*);
65: SLEPC_EXTERN PetscErrorCode STSetUp(ST);
66: SLEPC_EXTERN PetscErrorCode STSetFromOptions(ST);
67: SLEPC_EXTERN PetscErrorCode STView(ST,PetscViewer);
68: SLEPC_EXTERN PetscErrorCode STViewFromOptions(ST,PetscObject,const char[]);
70: PETSC_DEPRECATED_FUNCTION(3, 15, 0, "STSetMatrices()", ) static inline PetscErrorCode STSetOperators(ST st,PetscInt n,Mat *A) {return STSetMatrices(st,n,A);}
71: PETSC_DEPRECATED_FUNCTION(3, 15, 0, "STGetMatrix()", ) static inline PetscErrorCode STGetOperators(ST st,PetscInt k,Mat *A) {return STGetMatrix(st,k,A);}
72: PETSC_DEPRECATED_FUNCTION(3, 15, 0, "STGetMatrixTransformed()", ) static inline PetscErrorCode STGetTOperators(ST st,PetscInt k,Mat *A) {return STGetMatrixTransformed(st,k,A);}
73: PETSC_DEPRECATED_FUNCTION(3, 15, 0, "STGetOperator() followed by MatComputeOperator()", ) static inline PetscErrorCode STComputeExplicitOperator(ST st,Mat *A)
74: {
75: Mat Op;
77: PetscFunctionBegin;
78: PetscCall(STGetOperator(st,&Op));
79: PetscCall(MatComputeOperator(Op,MATAIJ,A));
80: PetscCall(STRestoreOperator(st,&Op));
81: PetscFunctionReturn(PETSC_SUCCESS);
82: }
84: SLEPC_EXTERN PetscErrorCode STApply(ST,Vec,Vec);
85: SLEPC_EXTERN PetscErrorCode STApplyMat(ST,Mat,Mat);
86: SLEPC_EXTERN PetscErrorCode STApplyTranspose(ST,Vec,Vec);
87: SLEPC_EXTERN PetscErrorCode STApplyHermitianTranspose(ST,Vec,Vec);
88: SLEPC_EXTERN PetscErrorCode STMatMult(ST,PetscInt,Vec,Vec);
89: SLEPC_EXTERN PetscErrorCode STMatMultTranspose(ST,PetscInt,Vec,Vec);
90: SLEPC_EXTERN PetscErrorCode STMatMultHermitianTranspose(ST,PetscInt,Vec,Vec);
91: SLEPC_EXTERN PetscErrorCode STMatSolve(ST,Vec,Vec);
92: SLEPC_EXTERN PetscErrorCode STMatSolveTranspose(ST,Vec,Vec);
93: SLEPC_EXTERN PetscErrorCode STMatSolveHermitianTranspose(ST,Vec,Vec);
94: SLEPC_EXTERN PetscErrorCode STMatMatSolve(ST,Mat,Mat);
95: SLEPC_EXTERN PetscErrorCode STGetBilinearForm(ST,Mat*);
96: SLEPC_EXTERN PetscErrorCode STMatSetUp(ST,PetscScalar,PetscScalar[]);
97: SLEPC_EXTERN PetscErrorCode STPostSolve(ST);
98: SLEPC_EXTERN PetscErrorCode STResetMatrixState(ST);
99: SLEPC_EXTERN PetscErrorCode STSetWorkVecs(ST,PetscInt);
101: SLEPC_EXTERN PetscErrorCode STSetKSP(ST,KSP);
102: SLEPC_EXTERN PetscErrorCode STGetKSP(ST,KSP*);
103: SLEPC_EXTERN PetscErrorCode STSetShift(ST,PetscScalar);
104: SLEPC_EXTERN PetscErrorCode STGetShift(ST,PetscScalar*);
105: SLEPC_EXTERN PetscErrorCode STSetDefaultShift(ST,PetscScalar);
106: SLEPC_EXTERN PetscErrorCode STScaleShift(ST,PetscScalar);
107: SLEPC_EXTERN PetscErrorCode STSetBalanceMatrix(ST,Vec);
108: SLEPC_EXTERN PetscErrorCode STGetBalanceMatrix(ST,Vec*);
109: SLEPC_EXTERN PetscErrorCode STSetTransform(ST,PetscBool);
110: SLEPC_EXTERN PetscErrorCode STGetTransform(ST,PetscBool*);
111: SLEPC_EXTERN PetscErrorCode STSetStructured(ST,PetscBool);
112: SLEPC_EXTERN PetscErrorCode STGetStructured(ST,PetscBool*);
114: SLEPC_EXTERN PetscErrorCode STSetOptionsPrefix(ST,const char[]);
115: SLEPC_EXTERN PetscErrorCode STAppendOptionsPrefix(ST,const char[]);
116: SLEPC_EXTERN PetscErrorCode STGetOptionsPrefix(ST,const char*[]);
118: SLEPC_EXTERN PetscErrorCode STBackTransform(ST,PetscInt,PetscScalar[],PetscScalar[]);
119: SLEPC_EXTERN PetscErrorCode STIsInjective(ST,PetscBool*);
121: SLEPC_EXTERN PetscErrorCode STCheckNullSpace(ST,BV);
123: SLEPC_EXTERN PetscErrorCode STSetPreconditionerMat(ST,Mat);
124: SLEPC_EXTERN PetscErrorCode STGetPreconditionerMat(ST,Mat*);
125: SLEPC_EXTERN PetscErrorCode STSetSplitPreconditioner(ST,PetscInt,Mat[],MatStructure);
126: SLEPC_EXTERN PetscErrorCode STGetSplitPreconditionerTerm(ST,PetscInt,Mat*);
127: SLEPC_EXTERN PetscErrorCode STGetSplitPreconditionerInfo(ST,PetscInt*,MatStructure*);
129: SLEPC_EXTERN PetscErrorCode STMatCreateVecs(ST,Vec*,Vec*);
130: SLEPC_EXTERN PetscErrorCode STMatCreateVecsEmpty(ST,Vec*,Vec*);
131: SLEPC_EXTERN PetscErrorCode STMatGetSize(ST,PetscInt*,PetscInt*);
132: SLEPC_EXTERN PetscErrorCode STMatGetLocalSize(ST,PetscInt*,PetscInt*);
134: /*E
135: STMatMode - Determines how to handle the coefficient matrix of the linear
136: system associated with the spectral transformation.
138: Values:
139: + `ST_MATMODE_COPY` - the coefficient matrix is built explicitly on a copy of $A$
140: . `ST_MATMODE_INPLACE` - the coefficient matrix is built explicitly overwriting $A$
141: - `ST_MATMODE_SHELL` - the coefficient matrix is handled implicitly
143: Level: intermediate
145: .seealso: [](ch:st), `STSetMatMode()`, `STGetMatMode()`
146: E*/
147: typedef enum { ST_MATMODE_COPY,
148: ST_MATMODE_INPLACE,
149: ST_MATMODE_SHELL } STMatMode;
150: SLEPC_EXTERN const char *STMatModes[];
152: /*MC
153: ST_MATMODE_COPY - The coefficient matrix of the linear system, $A-\sigma B$, is
154: built explicitly on a copy of $A$.
156: Note:
157: If memory is an issue, one may prefer one of the other two modes.
159: Level: intermediate
161: .seealso: [](ch:st), `STMatMode`, `STSetMatMode()`, `ST_MATMODE_INPLACE`, `ST_MATMODE_SHELL`
162: M*/
164: /*MC
165: ST_MATMODE_INPLACE - The coefficient matrix of the linear system, $A-\sigma B$, is
166: built explicitly overwritting $A$.
168: Note:
169: This mode uses less memory than `ST_MATMODE_COPY`, but it modifies $A$. This
170: alteration of $A$ is reversed after the eigensolution process has finished,
171: but due to roundoff, the result might not be exactly equal to the original $A$.
173: Level: intermediate
175: .seealso: [](ch:st), `STMatMode`, `STSetMatMode()`, `ST_MATMODE_COPY`, `ST_MATMODE_SHELL`
176: M*/
178: /*MC
179: ST_MATMODE_SHELL - The coefficient matrix of the linear system, $A-\sigma B$, is
180: not built explicitly, and instead it is handled implicitly via a `MATSHELL`.
182: Note:
183: This mode severely restricts the number of possibilities available for solving
184: the linear system via `KSP`.
186: Level: intermediate
188: .seealso: [](ch:st), `STMatMode`, `STSetMatMode()`, `ST_MATMODE_COPY`, `ST_MATMODE_INPLACE`
189: M*/
191: SLEPC_EXTERN PetscErrorCode STSetMatMode(ST,STMatMode);
192: SLEPC_EXTERN PetscErrorCode STGetMatMode(ST,STMatMode*);
193: SLEPC_EXTERN PetscErrorCode STSetMatStructure(ST,MatStructure);
194: SLEPC_EXTERN PetscErrorCode STGetMatStructure(ST,MatStructure*);
196: SLEPC_EXTERN PetscFunctionList STList;
197: SLEPC_EXTERN PetscErrorCode STRegister(const char[],PetscErrorCode(*)(ST));
199: /* --------- options specific to particular spectral transformations-------- */
201: /*S
202: STShellDestroyFn - A prototype of a function for destroying the context
203: set with `STShellSetContext()` in an `STSHELL`.
205: Calling Sequence:
206: . st - the spectral transformation context
208: Level: advanced
210: .seealso: [](ch:st), `STShellSetDestroy()`, `STShellSetContext()`
211: S*/
212: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode STShellDestroyFn(ST st);
214: /*S
215: STShellApplyFn - A prototype of a function for the user-defined `STApply()`
216: operation in an `STSHELL`.
218: Calling Sequence:
219: + st - the spectral transformation context
220: . xin - input vector
221: - xout - output vector
223: Level: advanced
225: .seealso: [](ch:st), `STShellSetApply()`, `STApply()`
226: S*/
227: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode STShellApplyFn(ST st,Vec xin,Vec xout);
229: /*S
230: STShellApplyTransposeFn - A prototype of a function for the user-defined `STApplyTranspose()`
231: operation in an `STSHELL`.
233: Calling Sequence:
234: + st - the spectral transformation context
235: . xin - input vector
236: - xout - output vector
238: Level: advanced
240: .seealso: [](ch:st), `STShellSetApplyTranspose()`, `STApplyTranspose()`
241: S*/
242: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode STShellApplyTransposeFn(ST st,Vec xin,Vec xout);
244: /*S
245: STShellApplyHermitianTransposeFn - A prototype of a function for the user-defined
246: `STApplyHermitianTranspose()` operation in an `STSHELL`.
248: Calling Sequence:
249: + st - the spectral transformation context
250: . xin - input vector
251: - xout - output vector
253: Level: advanced
255: .seealso: [](ch:st), `STShellSetApplyHermitianTranspose()`, `STApplyHermitianTranspose()`
256: S*/
257: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode STShellApplyHermitianTransposeFn(ST st,Vec xin,Vec xout);
259: /*S
260: STShellBackTransformFn - A prototype of a function for the user-defined `STBackTransform()`
261: operation in an `STSHELL`.
263: Calling Sequence:
264: + st - the spectral transformation context
265: . n - number of eigenvalues to be backtransformed
266: . eigr - pointer to the real parts of the eigenvalues to transform back
267: - eigi - pointer to the imaginary parts
269: Level: advanced
271: .seealso: [](ch:st), `STShellSetBackTransform()`, `STBackTransform()`
272: S*/
273: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode STShellBackTransformFn(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi);
275: SLEPC_EXTERN PetscErrorCode STShellGetContext(ST,PetscCtxRt);
276: SLEPC_EXTERN PetscErrorCode STShellSetContext(ST,PetscCtx);
277: SLEPC_EXTERN PetscErrorCode STShellSetDestroy(ST,STShellDestroyFn*);
278: SLEPC_EXTERN PetscErrorCode STShellSetApply(ST,STShellApplyFn*);
279: SLEPC_EXTERN PetscErrorCode STShellSetApplyTranspose(ST,STShellApplyTransposeFn*);
280: SLEPC_EXTERN PetscErrorCode STShellSetApplyHermitianTranspose(ST,STShellApplyHermitianTransposeFn*);
281: SLEPC_EXTERN PetscErrorCode STShellSetBackTransform(ST,STShellBackTransformFn*);
283: SLEPC_EXTERN PetscErrorCode STCayleyGetAntishift(ST,PetscScalar*);
284: SLEPC_EXTERN PetscErrorCode STCayleySetAntishift(ST,PetscScalar);
286: PETSC_DEPRECATED_FUNCTION(3, 15, 0, "STGetPreconditionerMat()", ) static inline PetscErrorCode STPrecondGetMatForPC(ST st,Mat *A) {return STGetPreconditionerMat(st,A);}
287: PETSC_DEPRECATED_FUNCTION(3, 15, 0, "STSetPreconditionerMat()", ) static inline PetscErrorCode STPrecondSetMatForPC(ST st,Mat A) {return STSetPreconditionerMat(st,A);}
288: SLEPC_EXTERN PetscErrorCode STPrecondGetKSPHasMat(ST,PetscBool*);
289: SLEPC_EXTERN PetscErrorCode STPrecondSetKSPHasMat(ST,PetscBool);
291: /*E
292: STFilterType - Selects the method used to build the filter.
294: Values:
295: + `ST_FILTER_FILTLAN` - the filter is built with FILTLAN (Filtered Lanczos)
296: - `ST_FILTER_CHEBYSHEV` - the filter is built with a Chebyshev series
298: Level: intermediate
300: .seealso: [](ch:st), `STFilterSetType()`, `STGetFilterType()`
301: E*/
302: typedef enum { ST_FILTER_FILTLAN = 1,
303: ST_FILTER_CHEBYSHEV = 2 } STFilterType;
304: SLEPC_EXTERN const char *STFilterTypes[];
306: /*MC
307: ST_FILTER_FILTLAN - The polynomial filter is built with FILTLAN (Filtered Lanczos).
309: Note:
310: This filter implements the Filtered Lanczos method {cite:p}`Fan12`. In fact,
311: the implementation adapts files from the FILTLAN package by the same authors,
312: which have been converted for a native integration in SLEPc.
314: Level: intermediate
316: .seealso: [](ch:st), `STFilterType`, `STFilterSetType()`, `ST_FILTER_CHEBYSHEV`
317: M*/
319: /*MC
320: ST_FILTER_CHEBYSHEV - The filter is based on a truncated Chebyshev series.
322: Note:
323: The application of the filter is implemented by means of a Clenshaw
324: algorithm. This filter also supports using damping, see `STFilterSetDamping()`.
326: Level: intermediate
328: .seealso: [](ch:st), `STFilterType`, `STSetFilterType()`, `STFilterSetDamping()`, `ST_FILTER_FILTLAN`
329: M*/
331: /*E
332: STFilterDamping - The damping type used to build the filter.
334: Values:
335: + `ST_FILTER_DAMPING_NONE` - no damping
336: . `ST_FILTER_DAMPING_JACKSON` - Jackson damping
337: . `ST_FILTER_DAMPING_LANCZOS` - Lanczos damping
338: - `ST_FILTER_DAMPING_FEJER` - Fejer damping
340: Notes:
341: The default is no damping.
343: For the definition of the damping coefficients, see for instance {cite:p}`Pie16`
345: Level: advanced
347: .seealso: [](ch:st), `STSetFilterDamping()`, `STGetFilterDamping()`
348: E*/
349: typedef enum { ST_FILTER_DAMPING_NONE,
350: ST_FILTER_DAMPING_JACKSON,
351: ST_FILTER_DAMPING_LANCZOS,
352: ST_FILTER_DAMPING_FEJER } STFilterDamping;
353: SLEPC_EXTERN const char *STFilterDampings[];
355: /*MC
356: ST_MATMODE_COPY - The coefficient matrix of the linear system, $A-\sigma B$, is
357: built explicitly on a copy of $A$.
359: Note:
360: If memory is an issue, one may prefer one of the other two modes.
362: Level: intermediate
364: .seealso: [](ch:st), `STMatMode`, `STSetMatMode()`, `ST_MATMODE_INPLACE`, `ST_MATMODE_SHELL`
365: M*/
367: SLEPC_EXTERN PetscErrorCode STFilterSetType(ST,STFilterType);
368: SLEPC_EXTERN PetscErrorCode STFilterGetType(ST,STFilterType*);
369: SLEPC_EXTERN PetscErrorCode STFilterSetInterval(ST,PetscReal,PetscReal);
370: SLEPC_EXTERN PetscErrorCode STFilterGetInterval(ST,PetscReal*,PetscReal*);
371: SLEPC_EXTERN PetscErrorCode STFilterSetRange(ST,PetscReal,PetscReal);
372: SLEPC_EXTERN PetscErrorCode STFilterGetRange(ST,PetscReal*,PetscReal*);
373: SLEPC_EXTERN PetscErrorCode STFilterSetDegree(ST,PetscInt);
374: SLEPC_EXTERN PetscErrorCode STFilterGetDegree(ST,PetscInt*);
375: SLEPC_EXTERN PetscErrorCode STFilterGetThreshold(ST,PetscReal*);
376: SLEPC_EXTERN PetscErrorCode STFilterSetDamping(ST,STFilterDamping);
377: SLEPC_EXTERN PetscErrorCode STFilterGetDamping(ST,STFilterDamping*);