LCOV - code coverage report
Current view: top level - eps/impls/davidson - dvdupdatev.c (source / functions) Hit Total Coverage
Test: coverage.info Lines: 177 186 95.2 %
Date: 2019-10-20 05:41:02 Functions: 9 9 100.0 %

          Line data    Source code
       1             : /*
       2             :    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
       3             :    SLEPc - Scalable Library for Eigenvalue Problem Computations
       4             :    Copyright (c) 2002-2019, 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 eigensolver: "davidson"
      12             : 
      13             :    Step: test for restarting, updateV, restartV
      14             : */
      15             : 
      16             : #include "davidson.h"
      17             : 
      18             : typedef struct {
      19             :   PetscInt          min_size_V;        /* restart with this number of eigenvectors */
      20             :   PetscInt          plusk;             /* at restart, save plusk vectors from last iteration */
      21             :   PetscInt          mpd;               /* max size of the searching subspace */
      22             :   void              *old_updateV_data; /* old updateV data */
      23             :   PetscErrorCode    (*old_isRestarting)(dvdDashboard*,PetscBool*);  /* old isRestarting */
      24             :   Mat               oldU;              /* previous projected right igenvectors */
      25             :   Mat               oldV;              /* previous projected left eigenvectors */
      26             :   PetscInt          size_oldU;         /* size of oldU */
      27             :   PetscBool         allResiduals;      /* if computing all the residuals */
      28             : } dvdManagV_basic;
      29             : 
      30          90 : static PetscErrorCode dvd_updateV_start(dvdDashboard *d)
      31             : {
      32          90 :   dvdManagV_basic *data = (dvdManagV_basic*)d->updateV_data;
      33             :   PetscInt        i;
      34             : 
      35          90 :   PetscFunctionBegin;
      36        1743 :   for (i=0;i<d->eps->ncv;i++) d->eigi[i] = 0.0;
      37          90 :   d->nR = d->real_nR;
      38          90 :   for (i=0;i<d->eps->ncv;i++) d->nR[i] = 1.0;
      39          90 :   d->nX = d->real_nX;
      40          90 :   for (i=0;i<d->eps->ncv;i++) d->errest[i] = 1.0;
      41          90 :   data->size_oldU = 0;
      42          90 :   d->nconv = 0;
      43          90 :   d->npreconv = 0;
      44          90 :   d->V_tra_s = d->V_tra_e = d->V_new_s = d->V_new_e = 0;
      45          90 :   d->size_D = 0;
      46          90 :   PetscFunctionReturn(0);
      47             : }
      48             : 
      49        8873 : static PetscErrorCode dvd_isrestarting_fullV(dvdDashboard *d,PetscBool *r)
      50             : {
      51             :   PetscErrorCode  ierr;
      52             :   PetscInt        l,k;
      53             :   PetscBool       restart;
      54        8873 :   dvdManagV_basic *data = (dvdManagV_basic*)d->updateV_data;
      55             : 
      56        8873 :   PetscFunctionBegin;
      57        8873 :   ierr = BVGetActiveColumns(d->eps->V,&l,&k);CHKERRQ(ierr);
      58        8873 :   restart = (k+2 > d->eps->ncv)? PETSC_TRUE: PETSC_FALSE;
      59             : 
      60             :   /* Check old isRestarting function */
      61        8873 :   if (!restart && data->old_isRestarting) {
      62           0 :     ierr = data->old_isRestarting(d,&restart);CHKERRQ(ierr);
      63             :   }
      64        8873 :   *r = restart;
      65        8873 :   PetscFunctionReturn(0);
      66             : }
      67             : 
      68          90 : static PetscErrorCode dvd_managementV_basic_d(dvdDashboard *d)
      69             : {
      70             :   PetscErrorCode  ierr;
      71          90 :   dvdManagV_basic *data = (dvdManagV_basic*)d->updateV_data;
      72             : 
      73          90 :   PetscFunctionBegin;
      74             :   /* Restore changes in dvdDashboard */
      75          90 :   d->updateV_data = data->old_updateV_data;
      76             : 
      77             :   /* Free local data */
      78          90 :   ierr = MatDestroy(&data->oldU);CHKERRQ(ierr);
      79          90 :   ierr = MatDestroy(&data->oldV);CHKERRQ(ierr);
      80          90 :   ierr = PetscFree(d->real_nR);CHKERRQ(ierr);
      81          90 :   ierr = PetscFree(d->real_nX);CHKERRQ(ierr);
      82          90 :   ierr = PetscFree(data);CHKERRQ(ierr);
      83          90 :   PetscFunctionReturn(0);
      84             : }
      85             : 
      86         299 : static PetscErrorCode dvd_updateV_conv_gen(dvdDashboard *d)
      87             : {
      88         299 :   dvdManagV_basic *data = (dvdManagV_basic*)d->updateV_data;
      89             :   PetscInt        npreconv,cMT,cMTX,lV,kV,nV;
      90             :   PetscErrorCode  ierr;
      91             :   Mat             Z;
      92             :   PetscBool       t;
      93             : #if !defined(PETSC_USE_COMPLEX)
      94             :   PetscInt        i;
      95             : #endif
      96             : 
      97         299 :   PetscFunctionBegin;
      98         299 :   npreconv = d->npreconv;
      99             :   /* Constrains the converged pairs to nev */
     100             : #if !defined(PETSC_USE_COMPLEX)
     101             :   /* Tries to maintain together conjugate eigenpairs */
     102         299 :   for (i=0; (i + (d->eigi[i]!=0.0?1:0) < npreconv) && (d->nconv + i < d->nev); i+= (d->eigi[i]!=0.0?2:1));
     103         299 :   npreconv = i;
     104             : #else
     105             :   npreconv = PetscMax(PetscMin(d->nev-d->nconv,npreconv),0);
     106             : #endif
     107             :   /* For GHEP without B-ortho, converge all of the requested pairs at once */
     108         299 :   ierr = PetscObjectTypeCompare((PetscObject)d->eps->ds,DSGHEP,&t);CHKERRQ(ierr);
     109         299 :   if (t && d->nconv+npreconv<d->nev) npreconv = 0;
     110             :   /* Quick exit */
     111         299 :   if (npreconv == 0) PetscFunctionReturn(0);
     112             : 
     113         299 :   ierr = BVGetActiveColumns(d->eps->V,&lV,&kV);CHKERRQ(ierr);
     114         299 :   nV  = kV - lV;
     115         299 :   cMT = nV - npreconv;
     116             :   /* Harmonics restarts wiht right eigenvectors, and other with the left ones.
     117             :      If the problem is standard or hermitian, left and right vectors are the same */
     118         299 :   if (!(d->W||DVD_IS(d->sEP,DVD_EP_STD)||DVD_IS(d->sEP,DVD_EP_HERMITIAN))) {
     119             :     /* ps.Q <- [ps.Q(0:npreconv-1) ps.Z(npreconv:size_H-1)] */
     120           0 :     ierr = DSGetMat(d->eps->ds,DS_MAT_Z,&Z);CHKERRQ(ierr);
     121           0 :     ierr = DSCopyMat(d->eps->ds,DS_MAT_Q,0,npreconv,Z,0,npreconv,nV,cMT,PETSC_FALSE);CHKERRQ(ierr);
     122           0 :     ierr = MatDestroy(&Z);CHKERRQ(ierr);
     123             :   }
     124         299 :   if (DVD_IS(d->sEP,DVD_EP_INDEFINITE)) {
     125           0 :     ierr = DSPseudoOrthogonalize(d->eps->ds,DS_MAT_Q,nV,d->nBds,&cMTX,d->nBds);CHKERRQ(ierr);
     126             :   } else {
     127         299 :     ierr = DSOrthogonalize(d->eps->ds,DS_MAT_Q,nV,&cMTX);CHKERRQ(ierr);
     128             :   }
     129         299 :   cMT = cMTX - npreconv;
     130             : 
     131         299 :   if (d->W) {
     132          57 :     ierr = DSOrthogonalize(d->eps->ds,DS_MAT_Z,nV,&cMTX);CHKERRQ(ierr);
     133          57 :     cMT = PetscMin(cMT,cMTX - npreconv);
     134             :   }
     135             : 
     136             :   /* Lock the converged pairs */
     137         299 :   d->eigr+= npreconv;
     138             : #if !defined(PETSC_USE_COMPLEX)
     139         299 :   if (d->eigi) d->eigi+= npreconv;
     140             : #endif
     141         299 :   d->nconv+= npreconv;
     142         299 :   d->errest+= npreconv;
     143             :   /* Notify the changes in V and update the other subspaces */
     144         299 :   d->V_tra_s = npreconv;          d->V_tra_e = nV;
     145         299 :   d->V_new_s = cMT;               d->V_new_e = d->V_new_s;
     146             :   /* Remove oldU */
     147         299 :   data->size_oldU = 0;
     148             : 
     149         299 :   d->npreconv-= npreconv;
     150         299 :   PetscFunctionReturn(0);
     151             : }
     152             : 
     153        1333 : static PetscErrorCode dvd_updateV_restart_gen(dvdDashboard *d)
     154             : {
     155        1333 :   dvdManagV_basic *data = (dvdManagV_basic*)d->updateV_data;
     156             :   PetscInt        lV,kV,nV,size_plusk,size_X,cMTX,cMTY,max_restart_size;
     157             :   Mat             Z;
     158             :   PetscErrorCode  ierr;
     159             : 
     160        1333 :   PetscFunctionBegin;
     161             :   /* Select size_X desired pairs from V */
     162             :   /* The restarted basis should:
     163             :      - have at least one spot to add a new direction;
     164             :      - keep converged vectors, npreconv;
     165             :      - keep at least 1 oldU direction if possible.
     166             :   */
     167        1333 :   ierr = BVGetActiveColumns(d->eps->V,&lV,&kV);CHKERRQ(ierr);
     168        1333 :   nV = kV - lV;
     169        1333 :   max_restart_size = PetscMax(0,PetscMin(d->eps->mpd - 1,d->eps->ncv - lV - 2));
     170        1333 :   size_X = PetscMin(PetscMin(data->min_size_V+d->npreconv,max_restart_size - (max_restart_size - d->npreconv > 1 && data->plusk > 0 && data->size_oldU > 0 ? 1 : 0)), nV);
     171             : 
     172             :   /* Add plusk eigenvectors from the previous iteration */
     173        1333 :   size_plusk = PetscMax(0,PetscMin(PetscMin(PetscMin(data->plusk,data->size_oldU),max_restart_size - size_X),nV - size_X));
     174             : 
     175        1333 :   d->size_MT = nV;
     176             :   /* ps.Q <- orth([pX(0:size_X-1) [oldU(0:size_plusk-1); 0] ]) */
     177             :   /* Harmonics restarts wiht right eigenvectors, and other with the left ones.
     178             :      If the problem is standard or hermitian, left and right vectors are the same */
     179        1333 :   if (!(d->W||DVD_IS(d->sEP,DVD_EP_STD)||DVD_IS(d->sEP,DVD_EP_HERMITIAN))) {
     180           0 :     ierr = DSGetMat(d->eps->ds,DS_MAT_Z,&Z);CHKERRQ(ierr);
     181           0 :     ierr = DSCopyMat(d->eps->ds,DS_MAT_Q,0,0,Z,0,0,nV,size_X,PETSC_FALSE);CHKERRQ(ierr);
     182           0 :     ierr = MatDestroy(&Z);CHKERRQ(ierr);
     183             :   }
     184        1333 :   if (size_plusk > 0 && DVD_IS(d->sEP,DVD_EP_INDEFINITE)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Unsupported plusk>0 in indefinite eigenvalue problems");
     185        1333 :   if (size_plusk > 0) {
     186        1280 :     ierr = DSCopyMat(d->eps->ds,DS_MAT_Q,0,size_X,data->oldU,0,0,nV,size_plusk,PETSC_FALSE);CHKERRQ(ierr);
     187             :   }
     188        1333 :   if (DVD_IS(d->sEP,DVD_EP_INDEFINITE)) {
     189           0 :     ierr = DSPseudoOrthogonalize(d->eps->ds,DS_MAT_Q,size_X,d->nBds,&cMTX,d->nBds);CHKERRQ(ierr);
     190             :   } else {
     191        1333 :     ierr = DSOrthogonalize(d->eps->ds,DS_MAT_Q,size_X+size_plusk,&cMTX);CHKERRQ(ierr);
     192             :   }
     193             : 
     194        1333 :   if (d->W && size_plusk > 0) {
     195             :     /* ps.Z <- orth([ps.Z(0:size_X-1) [oldV(0:size_plusk-1); 0] ]) */
     196         123 :     ierr = DSCopyMat(d->eps->ds,DS_MAT_Z,0,size_X,data->oldV,0,0,nV,size_plusk,PETSC_FALSE);CHKERRQ(ierr);
     197         123 :     ierr = DSOrthogonalize(d->eps->ds,DS_MAT_Z,size_X+size_plusk,&cMTY);CHKERRQ(ierr);
     198         123 :     cMTX = PetscMin(cMTX, cMTY);
     199             :   }
     200        1333 :   if (cMTX > size_X+size_plusk) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Invalid number of columns to restart");
     201             : 
     202             :   /* Notify the changes in V and update the other subspaces */
     203        1333 :   d->V_tra_s = 0;                     d->V_tra_e = cMTX;
     204        1333 :   d->V_new_s = d->V_tra_e;            d->V_new_e = d->V_new_s;
     205             : 
     206             :   /* Remove oldU */
     207        1333 :   data->size_oldU = 0;
     208             : 
     209             :   /* Remove npreconv */
     210        1333 :   d->npreconv = 0;
     211        1333 :   PetscFunctionReturn(0);
     212             : }
     213             : 
     214       14781 : static PetscErrorCode dvd_updateV_testConv(dvdDashboard *d,PetscInt s,PetscInt pre,PetscInt e,PetscInt *nConv)
     215             : {
     216             :   PetscInt        i,j,b;
     217             :   PetscReal       norm;
     218             :   PetscErrorCode  ierr;
     219             :   PetscBool       conv, c;
     220       14781 :   dvdManagV_basic *data = (dvdManagV_basic*)d->updateV_data;
     221             : 
     222       14781 :   PetscFunctionBegin;
     223       14781 :   if (nConv) *nConv = s;
     224       25096 :   for (i=s,conv=PETSC_TRUE;(conv || data->allResiduals) && (i < e);i+=b) {
     225             : #if !defined(PETSC_USE_COMPLEX)
     226       10315 :     b = d->eigi[i]!=0.0?2:1;
     227             : #else
     228             :     b = 1;
     229             : #endif
     230       10315 :     if (i+b-1 >= pre) {
     231        2778 :       ierr = d->calcpairs_residual(d,i,i+b);CHKERRQ(ierr);
     232             :     }
     233             :     /* Test the Schur vector */
     234       10317 :     for (j=0,c=PETSC_TRUE;j<b && c;j++) {
     235       10317 :       norm = d->nR[i+j]/d->nX[i+j];
     236       10317 :       c = d->testConv(d,d->eigr[i+j],d->eigi[i+j],norm,&d->errest[i+j]);
     237             :     }
     238       10315 :     if (conv && c) { if (nConv) *nConv = i+b; }
     239             :     else conv = PETSC_FALSE;
     240             :   }
     241       14781 :   pre = PetscMax(pre,i);
     242             : 
     243             : #if !defined(PETSC_USE_COMPLEX)
     244             :   /* Enforce converged conjugate complex eigenpairs */
     245       14781 :   if (nConv) {
     246         299 :     for (j=0;j<*nConv;j++) if (d->eigi[j] != 0.0) j++;
     247        7540 :     if (j>*nConv) (*nConv)--;
     248             :   }
     249             : #endif
     250         221 :   for (i=pre;i<e;i++) d->errest[i] = d->nR[i] = 1.0;
     251       14781 :   PetscFunctionReturn(0);
     252             : }
     253             : 
     254        7540 : static PetscErrorCode dvd_updateV_update_gen(dvdDashboard *d)
     255             : {
     256        7540 :   dvdManagV_basic *data = (dvdManagV_basic*)d->updateV_data;
     257             :   PetscInt        size_D,s,lV,kV,nV;
     258             :   PetscErrorCode  ierr;
     259             : 
     260        7540 :   PetscFunctionBegin;
     261             :   /* Select the desired pairs */
     262        7540 :   ierr = BVGetActiveColumns(d->eps->V,&lV,&kV);CHKERRQ(ierr);
     263        7540 :   nV = kV - lV;
     264        7540 :   size_D = PetscMin(PetscMin(PetscMin(d->bs,nV),d->eps->ncv-nV),nV);
     265        7540 :   if (size_D == 0) PetscFunctionReturn(0);
     266             : 
     267             :   /* Fill V with D */
     268        7540 :   ierr = d->improveX(d,d->npreconv,d->npreconv+size_D,&size_D);CHKERRQ(ierr);
     269             : 
     270             :   /* If D is empty, exit */
     271        7540 :   d->size_D = size_D;
     272        7540 :   if (size_D == 0) PetscFunctionReturn(0);
     273             : 
     274             :   /* Get the residual of all pairs */
     275             : #if !defined(PETSC_USE_COMPLEX)
     276        7241 :   s = (d->eigi[0]!=0.0)? 2: 1;
     277             : #else
     278             :   s = 1;
     279             : #endif
     280        7241 :   ierr = BVGetActiveColumns(d->eps->V,&lV,&kV);CHKERRQ(ierr);
     281        7241 :   nV = kV - lV;
     282        7241 :   ierr = dvd_updateV_testConv(d,s,s,data->allResiduals?nV:size_D,NULL);CHKERRQ(ierr);
     283             : 
     284             :   /* Notify the changes in V */
     285        7241 :   d->V_tra_s = 0;                 d->V_tra_e = 0;
     286        7241 :   d->V_new_s = nV;                d->V_new_e = nV+size_D;
     287             : 
     288             :   /* Save the projected eigenvectors */
     289        7241 :   if (data->plusk > 0) {
     290        7079 :     ierr = MatZeroEntries(data->oldU);CHKERRQ(ierr);
     291        7079 :     data->size_oldU = nV;
     292        7079 :     ierr = DSCopyMat(d->eps->ds,DS_MAT_Q,0,0,data->oldU,0,0,nV,nV,PETSC_TRUE);CHKERRQ(ierr);
     293        7079 :     if (d->W) {
     294         890 :       ierr = MatZeroEntries(data->oldV);CHKERRQ(ierr);
     295         890 :       ierr = DSCopyMat(d->eps->ds,DS_MAT_Z,0,0,data->oldV,0,0,nV,nV,PETSC_TRUE);CHKERRQ(ierr);
     296             :     }
     297             :   }
     298        7241 :   PetscFunctionReturn(0);
     299             : }
     300             : 
     301        8873 : static PetscErrorCode dvd_updateV_extrapol(dvdDashboard *d)
     302             : {
     303        8873 :   dvdManagV_basic *data = (dvdManagV_basic*)d->updateV_data;
     304             :   PetscInt        i;
     305             :   PetscBool       restart,t;
     306             :   PetscErrorCode  ierr;
     307             : 
     308        8873 :   PetscFunctionBegin;
     309             :   /* TODO: restrict select pairs to each case */
     310        8873 :   ierr = d->calcpairs_selectPairs(d, data->min_size_V+d->npreconv);CHKERRQ(ierr);
     311             : 
     312             :   /* If the subspaces doesn't need restart, add new vector */
     313        8873 :   ierr = d->isRestarting(d,&restart);CHKERRQ(ierr);
     314        8873 :   if (!restart) {
     315        7540 :     d->size_D = 0;
     316        7540 :     ierr = dvd_updateV_update_gen(d);CHKERRQ(ierr);
     317             : 
     318             :     /* If no vector were converged, exit */
     319             :     /* For GHEP without B-ortho, converge all of the requested pairs at once */
     320        7540 :     ierr = PetscObjectTypeCompareAny((PetscObject)d->eps->ds,&t,DSGHEP,"");CHKERRQ(ierr);
     321        7540 :     if (d->nconv+d->npreconv < d->nev && (t || d->npreconv == 0)) PetscFunctionReturn(0);
     322             :   }
     323             : 
     324             :   /* If some eigenpairs were converged, lock them  */
     325        1632 :   if (d->npreconv > 0) {
     326         299 :     i = d->npreconv;
     327         299 :     ierr = dvd_updateV_conv_gen(d);CHKERRQ(ierr);
     328             : 
     329             :     /* If some eigenpair was locked, exit */
     330         299 :     if (i > d->npreconv) PetscFunctionReturn(0);
     331             :   }
     332             : 
     333             :   /* Else, a restarting is performed */
     334        1333 :   ierr = dvd_updateV_restart_gen(d);CHKERRQ(ierr);
     335        1333 :   PetscFunctionReturn(0);
     336             : }
     337             : 
     338         270 : PetscErrorCode dvd_managementV_basic(dvdDashboard *d,dvdBlackboard *b,PetscInt bs,PetscInt mpd,PetscInt min_size_V,PetscInt plusk,PetscBool harm,PetscBool allResiduals)
     339             : {
     340             :   PetscErrorCode  ierr;
     341             :   dvdManagV_basic *data;
     342             : #if !defined(PETSC_USE_COMPLEX)
     343             :   PetscBool       her_probl,std_probl;
     344             : #endif
     345             : 
     346         270 :   PetscFunctionBegin;
     347             :   /* Setting configuration constrains */
     348             : #if !defined(PETSC_USE_COMPLEX)
     349             :   /* if the last converged eigenvalue is complex its conjugate pair is also
     350             :      converged */
     351         270 :   her_probl = DVD_IS(d->sEP,DVD_EP_HERMITIAN)? PETSC_TRUE: PETSC_FALSE;
     352         270 :   std_probl = DVD_IS(d->sEP,DVD_EP_STD)? PETSC_TRUE: PETSC_FALSE;
     353         270 :   b->max_size_X = PetscMax(b->max_size_X,bs+((her_probl && std_probl)?0:1));
     354             : #else
     355             :   b->max_size_X = PetscMax(b->max_size_X,bs);
     356             : #endif
     357             : 
     358         270 :   b->max_size_V = PetscMax(b->max_size_V,mpd);
     359         270 :   min_size_V = PetscMin(min_size_V,mpd-bs);
     360         270 :   b->size_V = PetscMax(b->size_V,b->max_size_V+b->max_size_P+b->max_nev);
     361         270 :   b->max_size_oldX = plusk;
     362             : 
     363             :   /* Setup the step */
     364         270 :   if (b->state >= DVD_STATE_CONF) {
     365          90 :     ierr = PetscNewLog(d->eps,&data);CHKERRQ(ierr);
     366          90 :     data->mpd = b->max_size_V;
     367          90 :     data->min_size_V = min_size_V;
     368          90 :     d->bs = bs;
     369          90 :     data->plusk = plusk;
     370          90 :     data->allResiduals = allResiduals;
     371             : 
     372          90 :     d->eigr = d->eps->eigr;
     373          90 :     d->eigi = d->eps->eigi;
     374          90 :     d->errest = d->eps->errest;
     375          90 :     ierr = PetscMalloc1(d->eps->ncv,&d->real_nR);CHKERRQ(ierr);
     376          90 :     ierr = PetscMalloc1(d->eps->ncv,&d->real_nX);CHKERRQ(ierr);
     377          90 :     if (plusk > 0) { ierr = MatCreateSeqDense(PETSC_COMM_SELF,d->eps->ncv,d->eps->ncv,NULL,&data->oldU);CHKERRQ(ierr); }
     378           2 :     else data->oldU = NULL;
     379          90 :     if (harm && plusk>0) { ierr = MatCreateSeqDense(PETSC_COMM_SELF,d->eps->ncv,d->eps->ncv,NULL,&data->oldV);CHKERRQ(ierr); }
     380          71 :     else data->oldV = NULL;
     381             : 
     382          90 :     data->old_updateV_data = d->updateV_data;
     383          90 :     d->updateV_data = data;
     384          90 :     data->old_isRestarting = d->isRestarting;
     385          90 :     d->isRestarting = dvd_isrestarting_fullV;
     386          90 :     d->updateV = dvd_updateV_extrapol;
     387          90 :     d->preTestConv = dvd_updateV_testConv;
     388          90 :     ierr = EPSDavidsonFLAdd(&d->startList,dvd_updateV_start);CHKERRQ(ierr);
     389          90 :     ierr = EPSDavidsonFLAdd(&d->destroyList,dvd_managementV_basic_d);CHKERRQ(ierr);
     390             :   }
     391         270 :   PetscFunctionReturn(0);
     392             : }
     393             : 

Generated by: LCOV version 1.13