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