Actual source code: slepclme.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: User interface for the SLEPc object for solving linear matrix equations
12: */
14: #pragma once
16: #include <slepcbv.h>
18: /* SUBMANSEC = LME */
20: SLEPC_EXTERN PetscErrorCode LMEInitializePackage(void);
21: SLEPC_EXTERN PetscErrorCode LMEFinalizePackage(void);
23: /*S
24: LME - SLEPc object that encapsulates functionality for linear matrix equations.
26: Level: beginner
28: .seealso: [](ch:lme), `LMECreate()`
29: S*/
30: typedef struct _p_LME* LME;
32: /*J
33: LMEType - String with the name of a method for solving linear matrix equations.
35: Level: beginner
37: .seealso: [](ch:lme), `LMESetType()`, `LME`
38: J*/
39: typedef const char *LMEType;
40: #define LMEKRYLOV "krylov"
42: /* Logging support */
43: SLEPC_EXTERN PetscClassId LME_CLASSID;
45: /*E
46: LMEProblemType - Determines the type of linear matrix equation.
48: Values:
49: + `LME_LYAPUNOV` - continuous-time Lyapunov equation
50: . `LME_SYLVESTER` - Sylvester equation
51: . `LME_GEN_LYAPUNOV` - generalized Lyapunov equation
52: . `LME_GEN_SYLVESTER` - generalized Sylvester equation
53: . `LME_DT_LYAPUNOV` - discrete-time Lyapunov equation
54: - `LME_STEIN` - Stein equation
56: Level: beginner
58: .seealso: [](ch:lme), `LMESetProblemType()`, `LMEGetProblemType()`
59: E*/
60: typedef enum { LME_LYAPUNOV,
61: LME_SYLVESTER,
62: LME_GEN_LYAPUNOV,
63: LME_GEN_SYLVESTER,
64: LME_DT_LYAPUNOV ,
65: LME_STEIN} LMEProblemType;
66: SLEPC_EXTERN const char *LMEProblemTypes[];
68: /*MC
69: LME_LYAPUNOV - The continuous-time Lyapunov equation, $AX+XA^*=-C$.
71: Level: beginner
73: .seealso: [](ch:lme), `LMEProblemType`, `LMESetProblemType()`, `LME_SYLVESTER`, `LME_GEN_LYAPUNOV`, `LME_GEN_SYLVESTER`, `LME_DT_LYAPUNOV`, `LME_STEIN`
74: M*/
76: /*MC
77: LME_SYLVESTER - The Sylvester equation, $AX+XB=C$.
79: Level: beginner
81: .seealso: [](ch:lme), `LMEProblemType`, `LMESetProblemType()`, `LME_LYAPUNOV`, `LME_GEN_LYAPUNOV`, `LME_GEN_SYLVESTER`, `LME_DT_LYAPUNOV`, `LME_STEIN`
82: M*/
84: /*MC
85: LME_GEN_LYAPUNOV - The generalized Lyapunov equation, $AXD^*+DXA^*=-C$.
87: Level: beginner
89: .seealso: [](ch:lme), `LMEProblemType`, `LMESetProblemType()`, `LME_LYAPUNOV`, `LME_GEN_SYLVESTER`, `LME_SYLVESTER`, `LME_DT_LYAPUNOV`, `LME_STEIN`
90: M*/
92: /*MC
93: LME_GEN_SYLVESTER - The generalized Sylvester equation, $AXE+DXB=C$.
95: Level: beginner
97: .seealso: [](ch:lme), `LMEProblemType`, `LMESetProblemType()`, `LME_LYAPUNOV`, `LME_DT_LYAPUNOV`, `LME_SYLVESTER`, `LME_GEN_LYAPUNOV`, `LME_STEIN`
98: M*/
100: /*MC
101: LME_DT_LYAPUNOV - The discrete-time Lyapunov equation, $AXA^*-X=-C$.
103: Level: beginner
105: .seealso: [](ch:lme), `LMEProblemType`, `LMESetProblemType()`, `LME_LYAPUNOV`, `LME_SYLVESTER`, `LME_GEN_LYAPUNOV`, `LME_GEN_SYLVESTER`, `LME_STEIN`
106: M*/
108: /*MC
109: LME_STEIN - The Stein equation, $AXE-X=-C$.
111: Level: beginner
113: .seealso: [](ch:lme), `LMEProblemType`, `LMESetProblemType()`, `LME_LYAPUNOV`, `LME_SYLVESTER`, `LME_GEN_LYAPUNOV`, `LME_GEN_SYLVESTER`, `LME_DT_LYAPUNOV`
114: M*/
116: SLEPC_EXTERN PetscErrorCode LMECreate(MPI_Comm,LME*);
117: SLEPC_EXTERN PetscErrorCode LMEDestroy(LME*);
118: SLEPC_EXTERN PetscErrorCode LMEReset(LME);
119: SLEPC_EXTERN PetscErrorCode LMESetType(LME,LMEType);
120: SLEPC_EXTERN PetscErrorCode LMEGetType(LME,LMEType*);
121: SLEPC_EXTERN PetscErrorCode LMESetProblemType(LME,LMEProblemType);
122: SLEPC_EXTERN PetscErrorCode LMEGetProblemType(LME,LMEProblemType*);
123: SLEPC_EXTERN PetscErrorCode LMESetCoefficients(LME,Mat,Mat,Mat,Mat);
124: SLEPC_EXTERN PetscErrorCode LMEGetCoefficients(LME,Mat*,Mat*,Mat*,Mat*);
125: SLEPC_EXTERN PetscErrorCode LMESetRHS(LME,Mat);
126: SLEPC_EXTERN PetscErrorCode LMEGetRHS(LME,Mat*);
127: SLEPC_EXTERN PetscErrorCode LMESetSolution(LME,Mat);
128: SLEPC_EXTERN PetscErrorCode LMEGetSolution(LME,Mat*);
129: SLEPC_EXTERN PetscErrorCode LMESetFromOptions(LME);
130: SLEPC_EXTERN PetscErrorCode LMESetUp(LME);
131: SLEPC_EXTERN PetscErrorCode LMESolve(LME);
132: SLEPC_EXTERN PetscErrorCode LMEView(LME,PetscViewer);
133: SLEPC_EXTERN PetscErrorCode LMEViewFromOptions(LME,PetscObject,const char[]);
134: SLEPC_EXTERN PetscErrorCode LMEConvergedReasonView(LME,PetscViewer);
135: SLEPC_EXTERN PetscErrorCode LMEConvergedReasonViewFromOptions(LME);
136: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "LMEConvergedReasonView()", ) static inline PetscErrorCode LMEReasonView(LME lme,PetscViewer v) {return LMEConvergedReasonView(lme,v);}
137: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "LMEConvergedReasonViewFromOptions()", ) static inline PetscErrorCode LMEReasonViewFromOptions(LME lme) {return LMEConvergedReasonViewFromOptions(lme);}
139: SLEPC_EXTERN PetscErrorCode LMESetBV(LME,BV);
140: SLEPC_EXTERN PetscErrorCode LMEGetBV(LME,BV*);
141: SLEPC_EXTERN PetscErrorCode LMESetTolerances(LME,PetscReal,PetscInt);
142: SLEPC_EXTERN PetscErrorCode LMEGetTolerances(LME,PetscReal*,PetscInt*);
143: SLEPC_EXTERN PetscErrorCode LMESetDimensions(LME,PetscInt);
144: SLEPC_EXTERN PetscErrorCode LMEGetDimensions(LME,PetscInt*);
145: SLEPC_EXTERN PetscErrorCode LMEGetIterationNumber(LME,PetscInt*);
147: SLEPC_EXTERN PetscErrorCode LMEGetErrorEstimate(LME,PetscReal*);
148: SLEPC_EXTERN PetscErrorCode LMEComputeError(LME,PetscReal*);
149: SLEPC_EXTERN PetscErrorCode LMESetErrorIfNotConverged(LME,PetscBool);
150: SLEPC_EXTERN PetscErrorCode LMEGetErrorIfNotConverged(LME,PetscBool*);
152: SLEPC_EXTERN PetscErrorCode LMEDenseLyapunov(LME,PetscInt,PetscScalar[],PetscInt,PetscScalar[],PetscInt,PetscScalar[],PetscInt);
153: SLEPC_EXTERN PetscErrorCode LMEDenseHessLyapunovChol(LME,PetscInt,PetscScalar[],PetscInt,PetscInt,PetscScalar[],PetscInt,PetscScalar[],PetscInt,PetscReal*);
154: PETSC_DEPRECATED_FUNCTION(3, 8, 0, "LMEDenseHessLyapunovChol()", ) static inline PetscErrorCode LMEDenseLyapunovChol(LME lme,PetscScalar *H,PetscInt m,PetscInt ldh,PetscScalar *r,PetscScalar *L,PetscInt ldl,PetscReal *res) {return LMEDenseHessLyapunovChol(lme,m,H,ldh,1,r,m,L,ldl,res);}
156: /*S
157: LMEMonitorFn - A function prototype for functions provided to `LMEMonitorSet()`.
159: Calling Sequence:
160: + lme - the linear matrix equation solver context
161: . its - iteration number
162: . errest - error estimate
163: - ctx - optional monitoring context, as provided with `LMEMonitorSet()`
165: Level: intermediate
167: .seealso: [](ch:lme), `LMEMonitorSet()`
168: S*/
169: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode LMEMonitorFn(LME lme,PetscInt its,PetscReal errest,void *ctx);
171: /*S
172: LMEMonitorRegisterFn - A function prototype for functions provided to `LMEMonitorRegister()`.
174: Calling Sequence:
175: + lme - the linear matrix equation solver context
176: . its - iteration number
177: . errest - error estimate
178: - ctx - `PetscViewerAndFormat` object
180: Level: advanced
182: Note:
183: This is an `LMEMonitorFn` specialized for a context of `PetscViewerAndFormat`.
185: .seealso: [](ch:lme), `LMEMonitorSet()`, `LMEMonitorRegister()`, `LMEMonitorFn`, `LMEMonitorRegisterCreateFn`, `LMEMonitorRegisterDestroyFn`
186: S*/
187: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode LMEMonitorRegisterFn(LME lme,PetscInt its,PetscReal errest,PetscViewerAndFormat *ctx);
189: /*S
190: LMEMonitorRegisterCreateFn - A function prototype for functions that do the
191: creation when provided to `LMEMonitorRegister()`.
193: Calling Sequence:
194: + viewer - the viewer to be used with the `LMEMonitorRegisterFn`
195: . format - the format of the viewer
196: . ctx - a context for the monitor
197: - result - a `PetscViewerAndFormat` object
199: Level: advanced
201: .seealso: [](ch:lme), `LMEMonitorRegisterFn`, `LMEMonitorSet()`, `LMEMonitorRegister()`, `LMEMonitorFn`, `LMEMonitorRegisterDestroyFn`
202: S*/
203: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode LMEMonitorRegisterCreateFn(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **result);
205: /*S
206: LMEMonitorRegisterDestroyFn - A function prototype for functions that do the after
207: use destruction when provided to `LMEMonitorRegister()`.
209: Calling Sequence:
210: . vf - a `PetscViewerAndFormat` object to be destroyed, including any context
212: Level: advanced
214: .seealso: [](ch:lme), `LMEMonitorRegisterFn`, `LMEMonitorSet()`, `LMEMonitorRegister()`, `LMEMonitorFn`, `LMEMonitorRegisterCreateFn`
215: S*/
216: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode LMEMonitorRegisterDestroyFn(PetscViewerAndFormat **result);
218: SLEPC_EXTERN PetscErrorCode LMEMonitor(LME,PetscInt,PetscReal);
219: SLEPC_EXTERN PetscErrorCode LMEMonitorSet(LME,LMEMonitorFn,void*,PetscCtxDestroyFn*);
220: SLEPC_EXTERN PetscErrorCode LMEMonitorCancel(LME);
221: SLEPC_EXTERN PetscErrorCode LMEGetMonitorContext(LME,void*);
223: SLEPC_EXTERN PetscErrorCode LMEMonitorSetFromOptions(LME,const char[],const char[],void*);
224: SLEPC_EXTERN LMEMonitorRegisterFn LMEMonitorDefault;
225: SLEPC_EXTERN LMEMonitorRegisterFn LMEMonitorDefaultDrawLG;
226: SLEPC_EXTERN LMEMonitorRegisterCreateFn LMEMonitorDefaultDrawLGCreate;
228: SLEPC_EXTERN PetscErrorCode LMESetOptionsPrefix(LME,const char[]);
229: SLEPC_EXTERN PetscErrorCode LMEAppendOptionsPrefix(LME,const char[]);
230: SLEPC_EXTERN PetscErrorCode LMEGetOptionsPrefix(LME,const char*[]);
232: /*E
233: LMEConvergedReason - Reason a matrix function iteration was determined to
234: have converged or diverged.
236: Values:
237: + `LME_CONVERGED_TOL` - requested decrease in the residual
238: . `LME_DIVERGED_ITS` - requested number of iterations
239: . `LME_DIVERGED_BREAKDOWN` - breakdown in the solver
240: - `LME_CONVERGED_ITERATING` - the solver is still running
242: Level: intermediate
244: .seealso: [](ch:lme), `LMESolve()`, `LMEGetConvergedReason()`, `LMESetTolerances()`
245: E*/
246: typedef enum {/* converged */
247: LME_CONVERGED_TOL = 1,
248: /* diverged */
249: LME_DIVERGED_ITS = -1,
250: LME_DIVERGED_BREAKDOWN = -2,
251: LME_CONVERGED_ITERATING = 0} LMEConvergedReason;
252: SLEPC_EXTERN const char *const*LMEConvergedReasons;
254: /*MC
255: LME_CONVERGED_TOL - The computed residual of the matrix equation is below the
256: tolerance.
258: Level: intermediate
260: .seealso: [](ch:lme), `LMESolve()`, `LMEGetConvergedReason()`, `LMEConvergedReason`
261: M*/
263: /*MC
264: LME_DIVERGED_ITS - Ran out of iterations before the convergence criterion was
265: reached.
267: Level: intermediate
269: .seealso: [](ch:lme), `LMESolve()`, `LMEGetConvergedReason()`, `LMEConvergedReason`
270: M*/
272: /*MC
273: LME_DIVERGED_BREAKDOWN - A breakdown in the solver was detected so the
274: method could not continue.
276: Level: intermediate
278: .seealso: [](ch:lme), `LMESolve()`, `LMEGetConvergedReason()`, `LMEConvergedReason`
279: M*/
281: /*MC
282: LME_CONVERGED_ITERATING - This flag is returned if `LMEGetConvergedReason()` is called
283: while `LMESolve()` is still running..
285: Level: intermediate
287: .seealso: [](ch:lme), `LMESolve()`, `LMEGetConvergedReason()`, `LMEConvergedReason`
288: M*/
290: SLEPC_EXTERN PetscErrorCode LMEGetConvergedReason(LME,LMEConvergedReason*);
292: SLEPC_EXTERN PetscFunctionList LMEList;
293: SLEPC_EXTERN PetscFunctionList LMEMonitorList;
294: SLEPC_EXTERN PetscFunctionList LMEMonitorCreateList;
295: SLEPC_EXTERN PetscFunctionList LMEMonitorDestroyList;
296: SLEPC_EXTERN PetscErrorCode LMERegister(const char[],PetscErrorCode(*)(LME));
297: SLEPC_EXTERN PetscErrorCode LMEMonitorRegister(const char[],PetscViewerType,PetscViewerFormat,LMEMonitorRegisterFn*,LMEMonitorRegisterCreateFn*,LMEMonitorRegisterDestroyFn*);
299: SLEPC_EXTERN PetscErrorCode LMEAllocateSolution(LME,PetscInt);