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

Generated by: LCOV version 1.14