LCOV - code coverage report
Current view: top level - nep/interface - nepopts.c (source / functions) Hit Total Coverage
Test: SLEPc Lines: 368 393 93.6 %
Date: 2021-08-02 00:35:43 Functions: 25 26 96.2 %
Legend: Lines: hit not hit

          Line data    Source code
       1             : /*
       2             :    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
       3             :    SLEPc - Scalable Library for Eigenvalue Problem Computations
       4             :    Copyright (c) 2002-2021, Universitat Politecnica de Valencia, Spain
       5             : 
       6             :    This file is part of SLEPc.
       7             :    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
       8             :    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
       9             : */
      10             : /*
      11             :    NEP routines related to options that can be set via the command-line
      12             :    or procedurally
      13             : */
      14             : 
      15             : #include <slepc/private/nepimpl.h>       /*I "slepcnep.h" I*/
      16             : #include <petscdraw.h>
      17             : 
      18             : /*@C
      19             :    NEPMonitorSetFromOptions - Sets a monitor function and viewer appropriate for the type
      20             :    indicated by the user.
      21             : 
      22             :    Collective on nep
      23             : 
      24             :    Input Parameters:
      25             : +  nep      - the nonlinear eigensolver context
      26             : .  opt      - the command line option for this monitor
      27             : .  name     - the monitor type one is seeking
      28             : .  ctx      - an optional user context for the monitor, or NULL
      29             : -  trackall - whether this monitor tracks all eigenvalues or not
      30             : 
      31             :    Level: developer
      32             : 
      33             : .seealso: NEPMonitorSet(), NEPSetTrackAll()
      34             : @*/
      35         330 : PetscErrorCode NEPMonitorSetFromOptions(NEP nep,const char opt[],const char name[],void *ctx,PetscBool trackall)
      36             : {
      37         330 :   PetscErrorCode       (*mfunc)(NEP,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
      38         330 :   PetscErrorCode       (*cfunc)(PetscViewer,PetscViewerFormat,void*,PetscViewerAndFormat**);
      39         330 :   PetscErrorCode       (*dfunc)(PetscViewerAndFormat**);
      40         330 :   PetscViewerAndFormat *vf;
      41         330 :   PetscViewer          viewer;
      42         330 :   PetscViewerFormat    format;
      43         330 :   PetscViewerType      vtype;
      44         330 :   char                 key[PETSC_MAX_PATH_LEN];
      45         330 :   PetscBool            flg;
      46         330 :   PetscErrorCode       ierr;
      47             : 
      48         330 :   PetscFunctionBegin;
      49         330 :   ierr = PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,opt,&viewer,&format,&flg);CHKERRQ(ierr);
      50         330 :   if (!flg) PetscFunctionReturn(0);
      51             : 
      52           6 :   ierr = PetscViewerGetType(viewer,&vtype);CHKERRQ(ierr);
      53           6 :   ierr = SlepcMonitorMakeKey_Internal(name,vtype,format,key);CHKERRQ(ierr);
      54           6 :   ierr = PetscFunctionListFind(NEPMonitorList,key,&mfunc);CHKERRQ(ierr);
      55           6 :   ierr = PetscFunctionListFind(NEPMonitorCreateList,key,&cfunc);CHKERRQ(ierr);
      56           6 :   ierr = PetscFunctionListFind(NEPMonitorDestroyList,key,&dfunc);CHKERRQ(ierr);
      57           6 :   if (!cfunc) cfunc = PetscViewerAndFormatCreate_Internal;
      58           6 :   if (!dfunc) dfunc = PetscViewerAndFormatDestroy;
      59             : 
      60           6 :   ierr = (*cfunc)(viewer,format,ctx,&vf);CHKERRQ(ierr);
      61           6 :   ierr = PetscObjectDereference((PetscObject)viewer);CHKERRQ(ierr);
      62           6 :   ierr = NEPMonitorSet(nep,mfunc,vf,(PetscErrorCode(*)(void **))dfunc);CHKERRQ(ierr);
      63           6 :   if (trackall) {
      64           2 :     ierr = NEPSetTrackAll(nep,PETSC_TRUE);CHKERRQ(ierr);
      65             :   }
      66           6 :   PetscFunctionReturn(0);
      67             : }
      68             : 
      69             : /*@
      70             :    NEPSetFromOptions - Sets NEP options from the options database.
      71             :    This routine must be called before NEPSetUp() if the user is to be
      72             :    allowed to set the solver type.
      73             : 
      74             :    Collective on nep
      75             : 
      76             :    Input Parameters:
      77             : .  nep - the nonlinear eigensolver context
      78             : 
      79             :    Notes:
      80             :    To see all options, run your program with the -help option.
      81             : 
      82             :    Level: beginner
      83             : @*/
      84         110 : PetscErrorCode NEPSetFromOptions(NEP nep)
      85             : {
      86         110 :   PetscErrorCode  ierr;
      87         110 :   char            type[256];
      88         110 :   PetscBool       set,flg,flg1,flg2,flg3,flg4,flg5,bval;
      89         110 :   PetscReal       r;
      90         110 :   PetscScalar     s;
      91         110 :   PetscInt        i,j,k;
      92         110 :   NEPRefine       refine;
      93         110 :   NEPRefineScheme scheme;
      94             : 
      95         110 :   PetscFunctionBegin;
      96         110 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
      97         110 :   ierr = NEPRegisterAll();CHKERRQ(ierr);
      98         330 :   ierr = PetscObjectOptionsBegin((PetscObject)nep);CHKERRQ(ierr);
      99         151 :     ierr = PetscOptionsFList("-nep_type","Nonlinear eigensolver method","NEPSetType",NEPList,(char*)(((PetscObject)nep)->type_name?((PetscObject)nep)->type_name:NEPRII),type,sizeof(type),&flg);CHKERRQ(ierr);
     100         110 :     if (flg) {
     101          78 :       ierr = NEPSetType(nep,type);CHKERRQ(ierr);
     102          32 :     } else if (!((PetscObject)nep)->type_name) {
     103           2 :       ierr = NEPSetType(nep,NEPRII);CHKERRQ(ierr);
     104             :     }
     105             : 
     106         110 :     ierr = PetscOptionsBoolGroupBegin("-nep_general","General nonlinear eigenvalue problem","NEPSetProblemType",&flg);CHKERRQ(ierr);
     107         110 :     if (flg) { ierr = NEPSetProblemType(nep,NEP_GENERAL);CHKERRQ(ierr); }
     108         110 :     ierr = PetscOptionsBoolGroupEnd("-nep_rational","Rational eigenvalue problem","NEPSetProblemType",&flg);CHKERRQ(ierr);
     109         110 :     if (flg) { ierr = NEPSetProblemType(nep,NEP_RATIONAL);CHKERRQ(ierr); }
     110             : 
     111         110 :     refine = nep->refine;
     112         110 :     ierr = PetscOptionsEnum("-nep_refine","Iterative refinement method","NEPSetRefine",NEPRefineTypes,(PetscEnum)refine,(PetscEnum*)&refine,&flg1);CHKERRQ(ierr);
     113         110 :     i = nep->npart;
     114         110 :     ierr = PetscOptionsInt("-nep_refine_partitions","Number of partitions of the communicator for iterative refinement","NEPSetRefine",nep->npart,&i,&flg2);CHKERRQ(ierr);
     115         110 :     r = nep->rtol;
     116         111 :     ierr = PetscOptionsReal("-nep_refine_tol","Tolerance for iterative refinement","NEPSetRefine",nep->rtol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL/1000:nep->rtol,&r,&flg3);CHKERRQ(ierr);
     117         110 :     j = nep->rits;
     118         110 :     ierr = PetscOptionsInt("-nep_refine_its","Maximum number of iterations for iterative refinement","NEPSetRefine",nep->rits,&j,&flg4);CHKERRQ(ierr);
     119         110 :     scheme = nep->scheme;
     120         110 :     ierr = PetscOptionsEnum("-nep_refine_scheme","Scheme used for linear systems within iterative refinement","NEPSetRefine",NEPRefineSchemes,(PetscEnum)scheme,(PetscEnum*)&scheme,&flg5);CHKERRQ(ierr);
     121         110 :     if (flg1 || flg2 || flg3 || flg4 || flg5) { ierr = NEPSetRefine(nep,refine,i,r,j,scheme);CHKERRQ(ierr); }
     122             : 
     123         110 :     i = nep->max_it;
     124         110 :     ierr = PetscOptionsInt("-nep_max_it","Maximum number of iterations","NEPSetTolerances",nep->max_it,&i,&flg1);CHKERRQ(ierr);
     125         110 :     r = nep->tol;
     126         160 :     ierr = PetscOptionsReal("-nep_tol","Tolerance","NEPSetTolerances",SlepcDefaultTol(nep->tol),&r,&flg2);CHKERRQ(ierr);
     127         110 :     if (flg1 || flg2) { ierr = NEPSetTolerances(nep,r,i);CHKERRQ(ierr); }
     128             : 
     129         110 :     ierr = PetscOptionsBoolGroupBegin("-nep_conv_rel","Relative error convergence test","NEPSetConvergenceTest",&flg);CHKERRQ(ierr);
     130         110 :     if (flg) { ierr = NEPSetConvergenceTest(nep,NEP_CONV_REL);CHKERRQ(ierr); }
     131         110 :     ierr = PetscOptionsBoolGroup("-nep_conv_norm","Convergence test relative to the matrix norms","NEPSetConvergenceTest",&flg);CHKERRQ(ierr);
     132         110 :     if (flg) { ierr = NEPSetConvergenceTest(nep,NEP_CONV_NORM);CHKERRQ(ierr); }
     133         110 :     ierr = PetscOptionsBoolGroup("-nep_conv_abs","Absolute error convergence test","NEPSetConvergenceTest",&flg);CHKERRQ(ierr);
     134         110 :     if (flg) { ierr = NEPSetConvergenceTest(nep,NEP_CONV_ABS);CHKERRQ(ierr); }
     135         110 :     ierr = PetscOptionsBoolGroupEnd("-nep_conv_user","User-defined convergence test","NEPSetConvergenceTest",&flg);CHKERRQ(ierr);
     136         110 :     if (flg) { ierr = NEPSetConvergenceTest(nep,NEP_CONV_USER);CHKERRQ(ierr); }
     137             : 
     138         110 :     ierr = PetscOptionsBoolGroupBegin("-nep_stop_basic","Stop iteration if all eigenvalues converged or max_it reached","NEPSetStoppingTest",&flg);CHKERRQ(ierr);
     139         110 :     if (flg) { ierr = NEPSetStoppingTest(nep,NEP_STOP_BASIC);CHKERRQ(ierr); }
     140         110 :     ierr = PetscOptionsBoolGroupEnd("-nep_stop_user","User-defined stopping test","NEPSetStoppingTest",&flg);CHKERRQ(ierr);
     141         110 :     if (flg) { ierr = NEPSetStoppingTest(nep,NEP_STOP_USER);CHKERRQ(ierr); }
     142             : 
     143         110 :     i = nep->nev;
     144         110 :     ierr = PetscOptionsInt("-nep_nev","Number of eigenvalues to compute","NEPSetDimensions",nep->nev,&i,&flg1);CHKERRQ(ierr);
     145         110 :     j = nep->ncv;
     146         110 :     ierr = PetscOptionsInt("-nep_ncv","Number of basis vectors","NEPSetDimensions",nep->ncv,&j,&flg2);CHKERRQ(ierr);
     147         110 :     k = nep->mpd;
     148         110 :     ierr = PetscOptionsInt("-nep_mpd","Maximum dimension of projected problem","NEPSetDimensions",nep->mpd,&k,&flg3);CHKERRQ(ierr);
     149         110 :     if (flg1 || flg2 || flg3) {
     150          57 :       ierr = NEPSetDimensions(nep,i,j,k);CHKERRQ(ierr);
     151             :     }
     152             : 
     153         110 :     ierr = PetscOptionsBoolGroupBegin("-nep_largest_magnitude","Compute largest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     154         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_LARGEST_MAGNITUDE);CHKERRQ(ierr); }
     155         110 :     ierr = PetscOptionsBoolGroup("-nep_smallest_magnitude","Compute smallest eigenvalues in magnitude","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     156         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_SMALLEST_MAGNITUDE);CHKERRQ(ierr); }
     157         110 :     ierr = PetscOptionsBoolGroup("-nep_largest_real","Compute eigenvalues with largest real parts","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     158         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_LARGEST_REAL);CHKERRQ(ierr); }
     159         110 :     ierr = PetscOptionsBoolGroup("-nep_smallest_real","Compute eigenvalues with smallest real parts","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     160         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_SMALLEST_REAL);CHKERRQ(ierr); }
     161         110 :     ierr = PetscOptionsBoolGroup("-nep_largest_imaginary","Compute eigenvalues with largest imaginary parts","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     162         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_LARGEST_IMAGINARY);CHKERRQ(ierr); }
     163         110 :     ierr = PetscOptionsBoolGroup("-nep_smallest_imaginary","Compute eigenvalues with smallest imaginary parts","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     164         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_SMALLEST_IMAGINARY);CHKERRQ(ierr); }
     165         110 :     ierr = PetscOptionsBoolGroup("-nep_target_magnitude","Compute eigenvalues closest to target","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     166         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE);CHKERRQ(ierr); }
     167         110 :     ierr = PetscOptionsBoolGroup("-nep_target_real","Compute eigenvalues with real parts closest to target","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     168         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_TARGET_REAL);CHKERRQ(ierr); }
     169         110 :     ierr = PetscOptionsBoolGroup("-nep_target_imaginary","Compute eigenvalues with imaginary parts closest to target","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     170         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_TARGET_IMAGINARY);CHKERRQ(ierr); }
     171         110 :     ierr = PetscOptionsBoolGroupEnd("-nep_all","Compute all eigenvalues in a region","NEPSetWhichEigenpairs",&flg);CHKERRQ(ierr);
     172         110 :     if (flg) { ierr = NEPSetWhichEigenpairs(nep,NEP_ALL);CHKERRQ(ierr); }
     173             : 
     174         110 :     ierr = PetscOptionsScalar("-nep_target","Value of the target","NEPSetTarget",nep->target,&s,&flg);CHKERRQ(ierr);
     175         110 :     if (flg) {
     176          62 :       if (nep->which!=NEP_TARGET_REAL && nep->which!=NEP_TARGET_IMAGINARY) {
     177          61 :         ierr = NEPSetWhichEigenpairs(nep,NEP_TARGET_MAGNITUDE);CHKERRQ(ierr);
     178             :       }
     179          62 :       ierr = NEPSetTarget(nep,s);CHKERRQ(ierr);
     180             :     }
     181             : 
     182         110 :     ierr = PetscOptionsBool("-nep_two_sided","Use two-sided variant (to compute left eigenvectors)","NEPSetTwoSided",nep->twosided,&bval,&flg);CHKERRQ(ierr);
     183         110 :     if (flg) { ierr = NEPSetTwoSided(nep,bval);CHKERRQ(ierr); }
     184             : 
     185             :     /* -----------------------------------------------------------------------*/
     186             :     /*
     187             :       Cancels all monitors hardwired into code before call to NEPSetFromOptions()
     188             :     */
     189         110 :     ierr = PetscOptionsBool("-nep_monitor_cancel","Remove any hardwired monitor routines","NEPMonitorCancel",PETSC_FALSE,&flg,&set);CHKERRQ(ierr);
     190         110 :     if (set && flg) { ierr = NEPMonitorCancel(nep);CHKERRQ(ierr); }
     191         110 :     ierr = NEPMonitorSetFromOptions(nep,"-nep_monitor","first_approximation",NULL,PETSC_FALSE);CHKERRQ(ierr);
     192         110 :     ierr = NEPMonitorSetFromOptions(nep,"-nep_monitor_all","all_approximations",NULL,PETSC_TRUE);CHKERRQ(ierr);
     193         110 :     ierr = NEPMonitorSetFromOptions(nep,"-nep_monitor_conv","convergence_history",NULL,PETSC_FALSE);CHKERRQ(ierr);
     194             : 
     195             :     /* -----------------------------------------------------------------------*/
     196         110 :     ierr = PetscOptionsName("-nep_view","Print detailed information on solver used","NEPView",NULL);CHKERRQ(ierr);
     197         110 :     ierr = PetscOptionsName("-nep_view_vectors","View computed eigenvectors","NEPVectorsView",NULL);CHKERRQ(ierr);
     198         110 :     ierr = PetscOptionsName("-nep_view_values","View computed eigenvalues","NEPValuesView",NULL);CHKERRQ(ierr);
     199         110 :     ierr = PetscOptionsName("-nep_converged_reason","Print reason for convergence, and number of iterations","NEPConvergedReasonView",NULL);CHKERRQ(ierr);
     200         110 :     ierr = PetscOptionsName("-nep_error_absolute","Print absolute errors of each eigenpair","NEPErrorView",NULL);CHKERRQ(ierr);
     201         110 :     ierr = PetscOptionsName("-nep_error_relative","Print relative errors of each eigenpair","NEPErrorView",NULL);CHKERRQ(ierr);
     202             : 
     203         110 :     if (nep->ops->setfromoptions) {
     204         110 :       ierr = (*nep->ops->setfromoptions)(PetscOptionsObject,nep);CHKERRQ(ierr);
     205             :     }
     206         110 :     ierr = PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject)nep);CHKERRQ(ierr);
     207         220 :   ierr = PetscOptionsEnd();CHKERRQ(ierr);
     208             : 
     209         110 :   if (!nep->V) { ierr = NEPGetBV(nep,&nep->V);CHKERRQ(ierr); }
     210         110 :   ierr = BVSetFromOptions(nep->V);CHKERRQ(ierr);
     211         110 :   if (!nep->rg) { ierr = NEPGetRG(nep,&nep->rg);CHKERRQ(ierr); }
     212         110 :   ierr = RGSetFromOptions(nep->rg);CHKERRQ(ierr);
     213         110 :   if (nep->useds) {
     214          86 :     if (!nep->ds) { ierr = NEPGetDS(nep,&nep->ds);CHKERRQ(ierr); }
     215          86 :     ierr = DSSetFromOptions(nep->ds);CHKERRQ(ierr);
     216             :   }
     217         110 :   if (!nep->refineksp) { ierr = NEPRefineGetKSP(nep,&nep->refineksp);CHKERRQ(ierr); }
     218         110 :   ierr = KSPSetFromOptions(nep->refineksp);CHKERRQ(ierr);
     219         372 :   if (nep->fui==NEP_USER_INTERFACE_SPLIT) for (i=0;i<nep->nt;i++) {ierr = FNSetFromOptions(nep->f[i]);CHKERRQ(ierr);}
     220         110 :   PetscFunctionReturn(0);
     221             : }
     222             : 
     223             : /*@C
     224             :    NEPGetTolerances - Gets the tolerance and maximum iteration count used
     225             :    by the NEP convergence tests.
     226             : 
     227             :    Not Collective
     228             : 
     229             :    Input Parameter:
     230             : .  nep - the nonlinear eigensolver context
     231             : 
     232             :    Output Parameters:
     233             : +  tol - the convergence tolerance
     234             : -  maxits - maximum number of iterations
     235             : 
     236             :    Notes:
     237             :    The user can specify NULL for any parameter that is not needed.
     238             : 
     239             :    Level: intermediate
     240             : 
     241             : .seealso: NEPSetTolerances()
     242             : @*/
     243           7 : PetscErrorCode NEPGetTolerances(NEP nep,PetscReal *tol,PetscInt *maxits)
     244             : {
     245           7 :   PetscFunctionBegin;
     246           7 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     247           7 :   if (tol)    *tol    = nep->tol;
     248           7 :   if (maxits) *maxits = nep->max_it;
     249           7 :   PetscFunctionReturn(0);
     250             : }
     251             : 
     252             : /*@
     253             :    NEPSetTolerances - Sets the tolerance and maximum iteration count used
     254             :    by the NEP convergence tests.
     255             : 
     256             :    Logically Collective on nep
     257             : 
     258             :    Input Parameters:
     259             : +  nep    - the nonlinear eigensolver context
     260             : .  tol    - the convergence tolerance
     261             : -  maxits - maximum number of iterations to use
     262             : 
     263             :    Options Database Keys:
     264             : +  -nep_tol <tol> - Sets the convergence tolerance
     265             : -  -nep_max_it <maxits> - Sets the maximum number of iterations allowed
     266             : 
     267             :    Notes:
     268             :    Use PETSC_DEFAULT for either argument to assign a reasonably good value.
     269             : 
     270             :    Level: intermediate
     271             : 
     272             : .seealso: NEPGetTolerances()
     273             : @*/
     274          61 : PetscErrorCode NEPSetTolerances(NEP nep,PetscReal tol,PetscInt maxits)
     275             : {
     276          61 :   PetscFunctionBegin;
     277          61 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     278         122 :   PetscValidLogicalCollectiveReal(nep,tol,2);
     279         122 :   PetscValidLogicalCollectiveInt(nep,maxits,3);
     280          61 :   if (tol == PETSC_DEFAULT) {
     281           1 :     nep->tol   = PETSC_DEFAULT;
     282           1 :     nep->state = NEP_STATE_INITIAL;
     283             :   } else {
     284          60 :     if (tol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
     285          60 :     nep->tol = tol;
     286             :   }
     287          61 :   if (maxits == PETSC_DEFAULT || maxits == PETSC_DECIDE) {
     288          59 :     nep->max_it = PETSC_DEFAULT;
     289          59 :     nep->state  = NEP_STATE_INITIAL;
     290             :   } else {
     291           2 :     if (maxits <= 0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of maxits. Must be > 0");
     292           2 :     nep->max_it = maxits;
     293             :   }
     294          61 :   PetscFunctionReturn(0);
     295             : }
     296             : 
     297             : /*@C
     298             :    NEPGetDimensions - Gets the number of eigenvalues to compute
     299             :    and the dimension of the subspace.
     300             : 
     301             :    Not Collective
     302             : 
     303             :    Input Parameter:
     304             : .  nep - the nonlinear eigensolver context
     305             : 
     306             :    Output Parameters:
     307             : +  nev - number of eigenvalues to compute
     308             : .  ncv - the maximum dimension of the subspace to be used by the solver
     309             : -  mpd - the maximum dimension allowed for the projected problem
     310             : 
     311             :    Notes:
     312             :    The user can specify NULL for any parameter that is not needed.
     313             : 
     314             :    Level: intermediate
     315             : 
     316             : .seealso: NEPSetDimensions()
     317             : @*/
     318          45 : PetscErrorCode NEPGetDimensions(NEP nep,PetscInt *nev,PetscInt *ncv,PetscInt *mpd)
     319             : {
     320          45 :   PetscFunctionBegin;
     321          45 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     322          45 :   if (nev) *nev = nep->nev;
     323          45 :   if (ncv) *ncv = nep->ncv;
     324          45 :   if (mpd) *mpd = nep->mpd;
     325          45 :   PetscFunctionReturn(0);
     326             : }
     327             : 
     328             : /*@
     329             :    NEPSetDimensions - Sets the number of eigenvalues to compute
     330             :    and the dimension of the subspace.
     331             : 
     332             :    Logically Collective on nep
     333             : 
     334             :    Input Parameters:
     335             : +  nep - the nonlinear eigensolver context
     336             : .  nev - number of eigenvalues to compute
     337             : .  ncv - the maximum dimension of the subspace to be used by the solver
     338             : -  mpd - the maximum dimension allowed for the projected problem
     339             : 
     340             :    Options Database Keys:
     341             : +  -nep_nev <nev> - Sets the number of eigenvalues
     342             : .  -nep_ncv <ncv> - Sets the dimension of the subspace
     343             : -  -nep_mpd <mpd> - Sets the maximum projected dimension
     344             : 
     345             :    Notes:
     346             :    Use PETSC_DEFAULT for ncv and mpd to assign a reasonably good value, which is
     347             :    dependent on the solution method.
     348             : 
     349             :    The parameters ncv and mpd are intimately related, so that the user is advised
     350             :    to set one of them at most. Normal usage is that
     351             :    (a) in cases where nev is small, the user sets ncv (a reasonable default is 2*nev); and
     352             :    (b) in cases where nev is large, the user sets mpd.
     353             : 
     354             :    The value of ncv should always be between nev and (nev+mpd), typically
     355             :    ncv=nev+mpd. If nev is not too large, mpd=nev is a reasonable choice, otherwise
     356             :    a smaller value should be used.
     357             : 
     358             :    Level: intermediate
     359             : 
     360             : .seealso: NEPGetDimensions()
     361             : @*/
     362          95 : PetscErrorCode NEPSetDimensions(NEP nep,PetscInt nev,PetscInt ncv,PetscInt mpd)
     363             : {
     364          95 :   PetscFunctionBegin;
     365          95 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     366         190 :   PetscValidLogicalCollectiveInt(nep,nev,2);
     367         190 :   PetscValidLogicalCollectiveInt(nep,ncv,3);
     368         190 :   PetscValidLogicalCollectiveInt(nep,mpd,4);
     369          95 :   if (nev<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of nev. Must be > 0");
     370          95 :   nep->nev = nev;
     371          95 :   if (ncv == PETSC_DECIDE || ncv == PETSC_DEFAULT) {
     372          74 :     nep->ncv = PETSC_DEFAULT;
     373             :   } else {
     374          21 :     if (ncv<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of ncv. Must be > 0");
     375          21 :     nep->ncv = ncv;
     376             :   }
     377          95 :   if (mpd == PETSC_DECIDE || mpd == PETSC_DEFAULT) {
     378          94 :     nep->mpd = PETSC_DEFAULT;
     379             :   } else {
     380           1 :     if (mpd<1) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of mpd. Must be > 0");
     381           1 :     nep->mpd = mpd;
     382             :   }
     383          95 :   nep->state = NEP_STATE_INITIAL;
     384          95 :   PetscFunctionReturn(0);
     385             : }
     386             : 
     387             : /*@
     388             :     NEPSetWhichEigenpairs - Specifies which portion of the spectrum is
     389             :     to be sought.
     390             : 
     391             :     Logically Collective on nep
     392             : 
     393             :     Input Parameters:
     394             : +   nep   - eigensolver context obtained from NEPCreate()
     395             : -   which - the portion of the spectrum to be sought
     396             : 
     397             :     Possible values:
     398             :     The parameter 'which' can have one of these values
     399             : 
     400             : +     NEP_LARGEST_MAGNITUDE - largest eigenvalues in magnitude (default)
     401             : .     NEP_SMALLEST_MAGNITUDE - smallest eigenvalues in magnitude
     402             : .     NEP_LARGEST_REAL - largest real parts
     403             : .     NEP_SMALLEST_REAL - smallest real parts
     404             : .     NEP_LARGEST_IMAGINARY - largest imaginary parts
     405             : .     NEP_SMALLEST_IMAGINARY - smallest imaginary parts
     406             : .     NEP_TARGET_MAGNITUDE - eigenvalues closest to the target (in magnitude)
     407             : .     NEP_TARGET_REAL - eigenvalues with real part closest to target
     408             : .     NEP_TARGET_IMAGINARY - eigenvalues with imaginary part closest to target
     409             : .     NEP_ALL - all eigenvalues contained in a given region
     410             : -     NEP_WHICH_USER - user defined ordering set with NEPSetEigenvalueComparison()
     411             : 
     412             :     Options Database Keys:
     413             : +   -nep_largest_magnitude - Sets largest eigenvalues in magnitude
     414             : .   -nep_smallest_magnitude - Sets smallest eigenvalues in magnitude
     415             : .   -nep_largest_real - Sets largest real parts
     416             : .   -nep_smallest_real - Sets smallest real parts
     417             : .   -nep_largest_imaginary - Sets largest imaginary parts
     418             : .   -nep_smallest_imaginary - Sets smallest imaginary parts
     419             : .   -nep_target_magnitude - Sets eigenvalues closest to target
     420             : .   -nep_target_real - Sets real parts closest to target
     421             : .   -nep_target_imaginary - Sets imaginary parts closest to target
     422             : -   -nep_all - Sets all eigenvalues in a region
     423             : 
     424             :     Notes:
     425             :     Not all eigensolvers implemented in NEP account for all the possible values
     426             :     stated above. If SLEPc is compiled for real numbers NEP_LARGEST_IMAGINARY
     427             :     and NEP_SMALLEST_IMAGINARY use the absolute value of the imaginary part
     428             :     for eigenvalue selection.
     429             : 
     430             :     The target is a scalar value provided with NEPSetTarget().
     431             : 
     432             :     NEP_ALL is intended for use in the context of the CISS solver for
     433             :     computing all eigenvalues in a region.
     434             : 
     435             :     Level: intermediate
     436             : 
     437             : .seealso: NEPGetWhichEigenpairs(), NEPSetTarget(), NEPSetEigenvalueComparison(), NEPWhich
     438             : @*/
     439          68 : PetscErrorCode NEPSetWhichEigenpairs(NEP nep,NEPWhich which)
     440             : {
     441          68 :   PetscFunctionBegin;
     442          68 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     443         136 :   PetscValidLogicalCollectiveEnum(nep,which,2);
     444          68 :   switch (which) {
     445          68 :     case NEP_LARGEST_MAGNITUDE:
     446             :     case NEP_SMALLEST_MAGNITUDE:
     447             :     case NEP_LARGEST_REAL:
     448             :     case NEP_SMALLEST_REAL:
     449             :     case NEP_LARGEST_IMAGINARY:
     450             :     case NEP_SMALLEST_IMAGINARY:
     451             :     case NEP_TARGET_MAGNITUDE:
     452             :     case NEP_TARGET_REAL:
     453             : #if defined(PETSC_USE_COMPLEX)
     454             :     case NEP_TARGET_IMAGINARY:
     455             : #endif
     456             :     case NEP_ALL:
     457             :     case NEP_WHICH_USER:
     458          68 :       if (nep->which != which) {
     459          66 :         nep->state = NEP_STATE_INITIAL;
     460          66 :         nep->which = which;
     461             :       }
     462          68 :       break;
     463             : #if !defined(PETSC_USE_COMPLEX)
     464             :     case NEP_TARGET_IMAGINARY:
     465             :       SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"NEP_TARGET_IMAGINARY can be used only with complex scalars");
     466             : #endif
     467           0 :     default:
     468           0 :       SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'which' value");
     469             :   }
     470          68 :   PetscFunctionReturn(0);
     471             : }
     472             : 
     473             : /*@
     474             :     NEPGetWhichEigenpairs - Returns which portion of the spectrum is to be
     475             :     sought.
     476             : 
     477             :     Not Collective
     478             : 
     479             :     Input Parameter:
     480             : .   nep - eigensolver context obtained from NEPCreate()
     481             : 
     482             :     Output Parameter:
     483             : .   which - the portion of the spectrum to be sought
     484             : 
     485             :     Notes:
     486             :     See NEPSetWhichEigenpairs() for possible values of 'which'.
     487             : 
     488             :     Level: intermediate
     489             : 
     490             : .seealso: NEPSetWhichEigenpairs(), NEPWhich
     491             : @*/
     492           1 : PetscErrorCode NEPGetWhichEigenpairs(NEP nep,NEPWhich *which)
     493             : {
     494           1 :   PetscFunctionBegin;
     495           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     496           1 :   PetscValidPointer(which,2);
     497           1 :   *which = nep->which;
     498           1 :   PetscFunctionReturn(0);
     499             : }
     500             : 
     501             : /*@C
     502             :    NEPSetEigenvalueComparison - Specifies the eigenvalue comparison function
     503             :    when NEPSetWhichEigenpairs() is set to NEP_WHICH_USER.
     504             : 
     505             :    Logically Collective on nep
     506             : 
     507             :    Input Parameters:
     508             : +  nep  - eigensolver context obtained from NEPCreate()
     509             : .  func - a pointer to the comparison function
     510             : -  ctx  - a context pointer (the last parameter to the comparison function)
     511             : 
     512             :    Calling Sequence of func:
     513             : $   func(PetscScalar ar,PetscScalar ai,PetscScalar br,PetscScalar bi,PetscInt *res,void *ctx)
     514             : 
     515             : +   ar     - real part of the 1st eigenvalue
     516             : .   ai     - imaginary part of the 1st eigenvalue
     517             : .   br     - real part of the 2nd eigenvalue
     518             : .   bi     - imaginary part of the 2nd eigenvalue
     519             : .   res    - result of comparison
     520             : -   ctx    - optional context, as set by NEPSetEigenvalueComparison()
     521             : 
     522             :    Note:
     523             :    The returning parameter 'res' can be
     524             : +  negative - if the 1st eigenvalue is preferred to the 2st one
     525             : .  zero     - if both eigenvalues are equally preferred
     526             : -  positive - if the 2st eigenvalue is preferred to the 1st one
     527             : 
     528             :    Level: advanced
     529             : 
     530             : .seealso: NEPSetWhichEigenpairs(), NEPWhich
     531             : @*/
     532           0 : PetscErrorCode NEPSetEigenvalueComparison(NEP nep,PetscErrorCode (*func)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*),void* ctx)
     533             : {
     534           0 :   PetscFunctionBegin;
     535           0 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     536           0 :   nep->sc->comparison    = func;
     537           0 :   nep->sc->comparisonctx = ctx;
     538           0 :   nep->which             = NEP_WHICH_USER;
     539           0 :   PetscFunctionReturn(0);
     540             : }
     541             : 
     542             : /*@
     543             :    NEPSetProblemType - Specifies the type of the nonlinear eigenvalue problem.
     544             : 
     545             :    Logically Collective on nep
     546             : 
     547             :    Input Parameters:
     548             : +  nep  - the nonlinear eigensolver context
     549             : -  type - a known type of nonlinear eigenvalue problem
     550             : 
     551             :    Options Database Keys:
     552             : +  -nep_general - general problem with no particular structure
     553             : -  -nep_rational - a rational eigenvalue problem defined in split form with all f_i rational
     554             : 
     555             :    Notes:
     556             :    Allowed values for the problem type are: general (NEP_GENERAL), and rational
     557             :    (NEP_RATIONAL).
     558             : 
     559             :    This function is used to provide a hint to the NEP solver to exploit certain
     560             :    properties of the nonlinear eigenproblem. This hint may be used or not,
     561             :    depending on the solver. By default, no particular structure is assumed.
     562             : 
     563             :    Level: intermediate
     564             : 
     565             : .seealso: NEPSetType(), NEPGetProblemType(), NEPProblemType
     566             : @*/
     567         111 : PetscErrorCode NEPSetProblemType(NEP nep,NEPProblemType type)
     568             : {
     569         111 :   PetscFunctionBegin;
     570         111 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     571         222 :   PetscValidLogicalCollectiveEnum(nep,type,2);
     572         111 :   if (type!=NEP_GENERAL && type!=NEP_RATIONAL) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_WRONG,"Unknown eigenvalue problem type");
     573         111 :   if (type != nep->problem_type) {
     574         111 :     nep->problem_type = type;
     575         111 :     nep->state = NEP_STATE_INITIAL;
     576             :   }
     577         111 :   PetscFunctionReturn(0);
     578             : }
     579             : 
     580             : /*@
     581             :    NEPGetProblemType - Gets the problem type from the NEP object.
     582             : 
     583             :    Not Collective
     584             : 
     585             :    Input Parameter:
     586             : .  nep - the nonlinear eigensolver context
     587             : 
     588             :    Output Parameter:
     589             : .  type - the problem type
     590             : 
     591             :    Level: intermediate
     592             : 
     593             : .seealso: NEPSetProblemType(), NEPProblemType
     594             : @*/
     595           2 : PetscErrorCode NEPGetProblemType(NEP nep,NEPProblemType *type)
     596             : {
     597           2 :   PetscFunctionBegin;
     598           2 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     599           2 :   PetscValidPointer(type,2);
     600           2 :   *type = nep->problem_type;
     601           2 :   PetscFunctionReturn(0);
     602             : }
     603             : 
     604             : /*@
     605             :    NEPSetTwoSided - Sets the solver to use a two-sided variant so that left
     606             :    eigenvectors are also computed.
     607             : 
     608             :    Logically Collective on nep
     609             : 
     610             :    Input Parameters:
     611             : +  nep      - the eigensolver context
     612             : -  twosided - whether the two-sided variant is to be used or not
     613             : 
     614             :    Options Database Keys:
     615             : .  -nep_two_sided <boolean> - Sets/resets the twosided flag
     616             : 
     617             :    Notes:
     618             :    If the user sets twosided=PETSC_TRUE then the solver uses a variant of
     619             :    the algorithm that computes both right and left eigenvectors. This is
     620             :    usually much more costly. This option is not available in all solvers.
     621             : 
     622             :    When using two-sided solvers, the problem matrices must have both the
     623             :    MatMult and MatMultTranspose operations defined.
     624             : 
     625             :    Level: advanced
     626             : 
     627             : .seealso: NEPGetTwoSided(), NEPGetLeftEigenvector()
     628             : @*/
     629           6 : PetscErrorCode NEPSetTwoSided(NEP nep,PetscBool twosided)
     630             : {
     631           6 :   PetscFunctionBegin;
     632           6 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     633          12 :   PetscValidLogicalCollectiveBool(nep,twosided,2);
     634           6 :   if (twosided!=nep->twosided) {
     635           5 :     nep->twosided = twosided;
     636           5 :     nep->state    = NEP_STATE_INITIAL;
     637             :   }
     638           6 :   PetscFunctionReturn(0);
     639             : }
     640             : 
     641             : /*@
     642             :    NEPGetTwoSided - Returns the flag indicating whether a two-sided variant
     643             :    of the algorithm is being used or not.
     644             : 
     645             :    Not Collective
     646             : 
     647             :    Input Parameter:
     648             : .  nep - the eigensolver context
     649             : 
     650             :    Output Parameter:
     651             : .  twosided - the returned flag
     652             : 
     653             :    Level: advanced
     654             : 
     655             : .seealso: NEPSetTwoSided()
     656             : @*/
     657           3 : PetscErrorCode NEPGetTwoSided(NEP nep,PetscBool *twosided)
     658             : {
     659           3 :   PetscFunctionBegin;
     660           3 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     661           3 :   PetscValidBoolPointer(twosided,2);
     662           3 :   *twosided = nep->twosided;
     663           3 :   PetscFunctionReturn(0);
     664             : }
     665             : 
     666             : /*@C
     667             :    NEPSetConvergenceTestFunction - Sets a function to compute the error estimate
     668             :    used in the convergence test.
     669             : 
     670             :    Logically Collective on nep
     671             : 
     672             :    Input Parameters:
     673             : +  nep     - nonlinear eigensolver context obtained from NEPCreate()
     674             : .  func    - a pointer to the convergence test function
     675             : .  ctx     - context for private data for the convergence routine (may be null)
     676             : -  destroy - a routine for destroying the context (may be null)
     677             : 
     678             :    Calling Sequence of func:
     679             : $   func(NEP nep,PetscScalar eigr,PetscScalar eigi,PetscReal res,PetscReal *errest,void *ctx)
     680             : 
     681             : +   nep    - nonlinear eigensolver context obtained from NEPCreate()
     682             : .   eigr   - real part of the eigenvalue
     683             : .   eigi   - imaginary part of the eigenvalue
     684             : .   res    - residual norm associated to the eigenpair
     685             : .   errest - (output) computed error estimate
     686             : -   ctx    - optional context, as set by NEPSetConvergenceTestFunction()
     687             : 
     688             :    Note:
     689             :    If the error estimate returned by the convergence test function is less than
     690             :    the tolerance, then the eigenvalue is accepted as converged.
     691             : 
     692             :    Level: advanced
     693             : 
     694             : .seealso: NEPSetConvergenceTest(), NEPSetTolerances()
     695             : @*/
     696           1 : PetscErrorCode NEPSetConvergenceTestFunction(NEP nep,PetscErrorCode (*func)(NEP,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
     697             : {
     698           1 :   PetscErrorCode ierr;
     699             : 
     700           1 :   PetscFunctionBegin;
     701           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     702           1 :   if (nep->convergeddestroy) {
     703           0 :     ierr = (*nep->convergeddestroy)(nep->convergedctx);CHKERRQ(ierr);
     704             :   }
     705           1 :   nep->convergeduser    = func;
     706           1 :   nep->convergeddestroy = destroy;
     707           1 :   nep->convergedctx     = ctx;
     708           1 :   if (func == NEPConvergedRelative) nep->conv = NEP_CONV_REL;
     709           1 :   else if (func == NEPConvergedNorm) nep->conv = NEP_CONV_NORM;
     710           1 :   else if (func == NEPConvergedAbsolute) nep->conv = NEP_CONV_ABS;
     711             :   else {
     712           1 :     nep->conv      = NEP_CONV_USER;
     713           1 :     nep->converged = nep->convergeduser;
     714             :   }
     715           1 :   PetscFunctionReturn(0);
     716             : }
     717             : 
     718             : /*@
     719             :    NEPSetConvergenceTest - Specifies how to compute the error estimate
     720             :    used in the convergence test.
     721             : 
     722             :    Logically Collective on nep
     723             : 
     724             :    Input Parameters:
     725             : +  nep  - nonlinear eigensolver context obtained from NEPCreate()
     726             : -  conv - the type of convergence test
     727             : 
     728             :    Options Database Keys:
     729             : +  -nep_conv_abs  - Sets the absolute convergence test
     730             : .  -nep_conv_rel  - Sets the convergence test relative to the eigenvalue
     731             : -  -nep_conv_user - Selects the user-defined convergence test
     732             : 
     733             :    Note:
     734             :    The parameter 'conv' can have one of these values
     735             : +     NEP_CONV_ABS  - absolute error ||r||
     736             : .     NEP_CONV_REL  - error relative to the eigenvalue l, ||r||/|l|
     737             : .     NEP_CONV_NORM - error relative matrix norms, ||r||/sum_i(|f_i(l)|*||A_i||)
     738             : -     NEP_CONV_USER - function set by NEPSetConvergenceTestFunction()
     739             : 
     740             :    Level: intermediate
     741             : 
     742             : .seealso: NEPGetConvergenceTest(), NEPSetConvergenceTestFunction(), NEPSetStoppingTest(), NEPConv
     743             : @*/
     744           2 : PetscErrorCode NEPSetConvergenceTest(NEP nep,NEPConv conv)
     745             : {
     746           2 :   PetscFunctionBegin;
     747           2 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     748           4 :   PetscValidLogicalCollectiveEnum(nep,conv,2);
     749           2 :   switch (conv) {
     750           1 :     case NEP_CONV_ABS:  nep->converged = NEPConvergedAbsolute; break;
     751           0 :     case NEP_CONV_REL:  nep->converged = NEPConvergedRelative; break;
     752           1 :     case NEP_CONV_NORM: nep->converged = NEPConvergedNorm; break;
     753           0 :     case NEP_CONV_USER:
     754           0 :       if (!nep->convergeduser) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ORDER,"Must call NEPSetConvergenceTestFunction() first");
     755           0 :       nep->converged = nep->convergeduser;
     756           0 :       break;
     757           0 :     default:
     758           0 :       SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'conv' value");
     759             :   }
     760           2 :   nep->conv = conv;
     761           2 :   PetscFunctionReturn(0);
     762             : }
     763             : 
     764             : /*@
     765             :    NEPGetConvergenceTest - Gets the method used to compute the error estimate
     766             :    used in the convergence test.
     767             : 
     768             :    Not Collective
     769             : 
     770             :    Input Parameters:
     771             : .  nep   - nonlinear eigensolver context obtained from NEPCreate()
     772             : 
     773             :    Output Parameters:
     774             : .  conv  - the type of convergence test
     775             : 
     776             :    Level: intermediate
     777             : 
     778             : .seealso: NEPSetConvergenceTest(), NEPConv
     779             : @*/
     780           1 : PetscErrorCode NEPGetConvergenceTest(NEP nep,NEPConv *conv)
     781             : {
     782           1 :   PetscFunctionBegin;
     783           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     784           1 :   PetscValidPointer(conv,2);
     785           1 :   *conv = nep->conv;
     786           1 :   PetscFunctionReturn(0);
     787             : }
     788             : 
     789             : /*@C
     790             :    NEPSetStoppingTestFunction - Sets a function to decide when to stop the outer
     791             :    iteration of the eigensolver.
     792             : 
     793             :    Logically Collective on nep
     794             : 
     795             :    Input Parameters:
     796             : +  nep     - nonlinear eigensolver context obtained from NEPCreate()
     797             : .  func    - pointer to the stopping test function
     798             : .  ctx     - context for private data for the stopping routine (may be null)
     799             : -  destroy - a routine for destroying the context (may be null)
     800             : 
     801             :    Calling Sequence of func:
     802             : $   func(NEP nep,PetscInt its,PetscInt max_it,PetscInt nconv,PetscInt nev,NEPConvergedReason *reason,void *ctx)
     803             : 
     804             : +   nep    - nonlinear eigensolver context obtained from NEPCreate()
     805             : .   its    - current number of iterations
     806             : .   max_it - maximum number of iterations
     807             : .   nconv  - number of currently converged eigenpairs
     808             : .   nev    - number of requested eigenpairs
     809             : .   reason - (output) result of the stopping test
     810             : -   ctx    - optional context, as set by NEPSetStoppingTestFunction()
     811             : 
     812             :    Note:
     813             :    Normal usage is to first call the default routine NEPStoppingBasic() and then
     814             :    set reason to NEP_CONVERGED_USER if some user-defined conditions have been
     815             :    met. To let the eigensolver continue iterating, the result must be left as
     816             :    NEP_CONVERGED_ITERATING.
     817             : 
     818             :    Level: advanced
     819             : 
     820             : .seealso: NEPSetStoppingTest(), NEPStoppingBasic()
     821             : @*/
     822           1 : PetscErrorCode NEPSetStoppingTestFunction(NEP nep,PetscErrorCode (*func)(NEP,PetscInt,PetscInt,PetscInt,PetscInt,NEPConvergedReason*,void*),void* ctx,PetscErrorCode (*destroy)(void*))
     823             : {
     824           1 :   PetscErrorCode ierr;
     825             : 
     826           1 :   PetscFunctionBegin;
     827           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     828           1 :   if (nep->stoppingdestroy) {
     829           0 :     ierr = (*nep->stoppingdestroy)(nep->stoppingctx);CHKERRQ(ierr);
     830             :   }
     831           1 :   nep->stoppinguser    = func;
     832           1 :   nep->stoppingdestroy = destroy;
     833           1 :   nep->stoppingctx     = ctx;
     834           1 :   if (func == NEPStoppingBasic) nep->stop = NEP_STOP_BASIC;
     835             :   else {
     836           1 :     nep->stop     = NEP_STOP_USER;
     837           1 :     nep->stopping = nep->stoppinguser;
     838             :   }
     839           1 :   PetscFunctionReturn(0);
     840             : }
     841             : 
     842             : /*@
     843             :    NEPSetStoppingTest - Specifies how to decide the termination of the outer
     844             :    loop of the eigensolver.
     845             : 
     846             :    Logically Collective on nep
     847             : 
     848             :    Input Parameters:
     849             : +  nep  - nonlinear eigensolver context obtained from NEPCreate()
     850             : -  stop - the type of stopping test
     851             : 
     852             :    Options Database Keys:
     853             : +  -nep_stop_basic - Sets the default stopping test
     854             : -  -nep_stop_user  - Selects the user-defined stopping test
     855             : 
     856             :    Note:
     857             :    The parameter 'stop' can have one of these values
     858             : +     NEP_STOP_BASIC - default stopping test
     859             : -     NEP_STOP_USER  - function set by NEPSetStoppingTestFunction()
     860             : 
     861             :    Level: advanced
     862             : 
     863             : .seealso: NEPGetStoppingTest(), NEPSetStoppingTestFunction(), NEPSetConvergenceTest(), NEPStop
     864             : @*/
     865           1 : PetscErrorCode NEPSetStoppingTest(NEP nep,NEPStop stop)
     866             : {
     867           1 :   PetscFunctionBegin;
     868           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     869           2 :   PetscValidLogicalCollectiveEnum(nep,stop,2);
     870           1 :   switch (stop) {
     871           1 :     case NEP_STOP_BASIC: nep->stopping = NEPStoppingBasic; break;
     872           0 :     case NEP_STOP_USER:
     873           0 :       if (!nep->stoppinguser) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ORDER,"Must call NEPSetStoppingTestFunction() first");
     874           0 :       nep->stopping = nep->stoppinguser;
     875           0 :       break;
     876           0 :     default:
     877           0 :       SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Invalid 'stop' value");
     878             :   }
     879           1 :   nep->stop = stop;
     880           1 :   PetscFunctionReturn(0);
     881             : }
     882             : 
     883             : /*@
     884             :    NEPGetStoppingTest - Gets the method used to decide the termination of the outer
     885             :    loop of the eigensolver.
     886             : 
     887             :    Not Collective
     888             : 
     889             :    Input Parameters:
     890             : .  nep   - nonlinear eigensolver context obtained from NEPCreate()
     891             : 
     892             :    Output Parameters:
     893             : .  stop  - the type of stopping test
     894             : 
     895             :    Level: advanced
     896             : 
     897             : .seealso: NEPSetStoppingTest(), NEPStop
     898             : @*/
     899           1 : PetscErrorCode NEPGetStoppingTest(NEP nep,NEPStop *stop)
     900             : {
     901           1 :   PetscFunctionBegin;
     902           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     903           1 :   PetscValidPointer(stop,2);
     904           1 :   *stop = nep->stop;
     905           1 :   PetscFunctionReturn(0);
     906             : }
     907             : 
     908             : /*@
     909             :    NEPSetTrackAll - Specifies if the solver must compute the residual of all
     910             :    approximate eigenpairs or not.
     911             : 
     912             :    Logically Collective on nep
     913             : 
     914             :    Input Parameters:
     915             : +  nep      - the eigensolver context
     916             : -  trackall - whether compute all residuals or not
     917             : 
     918             :    Notes:
     919             :    If the user sets trackall=PETSC_TRUE then the solver explicitly computes
     920             :    the residual for each eigenpair approximation. Computing the residual is
     921             :    usually an expensive operation and solvers commonly compute the associated
     922             :    residual to the first unconverged eigenpair.
     923             : 
     924             :    The option '-nep_monitor_all' automatically activates this option.
     925             : 
     926             :    Level: developer
     927             : 
     928             : .seealso: NEPGetTrackAll()
     929             : @*/
     930           2 : PetscErrorCode NEPSetTrackAll(NEP nep,PetscBool trackall)
     931             : {
     932           2 :   PetscFunctionBegin;
     933           2 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     934           4 :   PetscValidLogicalCollectiveBool(nep,trackall,2);
     935           2 :   nep->trackall = trackall;
     936           2 :   PetscFunctionReturn(0);
     937             : }
     938             : 
     939             : /*@
     940             :    NEPGetTrackAll - Returns the flag indicating whether all residual norms must
     941             :    be computed or not.
     942             : 
     943             :    Not Collective
     944             : 
     945             :    Input Parameter:
     946             : .  nep - the eigensolver context
     947             : 
     948             :    Output Parameter:
     949             : .  trackall - the returned flag
     950             : 
     951             :    Level: developer
     952             : 
     953             : .seealso: NEPSetTrackAll()
     954             : @*/
     955          27 : PetscErrorCode NEPGetTrackAll(NEP nep,PetscBool *trackall)
     956             : {
     957          27 :   PetscFunctionBegin;
     958          27 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     959          27 :   PetscValidBoolPointer(trackall,2);
     960          27 :   *trackall = nep->trackall;
     961          27 :   PetscFunctionReturn(0);
     962             : }
     963             : 
     964             : /*@
     965             :    NEPSetRefine - Specifies the refinement type (and options) to be used
     966             :    after the solve.
     967             : 
     968             :    Logically Collective on nep
     969             : 
     970             :    Input Parameters:
     971             : +  nep    - the nonlinear eigensolver context
     972             : .  refine - refinement type
     973             : .  npart  - number of partitions of the communicator
     974             : .  tol    - the convergence tolerance
     975             : .  its    - maximum number of refinement iterations
     976             : -  scheme - which scheme to be used for solving the involved linear systems
     977             : 
     978             :    Options Database Keys:
     979             : +  -nep_refine <type> - refinement type, one of <none,simple,multiple>
     980             : .  -nep_refine_partitions <n> - the number of partitions
     981             : .  -nep_refine_tol <tol> - the tolerance
     982             : .  -nep_refine_its <its> - number of iterations
     983             : -  -nep_refine_scheme - to set the scheme for the linear solves
     984             : 
     985             :    Notes:
     986             :    By default, iterative refinement is disabled, since it may be very
     987             :    costly. There are two possible refinement strategies: simple and multiple.
     988             :    The simple approach performs iterative refinement on each of the
     989             :    converged eigenpairs individually, whereas the multiple strategy works
     990             :    with the invariant pair as a whole, refining all eigenpairs simultaneously.
     991             :    The latter may be required for the case of multiple eigenvalues.
     992             : 
     993             :    In some cases, especially when using direct solvers within the
     994             :    iterative refinement method, it may be helpful for improved scalability
     995             :    to split the communicator in several partitions. The npart parameter
     996             :    indicates how many partitions to use (defaults to 1).
     997             : 
     998             :    The tol and its parameters specify the stopping criterion. In the simple
     999             :    method, refinement continues until the residual of each eigenpair is
    1000             :    below the tolerance (tol defaults to the NEP tol, but may be set to a
    1001             :    different value). In contrast, the multiple method simply performs its
    1002             :    refinement iterations (just one by default).
    1003             : 
    1004             :    The scheme argument is used to change the way in which linear systems are
    1005             :    solved. Possible choices are: explicit, mixed block elimination (MBE),
    1006             :    and Schur complement.
    1007             : 
    1008             :    Level: intermediate
    1009             : 
    1010             : .seealso: NEPGetRefine()
    1011             : @*/
    1012          13 : PetscErrorCode NEPSetRefine(NEP nep,NEPRefine refine,PetscInt npart,PetscReal tol,PetscInt its,NEPRefineScheme scheme)
    1013             : {
    1014          13 :   PetscErrorCode ierr;
    1015          13 :   PetscMPIInt    size;
    1016             : 
    1017          13 :   PetscFunctionBegin;
    1018          13 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
    1019          26 :   PetscValidLogicalCollectiveEnum(nep,refine,2);
    1020          26 :   PetscValidLogicalCollectiveInt(nep,npart,3);
    1021          26 :   PetscValidLogicalCollectiveReal(nep,tol,4);
    1022          26 :   PetscValidLogicalCollectiveInt(nep,its,5);
    1023          26 :   PetscValidLogicalCollectiveEnum(nep,scheme,6);
    1024          13 :   nep->refine = refine;
    1025          13 :   if (refine) {  /* process parameters only if not REFINE_NONE */
    1026          13 :     if (npart!=nep->npart) {
    1027           8 :       ierr = PetscSubcommDestroy(&nep->refinesubc);CHKERRQ(ierr);
    1028           8 :       ierr = KSPDestroy(&nep->refineksp);CHKERRQ(ierr);
    1029             :     }
    1030          13 :     if (npart == PETSC_DEFAULT || npart == PETSC_DECIDE) {
    1031           0 :       nep->npart = 1;
    1032             :     } else {
    1033          13 :       ierr = MPI_Comm_size(PetscObjectComm((PetscObject)nep),&size);CHKERRMPI(ierr);
    1034          13 :       if (npart<1 || npart>size) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of npart");
    1035          13 :       nep->npart = npart;
    1036             :     }
    1037          13 :     if (tol == PETSC_DEFAULT || tol == PETSC_DECIDE) {
    1038          12 :       nep->rtol = PETSC_DEFAULT;
    1039             :     } else {
    1040           1 :       if (tol<=0.0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of tol. Must be > 0");
    1041           1 :       nep->rtol = tol;
    1042             :     }
    1043          13 :     if (its==PETSC_DECIDE || its==PETSC_DEFAULT) {
    1044          11 :       nep->rits = PETSC_DEFAULT;
    1045             :     } else {
    1046           2 :       if (its<0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of its. Must be >= 0");
    1047           2 :       nep->rits = its;
    1048             :     }
    1049          13 :     nep->scheme = scheme;
    1050             :   }
    1051          13 :   nep->state = NEP_STATE_INITIAL;
    1052          13 :   PetscFunctionReturn(0);
    1053             : }
    1054             : 
    1055             : /*@C
    1056             :    NEPGetRefine - Gets the refinement strategy used by the NEP object, and the
    1057             :    associated parameters.
    1058             : 
    1059             :    Not Collective
    1060             : 
    1061             :    Input Parameter:
    1062             : .  nep - the nonlinear eigensolver context
    1063             : 
    1064             :    Output Parameters:
    1065             : +  refine - refinement type
    1066             : .  npart  - number of partitions of the communicator
    1067             : .  tol    - the convergence tolerance
    1068             : .  its    - maximum number of refinement iterations
    1069             : -  scheme - the scheme used for solving linear systems
    1070             : 
    1071             :    Level: intermediate
    1072             : 
    1073             :    Note:
    1074             :    The user can specify NULL for any parameter that is not needed.
    1075             : 
    1076             : .seealso: NEPSetRefine()
    1077             : @*/
    1078           1 : PetscErrorCode NEPGetRefine(NEP nep,NEPRefine *refine,PetscInt *npart,PetscReal *tol,PetscInt *its,NEPRefineScheme *scheme)
    1079             : {
    1080           1 :   PetscFunctionBegin;
    1081           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
    1082           1 :   if (refine) *refine = nep->refine;
    1083           1 :   if (npart)  *npart  = nep->npart;
    1084           1 :   if (tol)    *tol    = nep->rtol;
    1085           1 :   if (its)    *its    = nep->rits;
    1086           1 :   if (scheme) *scheme = nep->scheme;
    1087           1 :   PetscFunctionReturn(0);
    1088             : }
    1089             : 
    1090             : /*@C
    1091             :    NEPSetOptionsPrefix - Sets the prefix used for searching for all
    1092             :    NEP options in the database.
    1093             : 
    1094             :    Logically Collective on nep
    1095             : 
    1096             :    Input Parameters:
    1097             : +  nep - the nonlinear eigensolver context
    1098             : -  prefix - the prefix string to prepend to all NEP option requests
    1099             : 
    1100             :    Notes:
    1101             :    A hyphen (-) must NOT be given at the beginning of the prefix name.
    1102             :    The first character of all runtime options is AUTOMATICALLY the
    1103             :    hyphen.
    1104             : 
    1105             :    For example, to distinguish between the runtime options for two
    1106             :    different NEP contexts, one could call
    1107             : .vb
    1108             :       NEPSetOptionsPrefix(nep1,"neig1_")
    1109             :       NEPSetOptionsPrefix(nep2,"neig2_")
    1110             : .ve
    1111             : 
    1112             :    Level: advanced
    1113             : 
    1114             : .seealso: NEPAppendOptionsPrefix(), NEPGetOptionsPrefix()
    1115             : @*/
    1116           1 : PetscErrorCode NEPSetOptionsPrefix(NEP nep,const char *prefix)
    1117             : {
    1118           1 :   PetscErrorCode ierr;
    1119             : 
    1120           1 :   PetscFunctionBegin;
    1121           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
    1122           1 :   if (!nep->V) { ierr = NEPGetBV(nep,&nep->V);CHKERRQ(ierr); }
    1123           1 :   ierr = BVSetOptionsPrefix(nep->V,prefix);CHKERRQ(ierr);
    1124           1 :   if (!nep->ds) { ierr = NEPGetDS(nep,&nep->ds);CHKERRQ(ierr); }
    1125           1 :   ierr = DSSetOptionsPrefix(nep->ds,prefix);CHKERRQ(ierr);
    1126           1 :   if (!nep->rg) { ierr = NEPGetRG(nep,&nep->rg);CHKERRQ(ierr); }
    1127           1 :   ierr = RGSetOptionsPrefix(nep->rg,prefix);CHKERRQ(ierr);
    1128           1 :   ierr = PetscObjectSetOptionsPrefix((PetscObject)nep,prefix);CHKERRQ(ierr);
    1129           1 :   PetscFunctionReturn(0);
    1130             : }
    1131             : 
    1132             : /*@C
    1133             :    NEPAppendOptionsPrefix - Appends to the prefix used for searching for all
    1134             :    NEP options in the database.
    1135             : 
    1136             :    Logically Collective on nep
    1137             : 
    1138             :    Input Parameters:
    1139             : +  nep - the nonlinear eigensolver context
    1140             : -  prefix - the prefix string to prepend to all NEP option requests
    1141             : 
    1142             :    Notes:
    1143             :    A hyphen (-) must NOT be given at the beginning of the prefix name.
    1144             :    The first character of all runtime options is AUTOMATICALLY the hyphen.
    1145             : 
    1146             :    Level: advanced
    1147             : 
    1148             : .seealso: NEPSetOptionsPrefix(), NEPGetOptionsPrefix()
    1149             : @*/
    1150           1 : PetscErrorCode NEPAppendOptionsPrefix(NEP nep,const char *prefix)
    1151             : {
    1152           1 :   PetscErrorCode ierr;
    1153             : 
    1154           1 :   PetscFunctionBegin;
    1155           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
    1156           1 :   if (!nep->V) { ierr = NEPGetBV(nep,&nep->V);CHKERRQ(ierr); }
    1157           1 :   ierr = BVAppendOptionsPrefix(nep->V,prefix);CHKERRQ(ierr);
    1158           1 :   if (!nep->ds) { ierr = NEPGetDS(nep,&nep->ds);CHKERRQ(ierr); }
    1159           1 :   ierr = DSAppendOptionsPrefix(nep->ds,prefix);CHKERRQ(ierr);
    1160           1 :   if (!nep->rg) { ierr = NEPGetRG(nep,&nep->rg);CHKERRQ(ierr); }
    1161           1 :   ierr = RGAppendOptionsPrefix(nep->rg,prefix);CHKERRQ(ierr);
    1162           1 :   ierr = PetscObjectAppendOptionsPrefix((PetscObject)nep,prefix);CHKERRQ(ierr);
    1163           1 :   PetscFunctionReturn(0);
    1164             : }
    1165             : 
    1166             : /*@C
    1167             :    NEPGetOptionsPrefix - Gets the prefix used for searching for all
    1168             :    NEP options in the database.
    1169             : 
    1170             :    Not Collective
    1171             : 
    1172             :    Input Parameters:
    1173             : .  nep - the nonlinear eigensolver context
    1174             : 
    1175             :    Output Parameters:
    1176             : .  prefix - pointer to the prefix string used is returned
    1177             : 
    1178             :    Note:
    1179             :    On the Fortran side, the user should pass in a string 'prefix' of
    1180             :    sufficient length to hold the prefix.
    1181             : 
    1182             :    Level: advanced
    1183             : 
    1184             : .seealso: NEPSetOptionsPrefix(), NEPAppendOptionsPrefix()
    1185             : @*/
    1186           1 : PetscErrorCode NEPGetOptionsPrefix(NEP nep,const char *prefix[])
    1187             : {
    1188           1 :   PetscErrorCode ierr;
    1189             : 
    1190           1 :   PetscFunctionBegin;
    1191           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
    1192           1 :   PetscValidPointer(prefix,2);
    1193           1 :   ierr = PetscObjectGetOptionsPrefix((PetscObject)nep,prefix);CHKERRQ(ierr);
    1194           1 :   PetscFunctionReturn(0);
    1195             : }
    1196             : 

Generated by: LCOV version 1.14