Actual source code: epsview.c
slepc-3.21.0 2024-03-30
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: /*@C
19: EPSView - Prints the EPS data structure.
21: Collective
23: Input Parameters:
24: + eps - the eigenproblem solver context
25: - viewer - optional visualization context
27: Options Database Key:
28: . -eps_view - Calls EPSView() at end of EPSSolve()
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.
38: The user can open an alternative visualization context with
39: PetscViewerASCIIOpen() - output to a specified file.
41: Level: beginner
43: .seealso: STView()
44: @*/
45: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
46: {
47: const char *type=NULL,*extr=NULL,*bal=NULL;
48: char str[50];
49: PetscBool isascii,isexternal,istrivial;
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: }
72: } else type = "not yet set";
73: PetscCall(PetscViewerASCIIPrintf(viewer," problem type: %s\n",type));
74: if (eps->extraction) {
75: switch (eps->extraction) {
76: case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
77: case EPS_HARMONIC: extr = "harmonic Ritz"; break;
78: case EPS_HARMONIC_RELATIVE: extr = "relative harmonic Ritz"; break;
79: case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
80: case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
81: case EPS_REFINED: extr = "refined Ritz"; break;
82: case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
83: }
84: PetscCall(PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr));
85: }
86: if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
87: switch (eps->balance) {
88: case EPS_BALANCE_NONE: break;
89: case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
90: case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
91: case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
92: }
93: PetscCall(PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal));
94: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
95: if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) PetscCall(PetscViewerASCIIPrintf(viewer,", with its=%" PetscInt_FMT,eps->balance_its));
96: if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff));
97: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
98: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
99: }
100: PetscCall(PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: "));
101: PetscCall(SlepcSNPrintfScalar(str,sizeof(str),eps->target,PETSC_FALSE));
102: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
103: if (!eps->which) PetscCall(PetscViewerASCIIPrintf(viewer,"not yet set\n"));
104: else switch (eps->which) {
105: case EPS_WHICH_USER:
106: PetscCall(PetscViewerASCIIPrintf(viewer,"user defined\n"));
107: break;
108: case EPS_TARGET_MAGNITUDE:
109: PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str));
110: break;
111: case EPS_TARGET_REAL:
112: PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str));
113: break;
114: case EPS_TARGET_IMAGINARY:
115: PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str));
116: break;
117: case EPS_LARGEST_MAGNITUDE:
118: PetscCall(PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n"));
119: break;
120: case EPS_SMALLEST_MAGNITUDE:
121: PetscCall(PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n"));
122: break;
123: case EPS_LARGEST_REAL:
124: PetscCall(PetscViewerASCIIPrintf(viewer,"largest real parts\n"));
125: break;
126: case EPS_SMALLEST_REAL:
127: PetscCall(PetscViewerASCIIPrintf(viewer,"smallest real parts\n"));
128: break;
129: case EPS_LARGEST_IMAGINARY:
130: PetscCall(PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n"));
131: break;
132: case EPS_SMALLEST_IMAGINARY:
133: PetscCall(PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n"));
134: break;
135: case EPS_ALL:
136: if (eps->inta || eps->intb) PetscCall(PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb));
137: else PetscCall(PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n"));
138: break;
139: }
140: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
141: 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: if (eps->purify) PetscCall(PetscViewerASCIIPrintf(viewer," postprocessing eigenvectors with purification\n"));
143: if (eps->trueres) PetscCall(PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n"));
144: if (eps->trackall) PetscCall(PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n"));
145: PetscCall(PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %" PetscInt_FMT "\n",eps->nev));
146: PetscCall(PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",eps->ncv));
147: PetscCall(PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",eps->mpd));
148: PetscCall(PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",eps->max_it));
149: PetscCall(PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)eps->tol));
150: PetscCall(PetscViewerASCIIPrintf(viewer," convergence test: "));
151: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
152: switch (eps->conv) {
153: case EPS_CONV_ABS:
154: PetscCall(PetscViewerASCIIPrintf(viewer,"absolute\n"));break;
155: case EPS_CONV_REL:
156: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n"));break;
157: case EPS_CONV_NORM:
158: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n"));
159: PetscCall(PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)eps->nrma));
160: if (eps->isgeneralized) PetscCall(PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb));
161: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
162: break;
163: case EPS_CONV_USER:
164: PetscCall(PetscViewerASCIIPrintf(viewer,"user-defined\n"));break;
165: }
166: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
167: if (eps->nini) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(eps->nini)));
168: if (eps->ninil) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided left initial space: %" PetscInt_FMT "\n",PetscAbs(eps->ninil)));
169: if (eps->nds) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %" PetscInt_FMT "\n",PetscAbs(eps->nds)));
170: } else PetscTryTypeMethod(eps,view,viewer);
171: PetscCall(PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSELEMENTAL,EPSFEAST,EPSPRIMME,EPSSCALAPACK,EPSELPA,EPSEVSL,EPSTRLAN,""));
172: if (!isexternal) {
173: PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
174: if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
175: PetscCall(BVView(eps->V,viewer));
176: if (eps->rg) {
177: PetscCall(RGIsTrivial(eps->rg,&istrivial));
178: if (!istrivial) PetscCall(RGView(eps->rg,viewer));
179: }
180: if (eps->useds) {
181: if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
182: PetscCall(DSView(eps->ds,viewer));
183: }
184: PetscCall(PetscViewerPopFormat(viewer));
185: }
186: if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
187: PetscCall(STView(eps->st,viewer));
188: PetscFunctionReturn(PETSC_SUCCESS);
189: }
191: /*@C
192: EPSViewFromOptions - View from options
194: Collective
196: Input Parameters:
197: + eps - the eigensolver context
198: . obj - optional object
199: - name - command line option
201: Level: intermediate
203: .seealso: EPSView(), EPSCreate()
204: @*/
205: PetscErrorCode EPSViewFromOptions(EPS eps,PetscObject obj,const char name[])
206: {
207: PetscFunctionBegin;
209: PetscCall(PetscObjectViewFromOptions((PetscObject)eps,obj,name));
210: PetscFunctionReturn(PETSC_SUCCESS);
211: }
213: /*@C
214: EPSConvergedReasonView - Displays the reason an EPS solve converged or diverged.
216: Collective
218: Input Parameters:
219: + eps - the eigensolver context
220: - viewer - the viewer to display the reason
222: Options Database Keys:
223: . -eps_converged_reason - print reason for convergence, and number of iterations
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
231: Level: intermediate
233: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber(), EPSConvergedReasonViewFromOptions()
234: @*/
235: PetscErrorCode EPSConvergedReasonView(EPS eps,PetscViewer viewer)
236: {
237: PetscBool isAscii;
238: PetscViewerFormat format;
240: PetscFunctionBegin;
241: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
242: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
243: if (isAscii) {
244: PetscCall(PetscViewerGetFormat(viewer,&format));
245: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
246: 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: 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: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
249: }
250: PetscFunctionReturn(PETSC_SUCCESS);
251: }
253: /*@
254: EPSConvergedReasonViewFromOptions - Processes command line options to determine if/how
255: the EPS converged reason is to be viewed.
257: Collective
259: Input Parameter:
260: . eps - the eigensolver context
262: Level: developer
264: .seealso: EPSConvergedReasonView()
265: @*/
266: PetscErrorCode EPSConvergedReasonViewFromOptions(EPS eps)
267: {
268: PetscViewer viewer;
269: PetscBool flg;
270: static PetscBool incall = PETSC_FALSE;
271: PetscViewerFormat format;
273: PetscFunctionBegin;
274: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
275: incall = PETSC_TRUE;
276: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg));
277: if (flg) {
278: PetscCall(PetscViewerPushFormat(viewer,format));
279: PetscCall(EPSConvergedReasonView(eps,viewer));
280: PetscCall(PetscViewerPopFormat(viewer));
281: PetscCall(PetscOptionsRestoreViewer(&viewer));
282: }
283: incall = PETSC_FALSE;
284: PetscFunctionReturn(PETSC_SUCCESS);
285: }
287: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
288: {
289: PetscReal error;
290: PetscInt i,j,k,nvals;
292: PetscFunctionBegin;
293: nvals = (eps->which==EPS_ALL)? eps->nconv: eps->nev;
294: if (eps->which!=EPS_ALL && eps->nconv<nvals) {
295: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " eigenvalues converged\n\n",eps->nev));
296: PetscFunctionReturn(PETSC_SUCCESS);
297: }
298: if (eps->which==EPS_ALL && !nvals) {
299: PetscCall(PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n"));
300: PetscFunctionReturn(PETSC_SUCCESS);
301: }
302: for (i=0;i<nvals;i++) {
303: PetscCall(EPSComputeError(eps,i,etype,&error));
304: if (error>=5.0*eps->tol) {
305: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",nvals));
306: PetscFunctionReturn(PETSC_SUCCESS);
307: }
308: }
309: if (eps->which==EPS_ALL) PetscCall(PetscViewerASCIIPrintf(viewer," Found %" PetscInt_FMT " eigenvalues, all of them computed up to the required tolerance:",nvals));
310: else PetscCall(PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:"));
311: for (i=0;i<=(nvals-1)/8;i++) {
312: PetscCall(PetscViewerASCIIPrintf(viewer,"\n "));
313: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
314: k = eps->perm[8*i+j];
315: PetscCall(SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]));
316: if (8*i+j+1<nvals) PetscCall(PetscViewerASCIIPrintf(viewer,", "));
317: }
318: }
319: PetscCall(PetscViewerASCIIPrintf(viewer,"\n\n"));
320: PetscFunctionReturn(PETSC_SUCCESS);
321: }
323: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
324: {
325: PetscReal error,re,im;
326: PetscScalar kr,ki;
327: PetscInt i;
328: char ex[30],sep[]=" ---------------------- --------------------\n";
330: PetscFunctionBegin;
331: if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
332: switch (etype) {
333: case EPS_ERROR_ABSOLUTE:
334: PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||Ax-k%sx||",eps->isgeneralized?"B":""));
335: break;
336: case EPS_ERROR_RELATIVE:
337: PetscCall(PetscSNPrintf(ex,sizeof(ex),"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":""));
338: break;
339: case EPS_ERROR_BACKWARD:
340: PetscCall(PetscSNPrintf(ex,sizeof(ex)," eta(x,k)"));
341: break;
342: }
343: PetscCall(PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep));
344: for (i=0;i<eps->nconv;i++) {
345: PetscCall(EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL));
346: PetscCall(EPSComputeError(eps,i,etype,&error));
347: #if defined(PETSC_USE_COMPLEX)
348: re = PetscRealPart(kr);
349: im = PetscImaginaryPart(kr);
350: #else
351: re = kr;
352: im = ki;
353: #endif
354: if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error));
355: else PetscCall(PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error));
356: }
357: PetscCall(PetscViewerASCIIPrintf(viewer,"%s",sep));
358: PetscFunctionReturn(PETSC_SUCCESS);
359: }
361: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
362: {
363: PetscReal error;
364: PetscInt i;
365: const char *name;
367: PetscFunctionBegin;
368: PetscCall(PetscObjectGetName((PetscObject)eps,&name));
369: PetscCall(PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name));
370: for (i=0;i<eps->nconv;i++) {
371: PetscCall(EPSComputeError(eps,i,etype,&error));
372: PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error));
373: }
374: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
375: PetscFunctionReturn(PETSC_SUCCESS);
376: }
378: /*@C
379: EPSErrorView - Displays the errors associated with the computed solution
380: (as well as the eigenvalues).
382: Collective
384: Input Parameters:
385: + eps - the eigensolver context
386: . etype - error type
387: - viewer - optional visualization context
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
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.
400: Level: intermediate
402: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
403: @*/
404: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
405: {
406: PetscBool isascii;
407: PetscViewerFormat format;
409: PetscFunctionBegin;
411: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
413: PetscCheckSameComm(eps,1,viewer,3);
414: EPSCheckSolved(eps,1);
415: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
416: if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);
418: PetscCall(PetscViewerGetFormat(viewer,&format));
419: switch (format) {
420: case PETSC_VIEWER_DEFAULT:
421: case PETSC_VIEWER_ASCII_INFO:
422: PetscCall(EPSErrorView_ASCII(eps,etype,viewer));
423: break;
424: case PETSC_VIEWER_ASCII_INFO_DETAIL:
425: PetscCall(EPSErrorView_DETAIL(eps,etype,viewer));
426: break;
427: case PETSC_VIEWER_ASCII_MATLAB:
428: PetscCall(EPSErrorView_MATLAB(eps,etype,viewer));
429: break;
430: default:
431: PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
432: }
433: PetscFunctionReturn(PETSC_SUCCESS);
434: }
436: /*@
437: EPSErrorViewFromOptions - Processes command line options to determine if/how
438: the errors of the computed solution are to be viewed.
440: Collective
442: Input Parameter:
443: . eps - the eigensolver context
445: Level: developer
447: .seealso: EPSErrorView()
448: @*/
449: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
450: {
451: PetscViewer viewer;
452: PetscBool flg;
453: static PetscBool incall = PETSC_FALSE;
454: PetscViewerFormat format;
456: PetscFunctionBegin;
457: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
458: incall = PETSC_TRUE;
459: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg));
460: if (flg) {
461: PetscCall(PetscViewerPushFormat(viewer,format));
462: PetscCall(EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer));
463: PetscCall(PetscViewerPopFormat(viewer));
464: PetscCall(PetscOptionsRestoreViewer(&viewer));
465: }
466: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg));
467: if (flg) {
468: PetscCall(PetscViewerPushFormat(viewer,format));
469: PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer));
470: PetscCall(PetscViewerPopFormat(viewer));
471: PetscCall(PetscOptionsRestoreViewer(&viewer));
472: }
473: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg));
474: if (flg) {
475: PetscCall(PetscViewerPushFormat(viewer,format));
476: PetscCall(EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer));
477: PetscCall(PetscViewerPopFormat(viewer));
478: PetscCall(PetscOptionsRestoreViewer(&viewer));
479: }
480: incall = PETSC_FALSE;
481: PetscFunctionReturn(PETSC_SUCCESS);
482: }
484: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
485: {
486: PetscDraw draw;
487: PetscDrawSP drawsp;
488: PetscReal re,im;
489: PetscInt i,k;
491: PetscFunctionBegin;
492: if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
493: PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
494: PetscCall(PetscDrawSetTitle(draw,"Computed Eigenvalues"));
495: PetscCall(PetscDrawSPCreate(draw,1,&drawsp));
496: for (i=0;i<eps->nconv;i++) {
497: 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: re = eps->eigr[k];
503: im = eps->eigi[k];
504: #endif
505: PetscCall(PetscDrawSPAddPoint(drawsp,&re,&im));
506: }
507: PetscCall(PetscDrawSPDraw(drawsp,PETSC_TRUE));
508: PetscCall(PetscDrawSPSave(drawsp));
509: PetscCall(PetscDrawSPDestroy(&drawsp));
510: PetscFunctionReturn(PETSC_SUCCESS);
511: }
513: static PetscErrorCode EPSValuesView_BINARY(EPS eps,PetscViewer viewer)
514: {
515: #if defined(PETSC_HAVE_COMPLEX)
516: PetscInt i,k;
517: PetscComplex *ev;
518: #endif
520: PetscFunctionBegin;
521: #if defined(PETSC_HAVE_COMPLEX)
522: PetscCall(PetscMalloc1(eps->nconv,&ev));
523: for (i=0;i<eps->nconv;i++) {
524: k = eps->perm[i];
525: #if defined(PETSC_USE_COMPLEX)
526: ev[i] = eps->eigr[k];
527: #else
528: ev[i] = PetscCMPLX(eps->eigr[k],eps->eigi[k]);
529: #endif
530: }
531: PetscCall(PetscViewerBinaryWrite(viewer,ev,eps->nconv,PETSC_COMPLEX));
532: PetscCall(PetscFree(ev));
533: #endif
534: PetscFunctionReturn(PETSC_SUCCESS);
535: }
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;
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
583: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
584: {
585: PetscInt i,k;
587: PetscFunctionBegin;
588: PetscCall(PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n"));
589: for (i=0;i<eps->nconv;i++) {
590: k = eps->perm[i];
591: PetscCall(PetscViewerASCIIPrintf(viewer," "));
592: PetscCall(SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]));
593: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
594: }
595: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
596: PetscFunctionReturn(PETSC_SUCCESS);
597: }
599: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
600: {
601: PetscInt i,k;
602: PetscReal re,im;
603: const char *name;
605: PetscFunctionBegin;
606: PetscCall(PetscObjectGetName((PetscObject)eps,&name));
607: PetscCall(PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name));
608: for (i=0;i<eps->nconv;i++) {
609: 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: re = eps->eigr[k];
615: im = eps->eigi[k];
616: #endif
617: if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im));
618: else PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re));
619: }
620: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
621: PetscFunctionReturn(PETSC_SUCCESS);
622: }
624: /*@C
625: EPSValuesView - Displays the computed eigenvalues in a viewer.
627: Collective
629: Input Parameters:
630: + eps - the eigensolver context
631: - viewer - the viewer
633: Options Database Key:
634: . -eps_view_values - print computed eigenvalues
636: Level: intermediate
638: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
639: @*/
640: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
641: {
642: PetscBool isascii,isdraw,isbinary;
643: PetscViewerFormat format;
644: #if defined(PETSC_HAVE_HDF5)
645: PetscBool ishdf5;
646: #endif
648: PetscFunctionBegin;
650: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
652: PetscCheckSameComm(eps,1,viewer,2);
653: EPSCheckSolved(eps,1);
654: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
655: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
656: #if defined(PETSC_HAVE_HDF5)
657: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
658: #endif
659: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
660: if (isdraw) PetscCall(EPSValuesView_DRAW(eps,viewer));
661: 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: else if (isascii) {
666: PetscCall(PetscViewerGetFormat(viewer,&format));
667: switch (format) {
668: case PETSC_VIEWER_DEFAULT:
669: case PETSC_VIEWER_ASCII_INFO:
670: case PETSC_VIEWER_ASCII_INFO_DETAIL:
671: PetscCall(EPSValuesView_ASCII(eps,viewer));
672: break;
673: case PETSC_VIEWER_ASCII_MATLAB:
674: PetscCall(EPSValuesView_MATLAB(eps,viewer));
675: break;
676: default:
677: PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
678: }
679: }
680: PetscFunctionReturn(PETSC_SUCCESS);
681: }
683: /*@
684: EPSValuesViewFromOptions - Processes command line options to determine if/how
685: the computed eigenvalues are to be viewed.
687: Collective
689: Input Parameters:
690: . eps - the eigensolver context
692: Level: developer
694: .seealso: EPSValuesView()
695: @*/
696: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
697: {
698: PetscViewer viewer;
699: PetscBool flg;
700: static PetscBool incall = PETSC_FALSE;
701: PetscViewerFormat format;
703: PetscFunctionBegin;
704: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
705: incall = PETSC_TRUE;
706: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg));
707: if (flg) {
708: PetscCall(PetscViewerPushFormat(viewer,format));
709: PetscCall(EPSValuesView(eps,viewer));
710: PetscCall(PetscViewerPopFormat(viewer));
711: PetscCall(PetscOptionsRestoreViewer(&viewer));
712: }
713: incall = PETSC_FALSE;
714: PetscFunctionReturn(PETSC_SUCCESS);
715: }
717: /*@C
718: EPSVectorsView - Outputs computed eigenvectors to a viewer.
720: Collective
722: Input Parameters:
723: + eps - the eigensolver context
724: - viewer - the viewer
726: Options Database Key:
727: . -eps_view_vectors - output eigenvectors.
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.
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, ...
738: Level: intermediate
740: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
741: @*/
742: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
743: {
744: PetscInt i,k;
745: Vec xr,xi=NULL;
747: PetscFunctionBegin;
749: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
751: PetscCheckSameComm(eps,1,viewer,2);
752: EPSCheckSolved(eps,1);
753: if (eps->nconv) {
754: PetscCall(EPSComputeVectors(eps));
755: PetscCall(BVCreateVec(eps->V,&xr));
756: #if !defined(PETSC_USE_COMPLEX)
757: PetscCall(BVCreateVec(eps->V,&xi));
758: #endif
759: for (i=0;i<eps->nconv;i++) {
760: k = eps->perm[i];
761: PetscCall(BV_GetEigenvector(eps->V,k,eps->eigi[k],xr,xi));
762: PetscCall(SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)eps));
763: if (eps->twosided) {
764: PetscCall(BV_GetEigenvector(eps->W,k,eps->eigi[k],xr,xi));
765: PetscCall(SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)eps));
766: }
767: }
768: PetscCall(VecDestroy(&xr));
769: #if !defined(PETSC_USE_COMPLEX)
770: PetscCall(VecDestroy(&xi));
771: #endif
772: }
773: PetscFunctionReturn(PETSC_SUCCESS);
774: }
776: /*@
777: EPSVectorsViewFromOptions - Processes command line options to determine if/how
778: the computed eigenvectors are to be viewed.
780: Collective
782: Input Parameter:
783: . eps - the eigensolver context
785: Level: developer
787: .seealso: EPSVectorsView()
788: @*/
789: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
790: {
791: PetscViewer viewer;
792: PetscBool flg = PETSC_FALSE;
793: static PetscBool incall = PETSC_FALSE;
794: PetscViewerFormat format;
796: PetscFunctionBegin;
797: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
798: incall = PETSC_TRUE;
799: PetscCall(PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg));
800: if (flg) {
801: PetscCall(PetscViewerPushFormat(viewer,format));
802: PetscCall(EPSVectorsView(eps,viewer));
803: PetscCall(PetscViewerPopFormat(viewer));
804: PetscCall(PetscOptionsRestoreViewer(&viewer));
805: }
806: incall = PETSC_FALSE;
807: PetscFunctionReturn(PETSC_SUCCESS);
808: }