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