Actual source code: slepcds.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 direct solver object in SLEPc
12: */
14: #pragma once
16: #include <slepcsc.h>
17: #include <slepcfn.h>
18: #include <slepcrg.h>
20: /* SUBMANSEC = DS */
22: #define DS_MAX_SOLVE 6
24: SLEPC_EXTERN PetscErrorCode DSInitializePackage(void);
25: SLEPC_EXTERN PetscErrorCode DSFinalizePackage(void);
27: /*S
28: DS - Direct solver (or dense system), to represent low-dimensional
29: eigenproblems that must be solved within iterative solvers. This is an
30: auxiliary object and is not normally needed by application programmers.
32: Level: beginner
34: .seealso: [](sec:ds), `DSCreate()`
35: S*/
36: typedef struct _p_DS* DS;
38: /*J
39: DSType - String with the name of the type of direct solver. Roughly,
40: there are as many types as problem types are available within SLEPc.
42: Level: beginner
44: .seealso: [](sec:ds), `DSSetType()`, `DS`
45: J*/
46: typedef const char *DSType;
47: #define DSHEP "hep"
48: #define DSNHEP "nhep"
49: #define DSGHEP "ghep"
50: #define DSGHIEP "ghiep"
51: #define DSGNHEP "gnhep"
52: #define DSNHEPTS "nhepts"
53: #define DSSVD "svd"
54: #define DSHSVD "hsvd"
55: #define DSGSVD "gsvd"
56: #define DSPEP "pep"
57: #define DSNEP "nep"
59: /* Logging support */
60: SLEPC_EXTERN PetscClassId DS_CLASSID;
62: /*E
63: DSStateType - Indicates in which state the direct solver is.
65: Values:
66: + `DS_STATE_RAW` - initial state, the matrices have not been modified yet
67: . `DS_STATE_INTERMEDIATE` - matrices have been reduced to intermediate form
68: . `DS_STATE_CONDENSED` - problem solved, matrices in condensed form
69: - `DS_STATE_TRUNCATED` - problem solved and in addition the dimension has been truncated
71: Level: advanced
73: .seealso: [](sec:ds), `DSSetState()`
74: E*/
75: typedef enum { DS_STATE_RAW,
76: DS_STATE_INTERMEDIATE,
77: DS_STATE_CONDENSED,
78: DS_STATE_TRUNCATED } DSStateType;
79: SLEPC_EXTERN const char *DSStateTypes[];
81: /*MC
82: DS_STATE_RAW - The `DS` object is in the initial state, where the matrices
83: have not been modified yet, and have no particular structure.
85: Level: advanced
87: .seealso: [](sec:ds), `DSStateType`, `DSSetState()`, `DS_STATE_INTERMEDIATE`, `DS_STATE_CONDENSED`, `DS_STATE_TRUNCATED`
88: M*/
90: /*MC
91: DS_STATE_INTERMEDIATE - The `DS` object is a intermediate state, where the
92: matrices have been reduced to an intermediate form, but the solve is not
93: finished completely.
95: Note:
96: In some cases the solve starts at the intermediate stage, e.g., in Lanczos
97: methods the projected matrix is already in tridiagonal form.
99: Level: advanced
101: .seealso: [](sec:ds), `DSStateType`, `DSSetState()`, `DS_STATE_RAW`, `DS_STATE_CONDENSED`, `DS_STATE_TRUNCATED`
102: M*/
104: /*MC
105: DS_STATE_CONDENSED - The `DS` object is in a state where the problem has
106: been solved, and matrices have been reduced to condensed form.
108: Note:
109: This state is reached after a call to `DSSolve()`.
111: Level: advanced
113: .seealso: [](sec:ds), `DSStateType`, `DSSetState()`, `DS_STATE_RAW`, `DS_STATE_INTERMEDIATE`, `DS_STATE_TRUNCATED`
114: M*/
116: /*MC
117: DS_STATE_TRUNCATED - The `DS` object is in a state where the problem is
118: solved and in addition the dimension has been truncated.
120: Note:
121: This state is reached after a call to `DSTruncate()`. The truncated size
122: can be obtained with `DSGetDimensions()`.
124: Level: advanced
126: .seealso: [](sec:ds), `DSStateType`, `DSSetState()`, `DSTruncate()`, `DSGetDimensions()`, `DS_STATE_RAW`, `DS_STATE_INTERMEDIATE`, `DS_STATE_CONDENSED`
127: M*/
129: /*E
130: DSMatType - Used to refer to one of the matrices stored internally in `DS`.
132: Values:
133: + `DS_MAT_A` - first matrix of eigenproblem/singular value problem
134: . `DS_MAT_B` - second matrix of a generalized eigenproblem
135: . `DS_MAT_C` - third matrix of a quadratic eigenproblem (deprecated)
136: . `DS_MAT_T` - tridiagonal matrix
137: . `DS_MAT_D` - diagonal matrix
138: . `DS_MAT_Q` - orthogonal matrix of (right) Schur vectors
139: . `DS_MAT_Z` - orthogonal matrix of left Schur vectors
140: . `DS_MAT_X` - right eigenvectors
141: . `DS_MAT_Y` - left eigenvectors
142: . `DS_MAT_U` - left singular vectors
143: . `DS_MAT_V` - right singular vectors
144: . `DS_MAT_W` - workspace matrix
145: - `DS_MAT_E0` to `DS_MAT_E9` - extra matrices, used in `DSPEP` and `DSNEP`
147: Notes:
148: The matrices preferentially refer to the description above, but they
149: may be used for a different purpose depending on the `DSType`.
151: All matrices can have space to hold `ld x ld` elements, except for
152: `DS_MAT_T` that has space for `3 x ld` elements (`ld` = leading dimension)
153: and `DS_MAT_D` that has space for just `ld` elements.
155: In `DSPEP` problems, matrices `A`, `B`, `W` can have space for `d*ld x d*ld`,
156: where `d` is the polynomial degree, and `X` can have `ld x d*ld`.
157: Also `DSNEP` has exceptions. Check the manual page of each `DS` type
158: for details.
160: Level: advanced
162: .seealso: [](sec:ds), `DSAllocate()`, `DSGetArray()`, `DSGetArrayReal()`, `DSVectors()`, `DSGetLeadingDimension()`
163: E*/
164: typedef enum { DS_MAT_A,
165: DS_MAT_B,
166: DS_MAT_C,
167: DS_MAT_T,
168: DS_MAT_D,
169: DS_MAT_Q,
170: DS_MAT_Z,
171: DS_MAT_X,
172: DS_MAT_Y,
173: DS_MAT_U,
174: DS_MAT_V,
175: DS_MAT_W,
176: DS_MAT_E0,
177: DS_MAT_E1,
178: DS_MAT_E2,
179: DS_MAT_E3,
180: DS_MAT_E4,
181: DS_MAT_E5,
182: DS_MAT_E6,
183: DS_MAT_E7,
184: DS_MAT_E8,
185: DS_MAT_E9,
186: DS_NUM_MAT } DSMatType;
188: /* Convenience for indexing extra matrices */
189: SLEPC_EXTERN DSMatType DSMatExtra[];
190: #define DS_NUM_EXTRA 10
192: /*E
193: DSParallelType - Indicates the parallel mode that the direct solver will use.
195: Values:
196: + `DS_PARALLEL_REDUNDANT` - redundant computation
197: . `DS_PARALLEL_SYNCHRONIZED` - only one process computes the solution
198: - `DS_PARALLEL_DISTRIBUTED` - all processes participate in the solution
200: Level: advanced
202: .seealso: [](sec:ds), `DSSetParallel()`
203: E*/
204: typedef enum { DS_PARALLEL_REDUNDANT,
205: DS_PARALLEL_SYNCHRONIZED,
206: DS_PARALLEL_DISTRIBUTED } DSParallelType;
207: SLEPC_EXTERN const char *DSParallelTypes[];
209: /*MC
210: DS_PARALLEL_REDUNDANT - In this parallel mode, all processes will do
211: the computation redundantly, starting from the same data, and producing
212: the same result.
214: Note:
215: The result may be slightly different in the different processes if using a
216: multithreaded BLAS library, which may cause issues in ill-conditioned problems.
218: Level: advanced
220: .seealso: [](sec:ds), `DSParallelType`, `DSSetParallel()`, `DS_PARALLEL_SYNCHRONIZED`, `DS_PARALLEL_DISTRIBUTED`
221: M*/
223: /*MC
224: DS_PARALLEL_SYNCHRONIZED - In this parallel mode, only the first MPI process
225: performs the computation and then the computed quantities are broadcast to the
226: other processes in the communicator.
228: Note:
229: The communication is not done automatically, an explicit call to `DSSynchronize()`
230: is required.
232: Level: advanced
234: .seealso: [](sec:ds), `DSParallelType`, `DSSetParallel()`, `DSSynchronize()`, `DS_PARALLEL_REDUNDANT`, `DS_PARALLEL_DISTRIBUTED`
235: M*/
237: /*MC
238: DS_PARALLEL_DISTRIBUTED - In this parallel mode, every MPI process will be
239: in charge of part of the computation.
241: Note:
242: This parallel mode can be used in some `DS` types only, such as the contour
243: integral method of `DSNEP`.
245: Level: advanced
247: .seealso: [](sec:ds), `DSParallelType`, `DSSetParallel()`, `DS_PARALLEL_REDUNDANT`, `DS_PARALLEL_SYNCHRONIZED`
248: M*/
250: SLEPC_EXTERN PetscErrorCode DSCreate(MPI_Comm,DS*);
251: SLEPC_EXTERN PetscErrorCode DSSetType(DS,DSType);
252: SLEPC_EXTERN PetscErrorCode DSGetType(DS,DSType*);
253: SLEPC_EXTERN PetscErrorCode DSSetOptionsPrefix(DS,const char[]);
254: SLEPC_EXTERN PetscErrorCode DSAppendOptionsPrefix(DS,const char[]);
255: SLEPC_EXTERN PetscErrorCode DSGetOptionsPrefix(DS,const char*[]);
256: SLEPC_EXTERN PetscErrorCode DSSetFromOptions(DS);
257: SLEPC_EXTERN PetscErrorCode DSView(DS,PetscViewer);
258: SLEPC_EXTERN PetscErrorCode DSViewFromOptions(DS,PetscObject,const char[]);
259: SLEPC_EXTERN PetscErrorCode DSViewMat(DS,PetscViewer,DSMatType);
260: SLEPC_EXTERN PetscErrorCode DSDestroy(DS*);
261: SLEPC_EXTERN PetscErrorCode DSReset(DS);
262: SLEPC_EXTERN PetscErrorCode DSDuplicate(DS,DS*);
264: SLEPC_EXTERN PetscErrorCode DSAllocate(DS,PetscInt);
265: SLEPC_EXTERN PetscErrorCode DSReallocate(DS,PetscInt);
266: SLEPC_EXTERN PetscErrorCode DSGetLeadingDimension(DS,PetscInt*);
267: SLEPC_EXTERN PetscErrorCode DSSetState(DS,DSStateType);
268: SLEPC_EXTERN PetscErrorCode DSGetState(DS,DSStateType*);
269: SLEPC_EXTERN PetscErrorCode DSSetDimensions(DS,PetscInt,PetscInt,PetscInt);
270: SLEPC_EXTERN PetscErrorCode DSGetDimensions(DS,PetscInt*,PetscInt*,PetscInt*,PetscInt*);
271: SLEPC_EXTERN PetscErrorCode DSSetBlockSize(DS,PetscInt);
272: SLEPC_EXTERN PetscErrorCode DSGetBlockSize(DS,PetscInt*);
273: SLEPC_EXTERN PetscErrorCode DSGetTruncateSize(DS,PetscInt,PetscInt,PetscInt*);
274: SLEPC_EXTERN PetscErrorCode DSTruncate(DS,PetscInt,PetscBool);
275: SLEPC_EXTERN PetscErrorCode DSSetIdentity(DS,DSMatType);
276: SLEPC_EXTERN PetscErrorCode DSSetMethod(DS,PetscInt);
277: SLEPC_EXTERN PetscErrorCode DSGetMethod(DS,PetscInt*);
278: SLEPC_EXTERN PetscErrorCode DSSetParallel(DS,DSParallelType);
279: SLEPC_EXTERN PetscErrorCode DSGetParallel(DS,DSParallelType*);
280: SLEPC_EXTERN PetscErrorCode DSSetCompact(DS,PetscBool);
281: SLEPC_EXTERN PetscErrorCode DSGetCompact(DS,PetscBool*);
282: SLEPC_EXTERN PetscErrorCode DSSetExtraRow(DS,PetscBool);
283: SLEPC_EXTERN PetscErrorCode DSGetExtraRow(DS,PetscBool*);
284: SLEPC_EXTERN PetscErrorCode DSSetRefined(DS,PetscBool);
285: SLEPC_EXTERN PetscErrorCode DSGetRefined(DS,PetscBool*);
286: SLEPC_EXTERN PetscErrorCode DSGetMat(DS,DSMatType,Mat*);
287: SLEPC_EXTERN PetscErrorCode DSRestoreMat(DS,DSMatType,Mat*);
288: SLEPC_EXTERN PetscErrorCode DSGetMatAndColumn(DS,DSMatType,PetscInt,Mat*,Vec*);
289: SLEPC_EXTERN PetscErrorCode DSRestoreMatAndColumn(DS,DSMatType,PetscInt,Mat*,Vec*);
290: SLEPC_EXTERN PetscErrorCode DSGetArray(DS,DSMatType,PetscScalar*[]);
291: SLEPC_EXTERN PetscErrorCode DSRestoreArray(DS,DSMatType,PetscScalar*[]);
292: SLEPC_EXTERN PetscErrorCode DSGetArrayReal(DS,DSMatType,PetscReal*[]);
293: SLEPC_EXTERN PetscErrorCode DSRestoreArrayReal(DS,DSMatType,PetscReal*[]);
294: SLEPC_EXTERN PetscErrorCode DSVectors(DS,DSMatType,PetscInt*,PetscReal*);
295: SLEPC_EXTERN PetscErrorCode DSSolve(DS,PetscScalar[],PetscScalar[]);
296: SLEPC_EXTERN PetscErrorCode DSSort(DS,PetscScalar[],PetscScalar[],PetscScalar[],PetscScalar[],PetscInt*);
297: SLEPC_EXTERN PetscErrorCode DSSortWithPermutation(DS,PetscInt[],PetscScalar[],PetscScalar[]);
298: SLEPC_EXTERN PetscErrorCode DSSynchronize(DS,PetscScalar[],PetscScalar[]);
299: PETSC_DEPRECATED_FUNCTION(3, 18, 0, "DSGetMat()+MatDenseGetSubMatrix()+MatCopy()", ) static inline PetscErrorCode DSCopyMat(DS ds,DSMatType m,PetscInt mr,PetscInt mc,Mat A,PetscInt Ar,PetscInt Ac,PetscInt rows,PetscInt cols,PetscBool out)
300: {
301: Mat M,M0,A0;
303: PetscFunctionBegin;
304: PetscCall(DSGetMat(ds,m,&M));
305: PetscCall(MatDenseGetSubMatrix(M,mr,mr+rows,mc,mc+cols,&M0));
306: PetscCall(MatDenseGetSubMatrix(A,Ar,Ar+rows,Ac,Ac+cols,&A0));
307: if (out) PetscCall(MatCopy(M0,A0,SAME_NONZERO_PATTERN));
308: else PetscCall(MatCopy(A0,M0,SAME_NONZERO_PATTERN));
309: PetscCall(MatDenseRestoreSubMatrix(M,&M0));
310: PetscCall(MatDenseRestoreSubMatrix(A,&A0));
311: PetscCall(DSRestoreMat(ds,m,&M));
312: PetscFunctionReturn(PETSC_SUCCESS);
313: }
314: SLEPC_EXTERN PetscErrorCode DSMatGetSize(DS,DSMatType,PetscInt*,PetscInt*);
315: SLEPC_EXTERN PetscErrorCode DSMatIsHermitian(DS,DSMatType,PetscBool*);
316: SLEPC_EXTERN PetscErrorCode DSSetSlepcSC(DS,SlepcSC);
317: SLEPC_EXTERN PetscErrorCode DSGetSlepcSC(DS,SlepcSC*);
318: SLEPC_EXTERN PetscErrorCode DSUpdateExtraRow(DS);
319: SLEPC_EXTERN PetscErrorCode DSCond(DS,PetscReal*);
320: SLEPC_EXTERN PetscErrorCode DSTranslateHarmonic(DS,PetscScalar,PetscReal,PetscBool,PetscScalar[],PetscReal*);
321: SLEPC_EXTERN PetscErrorCode DSTranslateRKS(DS,PetscScalar);
322: SLEPC_EXTERN PetscErrorCode DSOrthogonalize(DS,DSMatType,PetscInt,PetscInt*);
323: SLEPC_EXTERN PetscErrorCode DSPseudoOrthogonalize(DS,DSMatType,PetscInt,PetscReal[],PetscInt*,PetscReal[]);
325: /* --------- options specific to particular solvers -------- */
327: SLEPC_EXTERN PetscErrorCode DSSVDSetDimensions(DS,PetscInt);
328: SLEPC_EXTERN PetscErrorCode DSSVDGetDimensions(DS,PetscInt*);
329: SLEPC_EXTERN PetscErrorCode DSGSVDSetDimensions(DS,PetscInt,PetscInt);
330: SLEPC_EXTERN PetscErrorCode DSGSVDGetDimensions(DS,PetscInt*,PetscInt*);
331: SLEPC_EXTERN PetscErrorCode DSHSVDSetDimensions(DS,PetscInt);
332: SLEPC_EXTERN PetscErrorCode DSHSVDGetDimensions(DS,PetscInt*);
333: SLEPC_EXTERN PetscErrorCode DSHSVDSetReorthogonalize(DS,PetscBool);
334: SLEPC_EXTERN PetscErrorCode DSHSVDGetReorthogonalize(DS,PetscBool*);
336: SLEPC_EXTERN PetscErrorCode DSPEPSetDegree(DS,PetscInt);
337: SLEPC_EXTERN PetscErrorCode DSPEPGetDegree(DS,PetscInt*);
338: SLEPC_EXTERN PetscErrorCode DSPEPSetCoefficients(DS,PetscReal[]);
339: SLEPC_EXTERN PetscErrorCode DSPEPGetCoefficients(DS,PetscReal*[]);
341: SLEPC_EXTERN PetscErrorCode DSNEPSetFN(DS,PetscInt,FN[]);
342: SLEPC_EXTERN PetscErrorCode DSNEPGetFN(DS,PetscInt,FN*);
343: SLEPC_EXTERN PetscErrorCode DSNEPGetNumFN(DS,PetscInt*);
344: SLEPC_EXTERN PetscErrorCode DSNEPSetMinimality(DS,PetscInt);
345: SLEPC_EXTERN PetscErrorCode DSNEPGetMinimality(DS,PetscInt*);
346: SLEPC_EXTERN PetscErrorCode DSNEPSetRefine(DS,PetscReal,PetscInt);
347: SLEPC_EXTERN PetscErrorCode DSNEPGetRefine(DS,PetscReal*,PetscInt*);
348: SLEPC_EXTERN PetscErrorCode DSNEPSetIntegrationPoints(DS,PetscInt);
349: SLEPC_EXTERN PetscErrorCode DSNEPGetIntegrationPoints(DS,PetscInt*);
350: SLEPC_EXTERN PetscErrorCode DSNEPSetSamplingSize(DS,PetscInt);
351: SLEPC_EXTERN PetscErrorCode DSNEPGetSamplingSize(DS,PetscInt*);
352: SLEPC_EXTERN PetscErrorCode DSNEPSetRG(DS,RG);
353: SLEPC_EXTERN PetscErrorCode DSNEPGetRG(DS,RG*);
355: /*S
356: DSNEPMatrixFunctionFn - A prototype of a `DSNEP` compute matrix function that
357: would be passed to `DSNEPSetComputeMatrixFunction()`.
359: Calling Sequence:
360: + ds - the direct solver object
361: . lambda - point where $T(\lambda)$ or $T'(\lambda)$ must be evaluated
362: . deriv - if true compute $T'(\lambda)$, otherwise compute $T(\lambda)$
363: . mat - the `DS` matrix where the result must be stored
364: - ctx - [optional] user-defined context for private data for the
365: matrix evaluation routine (may be `NULL`)
367: Level: developer
369: .seealso: [](sec:ds), `DSNEPSetComputeMatrixFunction()`
370: S*/
371: PETSC_EXTERN_TYPEDEF typedef PetscErrorCode DSNEPMatrixFunctionFn(DS ds,PetscScalar lambda,PetscBool deriv,DSMatType mat,void *ctx);
373: SLEPC_EXTERN PetscErrorCode DSNEPSetComputeMatrixFunction(DS,DSNEPMatrixFunctionFn*,void*);
374: SLEPC_EXTERN PetscErrorCode DSNEPGetComputeMatrixFunction(DS,DSNEPMatrixFunctionFn**,void*);
376: SLEPC_EXTERN PetscFunctionList DSList;
377: SLEPC_EXTERN PetscErrorCode DSRegister(const char[],PetscErrorCode(*)(DS));