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