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 : EPS routines related to various viewers
12 : */
13 :
14 : #include <slepc/private/epsimpl.h> /*I "slepceps.h" I*/
15 : #include <slepc/private/bvimpl.h>
16 : #include <petscdraw.h>
17 :
18 : /*@C
19 : EPSView - Prints the EPS data structure.
20 :
21 : Collective
22 :
23 : Input Parameters:
24 : + eps - the eigenproblem solver context
25 : - viewer - optional visualization context
26 :
27 : Options Database Key:
28 : . -eps_view - Calls EPSView() at end of EPSSolve()
29 :
30 : Note:
31 : The available visualization contexts include
32 : + PETSC_VIEWER_STDOUT_SELF - standard output (default)
33 : - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
34 : output where only the first processor opens
35 : the file. All other processors send their
36 : data to the first processor to print.
37 :
38 : The user can open an alternative visualization context with
39 : PetscViewerASCIIOpen() - output to a specified file.
40 :
41 : Level: beginner
42 :
43 : .seealso: STView()
44 : @*/
45 9 : PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
46 : {
47 9 : const char *type=NULL,*extr=NULL,*bal=NULL;
48 9 : char str[50];
49 9 : PetscBool isascii,isexternal,istrivial;
50 :
51 9 : PetscFunctionBegin;
52 9 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
53 9 : if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
54 9 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2);
55 9 : PetscCheckSameComm(eps,1,viewer,2);
56 :
57 9 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
58 9 : if (isascii) {
59 9 : PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer));
60 9 : PetscCall(PetscViewerASCIIPushTab(viewer));
61 9 : PetscTryTypeMethod(eps,view,viewer);
62 9 : PetscCall(PetscViewerASCIIPopTab(viewer));
63 9 : if (eps->problem_type) {
64 9 : switch (eps->problem_type) {
65 6 : case EPS_HEP: type = SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
66 1 : case EPS_GHEP: type = "generalized " SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
67 2 : case EPS_NHEP: type = "non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
68 0 : case EPS_GNHEP: type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
69 0 : case EPS_PGNHEP: type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem with " SLEPC_STRING_HERMITIAN " positive definite B"; break;
70 0 : case EPS_GHIEP: type = "generalized " SLEPC_STRING_HERMITIAN "-indefinite eigenvalue problem"; break;
71 : }
72 : } else type = "not yet set";
73 9 : PetscCall(PetscViewerASCIIPrintf(viewer," problem type: %s\n",type));
74 9 : if (eps->extraction) {
75 2 : switch (eps->extraction) {
76 0 : case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
77 2 : case EPS_HARMONIC: extr = "harmonic Ritz"; break;
78 0 : case EPS_HARMONIC_RELATIVE: extr = "relative harmonic Ritz"; break;
79 0 : case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
80 0 : case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
81 0 : case EPS_REFINED: extr = "refined Ritz"; break;
82 0 : case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
83 : }
84 2 : PetscCall(PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr));
85 : }
86 9 : if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
87 1 : switch (eps->balance) {
88 : case EPS_BALANCE_NONE: break;
89 1 : case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
90 0 : case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
91 0 : case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
92 : }
93 1 : PetscCall(PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal));
94 1 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
95 1 : if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) PetscCall(PetscViewerASCIIPrintf(viewer,", with its=%" PetscInt_FMT,eps->balance_its));
96 1 : if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff));
97 1 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
98 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
99 : }
100 9 : PetscCall(PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: "));
101 9 : PetscCall(SlepcSNPrintfScalar(str,sizeof(str),eps->target,PETSC_FALSE));
102 9 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
103 9 : if (!eps->which) PetscCall(PetscViewerASCIIPrintf(viewer,"not yet set\n"));
104 9 : else switch (eps->which) {
105 0 : case EPS_WHICH_USER:
106 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"user defined\n"));
107 : break;
108 3 : case EPS_TARGET_MAGNITUDE:
109 3 : PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str));
110 : break;
111 0 : case EPS_TARGET_REAL:
112 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str));
113 : break;
114 0 : case EPS_TARGET_IMAGINARY:
115 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str));
116 : break;
117 1 : case EPS_LARGEST_MAGNITUDE:
118 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n"));
119 : break;
120 0 : case EPS_SMALLEST_MAGNITUDE:
121 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n"));
122 : break;
123 3 : case EPS_LARGEST_REAL:
124 3 : PetscCall(PetscViewerASCIIPrintf(viewer,"largest real parts\n"));
125 : break;
126 1 : case EPS_SMALLEST_REAL:
127 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"smallest real parts\n"));
128 : break;
129 0 : case EPS_LARGEST_IMAGINARY:
130 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n"));
131 : break;
132 0 : case EPS_SMALLEST_IMAGINARY:
133 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n"));
134 : break;
135 1 : case EPS_ALL:
136 1 : if (eps->inta || eps->intb) PetscCall(PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb));
137 0 : else PetscCall(PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n"));
138 : break;
139 : }
140 9 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
141 9 : if (eps->twosided && eps->problem_type!=EPS_HEP && eps->problem_type!=EPS_GHEP) PetscCall(PetscViewerASCIIPrintf(viewer," using two-sided variant (for left eigenvectors)\n"));
142 9 : if (eps->purify) PetscCall(PetscViewerASCIIPrintf(viewer," postprocessing eigenvectors with purification\n"));
143 9 : if (eps->trueres) PetscCall(PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n"));
144 9 : if (eps->trackall) PetscCall(PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n"));
145 9 : PetscCall(PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %" PetscInt_FMT "\n",eps->nev));
146 9 : PetscCall(PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",eps->ncv));
147 9 : PetscCall(PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",eps->mpd));
148 9 : PetscCall(PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",eps->max_it));
149 9 : PetscCall(PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)eps->tol));
150 9 : PetscCall(PetscViewerASCIIPrintf(viewer," convergence test: "));
151 9 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
152 9 : switch (eps->conv) {
153 0 : case EPS_CONV_ABS:
154 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"absolute\n"));break;
155 7 : case EPS_CONV_REL:
156 7 : PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n"));break;
157 2 : case EPS_CONV_NORM:
158 2 : PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n"));
159 2 : PetscCall(PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)eps->nrma));
160 2 : if (eps->isgeneralized) PetscCall(PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb));
161 2 : PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
162 : break;
163 0 : case EPS_CONV_USER:
164 0 : PetscCall(PetscViewerASCIIPrintf(viewer,"user-defined\n"));break;
165 : }
166 9 : PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
167 9 : if (eps->nini) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(eps->nini)));
168 9 : if (eps->ninil) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided left initial space: %" PetscInt_FMT "\n",PetscAbs(eps->ninil)));
169 9 : if (eps->nds) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %" PetscInt_FMT "\n",PetscAbs(eps->nds)));
170 0 : } else PetscTryTypeMethod(eps,view,viewer);
171 9 : PetscCall(PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSELEMENTAL,EPSFEAST,EPSPRIMME,EPSSCALAPACK,EPSELPA,EPSEVSL,EPSTRLAN,""));
172 9 : if (!isexternal) {
173 7 : PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
174 7 : if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
175 7 : PetscCall(BVView(eps->V,viewer));
176 7 : if (eps->rg) {
177 7 : PetscCall(RGIsTrivial(eps->rg,&istrivial));
178 7 : if (!istrivial) PetscCall(RGView(eps->rg,viewer));
179 : }
180 7 : if (eps->useds) {
181 7 : if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
182 7 : PetscCall(DSView(eps->ds,viewer));
183 : }
184 7 : PetscCall(PetscViewerPopFormat(viewer));
185 : }
186 9 : if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
187 9 : PetscCall(STView(eps->st,viewer));
188 9 : PetscFunctionReturn(PETSC_SUCCESS);
189 : }
190 :
191 : /*@C
192 : EPSViewFromOptions - View from options
193 :
194 : Collective
195 :
196 : Input Parameters:
197 : + eps - the eigensolver context
198 : . obj - optional object
199 : - name - command line option
200 :
201 : Level: intermediate
202 :
203 : .seealso: EPSView(), EPSCreate()
204 : @*/
205 1710 : PetscErrorCode EPSViewFromOptions(EPS eps,PetscObject obj,const char name[])
206 : {
207 1710 : PetscFunctionBegin;
208 1710 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
209 1710 : PetscCall(PetscObjectViewFromOptions((PetscObject)eps,obj,name));
210 1710 : PetscFunctionReturn(PETSC_SUCCESS);
211 : }
212 :
213 : /*@C
214 : EPSConvergedReasonView - Displays the reason an EPS solve converged or diverged.
215 :
216 : Collective
217 :
218 : Input Parameters:
219 : + eps - the eigensolver context
220 : - viewer - the viewer to display the reason
221 :
222 : Options Database Keys:
223 : . -eps_converged_reason - print reason for convergence, and number of iterations
224 :
225 : Note:
226 : To change the format of the output call PetscViewerPushFormat(viewer,format) before
227 : this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
228 : display a reason if it fails. The latter can be set in the command line with
229 : -eps_converged_reason ::failed
230 :
231 : Level: intermediate
232 :
233 : .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber(), EPSConvergedReasonViewFromOptions()
234 : @*/
235 2 : PetscErrorCode EPSConvergedReasonView(EPS eps,PetscViewer viewer)
236 : {
237 2 : PetscBool isAscii;
238 2 : PetscViewerFormat format;
239 :
240 2 : PetscFunctionBegin;
241 2 : if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
242 2 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
243 2 : if (isAscii) {
244 2 : PetscCall(PetscViewerGetFormat(viewer,&format));
245 2 : PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
246 2 : if (eps->reason > 0 && format != PETSC_VIEWER_FAILED) PetscCall(PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve converged (%" PetscInt_FMT " eigenpair%s) due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",eps->nconv,(eps->nconv>1)?"s":"",EPSConvergedReasons[eps->reason],eps->its));
247 0 : else if (eps->reason <= 0) PetscCall(PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve did not converge due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",EPSConvergedReasons[eps->reason],eps->its));
248 2 : PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
249 : }
250 2 : PetscFunctionReturn(PETSC_SUCCESS);
251 : }
252 :
253 : /*@
254 : EPSConvergedReasonViewFromOptions - Processes command line options to determine if/how
255 : the EPS converged reason is to be viewed.
256 :
257 : Collective
258 :
259 : Input Parameter:
260 : . eps - the eigensolver context
261 :
262 : Level: developer
263 :
264 : .seealso: EPSConvergedReasonView()
265 : @*/
266 855 : PetscErrorCode EPSConvergedReasonViewFromOptions(EPS eps)
267 : {
268 855 : PetscViewer viewer;
269 855 : PetscBool flg;
270 855 : static PetscBool incall = PETSC_FALSE;
271 855 : PetscViewerFormat format;
272 :
273 855 : PetscFunctionBegin;
274 855 : if (incall) PetscFunctionReturn(PETSC_SUCCESS);
275 855 : incall = PETSC_TRUE;
276 855 : PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg));
277 855 : if (flg) {
278 1 : PetscCall(PetscViewerPushFormat(viewer,format));
279 1 : PetscCall(EPSConvergedReasonView(eps,viewer));
280 1 : PetscCall(PetscViewerPopFormat(viewer));
281 1 : PetscCall(PetscOptionsRestoreViewer(&viewer));
282 : }
283 855 : incall = PETSC_FALSE;
284 855 : PetscFunctionReturn(PETSC_SUCCESS);
285 : }
286 :
287 481 : static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
288 : {
289 481 : PetscReal error;
290 481 : PetscInt i,j,k,nvals;
291 :
292 481 : PetscFunctionBegin;
293 481 : nvals = (eps->which==EPS_ALL)? eps->nconv: eps->nev;
294 481 : if (eps->which!=EPS_ALL && eps->nconv<nvals) {
295 0 : PetscCall(PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " eigenvalues converged\n\n",eps->nev));
296 0 : PetscFunctionReturn(PETSC_SUCCESS);
297 : }
298 481 : if (eps->which==EPS_ALL && !nvals) {
299 0 : PetscCall(PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n"));
300 0 : PetscFunctionReturn(PETSC_SUCCESS);
301 : }
302 4132 : for (i=0;i<nvals;i++) {
303 3651 : PetscCall(EPSComputeError(eps,i,etype,&error));
304 3651 : if (error>=5.0*eps->tol) {
305 0 : PetscCall(PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",nvals));
306 3651 : PetscFunctionReturn(PETSC_SUCCESS);
307 : }
308 : }
309 481 : if (eps->which==EPS_ALL) PetscCall(PetscViewerASCIIPrintf(viewer," Found %" PetscInt_FMT " eigenvalues, all of them computed up to the required tolerance:",nvals));
310 415 : else PetscCall(PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:"));
311 1189 : for (i=0;i<=(nvals-1)/8;i++) {
312 708 : PetscCall(PetscViewerASCIIPrintf(viewer,"\n "));
313 4359 : for (j=0;j<PetscMin(8,nvals-8*i);j++) {
314 3651 : k = eps->perm[8*i+j];
315 3651 : PetscCall(SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]));
316 3651 : if (8*i+j+1<nvals) PetscCall(PetscViewerASCIIPrintf(viewer,", "));
317 : }
318 : }
319 481 : PetscCall(PetscViewerASCIIPrintf(viewer,"\n\n"));
320 481 : PetscFunctionReturn(PETSC_SUCCESS);
321 : }
322 :
323 1 : static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
324 : {
325 1 : PetscReal error,re,im;
326 1 : PetscScalar kr,ki;
327 1 : PetscInt i;
328 1 : char ex[30],sep[]=" ---------------------- --------------------\n";
329 :
330 1 : PetscFunctionBegin;
331 1 : if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
332 1 : switch (etype) {
333 0 : case EPS_ERROR_ABSOLUTE:
334 0 : PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||Ax-k%sx||",eps->isgeneralized?"B":""));
335 : break;
336 0 : case EPS_ERROR_RELATIVE:
337 0 : PetscCall(PetscSNPrintf(ex,sizeof(ex),"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":""));
338 : break;
339 1 : case EPS_ERROR_BACKWARD:
340 1 : PetscCall(PetscSNPrintf(ex,sizeof(ex)," eta(x,k)"));
341 : break;
342 : }
343 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep));
344 7 : for (i=0;i<eps->nconv;i++) {
345 6 : PetscCall(EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL));
346 6 : PetscCall(EPSComputeError(eps,i,etype,&error));
347 : #if defined(PETSC_USE_COMPLEX)
348 : re = PetscRealPart(kr);
349 : im = PetscImaginaryPart(kr);
350 : #else
351 6 : re = kr;
352 6 : im = ki;
353 : #endif
354 6 : if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error));
355 6 : else PetscCall(PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error));
356 : }
357 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"%s",sep));
358 1 : PetscFunctionReturn(PETSC_SUCCESS);
359 : }
360 :
361 1 : static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
362 : {
363 1 : PetscReal error;
364 1 : PetscInt i;
365 1 : const char *name;
366 :
367 1 : PetscFunctionBegin;
368 1 : PetscCall(PetscObjectGetName((PetscObject)eps,&name));
369 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name));
370 7 : for (i=0;i<eps->nconv;i++) {
371 6 : PetscCall(EPSComputeError(eps,i,etype,&error));
372 6 : PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error));
373 : }
374 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
375 1 : PetscFunctionReturn(PETSC_SUCCESS);
376 : }
377 :
378 : /*@C
379 : EPSErrorView - Displays the errors associated with the computed solution
380 : (as well as the eigenvalues).
381 :
382 : Collective
383 :
384 : Input Parameters:
385 : + eps - the eigensolver context
386 : . etype - error type
387 : - viewer - optional visualization context
388 :
389 : Options Database Keys:
390 : + -eps_error_absolute - print absolute errors of each eigenpair
391 : . -eps_error_relative - print relative errors of each eigenpair
392 : - -eps_error_backward - print backward errors of each eigenpair
393 :
394 : Notes:
395 : By default, this function checks the error of all eigenpairs and prints
396 : the eigenvalues if all of them are below the requested tolerance.
397 : If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
398 : eigenvalues and corresponding errors is printed.
399 :
400 : Level: intermediate
401 :
402 : .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
403 : @*/
404 483 : PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
405 : {
406 483 : PetscBool isascii;
407 483 : PetscViewerFormat format;
408 :
409 483 : PetscFunctionBegin;
410 483 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
411 483 : if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
412 483 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,3);
413 483 : PetscCheckSameComm(eps,1,viewer,3);
414 483 : EPSCheckSolved(eps,1);
415 483 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
416 483 : if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);
417 :
418 483 : PetscCall(PetscViewerGetFormat(viewer,&format));
419 483 : switch (format) {
420 481 : case PETSC_VIEWER_DEFAULT:
421 : case PETSC_VIEWER_ASCII_INFO:
422 481 : PetscCall(EPSErrorView_ASCII(eps,etype,viewer));
423 : break;
424 1 : case PETSC_VIEWER_ASCII_INFO_DETAIL:
425 1 : PetscCall(EPSErrorView_DETAIL(eps,etype,viewer));
426 : break;
427 1 : case PETSC_VIEWER_ASCII_MATLAB:
428 1 : PetscCall(EPSErrorView_MATLAB(eps,etype,viewer));
429 : break;
430 0 : default:
431 0 : PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
432 : }
433 483 : PetscFunctionReturn(PETSC_SUCCESS);
434 : }
435 :
436 : /*@
437 : EPSErrorViewFromOptions - Processes command line options to determine if/how
438 : the errors of the computed solution are to be viewed.
439 :
440 : Collective
441 :
442 : Input Parameter:
443 : . eps - the eigensolver context
444 :
445 : Level: developer
446 :
447 : .seealso: EPSErrorView()
448 : @*/
449 855 : PetscErrorCode EPSErrorViewFromOptions(EPS eps)
450 : {
451 855 : PetscViewer viewer;
452 855 : PetscBool flg;
453 855 : static PetscBool incall = PETSC_FALSE;
454 855 : PetscViewerFormat format;
455 :
456 855 : PetscFunctionBegin;
457 855 : if (incall) PetscFunctionReturn(PETSC_SUCCESS);
458 855 : incall = PETSC_TRUE;
459 855 : PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg));
460 855 : if (flg) {
461 1 : PetscCall(PetscViewerPushFormat(viewer,format));
462 1 : PetscCall(EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer));
463 1 : PetscCall(PetscViewerPopFormat(viewer));
464 1 : PetscCall(PetscOptionsRestoreViewer(&viewer));
465 : }
466 855 : PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg));
467 855 : if (flg) {
468 2 : PetscCall(PetscViewerPushFormat(viewer,format));
469 2 : PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer));
470 2 : PetscCall(PetscViewerPopFormat(viewer));
471 2 : PetscCall(PetscOptionsRestoreViewer(&viewer));
472 : }
473 855 : PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg));
474 855 : if (flg) {
475 1 : PetscCall(PetscViewerPushFormat(viewer,format));
476 1 : PetscCall(EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer));
477 1 : PetscCall(PetscViewerPopFormat(viewer));
478 1 : PetscCall(PetscOptionsRestoreViewer(&viewer));
479 : }
480 855 : incall = PETSC_FALSE;
481 855 : PetscFunctionReturn(PETSC_SUCCESS);
482 : }
483 :
484 0 : static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
485 : {
486 0 : PetscDraw draw;
487 0 : PetscDrawSP drawsp;
488 0 : PetscReal re,im;
489 0 : PetscInt i,k;
490 :
491 0 : PetscFunctionBegin;
492 0 : if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
493 0 : PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
494 0 : PetscCall(PetscDrawSetTitle(draw,"Computed Eigenvalues"));
495 0 : PetscCall(PetscDrawSPCreate(draw,1,&drawsp));
496 0 : for (i=0;i<eps->nconv;i++) {
497 0 : k = eps->perm[i];
498 : #if defined(PETSC_USE_COMPLEX)
499 : re = PetscRealPart(eps->eigr[k]);
500 : im = PetscImaginaryPart(eps->eigr[k]);
501 : #else
502 0 : re = eps->eigr[k];
503 0 : im = eps->eigi[k];
504 : #endif
505 0 : PetscCall(PetscDrawSPAddPoint(drawsp,&re,&im));
506 : }
507 0 : PetscCall(PetscDrawSPDraw(drawsp,PETSC_TRUE));
508 0 : PetscCall(PetscDrawSPSave(drawsp));
509 0 : PetscCall(PetscDrawSPDestroy(&drawsp));
510 0 : PetscFunctionReturn(PETSC_SUCCESS);
511 : }
512 :
513 1 : static PetscErrorCode EPSValuesView_BINARY(EPS eps,PetscViewer viewer)
514 : {
515 : #if defined(PETSC_HAVE_COMPLEX)
516 1 : PetscInt i,k;
517 1 : PetscComplex *ev;
518 : #endif
519 :
520 1 : PetscFunctionBegin;
521 : #if defined(PETSC_HAVE_COMPLEX)
522 1 : PetscCall(PetscMalloc1(eps->nconv,&ev));
523 5 : for (i=0;i<eps->nconv;i++) {
524 4 : k = eps->perm[i];
525 : #if defined(PETSC_USE_COMPLEX)
526 : ev[i] = eps->eigr[k];
527 : #else
528 4 : ev[i] = PetscCMPLX(eps->eigr[k],eps->eigi[k]);
529 : #endif
530 : }
531 1 : PetscCall(PetscViewerBinaryWrite(viewer,ev,eps->nconv,PETSC_COMPLEX));
532 1 : PetscCall(PetscFree(ev));
533 : #endif
534 1 : PetscFunctionReturn(PETSC_SUCCESS);
535 : }
536 :
537 : #if defined(PETSC_HAVE_HDF5)
538 : static PetscErrorCode EPSValuesView_HDF5(EPS eps,PetscViewer viewer)
539 : {
540 : PetscInt i,k,n,N;
541 : PetscMPIInt rank;
542 : Vec v;
543 : char vname[30];
544 : const char *ename;
545 :
546 : PetscFunctionBegin;
547 : PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank));
548 : N = eps->nconv;
549 : n = rank? 0: N;
550 : /* create a vector containing the eigenvalues */
551 : PetscCall(VecCreateMPI(PetscObjectComm((PetscObject)eps),n,N,&v));
552 : PetscCall(PetscObjectGetName((PetscObject)eps,&ename));
553 : PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename));
554 : PetscCall(PetscObjectSetName((PetscObject)v,vname));
555 : if (!rank) {
556 : for (i=0;i<eps->nconv;i++) {
557 : k = eps->perm[i];
558 : PetscCall(VecSetValue(v,i,eps->eigr[k],INSERT_VALUES));
559 : }
560 : }
561 : PetscCall(VecAssemblyBegin(v));
562 : PetscCall(VecAssemblyEnd(v));
563 : PetscCall(VecView(v,viewer));
564 : #if !defined(PETSC_USE_COMPLEX)
565 : /* in real scalars write the imaginary part as a separate vector */
566 : PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename));
567 : PetscCall(PetscObjectSetName((PetscObject)v,vname));
568 : if (!rank) {
569 : for (i=0;i<eps->nconv;i++) {
570 : k = eps->perm[i];
571 : PetscCall(VecSetValue(v,i,eps->eigi[k],INSERT_VALUES));
572 : }
573 : }
574 : PetscCall(VecAssemblyBegin(v));
575 : PetscCall(VecAssemblyEnd(v));
576 : PetscCall(VecView(v,viewer));
577 : #endif
578 : PetscCall(VecDestroy(&v));
579 : PetscFunctionReturn(PETSC_SUCCESS);
580 : }
581 : #endif
582 :
583 2 : static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
584 : {
585 2 : PetscInt i,k;
586 :
587 2 : PetscFunctionBegin;
588 2 : PetscCall(PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n"));
589 12 : for (i=0;i<eps->nconv;i++) {
590 10 : k = eps->perm[i];
591 10 : PetscCall(PetscViewerASCIIPrintf(viewer," "));
592 10 : PetscCall(SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]));
593 10 : PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
594 : }
595 2 : PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
596 2 : PetscFunctionReturn(PETSC_SUCCESS);
597 : }
598 :
599 1 : static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
600 : {
601 1 : PetscInt i,k;
602 1 : PetscReal re,im;
603 1 : const char *name;
604 :
605 1 : PetscFunctionBegin;
606 1 : PetscCall(PetscObjectGetName((PetscObject)eps,&name));
607 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name));
608 2 : for (i=0;i<eps->nconv;i++) {
609 1 : k = eps->perm[i];
610 : #if defined(PETSC_USE_COMPLEX)
611 : re = PetscRealPart(eps->eigr[k]);
612 : im = PetscImaginaryPart(eps->eigr[k]);
613 : #else
614 1 : re = eps->eigr[k];
615 1 : im = eps->eigi[k];
616 : #endif
617 1 : if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im));
618 1 : else PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re));
619 : }
620 1 : PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
621 1 : PetscFunctionReturn(PETSC_SUCCESS);
622 : }
623 :
624 : /*@C
625 : EPSValuesView - Displays the computed eigenvalues in a viewer.
626 :
627 : Collective
628 :
629 : Input Parameters:
630 : + eps - the eigensolver context
631 : - viewer - the viewer
632 :
633 : Options Database Key:
634 : . -eps_view_values - print computed eigenvalues
635 :
636 : Level: intermediate
637 :
638 : .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
639 : @*/
640 4 : PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
641 : {
642 4 : PetscBool isascii,isdraw,isbinary;
643 4 : PetscViewerFormat format;
644 : #if defined(PETSC_HAVE_HDF5)
645 : PetscBool ishdf5;
646 : #endif
647 :
648 4 : PetscFunctionBegin;
649 4 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
650 4 : if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
651 4 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2);
652 4 : PetscCheckSameComm(eps,1,viewer,2);
653 4 : EPSCheckSolved(eps,1);
654 4 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
655 4 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
656 : #if defined(PETSC_HAVE_HDF5)
657 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
658 : #endif
659 4 : PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
660 4 : if (isdraw) PetscCall(EPSValuesView_DRAW(eps,viewer));
661 4 : else if (isbinary) PetscCall(EPSValuesView_BINARY(eps,viewer));
662 : #if defined(PETSC_HAVE_HDF5)
663 : else if (ishdf5) PetscCall(EPSValuesView_HDF5(eps,viewer));
664 : #endif
665 3 : else if (isascii) {
666 3 : PetscCall(PetscViewerGetFormat(viewer,&format));
667 3 : switch (format) {
668 2 : case PETSC_VIEWER_DEFAULT:
669 : case PETSC_VIEWER_ASCII_INFO:
670 : case PETSC_VIEWER_ASCII_INFO_DETAIL:
671 2 : PetscCall(EPSValuesView_ASCII(eps,viewer));
672 : break;
673 1 : case PETSC_VIEWER_ASCII_MATLAB:
674 1 : PetscCall(EPSValuesView_MATLAB(eps,viewer));
675 : break;
676 0 : default:
677 0 : PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
678 : }
679 : }
680 4 : PetscFunctionReturn(PETSC_SUCCESS);
681 : }
682 :
683 : /*@
684 : EPSValuesViewFromOptions - Processes command line options to determine if/how
685 : the computed eigenvalues are to be viewed.
686 :
687 : Collective
688 :
689 : Input Parameters:
690 : . eps - the eigensolver context
691 :
692 : Level: developer
693 :
694 : .seealso: EPSValuesView()
695 : @*/
696 855 : PetscErrorCode EPSValuesViewFromOptions(EPS eps)
697 : {
698 855 : PetscViewer viewer;
699 855 : PetscBool flg;
700 855 : static PetscBool incall = PETSC_FALSE;
701 855 : PetscViewerFormat format;
702 :
703 855 : PetscFunctionBegin;
704 855 : if (incall) PetscFunctionReturn(PETSC_SUCCESS);
705 855 : incall = PETSC_TRUE;
706 855 : PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg));
707 855 : if (flg) {
708 4 : PetscCall(PetscViewerPushFormat(viewer,format));
709 4 : PetscCall(EPSValuesView(eps,viewer));
710 4 : PetscCall(PetscViewerPopFormat(viewer));
711 4 : PetscCall(PetscOptionsRestoreViewer(&viewer));
712 : }
713 855 : incall = PETSC_FALSE;
714 855 : PetscFunctionReturn(PETSC_SUCCESS);
715 : }
716 :
717 : /*@C
718 : EPSVectorsView - Outputs computed eigenvectors to a viewer.
719 :
720 : Collective
721 :
722 : Input Parameters:
723 : + eps - the eigensolver context
724 : - viewer - the viewer
725 :
726 : Options Database Key:
727 : . -eps_view_vectors - output eigenvectors.
728 :
729 : Notes:
730 : If PETSc was configured with real scalars, complex conjugate eigenvectors
731 : will be viewed as two separate real vectors, one containing the real part
732 : and another one containing the imaginary part.
733 :
734 : If left eigenvectors were computed with a two-sided eigensolver, the right
735 : and left eigenvectors are interleaved, that is, the vectors are output in
736 : the following order X0, Y0, X1, Y1, X2, Y2, ...
737 :
738 : Level: intermediate
739 :
740 : .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
741 : @*/
742 1 : PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
743 : {
744 1 : PetscInt i,k;
745 1 : Vec xr,xi=NULL;
746 :
747 1 : PetscFunctionBegin;
748 1 : PetscValidHeaderSpecific(eps,EPS_CLASSID,1);
749 1 : if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
750 1 : PetscValidHeaderSpecific(viewer,PETSC_VIEWER_CLASSID,2);
751 1 : PetscCheckSameComm(eps,1,viewer,2);
752 1 : EPSCheckSolved(eps,1);
753 1 : if (eps->nconv) {
754 1 : PetscCall(EPSComputeVectors(eps));
755 1 : PetscCall(BVCreateVec(eps->V,&xr));
756 : #if !defined(PETSC_USE_COMPLEX)
757 1 : PetscCall(BVCreateVec(eps->V,&xi));
758 : #endif
759 5 : for (i=0;i<eps->nconv;i++) {
760 4 : k = eps->perm[i];
761 4 : PetscCall(BV_GetEigenvector(eps->V,k,eps->eigi[k],xr,xi));
762 4 : PetscCall(SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)eps));
763 4 : if (eps->twosided) {
764 4 : PetscCall(BV_GetEigenvector(eps->W,k,eps->eigi[k],xr,xi));
765 4 : PetscCall(SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)eps));
766 : }
767 : }
768 1 : PetscCall(VecDestroy(&xr));
769 : #if !defined(PETSC_USE_COMPLEX)
770 1 : PetscCall(VecDestroy(&xi));
771 : #endif
772 : }
773 1 : PetscFunctionReturn(PETSC_SUCCESS);
774 : }
775 :
776 : /*@
777 : EPSVectorsViewFromOptions - Processes command line options to determine if/how
778 : the computed eigenvectors are to be viewed.
779 :
780 : Collective
781 :
782 : Input Parameter:
783 : . eps - the eigensolver context
784 :
785 : Level: developer
786 :
787 : .seealso: EPSVectorsView()
788 : @*/
789 855 : PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
790 : {
791 855 : PetscViewer viewer;
792 855 : PetscBool flg = PETSC_FALSE;
793 855 : static PetscBool incall = PETSC_FALSE;
794 855 : PetscViewerFormat format;
795 :
796 855 : PetscFunctionBegin;
797 855 : if (incall) PetscFunctionReturn(PETSC_SUCCESS);
798 855 : incall = PETSC_TRUE;
799 855 : PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg));
800 855 : if (flg) {
801 1 : PetscCall(PetscViewerPushFormat(viewer,format));
802 1 : PetscCall(EPSVectorsView(eps,viewer));
803 1 : PetscCall(PetscViewerPopFormat(viewer));
804 1 : PetscCall(PetscOptionsRestoreViewer(&viewer));
805 : }
806 855 : incall = PETSC_FALSE;
807 855 : PetscFunctionReturn(PETSC_SUCCESS);
808 : }
|