LCOV - code coverage report
Current view: top level - nep/impls/slp - slp.c (source / functions) Hit Total Coverage
Test: SLEPc Lines: 367 394 93.1 %
Date: 2021-08-02 00:32:28 Functions: 24 26 92.3 %
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             :    SLEPc nonlinear eigensolver: "slp"
      12             : 
      13             :    Method: Succesive linear problems
      14             : 
      15             :    Algorithm:
      16             : 
      17             :        Newton-type iteration based on first order Taylor approximation.
      18             : 
      19             :    References:
      20             : 
      21             :        [1] A. Ruhe, "Algorithms for the nonlinear eigenvalue problem", SIAM J.
      22             :            Numer. Anal. 10(4):674-689, 1973.
      23             : */
      24             : 
      25             : #include <slepc/private/nepimpl.h>         /*I "slepcnep.h" I*/
      26             : #include <../src/nep/impls/nepdefl.h>
      27             : #include "slp.h"
      28             : 
      29             : typedef struct {
      30             :   NEP_EXT_OP extop;
      31             :   Vec        w;
      32             : } NEP_SLP_MATSHELL;
      33             : 
      34          20 : PetscErrorCode NEPSetUp_SLP(NEP nep)
      35             : {
      36          20 :   PetscErrorCode ierr;
      37          20 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
      38          20 :   PetscBool      flg;
      39          20 :   ST             st;
      40             : 
      41          20 :   PetscFunctionBegin;
      42          20 :   if (nep->ncv!=PETSC_DEFAULT) { ierr = PetscInfo(nep,"Setting ncv = nev, ignoring user-provided value\n");CHKERRQ(ierr); }
      43          20 :   nep->ncv = nep->nev;
      44          20 :   if (nep->mpd!=PETSC_DEFAULT) { ierr = PetscInfo(nep,"Setting mpd = nev, ignoring user-provided value\n");CHKERRQ(ierr); }
      45          20 :   nep->mpd = nep->nev;
      46          20 :   if (nep->ncv>nep->nev+nep->mpd) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_USER_INPUT,"The value of ncv must not be larger than nev+mpd");
      47          20 :   if (nep->max_it==PETSC_DEFAULT) nep->max_it = PetscMax(5000,2*nep->n/nep->ncv);
      48          20 :   if (!nep->which) nep->which = NEP_TARGET_MAGNITUDE;
      49          20 :   if (nep->which!=NEP_TARGET_MAGNITUDE) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"This solver supports only target magnitude eigenvalues");
      50          20 :   NEPCheckUnsupported(nep,NEP_FEATURE_REGION);
      51             : 
      52          20 :   if (!ctx->eps) { ierr = NEPSLPGetEPS(nep,&ctx->eps);CHKERRQ(ierr); }
      53          20 :   ierr = EPSGetST(ctx->eps,&st);CHKERRQ(ierr);
      54          20 :   ierr = PetscObjectTypeCompareAny((PetscObject)st,&flg,STSINVERT,STCAYLEY,"");CHKERRQ(ierr);
      55          20 :   if (flg) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"SLP does not support spectral transformation");
      56          20 :   ierr = EPSSetDimensions(ctx->eps,1,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr);
      57          20 :   ierr = EPSSetWhichEigenpairs(ctx->eps,EPS_LARGEST_MAGNITUDE);CHKERRQ(ierr);
      58          31 :   ierr = EPSSetTolerances(ctx->eps,SlepcDefaultTol(nep->tol)/10.0,nep->max_it);CHKERRQ(ierr);
      59          20 :   if (nep->tol==PETSC_DEFAULT) nep->tol = SLEPC_DEFAULT_TOL;
      60          20 :   if (ctx->deftol==PETSC_DEFAULT) ctx->deftol = nep->tol;
      61             : 
      62          20 :   if (nep->twosided) {
      63           3 :     nep->ops->solve = NEPSolve_SLP_Twosided;
      64           3 :     nep->ops->computevectors = NULL;
      65           3 :     if (!ctx->epsts) { ierr = NEPSLPGetEPSLeft(nep,&ctx->epsts);CHKERRQ(ierr); }
      66           3 :     ierr = EPSGetST(ctx->epsts,&st);CHKERRQ(ierr);
      67           3 :     ierr = PetscObjectTypeCompareAny((PetscObject)st,&flg,STSINVERT,STCAYLEY,"");CHKERRQ(ierr);
      68           3 :     if (flg) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"SLP does not support spectral transformation");
      69           3 :     ierr = EPSSetDimensions(ctx->epsts,1,PETSC_DECIDE,PETSC_DECIDE);CHKERRQ(ierr);
      70           3 :     ierr = EPSSetWhichEigenpairs(ctx->epsts,EPS_LARGEST_MAGNITUDE);CHKERRQ(ierr);
      71           6 :     ierr = EPSSetTolerances(ctx->epsts,SlepcDefaultTol(nep->tol)/10.0,nep->max_it);CHKERRQ(ierr);
      72             :   } else {
      73          17 :     nep->ops->solve = NEPSolve_SLP;
      74          17 :     nep->ops->computevectors = NEPComputeVectors_Schur;
      75             :   }
      76          20 :   ierr = NEPAllocateSolution(nep,0);CHKERRQ(ierr);
      77          20 :   PetscFunctionReturn(0);
      78             : }
      79             : 
      80         894 : static PetscErrorCode MatMult_SLP(Mat M,Vec x,Vec y)
      81             : {
      82         894 :   PetscErrorCode   ierr;
      83         894 :   NEP_SLP_MATSHELL *ctx;
      84             : 
      85         894 :   PetscFunctionBegin;
      86         894 :   ierr = MatShellGetContext(M,(void**)&ctx);CHKERRQ(ierr);
      87         894 :   ierr = MatMult(ctx->extop->MJ,x,ctx->w);CHKERRQ(ierr);
      88         894 :   ierr = NEPDeflationFunctionSolve(ctx->extop,ctx->w,y);CHKERRQ(ierr);
      89         894 :   PetscFunctionReturn(0);
      90             : }
      91             : 
      92          17 : static PetscErrorCode MatDestroy_SLP(Mat M)
      93             : {
      94          17 :   PetscErrorCode   ierr;
      95          17 :   NEP_SLP_MATSHELL *ctx;
      96             : 
      97          17 :   PetscFunctionBegin;
      98          17 :   ierr = MatShellGetContext(M,(void**)&ctx);CHKERRQ(ierr);
      99          17 :   ierr = VecDestroy(&ctx->w);CHKERRQ(ierr);
     100          17 :   ierr = PetscFree(ctx);CHKERRQ(ierr);
     101          17 :   PetscFunctionReturn(0);
     102             : }
     103             : 
     104             : #if defined(PETSC_HAVE_CUDA)
     105             : static PetscErrorCode MatCreateVecs_SLP(Mat M,Vec *left,Vec *right)
     106             : {
     107             :   PetscErrorCode   ierr;
     108             :   NEP_SLP_MATSHELL *ctx;
     109             : 
     110             :   PetscFunctionBegin;
     111             :   ierr = MatShellGetContext(M,(void**)&ctx);CHKERRQ(ierr);
     112             :   if (right) {
     113             :     ierr = VecDuplicate(ctx->w,right);CHKERRQ(ierr);
     114             :   }
     115             :   if (left) {
     116             :     ierr = VecDuplicate(ctx->w,left);CHKERRQ(ierr);
     117             :   }
     118             :   PetscFunctionReturn(0);
     119             : }
     120             : #endif
     121             : 
     122          71 : static PetscErrorCode NEPSLPSetUpLinearEP(NEP nep,NEP_EXT_OP extop,PetscScalar lambda,Vec u,PetscBool ini)
     123             : {
     124          71 :   PetscErrorCode   ierr;
     125          71 :   NEP_SLP          *slpctx = (NEP_SLP*)nep->data;
     126          71 :   Mat              Mshell;
     127          71 :   PetscInt         nloc,mloc;
     128          71 :   NEP_SLP_MATSHELL *shellctx;
     129             : 
     130          71 :   PetscFunctionBegin;
     131          71 :   if (ini) {
     132             :     /* Create mat shell */
     133          17 :     ierr = PetscNew(&shellctx);CHKERRQ(ierr);
     134          17 :     shellctx->extop = extop;
     135          17 :     ierr = NEPDeflationCreateVec(extop,&shellctx->w);CHKERRQ(ierr);
     136          17 :     ierr = MatGetLocalSize(nep->function,&mloc,&nloc);CHKERRQ(ierr);
     137          17 :     nloc += extop->szd; mloc += extop->szd;
     138          17 :     ierr = MatCreateShell(PetscObjectComm((PetscObject)nep),nloc,mloc,PETSC_DETERMINE,PETSC_DETERMINE,shellctx,&Mshell);CHKERRQ(ierr);
     139          17 :     ierr = MatShellSetOperation(Mshell,MATOP_MULT,(void(*)(void))MatMult_SLP);CHKERRQ(ierr);
     140          17 :     ierr = MatShellSetOperation(Mshell,MATOP_DESTROY,(void(*)(void))MatDestroy_SLP);CHKERRQ(ierr);
     141             : #if defined(PETSC_HAVE_CUDA)
     142             :     ierr = MatShellSetOperation(Mshell,MATOP_CREATE_VECS,(void(*)(void))MatCreateVecs_SLP);CHKERRQ(ierr);
     143             : #endif
     144          17 :     ierr = EPSSetOperators(slpctx->eps,Mshell,NULL);CHKERRQ(ierr);
     145          17 :     ierr = MatDestroy(&Mshell);CHKERRQ(ierr);
     146             :   }
     147          71 :   ierr = NEPDeflationSolveSetUp(extop,lambda);CHKERRQ(ierr);
     148          71 :   ierr = NEPDeflationComputeJacobian(extop,lambda,NULL);CHKERRQ(ierr);
     149          71 :   ierr = EPSSetInitialSpace(slpctx->eps,1,&u);CHKERRQ(ierr);
     150          71 :   PetscFunctionReturn(0);
     151             : }
     152             : 
     153          17 : PetscErrorCode NEPSolve_SLP(NEP nep)
     154             : {
     155          17 :   PetscErrorCode    ierr;
     156          17 :   NEP_SLP           *ctx = (NEP_SLP*)nep->data;
     157          17 :   Mat               F,H;
     158          17 :   Vec               uu,u,r;
     159          17 :   PetscScalar       sigma,lambda,mu,im,*Ap;
     160          17 :   const PetscScalar *Hp;
     161          17 :   PetscReal         resnorm;
     162          17 :   PetscInt          nconv,ldh,ldds,i,j;
     163          17 :   PetscBool         skip=PETSC_FALSE,lock=PETSC_FALSE;
     164          17 :   NEP_EXT_OP        extop=NULL;    /* Extended operator for deflation */
     165             : 
     166          17 :   PetscFunctionBegin;
     167             :   /* get initial approximation of eigenvalue and eigenvector */
     168          17 :   ierr = NEPGetDefaultShift(nep,&sigma);CHKERRQ(ierr);
     169          17 :   if (!nep->nini) {
     170          17 :     ierr = BVSetRandomColumn(nep->V,0);CHKERRQ(ierr);
     171             :   }
     172          17 :   lambda = sigma;
     173          17 :   if (!ctx->ksp) { ierr = NEPSLPGetKSP(nep,&ctx->ksp);CHKERRQ(ierr); }
     174          17 :   ierr = NEPDeflationInitialize(nep,nep->V,ctx->ksp,PETSC_TRUE,nep->nev,&extop);CHKERRQ(ierr);
     175          17 :   ierr = NEPDeflationCreateVec(extop,&u);CHKERRQ(ierr);
     176          17 :   ierr = VecDuplicate(u,&r);CHKERRQ(ierr);
     177          17 :   ierr = BVGetColumn(nep->V,0,&uu);CHKERRQ(ierr);
     178          17 :   ierr = NEPDeflationCopyToExtendedVec(extop,uu,NULL,u,PETSC_FALSE);CHKERRQ(ierr);
     179          17 :   ierr = BVRestoreColumn(nep->V,0,&uu);CHKERRQ(ierr);
     180             : 
     181             :   /* Restart loop */
     182         122 :   while (nep->reason == NEP_CONVERGED_ITERATING) {
     183         105 :     nep->its++;
     184             : 
     185             :     /* form residual,  r = T(lambda)*u (used in convergence test only) */
     186         105 :     ierr = NEPDeflationComputeFunction(extop,lambda,&F);CHKERRQ(ierr);
     187         105 :     ierr = MatMult(F,u,r);CHKERRQ(ierr);
     188             : 
     189             :     /* convergence test */
     190         105 :     ierr = VecNorm(r,NORM_2,&resnorm);CHKERRQ(ierr);
     191         105 :     ierr = (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);CHKERRQ(ierr);
     192         105 :     nep->eigr[nep->nconv] = lambda;
     193         105 :     if (nep->errest[nep->nconv]<=nep->tol || nep->errest[nep->nconv]<=ctx->deftol) {
     194          34 :       if (nep->errest[nep->nconv]<=ctx->deftol && !extop->ref && nep->nconv) {
     195          17 :         ierr = NEPDeflationExtractEigenpair(extop,nep->nconv,u,lambda,nep->ds);CHKERRQ(ierr);
     196          17 :         ierr = NEPDeflationSetRefine(extop,PETSC_TRUE);CHKERRQ(ierr);
     197          17 :         ierr = MatMult(F,u,r);CHKERRQ(ierr);
     198          17 :         ierr = VecNorm(r,NORM_2,&resnorm);CHKERRQ(ierr);
     199          17 :         ierr = (*nep->converged)(nep,lambda,0,resnorm,&nep->errest[nep->nconv],nep->convergedctx);CHKERRQ(ierr);
     200          17 :         if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
     201          17 :       } else if (nep->errest[nep->nconv]<=nep->tol) lock = PETSC_TRUE;
     202             :     }
     203             : 
     204          71 :     if (lock) {
     205          34 :       ierr = NEPDeflationSetRefine(extop,PETSC_FALSE);CHKERRQ(ierr);
     206          34 :       nep->nconv = nep->nconv + 1;
     207          34 :       skip = PETSC_TRUE;
     208          34 :       lock = PETSC_FALSE;
     209          34 :       ierr = NEPDeflationLocking(extop,u,lambda);CHKERRQ(ierr);
     210             :     }
     211         105 :     ierr = (*nep->stopping)(nep,nep->its,nep->max_it,nep->nconv,nep->nev,&nep->reason,nep->stoppingctx);CHKERRQ(ierr);
     212         105 :     if (!skip || nep->reason>0) {
     213          88 :       ierr = NEPMonitor(nep,nep->its,nep->nconv,nep->eigr,nep->eigi,nep->errest,(nep->reason>0)?nep->nconv:nep->nconv+1);CHKERRQ(ierr);
     214             :     }
     215             : 
     216         105 :     if (nep->reason == NEP_CONVERGED_ITERATING) {
     217          88 :       if (!skip) {
     218             :         /* evaluate T(lambda) and T'(lambda) */
     219          71 :         ierr = NEPSLPSetUpLinearEP(nep,extop,lambda,u,nep->its==1?PETSC_TRUE:PETSC_FALSE);CHKERRQ(ierr);
     220             :         /* compute new eigenvalue correction mu and eigenvector approximation u */
     221          71 :         ierr = EPSSolve(ctx->eps);CHKERRQ(ierr);
     222          71 :         ierr = EPSGetConverged(ctx->eps,&nconv);CHKERRQ(ierr);
     223          71 :         if (!nconv) {
     224           0 :           ierr = PetscInfo1(nep,"iter=%D, inner iteration failed, stopping solve\n",nep->its);CHKERRQ(ierr);
     225           0 :           nep->reason = NEP_DIVERGED_LINEAR_SOLVE;
     226           0 :           break;
     227             :         }
     228          71 :         ierr = EPSGetEigenpair(ctx->eps,0,&mu,&im,u,NULL);CHKERRQ(ierr);
     229          71 :         mu = 1.0/mu;
     230          71 :         if (PetscAbsScalar(im)>PETSC_MACHINE_EPSILON) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_SUP,"Complex eigenvalue approximation - not implemented in real scalars");
     231             :       } else {
     232          17 :         nep->its--;  /* do not count this as a full iteration */
     233             :         /* use second eigenpair computed in previous iteration */
     234          17 :         ierr = EPSGetConverged(ctx->eps,&nconv);CHKERRQ(ierr);
     235          17 :         if (nconv>=2) {
     236          17 :           ierr = EPSGetEigenpair(ctx->eps,1,&mu,&im,u,NULL);CHKERRQ(ierr);
     237          17 :           mu = 1.0/mu;
     238             :         } else {
     239           0 :           ierr = NEPDeflationSetRandomVec(extop,u);CHKERRQ(ierr);
     240           0 :           mu = lambda-sigma;
     241             :         }
     242             :         skip = PETSC_FALSE;
     243             :       }
     244             :       /* correct eigenvalue */
     245          88 :       lambda = lambda - mu;
     246             :     }
     247             :   }
     248          17 :   ierr = NEPDeflationGetInvariantPair(extop,NULL,&H);CHKERRQ(ierr);
     249          17 :   ierr = MatGetSize(H,NULL,&ldh);CHKERRQ(ierr);
     250          17 :   ierr = DSReset(nep->ds);CHKERRQ(ierr);
     251          17 :   ierr = DSSetType(nep->ds,DSNHEP);CHKERRQ(ierr);
     252          17 :   ierr = DSAllocate(nep->ds,PetscMax(nep->nconv,1));CHKERRQ(ierr);
     253          17 :   ierr = DSGetLeadingDimension(nep->ds,&ldds);CHKERRQ(ierr);
     254          17 :   ierr = MatDenseGetArrayRead(H,&Hp);CHKERRQ(ierr);
     255          17 :   ierr = DSGetArray(nep->ds,DS_MAT_A,&Ap);CHKERRQ(ierr);
     256          51 :   for (j=0;j<nep->nconv;j++)
     257         124 :     for (i=0;i<nep->nconv;i++) Ap[j*ldds+i] = Hp[j*ldh+i];
     258          17 :   ierr = DSRestoreArray(nep->ds,DS_MAT_A,&Ap);CHKERRQ(ierr);
     259          17 :   ierr = MatDenseRestoreArrayRead(H,&Hp);CHKERRQ(ierr);
     260          17 :   ierr = MatDestroy(&H);CHKERRQ(ierr);
     261          17 :   ierr = DSSetDimensions(nep->ds,nep->nconv,0,nep->nconv);CHKERRQ(ierr);
     262          17 :   ierr = DSSolve(nep->ds,nep->eigr,nep->eigi);CHKERRQ(ierr);
     263          17 :   ierr = NEPDeflationReset(extop);CHKERRQ(ierr);
     264          17 :   ierr = VecDestroy(&u);CHKERRQ(ierr);
     265          17 :   ierr = VecDestroy(&r);CHKERRQ(ierr);
     266          17 :   PetscFunctionReturn(0);
     267             : }
     268             : 
     269          18 : PetscErrorCode NEPSetFromOptions_SLP(PetscOptionItems *PetscOptionsObject,NEP nep)
     270             : {
     271          18 :   PetscErrorCode ierr;
     272          18 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     273          18 :   PetscBool      flg;
     274          18 :   PetscReal      r;
     275             : 
     276          18 :   PetscFunctionBegin;
     277          18 :   ierr = PetscOptionsHead(PetscOptionsObject,"NEP SLP Options");CHKERRQ(ierr);
     278             : 
     279          18 :     r = 0.0;
     280          18 :     ierr = PetscOptionsReal("-nep_slp_deflation_threshold","Tolerance used as a threshold for including deflated eigenpairs","NEPSLPSetDeflationThreshold",ctx->deftol,&r,&flg);CHKERRQ(ierr);
     281          18 :     if (flg) { ierr = NEPSLPSetDeflationThreshold(nep,r);CHKERRQ(ierr); }
     282             : 
     283          18 :   ierr = PetscOptionsTail();CHKERRQ(ierr);
     284             : 
     285          18 :   if (!ctx->eps) { ierr = NEPSLPGetEPS(nep,&ctx->eps);CHKERRQ(ierr); }
     286          18 :   ierr = EPSSetFromOptions(ctx->eps);CHKERRQ(ierr);
     287          18 :   if (nep->twosided) {
     288           3 :     if (!ctx->epsts) { ierr = NEPSLPGetEPSLeft(nep,&ctx->epsts);CHKERRQ(ierr); }
     289           3 :     ierr = EPSSetFromOptions(ctx->epsts);CHKERRQ(ierr);
     290             :   }
     291          18 :   if (!ctx->ksp) { ierr = NEPSLPGetKSP(nep,&ctx->ksp);CHKERRQ(ierr); }
     292          18 :   ierr = KSPSetFromOptions(ctx->ksp);CHKERRQ(ierr);
     293          18 :   PetscFunctionReturn(0);
     294             : }
     295             : 
     296           3 : static PetscErrorCode NEPSLPSetDeflationThreshold_SLP(NEP nep,PetscReal deftol)
     297             : {
     298           3 :   NEP_SLP *ctx = (NEP_SLP*)nep->data;
     299             : 
     300           3 :   PetscFunctionBegin;
     301           3 :   if (deftol == PETSC_DEFAULT) {
     302           0 :     ctx->deftol = PETSC_DEFAULT;
     303           0 :     nep->state  = NEP_STATE_INITIAL;
     304             :   } else {
     305           3 :     if (deftol <= 0.0) SETERRQ(PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Illegal value of deftol. Must be > 0");
     306           3 :     ctx->deftol = deftol;
     307             :   }
     308           3 :   PetscFunctionReturn(0);
     309             : }
     310             : 
     311             : /*@
     312             :    NEPSLPSetDeflationThreshold - Sets the threshold value used to switch between
     313             :    deflated and non-deflated iteration.
     314             : 
     315             :    Logically Collective on nep
     316             : 
     317             :    Input Parameters:
     318             : +  nep    - nonlinear eigenvalue solver
     319             : -  deftol - the threshold value
     320             : 
     321             :    Options Database Keys:
     322             : .  -nep_slp_deflation_threshold <deftol> - set the threshold
     323             : 
     324             :    Notes:
     325             :    Normally, the solver iterates on the extended problem in order to deflate
     326             :    previously converged eigenpairs. If this threshold is set to a nonzero value,
     327             :    then once the residual error is below this threshold the solver will
     328             :    continue the iteration without deflation. The intention is to be able to
     329             :    improve the current eigenpair further, despite having previous eigenpairs
     330             :    with somewhat bad precision.
     331             : 
     332             :    Level: advanced
     333             : 
     334             : .seealso: NEPSLPGetDeflationThreshold()
     335             : @*/
     336           3 : PetscErrorCode NEPSLPSetDeflationThreshold(NEP nep,PetscReal deftol)
     337             : {
     338           3 :   PetscErrorCode ierr;
     339             : 
     340           3 :   PetscFunctionBegin;
     341           3 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     342           6 :   PetscValidLogicalCollectiveReal(nep,deftol,2);
     343           3 :   ierr = PetscTryMethod(nep,"NEPSLPSetDeflationThreshold_C",(NEP,PetscReal),(nep,deftol));CHKERRQ(ierr);
     344           3 :   PetscFunctionReturn(0);
     345             : }
     346             : 
     347           1 : static PetscErrorCode NEPSLPGetDeflationThreshold_SLP(NEP nep,PetscReal *deftol)
     348             : {
     349           1 :   NEP_SLP *ctx = (NEP_SLP*)nep->data;
     350             : 
     351           1 :   PetscFunctionBegin;
     352           1 :   *deftol = ctx->deftol;
     353           1 :   PetscFunctionReturn(0);
     354             : }
     355             : 
     356             : /*@
     357             :    NEPSLPGetDeflationThreshold - Returns the threshold value that controls deflation.
     358             : 
     359             :    Not Collective
     360             : 
     361             :    Input Parameter:
     362             : .  nep - nonlinear eigenvalue solver
     363             : 
     364             :    Output Parameter:
     365             : .  deftol - the threshold
     366             : 
     367             :    Level: advanced
     368             : 
     369             : .seealso: NEPSLPSetDeflationThreshold()
     370             : @*/
     371           1 : PetscErrorCode NEPSLPGetDeflationThreshold(NEP nep,PetscReal *deftol)
     372             : {
     373           1 :   PetscErrorCode ierr;
     374             : 
     375           1 :   PetscFunctionBegin;
     376           1 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     377           1 :   PetscValidRealPointer(deftol,2);
     378           1 :   ierr = PetscUseMethod(nep,"NEPSLPGetDeflationThreshold_C",(NEP,PetscReal*),(nep,deftol));CHKERRQ(ierr);
     379           1 :   PetscFunctionReturn(0);
     380             : }
     381             : 
     382           2 : static PetscErrorCode NEPSLPSetEPS_SLP(NEP nep,EPS eps)
     383             : {
     384           2 :   PetscErrorCode ierr;
     385           2 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     386             : 
     387           2 :   PetscFunctionBegin;
     388           2 :   ierr = PetscObjectReference((PetscObject)eps);CHKERRQ(ierr);
     389           2 :   ierr = EPSDestroy(&ctx->eps);CHKERRQ(ierr);
     390           2 :   ctx->eps = eps;
     391           2 :   ierr = PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->eps);CHKERRQ(ierr);
     392           2 :   nep->state = NEP_STATE_INITIAL;
     393           2 :   PetscFunctionReturn(0);
     394             : }
     395             : 
     396             : /*@
     397             :    NEPSLPSetEPS - Associate a linear eigensolver object (EPS) to the
     398             :    nonlinear eigenvalue solver.
     399             : 
     400             :    Collective on nep
     401             : 
     402             :    Input Parameters:
     403             : +  nep - nonlinear eigenvalue solver
     404             : -  eps - the eigensolver object
     405             : 
     406             :    Level: advanced
     407             : 
     408             : .seealso: NEPSLPGetEPS()
     409             : @*/
     410           2 : PetscErrorCode NEPSLPSetEPS(NEP nep,EPS eps)
     411             : {
     412           2 :   PetscErrorCode ierr;
     413             : 
     414           2 :   PetscFunctionBegin;
     415           2 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     416           2 :   PetscValidHeaderSpecific(eps,EPS_CLASSID,2);
     417           2 :   PetscCheckSameComm(nep,1,eps,2);
     418           2 :   ierr = PetscTryMethod(nep,"NEPSLPSetEPS_C",(NEP,EPS),(nep,eps));CHKERRQ(ierr);
     419           2 :   PetscFunctionReturn(0);
     420             : }
     421             : 
     422          16 : static PetscErrorCode NEPSLPGetEPS_SLP(NEP nep,EPS *eps)
     423             : {
     424          16 :   PetscErrorCode ierr;
     425          16 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     426             : 
     427          16 :   PetscFunctionBegin;
     428          16 :   if (!ctx->eps) {
     429          16 :     ierr = EPSCreate(PetscObjectComm((PetscObject)nep),&ctx->eps);CHKERRQ(ierr);
     430          16 :     ierr = PetscObjectIncrementTabLevel((PetscObject)ctx->eps,(PetscObject)nep,1);CHKERRQ(ierr);
     431          16 :     ierr = EPSSetOptionsPrefix(ctx->eps,((PetscObject)nep)->prefix);CHKERRQ(ierr);
     432          16 :     ierr = EPSAppendOptionsPrefix(ctx->eps,"nep_slp_");CHKERRQ(ierr);
     433          16 :     ierr = PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->eps);CHKERRQ(ierr);
     434          16 :     ierr = PetscObjectSetOptions((PetscObject)ctx->eps,((PetscObject)nep)->options);CHKERRQ(ierr);
     435             :   }
     436          16 :   *eps = ctx->eps;
     437          16 :   PetscFunctionReturn(0);
     438             : }
     439             : 
     440             : /*@
     441             :    NEPSLPGetEPS - Retrieve the linear eigensolver object (EPS) associated
     442             :    to the nonlinear eigenvalue solver.
     443             : 
     444             :    Not Collective
     445             : 
     446             :    Input Parameter:
     447             : .  nep - nonlinear eigenvalue solver
     448             : 
     449             :    Output Parameter:
     450             : .  eps - the eigensolver object
     451             : 
     452             :    Level: advanced
     453             : 
     454             : .seealso: NEPSLPSetEPS()
     455             : @*/
     456          16 : PetscErrorCode NEPSLPGetEPS(NEP nep,EPS *eps)
     457             : {
     458          16 :   PetscErrorCode ierr;
     459             : 
     460          16 :   PetscFunctionBegin;
     461          16 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     462          16 :   PetscValidPointer(eps,2);
     463          16 :   ierr = PetscUseMethod(nep,"NEPSLPGetEPS_C",(NEP,EPS*),(nep,eps));CHKERRQ(ierr);
     464          16 :   PetscFunctionReturn(0);
     465             : }
     466             : 
     467           0 : static PetscErrorCode NEPSLPSetEPSLeft_SLP(NEP nep,EPS eps)
     468             : {
     469           0 :   PetscErrorCode ierr;
     470           0 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     471             : 
     472           0 :   PetscFunctionBegin;
     473           0 :   ierr = PetscObjectReference((PetscObject)eps);CHKERRQ(ierr);
     474           0 :   ierr = EPSDestroy(&ctx->epsts);CHKERRQ(ierr);
     475           0 :   ctx->epsts = eps;
     476           0 :   ierr = PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->epsts);CHKERRQ(ierr);
     477           0 :   nep->state = NEP_STATE_INITIAL;
     478           0 :   PetscFunctionReturn(0);
     479             : }
     480             : 
     481             : /*@
     482             :    NEPSLPSetEPSLeft - Associate a linear eigensolver object (EPS) to the
     483             :    nonlinear eigenvalue solver, used to compute left eigenvectors in the
     484             :    two-sided variant of SLP.
     485             : 
     486             :    Collective on nep
     487             : 
     488             :    Input Parameters:
     489             : +  nep - nonlinear eigenvalue solver
     490             : -  eps - the eigensolver object
     491             : 
     492             :    Level: advanced
     493             : 
     494             : .seealso: NEPSLPGetEPSLeft(), NEPSetTwoSided()
     495             : @*/
     496           0 : PetscErrorCode NEPSLPSetEPSLeft(NEP nep,EPS eps)
     497             : {
     498           0 :   PetscErrorCode ierr;
     499             : 
     500           0 :   PetscFunctionBegin;
     501           0 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     502           0 :   PetscValidHeaderSpecific(eps,EPS_CLASSID,2);
     503           0 :   PetscCheckSameComm(nep,1,eps,2);
     504           0 :   ierr = PetscTryMethod(nep,"NEPSLPSetEPSLeft_C",(NEP,EPS),(nep,eps));CHKERRQ(ierr);
     505           0 :   PetscFunctionReturn(0);
     506             : }
     507             : 
     508           3 : static PetscErrorCode NEPSLPGetEPSLeft_SLP(NEP nep,EPS *eps)
     509             : {
     510           3 :   PetscErrorCode ierr;
     511           3 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     512             : 
     513           3 :   PetscFunctionBegin;
     514           3 :   if (!ctx->epsts) {
     515           3 :     ierr = EPSCreate(PetscObjectComm((PetscObject)nep),&ctx->epsts);CHKERRQ(ierr);
     516           3 :     ierr = PetscObjectIncrementTabLevel((PetscObject)ctx->epsts,(PetscObject)nep,1);CHKERRQ(ierr);
     517           3 :     ierr = EPSSetOptionsPrefix(ctx->epsts,((PetscObject)nep)->prefix);CHKERRQ(ierr);
     518           3 :     ierr = EPSAppendOptionsPrefix(ctx->epsts,"nep_slp_left_");CHKERRQ(ierr);
     519           3 :     ierr = PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->epsts);CHKERRQ(ierr);
     520           3 :     ierr = PetscObjectSetOptions((PetscObject)ctx->epsts,((PetscObject)nep)->options);CHKERRQ(ierr);
     521             :   }
     522           3 :   *eps = ctx->epsts;
     523           3 :   PetscFunctionReturn(0);
     524             : }
     525             : 
     526             : /*@
     527             :    NEPSLPGetEPSLeft - Retrieve the linear eigensolver object (EPS) associated
     528             :    to the nonlinear eigenvalue solver, used to compute left eigenvectors in the
     529             :    two-sided variant of SLP.
     530             : 
     531             :    Not Collective
     532             : 
     533             :    Input Parameter:
     534             : .  nep - nonlinear eigenvalue solver
     535             : 
     536             :    Output Parameter:
     537             : .  eps - the eigensolver object
     538             : 
     539             :    Level: advanced
     540             : 
     541             : .seealso: NEPSLPSetEPSLeft(), NEPSetTwoSided()
     542             : @*/
     543           3 : PetscErrorCode NEPSLPGetEPSLeft(NEP nep,EPS *eps)
     544             : {
     545           3 :   PetscErrorCode ierr;
     546             : 
     547           3 :   PetscFunctionBegin;
     548           3 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     549           3 :   PetscValidPointer(eps,2);
     550           3 :   ierr = PetscUseMethod(nep,"NEPSLPGetEPSLeft_C",(NEP,EPS*),(nep,eps));CHKERRQ(ierr);
     551           3 :   PetscFunctionReturn(0);
     552             : }
     553             : 
     554           2 : static PetscErrorCode NEPSLPSetKSP_SLP(NEP nep,KSP ksp)
     555             : {
     556           2 :   PetscErrorCode ierr;
     557           2 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     558             : 
     559           2 :   PetscFunctionBegin;
     560           2 :   ierr = PetscObjectReference((PetscObject)ksp);CHKERRQ(ierr);
     561           2 :   ierr = KSPDestroy(&ctx->ksp);CHKERRQ(ierr);
     562           2 :   ctx->ksp = ksp;
     563           2 :   ierr = PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);CHKERRQ(ierr);
     564           2 :   nep->state = NEP_STATE_INITIAL;
     565           2 :   PetscFunctionReturn(0);
     566             : }
     567             : 
     568             : /*@
     569             :    NEPSLPSetKSP - Associate a linear solver object (KSP) to the nonlinear
     570             :    eigenvalue solver.
     571             : 
     572             :    Collective on nep
     573             : 
     574             :    Input Parameters:
     575             : +  nep - eigenvalue solver
     576             : -  ksp - the linear solver object
     577             : 
     578             :    Level: advanced
     579             : 
     580             : .seealso: NEPSLPGetKSP()
     581             : @*/
     582           2 : PetscErrorCode NEPSLPSetKSP(NEP nep,KSP ksp)
     583             : {
     584           2 :   PetscErrorCode ierr;
     585             : 
     586           2 :   PetscFunctionBegin;
     587           2 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     588           2 :   PetscValidHeaderSpecific(ksp,KSP_CLASSID,2);
     589           2 :   PetscCheckSameComm(nep,1,ksp,2);
     590           2 :   ierr = PetscTryMethod(nep,"NEPSLPSetKSP_C",(NEP,KSP),(nep,ksp));CHKERRQ(ierr);
     591           2 :   PetscFunctionReturn(0);
     592             : }
     593             : 
     594          16 : static PetscErrorCode NEPSLPGetKSP_SLP(NEP nep,KSP *ksp)
     595             : {
     596          16 :   PetscErrorCode ierr;
     597          16 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     598             : 
     599          16 :   PetscFunctionBegin;
     600          16 :   if (!ctx->ksp) {
     601          16 :     ierr = KSPCreate(PetscObjectComm((PetscObject)nep),&ctx->ksp);CHKERRQ(ierr);
     602          16 :     ierr = PetscObjectIncrementTabLevel((PetscObject)ctx->ksp,(PetscObject)nep,1);CHKERRQ(ierr);
     603          16 :     ierr = KSPSetOptionsPrefix(ctx->ksp,((PetscObject)nep)->prefix);CHKERRQ(ierr);
     604          16 :     ierr = KSPAppendOptionsPrefix(ctx->ksp,"nep_slp_");CHKERRQ(ierr);
     605          16 :     ierr = PetscLogObjectParent((PetscObject)nep,(PetscObject)ctx->ksp);CHKERRQ(ierr);
     606          16 :     ierr = PetscObjectSetOptions((PetscObject)ctx->ksp,((PetscObject)nep)->options);CHKERRQ(ierr);
     607          16 :     ierr = KSPSetErrorIfNotConverged(ctx->ksp,PETSC_TRUE);CHKERRQ(ierr);
     608          25 :     ierr = KSPSetTolerances(ctx->ksp,SlepcDefaultTol(nep->tol),PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT);CHKERRQ(ierr);
     609             :   }
     610          16 :   *ksp = ctx->ksp;
     611          16 :   PetscFunctionReturn(0);
     612             : }
     613             : 
     614             : /*@
     615             :    NEPSLPGetKSP - Retrieve the linear solver object (KSP) associated with
     616             :    the nonlinear eigenvalue solver.
     617             : 
     618             :    Not Collective
     619             : 
     620             :    Input Parameter:
     621             : .  nep - nonlinear eigenvalue solver
     622             : 
     623             :    Output Parameter:
     624             : .  ksp - the linear solver object
     625             : 
     626             :    Level: advanced
     627             : 
     628             : .seealso: NEPSLPSetKSP()
     629             : @*/
     630          16 : PetscErrorCode NEPSLPGetKSP(NEP nep,KSP *ksp)
     631             : {
     632          16 :   PetscErrorCode ierr;
     633             : 
     634          16 :   PetscFunctionBegin;
     635          16 :   PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
     636          16 :   PetscValidPointer(ksp,2);
     637          16 :   ierr = PetscUseMethod(nep,"NEPSLPGetKSP_C",(NEP,KSP*),(nep,ksp));CHKERRQ(ierr);
     638          16 :   PetscFunctionReturn(0);
     639             : }
     640             : 
     641           1 : PetscErrorCode NEPView_SLP(NEP nep,PetscViewer viewer)
     642             : {
     643           1 :   PetscErrorCode ierr;
     644           1 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     645           1 :   PetscBool      isascii;
     646             : 
     647           1 :   PetscFunctionBegin;
     648           1 :   ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);CHKERRQ(ierr);
     649           1 :   if (isascii) {
     650           1 :     if (ctx->deftol) {
     651           1 :       ierr = PetscViewerASCIIPrintf(viewer,"  deflation threshold: %g\n",(double)ctx->deftol);CHKERRQ(ierr);
     652             :     }
     653           1 :     if (!ctx->eps) { ierr = NEPSLPGetEPS(nep,&ctx->eps);CHKERRQ(ierr); }
     654           1 :     ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
     655           1 :     ierr = EPSView(ctx->eps,viewer);CHKERRQ(ierr);
     656           1 :     if (nep->twosided) {
     657           0 :       if (!ctx->epsts) { ierr = NEPSLPGetEPSLeft(nep,&ctx->epsts);CHKERRQ(ierr); }
     658           0 :       ierr = EPSView(ctx->epsts,viewer);CHKERRQ(ierr);
     659             :     }
     660           1 :     if (!ctx->ksp) { ierr = NEPSLPGetKSP(nep,&ctx->ksp);CHKERRQ(ierr); }
     661           1 :     ierr = KSPView(ctx->ksp,viewer);CHKERRQ(ierr);
     662           1 :     ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
     663             :   }
     664           1 :   PetscFunctionReturn(0);
     665             : }
     666             : 
     667          20 : PetscErrorCode NEPReset_SLP(NEP nep)
     668             : {
     669          20 :   PetscErrorCode ierr;
     670          20 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     671             : 
     672          20 :   PetscFunctionBegin;
     673          20 :   ierr = EPSReset(ctx->eps);CHKERRQ(ierr);
     674          20 :   if (nep->twosided) { ierr = EPSReset(ctx->epsts);CHKERRQ(ierr); }
     675          20 :   ierr = KSPReset(ctx->ksp);CHKERRQ(ierr);
     676          20 :   PetscFunctionReturn(0);
     677             : }
     678             : 
     679          18 : PetscErrorCode NEPDestroy_SLP(NEP nep)
     680             : {
     681          18 :   PetscErrorCode ierr;
     682          18 :   NEP_SLP        *ctx = (NEP_SLP*)nep->data;
     683             : 
     684          18 :   PetscFunctionBegin;
     685          18 :   ierr = KSPDestroy(&ctx->ksp);CHKERRQ(ierr);
     686          18 :   ierr = EPSDestroy(&ctx->eps);CHKERRQ(ierr);
     687          18 :   ierr = EPSDestroy(&ctx->epsts);CHKERRQ(ierr);
     688          18 :   ierr = PetscFree(nep->data);CHKERRQ(ierr);
     689          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetDeflationThreshold_C",NULL);CHKERRQ(ierr);
     690          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetDeflationThreshold_C",NULL);CHKERRQ(ierr);
     691          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPS_C",NULL);CHKERRQ(ierr);
     692          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPS_C",NULL);CHKERRQ(ierr);
     693          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPSLeft_C",NULL);CHKERRQ(ierr);
     694          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPSLeft_C",NULL);CHKERRQ(ierr);
     695          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetKSP_C",NULL);CHKERRQ(ierr);
     696          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetKSP_C",NULL);CHKERRQ(ierr);
     697          18 :   PetscFunctionReturn(0);
     698             : }
     699             : 
     700          18 : SLEPC_EXTERN PetscErrorCode NEPCreate_SLP(NEP nep)
     701             : {
     702          18 :   PetscErrorCode ierr;
     703          18 :   NEP_SLP        *ctx;
     704             : 
     705          18 :   PetscFunctionBegin;
     706          18 :   ierr = PetscNewLog(nep,&ctx);CHKERRQ(ierr);
     707          18 :   nep->data = (void*)ctx;
     708             : 
     709          18 :   nep->useds  = PETSC_TRUE;
     710          18 :   ctx->deftol = PETSC_DEFAULT;
     711             : 
     712          18 :   nep->ops->solve          = NEPSolve_SLP;
     713          18 :   nep->ops->setup          = NEPSetUp_SLP;
     714          18 :   nep->ops->setfromoptions = NEPSetFromOptions_SLP;
     715          18 :   nep->ops->reset          = NEPReset_SLP;
     716          18 :   nep->ops->destroy        = NEPDestroy_SLP;
     717          18 :   nep->ops->view           = NEPView_SLP;
     718          18 :   nep->ops->computevectors = NEPComputeVectors_Schur;
     719             : 
     720          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetDeflationThreshold_C",NEPSLPSetDeflationThreshold_SLP);CHKERRQ(ierr);
     721          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetDeflationThreshold_C",NEPSLPGetDeflationThreshold_SLP);CHKERRQ(ierr);
     722          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPS_C",NEPSLPSetEPS_SLP);CHKERRQ(ierr);
     723          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPS_C",NEPSLPGetEPS_SLP);CHKERRQ(ierr);
     724          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetEPSLeft_C",NEPSLPSetEPSLeft_SLP);CHKERRQ(ierr);
     725          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetEPSLeft_C",NEPSLPGetEPSLeft_SLP);CHKERRQ(ierr);
     726          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPSetKSP_C",NEPSLPSetKSP_SLP);CHKERRQ(ierr);
     727          18 :   ierr = PetscObjectComposeFunction((PetscObject)nep,"NEPSLPGetKSP_C",NEPSLPGetKSP_SLP);CHKERRQ(ierr);
     728          18 :   PetscFunctionReturn(0);
     729             : }
     730             : 

Generated by: LCOV version 1.14