Actual source code: nepdefl.h

slepc-3.11.2 2019-07-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2017, 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 Deflation in NEP
 12: */
 13: #if !defined(SLEPC_NEPDEFL_H)
 14: #define SLEPC_NEPDEFL_H

 16: # define MAX_MINIDX 1

 18: typedef struct _n_nep_ext_op *NEP_EXT_OP;
 19: typedef struct _n_nep_def_fun_solve *NEP_DEF_FUN_SOLVE;
 20: typedef struct _n_nep_def_project *NEP_DEF_PROJECT;

 22: /* Structure characterizing a deflation context */
 23: struct _n_nep_ext_op {
 24:   NEP               nep;
 25:   PetscScalar       *H;     /* invariant pair (X,H) */
 26:   BV                X;      /* locked eigenvectors */
 27:   PetscScalar       *bc;    /* polinomial basis roots */
 28:   RG                rg;
 29:   PetscInt          midx;   /* minimality index */
 30:   PetscInt          max_midx;
 31:   PetscInt          szd;    /* maximum size for deflation */
 32:   PetscInt          n;      /* invariant pair size */
 33:   Mat               MF;     /* function shell matrix */
 34:   Mat               MJ;     /* Jacobian shell matrix */
 35:   PetscBool         simpU;  /* the way U is computed */
 36:   NEP_DEF_FUN_SOLVE solve;  /* MatSolve context for the operator */
 37:   NEP_DEF_PROJECT   proj;   /* context for the projected eigenproblem */
 38:   /* auxiliary computations */
 39:   BV                W;
 40:   PetscScalar       *Hj;    /* matrix containing the powers of the invariant pair matrix */
 41:   PetscScalar       *XpX;   /* X^*X */
 42: };

 44: struct _n_nep_def_fun_solve {
 45:   KSP          ksp;   /* */
 46:   PetscBool    sincf;
 47:   Mat          T;
 48:   PetscScalar  theta;
 49:   PetscInt     n;
 50:   PetscScalar  *M;
 51:   PetscScalar  *work;
 52:   Vec          w[2];
 53:   BV           T_1U;
 54:   NEP_EXT_OP   extop;
 55: };

 57: typedef struct {
 58:   NEP          nep;
 59:   Mat          T;
 60:   BV           U;
 61:   PetscScalar  *A;
 62:   PetscScalar  *B;  
 63:   PetscScalar  theta;
 64:   PetscInt     n;
 65:   NEP_EXT_OP   extop;
 66:   PetscBool    jacob;
 67:   Vec          w[2];
 68:   PetscScalar  *work;
 69:   PetscScalar  *hfj;
 70:   PetscScalar  *hfjp;
 71:   PetscBool    hfjset;
 72: } NEP_DEF_MATSHELL;

 74: struct _n_nep_def_project {
 75:   Mat          *V1pApX;
 76:   Mat          XpV1;
 77:   PetscScalar  *V2;
 78:   Vec          w;
 79:   BV           V1;
 80:   PetscInt     dim;
 81:   PetscScalar  *work;
 82:   PetscInt     lwork;
 83:   NEP_EXT_OP   extop;
 84: };

 86: #if 0
 87: typedef struct {
 88:   PC          pc;      /* basic preconditioner */
 89:   PetscScalar *M;
 90:   PetscScalar *ps;
 91:   PetscInt    ld;
 92:   Vec         *work;
 93:   BV          X;
 94:   PetscInt    n;
 95: } NEP_DEF_PCSHELL;
 96: #endif
 97: #endif

 99: SLEPC_INTERN PetscErrorCode NEPDeflationCopyToExtendedVec(NEP_EXT_OP,Vec,PetscScalar*,Vec,PetscBool);
100: SLEPC_INTERN PetscErrorCode NEPDeflationReset(NEP_EXT_OP);
101: SLEPC_INTERN PetscErrorCode NEPDeflationInitialize(NEP,BV,KSP,PetscBool,PetscInt,NEP_EXT_OP*);
102: SLEPC_INTERN PetscErrorCode NEPDeflationCreateVec(NEP_EXT_OP,Vec*);
103: SLEPC_INTERN PetscErrorCode NEPDeflationComputeFunction(NEP_EXT_OP,PetscScalar,Mat*);
104: SLEPC_INTERN PetscErrorCode NEPDeflationComputeJacobian(NEP_EXT_OP,PetscScalar,Mat*);
105: SLEPC_INTERN PetscErrorCode NEPDeflationSolveSetUp(NEP_EXT_OP,PetscScalar);
106: SLEPC_INTERN PetscErrorCode NEPDeflationFunctionSolve(NEP_EXT_OP,Vec,Vec);
107: SLEPC_INTERN PetscErrorCode NEPDeflationGetInvariantPair(NEP_EXT_OP,BV*,Mat*);
108: SLEPC_INTERN PetscErrorCode NEPDeflationLocking(NEP_EXT_OP,Vec,PetscScalar);
109: SLEPC_INTERN PetscErrorCode NEPDeflationSetRandomVec(NEP_EXT_OP,Vec);
110: SLEPC_INTERN PetscErrorCode NEPDeflationProjectOperator(NEP_EXT_OP,BV,DS,PetscInt,PetscInt);
111: SLEPC_INTERN PetscErrorCode NEPDeflationCreateBV(NEP_EXT_OP,PetscInt,BV*);