Actual source code: pjdp.h

slepc-3.7.3 2016-09-29
Report Typos and Errors
  1: /*
  2:    Private header for PEPJD.

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain

  8:    This file is part of SLEPc.

 10:    SLEPc is free software: you can redistribute it and/or modify it under  the
 11:    terms of version 3 of the GNU Lesser General Public License as published by
 12:    the Free Software Foundation.

 14:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY
 15:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS
 16:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for
 17:    more details.

 19:    You  should have received a copy of the GNU Lesser General  Public  License
 20:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 21:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 22: */


 27: typedef struct {
 28:   PetscReal   keep;          /* restart parameter */
 29:   BV          V;             /* work basis vectors to store the search space */
 30:   BV          W;             /* work basis vectors to store the test space */
 31:   BV          *TV;           /* work basis vectors to store T*V (each TV[i] is the coefficient for \lambda^i of T*V for the extended T) */
 32:   BV          *AX;           /* work basis vectors to store A_i*X for locked eigenvectors */
 33:   BV          X;             /* locked eigenvectors */
 34:   PetscScalar *T;            /* matrix of the invariant pair */
 35:   PetscScalar *Tj;           /* matrix containing the powers of the invariant pair matrix */
 36:   PetscScalar *XpX;          /* X^H*X */
 37:   PC          pcshell;       /* preconditioner including basic precond+projector */
 38:   Mat         Pshell;        /* auxiliary shell matrix */
 39:   PetscInt    nconv;         /* number of locked vectors in the invariant pair */
 40: } PEP_JD;

 42: typedef struct {
 43:   PC          pc;            /* basic preconditioner */
 44:   Vec         Bp;            /* preconditioned residual of derivative polynomial, B\p */
 45:   Vec         u;             /* Ritz vector */
 46:   PetscScalar gamma;         /* precomputed scalar u'*B\p */
 47:   PetscScalar *M;
 48:   PetscScalar *ps;
 49:   PetscInt    ld;
 50:   Vec         *work;
 51:   BV          X;
 52:   PetscInt    n;
 53: } PEP_JD_PCSHELL;

 55: typedef struct {
 56:   Mat         P;             /*  */
 57:   PEP         pep;
 58:   Vec         *work;
 59:   PetscScalar theta;
 60: } PEP_JD_MATSHELL;

 62: PETSC_INTERN PetscErrorCode PEPView_JD(PEP,PetscViewer);
 63: PETSC_INTERN PetscErrorCode PEPSetFromOptions_JD(PetscOptionItems*,PEP);
 64: PETSC_INTERN PetscErrorCode PEPJDSetRestart_JD(PEP,PetscReal);
 65: PETSC_INTERN PetscErrorCode PEPJDGetRestart_JD(PEP,PetscReal*);

 67: #endif