Actual source code: slepcst.h
1: /*
2: Spectral transformation module for eigenvalue problems.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2011, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
9:
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
26: #include "petscksp.h"
27: PETSC_EXTERN_CXX_BEGIN
29: extern PetscErrorCode STInitializePackage(const char[]);
31: /*S
32: ST - Abstract SLEPc object that manages spectral transformations.
33: This object is accessed only in advanced applications.
35: Level: beginner
37: .seealso: STCreate(), EPS
38: S*/
39: typedef struct _p_ST* ST;
41: /*E
42: STType - String with the name of a SLEPc spectral transformation
44: Level: beginner
46: .seealso: STSetType(), ST
47: E*/
48: #define STType char*
49: #define STSHELL "shell"
50: #define STSHIFT "shift"
51: #define STSINVERT "sinvert"
52: #define STCAYLEY "cayley"
53: #define STFOLD "fold"
54: #define STPRECOND "precond"
56: /* Logging support */
57: extern PetscClassId ST_CLASSID;
59: extern PetscErrorCode STCreate(MPI_Comm,ST*);
60: extern PetscErrorCode STDestroy(ST*);
61: extern PetscErrorCode STReset(ST);
62: extern PetscErrorCode STSetType(ST,const STType);
63: extern PetscErrorCode STGetType(ST,const STType*);
64: extern PetscErrorCode STSetOperators(ST,Mat,Mat);
65: extern PetscErrorCode STGetOperators(ST,Mat*,Mat*);
66: extern PetscErrorCode STSetUp(ST);
67: extern PetscErrorCode STSetFromOptions(ST);
68: extern PetscErrorCode STView(ST,PetscViewer);
70: extern PetscErrorCode STApply(ST,Vec,Vec);
71: extern PetscErrorCode STGetBilinearForm(ST,Mat*);
72: extern PetscErrorCode STApplyTranspose(ST,Vec,Vec);
73: extern PetscErrorCode STComputeExplicitOperator(ST,Mat*);
74: extern PetscErrorCode STPostSolve(ST);
76: extern PetscErrorCode STSetKSP(ST,KSP);
77: extern PetscErrorCode STGetKSP(ST,KSP*);
78: extern PetscErrorCode STSetShift(ST,PetscScalar);
79: extern PetscErrorCode STGetShift(ST,PetscScalar*);
80: extern PetscErrorCode STSetDefaultShift(ST,PetscScalar);
81: extern PetscErrorCode STSetBalanceMatrix(ST,Vec);
82: extern PetscErrorCode STGetBalanceMatrix(ST,Vec*);
84: extern PetscErrorCode STSetOptionsPrefix(ST,const char*);
85: extern PetscErrorCode STAppendOptionsPrefix(ST,const char*);
86: extern PetscErrorCode STGetOptionsPrefix(ST,const char*[]);
88: extern PetscErrorCode STBackTransform(ST,PetscInt,PetscScalar*,PetscScalar*);
90: extern PetscErrorCode STCheckNullSpace(ST,PetscInt,const Vec[]);
92: extern PetscErrorCode STGetOperationCounters(ST,PetscInt*,PetscInt*);
93: extern PetscErrorCode STResetOperationCounters(ST);
95: /*E
96: STMatMode - determines how to handle the coefficient matrix associated
97: to the spectral transformation
99: Level: intermediate
101: .seealso: STSetMatMode(), STGetMatMode()
102: E*/
103: typedef enum { ST_MATMODE_COPY,
104: ST_MATMODE_INPLACE,
105: ST_MATMODE_SHELL } STMatMode;
106: extern PetscErrorCode STSetMatMode(ST,STMatMode);
107: extern PetscErrorCode STGetMatMode(ST,STMatMode*);
108: extern PetscErrorCode STSetMatStructure(ST,MatStructure);
109: extern PetscErrorCode STGetMatStructure(ST,MatStructure*);
111: extern PetscFList STList;
112: extern PetscBool STRegisterAllCalled;
113: extern PetscErrorCode STRegisterAll(const char[]);
114: extern PetscErrorCode STRegisterDestroy(void);
115: extern PetscErrorCode STRegister(const char[],const char[],const char[],PetscErrorCode(*)(ST));
117: /*MC
118: STRegisterDynamic - Adds a method to the spectral transformation package.
120: Synopsis:
121: PetscErrorCode STRegisterDynamic(const char *name,const char *path,const char *name_create,PetscErrorCode (*routine_create)(ST))
123: Not collective
125: Input Parameters:
126: + name - name of a new user-defined transformation
127: . path - path (either absolute or relative) the library containing this solver
128: . name_create - name of routine to create method context
129: - routine_create - routine to create method context
131: Notes:
132: STRegisterDynamic() may be called multiple times to add several user-defined
133: spectral transformations.
135: If dynamic libraries are used, then the fourth input argument (routine_create)
136: is ignored.
138: Sample usage:
139: .vb
140: STRegisterDynamic("my_solver","/home/username/my_lib/lib/libO/solaris/mylib.a",
141: "MySolverCreate",MySolverCreate);
142: .ve
144: Then, your solver can be chosen with the procedural interface via
145: $ STSetType(st,"my_solver")
146: or at runtime via the option
147: $ -st_type my_solver
149: Level: advanced
151: .seealso: STRegisterDestroy(), STRegisterAll()
152: M*/
153: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
154: #define STRegisterDynamic(a,b,c,d) STRegister(a,b,c,0)
155: #else
156: #define STRegisterDynamic(a,b,c,d) STRegister(a,b,c,d)
157: #endif
159: /* --------- options specific to particular spectral transformations-------- */
161: extern PetscErrorCode STShellGetContext(ST st,void **ctx);
162: extern PetscErrorCode STShellSetContext(ST st,void *ctx);
163: extern PetscErrorCode STShellSetApply(ST st,PetscErrorCode (*apply)(ST,Vec,Vec));
164: extern PetscErrorCode STShellSetApplyTranspose(ST st,PetscErrorCode (*applytrans)(ST,Vec,Vec));
165: extern PetscErrorCode STShellSetBackTransform(ST st,PetscErrorCode (*backtr)(ST,PetscInt,PetscScalar*,PetscScalar*));
167: extern PetscErrorCode STCayleyGetAntishift(ST,PetscScalar*);
168: extern PetscErrorCode STCayleySetAntishift(ST,PetscScalar);
170: extern PetscErrorCode STPrecondGetMatForPC(ST st,Mat *mat);
171: extern PetscErrorCode STPrecondSetMatForPC(ST st,Mat mat);
172: extern PetscErrorCode STPrecondGetKSPHasMat(ST st,PetscBool *setmat);
173: extern PetscErrorCode STPrecondSetKSPHasMat(ST st,PetscBool setmat);
175: PETSC_EXTERN_CXX_END
176: #endif