Actual source code: slepcst.h
slepc-3.16.0 2021-09-30
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: Spectral transformation module for eigenvalue problems
12: */
14: #if !defined(SLEPCST_H)
15: #define SLEPCST_H
16: #include <slepcsys.h>
17: #include <slepcbv.h>
18: #include <petscksp.h>
20: SLEPC_EXTERN PetscErrorCode STInitializePackage(void);
22: /*S
23: ST - Abstract SLEPc object that manages spectral transformations.
24: This object is accessed only in advanced applications.
26: Level: beginner
28: .seealso: STCreate(), EPS
29: S*/
30: typedef struct _p_ST* ST;
32: /*J
33: STType - String with the name of a SLEPc spectral transformation
35: Level: beginner
37: .seealso: STSetType(), ST
38: J*/
39: typedef const char* STType;
40: #define STSHELL "shell"
41: #define STSHIFT "shift"
42: #define STSINVERT "sinvert"
43: #define STCAYLEY "cayley"
44: #define STPRECOND "precond"
45: #define STFILTER "filter"
47: /* Logging support */
48: SLEPC_EXTERN PetscClassId ST_CLASSID;
50: SLEPC_EXTERN PetscErrorCode STCreate(MPI_Comm,ST*);
51: SLEPC_EXTERN PetscErrorCode STDestroy(ST*);
52: SLEPC_EXTERN PetscErrorCode STReset(ST);
53: SLEPC_EXTERN PetscErrorCode STSetType(ST,STType);
54: SLEPC_EXTERN PetscErrorCode STGetType(ST,STType*);
55: SLEPC_EXTERN PetscErrorCode STSetMatrices(ST,PetscInt,Mat*);
56: SLEPC_EXTERN PetscErrorCode STGetMatrix(ST,PetscInt,Mat*);
57: SLEPC_EXTERN PetscErrorCode STGetMatrixTransformed(ST,PetscInt,Mat*);
58: SLEPC_EXTERN PetscErrorCode STGetNumMatrices(ST,PetscInt*);
59: SLEPC_EXTERN PetscErrorCode STGetOperator(ST,Mat*);
60: SLEPC_EXTERN PetscErrorCode STRestoreOperator(ST,Mat*);
61: SLEPC_EXTERN PetscErrorCode STSetUp(ST);
62: SLEPC_EXTERN PetscErrorCode STSetFromOptions(ST);
63: SLEPC_EXTERN PetscErrorCode STView(ST,PetscViewer);
64: SLEPC_EXTERN PetscErrorCode STViewFromOptions(ST,PetscObject,const char[]);
66: PETSC_DEPRECATED_FUNCTION("Use STSetMatrices()") PETSC_STATIC_INLINE PetscErrorCode STSetOperators(ST st,PetscInt n,Mat *A) {return STSetMatrices(st,n,A);}
67: PETSC_DEPRECATED_FUNCTION("Use STGetMatrix()") PETSC_STATIC_INLINE PetscErrorCode STGetOperators(ST st,PetscInt k,Mat *A) {return STGetMatrix(st,k,A);}
68: PETSC_DEPRECATED_FUNCTION("Use STGetMatrixTransformed()") PETSC_STATIC_INLINE PetscErrorCode STGetTOperators(ST st,PetscInt k,Mat *A) {return STGetMatrixTransformed(st,k,A);}
69: PETSC_DEPRECATED_FUNCTION("Use STGetOperator() followed by MatComputeOperator()") PETSC_STATIC_INLINE PetscErrorCode STComputeExplicitOperator(ST st,Mat *A)
70: {
72: STGetOperator(st,&Op);
73: MatComputeOperator(Op,MATAIJ,A);
74: STRestoreOperator(st,&Op);
75: return(0);
76: }
78: SLEPC_EXTERN PetscErrorCode STApply(ST,Vec,Vec);
79: SLEPC_EXTERN PetscErrorCode STApplyMat(ST,Mat,Mat);
80: SLEPC_EXTERN PetscErrorCode STApplyTranspose(ST,Vec,Vec);
81: SLEPC_EXTERN PetscErrorCode STApplyHermitianTranspose(ST,Vec,Vec);
82: SLEPC_EXTERN PetscErrorCode STMatMult(ST,PetscInt,Vec,Vec);
83: SLEPC_EXTERN PetscErrorCode STMatMultTranspose(ST,PetscInt,Vec,Vec);
84: SLEPC_EXTERN PetscErrorCode STMatSolve(ST,Vec,Vec);
85: SLEPC_EXTERN PetscErrorCode STMatSolveTranspose(ST,Vec,Vec);
86: SLEPC_EXTERN PetscErrorCode STMatMatSolve(ST,Mat,Mat);
87: SLEPC_EXTERN PetscErrorCode STGetBilinearForm(ST,Mat*);
88: SLEPC_EXTERN PetscErrorCode STMatSetUp(ST,PetscScalar,PetscScalar*);
89: SLEPC_EXTERN PetscErrorCode STPostSolve(ST);
90: SLEPC_EXTERN PetscErrorCode STResetMatrixState(ST);
91: SLEPC_EXTERN PetscErrorCode STSetWorkVecs(ST,PetscInt);
93: SLEPC_EXTERN PetscErrorCode STSetKSP(ST,KSP);
94: SLEPC_EXTERN PetscErrorCode STGetKSP(ST,KSP*);
95: SLEPC_EXTERN PetscErrorCode STSetShift(ST,PetscScalar);
96: SLEPC_EXTERN PetscErrorCode STGetShift(ST,PetscScalar*);
97: SLEPC_EXTERN PetscErrorCode STSetDefaultShift(ST,PetscScalar);
98: SLEPC_EXTERN PetscErrorCode STScaleShift(ST,PetscScalar);
99: SLEPC_EXTERN PetscErrorCode STSetBalanceMatrix(ST,Vec);
100: SLEPC_EXTERN PetscErrorCode STGetBalanceMatrix(ST,Vec*);
101: SLEPC_EXTERN PetscErrorCode STSetTransform(ST,PetscBool);
102: SLEPC_EXTERN PetscErrorCode STGetTransform(ST,PetscBool*);
104: SLEPC_EXTERN PetscErrorCode STSetOptionsPrefix(ST,const char*);
105: SLEPC_EXTERN PetscErrorCode STAppendOptionsPrefix(ST,const char*);
106: SLEPC_EXTERN PetscErrorCode STGetOptionsPrefix(ST,const char*[]);
108: SLEPC_EXTERN PetscErrorCode STBackTransform(ST,PetscInt,PetscScalar*,PetscScalar*);
109: SLEPC_EXTERN PetscErrorCode STIsInjective(ST,PetscBool*);
111: SLEPC_EXTERN PetscErrorCode STCheckNullSpace(ST,BV);
113: SLEPC_EXTERN PetscErrorCode STSetPreconditionerMat(ST,Mat);
114: SLEPC_EXTERN PetscErrorCode STGetPreconditionerMat(ST,Mat*);
116: SLEPC_EXTERN PetscErrorCode STMatCreateVecs(ST,Vec*,Vec*);
117: SLEPC_EXTERN PetscErrorCode STMatCreateVecsEmpty(ST,Vec*,Vec*);
118: SLEPC_EXTERN PetscErrorCode STMatGetSize(ST,PetscInt*,PetscInt*);
119: SLEPC_EXTERN PetscErrorCode STMatGetLocalSize(ST,PetscInt*,PetscInt*);
121: /*E
122: STMatMode - Determines how to handle the coefficient matrix associated
123: to the spectral transformation
125: Level: intermediate
127: .seealso: STSetMatMode(), STGetMatMode()
128: E*/
129: typedef enum { ST_MATMODE_COPY,
130: ST_MATMODE_INPLACE,
131: ST_MATMODE_SHELL } STMatMode;
132: SLEPC_EXTERN const char *STMatModes[];
134: SLEPC_EXTERN PetscErrorCode STSetMatMode(ST,STMatMode);
135: SLEPC_EXTERN PetscErrorCode STGetMatMode(ST,STMatMode*);
136: SLEPC_EXTERN PetscErrorCode STSetMatStructure(ST,MatStructure);
137: SLEPC_EXTERN PetscErrorCode STGetMatStructure(ST,MatStructure*);
139: SLEPC_EXTERN PetscFunctionList STList;
140: SLEPC_EXTERN PetscErrorCode STRegister(const char[],PetscErrorCode(*)(ST));
142: /* --------- options specific to particular spectral transformations-------- */
144: SLEPC_EXTERN PetscErrorCode STShellGetContext(ST,void*);
145: SLEPC_EXTERN PetscErrorCode STShellSetContext(ST,void*);
146: SLEPC_EXTERN PetscErrorCode STShellSetApply(ST,PetscErrorCode (*)(ST,Vec,Vec));
147: SLEPC_EXTERN PetscErrorCode STShellSetApplyTranspose(ST,PetscErrorCode (*)(ST,Vec,Vec));
148: SLEPC_EXTERN PetscErrorCode STShellSetBackTransform(ST,PetscErrorCode (*)(ST,PetscInt,PetscScalar*,PetscScalar*));
150: SLEPC_EXTERN PetscErrorCode STCayleyGetAntishift(ST,PetscScalar*);
151: SLEPC_EXTERN PetscErrorCode STCayleySetAntishift(ST,PetscScalar);
153: PETSC_DEPRECATED_FUNCTION("Use STGetPreconditionerMat()") PETSC_STATIC_INLINE PetscErrorCode STPrecondGetMatForPC(ST st,Mat *A) {return STGetPreconditionerMat(st,A);}
154: PETSC_DEPRECATED_FUNCTION("Use STSetPreconditionerMat()") PETSC_STATIC_INLINE PetscErrorCode STPrecondSetMatForPC(ST st,Mat A) {return STSetPreconditionerMat(st,A);}
155: SLEPC_EXTERN PetscErrorCode STPrecondGetKSPHasMat(ST,PetscBool*);
156: SLEPC_EXTERN PetscErrorCode STPrecondSetKSPHasMat(ST,PetscBool);
158: SLEPC_EXTERN PetscErrorCode STFilterSetInterval(ST,PetscReal,PetscReal);
159: SLEPC_EXTERN PetscErrorCode STFilterGetInterval(ST,PetscReal*,PetscReal*);
160: SLEPC_EXTERN PetscErrorCode STFilterSetRange(ST,PetscReal,PetscReal);
161: SLEPC_EXTERN PetscErrorCode STFilterGetRange(ST,PetscReal*,PetscReal*);
162: SLEPC_EXTERN PetscErrorCode STFilterSetDegree(ST,PetscInt);
163: SLEPC_EXTERN PetscErrorCode STFilterGetDegree(ST,PetscInt*);
164: SLEPC_EXTERN PetscErrorCode STFilterGetThreshold(ST,PetscReal*);
166: #endif