LCOV - code coverage report
Current view: top level - eps/impls/external/blopex - blopex.c (source / functions) Hit Total Coverage
Test: coverage.info Lines: 196 209 93.8 %
Date: 2020-05-28 03:14:05 Functions: 18 18 100.0 %

          Line data    Source code
       1             : /*
       2             :    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
       3             :    SLEPc - Scalable Library for Eigenvalue Problem Computations
       4             :    Copyright (c) 2002-2020, 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             :    This file implements a wrapper to the BLOPEX package
      12             : */
      13             : 
      14             : #include <slepc/private/epsimpl.h>                /*I "slepceps.h" I*/
      15             : #include "blopex.h"
      16             : #include <lobpcg.h>
      17             : #include <interpreter.h>
      18             : #include <multivector.h>
      19             : #include <temp_multivector.h>
      20             : 
      21             : PetscInt slepc_blopex_useconstr = -1;
      22             : 
      23             : typedef struct {
      24             :   lobpcg_Tolerance           tol;
      25             :   lobpcg_BLASLAPACKFunctions blap_fn;
      26             :   mv_InterfaceInterpreter    ii;
      27             :   ST                         st;
      28             :   Vec                        w;
      29             :   PetscInt                   bs;     /* block size */
      30             : } EPS_BLOPEX;
      31             : 
      32         760 : static void Precond_FnSingleVector(void *data,void *x,void *y)
      33             : {
      34             :   PetscErrorCode ierr;
      35         760 :   EPS_BLOPEX     *blopex = (EPS_BLOPEX*)data;
      36         760 :   MPI_Comm       comm = PetscObjectComm((PetscObject)blopex->st);
      37             :   KSP            ksp;
      38             : 
      39         760 :   PetscFunctionBegin;
      40         760 :   ierr = STGetKSP(blopex->st,&ksp);CHKERRABORT(comm,ierr);
      41         760 :   ierr = KSPSolve(ksp,(Vec)x,(Vec)y);CHKERRABORT(comm,ierr);
      42         760 :   PetscFunctionReturnVoid();
      43             : }
      44             : 
      45         264 : static void Precond_FnMultiVector(void *data,void *x,void *y)
      46             : {
      47         264 :   EPS_BLOPEX *blopex = (EPS_BLOPEX*)data;
      48             : 
      49         264 :   PetscFunctionBegin;
      50         264 :   blopex->ii.Eval(Precond_FnSingleVector,data,x,y);
      51         264 :   PetscFunctionReturnVoid();
      52             : }
      53             : 
      54         795 : static void OperatorASingleVector(void *data,void *x,void *y)
      55             : {
      56             :   PetscErrorCode ierr;
      57         795 :   EPS_BLOPEX     *blopex = (EPS_BLOPEX*)data;
      58         795 :   MPI_Comm       comm = PetscObjectComm((PetscObject)blopex->st);
      59             :   Mat            A,B;
      60             :   PetscScalar    sigma;
      61             :   PetscInt       nmat;
      62             : 
      63         795 :   PetscFunctionBegin;
      64         795 :   ierr = STGetNumMatrices(blopex->st,&nmat);CHKERRABORT(comm,ierr);
      65         795 :   ierr = STGetMatrix(blopex->st,0,&A);CHKERRABORT(comm,ierr);
      66         795 :   if (nmat>1) { ierr = STGetMatrix(blopex->st,1,&B);CHKERRABORT(comm,ierr); }
      67         795 :   ierr = MatMult(A,(Vec)x,(Vec)y);CHKERRABORT(comm,ierr);
      68         795 :   ierr = STGetShift(blopex->st,&sigma);CHKERRABORT(comm,ierr);
      69         795 :   if (sigma != 0.0) {
      70         217 :     if (nmat>1) {
      71          49 :       ierr = MatMult(B,(Vec)x,blopex->w);CHKERRABORT(comm,ierr);
      72             :     } else {
      73         168 :       ierr = VecCopy((Vec)x,blopex->w);CHKERRABORT(comm,ierr);
      74             :     }
      75         217 :     ierr = VecAXPY((Vec)y,-sigma,blopex->w);CHKERRABORT(comm,ierr);
      76             :   }
      77         795 :   PetscFunctionReturnVoid();
      78             : }
      79             : 
      80         273 : static void OperatorAMultiVector(void *data,void *x,void *y)
      81             : {
      82         273 :   EPS_BLOPEX *blopex = (EPS_BLOPEX*)data;
      83             : 
      84         273 :   PetscFunctionBegin;
      85         273 :   blopex->ii.Eval(OperatorASingleVector,data,x,y);
      86         273 :   PetscFunctionReturnVoid();
      87             : }
      88             : 
      89         478 : static void OperatorBSingleVector(void *data,void *x,void *y)
      90             : {
      91             :   PetscErrorCode ierr;
      92         478 :   EPS_BLOPEX     *blopex = (EPS_BLOPEX*)data;
      93         478 :   MPI_Comm       comm = PetscObjectComm((PetscObject)blopex->st);
      94             :   Mat            B;
      95             : 
      96         478 :   PetscFunctionBegin;
      97         478 :   ierr = STGetMatrix(blopex->st,1,&B);CHKERRABORT(comm,ierr);
      98         478 :   ierr = MatMult(B,(Vec)x,(Vec)y);CHKERRABORT(comm,ierr);
      99         478 :   PetscFunctionReturnVoid();
     100             : }
     101             : 
     102         171 : static void OperatorBMultiVector(void *data,void *x,void *y)
     103             : {
     104         171 :   EPS_BLOPEX *blopex = (EPS_BLOPEX*)data;
     105             : 
     106         171 :   PetscFunctionBegin;
     107         171 :   blopex->ii.Eval(OperatorBSingleVector,data,x,y);
     108         171 :   PetscFunctionReturnVoid();
     109             : }
     110             : 
     111           7 : PetscErrorCode EPSSetDimensions_BLOPEX(EPS eps,PetscInt nev,PetscInt *ncv,PetscInt *mpd)
     112             : {
     113             :   PetscErrorCode ierr;
     114           7 :   EPS_BLOPEX     *ctx = (EPS_BLOPEX*)eps->data;
     115             :   PetscInt       k;
     116             : 
     117           7 :   PetscFunctionBegin;
     118           7 :   k = ((eps->nev-1)/ctx->bs+1)*ctx->bs;
     119           7 :   if (*ncv!=PETSC_DEFAULT) { /* ncv set */
     120           1 :     if (*ncv<k) SETERRQ(PetscObjectComm((PetscObject)eps),1,"The value of ncv is not sufficiently large");
     121           6 :   } else *ncv = k;
     122           7 :   if (*mpd==PETSC_DEFAULT) *mpd = *ncv;
     123           1 :   else { ierr = PetscInfo(eps,"Warning: given value of mpd ignored\n");CHKERRQ(ierr); }
     124           7 :   PetscFunctionReturn(0);
     125             : }
     126             : 
     127           7 : PetscErrorCode EPSSetUp_BLOPEX(EPS eps)
     128             : {
     129             :   PetscErrorCode ierr;
     130           7 :   EPS_BLOPEX     *blopex = (EPS_BLOPEX*)eps->data;
     131             :   PetscBool      flg;
     132             :   KSP            ksp;
     133             : 
     134           7 :   PetscFunctionBegin;
     135           7 :   EPSCheckHermitianDefinite(eps);
     136           7 :   if (!blopex->bs) blopex->bs = PetscMin(16,eps->nev);
     137           7 :   ierr = EPSSetDimensions_BLOPEX(eps,eps->nev,&eps->ncv,&eps->mpd);CHKERRQ(ierr);
     138           7 :   if (eps->max_it==PETSC_DEFAULT) eps->max_it = PetscMax(100,2*eps->n/eps->ncv);
     139           7 :   if (!eps->which) eps->which = EPS_SMALLEST_REAL;
     140           7 :   if (eps->which!=EPS_SMALLEST_REAL) SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"This solver supports only smallest real eigenvalues");
     141           7 :   EPSCheckUnsupported(eps,EPS_FEATURE_ARBITRARY | EPS_FEATURE_REGION | EPS_FEATURE_STOPPING);
     142           7 :   EPSCheckIgnored(eps,EPS_FEATURE_BALANCE | EPS_FEATURE_EXTRACTION);
     143             : 
     144           7 :   blopex->st = eps->st;
     145             : 
     146           7 :   if (eps->converged == EPSConvergedRelative) {
     147           4 :     blopex->tol.absolute = 0.0;
     148           4 :     blopex->tol.relative = eps->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:eps->tol;
     149           3 :   } else if (eps->converged == EPSConvergedAbsolute) {
     150           3 :     blopex->tol.absolute = eps->tol==PETSC_DEFAULT?SLEPC_DEFAULT_TOL:eps->tol;
     151           3 :     blopex->tol.relative = 0.0;
     152           0 :   } else SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"Convergence test not supported in this solver");
     153             : 
     154           7 :   SLEPCSetupInterpreter(&blopex->ii);
     155             : 
     156           7 :   ierr = STGetKSP(eps->st,&ksp);CHKERRQ(ierr);
     157           7 :   ierr = PetscObjectTypeCompare((PetscObject)ksp,KSPPREONLY,&flg);CHKERRQ(ierr);
     158           7 :   if (!flg) { ierr = PetscInfo(eps,"Warning: ignoring KSP, should use KSPPREONLY\n");CHKERRQ(ierr); }
     159             : 
     160             :   /* allocate memory */
     161           7 :   if (!eps->V) { ierr = EPSGetBV(eps,&eps->V);CHKERRQ(ierr); }
     162           7 :   ierr = PetscObjectTypeCompareAny((PetscObject)eps->V,&flg,BVVECS,BVCONTIGUOUS,"");CHKERRQ(ierr);
     163           7 :   if (!flg) {  /* blopex only works with BVVECS or BVCONTIGUOUS */
     164           6 :     ierr = BVSetType(eps->V,BVCONTIGUOUS);CHKERRQ(ierr);
     165             :   }
     166           7 :   ierr = EPSAllocateSolution(eps,0);CHKERRQ(ierr);
     167           7 :   if (!blopex->w) {
     168           6 :     ierr = BVCreateVec(eps->V,&blopex->w);CHKERRQ(ierr);
     169           6 :     ierr = PetscLogObjectParent((PetscObject)eps,(PetscObject)blopex->w);CHKERRQ(ierr);
     170             :   }
     171             : 
     172             : #if defined(PETSC_USE_COMPLEX)
     173             :   blopex->blap_fn.zpotrf = PETSC_zpotrf_interface;
     174             :   blopex->blap_fn.zhegv = PETSC_zsygv_interface;
     175             : #else
     176           7 :   blopex->blap_fn.dpotrf = PETSC_dpotrf_interface;
     177           7 :   blopex->blap_fn.dsygv = PETSC_dsygv_interface;
     178             : #endif
     179           7 :   PetscFunctionReturn(0);
     180             : }
     181             : 
     182           7 : PetscErrorCode EPSSolve_BLOPEX(EPS eps)
     183             : {
     184           7 :   EPS_BLOPEX        *blopex = (EPS_BLOPEX*)eps->data;
     185           7 :   PetscScalar       sigma,*eigr=NULL;
     186           7 :   PetscReal         *errest=NULL;
     187             :   int               i,j,info,its,nconv;
     188           7 :   double            *residhist=NULL;
     189             :   PetscErrorCode    ierr;
     190             :   mv_MultiVectorPtr eigenvectors,constraints;
     191             : #if defined(PETSC_USE_COMPLEX)
     192             :   komplex           *lambda=NULL,*lambdahist=NULL;
     193             : #else
     194           7 :   double            *lambda=NULL,*lambdahist=NULL;
     195             : #endif
     196             : 
     197           7 :   PetscFunctionBegin;
     198           7 :   ierr = STGetShift(eps->st,&sigma);CHKERRQ(ierr);
     199           7 :   ierr = PetscMalloc1(blopex->bs,&lambda);CHKERRQ(ierr);
     200           7 :   if (eps->numbermonitors>0) {
     201           0 :     ierr = PetscMalloc4(blopex->bs*(eps->max_it+1),&lambdahist,eps->ncv,&eigr,blopex->bs*(eps->max_it+1),&residhist,eps->ncv,&errest);CHKERRQ(ierr);
     202             :   }
     203             : 
     204             :   /* Complete the initial basis with random vectors */
     205           1 :   for (i=0;i<eps->nini;i++) {  /* in case the initial vectors were also set with VecSetRandom */
     206           1 :     ierr = BVSetRandomColumn(eps->V,eps->nini);CHKERRQ(ierr);
     207             :   }
     208          34 :   for (i=eps->nini;i<eps->ncv;i++) {
     209          34 :     ierr = BVSetRandomColumn(eps->V,i);CHKERRQ(ierr);
     210             :   }
     211             : 
     212          16 :   while (eps->reason == EPS_CONVERGED_ITERATING) {
     213             : 
     214             :     /* Create multivector of constraints from leading columns of V */
     215           9 :     ierr = PetscObjectComposedDataSetInt((PetscObject)eps->V,slepc_blopex_useconstr,1);CHKERRQ(ierr);
     216           9 :     ierr = BVSetActiveColumns(eps->V,0,eps->nconv);CHKERRQ(ierr);
     217           9 :     constraints = mv_MultiVectorCreateFromSampleVector(&blopex->ii,eps->nds+eps->nconv,eps->V);
     218             : 
     219             :     /* Create multivector where eigenvectors of this run will be stored */
     220           9 :     ierr = PetscObjectComposedDataSetInt((PetscObject)eps->V,slepc_blopex_useconstr,0);CHKERRQ(ierr);
     221           9 :     ierr = BVSetActiveColumns(eps->V,eps->nconv,eps->nconv+blopex->bs);CHKERRQ(ierr);
     222           9 :     eigenvectors = mv_MultiVectorCreateFromSampleVector(&blopex->ii,blopex->bs,eps->V);
     223             : 
     224             : #if defined(PETSC_USE_COMPLEX)
     225             :     info = lobpcg_solve_complex(eigenvectors,blopex,OperatorAMultiVector,
     226             :           eps->isgeneralized?blopex:NULL,eps->isgeneralized?OperatorBMultiVector:NULL,
     227             :           blopex,Precond_FnMultiVector,constraints,
     228             :           blopex->blap_fn,blopex->tol,eps->max_it,0,&its,
     229             :           lambda,lambdahist,blopex->bs,eps->errest+eps->nconv,residhist,blopex->bs);
     230             : #else
     231          27 :     info = lobpcg_solve_double(eigenvectors,blopex,OperatorAMultiVector,
     232           9 :           eps->isgeneralized?blopex:NULL,eps->isgeneralized?OperatorBMultiVector:NULL,
     233             :           blopex,Precond_FnMultiVector,constraints,
     234             :           blopex->blap_fn,blopex->tol,eps->max_it,0,&its,
     235           9 :           lambda,lambdahist,blopex->bs,eps->errest+eps->nconv,residhist,blopex->bs);
     236             : #endif
     237           9 :     if (info>0) SETERRQ1(PetscObjectComm((PetscObject)eps),PETSC_ERR_LIB,"BLOPEX failed with exit code=%d",info);
     238           9 :     mv_MultiVectorDestroy(constraints);
     239           9 :     mv_MultiVectorDestroy(eigenvectors);
     240             : 
     241          44 :     for (j=0;j<blopex->bs;j++) {
     242             : #if defined(PETSC_USE_COMPLEX)
     243             :       eps->eigr[eps->nconv+j] = PetscCMPLX(lambda[j].real,lambda[j].imag);
     244             : #else
     245          35 :       eps->eigr[eps->nconv+j] = lambda[j];
     246             : #endif
     247             :     }
     248             : 
     249           9 :     if (eps->numbermonitors>0) {
     250           0 :       for (i=0;i<its;i++) {
     251             :         nconv = 0;
     252           0 :         for (j=0;j<blopex->bs;j++) {
     253             : #if defined(PETSC_USE_COMPLEX)
     254             :           eigr[eps->nconv+j] = PetscCMPLX(lambdahist[j+i*blopex->bs].real,lambdahist[j+i*blopex->bs].imag);
     255             : #else
     256           0 :           eigr[eps->nconv+j] = lambdahist[j+i*blopex->bs];
     257             : #endif
     258           0 :           errest[eps->nconv+j] = residhist[j+i*blopex->bs];
     259           0 :           if (residhist[j+i*blopex->bs]<=eps->tol) nconv++;
     260             :         }
     261           0 :         ierr = EPSMonitor(eps,eps->its+i,eps->nconv+nconv,eigr,eps->eigi,errest,eps->nconv+blopex->bs);CHKERRQ(ierr);
     262             :       }
     263             :     }
     264             : 
     265           9 :     eps->its += its;
     266           9 :     if (info==-1) {
     267           0 :       eps->reason = EPS_DIVERGED_ITS;
     268           0 :       break;
     269             :     } else {
     270          35 :       for (i=0;i<blopex->bs;i++) {
     271          35 :         if (sigma != 0.0) eps->eigr[eps->nconv+i] += sigma;
     272             :       }
     273           9 :       eps->nconv += blopex->bs;
     274           9 :       if (eps->nconv>=eps->nev) eps->reason = EPS_CONVERGED_TOL;
     275             :     }
     276             :   }
     277             : 
     278           7 :   ierr = PetscFree(lambda);CHKERRQ(ierr);
     279           7 :   if (eps->numbermonitors>0) {
     280           0 :     ierr = PetscFree4(lambdahist,eigr,residhist,errest);CHKERRQ(ierr);
     281             :   }
     282           7 :   PetscFunctionReturn(0);
     283             : }
     284             : 
     285           2 : static PetscErrorCode EPSBLOPEXSetBlockSize_BLOPEX(EPS eps,PetscInt bs)
     286             : {
     287           2 :   EPS_BLOPEX *ctx = (EPS_BLOPEX*)eps->data;
     288             : 
     289           2 :   PetscFunctionBegin;
     290           2 :   if (bs==PETSC_DEFAULT) {
     291           0 :     ctx->bs    = 0;
     292           0 :     eps->state = EPS_STATE_INITIAL;
     293             :   } else {
     294           2 :     if (bs<=0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Block size must be >0");
     295           2 :     ctx->bs = bs;
     296             :   }
     297           2 :   PetscFunctionReturn(0);
     298             : }
     299             : 
     300             : /*@
     301             :    EPSBLOPEXSetBlockSize - Sets the block size of the BLOPEX solver.
     302             : 
     303             :    Logically Collective on eps
     304             : 
     305             :    Input Parameters:
     306             : +  eps - the eigenproblem solver context
     307             : -  bs  - the block size
     308             : 
     309             :    Options Database Key:
     310             : .  -eps_blopex_blocksize - Sets the block size
     311             : 
     312             :    Level: advanced
     313             : 
     314             : .seealso: EPSBLOPEXGetBlockSize()
     315             : @*/
     316           2 : PetscErrorCode EPSBLOPEXSetBlockSize(EPS eps,PetscInt bs)
     317             : {
     318             :   PetscErrorCode ierr;
     319             : 
     320           2 :   PetscFunctionBegin;
     321           2 :   PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
     322           4 :   PetscValidLogicalCollectiveInt(eps,bs,2);
     323           2 :   ierr = PetscTryMethod(eps,"EPSBLOPEXSetBlockSize_C",(EPS,PetscInt),(eps,bs));CHKERRQ(ierr);
     324           4 :   PetscFunctionReturn(0);
     325             : }
     326             : 
     327           1 : static PetscErrorCode EPSBLOPEXGetBlockSize_BLOPEX(EPS eps,PetscInt *bs)
     328             : {
     329           1 :   EPS_BLOPEX *ctx = (EPS_BLOPEX*)eps->data;
     330             : 
     331           1 :   PetscFunctionBegin;
     332           1 :   *bs = ctx->bs;
     333           1 :   PetscFunctionReturn(0);
     334             : }
     335             : 
     336             : /*@
     337             :    EPSBLOPEXGetBlockSize - Gets the block size used in the BLOPEX solver.
     338             : 
     339             :    Not Collective
     340             : 
     341             :    Input Parameter:
     342             : .  eps - the eigenproblem solver context
     343             : 
     344             :    Output Parameter:
     345             : .  bs - the block size
     346             : 
     347             :    Level: advanced
     348             : 
     349             : .seealso: EPSBLOPEXSetBlockSize()
     350             : @*/
     351           1 : PetscErrorCode EPSBLOPEXGetBlockSize(EPS eps,PetscInt *bs)
     352             : {
     353             :   PetscErrorCode ierr;
     354             : 
     355           1 :   PetscFunctionBegin;
     356           1 :   PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
     357           1 :   PetscValidIntPointer(bs,2);
     358           1 :   ierr = PetscUseMethod(eps,"EPSBLOPEXGetBlockSize_C",(EPS,PetscInt*),(eps,bs));CHKERRQ(ierr);
     359           2 :   PetscFunctionReturn(0);
     360             : }
     361             : 
     362           6 : PetscErrorCode EPSReset_BLOPEX(EPS eps)
     363             : {
     364             :   PetscErrorCode ierr;
     365           6 :   EPS_BLOPEX     *blopex = (EPS_BLOPEX*)eps->data;
     366             : 
     367           6 :   PetscFunctionBegin;
     368           6 :   ierr = VecDestroy(&blopex->w);CHKERRQ(ierr);
     369           6 :   PetscFunctionReturn(0);
     370             : }
     371             : 
     372           6 : PetscErrorCode EPSDestroy_BLOPEX(EPS eps)
     373             : {
     374             :   PetscErrorCode ierr;
     375             : 
     376           6 :   PetscFunctionBegin;
     377           6 :   LOBPCG_DestroyRandomContext();
     378           6 :   ierr = PetscFree(eps->data);CHKERRQ(ierr);
     379           6 :   ierr = PetscObjectComposeFunction((PetscObject)eps,"EPSBLOPEXSetBlockSize_C",NULL);CHKERRQ(ierr);
     380           6 :   ierr = PetscObjectComposeFunction((PetscObject)eps,"EPSBLOPEXGetBlockSize_C",NULL);CHKERRQ(ierr);
     381           6 :   PetscFunctionReturn(0);
     382             : }
     383             : 
     384           1 : PetscErrorCode EPSView_BLOPEX(EPS eps,PetscViewer viewer)
     385             : {
     386             :   PetscErrorCode ierr;
     387           1 :   EPS_BLOPEX     *ctx = (EPS_BLOPEX*)eps->data;
     388             :   PetscBool      isascii;
     389             : 
     390           1 :   PetscFunctionBegin;
     391           1 :   ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);CHKERRQ(ierr);
     392           1 :   if (isascii) {
     393           1 :     ierr = PetscViewerASCIIPrintf(viewer,"  block size %D\n",ctx->bs);CHKERRQ(ierr);
     394             :   }
     395           1 :   PetscFunctionReturn(0);
     396             : }
     397             : 
     398           6 : PetscErrorCode EPSSetFromOptions_BLOPEX(PetscOptionItems *PetscOptionsObject,EPS eps)
     399             : {
     400             :   PetscErrorCode ierr;
     401             :   PetscBool      flg;
     402             :   PetscInt       bs;
     403             : 
     404           6 :   PetscFunctionBegin;
     405           6 :   ierr = PetscOptionsHead(PetscOptionsObject,"EPS BLOPEX Options");CHKERRQ(ierr);
     406             : 
     407           6 :     ierr = PetscOptionsInt("-eps_blopex_blocksize","Block size","EPSBLOPEXSetBlockSize",20,&bs,&flg);CHKERRQ(ierr);
     408           6 :     if (flg) { ierr = EPSBLOPEXSetBlockSize(eps,bs);CHKERRQ(ierr); }
     409             : 
     410           6 :   ierr = PetscOptionsTail();CHKERRQ(ierr);
     411             : 
     412           6 :   LOBPCG_SetFromOptionsRandomContext();
     413           6 :   PetscFunctionReturn(0);
     414             : }
     415             : 
     416           6 : SLEPC_EXTERN PetscErrorCode EPSCreate_BLOPEX(EPS eps)
     417             : {
     418             :   EPS_BLOPEX     *ctx;
     419             :   PetscErrorCode ierr;
     420             : 
     421           6 :   PetscFunctionBegin;
     422           6 :   ierr = PetscNewLog(eps,&ctx);CHKERRQ(ierr);
     423           6 :   eps->data = (void*)ctx;
     424             : 
     425           6 :   eps->categ = EPS_CATEGORY_PRECOND;
     426             : 
     427           6 :   eps->ops->solve          = EPSSolve_BLOPEX;
     428           6 :   eps->ops->setup          = EPSSetUp_BLOPEX;
     429           6 :   eps->ops->setupsort      = EPSSetUpSort_Basic;
     430           6 :   eps->ops->setfromoptions = EPSSetFromOptions_BLOPEX;
     431           6 :   eps->ops->destroy        = EPSDestroy_BLOPEX;
     432           6 :   eps->ops->reset          = EPSReset_BLOPEX;
     433           6 :   eps->ops->view           = EPSView_BLOPEX;
     434           6 :   eps->ops->backtransform  = EPSBackTransform_Default;
     435           6 :   eps->ops->setdefaultst   = EPSSetDefaultST_GMRES;
     436             : 
     437           6 :   LOBPCG_InitRandomContext(PetscObjectComm((PetscObject)eps),NULL);
     438           6 :   ierr = PetscObjectComposeFunction((PetscObject)eps,"EPSBLOPEXSetBlockSize_C",EPSBLOPEXSetBlockSize_BLOPEX);CHKERRQ(ierr);
     439           6 :   ierr = PetscObjectComposeFunction((PetscObject)eps,"EPSBLOPEXGetBlockSize_C",EPSBLOPEXGetBlockSize_BLOPEX);CHKERRQ(ierr);
     440           6 :   if (slepc_blopex_useconstr < 0) { ierr = PetscObjectComposedDataRegister(&slepc_blopex_useconstr);CHKERRQ(ierr); }
     441           6 :   PetscFunctionReturn(0);
     442             : }
     443             : 

Generated by: LCOV version 1.13