Actual source code: pjdopt.c

slepc-3.7.3 2016-09-29
Report Typos and Errors
  1: /*
  2:    Options of polynomial JD solver.

  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: */

 24: #include <slepc/private/pepimpl.h>    /*I "slepcpep.h" I*/
 25:  #include pjdp.h

 29: PetscErrorCode PEPJDSetRestart_JD(PEP pep,PetscReal keep)
 30: {
 31:   PEP_JD *pjd = (PEP_JD*)pep->data;

 34:   if (keep==PETSC_DEFAULT) pjd->keep = 0.5;
 35:   else {
 36:     if (keep<0.1 || keep>0.9) SETERRQ(PetscObjectComm((PetscObject)pep),PETSC_ERR_ARG_OUTOFRANGE,"The keep argument must be in the range [0.1,0.9]");
 37:     pjd->keep = keep;
 38:   }
 39:   return(0);
 40: }

 44: /*@
 45:    PEPJDSetRestart - Sets the restart parameter for the Jacobi-Davidson
 46:    method, in particular the proportion of basis vectors that must be kept
 47:    after restart.

 49:    Logically Collective on PEP

 51:    Input Parameters:
 52: +  pep  - the eigenproblem solver context
 53: -  keep - the number of vectors to be kept at restart

 55:    Options Database Key:
 56: .  -pep_jd_restart - Sets the restart parameter

 58:    Notes:
 59:    Allowed values are in the range [0.1,0.9]. The default is 0.5.

 61:    Level: advanced

 63: .seealso: PEPJDGetRestart()
 64: @*/
 65: PetscErrorCode PEPJDSetRestart(PEP pep,PetscReal keep)
 66: {

 72:   PetscTryMethod(pep,"PEPJDSetRestart_C",(PEP,PetscReal),(pep,keep));
 73:   return(0);
 74: }

 78: PetscErrorCode PEPJDGetRestart_JD(PEP pep,PetscReal *keep)
 79: {
 80:   PEP_JD *pjd = (PEP_JD*)pep->data;

 83:   *keep = pjd->keep;
 84:   return(0);
 85: }

 89: /*@
 90:    PEPJDGetRestart - Gets the restart parameter used in the Jacobi-Davidson method.

 92:    Not Collective

 94:    Input Parameter:
 95: .  pep - the eigenproblem solver context

 97:    Output Parameter:
 98: .  keep - the restart parameter

100:    Level: advanced

102: .seealso: PEPJDSetRestart()
103: @*/
104: PetscErrorCode PEPJDGetRestart(PEP pep,PetscReal *keep)
105: {

111:   PetscUseMethod(pep,"PEPJDGetRestart_C",(PEP,PetscReal*),(pep,keep));
112:   return(0);
113: }

117: PetscErrorCode PEPSetFromOptions_JD(PetscOptionItems *PetscOptionsObject,PEP pep)
118: {
120:   PetscBool      flg;
121:   PetscReal      r1;
122:   KSP            ksp;

125:   PetscOptionsHead(PetscOptionsObject,"PEP JD Options");
126:   PetscOptionsReal("-pep_jd_restart","Proportion of vectors kept after restart","PEPJDSetRestart",0.5,&r1,&flg);
127:   if (flg) {
128:     PEPJDSetRestart(pep,r1);
129:   }
130:   /* Set STPRECOND as the default ST */
131:   if (!pep->st) { PEPGetST(pep,&pep->st); }
132:   if (!((PetscObject)pep->st)->type_name) {
133:     STSetType(pep->st,STPRECOND);
134:   }

136:   /* Set the default options of the KSP */
137:   STGetKSP(pep->st,&ksp);
138:   if (!((PetscObject)ksp)->type_name) {
139:     KSPSetType(ksp,KSPBCGSL);
140:     KSPSetTolerances(ksp,1e-5,PETSC_DEFAULT,PETSC_DEFAULT,100);
141:   }
142:   PetscOptionsTail();
143:   return(0);
144: }

148: PetscErrorCode PEPView_JD(PEP pep,PetscViewer viewer)
149: {
151:   PEP_JD         *pjd = (PEP_JD*)pep->data;
152:   PetscBool      isascii;

155:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
156:   if (isascii) {
157:     PetscViewerASCIIPrintf(viewer,"  JD: %d%% of basis vectors kept after restart\n",(int)(100*pjd->keep));
158:   }
159:   return(0);
160: }