Actual source code: slepcst.h
slepc-3.10.0 2018-09-18
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2018, 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: */
16: #include <slepcsys.h>
17: #include <slepcbv.h>
18: #include <petscksp.h>
20: PETSC_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: PETSC_EXTERN PetscClassId ST_CLASSID;
50: PETSC_EXTERN PetscErrorCode STCreate(MPI_Comm,ST*);
51: PETSC_EXTERN PetscErrorCode STDestroy(ST*);
52: PETSC_EXTERN PetscErrorCode STReset(ST);
53: PETSC_EXTERN PetscErrorCode STSetType(ST,STType);
54: PETSC_EXTERN PetscErrorCode STGetType(ST,STType*);
55: PETSC_EXTERN PetscErrorCode STSetMatrices(ST,PetscInt,Mat*);
56: PETSC_EXTERN PetscErrorCode STGetMatrix(ST,PetscInt,Mat*);
57: PETSC_EXTERN PetscErrorCode STGetMatrixTransformed(ST,PetscInt,Mat*);
58: PETSC_EXTERN PetscErrorCode STGetNumMatrices(ST,PetscInt*);
59: PETSC_EXTERN PetscErrorCode STGetOperator(ST,Mat*);
60: PETSC_EXTERN PetscErrorCode STSetUp(ST);
61: PETSC_EXTERN PetscErrorCode STSetFromOptions(ST);
62: PETSC_EXTERN PetscErrorCode STView(ST,PetscViewer);
64: PETSC_DEPRECATED("Use STSetMatrices()") PETSC_STATIC_INLINE PetscErrorCode STSetOperators(ST st,PetscInt n,Mat *A) {return STSetMatrices(st,n,A);}
65: PETSC_DEPRECATED("Use STGetMatrix()") PETSC_STATIC_INLINE PetscErrorCode STGetOperators(ST st,PetscInt k,Mat *A) {return STGetMatrix(st,k,A);}
66: PETSC_DEPRECATED("Use STGetMatrixTransformed()") PETSC_STATIC_INLINE PetscErrorCode STGetTOperators(ST st,PetscInt k,Mat *A) {return STGetMatrixTransformed(st,k,A);}
67: PETSC_DEPRECATED("Use STGetOperator() followed by MatComputeExplicitOperator()") PETSC_STATIC_INLINE PetscErrorCode STComputeExplicitOperator(ST st,Mat *A) {
69: STGetOperator(st,&Op);
70: MatComputeExplicitOperator(Op,A);
71: MatDestroy(&Op);
72: return(0);
73: }
75: PETSC_EXTERN PetscErrorCode STApply(ST,Vec,Vec);
76: PETSC_EXTERN PetscErrorCode STMatMult(ST,PetscInt,Vec,Vec);
77: PETSC_EXTERN PetscErrorCode STMatMultTranspose(ST,PetscInt,Vec,Vec);
78: PETSC_EXTERN PetscErrorCode STMatSolve(ST,Vec,Vec);
79: PETSC_EXTERN PetscErrorCode STMatSolveTranspose(ST,Vec,Vec);
80: PETSC_EXTERN PetscErrorCode STGetBilinearForm(ST,Mat*);
81: PETSC_EXTERN PetscErrorCode STApplyTranspose(ST,Vec,Vec);
82: PETSC_EXTERN PetscErrorCode STMatSetUp(ST,PetscScalar,PetscScalar*);
83: PETSC_EXTERN PetscErrorCode STPostSolve(ST);
84: PETSC_EXTERN PetscErrorCode STResetMatrixState(ST);
85: PETSC_EXTERN PetscErrorCode STSetWorkVecs(ST,PetscInt);
87: PETSC_EXTERN PetscErrorCode STSetKSP(ST,KSP);
88: PETSC_EXTERN PetscErrorCode STGetKSP(ST,KSP*);
89: PETSC_EXTERN PetscErrorCode STSetShift(ST,PetscScalar);
90: PETSC_EXTERN PetscErrorCode STGetShift(ST,PetscScalar*);
91: PETSC_EXTERN PetscErrorCode STSetDefaultShift(ST,PetscScalar);
92: PETSC_EXTERN PetscErrorCode STScaleShift(ST,PetscScalar);
93: PETSC_EXTERN PetscErrorCode STSetBalanceMatrix(ST,Vec);
94: PETSC_EXTERN PetscErrorCode STGetBalanceMatrix(ST,Vec*);
95: PETSC_EXTERN PetscErrorCode STSetTransform(ST,PetscBool);
96: PETSC_EXTERN PetscErrorCode STGetTransform(ST,PetscBool*);
98: PETSC_EXTERN PetscErrorCode STSetOptionsPrefix(ST,const char*);
99: PETSC_EXTERN PetscErrorCode STAppendOptionsPrefix(ST,const char*);
100: PETSC_EXTERN PetscErrorCode STGetOptionsPrefix(ST,const char*[]);
102: PETSC_EXTERN PetscErrorCode STBackTransform(ST,PetscInt,PetscScalar*,PetscScalar*);
103: PETSC_EXTERN PetscErrorCode STIsInjective(ST,PetscBool*);
105: PETSC_EXTERN PetscErrorCode STCheckNullSpace(ST,BV);
107: PETSC_EXTERN PetscErrorCode STMatCreateVecs(ST,Vec*,Vec*);
108: PETSC_EXTERN PetscErrorCode STMatCreateVecsEmpty(ST,Vec*,Vec*);
109: PETSC_EXTERN PetscErrorCode STMatGetSize(ST,PetscInt*,PetscInt*);
110: PETSC_EXTERN PetscErrorCode STMatGetLocalSize(ST,PetscInt*,PetscInt*);
112: /*E
113: STMatMode - Determines how to handle the coefficient matrix associated
114: to the spectral transformation
116: Level: intermediate
118: .seealso: STSetMatMode(), STGetMatMode()
119: E*/
120: typedef enum { ST_MATMODE_COPY,
121: ST_MATMODE_INPLACE,
122: ST_MATMODE_SHELL } STMatMode;
123: PETSC_EXTERN const char *STMatModes[];
125: PETSC_EXTERN PetscErrorCode STSetMatMode(ST,STMatMode);
126: PETSC_EXTERN PetscErrorCode STGetMatMode(ST,STMatMode*);
127: PETSC_EXTERN PetscErrorCode STSetMatStructure(ST,MatStructure);
128: PETSC_EXTERN PetscErrorCode STGetMatStructure(ST,MatStructure*);
130: PETSC_EXTERN PetscFunctionList STList;
131: PETSC_EXTERN PetscErrorCode STRegister(const char[],PetscErrorCode(*)(ST));
133: /* --------- options specific to particular spectral transformations-------- */
135: PETSC_EXTERN PetscErrorCode STShellGetContext(ST st,void **ctx);
136: PETSC_EXTERN PetscErrorCode STShellSetContext(ST st,void *ctx);
137: PETSC_EXTERN PetscErrorCode STShellSetApply(ST st,PetscErrorCode (*apply)(ST,Vec,Vec));
138: PETSC_EXTERN PetscErrorCode STShellSetApplyTranspose(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec));
139: PETSC_EXTERN PetscErrorCode STShellSetBackTransform(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*));
141: PETSC_EXTERN PetscErrorCode STCayleyGetAntishift(ST,PetscScalar*);
142: PETSC_EXTERN PetscErrorCode STCayleySetAntishift(ST,PetscScalar);
144: PETSC_EXTERN PetscErrorCode STPrecondGetMatForPC(ST,Mat*);
145: PETSC_EXTERN PetscErrorCode STPrecondSetMatForPC(ST,Mat);
146: PETSC_EXTERN PetscErrorCode STPrecondGetKSPHasMat(ST,PetscBool*);
147: PETSC_EXTERN PetscErrorCode STPrecondSetKSPHasMat(ST,PetscBool);
149: PETSC_EXTERN PetscErrorCode STFilterSetInterval(ST,PetscReal,PetscReal);
150: PETSC_EXTERN PetscErrorCode STFilterGetInterval(ST,PetscReal*,PetscReal*);
151: PETSC_EXTERN PetscErrorCode STFilterSetRange(ST,PetscReal,PetscReal);
152: PETSC_EXTERN PetscErrorCode STFilterGetRange(ST,PetscReal*,PetscReal*);
153: PETSC_EXTERN PetscErrorCode STFilterSetDegree(ST,PetscInt);
154: PETSC_EXTERN PetscErrorCode STFilterGetDegree(ST,PetscInt*);
155: PETSC_EXTERN PetscErrorCode STFilterGetThreshold(ST,PetscReal*);
157: #endif