Actual source code: slepcsvd.h
1: /*
2: User interface for the SLEPC singular value solvers.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2007, Universidad Politecnica de Valencia, Spain
8: This file is part of SLEPc. See the README file for conditions of use
9: and additional information.
10: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
11: */
15: #include slepc.h
16: #include slepceps.h
21: /*S
22: SVD - Abstract SLEPc object that manages all the singular value
23: problem solvers.
25: Level: beginner
27: .seealso: SVDCreate()
28: S*/
29: typedef struct _p_SVD* SVD;
31: /*E
32: SVDType - String with the name of a SLEPc singular value solver
34: Level: beginner
36: .seealso: SVDSetType(), SVD
37: E*/
38: #define SVDType const char*
39: #define SVDCROSS "cross"
40: #define SVDCYCLIC "cyclic"
41: #define SVDLAPACK "lapack"
42: #define SVDLANCZOS "lanczos"
43: #define SVDTRLANCZOS "trlanczos"
45: /*E
46: SVDTransposeMode - determines how to handle the transpose of the matrix
48: Level: advanced
50: .seealso: SVDSetTransposeMode(), SVDGetTransposeMode()
51: E*/
52: typedef enum { SVD_TRANSPOSE_EXPLICIT, SVD_TRANSPOSE_IMPLICIT } SVDTransposeMode;
54: /*E
55: SVDWhich - determines whether largest or smallest singular triplets
56: are to be computed
58: Level: intermediate
60: .seealso: SVDSetWhichSingularTriplets(), SVDGetWhichSingularTriplets()
61: E*/
62: typedef enum { SVD_LARGEST, SVD_SMALLEST } SVDWhich;
64: /*E
65: SVDConvergedReason - reason a singular value solver was said to
66: have converged or diverged
68: Level: beginner
70: .seealso: SVDSolve(), SVDGetConvergedReason(), SVDSetTolerances()
71: E*/
72: typedef enum {/* converged */
73: SVD_CONVERGED_TOL = 2,
74: /* diverged */
75: SVD_DIVERGED_ITS = -3,
76: SVD_DIVERGED_BREAKDOWN = -4,
77: SVD_CONVERGED_ITERATING = 0 } SVDConvergedReason;
79: EXTERN PetscErrorCode SVDCreate(MPI_Comm,SVD*);
80: EXTERN PetscErrorCode SVDSetIP(SVD,IP);
81: EXTERN PetscErrorCode SVDGetIP(SVD,IP*);
82: EXTERN PetscErrorCode SVDSetType(SVD,SVDType);
83: EXTERN PetscErrorCode SVDGetType(SVD,SVDType*);
84: EXTERN PetscErrorCode SVDSetOperator(SVD,Mat);
85: EXTERN PetscErrorCode SVDGetOperator(SVD,Mat*);
86: EXTERN PetscErrorCode SVDSetInitialVector(SVD,Vec);
87: EXTERN PetscErrorCode SVDGetInitialVector(SVD,Vec*);
88: EXTERN PetscErrorCode SVDSetTransposeMode(SVD,SVDTransposeMode);
89: EXTERN PetscErrorCode SVDGetTransposeMode(SVD,SVDTransposeMode*);
90: EXTERN PetscErrorCode SVDSetDimensions(SVD,int,int);
91: EXTERN PetscErrorCode SVDGetDimensions(SVD,int*,int*);
92: EXTERN PetscErrorCode SVDSetTolerances(SVD,PetscReal,int);
93: EXTERN PetscErrorCode SVDGetTolerances(SVD,PetscReal*,int*);
94: EXTERN PetscErrorCode SVDSetWhichSingularTriplets(SVD,SVDWhich);
95: EXTERN PetscErrorCode SVDGetWhichSingularTriplets(SVD,SVDWhich*);
96: EXTERN PetscErrorCode SVDSetFromOptions(SVD);
97: EXTERN PetscErrorCode SVDSetOptionsPrefix(SVD,const char*);
98: EXTERN PetscErrorCode SVDAppendOptionsPrefix(SVD,const char*);
99: EXTERN PetscErrorCode SVDGetOptionsPrefix(SVD,const char*[]);
100: EXTERN PetscErrorCode SVDSetUp(SVD);
101: EXTERN PetscErrorCode SVDSolve(SVD);
102: EXTERN PetscErrorCode SVDGetIterationNumber(SVD,int*);
103: EXTERN PetscErrorCode SVDGetConvergedReason(SVD,SVDConvergedReason*);
104: EXTERN PetscErrorCode SVDGetConverged(SVD,int*);
105: EXTERN PetscErrorCode SVDGetSingularTriplet(SVD,int,PetscReal*,Vec,Vec);
106: EXTERN PetscErrorCode SVDComputeResidualNorms(SVD,int,PetscReal*,PetscReal*);
107: EXTERN PetscErrorCode SVDComputeRelativeError(SVD,int,PetscReal*);
108: EXTERN PetscErrorCode SVDGetOperationCounters(SVD,int*,int*);
109: EXTERN PetscErrorCode SVDView(SVD,PetscViewer);
110: EXTERN PetscErrorCode SVDDestroy(SVD);
111: EXTERN PetscErrorCode SVDInitializePackage(char*);
113: EXTERN PetscErrorCode SVDMonitorSet(SVD,PetscErrorCode (*)(SVD,int,int,PetscReal*,PetscReal*,int,void*),
114: void*,PetscErrorCode (*monitordestroy)(void*));
115: EXTERN PetscErrorCode SVDMonitorCancel(SVD);
116: EXTERN PetscErrorCode SVDGetMonitorContext(SVD,void **);
117: EXTERN PetscErrorCode SVDMonitorDefault(SVD,int,int,PetscReal*,PetscReal*,int,void*);
118: EXTERN PetscErrorCode SVDMonitorLG(SVD,int,int,PetscReal*,PetscReal*,int,void*);
120: EXTERN PetscErrorCode SVDDense(int,int,PetscScalar*,PetscReal*,PetscScalar*,PetscScalar*);
122: EXTERN PetscErrorCode SVDCrossSetEPS(SVD,EPS);
123: EXTERN PetscErrorCode SVDCrossGetEPS(SVD,EPS*);
125: EXTERN PetscErrorCode SVDCyclicSetExplicitMatrix(SVD,PetscTruth);
126: EXTERN PetscErrorCode SVDCyclicGetExplicitMatrix(SVD,PetscTruth*);
127: EXTERN PetscErrorCode SVDCyclicSetEPS(SVD,EPS);
128: EXTERN PetscErrorCode SVDCyclicGetEPS(SVD,EPS*);
130: EXTERN PetscErrorCode SVDLanczosSetOneSide(SVD,PetscTruth);
132: EXTERN PetscErrorCode SVDTRLanczosSetOneSide(SVD,PetscTruth);
134: EXTERN PetscErrorCode SVDRegister(const char*,const char*,const char*,int(*)(SVD));
135: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
136: #define SVDRegisterDynamic(a,b,c,d) SVDRegister(a,b,c,0)
137: #else
138: #define SVDRegisterDynamic(a,b,c,d) SVDRegister(a,b,c,d)
139: #endif
140: EXTERN PetscErrorCode SVDRegisterDestroy(void);
142: #endif