Line data Source code
1 : /*
2 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3 : SLEPc - Scalable Library for Eigenvalue Problem Computations
4 : Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5 :
6 : This file is part of SLEPc.
7 : SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9 : */
10 : /*
11 : NEP routines related to monitors
12 : */
13 :
14 : #include <slepc/private/nepimpl.h> /*I "slepcnep.h" I*/
15 : #include <petscdraw.h>
16 :
17 3 : PetscErrorCode NEPMonitorLGCreate(MPI_Comm comm,const char host[],const char label[],const char metric[],PetscInt l,const char *names[],int x,int y,int m,int n,PetscDrawLG *lgctx)
18 : {
19 3 : PetscDraw draw;
20 3 : PetscDrawAxis axis;
21 3 : PetscDrawLG lg;
22 :
23 3 : PetscFunctionBegin;
24 3 : PetscCall(PetscDrawCreate(comm,host,label,x,y,m,n,&draw));
25 3 : PetscCall(PetscDrawSetFromOptions(draw));
26 3 : PetscCall(PetscDrawLGCreate(draw,l,&lg));
27 3 : if (names) PetscCall(PetscDrawLGSetLegend(lg,names));
28 3 : PetscCall(PetscDrawLGSetFromOptions(lg));
29 3 : PetscCall(PetscDrawLGGetAxis(lg,&axis));
30 3 : PetscCall(PetscDrawAxisSetLabels(axis,"Convergence","Iteration",metric));
31 3 : PetscCall(PetscDrawDestroy(&draw));
32 3 : *lgctx = lg;
33 3 : PetscFunctionReturn(PETSC_SUCCESS);
34 : }
35 :
36 : /*
37 : Runs the user provided monitor routines, if any.
38 : */
39 1163 : PetscErrorCode NEPMonitor(NEP nep,PetscInt it,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest)
40 : {
41 1163 : PetscInt i,n = nep->numbermonitors;
42 :
43 1163 : PetscFunctionBegin;
44 1208 : for (i=0;i<n;i++) PetscCall((*nep->monitor[i])(nep,it,nconv,eigr,eigi,errest,nest,nep->monitorcontext[i]));
45 1163 : PetscFunctionReturn(PETSC_SUCCESS);
46 : }
47 :
48 : /*@C
49 : NEPMonitorSet - Sets an ADDITIONAL function to be called at every
50 : iteration to monitor the error estimates for each requested eigenpair.
51 :
52 : Logically Collective
53 :
54 : Input Parameters:
55 : + nep - eigensolver context obtained from NEPCreate()
56 : . monitor - pointer to function (if this is NULL, it turns off monitoring)
57 : . mctx - [optional] context for private data for the
58 : monitor routine (use NULL if no context is desired)
59 : - monitordestroy - [optional] routine that frees monitor context (may be NULL)
60 :
61 : Calling sequence of monitor:
62 : $ PetscErrorCode monitor(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx)
63 : + nep - nonlinear eigensolver context obtained from NEPCreate()
64 : . its - iteration number
65 : . nconv - number of converged eigenpairs
66 : . eigr - real part of the eigenvalues
67 : . eigi - imaginary part of the eigenvalues
68 : . errest - error estimates for each eigenpair
69 : . nest - number of error estimates
70 : - mctx - optional monitoring context, as set by NEPMonitorSet()
71 :
72 : Options Database Keys:
73 : + -nep_monitor - print only the first error estimate
74 : . -nep_monitor_all - print error estimates at each iteration
75 : . -nep_monitor_conv - print the eigenvalue approximations only when
76 : convergence has been reached
77 : . -nep_monitor draw::draw_lg - sets line graph monitor for the first unconverged
78 : approximate eigenvalue
79 : . -nep_monitor_all draw::draw_lg - sets line graph monitor for all unconverged
80 : approximate eigenvalues
81 : . -nep_monitor_conv draw::draw_lg - sets line graph monitor for convergence history
82 : - -nep_monitor_cancel - cancels all monitors that have been hardwired into
83 : a code by calls to NEPMonitorSet(), but does not cancel those set via
84 : the options database.
85 :
86 : Notes:
87 : Several different monitoring routines may be set by calling
88 : NEPMonitorSet() multiple times; all will be called in the
89 : order in which they were set.
90 :
91 : Level: intermediate
92 :
93 : .seealso: NEPMonitorFirst(), NEPMonitorAll(), NEPMonitorCancel()
94 : @*/
95 7 : PetscErrorCode NEPMonitorSet(NEP nep,PetscErrorCode (*monitor)(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,void *mctx),void *mctx,PetscErrorCode (*monitordestroy)(void**))
96 : {
97 7 : PetscFunctionBegin;
98 7 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
99 7 : PetscCheck(nep->numbermonitors<MAXNEPMONITORS,PetscObjectComm((PetscObject)nep),PETSC_ERR_ARG_OUTOFRANGE,"Too many NEP monitors set");
100 7 : nep->monitor[nep->numbermonitors] = monitor;
101 7 : nep->monitorcontext[nep->numbermonitors] = (void*)mctx;
102 7 : nep->monitordestroy[nep->numbermonitors++] = monitordestroy;
103 7 : PetscFunctionReturn(PETSC_SUCCESS);
104 : }
105 :
106 : /*@
107 : NEPMonitorCancel - Clears all monitors for a NEP object.
108 :
109 : Logically Collective
110 :
111 : Input Parameters:
112 : . nep - eigensolver context obtained from NEPCreate()
113 :
114 : Options Database Key:
115 : . -nep_monitor_cancel - Cancels all monitors that have been hardwired
116 : into a code by calls to NEPMonitorSet(),
117 : but does not cancel those set via the options database.
118 :
119 : Level: intermediate
120 :
121 : .seealso: NEPMonitorSet()
122 : @*/
123 144 : PetscErrorCode NEPMonitorCancel(NEP nep)
124 : {
125 144 : PetscInt i;
126 :
127 144 : PetscFunctionBegin;
128 144 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
129 151 : for (i=0; i<nep->numbermonitors; i++) {
130 7 : if (nep->monitordestroy[i]) PetscCall((*nep->monitordestroy[i])(&nep->monitorcontext[i]));
131 : }
132 144 : nep->numbermonitors = 0;
133 144 : PetscFunctionReturn(PETSC_SUCCESS);
134 : }
135 :
136 : /*@C
137 : NEPGetMonitorContext - Gets the monitor context, as set by
138 : NEPMonitorSet() for the FIRST monitor only.
139 :
140 : Not Collective
141 :
142 : Input Parameter:
143 : . nep - eigensolver context obtained from NEPCreate()
144 :
145 : Output Parameter:
146 : . ctx - monitor context
147 :
148 : Level: intermediate
149 :
150 : .seealso: NEPMonitorSet(), NEPDefaultMonitor()
151 : @*/
152 0 : PetscErrorCode NEPGetMonitorContext(NEP nep,void *ctx)
153 : {
154 0 : PetscFunctionBegin;
155 0 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
156 0 : *(void**)ctx = nep->monitorcontext[0];
157 0 : PetscFunctionReturn(PETSC_SUCCESS);
158 : }
159 :
160 : /*@C
161 : NEPMonitorFirst - Print the first unconverged approximate value and
162 : error estimate at each iteration of the nonlinear eigensolver.
163 :
164 : Collective
165 :
166 : Input Parameters:
167 : + nep - nonlinear eigensolver context
168 : . its - iteration number
169 : . nconv - number of converged eigenpairs so far
170 : . eigr - real part of the eigenvalues
171 : . eigi - imaginary part of the eigenvalues
172 : . errest - error estimates
173 : . nest - number of error estimates to display
174 : - vf - viewer and format for monitoring
175 :
176 : Options Database Key:
177 : . -nep_monitor - activates NEPMonitorFirst()
178 :
179 : Level: intermediate
180 :
181 : .seealso: NEPMonitorSet(), NEPMonitorAll(), NEPMonitorConverged()
182 : @*/
183 9 : PetscErrorCode NEPMonitorFirst(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
184 : {
185 9 : PetscViewer viewer = vf->viewer;
186 :
187 9 : PetscFunctionBegin;
188 9 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
189 9 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,8);
190 9 : if (its==1 && ((PetscObject)nep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix));
191 9 : if (nconv<nest) {
192 8 : PetscCall(PetscViewerPushFormat(viewer,vf->format));
193 8 : PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
194 8 : PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP nconv=%" PetscInt_FMT " first unconverged value (error)",its,nconv));
195 8 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
196 : #if defined(PETSC_USE_COMPLEX)
197 8 : PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[nconv]),(double)PetscImaginaryPart(eigr[nconv])));
198 : #else
199 : PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)eigr[nconv]));
200 : if (eigi[nconv]!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[nconv]));
201 : #endif
202 8 : PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[nconv]));
203 8 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
204 8 : PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
205 8 : PetscCall(PetscViewerPopFormat(viewer));
206 : }
207 9 : PetscFunctionReturn(PETSC_SUCCESS);
208 : }
209 :
210 : /*@C
211 : NEPMonitorAll - Print the current approximate values and
212 : error estimates at each iteration of the nonlinear eigensolver.
213 :
214 : Collective
215 :
216 : Input Parameters:
217 : + nep - nonlinear eigensolver context
218 : . its - iteration number
219 : . nconv - number of converged eigenpairs so far
220 : . eigr - real part of the eigenvalues
221 : . eigi - imaginary part of the eigenvalues
222 : . errest - error estimates
223 : . nest - number of error estimates to display
224 : - vf - viewer and format for monitoring
225 :
226 : Options Database Key:
227 : . -nep_monitor_all - activates NEPMonitorAll()
228 :
229 : Level: intermediate
230 :
231 : .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorConverged()
232 : @*/
233 5 : PetscErrorCode NEPMonitorAll(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
234 : {
235 5 : PetscInt i;
236 5 : PetscViewer viewer = vf->viewer;
237 :
238 5 : PetscFunctionBegin;
239 5 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
240 5 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,8);
241 5 : PetscCall(PetscViewerPushFormat(viewer,vf->format));
242 5 : PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
243 5 : if (its==1 && ((PetscObject)nep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Eigenvalue approximations and residual norms for %s solve.\n",((PetscObject)nep)->prefix));
244 5 : PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP nconv=%" PetscInt_FMT " Values (Errors)",its,nconv));
245 5 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
246 10 : for (i=0;i<nest;i++) {
247 : #if defined(PETSC_USE_COMPLEX)
248 5 : PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i])));
249 : #else
250 : PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]));
251 : if (eigi[i]!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]));
252 : #endif
253 5 : PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)",(double)errest[i]));
254 : }
255 5 : PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
256 5 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
257 5 : PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
258 5 : PetscCall(PetscViewerPopFormat(viewer));
259 5 : PetscFunctionReturn(PETSC_SUCCESS);
260 : }
261 :
262 : /*@C
263 : NEPMonitorConverged - Print the approximate values and
264 : error estimates as they converge.
265 :
266 : Collective
267 :
268 : Input Parameters:
269 : + nep - nonlinear eigensolver context
270 : . its - iteration number
271 : . nconv - number of converged eigenpairs so far
272 : . eigr - real part of the eigenvalues
273 : . eigi - imaginary part of the eigenvalues
274 : . errest - error estimates
275 : . nest - number of error estimates to display
276 : - vf - viewer and format for monitoring
277 :
278 : Options Database Key:
279 : . -nep_monitor_conv - activates NEPMonitorConverged()
280 :
281 : Level: intermediate
282 :
283 : .seealso: NEPMonitorSet(), NEPMonitorFirst(), NEPMonitorAll()
284 : @*/
285 6 : PetscErrorCode NEPMonitorConverged(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
286 : {
287 6 : PetscInt i;
288 6 : PetscViewer viewer = vf->viewer;
289 6 : SlepcConvMon ctx;
290 :
291 6 : PetscFunctionBegin;
292 6 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
293 6 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,8);
294 6 : ctx = (SlepcConvMon)vf->data;
295 6 : if (its==1 && ((PetscObject)nep)->prefix) PetscCall(PetscViewerASCIIPrintf(viewer," Convergence history for %s solve.\n",((PetscObject)nep)->prefix));
296 1 : if (its==1) ctx->oldnconv = 0;
297 6 : if (ctx->oldnconv!=nconv) {
298 1 : PetscCall(PetscViewerPushFormat(viewer,vf->format));
299 1 : PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
300 2 : for (i=ctx->oldnconv;i<nconv;i++) {
301 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"%3" PetscInt_FMT " NEP converged value (error) #%" PetscInt_FMT,its,i));
302 1 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
303 : #if defined(PETSC_USE_COMPLEX)
304 1 : PetscCall(PetscViewerASCIIPrintf(viewer," %g%+gi",(double)PetscRealPart(eigr[i]),(double)PetscImaginaryPart(eigr[i])));
305 : #else
306 : PetscCall(PetscViewerASCIIPrintf(viewer," %g",(double)eigr[i]));
307 : if (eigi[i]!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%+gi",(double)eigi[i]));
308 : #endif
309 1 : PetscCall(PetscViewerASCIIPrintf(viewer," (%10.8e)\n",(double)errest[i]));
310 1 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
311 : }
312 1 : PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
313 1 : PetscCall(PetscViewerPopFormat(viewer));
314 1 : ctx->oldnconv = nconv;
315 : }
316 6 : PetscFunctionReturn(PETSC_SUCCESS);
317 : }
318 :
319 1 : PetscErrorCode NEPMonitorConvergedCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
320 : {
321 1 : SlepcConvMon mctx;
322 :
323 1 : PetscFunctionBegin;
324 1 : PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
325 1 : PetscCall(PetscNew(&mctx));
326 1 : mctx->ctx = ctx;
327 1 : (*vf)->data = (void*)mctx;
328 1 : PetscFunctionReturn(PETSC_SUCCESS);
329 : }
330 :
331 1 : PetscErrorCode NEPMonitorConvergedDestroy(PetscViewerAndFormat **vf)
332 : {
333 1 : PetscFunctionBegin;
334 1 : if (!*vf) PetscFunctionReturn(PETSC_SUCCESS);
335 1 : PetscCall(PetscFree((*vf)->data));
336 1 : PetscCall(PetscOptionsRestoreViewer(&(*vf)->viewer));
337 1 : PetscCall(PetscDrawLGDestroy(&(*vf)->lg));
338 1 : PetscCall(PetscFree(*vf));
339 1 : PetscFunctionReturn(PETSC_SUCCESS);
340 : }
341 :
342 : /*@C
343 : NEPMonitorFirstDrawLG - Plots the error estimate of the first unconverged
344 : approximation at each iteration of the nonlinear eigensolver.
345 :
346 : Collective
347 :
348 : Input Parameters:
349 : + nep - nonlinear eigensolver context
350 : . its - iteration number
351 : . nconv - number of converged eigenpairs so far
352 : . eigr - real part of the eigenvalues
353 : . eigi - imaginary part of the eigenvalues
354 : . errest - error estimates
355 : . nest - number of error estimates to display
356 : - vf - viewer and format for monitoring
357 :
358 : Options Database Key:
359 : . -nep_monitor draw::draw_lg - activates NEPMonitorFirstDrawLG()
360 :
361 : Level: intermediate
362 :
363 : .seealso: NEPMonitorSet()
364 : @*/
365 13 : PetscErrorCode NEPMonitorFirstDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
366 : {
367 13 : PetscViewer viewer = vf->viewer;
368 13 : PetscDrawLG lg = vf->lg;
369 13 : PetscReal x,y;
370 :
371 13 : PetscFunctionBegin;
372 13 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
373 13 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,8);
374 13 : PetscValidHeaderSpecific(lg,PETSC_DRAWLG_CLASSID,8);
375 13 : PetscCall(PetscViewerPushFormat(viewer,vf->format));
376 13 : if (its==1) {
377 2 : PetscCall(PetscDrawLGReset(lg));
378 2 : PetscCall(PetscDrawLGSetDimension(lg,1));
379 2 : PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0));
380 : }
381 13 : if (nconv<nest) {
382 12 : x = (PetscReal)its;
383 12 : if (errest[nconv] > 0.0) y = PetscLog10Real(errest[nconv]);
384 0 : else y = 0.0;
385 12 : PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
386 12 : if (its <= 20 || !(its % 5) || nep->reason) {
387 12 : PetscCall(PetscDrawLGDraw(lg));
388 12 : PetscCall(PetscDrawLGSave(lg));
389 : }
390 : }
391 13 : PetscCall(PetscViewerPopFormat(viewer));
392 13 : PetscFunctionReturn(PETSC_SUCCESS);
393 : }
394 :
395 : /*@C
396 : NEPMonitorFirstDrawLGCreate - Creates the plotter for the first error estimate.
397 :
398 : Collective
399 :
400 : Input Parameters:
401 : + viewer - the viewer
402 : . format - the viewer format
403 : - ctx - an optional user context
404 :
405 : Output Parameter:
406 : . vf - the viewer and format context
407 :
408 : Level: intermediate
409 :
410 : .seealso: NEPMonitorSet()
411 : @*/
412 2 : PetscErrorCode NEPMonitorFirstDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
413 : {
414 2 : PetscFunctionBegin;
415 2 : PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
416 2 : (*vf)->data = ctx;
417 2 : PetscCall(NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"First Error Estimate","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
418 2 : PetscFunctionReturn(PETSC_SUCCESS);
419 : }
420 :
421 : /*@C
422 : NEPMonitorAllDrawLG - Plots the error estimate of all unconverged
423 : approximations at each iteration of the nonlinear eigensolver.
424 :
425 : Collective
426 :
427 : Input Parameters:
428 : + nep - nonlinear eigensolver context
429 : . its - iteration number
430 : . nconv - number of converged eigenpairs so far
431 : . eigr - real part of the eigenvalues
432 : . eigi - imaginary part of the eigenvalues
433 : . errest - error estimates
434 : . nest - number of error estimates to display
435 : - vf - viewer and format for monitoring
436 :
437 : Options Database Key:
438 : . -nep_monitor_all draw::draw_lg - activates NEPMonitorAllDrawLG()
439 :
440 : Level: intermediate
441 :
442 : .seealso: NEPMonitorSet()
443 : @*/
444 12 : PetscErrorCode NEPMonitorAllDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
445 : {
446 12 : PetscViewer viewer = vf->viewer;
447 12 : PetscDrawLG lg = vf->lg;
448 12 : PetscInt i,n = PetscMin(nep->nev,255);
449 12 : PetscReal *x,*y;
450 :
451 12 : PetscFunctionBegin;
452 12 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
453 12 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,8);
454 12 : PetscValidHeaderSpecific(lg,PETSC_DRAWLG_CLASSID,8);
455 12 : PetscCall(PetscViewerPushFormat(viewer,vf->format));
456 12 : if (its==1) {
457 1 : PetscCall(PetscDrawLGReset(lg));
458 1 : PetscCall(PetscDrawLGSetDimension(lg,n));
459 1 : PetscCall(PetscDrawLGSetLimits(lg,1,1,PetscLog10Real(nep->tol)-2,0.0));
460 : }
461 12 : PetscCall(PetscMalloc2(n,&x,n,&y));
462 60 : for (i=0;i<n;i++) {
463 48 : x[i] = (PetscReal)its;
464 48 : if (i < nest && errest[i] > 0.0) y[i] = PetscLog10Real(errest[i]);
465 18 : else y[i] = 0.0;
466 : }
467 12 : PetscCall(PetscDrawLGAddPoint(lg,x,y));
468 12 : if (its <= 20 || !(its % 5) || nep->reason) {
469 12 : PetscCall(PetscDrawLGDraw(lg));
470 12 : PetscCall(PetscDrawLGSave(lg));
471 : }
472 12 : PetscCall(PetscFree2(x,y));
473 12 : PetscCall(PetscViewerPopFormat(viewer));
474 12 : PetscFunctionReturn(PETSC_SUCCESS);
475 : }
476 :
477 : /*@C
478 : NEPMonitorAllDrawLGCreate - Creates the plotter for all the error estimates.
479 :
480 : Collective
481 :
482 : Input Parameters:
483 : + viewer - the viewer
484 : . format - the viewer format
485 : - ctx - an optional user context
486 :
487 : Output Parameter:
488 : . vf - the viewer and format context
489 :
490 : Level: intermediate
491 :
492 : .seealso: NEPMonitorSet()
493 : @*/
494 1 : PetscErrorCode NEPMonitorAllDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
495 : {
496 1 : PetscFunctionBegin;
497 1 : PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
498 1 : (*vf)->data = ctx;
499 1 : PetscCall(NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"All Error Estimates","Log Error Estimate",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
500 1 : PetscFunctionReturn(PETSC_SUCCESS);
501 : }
502 :
503 : /*@C
504 : NEPMonitorConvergedDrawLG - Plots the number of converged eigenvalues
505 : at each iteration of the nonlinear eigensolver.
506 :
507 : Collective
508 :
509 : Input Parameters:
510 : + nep - nonlinear eigensolver context
511 : . its - iteration number
512 : . nconv - number of converged eigenpairs so far
513 : . eigr - real part of the eigenvalues
514 : . eigi - imaginary part of the eigenvalues
515 : . errest - error estimates
516 : . nest - number of error estimates to display
517 : - vf - viewer and format for monitoring
518 :
519 : Options Database Key:
520 : . -nep_monitor_conv draw::draw_lg - activates NEPMonitorConvergedDrawLG()
521 :
522 : Level: intermediate
523 :
524 : .seealso: NEPMonitorSet()
525 : @*/
526 0 : PetscErrorCode NEPMonitorConvergedDrawLG(NEP nep,PetscInt its,PetscInt nconv,PetscScalar *eigr,PetscScalar *eigi,PetscReal *errest,PetscInt nest,PetscViewerAndFormat *vf)
527 : {
528 0 : PetscViewer viewer = vf->viewer;
529 0 : PetscDrawLG lg = vf->lg;
530 0 : PetscReal x,y;
531 :
532 0 : PetscFunctionBegin;
533 0 : PetscValidHeaderSpecific(nep,NEP_CLASSID,1);
534 0 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,8);
535 0 : PetscValidHeaderSpecific(lg,PETSC_DRAWLG_CLASSID,8);
536 0 : PetscCall(PetscViewerPushFormat(viewer,vf->format));
537 0 : if (its==1) {
538 0 : PetscCall(PetscDrawLGReset(lg));
539 0 : PetscCall(PetscDrawLGSetDimension(lg,1));
540 0 : PetscCall(PetscDrawLGSetLimits(lg,1,1,0,nep->nev));
541 : }
542 0 : x = (PetscReal)its;
543 0 : y = (PetscReal)nep->nconv;
544 0 : PetscCall(PetscDrawLGAddPoint(lg,&x,&y));
545 0 : if (its <= 20 || !(its % 5) || nep->reason) {
546 0 : PetscCall(PetscDrawLGDraw(lg));
547 0 : PetscCall(PetscDrawLGSave(lg));
548 : }
549 0 : PetscCall(PetscViewerPopFormat(viewer));
550 0 : PetscFunctionReturn(PETSC_SUCCESS);
551 : }
552 :
553 : /*@C
554 : NEPMonitorConvergedDrawLGCreate - Creates the plotter for the convergence history.
555 :
556 : Collective
557 :
558 : Input Parameters:
559 : + viewer - the viewer
560 : . format - the viewer format
561 : - ctx - an optional user context
562 :
563 : Output Parameter:
564 : . vf - the viewer and format context
565 :
566 : Level: intermediate
567 :
568 : .seealso: NEPMonitorSet()
569 : @*/
570 0 : PetscErrorCode NEPMonitorConvergedDrawLGCreate(PetscViewer viewer,PetscViewerFormat format,void *ctx,PetscViewerAndFormat **vf)
571 : {
572 0 : SlepcConvMon mctx;
573 :
574 0 : PetscFunctionBegin;
575 0 : PetscCall(PetscViewerAndFormatCreate(viewer,format,vf));
576 0 : PetscCall(PetscNew(&mctx));
577 0 : mctx->ctx = ctx;
578 0 : (*vf)->data = (void*)mctx;
579 0 : PetscCall(NEPMonitorLGCreate(PetscObjectComm((PetscObject)viewer),NULL,"Convergence History","Number Converged",1,NULL,PETSC_DECIDE,PETSC_DECIDE,400,300,&(*vf)->lg));
580 0 : PetscFunctionReturn(PETSC_SUCCESS);
581 : }
|