Actual source code: nepdefault.c

slepc-3.6.1 2015-09-03
Report Typos and Errors
  1: /*
  2:      This file contains some simple default routines for common NEP operations.

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2015, 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: */

 24: #include <slepc/private/nepimpl.h>     /*I "slepcnep.h" I*/

 28: /*@
 29:    NEPSetWorkVecs - Sets a number of work vectors into a NEP object

 31:    Collective on NEP

 33:    Input Parameters:
 34: +  nep - nonlinear eigensolver context
 35: -  nw  - number of work vectors to allocate

 37:    Developers Note:
 38:    This is PETSC_EXTERN because it may be required by user plugin NEP
 39:    implementations.

 41:    Level: developer
 42: @*/
 43: PetscErrorCode NEPSetWorkVecs(NEP nep,PetscInt nw)
 44: {
 46:   Vec            t;

 49:   if (nep->nwork < nw) {
 50:     VecDestroyVecs(nep->nwork,&nep->work);
 51:     nep->nwork = nw;
 52:     BVGetColumn(nep->V,0,&t);
 53:     VecDuplicateVecs(t,nw,&nep->work);
 54:     BVRestoreColumn(nep->V,0,&t);
 55:     PetscLogObjectParents(nep,nw,nep->work);
 56:   }
 57:   return(0);
 58: }

 62: /*
 63:   NEPGetDefaultShift - Return the value of sigma to start the nonlinear iteration.
 64:  */
 65: PetscErrorCode NEPGetDefaultShift(NEP nep,PetscScalar *sigma)
 66: {
 69:   switch (nep->which) {
 70:     case NEP_LARGEST_MAGNITUDE:
 71:     case NEP_LARGEST_IMAGINARY:
 72:       *sigma = 1.0;   /* arbitrary value */
 73:       break;
 74:     case NEP_SMALLEST_MAGNITUDE:
 75:     case NEP_SMALLEST_IMAGINARY:
 76:       *sigma = 0.0;
 77:       break;
 78:     case NEP_LARGEST_REAL:
 79:       *sigma = PETSC_MAX_REAL;
 80:       break;
 81:     case NEP_SMALLEST_REAL:
 82:       *sigma = PETSC_MIN_REAL;
 83:       break;
 84:     case NEP_TARGET_MAGNITUDE:
 85:     case NEP_TARGET_REAL:
 86:     case NEP_TARGET_IMAGINARY:
 87:       *sigma = nep->target;
 88:       break;
 89:   }
 90:   return(0);
 91: }

 95: /*
 96:   NEPConvergedDefault - Checks convergence of the nonlinear eigensolver.
 97: */
 98: PetscErrorCode NEPConvergedDefault(NEP nep,PetscInt it,PetscReal xnorm,PetscReal snorm,PetscReal fnorm,NEPConvergedReason *reason,void *ctx)
 99: {


106:   *reason = NEP_CONVERGED_ITERATING;

108:   if (!it) {
109:     /* set parameter for default relative tolerance convergence test */
110:     nep->ttol = fnorm*nep->rtol;
111:   }
112:   if (PetscIsInfOrNanReal(fnorm)) {
113:     PetscInfo(nep,"Failed to converged, function norm is NaN\n");
114:     *reason = NEP_DIVERGED_FNORM_NAN;
115:   } else if (fnorm < nep->abstol) {
116:     PetscInfo2(nep,"Converged due to function norm %14.12e < %14.12e\n",(double)fnorm,(double)nep->abstol);
117:     *reason = NEP_CONVERGED_FNORM_ABS;
118:   } else if (nep->nfuncs >= nep->max_funcs) {
119:     PetscInfo2(nep,"Exceeded maximum number of function evaluations: %D > %D\n",nep->nfuncs,nep->max_funcs);
120:     *reason = NEP_DIVERGED_FUNCTION_COUNT;
121:   }

123:   if (it && !*reason) {
124:     if (fnorm <= nep->ttol) {
125:       PetscInfo2(nep,"Converged due to function norm %14.12e < %14.12e (relative tolerance)\n",(double)fnorm,(double)nep->ttol);
126:       *reason = NEP_CONVERGED_FNORM_RELATIVE;
127:     } else if (snorm < nep->stol*xnorm) {
128:       PetscInfo3(nep,"Converged due to small update length: %14.12e < %14.12e * %14.12e\n",(double)snorm,(double)nep->stol,(double)xnorm);
129:       *reason = NEP_CONVERGED_SNORM_RELATIVE;
130:     }
131:   }
132:   return(0);
133: }