Actual source code: mfnimpl.h

slepc-3.7.0 2016-05-16
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.

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

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

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

 22: #if !defined(_MFNIMPL)
 23: #define _MFNIMPL

 25: #include <slepcmfn.h>
 26: #include <slepc/private/slepcimpl.h>

 28: PETSC_EXTERN PetscBool MFNRegisterAllCalled;
 29: PETSC_EXTERN PetscErrorCode MFNRegisterAll(void);
 30: PETSC_EXTERN PetscLogEvent MFN_SetUp, MFN_Solve;

 32: typedef struct _MFNOps *MFNOps;

 34: struct _MFNOps {
 35:   PetscErrorCode (*solve)(MFN,Vec,Vec);
 36:   PetscErrorCode (*setup)(MFN);
 37:   PetscErrorCode (*setfromoptions)(PetscOptionItems*,MFN);
 38:   PetscErrorCode (*publishoptions)(MFN);
 39:   PetscErrorCode (*destroy)(MFN);
 40:   PetscErrorCode (*reset)(MFN);
 41:   PetscErrorCode (*view)(MFN,PetscViewer);
 42: };

 44: /*
 45:      Maximum number of monitors you can run with a single MFN
 46: */
 47: #define MAXMFNMONITORS 5

 49: /*
 50:    Defines the MFN data structure.
 51: */
 52: struct _p_MFN {
 53:   PETSCHEADER(struct _MFNOps);
 54:   /*------------------------- User parameters ---------------------------*/
 55:   Mat            A;              /* the problem matrix */
 56:   FN             fn;             /* which function to compute */
 57:   PetscInt       max_it;         /* maximum number of iterations */
 58:   PetscInt       ncv;            /* number of basis vectors */
 59:   PetscReal      tol;            /* tolerance */
 60:   PetscBool      errorifnotconverged;    /* error out if MFNSolve() does not converge */

 62:   /*-------------- User-provided functions and contexts -----------------*/
 63:   PetscErrorCode (*monitor[MAXMFNMONITORS])(MFN,PetscInt,PetscReal,void*);
 64:   PetscErrorCode (*monitordestroy[MAXMFNMONITORS])(void**);
 65:   void           *monitorcontext[MAXMFNMONITORS];
 66:   PetscInt       numbermonitors;

 68:   /*----------------- Child objects and working data -------------------*/
 69:   BV             V;              /* set of basis vectors */
 70:   PetscInt       nwork;          /* number of work vectors */
 71:   Vec            *work;          /* work vectors */
 72:   void           *data;          /* placeholder for solver-specific stuff */

 74:   /* ----------------------- Status variables -------------------------- */
 75:   PetscInt       its;            /* number of iterations so far computed */
 76:   PetscInt       nv;             /* size of current Schur decomposition */
 77:   PetscReal      errest;         /* error estimate */
 78:   PetscReal      bnorm;          /* computed norm of right-hand side in current solve */
 79:   PetscInt       setupcalled;
 80:   MFNConvergedReason reason;
 81: };

 85: /*
 86:    MFN_CreateDenseMat - Creates a dense Mat of size k unless it already has that size
 87: */
 88: PETSC_STATIC_INLINE PetscErrorCode MFN_CreateDenseMat(PetscInt k,Mat *A)
 89: {
 91:   PetscBool      create=PETSC_FALSE;
 92:   PetscInt       m,n;

 95:   if (!*A) create=PETSC_TRUE;
 96:   else {
 97:     MatGetSize(*A,&m,&n);
 98:     if (m!=k || n!=k) {
 99:       MatDestroy(A);
100:       create=PETSC_TRUE;
101:     }
102:   }
103:   if (create) {
104:     MatCreateSeqDense(PETSC_COMM_SELF,k,k,NULL,A);
105:   }
106:   return(0);
107: }

111: /*
112:    MFN_CreateVec - Creates a Vec of size k unless it already has that size
113: */
114: PETSC_STATIC_INLINE PetscErrorCode MFN_CreateVec(PetscInt k,Vec *v)
115: {
117:   PetscBool      create=PETSC_FALSE;
118:   PetscInt       n;

121:   if (!*v) create=PETSC_TRUE;
122:   else {
123:     VecGetSize(*v,&n);
124:     if (n!=k) {
125:       VecDestroy(v);
126:       create=PETSC_TRUE;
127:     }
128:   }
129:   if (create) {
130:     VecCreateSeq(PETSC_COMM_SELF,k,v);
131:   }
132:   return(0);
133: }

135: PETSC_INTERN PetscErrorCode MFNBasicArnoldi(MFN,PetscScalar*,PetscInt,PetscInt,PetscInt*,PetscReal*,PetscBool*);

137: #endif