Actual source code: mfnimpl.h
slepc-3.19.0 2023-03-31
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: */
11: #if !defined(SLEPCMFNIMPL_H)
12: #define SLEPCMFNIMPL_H
14: #include <slepcmfn.h>
15: #include <slepc/private/slepcimpl.h>
17: /* SUBMANSEC = MFN */
19: SLEPC_EXTERN PetscBool MFNRegisterAllCalled;
20: SLEPC_EXTERN PetscBool MFNMonitorRegisterAllCalled;
21: SLEPC_EXTERN PetscErrorCode MFNRegisterAll(void);
22: SLEPC_EXTERN PetscErrorCode MFNMonitorRegisterAll(void);
23: SLEPC_EXTERN PetscLogEvent MFN_SetUp,MFN_Solve;
25: typedef struct _MFNOps *MFNOps;
27: struct _MFNOps {
28: PetscErrorCode (*solve)(MFN,Vec,Vec);
29: PetscErrorCode (*setup)(MFN);
30: PetscErrorCode (*setfromoptions)(MFN,PetscOptionItems*);
31: PetscErrorCode (*publishoptions)(MFN);
32: PetscErrorCode (*destroy)(MFN);
33: PetscErrorCode (*reset)(MFN);
34: PetscErrorCode (*view)(MFN,PetscViewer);
35: };
37: /*
38: Maximum number of monitors you can run with a single MFN
39: */
40: #define MAXMFNMONITORS 5
42: /*
43: Defines the MFN data structure.
44: */
45: struct _p_MFN {
46: PETSCHEADER(struct _MFNOps);
47: /*------------------------- User parameters ---------------------------*/
48: Mat A; /* the problem matrix */
49: FN fn; /* which function to compute */
50: PetscInt max_it; /* maximum number of iterations */
51: PetscInt ncv; /* number of basis vectors */
52: PetscReal tol; /* tolerance */
53: PetscBool errorifnotconverged; /* error out if MFNSolve() does not converge */
55: /*-------------- User-provided functions and contexts -----------------*/
56: PetscErrorCode (*monitor[MAXMFNMONITORS])(MFN,PetscInt,PetscReal,void*);
57: PetscErrorCode (*monitordestroy[MAXMFNMONITORS])(void**);
58: void *monitorcontext[MAXMFNMONITORS];
59: PetscInt numbermonitors;
61: /*----------------- Child objects and working data -------------------*/
62: BV V; /* set of basis vectors */
63: Mat AT; /* the transpose of A, used in MFNSolveTranspose */
64: PetscInt nwork; /* number of work vectors */
65: Vec *work; /* work vectors */
66: void *data; /* placeholder for solver-specific stuff */
68: /* ----------------------- Status variables -------------------------- */
69: PetscInt its; /* number of iterations so far computed */
70: PetscInt nv; /* size of current Schur decomposition */
71: PetscReal errest; /* error estimate */
72: PetscReal bnorm; /* computed norm of right-hand side in current solve */
73: PetscBool transpose_solve; /* solve transpose system instead */
74: PetscInt setupcalled;
75: MFNConvergedReason reason;
76: };
78: /*
79: MFN_CreateDenseMat - Creates a dense Mat of size k unless it already has that size
80: */
81: static inline PetscErrorCode MFN_CreateDenseMat(PetscInt k,Mat *A)
82: {
83: PetscBool create=PETSC_FALSE;
84: PetscInt m,n;
86: PetscFunctionBegin;
87: if (!*A) create=PETSC_TRUE;
88: else {
89: PetscCall(MatGetSize(*A,&m,&n));
90: if (m!=k || n!=k) {
91: PetscCall(MatDestroy(A));
92: create=PETSC_TRUE;
93: }
94: }
95: if (create) PetscCall(MatCreateSeqDense(PETSC_COMM_SELF,k,k,NULL,A));
96: PetscFunctionReturn(PETSC_SUCCESS);
97: }
99: /*
100: MFN_CreateVec - Creates a Vec of size k unless it already has that size
101: */
102: static inline PetscErrorCode MFN_CreateVec(PetscInt k,Vec *v)
103: {
104: PetscBool create=PETSC_FALSE;
105: PetscInt n;
107: PetscFunctionBegin;
108: if (!*v) create=PETSC_TRUE;
109: else {
110: PetscCall(VecGetSize(*v,&n));
111: if (n!=k) {
112: PetscCall(VecDestroy(v));
113: create=PETSC_TRUE;
114: }
115: }
116: if (create) PetscCall(VecCreateSeq(PETSC_COMM_SELF,k,v));
117: PetscFunctionReturn(PETSC_SUCCESS);
118: }
120: #endif