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:    STShellApplyFn - A prototype of a function for the user-defined `STApply()`
203:    operation in an `STSHELL`.

205:    Calling Sequence:
206: +  st   - the spectral transformation context
207: .  xin  - input vector
208: -  xout - output vector

210:    Level: advanced

212: .seealso: [](ch:st), `STShellSetApply()`, `STApply()`
213: S*/
214: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode STShellApplyFn(ST st,Vec xin,Vec xout);

216: /*S
217:    STShellApplyTransposeFn - A prototype of a function for the user-defined `STApplyTranspose()`
218:    operation in an `STSHELL`.

220:    Calling Sequence:
221: +  st   - the spectral transformation context
222: .  xin  - input vector
223: -  xout - output vector

225:    Level: advanced

227: .seealso: [](ch:st), `STShellSetApplyTranspose()`, `STApplyTranspose()`
228: S*/
229: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode STShellApplyTransposeFn(ST st,Vec xin,Vec xout);

231: /*S
232:    STShellApplyHermitianTransposeFn - A prototype of a function for the user-defined
233:    `STApplyHermitianTranspose()` operation in an `STSHELL`.

235:    Calling Sequence:
236: +  st   - the spectral transformation context
237: .  xin  - input vector
238: -  xout - output vector

240:    Level: advanced

242: .seealso: [](ch:st), `STShellSetApplyHermitianTranspose()`, `STApplyHermitianTranspose()`
243: S*/
244: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode STShellApplyHermitianTransposeFn(ST st,Vec xin,Vec xout);

246: /*S
247:    STShellBackTransformFn - A prototype of a function for the user-defined `STBackTransform()`
248:    operation in an `STSHELL`.

250:    Calling Sequence:
251: +  st   - the spectral transformation context
252: .  n    - number of eigenvalues to be backtransformed
253: .  eigr - pointer to the real parts of the eigenvalues to transform back
254: -  eigi - pointer to the imaginary parts

256:    Level: advanced

258: .seealso: [](ch:st), `STShellSetBackTransform()`, `STBackTransform()`
259: S*/
260: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode STShellBackTransformFn(ST st,PetscInt n,PetscScalar *eigr,PetscScalar *eigi);

262: SLEPC_EXTERN PetscErrorCode STShellGetContext(ST,void*);
263: SLEPC_EXTERN PetscErrorCode STShellSetContext(ST,void*);
264: SLEPC_EXTERN PetscErrorCode STShellSetApply(ST,STShellApplyFn*);
265: SLEPC_EXTERN PetscErrorCode STShellSetApplyTranspose(ST,STShellApplyTransposeFn*);
266: SLEPC_EXTERN PetscErrorCode STShellSetApplyHermitianTranspose(ST,STShellApplyHermitianTransposeFn*);
267: SLEPC_EXTERN PetscErrorCode STShellSetBackTransform(ST,STShellBackTransformFn*);

269: SLEPC_EXTERN PetscErrorCode STCayleyGetAntishift(ST,PetscScalar*);
270: SLEPC_EXTERN PetscErrorCode STCayleySetAntishift(ST,PetscScalar);

272: PETSC_DEPRECATED_FUNCTION(3, 15, 0, "STGetPreconditionerMat()", ) static inline PetscErrorCode STPrecondGetMatForPC(ST st,Mat *A) {return STGetPreconditionerMat(st,A);}
273: PETSC_DEPRECATED_FUNCTION(3, 15, 0, "STSetPreconditionerMat()", ) static inline PetscErrorCode STPrecondSetMatForPC(ST st,Mat A) {return STSetPreconditionerMat(st,A);}
274: SLEPC_EXTERN PetscErrorCode STPrecondGetKSPHasMat(ST,PetscBool*);
275: SLEPC_EXTERN PetscErrorCode STPrecondSetKSPHasMat(ST,PetscBool);

277: /*E
278:    STFilterType - Selects the method used to build the filter.

280:    Values:
281: +  `ST_FILTER_FILTLAN`   - the filter is built with FILTLAN (Filtered Lanczos)
282: -  `ST_FILTER_CHEBYSHEV` - the filter is built with a Chebyshev series

284:    Level: intermediate

286: .seealso: [](ch:st), `STFilterSetType()`, `STGetFilterType()`
287: E*/
288: typedef enum { ST_FILTER_FILTLAN   = 1,
289:                ST_FILTER_CHEBYSHEV = 2 } STFilterType;
290: SLEPC_EXTERN const char *STFilterTypes[];

292: /*MC
293:    ST_FILTER_FILTLAN - The polynomial filter is built with FILTLAN (Filtered Lanczos).

295:    Note:
296:    This filter implements the Filtered Lanczos method {cite:p}`Fan12`. In fact,
297:    the implementation adapts files from the FILTLAN package by the same authors,
298:    which have been converted for a native integration in SLEPc.

300:    Level: intermediate

302: .seealso: [](ch:st), `STFilterType`, `STFilterSetType()`, `ST_FILTER_CHEBYSHEV`
303: M*/

305: /*MC
306:    ST_FILTER_CHEBYSHEV - The filter is based on a truncated Chebyshev series.

308:    Note:
309:    The application of the filter is implemented by means of a Clenshaw
310:    algorithm. This filter also supports using damping, see `STFilterSetDamping()`.

312:    Level: intermediate

314: .seealso: [](ch:st), `STFilterType`, `STSetFilterType()`, `STFilterSetDamping()`, `ST_FILTER_FILTLAN`
315: M*/

317: /*E
318:    STFilterDamping - The damping type used to build the filter.

320:    Values:
321: +  `ST_FILTER_DAMPING_NONE`    - no damping
322: .  `ST_FILTER_DAMPING_JACKSON` - Jackson damping
323: .  `ST_FILTER_DAMPING_LANCZOS` - Lanczos damping
324: -  `ST_FILTER_DAMPING_FEJER`   - Fejer damping

326:    Notes:
327:    The default is no damping.

329:    For the definition of the damping coefficients, see for instance {cite:p}`Pie16`

331:    Level: advanced

333: .seealso: [](ch:st), `STSetFilterDamping()`, `STGetFilterDamping()`
334: E*/
335: typedef enum { ST_FILTER_DAMPING_NONE,
336:                ST_FILTER_DAMPING_JACKSON,
337:                ST_FILTER_DAMPING_LANCZOS,
338:                ST_FILTER_DAMPING_FEJER } STFilterDamping;
339: SLEPC_EXTERN const char *STFilterDampings[];

341: /*MC
342:    ST_MATMODE_COPY - The coefficient matrix of the linear system, $A-\sigma B$, is
343:    built explicitly on a copy of $A$.

345:    Note:
346:    If memory is an issue, one may prefer one of the other two modes.

348:    Level: intermediate

350: .seealso: [](ch:st), `STMatMode`, `STSetMatMode()`, `ST_MATMODE_INPLACE`, `ST_MATMODE_SHELL`
351: M*/

353: SLEPC_EXTERN PetscErrorCode STFilterSetType(ST,STFilterType);
354: SLEPC_EXTERN PetscErrorCode STFilterGetType(ST,STFilterType*);
355: SLEPC_EXTERN PetscErrorCode STFilterSetInterval(ST,PetscReal,PetscReal);
356: SLEPC_EXTERN PetscErrorCode STFilterGetInterval(ST,PetscReal*,PetscReal*);
357: SLEPC_EXTERN PetscErrorCode STFilterSetRange(ST,PetscReal,PetscReal);
358: SLEPC_EXTERN PetscErrorCode STFilterGetRange(ST,PetscReal*,PetscReal*);
359: SLEPC_EXTERN PetscErrorCode STFilterSetDegree(ST,PetscInt);
360: SLEPC_EXTERN PetscErrorCode STFilterGetDegree(ST,PetscInt*);
361: SLEPC_EXTERN PetscErrorCode STFilterGetThreshold(ST,PetscReal*);
362: SLEPC_EXTERN PetscErrorCode STFilterSetDamping(ST,STFilterDamping);
363: SLEPC_EXTERN PetscErrorCode STFilterGetDamping(ST,STFilterDamping*);