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 : This file implements a wrapper to eigensolvers in EVSL.
12 : */
13 :
14 : #include <slepc/private/epsimpl.h> /*I "slepceps.h" I*/
15 : #include <evsl.h>
16 :
17 : #define PetscCallEVSL(func, ...) do { \
18 : PetscStackPushExternal(PetscStringize(func)); \
19 : int evsl_ierr_ = func(__VA_ARGS__); \
20 : PetscStackPop; \
21 : PetscCheck(!evsl_ierr_,PETSC_COMM_SELF,PETSC_ERR_LIB,"Error calling %s: error code %d",PetscStringize(func(__VA_ARGS__)),evsl_ierr_); \
22 : } while (0)
23 :
24 : typedef struct {
25 : PetscBool initialized;
26 : Mat A; /* problem matrix */
27 : Vec x,y; /* auxiliary vectors */
28 : PetscReal *sli; /* slice bounds */
29 : PetscInt nev; /* approximate number of wanted eigenvalues in each slice */
30 : PetscLayout map; /* used to distribute slices among MPI processes */
31 : PetscBool estimrange; /* the filter range was not set by the user */
32 : /* user parameters */
33 : PetscInt nslices; /* number of slices */
34 : PetscReal lmin,lmax; /* numerical range (min and max eigenvalue) */
35 : EPSEVSLDOSMethod dos; /* DOS method, either KPM or Lanczos */
36 : PetscInt nvec; /* number of sample vectors used for DOS */
37 : PetscInt deg; /* polynomial degree used for DOS (KPM only) */
38 : PetscInt steps; /* number of Lanczos steps used for DOS (Lanczos only) */
39 : PetscInt npoints; /* number of sample points used for DOS (Lanczos only) */
40 : PetscInt max_deg; /* maximum degree allowed for the polynomial */
41 : PetscReal thresh; /* threshold for accepting polynomial */
42 : EPSEVSLDamping damping; /* type of damping (for polynomial and for DOS-KPM) */
43 : } EPS_EVSL;
44 :
45 497050 : static void AMatvec_EVSL(double *xa,double *ya,void *data)
46 : {
47 497050 : EPS_EVSL *ctx = (EPS_EVSL*)data;
48 497050 : Vec x = ctx->x,y = ctx->y;
49 497050 : Mat A = ctx->A;
50 :
51 497050 : PetscFunctionBegin;
52 497050 : PetscCallAbort(PetscObjectComm((PetscObject)A),VecPlaceArray(x,(PetscScalar*)xa));
53 497050 : PetscCallAbort(PetscObjectComm((PetscObject)A),VecPlaceArray(y,(PetscScalar*)ya));
54 497050 : PetscCallAbort(PetscObjectComm((PetscObject)A),MatMult(A,x,y));
55 497050 : PetscCallAbort(PetscObjectComm((PetscObject)A),VecResetArray(x));
56 497050 : PetscCallAbort(PetscObjectComm((PetscObject)A),VecResetArray(y));
57 497050 : PetscFunctionReturnVoid();
58 : }
59 :
60 15 : static PetscErrorCode EPSSetUp_EVSL(EPS eps)
61 : {
62 15 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
63 15 : PetscMPIInt size,rank;
64 15 : PetscBool isshift;
65 15 : PetscScalar *vinit;
66 15 : PetscReal *mu,ecount,xintv[4],*xdos,*ydos;
67 15 : Vec v0;
68 15 : Mat A;
69 15 : PetscRandom rnd;
70 :
71 15 : PetscFunctionBegin;
72 15 : EPSCheckStandard(eps);
73 15 : EPSCheckHermitian(eps);
74 15 : PetscCall(PetscObjectTypeCompare((PetscObject)eps->st,STSHIFT,&isshift));
75 15 : PetscCheck(isshift,PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"This solver does not support spectral transformations");
76 :
77 15 : if (ctx->initialized) EVSLFinish();
78 15 : EVSLStart();
79 15 : ctx->initialized=PETSC_TRUE;
80 :
81 : /* get number of slices per process */
82 15 : PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)eps),&size));
83 15 : PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank));
84 15 : if (!ctx->nslices) ctx->nslices = size;
85 15 : PetscCall(PetscLayoutDestroy(&ctx->map));
86 15 : PetscCall(PetscLayoutCreateFromSizes(PetscObjectComm((PetscObject)eps),PETSC_DECIDE,ctx->nslices,1,&ctx->map));
87 :
88 : /* get matrix and prepare auxiliary vectors */
89 15 : PetscCall(MatDestroy(&ctx->A));
90 15 : PetscCall(STGetMatrix(eps->st,0,&A));
91 15 : if (size==1) {
92 7 : PetscCall(PetscObjectReference((PetscObject)A));
93 7 : ctx->A = A;
94 8 : } else PetscCall(MatCreateRedundantMatrix(A,0,PETSC_COMM_SELF,MAT_INITIAL_MATRIX,&ctx->A));
95 15 : SetAMatvec(eps->n,&AMatvec_EVSL,(void*)ctx);
96 15 : if (!ctx->x) PetscCall(MatCreateVecsEmpty(ctx->A,&ctx->x,&ctx->y));
97 15 : EPSCheckUnsupported(eps,EPS_FEATURE_ARBITRARY | EPS_FEATURE_REGION | EPS_FEATURE_STOPPING);
98 15 : EPSCheckIgnored(eps,EPS_FEATURE_EXTRACTION | EPS_FEATURE_CONVERGENCE);
99 :
100 15 : if (!eps->which) eps->which=EPS_ALL;
101 15 : PetscCheck(eps->which==EPS_ALL && eps->inta!=eps->intb,PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"This solver requires setting an interval with EPSSetInterval()");
102 :
103 : /* estimate numerical range */
104 15 : if (ctx->estimrange || ctx->lmin == PETSC_MIN_REAL || ctx->lmax == PETSC_MAX_REAL) {
105 14 : PetscCall(MatCreateVecs(ctx->A,&v0,NULL));
106 14 : if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
107 14 : PetscCall(BVGetRandomContext(eps->V,&rnd));
108 14 : PetscCall(VecSetRandom(v0,rnd));
109 14 : PetscCall(VecGetArray(v0,&vinit));
110 14 : PetscCallEVSL(LanTrbounds,50,200,eps->tol,vinit,1,&ctx->lmin,&ctx->lmax,NULL);
111 14 : PetscCall(VecRestoreArray(v0,&vinit));
112 14 : PetscCall(VecDestroy(&v0));
113 14 : ctx->estimrange = PETSC_TRUE; /* estimate if called again with another matrix */
114 : }
115 15 : PetscCheck(ctx->lmin<=eps->inta && ctx->lmax>=eps->intb,PetscObjectComm((PetscObject)eps),PETSC_ERR_SUP,"The requested interval [%g,%g] must be contained in the numerical range [%g,%g]",(double)eps->inta,(double)eps->intb,(double)ctx->lmin,(double)ctx->lmax);
116 15 : xintv[0] = eps->inta;
117 15 : xintv[1] = eps->intb;
118 15 : xintv[2] = ctx->lmin;
119 15 : xintv[3] = ctx->lmax;
120 :
121 : /* estimate number of eigenvalues in the interval */
122 15 : switch (ctx->dos) {
123 9 : case EPS_EVSL_DOS_KPM:
124 9 : PetscCall(PetscMalloc1(ctx->deg+1,&mu));
125 9 : if (!rank) PetscCallEVSL(kpmdos,ctx->deg,(int)ctx->damping,ctx->nvec,xintv,mu,&ecount);
126 18 : PetscCallMPI(MPI_Bcast(mu,ctx->deg+1,MPIU_REAL,0,PetscObjectComm((PetscObject)eps)));
127 : break;
128 6 : case EPS_EVSL_DOS_LANCZOS:
129 6 : PetscCall(PetscMalloc2(ctx->npoints,&xdos,ctx->npoints,&ydos));
130 6 : if (!rank) PetscCallEVSL(LanDos,ctx->nvec,PetscMin(ctx->steps,eps->n/2),ctx->npoints,xdos,ydos,&ecount,xintv);
131 12 : PetscCallMPI(MPI_Bcast(xdos,ctx->npoints,MPIU_REAL,0,PetscObjectComm((PetscObject)eps)));
132 12 : PetscCallMPI(MPI_Bcast(ydos,ctx->npoints,MPIU_REAL,0,PetscObjectComm((PetscObject)eps)));
133 : break;
134 0 : default:
135 0 : SETERRQ(PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Invalid DOS method");
136 : }
137 30 : PetscCallMPI(MPI_Bcast(&ecount,1,MPIU_REAL,0,PetscObjectComm((PetscObject)eps)));
138 :
139 15 : PetscCall(PetscInfo(eps,"Estimated eigenvalue count in the interval: %g\n",ecount));
140 15 : eps->ncv = (PetscInt)PetscCeilReal(1.5*ecount);
141 :
142 : /* slice the spectrum */
143 15 : PetscCall(PetscFree(ctx->sli));
144 15 : PetscCall(PetscMalloc1(ctx->nslices+1,&ctx->sli));
145 15 : if (ctx->dos == EPS_EVSL_DOS_KPM) {
146 9 : PetscCallEVSL(spslicer,ctx->sli,mu,ctx->deg,xintv,ctx->nslices,10*(PetscInt)ecount);
147 9 : PetscCall(PetscFree(mu));
148 6 : } else if (ctx->dos == EPS_EVSL_DOS_LANCZOS) {
149 6 : spslicer2(xdos,ydos,ctx->nslices,ctx->npoints,ctx->sli);
150 6 : PetscCall(PetscFree2(xdos,ydos));
151 : }
152 :
153 : /* approximate number of eigenvalues wanted in each slice */
154 15 : ctx->nev = (PetscInt)(1.0 + ecount/(PetscReal)ctx->nslices) + 2;
155 :
156 15 : if (eps->mpd!=PETSC_DEFAULT) PetscCall(PetscInfo(eps,"Warning: parameter mpd ignored\n"));
157 15 : if (eps->max_it==PETSC_DEFAULT) eps->max_it = 1;
158 15 : PetscCall(EPSAllocateSolution(eps,0));
159 15 : PetscFunctionReturn(PETSC_SUCCESS);
160 : }
161 :
162 15 : static PetscErrorCode EPSSolve_EVSL(EPS eps)
163 : {
164 15 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
165 15 : PetscInt i,j,k=0,sl,mlan,nevout,*ind,nevmax,rstart,rend,*nevloc,*disp,N;
166 15 : PetscReal *res,xintv[4],*errest;
167 15 : PetscScalar *lam,*X,*Y,*vinit,*eigr;
168 15 : PetscMPIInt size,rank;
169 15 : PetscRandom rnd;
170 15 : Vec v,w,v0,x;
171 15 : VecScatter vs;
172 15 : IS is;
173 15 : polparams pol;
174 :
175 15 : PetscFunctionBegin;
176 15 : PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)eps),&size));
177 15 : PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank));
178 15 : PetscCall(PetscLayoutGetRange(ctx->map,&rstart,&rend));
179 15 : nevmax = (rend-rstart)*ctx->nev;
180 15 : PetscCall(MatCreateVecs(ctx->A,&v0,NULL));
181 15 : PetscCall(BVGetRandomContext(eps->V,&rnd));
182 15 : PetscCall(VecSetRandom(v0,rnd));
183 15 : PetscCall(VecGetArray(v0,&vinit));
184 15 : PetscCall(PetscMalloc5(size,&nevloc,size+1,&disp,nevmax,&eigr,nevmax,&errest,nevmax*eps->n,&X));
185 15 : mlan = PetscMin(PetscMax(5*ctx->nev,300),eps->n);
186 37 : for (sl=rstart; sl<rend; sl++) {
187 22 : xintv[0] = ctx->sli[sl];
188 22 : xintv[1] = ctx->sli[sl+1];
189 22 : xintv[2] = ctx->lmin;
190 22 : xintv[3] = ctx->lmax;
191 22 : PetscCall(PetscInfo(ctx->A,"Subinterval %" PetscInt_FMT ": [%.4e, %.4e]\n",sl+1,xintv[0],xintv[1]));
192 22 : set_pol_def(&pol);
193 22 : pol.max_deg = ctx->max_deg;
194 22 : pol.damping = (int)ctx->damping;
195 22 : pol.thresh_int = ctx->thresh;
196 22 : find_pol(xintv,&pol);
197 22 : PetscCall(PetscInfo(ctx->A,"Polynomial [type = %" PetscInt_FMT "], deg %" PetscInt_FMT ", bar %e gam %e\n",pol.type,pol.deg,pol.bar,pol.gam));
198 22 : PetscCallEVSL(ChebLanNr,xintv,mlan,eps->tol,vinit,&pol,&nevout,&lam,&Y,&res,NULL);
199 22 : PetscCheck(k+nevout<=nevmax,PetscObjectComm((PetscObject)eps),PETSC_ERR_LIB,"Too low estimation of eigenvalue count, try modifying the sampling parameters");
200 22 : free_pol(&pol);
201 22 : PetscCall(PetscInfo(ctx->A,"Computed %" PetscInt_FMT " eigenvalues\n",nevout));
202 22 : PetscCall(PetscMalloc1(nevout,&ind));
203 22 : sort_double(nevout,lam,ind);
204 170 : for (i=0;i<nevout;i++) {
205 126 : eigr[i+k] = lam[i];
206 126 : errest[i+k] = res[ind[i]];
207 126 : PetscCall(PetscArraycpy(X+(i+k)*eps->n,Y+ind[i]*eps->n,eps->n));
208 : }
209 22 : k += nevout;
210 22 : if (lam) evsl_Free(lam);
211 22 : if (Y) evsl_Free_device(Y);
212 22 : if (res) evsl_Free(res);
213 22 : PetscCall(PetscFree(ind));
214 : }
215 15 : PetscCall(VecRestoreArray(v0,&vinit));
216 15 : PetscCall(VecDestroy(&v0));
217 :
218 : /* gather eigenvalues computed by each MPI process */
219 30 : PetscCallMPI(MPI_Allgather(&k,1,MPIU_INT,nevloc,1,MPIU_INT,PetscObjectComm((PetscObject)eps)));
220 15 : eps->nev = nevloc[0];
221 15 : disp[0] = 0;
222 23 : for (i=1;i<size;i++) {
223 8 : eps->nev += nevloc[i];
224 8 : disp[i] = disp[i-1]+nevloc[i-1];
225 : }
226 15 : disp[size] = disp[size-1]+nevloc[size-1];
227 15 : PetscCheck(eps->nev<=eps->ncv,PetscObjectComm((PetscObject)eps),PETSC_ERR_LIB,"Too low estimation of eigenvalue count, try modifying the sampling parameters");
228 30 : PetscCallMPI(MPI_Allgatherv(eigr,k,MPIU_SCALAR,eps->eigr,nevloc,disp,MPIU_SCALAR,PetscObjectComm((PetscObject)eps)));
229 30 : PetscCallMPI(MPI_Allgatherv(errest,k,MPIU_REAL,eps->errest,nevloc,disp,MPIU_REAL,PetscObjectComm((PetscObject)eps)));
230 15 : eps->nconv = eps->nev;
231 15 : eps->its = 1;
232 15 : eps->reason = EPS_CONVERGED_TOL;
233 :
234 : /* scatter computed eigenvectors and store them in eps->V */
235 15 : PetscCall(BVCreateVec(eps->V,&w));
236 38 : for (i=0;i<size;i++) {
237 23 : N = (rank==i)? eps->n: 0;
238 23 : PetscCall(VecCreateSeq(PETSC_COMM_SELF,N,&x));
239 23 : PetscCall(VecSetFromOptions(x));
240 23 : PetscCall(ISCreateStride(PETSC_COMM_SELF,N,0,1,&is));
241 23 : PetscCall(VecScatterCreate(x,is,w,is,&vs));
242 23 : PetscCall(ISDestroy(&is));
243 176 : for (j=disp[i];j<disp[i+1];j++) {
244 153 : PetscCall(BVGetColumn(eps->V,j,&v));
245 153 : if (rank==i) PetscCall(VecPlaceArray(x,X+(j-disp[i])*eps->n));
246 153 : PetscCall(VecScatterBegin(vs,x,v,INSERT_VALUES,SCATTER_FORWARD));
247 153 : PetscCall(VecScatterEnd(vs,x,v,INSERT_VALUES,SCATTER_FORWARD));
248 153 : if (rank==i) PetscCall(VecResetArray(x));
249 153 : PetscCall(BVRestoreColumn(eps->V,j,&v));
250 : }
251 23 : PetscCall(VecScatterDestroy(&vs));
252 23 : PetscCall(VecDestroy(&x));
253 : }
254 15 : PetscCall(VecDestroy(&w));
255 15 : PetscCall(PetscFree5(nevloc,disp,eigr,errest,X));
256 15 : PetscFunctionReturn(PETSC_SUCCESS);
257 : }
258 :
259 2 : static PetscErrorCode EPSEVSLSetSlices_EVSL(EPS eps,PetscInt nslices)
260 : {
261 2 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
262 :
263 2 : PetscFunctionBegin;
264 2 : if (nslices == PETSC_DECIDE || nslices == PETSC_DEFAULT) nslices = 0;
265 2 : else PetscCheck(nslices>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"Number of slices must be 1 at least");
266 2 : if (ctx->nslices != nslices) {
267 2 : ctx->nslices = nslices;
268 2 : eps->state = EPS_STATE_INITIAL;
269 : }
270 2 : PetscFunctionReturn(PETSC_SUCCESS);
271 : }
272 :
273 : /*@
274 : EPSEVSLSetSlices - Set the number of slices in which the interval must be
275 : subdivided.
276 :
277 : Logically Collective
278 :
279 : Input Parameters:
280 : + eps - the eigensolver context
281 : - nslices - the number of slices
282 :
283 : Options Database Key:
284 : . -eps_evsl_slices <n> - set the number of slices to n
285 :
286 : Notes:
287 : By default, one slice per MPI process is used. Depending on the number of
288 : eigenvalues, using more slices may be beneficial, but very narrow subintervals
289 : imply higher polynomial degree.
290 :
291 : Level: intermediate
292 :
293 : .seealso: EPSEVSLGetSlices()
294 : @*/
295 2 : PetscErrorCode EPSEVSLSetSlices(EPS eps,PetscInt nslices)
296 : {
297 2 : PetscFunctionBegin;
298 2 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
299 8 : PetscValidLogicalCollectiveInt(eps,nslices,2);
300 2 : PetscTryMethod(eps,"EPSEVSLSetSlices_C",(EPS,PetscInt),(eps,nslices));
301 2 : PetscFunctionReturn(PETSC_SUCCESS);
302 : }
303 :
304 1 : static PetscErrorCode EPSEVSLGetSlices_EVSL(EPS eps,PetscInt *nslices)
305 : {
306 1 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
307 :
308 1 : PetscFunctionBegin;
309 1 : *nslices = ctx->nslices;
310 1 : PetscFunctionReturn(PETSC_SUCCESS);
311 : }
312 :
313 : /*@
314 : EPSEVSLGetSlices - Gets the number of slices in which the interval must be
315 : subdivided.
316 :
317 : Not Collective
318 :
319 : Input Parameter:
320 : . eps - the eigensolver context
321 :
322 : Output Parameter:
323 : . nslices - the number of slices
324 :
325 : Level: intermediate
326 :
327 : .seealso: EPSEVSLSetSlices()
328 : @*/
329 1 : PetscErrorCode EPSEVSLGetSlices(EPS eps,PetscInt *nslices)
330 : {
331 1 : PetscFunctionBegin;
332 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
333 1 : PetscAssertPointer(nslices,2);
334 1 : PetscUseMethod(eps,"EPSEVSLGetSlices_C",(EPS,PetscInt*),(eps,nslices));
335 1 : PetscFunctionReturn(PETSC_SUCCESS);
336 : }
337 :
338 1 : static PetscErrorCode EPSEVSLSetRange_EVSL(EPS eps,PetscReal lmin,PetscReal lmax)
339 : {
340 1 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
341 :
342 1 : PetscFunctionBegin;
343 1 : PetscCheck(lmin<lmax,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_WRONG,"Badly defined interval, must be lmin<lmax");
344 1 : if (ctx->lmin != lmin || ctx->lmax != lmax) {
345 1 : ctx->lmin = lmin;
346 1 : ctx->lmax = lmax;
347 1 : eps->state = EPS_STATE_INITIAL;
348 : }
349 1 : PetscFunctionReturn(PETSC_SUCCESS);
350 : }
351 :
352 : /*@
353 : EPSEVSLSetRange - Defines the numerical range (or field of values) of the problem,
354 : that is, the interval containing all eigenvalues.
355 :
356 : Logically Collective
357 :
358 : Input Parameters:
359 : + eps - the eigensolver context
360 : . lmin - left end of the interval
361 : - lmax - right end of the interval
362 :
363 : Options Database Key:
364 : . -eps_evsl_range <a,b> - set [a,b] as the numerical range
365 :
366 : Notes:
367 : The filter will be most effective if the numerical range is tight, that is, lmin
368 : and lmax are good approximations to the leftmost and rightmost eigenvalues,
369 : respectively. If not set by the user, an approximation is computed internally.
370 :
371 : The wanted computational interval specified via EPSSetInterval() must be
372 : contained in the numerical range.
373 :
374 : Level: intermediate
375 :
376 : .seealso: EPSEVSLGetRange(), EPSSetInterval()
377 : @*/
378 1 : PetscErrorCode EPSEVSLSetRange(EPS eps,PetscReal lmin,PetscReal lmax)
379 : {
380 1 : PetscFunctionBegin;
381 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
382 4 : PetscValidLogicalCollectiveReal(eps,lmin,2);
383 4 : PetscValidLogicalCollectiveReal(eps,lmax,3);
384 1 : PetscTryMethod(eps,"EPSEVSLSetRange_C",(EPS,PetscReal,PetscReal),(eps,lmin,lmax));
385 1 : PetscFunctionReturn(PETSC_SUCCESS);
386 : }
387 :
388 1 : static PetscErrorCode EPSEVSLGetRange_EVSL(EPS eps,PetscReal *lmin,PetscReal *lmax)
389 : {
390 1 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
391 :
392 1 : PetscFunctionBegin;
393 1 : if (lmin) *lmin = ctx->lmin;
394 1 : if (lmax) *lmax = ctx->lmax;
395 1 : PetscFunctionReturn(PETSC_SUCCESS);
396 : }
397 :
398 : /*@
399 : EPSEVSLGetRange - Gets the interval containing all eigenvalues.
400 :
401 : Not Collective
402 :
403 : Input Parameter:
404 : . eps - the eigensolver context
405 :
406 : Output Parameters:
407 : + lmin - left end of the interval
408 : - lmax - right end of the interval
409 :
410 : Level: intermediate
411 :
412 : .seealso: EPSEVSLSetRange()
413 : @*/
414 1 : PetscErrorCode EPSEVSLGetRange(EPS eps,PetscReal *lmin,PetscReal *lmax)
415 : {
416 1 : PetscFunctionBegin;
417 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
418 1 : PetscUseMethod(eps,"EPSEVSLGetRange_C",(EPS,PetscReal*,PetscReal*),(eps,lmin,lmax));
419 1 : PetscFunctionReturn(PETSC_SUCCESS);
420 : }
421 :
422 4 : static PetscErrorCode EPSEVSLSetDOSParameters_EVSL(EPS eps,EPSEVSLDOSMethod dos,PetscInt nvec,PetscInt deg,PetscInt steps,PetscInt npoints)
423 : {
424 4 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
425 :
426 4 : PetscFunctionBegin;
427 4 : ctx->dos = dos;
428 4 : if (nvec == PETSC_DECIDE || nvec == PETSC_DEFAULT) ctx->nvec = 80;
429 : else {
430 4 : PetscCheck(nvec>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"The nvec argument must be > 0");
431 4 : ctx->nvec = nvec;
432 : }
433 4 : switch (dos) {
434 1 : case EPS_EVSL_DOS_KPM:
435 1 : if (deg == PETSC_DECIDE || deg == PETSC_DEFAULT) ctx->deg = 300;
436 : else {
437 1 : PetscCheck(deg>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"The deg argument must be > 0");
438 1 : ctx->deg = deg;
439 : }
440 : break;
441 3 : case EPS_EVSL_DOS_LANCZOS:
442 3 : if (steps == PETSC_DECIDE || steps == PETSC_DEFAULT) ctx->steps = 40;
443 : else {
444 3 : PetscCheck(steps>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"The steps argument must be > 0");
445 3 : ctx->steps = steps;
446 : }
447 3 : if (npoints == PETSC_DECIDE || npoints == PETSC_DEFAULT) ctx->npoints = 200;
448 : else {
449 3 : PetscCheck(npoints>0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"The npoints argument must be > 0");
450 3 : ctx->npoints = npoints;
451 : }
452 : break;
453 : }
454 4 : eps->state = EPS_STATE_INITIAL;
455 4 : PetscFunctionReturn(PETSC_SUCCESS);
456 : }
457 :
458 : /*@
459 : EPSEVSLSetDOSParameters - Defines the parameters used for computing the
460 : density of states (DOS) in the EVSL solver.
461 :
462 : Logically Collective
463 :
464 : Input Parameters:
465 : + eps - the eigensolver context
466 : . dos - DOS method, either KPM or Lanczos
467 : . nvec - number of sample vectors
468 : . deg - polynomial degree (KPM only)
469 : . steps - number of Lanczos steps (Lanczos only)
470 : - npoints - number of sample points (Lanczos only)
471 :
472 : Options Database Keys:
473 : + -eps_evsl_dos_method <dos> - set the DOS method, either kpm or lanczos
474 : . -eps_evsl_dos_nvec <n> - set the number of sample vectors
475 : . -eps_evsl_dos_degree <n> - set the polynomial degree
476 : . -eps_evsl_dos_steps <n> - set the number of Lanczos steps
477 : - -eps_evsl_dos_npoints <n> - set the number of sample points
478 :
479 : Notes:
480 : The density of states (or spectral density) can be approximated with two
481 : methods, kernel polynomial method (KPM) or Lanczos. Some parameters for
482 : these methods can be set by the user with this function, with some of
483 : them being relevant for one of the methods only.
484 :
485 : Level: intermediate
486 :
487 : .seealso: EPSEVSLGetDOSParameters()
488 : @*/
489 4 : PetscErrorCode EPSEVSLSetDOSParameters(EPS eps,EPSEVSLDOSMethod dos,PetscInt nvec,PetscInt deg,PetscInt steps,PetscInt npoints)
490 : {
491 4 : PetscFunctionBegin;
492 4 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
493 16 : PetscValidLogicalCollectiveEnum(eps,dos,2);
494 16 : PetscValidLogicalCollectiveInt(eps,nvec,3);
495 16 : PetscValidLogicalCollectiveInt(eps,deg,4);
496 16 : PetscValidLogicalCollectiveInt(eps,steps,5);
497 16 : PetscValidLogicalCollectiveInt(eps,npoints,6);
498 4 : PetscTryMethod(eps,"EPSEVSLSetDOSParameters_C",(EPS,EPSEVSLDOSMethod,PetscInt,PetscInt,PetscInt,PetscInt),(eps,dos,nvec,deg,steps,npoints));
499 4 : PetscFunctionReturn(PETSC_SUCCESS);
500 : }
501 :
502 9 : static PetscErrorCode EPSEVSLGetDOSParameters_EVSL(EPS eps,EPSEVSLDOSMethod *dos,PetscInt *nvec,PetscInt *deg,PetscInt *steps,PetscInt *npoints)
503 : {
504 9 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
505 :
506 9 : PetscFunctionBegin;
507 9 : if (dos) *dos = ctx->dos;
508 9 : if (nvec) *nvec = ctx->nvec;
509 9 : if (deg) *deg = ctx->deg;
510 9 : if (steps) *steps = ctx->steps;
511 9 : if (npoints) *npoints = ctx->npoints;
512 9 : PetscFunctionReturn(PETSC_SUCCESS);
513 : }
514 :
515 : /*@
516 : EPSEVSLGetDOSParameters - Gets the parameters used for computing the
517 : density of states (DOS) in the EVSL solver.
518 :
519 : Not Collective
520 :
521 : Input Parameter:
522 : . eps - the eigensolver context
523 :
524 : Output Parameters:
525 : + dos - DOS method, either KPM or Lanczos
526 : . nvec - number of sample vectors
527 : . deg - polynomial degree (KPM only)
528 : . steps - number of Lanczos steps (Lanczos only)
529 : - npoints - number of sample points (Lanczos only)
530 :
531 : Level: intermediate
532 :
533 : .seealso: EPSEVSLSetDOSParameters()
534 : @*/
535 9 : PetscErrorCode EPSEVSLGetDOSParameters(EPS eps,EPSEVSLDOSMethod *dos,PetscInt *nvec,PetscInt *deg,PetscInt *steps,PetscInt *npoints)
536 : {
537 9 : PetscFunctionBegin;
538 9 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
539 9 : PetscUseMethod(eps,"EPSEVSLGetDOSParameters_C",(EPS,EPSEVSLDOSMethod*,PetscInt*,PetscInt*,PetscInt*,PetscInt*),(eps,dos,nvec,deg,steps,npoints));
540 9 : PetscFunctionReturn(PETSC_SUCCESS);
541 : }
542 :
543 1 : static PetscErrorCode EPSEVSLSetPolParameters_EVSL(EPS eps,PetscInt max_deg,PetscReal thresh)
544 : {
545 1 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
546 :
547 1 : PetscFunctionBegin;
548 1 : if (max_deg == PETSC_DECIDE || max_deg == PETSC_DEFAULT) ctx->max_deg = 10000;
549 : else {
550 1 : PetscCheck(max_deg>2,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"The max_deg argument must be > 2");
551 1 : ctx->max_deg = max_deg;
552 : }
553 1 : if (thresh == (PetscReal)PETSC_DECIDE || thresh == (PetscReal)PETSC_DEFAULT) ctx->thresh = 0.8;
554 : else {
555 1 : PetscCheck(thresh>0.0,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_OUTOFRANGE,"The thresh argument must be > 0.0");
556 1 : ctx->thresh = thresh;
557 : }
558 1 : eps->state = EPS_STATE_INITIAL;
559 1 : PetscFunctionReturn(PETSC_SUCCESS);
560 : }
561 :
562 : /*@
563 : EPSEVSLSetPolParameters - Defines the parameters used for building the
564 : building the polynomial in the EVSL solver.
565 :
566 : Logically Collective
567 :
568 : Input Parameters:
569 : + eps - the eigensolver context
570 : . max_deg - maximum degree allowed for the polynomial
571 : - thresh - threshold for accepting polynomial
572 :
573 : Options Database Keys:
574 : + -eps_evsl_pol_max_deg <d> - set maximum polynomial degree
575 : - -eps_evsl_pol_thresh <t> - set the threshold
576 :
577 : Level: intermediate
578 :
579 : .seealso: EPSEVSLGetPolParameters()
580 : @*/
581 1 : PetscErrorCode EPSEVSLSetPolParameters(EPS eps,PetscInt max_deg,PetscReal thresh)
582 : {
583 1 : PetscFunctionBegin;
584 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
585 4 : PetscValidLogicalCollectiveInt(eps,max_deg,2);
586 4 : PetscValidLogicalCollectiveReal(eps,thresh,3);
587 1 : PetscTryMethod(eps,"EPSEVSLSetPolParameters_C",(EPS,PetscInt,PetscReal),(eps,max_deg,thresh));
588 1 : PetscFunctionReturn(PETSC_SUCCESS);
589 : }
590 :
591 9 : static PetscErrorCode EPSEVSLGetPolParameters_EVSL(EPS eps,PetscInt *max_deg,PetscReal *thresh)
592 : {
593 9 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
594 :
595 9 : PetscFunctionBegin;
596 9 : if (max_deg) *max_deg = ctx->max_deg;
597 9 : if (thresh) *thresh = ctx->thresh;
598 9 : PetscFunctionReturn(PETSC_SUCCESS);
599 : }
600 :
601 : /*@
602 : EPSEVSLGetPolParameters - Gets the parameters used for building the
603 : polynomial in the EVSL solver.
604 :
605 : Not Collective
606 :
607 : Input Parameter:
608 : . eps - the eigensolver context
609 :
610 : Output Parameters:
611 : + max_deg - the maximum degree of the polynomial
612 : - thresh - the threshold
613 :
614 : Level: intermediate
615 :
616 : .seealso: EPSEVSLSetPolParameters()
617 : @*/
618 9 : PetscErrorCode EPSEVSLGetPolParameters(EPS eps,PetscInt *max_deg,PetscReal *thresh)
619 : {
620 9 : PetscFunctionBegin;
621 9 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
622 9 : PetscUseMethod(eps,"EPSEVSLGetPolParameters_C",(EPS,PetscInt*,PetscReal*),(eps,max_deg,thresh));
623 9 : PetscFunctionReturn(PETSC_SUCCESS);
624 : }
625 :
626 2 : static PetscErrorCode EPSEVSLSetDamping_EVSL(EPS eps,EPSEVSLDamping damping)
627 : {
628 2 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
629 :
630 2 : PetscFunctionBegin;
631 2 : if (ctx->damping != damping) {
632 1 : ctx->damping = damping;
633 1 : eps->state = EPS_STATE_INITIAL;
634 : }
635 2 : PetscFunctionReturn(PETSC_SUCCESS);
636 : }
637 :
638 : /*@
639 : EPSEVSLSetDamping - Set the type of damping to be used in EVSL.
640 :
641 : Logically Collective
642 :
643 : Input Parameters:
644 : + eps - the eigensolver context
645 : - damping - the type of damping
646 :
647 : Options Database Key:
648 : . -eps_evsl_damping <n> - set the type of damping
649 :
650 : Notes:
651 : Damping is applied when building the polynomial to be used when solving the
652 : eigenproblem, and also during estimation of DOS with the KPM method.
653 :
654 : Level: intermediate
655 :
656 : .seealso: EPSEVSLGetDamping(), EPSEVSLSetDOSParameters()
657 : @*/
658 2 : PetscErrorCode EPSEVSLSetDamping(EPS eps,EPSEVSLDamping damping)
659 : {
660 2 : PetscFunctionBegin;
661 2 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
662 8 : PetscValidLogicalCollectiveEnum(eps,damping,2);
663 2 : PetscTryMethod(eps,"EPSEVSLSetDamping_C",(EPS,EPSEVSLDamping),(eps,damping));
664 2 : PetscFunctionReturn(PETSC_SUCCESS);
665 : }
666 :
667 1 : static PetscErrorCode EPSEVSLGetDamping_EVSL(EPS eps,EPSEVSLDamping *damping)
668 : {
669 1 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
670 :
671 1 : PetscFunctionBegin;
672 1 : *damping = ctx->damping;
673 1 : PetscFunctionReturn(PETSC_SUCCESS);
674 : }
675 :
676 : /*@
677 : EPSEVSLGetDamping - Gets the type of damping.
678 :
679 : Not Collective
680 :
681 : Input Parameter:
682 : . eps - the eigensolver context
683 :
684 : Output Parameter:
685 : . damping - the type of damping
686 :
687 : Level: intermediate
688 :
689 : .seealso: EPSEVSLSetDamping()
690 : @*/
691 1 : PetscErrorCode EPSEVSLGetDamping(EPS eps,EPSEVSLDamping *damping)
692 : {
693 1 : PetscFunctionBegin;
694 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
695 1 : PetscAssertPointer(damping,2);
696 1 : PetscUseMethod(eps,"EPSEVSLGetDamping_C",(EPS,EPSEVSLDamping*),(eps,damping));
697 1 : PetscFunctionReturn(PETSC_SUCCESS);
698 : }
699 :
700 1 : static PetscErrorCode EPSView_EVSL(EPS eps,PetscViewer viewer)
701 : {
702 1 : PetscBool isascii;
703 1 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
704 :
705 1 : PetscFunctionBegin;
706 1 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
707 1 : if (isascii) {
708 1 : PetscCall(PetscViewerASCIIPrintf(viewer," numerical range = [%g,%g]\n",(double)ctx->lmin,(double)ctx->lmax));
709 1 : PetscCall(PetscViewerASCIIPrintf(viewer," number of slices = %" PetscInt_FMT "\n",ctx->nslices));
710 1 : PetscCall(PetscViewerASCIIPrintf(viewer," type of damping = %s\n",EPSEVSLDampings[ctx->damping]));
711 1 : PetscCall(PetscViewerASCIIPrintf(viewer," computing DOS with %s: nvec=%" PetscInt_FMT ", ",EPSEVSLDOSMethods[ctx->dos],ctx->nvec));
712 1 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
713 1 : switch (ctx->dos) {
714 1 : case EPS_EVSL_DOS_KPM:
715 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"degree=%" PetscInt_FMT "\n",ctx->deg));
716 : break;
717 0 : case EPS_EVSL_DOS_LANCZOS:
718 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"steps=%" PetscInt_FMT ", npoints=%" PetscInt_FMT "\n",ctx->steps,ctx->npoints));
719 : break;
720 : }
721 1 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
722 1 : PetscCall(PetscViewerASCIIPrintf(viewer," polynomial parameters: max degree = %" PetscInt_FMT ", threshold = %g\n",ctx->max_deg,(double)ctx->thresh));
723 : }
724 1 : PetscFunctionReturn(PETSC_SUCCESS);
725 : }
726 :
727 9 : static PetscErrorCode EPSSetFromOptions_EVSL(EPS eps,PetscOptionItems *PetscOptionsObject)
728 : {
729 9 : PetscReal array[2]={0,0},th;
730 9 : PetscInt k,i1,i2,i3,i4;
731 9 : PetscBool flg,flg1;
732 9 : EPSEVSLDOSMethod dos;
733 9 : EPSEVSLDamping damping;
734 9 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
735 :
736 9 : PetscFunctionBegin;
737 9 : PetscOptionsHeadBegin(PetscOptionsObject,"EPS EVSL Options");
738 :
739 9 : k = 2;
740 9 : PetscCall(PetscOptionsRealArray("-eps_evsl_range","Interval containing all eigenvalues (two real values separated with a comma without spaces)","EPSEVSLSetRange",array,&k,&flg));
741 9 : if (flg) {
742 0 : PetscCheck(k>1,PetscObjectComm((PetscObject)eps),PETSC_ERR_ARG_SIZ,"Must pass two values in -eps_evsl_range (comma-separated without spaces)");
743 0 : PetscCall(EPSEVSLSetRange(eps,array[0],array[1]));
744 : }
745 :
746 9 : PetscCall(PetscOptionsInt("-eps_evsl_slices","Number of slices","EPSEVSLSetSlices",ctx->nslices,&i1,&flg));
747 9 : if (flg) PetscCall(EPSEVSLSetSlices(eps,i1));
748 :
749 9 : PetscCall(PetscOptionsEnum("-eps_evsl_damping","Type of damping","EPSEVSLSetDamping",EPSEVSLDampings,(PetscEnum)ctx->damping,(PetscEnum*)&damping,&flg));
750 9 : if (flg) PetscCall(EPSEVSLSetDamping(eps,damping));
751 :
752 9 : PetscCall(EPSEVSLGetDOSParameters(eps,&dos,&i1,&i2,&i3,&i4));
753 9 : PetscCall(PetscOptionsEnum("-eps_evsl_dos_method","Method to compute the DOS","EPSEVSLSetDOSParameters",EPSEVSLDOSMethods,(PetscEnum)ctx->dos,(PetscEnum*)&dos,&flg));
754 9 : PetscCall(PetscOptionsInt("-eps_evsl_dos_nvec","Number of sample vectors for DOS","EPSEVSLSetDOSParameters",i1,&i1,&flg1));
755 9 : if (flg1) flg = PETSC_TRUE;
756 9 : PetscCall(PetscOptionsInt("-eps_evsl_dos_degree","Polynomial degree used for DOS","EPSEVSLSetDOSParameters",i2,&i2,&flg1));
757 9 : if (flg1) flg = PETSC_TRUE;
758 9 : PetscCall(PetscOptionsInt("-eps_evsl_dos_steps","Number of Lanczos steps in DOS","EPSEVSLSetDOSParameters",i3,&i3,&flg1));
759 9 : if (flg1) flg = PETSC_TRUE;
760 9 : PetscCall(PetscOptionsInt("-eps_evsl_dos_npoints","Number of sample points used for DOS","EPSEVSLSetDOSParameters",i4,&i4,&flg1));
761 9 : if (flg || flg1) PetscCall(EPSEVSLSetDOSParameters(eps,dos,i1,i2,i3,i4));
762 :
763 9 : PetscCall(EPSEVSLGetPolParameters(eps,&i1,&th));
764 9 : PetscCall(PetscOptionsInt("-eps_evsl_pol_max_deg","Maximum degree allowed for the polynomial","EPSEVSLSetPolParameters",i1,&i1,&flg));
765 9 : PetscCall(PetscOptionsReal("-eps_evsl_pol_threshold","Threshold for accepting polynomial","EPSEVSLSetPolParameters",th,&th,&flg1));
766 9 : if (flg || flg1) PetscCall(EPSEVSLSetPolParameters(eps,i1,th));
767 :
768 9 : PetscOptionsHeadEnd();
769 9 : PetscFunctionReturn(PETSC_SUCCESS);
770 : }
771 :
772 9 : static PetscErrorCode EPSDestroy_EVSL(EPS eps)
773 : {
774 9 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
775 :
776 9 : PetscFunctionBegin;
777 9 : if (ctx->initialized) EVSLFinish();
778 9 : PetscCall(PetscLayoutDestroy(&ctx->map));
779 9 : PetscCall(PetscFree(ctx->sli));
780 9 : PetscCall(PetscFree(eps->data));
781 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetRange_C",NULL));
782 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetRange_C",NULL));
783 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetSlices_C",NULL));
784 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetSlices_C",NULL));
785 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetDOSParameters_C",NULL));
786 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetDOSParameters_C",NULL));
787 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetPolParameters_C",NULL));
788 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetPolParameters_C",NULL));
789 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetDamping_C",NULL));
790 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetDamping_C",NULL));
791 9 : PetscFunctionReturn(PETSC_SUCCESS);
792 : }
793 :
794 12 : static PetscErrorCode EPSReset_EVSL(EPS eps)
795 : {
796 12 : EPS_EVSL *ctx = (EPS_EVSL*)eps->data;
797 :
798 12 : PetscFunctionBegin;
799 12 : PetscCall(MatDestroy(&ctx->A));
800 12 : PetscCall(VecDestroy(&ctx->x));
801 12 : PetscCall(VecDestroy(&ctx->y));
802 12 : PetscFunctionReturn(PETSC_SUCCESS);
803 : }
804 :
805 9 : SLEPC_EXTERN PetscErrorCode EPSCreate_EVSL(EPS eps)
806 : {
807 9 : EPS_EVSL *ctx;
808 :
809 9 : PetscFunctionBegin;
810 9 : PetscCall(PetscNew(&ctx));
811 9 : eps->data = (void*)ctx;
812 :
813 9 : ctx->nslices = 0;
814 9 : ctx->lmin = PETSC_MIN_REAL;
815 9 : ctx->lmax = PETSC_MAX_REAL;
816 9 : ctx->dos = EPS_EVSL_DOS_KPM;
817 9 : ctx->nvec = 80;
818 9 : ctx->deg = 300;
819 9 : ctx->steps = 40;
820 9 : ctx->npoints = 200;
821 9 : ctx->max_deg = 10000;
822 9 : ctx->thresh = 0.8;
823 9 : ctx->damping = EPS_EVSL_DAMPING_SIGMA;
824 :
825 9 : eps->categ = EPS_CATEGORY_OTHER;
826 :
827 9 : eps->ops->solve = EPSSolve_EVSL;
828 9 : eps->ops->setup = EPSSetUp_EVSL;
829 9 : eps->ops->setupsort = EPSSetUpSort_Basic;
830 9 : eps->ops->setfromoptions = EPSSetFromOptions_EVSL;
831 9 : eps->ops->destroy = EPSDestroy_EVSL;
832 9 : eps->ops->reset = EPSReset_EVSL;
833 9 : eps->ops->view = EPSView_EVSL;
834 9 : eps->ops->backtransform = EPSBackTransform_Default;
835 9 : eps->ops->setdefaultst = EPSSetDefaultST_NoFactor;
836 :
837 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetRange_C",EPSEVSLSetRange_EVSL));
838 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetRange_C",EPSEVSLGetRange_EVSL));
839 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetSlices_C",EPSEVSLSetSlices_EVSL));
840 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetSlices_C",EPSEVSLGetSlices_EVSL));
841 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetDOSParameters_C",EPSEVSLSetDOSParameters_EVSL));
842 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetDOSParameters_C",EPSEVSLGetDOSParameters_EVSL));
843 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetPolParameters_C",EPSEVSLSetPolParameters_EVSL));
844 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetPolParameters_C",EPSEVSLGetPolParameters_EVSL));
845 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLSetDamping_C",EPSEVSLSetDamping_EVSL));
846 9 : PetscCall(PetscObjectComposeFunction((PetscObject)eps,"EPSEVSLGetDamping_C",EPSEVSLGetDamping_EVSL));
847 9 : PetscFunctionReturn(PETSC_SUCCESS);
848 : }
|