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:  `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: `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:     Level: beginner

 50: .seealso: `LMESetProblemType()`, `LMEGetProblemType()`
 51: E*/
 52: typedef enum { LME_LYAPUNOV,
 53:                LME_SYLVESTER,
 54:                LME_GEN_LYAPUNOV,
 55:                LME_GEN_SYLVESTER,
 56:                LME_DT_LYAPUNOV ,
 57:                LME_STEIN} LMEProblemType;
 58: SLEPC_EXTERN const char *LMEProblemTypes[];

 60: SLEPC_EXTERN PetscErrorCode LMECreate(MPI_Comm,LME *);
 61: SLEPC_EXTERN PetscErrorCode LMEDestroy(LME*);
 62: SLEPC_EXTERN PetscErrorCode LMEReset(LME);
 63: SLEPC_EXTERN PetscErrorCode LMESetType(LME,LMEType);
 64: SLEPC_EXTERN PetscErrorCode LMEGetType(LME,LMEType*);
 65: SLEPC_EXTERN PetscErrorCode LMESetProblemType(LME,LMEProblemType);
 66: SLEPC_EXTERN PetscErrorCode LMEGetProblemType(LME,LMEProblemType*);
 67: SLEPC_EXTERN PetscErrorCode LMESetCoefficients(LME,Mat,Mat,Mat,Mat);
 68: SLEPC_EXTERN PetscErrorCode LMEGetCoefficients(LME,Mat*,Mat*,Mat*,Mat*);
 69: SLEPC_EXTERN PetscErrorCode LMESetRHS(LME,Mat);
 70: SLEPC_EXTERN PetscErrorCode LMEGetRHS(LME,Mat*);
 71: SLEPC_EXTERN PetscErrorCode LMESetSolution(LME,Mat);
 72: SLEPC_EXTERN PetscErrorCode LMEGetSolution(LME,Mat*);
 73: SLEPC_EXTERN PetscErrorCode LMESetFromOptions(LME);
 74: SLEPC_EXTERN PetscErrorCode LMESetUp(LME);
 75: SLEPC_EXTERN PetscErrorCode LMESolve(LME);
 76: SLEPC_EXTERN PetscErrorCode LMEView(LME,PetscViewer);
 77: SLEPC_EXTERN PetscErrorCode LMEViewFromOptions(LME,PetscObject,const char[]);
 78: SLEPC_EXTERN PetscErrorCode LMEConvergedReasonView(LME,PetscViewer);
 79: SLEPC_EXTERN PetscErrorCode LMEConvergedReasonViewFromOptions(LME);
 80: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "LMEConvergedReasonView()", ) static inline PetscErrorCode LMEReasonView(LME lme,PetscViewer v) {return LMEConvergedReasonView(lme,v);}
 81: PETSC_DEPRECATED_FUNCTION(3, 14, 0, "LMEConvergedReasonViewFromOptions()", ) static inline PetscErrorCode LMEReasonViewFromOptions(LME lme) {return LMEConvergedReasonViewFromOptions(lme);}

 83: SLEPC_EXTERN PetscErrorCode LMESetBV(LME,BV);
 84: SLEPC_EXTERN PetscErrorCode LMEGetBV(LME,BV*);
 85: SLEPC_EXTERN PetscErrorCode LMESetTolerances(LME,PetscReal,PetscInt);
 86: SLEPC_EXTERN PetscErrorCode LMEGetTolerances(LME,PetscReal*,PetscInt*);
 87: SLEPC_EXTERN PetscErrorCode LMESetDimensions(LME,PetscInt);
 88: SLEPC_EXTERN PetscErrorCode LMEGetDimensions(LME,PetscInt*);
 89: SLEPC_EXTERN PetscErrorCode LMEGetIterationNumber(LME,PetscInt*);

 91: SLEPC_EXTERN PetscErrorCode LMEGetErrorEstimate(LME,PetscReal*);
 92: SLEPC_EXTERN PetscErrorCode LMEComputeError(LME,PetscReal*);
 93: SLEPC_EXTERN PetscErrorCode LMESetErrorIfNotConverged(LME,PetscBool);
 94: SLEPC_EXTERN PetscErrorCode LMEGetErrorIfNotConverged(LME,PetscBool*);

 96: SLEPC_EXTERN PetscErrorCode LMEDenseLyapunov(LME,PetscInt,PetscScalar*,PetscInt,PetscScalar*,PetscInt,PetscScalar*,PetscInt);
 97: SLEPC_EXTERN PetscErrorCode LMEDenseHessLyapunovChol(LME,PetscInt,PetscScalar*,PetscInt,PetscInt,PetscScalar*,PetscInt,PetscScalar*,PetscInt,PetscReal*);
 98: 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);}

100: /*S
101:   LMEMonitorFn - A function prototype for functions provided to LMEMonitorSet()

103:   Calling Sequence:
104: +   lme    - linear matrix equation solver context obtained from LMECreate()
105: .   its    - iteration number
106: .   errest - error estimate
107: -   ctx    - optional monitoring context, as provided with LMEMonitorSet()

109:   Level: beginner

111: .seealso: `LMEMonitorSet()`
112: S*/
113: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode LMEMonitorFn(LME lme,PetscInt its,PetscReal errest,void *ctx);

115: /*S
116:   LMEMonitorRegisterFn - A function prototype for functions provided to LMEMonitorRegister()

118:   Calling Sequence:
119: +   lme    - linear matrix equation solver context obtained from LMECreate()
120: .   its    - iteration number
121: .   errest - error estimate
122: -   ctx    - PetscViewerAndFormat object

124:   Level: beginner

126:   Note:
127:   This is an LMEMonitorFn specialized for a context of PetscViewerAndFormat.

129: .seealso: `LMEMonitorSet()`, `LMEMonitorRegister()`, `LMEMonitorFn`, `LMEMonitorRegisterCreateFn`, `LMEMonitorRegisterDestroyFn`
130: S*/
131: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode LMEMonitorRegisterFn(LME lme,PetscInt its,PetscReal errest,PetscViewerAndFormat *ctx);

133: /*S
134:   LMEMonitorRegisterCreateFn - A function prototype for functions that do the creation when provided to LMEMonitorRegister()

136:   Calling Sequence:
137: +   viewer - the viewer to be used with the LMEMonitorRegisterFn
138: .   format - the format of the viewer
139: .   ctx    - a context for the monitor
140: -   result - a PetscViewerAndFormat object

142:   Level: beginner

144: .seealso: `LMEMonitorRegisterFn`, `LMEMonitorSet()`, `LMEMonitorRegister()`, `LMEMonitorFn`, `LMEMonitorRegisterDestroyFn`
145: S*/
146: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode LMEMonitorRegisterCreateFn(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **result);

148: /*S
149:   LMEMonitorRegisterDestroyFn - A function prototype for functions that do the after use destruction when provided to LMEMonitorRegister()

151:   Calling Sequence:
152: .   vf - a PetscViewerAndFormat object to be destroyed, including any context

154:   Level: beginner

156: .seealso: `LMEMonitorRegisterFn`, `LMEMonitorSet()`, `LMEMonitorRegister()`, `LMEMonitorFn`, `LMEMonitorRegisterCreateFn`
157: S*/
158: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode LMEMonitorRegisterDestroyFn(PetscViewerAndFormat **result);

160: SLEPC_EXTERN PetscErrorCode LMEMonitor(LME,PetscInt,PetscReal);
161: SLEPC_EXTERN PetscErrorCode LMEMonitorSet(LME,LMEMonitorFn,void*,PetscCtxDestroyFn*);
162: SLEPC_EXTERN PetscErrorCode LMEMonitorCancel(LME);
163: SLEPC_EXTERN PetscErrorCode LMEGetMonitorContext(LME,void*);

165: SLEPC_EXTERN PetscErrorCode LMEMonitorSetFromOptions(LME,const char[],const char[],void*);
166: SLEPC_EXTERN LMEMonitorRegisterFn       LMEMonitorDefault;
167: SLEPC_EXTERN LMEMonitorRegisterFn       LMEMonitorDefaultDrawLG;
168: SLEPC_EXTERN LMEMonitorRegisterCreateFn LMEMonitorDefaultDrawLGCreate;

170: SLEPC_EXTERN PetscErrorCode LMESetOptionsPrefix(LME,const char*);
171: SLEPC_EXTERN PetscErrorCode LMEAppendOptionsPrefix(LME,const char*);
172: SLEPC_EXTERN PetscErrorCode LMEGetOptionsPrefix(LME,const char*[]);

174: /*E
175:     LMEConvergedReason - reason a matrix function iteration was said to
176:          have converged or diverged

178:     Level: intermediate

180: .seealso: `LMESolve()`, `LMEGetConvergedReason()`, `LMESetTolerances()`
181: E*/
182: typedef enum {/* converged */
183:               LME_CONVERGED_TOL                =  1,
184:               /* diverged */
185:               LME_DIVERGED_ITS                 = -1,
186:               LME_DIVERGED_BREAKDOWN           = -2,
187:               LME_CONVERGED_ITERATING          =  0} LMEConvergedReason;
188: SLEPC_EXTERN const char *const*LMEConvergedReasons;

190: SLEPC_EXTERN PetscErrorCode LMEGetConvergedReason(LME,LMEConvergedReason *);

192: SLEPC_EXTERN PetscFunctionList LMEList;
193: SLEPC_EXTERN PetscFunctionList LMEMonitorList;
194: SLEPC_EXTERN PetscFunctionList LMEMonitorCreateList;
195: SLEPC_EXTERN PetscFunctionList LMEMonitorDestroyList;
196: SLEPC_EXTERN PetscErrorCode LMERegister(const char[],PetscErrorCode(*)(LME));
197: SLEPC_EXTERN PetscErrorCode LMEMonitorRegister(const char[],PetscViewerType,PetscViewerFormat,LMEMonitorRegisterFn*,LMEMonitorRegisterCreateFn*,LMEMonitorRegisterDestroyFn*);

199: SLEPC_EXTERN PetscErrorCode LMEAllocateSolution(LME,PetscInt);