1: /*
  2:    Private data structure used by the EKSM method
  3: */
  4: #pragma once

  6: #include <petsc/private/kspimpl.h>
  7: #include <../src/mat/impls/nest/matnestimpl.h>

  9: typedef struct {
 10:   /* User data */
 11:   PetscScalar shift;
 12:   PetscBool   shift_set;
 13:   KSP         ksps; /* linear solver for K + shift M */
 14:   KSP         kspm; /* linear solver for M */

 16:   /* Hessenberg matrix and orthogonalization information. */
 17:   PetscScalar *hh_origin;     /* holds the upper Hessenberg matrix */
 18:   PetscScalar *kk_origin;     /* holds the second matrix of the Krylov relation */
 19:   PetscScalar *tt_origin;     /* temporary matrix with the coefficient of the projected problem */
 20:   PetscScalar *yy_origin;     /* solutions of the projected problem */
 21:   PetscScalar *rs_origin;     /* holds the right-hand side of the Hessenberg system */
 22:   PetscInt     ldh, ldk, ldt; /* leading dimensions */
 23:   PetscInt     factor;        /* defaults to 1, will be 2 in case of complex-conjugate pairs */
 24:   PetscScalar *work;          /* LAPACK workspace */
 25:   PetscInt     lwork;         /* length of workspace */

 27:   PetscReal haptol; /* tolerance for happy breakdown */
 28:   PetscReal v0norm; /* norm of first basis vector (M^{-1}b) */

 30:   Vec      *vecs;           /* the work vectors */
 31:   PetscInt  delta_allocate; /* number of vectors to preallocate in each block if not preallocated */
 32:   PetscInt  vv_allocated;   /* number of allocated Krylov vectors */
 33:   PetscInt  vecs_allocated; /* total number of vecs available */
 34:   Vec     **user_work;
 35:   PetscInt *mwork_alloc; /* number of work vectors allocated as part of a work-vector chunk */
 36:   PetscInt  nwork_alloc; /* number of work vector chunks allocated */
 37: } KSP_EKSM;

 39: #define HH(a, b) (eksm->hh_origin + (b) * eksm->ldh + (a))
 40: #define KK(a, b) (eksm->kk_origin + (b) * eksm->ldk + (a))
 41: #define TT(a, b) (eksm->tt_origin + (b) * eksm->ldt + (a))
 42: #define YY(a, b) (eksm->yy_origin + (b) * eksm->ldk + (a))

 44: #define SHIFT_IS_COMPLEX(cmplx, i) (!PetscDefined(USE_COMPLEX) && (cmplx)[(i)])

 46: /* vector names */
 47: #define VEC_OFFSET 1
 48: #define VEC_TEMP   eksm->vecs[0]
 49: #define VEC_VV(i)  eksm->vecs[VEC_OFFSET + (i)]