Actual source code: nleigs.h
slepc-3.19.2 2023-09-05
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: Private header for NLEIGS
12: */
14: #if !defined(SLEPC_NLEIGS_H)
15: #define SLEPC_NLEIGS_H
17: #define LBPOINTS 100 /* default value of the maximum number of Leja-Bagby points */
18: #define NDPOINTS 1e4 /* number of discretization points */
20: typedef struct {
21: BV V; /* tensor vector basis for the linearization */
22: BV W; /* tensor vector basis for the linearization */
23: PetscInt nmat; /* number of interpolation points */
24: PetscScalar *s,*xi; /* Leja-Bagby points */
25: PetscScalar *beta; /* scaling factors */
26: Mat *D; /* divided difference matrices */
27: PetscScalar *coeffD; /* coefficients for divided differences in split form */
28: PetscInt nshifts; /* provided number of shifts */
29: PetscScalar *shifts; /* user-provided shifts for the Rational Krylov variant */
30: PetscInt nshiftsw; /* actual number of shifts (1 if Krylov-Schur) */
31: PetscReal ddtol; /* tolerance for divided difference convergence */
32: PetscInt ddmaxit; /* maximum number of divided difference terms */
33: PetscReal keep; /* restart parameter */
34: PetscBool lock; /* locking/non-locking variant */
35: PetscInt idxrk; /* index of next shift to use */
36: KSP *ksp; /* ksp array for storing shift factorizations */
37: Vec vrn; /* random vector with normally distributed value */
38: PetscBool fullbasis; /* use full Krylov basis instead of TOAR basis */
39: EPS eps; /* eigensolver used in the full basis variant */
40: Mat A; /* shell matrix used for the eps in full basis */
41: Vec w[6]; /* work vectors */
42: void *singularitiesctx;
43: PetscErrorCode (*computesingularities)(NEP,PetscInt*,PetscScalar*,void*);
44: } NEP_NLEIGS;
46: typedef struct {
47: PetscInt nmat,maxnmat;
48: PetscScalar *coeff;
49: Mat *A;
50: Vec t;
51: } NEP_NLEIGS_MATSHELL;
53: static inline PetscErrorCode NEPNLEIGSSetShifts(NEP nep,PetscInt *nshiftsw)
54: {
55: NEP_NLEIGS *ctx = (NEP_NLEIGS*)nep->data;
57: PetscFunctionBegin;
58: if (!ctx->nshifts) {
59: ctx->shifts = &nep->target;
60: *nshiftsw = 1;
61: } else *nshiftsw = ctx->nshifts;
62: PetscFunctionReturn(PETSC_SUCCESS);
63: }
65: SLEPC_INTERN PetscErrorCode NEPSetUp_NLEIGS_FullBasis(NEP);
66: SLEPC_INTERN PetscErrorCode NEPNLEIGSSetEPS_NLEIGS(NEP,EPS);
67: SLEPC_INTERN PetscErrorCode NEPNLEIGSGetEPS_NLEIGS(NEP,EPS*);
68: SLEPC_INTERN PetscErrorCode NEPNLEIGSBackTransform(PetscObject,PetscInt,PetscScalar*,PetscScalar *vali);
69: SLEPC_INTERN PetscErrorCode NEPNLEIGSEvalNRTFunct(NEP,PetscInt,PetscScalar,PetscScalar*);
70: SLEPC_INTERN PetscErrorCode NEPSolve_NLEIGS_FullBasis(NEP);
71: SLEPC_INTERN PetscErrorCode NEPSolve_NLEIGS(NEP);
73: #endif