Actual source code: epsview.c
slepc-3.22.1 2024-10-28
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: PetscCall(PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %" PetscInt_FMT "\n",eps->nev));
154: PetscCall(PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",eps->ncv));
155: PetscCall(PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",eps->mpd));
156: PetscCall(PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",eps->max_it));
157: PetscCall(PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)eps->tol));
158: PetscCall(PetscViewerASCIIPrintf(viewer," convergence test: "));
159: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
160: switch (eps->conv) {
161: case EPS_CONV_ABS:
162: PetscCall(PetscViewerASCIIPrintf(viewer,"absolute\n"));break;
163: case EPS_CONV_REL:
164: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n"));break;
165: case EPS_CONV_NORM:
166: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n"));
167: PetscCall(PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)eps->nrma));
168: if (eps->isgeneralized) PetscCall(PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb));
169: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
170: break;
171: case EPS_CONV_USER:
172: PetscCall(PetscViewerASCIIPrintf(viewer,"user-defined\n"));break;
173: }
174: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
175: if (eps->nini) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(eps->nini)));
176: if (eps->ninil) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided left initial space: %" PetscInt_FMT "\n",PetscAbs(eps->ninil)));
177: if (eps->nds) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %" PetscInt_FMT "\n",PetscAbs(eps->nds)));
178: } else PetscTryTypeMethod(eps,view,viewer);
179: PetscCall(PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSELEMENTAL,EPSFEAST,EPSPRIMME,EPSSCALAPACK,EPSELPA,EPSEVSL,EPSTRLAN,""));
180: if (!isexternal) {
181: PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
182: if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
183: PetscCall(BVView(eps->V,viewer));
184: if (eps->rg) {
185: PetscCall(RGIsTrivial(eps->rg,&istrivial));
186: if (!istrivial) PetscCall(RGView(eps->rg,viewer));
187: }
188: if (eps->useds) {
189: if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
190: PetscCall(DSView(eps->ds,viewer));
191: }
192: PetscCall(PetscViewerPopFormat(viewer));
193: }
194: if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
195: PetscCall(STView(eps->st,viewer));
196: PetscFunctionReturn(PETSC_SUCCESS);
197: }
199: /*@
200: EPSViewFromOptions - View from options
202: Collective
204: Input Parameters:
205: + eps - the eigensolver context
206: . obj - optional object
207: - name - command line option
209: Level: intermediate
211: .seealso: EPSView(), EPSCreate()
212: @*/
213: PetscErrorCode EPSViewFromOptions(EPS eps,PetscObject obj,const char name[])
214: {
215: PetscFunctionBegin;
217: PetscCall(PetscObjectViewFromOptions((PetscObject)eps,obj,name));
218: PetscFunctionReturn(PETSC_SUCCESS);
219: }
221: /*@
222: EPSConvergedReasonView - Displays the reason an EPS solve converged or diverged.
224: Collective
226: Input Parameters:
227: + eps - the eigensolver context
228: - viewer - the viewer to display the reason
230: Options Database Keys:
231: . -eps_converged_reason - print reason for convergence, and number of iterations
233: Note:
234: To change the format of the output call PetscViewerPushFormat(viewer,format) before
235: this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
236: display a reason if it fails. The latter can be set in the command line with
237: -eps_converged_reason ::failed
239: Level: intermediate
241: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber(), EPSConvergedReasonViewFromOptions()
242: @*/
243: PetscErrorCode EPSConvergedReasonView(EPS eps,PetscViewer viewer)
244: {
245: PetscBool isAscii;
246: PetscViewerFormat format;
247: PetscInt nconv;
249: PetscFunctionBegin;
250: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
251: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
252: if (isAscii) {
253: PetscCall(PetscViewerGetFormat(viewer,&format));
254: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
255: PetscCall(EPS_GetActualConverged(eps,&nconv));
256: 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));
257: 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));
258: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
259: }
260: PetscFunctionReturn(PETSC_SUCCESS);
261: }
263: /*@
264: EPSConvergedReasonViewFromOptions - Processes command line options to determine if/how
265: the EPS converged reason is to be viewed.
267: Collective
269: Input Parameter:
270: . eps - the eigensolver context
272: Level: developer
274: .seealso: EPSConvergedReasonView()
275: @*/
276: PetscErrorCode EPSConvergedReasonViewFromOptions(EPS eps)
277: {
278: PetscViewer viewer;
279: PetscBool flg;
280: static PetscBool incall = PETSC_FALSE;
281: PetscViewerFormat format;
283: PetscFunctionBegin;
284: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
285: incall = PETSC_TRUE;
286: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg));
287: if (flg) {
288: PetscCall(PetscViewerPushFormat(viewer,format));
289: PetscCall(EPSConvergedReasonView(eps,viewer));
290: PetscCall(PetscViewerPopFormat(viewer));
291: PetscCall(PetscViewerDestroy(&viewer));
292: }
293: incall = PETSC_FALSE;
294: PetscFunctionReturn(PETSC_SUCCESS);
295: }
297: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
298: {
299: PetscReal error;
300: PetscScalar kr,ki;
301: PetscInt i,j,nvals,nconv;
303: PetscFunctionBegin;
304: PetscCall(EPS_GetActualConverged(eps,&nconv));
305: nvals = (eps->which==EPS_ALL)? eps->nconv: eps->nev;
306: if (eps->which!=EPS_ALL && nconv<nvals) {
307: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " eigenvalues converged\n\n",eps->nev));
308: PetscFunctionReturn(PETSC_SUCCESS);
309: }
310: if (eps->which==EPS_ALL && !nvals) {
311: PetscCall(PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n"));
312: PetscFunctionReturn(PETSC_SUCCESS);
313: }
314: for (i=0;i<nvals;i++) {
315: PetscCall(EPSComputeError(eps,i,etype,&error));
316: if (error>=5.0*eps->tol) {
317: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",nvals));
318: PetscFunctionReturn(PETSC_SUCCESS);
319: }
320: }
321: if (eps->which==EPS_ALL) PetscCall(PetscViewerASCIIPrintf(viewer," Found %" PetscInt_FMT " eigenvalues, all of them computed up to the required tolerance:",nvals));
322: else PetscCall(PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:"));
323: for (i=0;i<=(nvals-1)/8;i++) {
324: PetscCall(PetscViewerASCIIPrintf(viewer,"\n "));
325: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
326: PetscCall(EPSGetEigenvalue(eps,8*i+j,&kr,&ki));
327: PetscCall(SlepcPrintEigenvalueASCII(viewer,kr,ki));
328: if (8*i+j+1<nvals) PetscCall(PetscViewerASCIIPrintf(viewer,", "));
329: }
330: }
331: PetscCall(PetscViewerASCIIPrintf(viewer,"\n\n"));
332: PetscFunctionReturn(PETSC_SUCCESS);
333: }
335: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
336: {
337: PetscReal error,re,im;
338: PetscScalar kr,ki;
339: PetscInt i,nconv;
340: char ex[30],sep[]=" ---------------------- --------------------\n";
342: PetscFunctionBegin;
343: if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
344: switch (etype) {
345: case EPS_ERROR_ABSOLUTE:
346: PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||Ax-k%sx||",eps->isgeneralized?"B":""));
347: break;
348: case EPS_ERROR_RELATIVE:
349: PetscCall(PetscSNPrintf(ex,sizeof(ex),"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":""));
350: break;
351: case EPS_ERROR_BACKWARD:
352: PetscCall(PetscSNPrintf(ex,sizeof(ex)," eta(x,k)"));
353: break;
354: }
355: PetscCall(PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep));
356: PetscCall(EPS_GetActualConverged(eps,&nconv));
357: for (i=0;i<nconv;i++) {
358: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
359: PetscCall(EPSComputeError(eps,i,etype,&error));
360: #if defined(PETSC_USE_COMPLEX)
361: re = PetscRealPart(kr);
362: im = PetscImaginaryPart(kr);
363: #else
364: re = kr;
365: im = ki;
366: #endif
367: if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error));
368: else PetscCall(PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error));
369: }
370: PetscCall(PetscViewerASCIIPrintf(viewer,"%s",sep));
371: PetscFunctionReturn(PETSC_SUCCESS);
372: }
374: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
375: {
376: PetscReal error;
377: PetscInt i,nconv;
378: const char *name;
380: PetscFunctionBegin;
381: PetscCall(PetscObjectGetName((PetscObject)eps,&name));
382: PetscCall(PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name));
383: PetscCall(EPS_GetActualConverged(eps,&nconv));
384: for (i=0;i<nconv;i++) {
385: PetscCall(EPSComputeError(eps,i,etype,&error));
386: PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error));
387: }
388: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
389: PetscFunctionReturn(PETSC_SUCCESS);
390: }
392: /*@
393: EPSErrorView - Displays the errors associated with the computed solution
394: (as well as the eigenvalues).
396: Collective
398: Input Parameters:
399: + eps - the eigensolver context
400: . etype - error type
401: - viewer - optional visualization context
403: Options Database Keys:
404: + -eps_error_absolute - print absolute errors of each eigenpair
405: . -eps_error_relative - print relative errors of each eigenpair
406: - -eps_error_backward - print backward errors of each eigenpair
408: Notes:
409: By default, this function checks the error of all eigenpairs and prints
410: the eigenvalues if all of them are below the requested tolerance.
411: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
412: eigenvalues and corresponding errors is printed.
414: Level: intermediate
416: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
417: @*/
418: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
419: {
420: PetscBool isascii;
421: PetscViewerFormat format;
423: PetscFunctionBegin;
425: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
427: PetscCheckSameComm(eps,1,viewer,3);
428: EPSCheckSolved(eps,1);
429: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
430: if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);
432: PetscCall(PetscViewerGetFormat(viewer,&format));
433: switch (format) {
434: case PETSC_VIEWER_DEFAULT:
435: case PETSC_VIEWER_ASCII_INFO:
436: PetscCall(EPSErrorView_ASCII(eps,etype,viewer));
437: break;
438: case PETSC_VIEWER_ASCII_INFO_DETAIL:
439: PetscCall(EPSErrorView_DETAIL(eps,etype,viewer));
440: break;
441: case PETSC_VIEWER_ASCII_MATLAB:
442: PetscCall(EPSErrorView_MATLAB(eps,etype,viewer));
443: break;
444: default:
445: PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
446: }
447: PetscFunctionReturn(PETSC_SUCCESS);
448: }
450: /*@
451: EPSErrorViewFromOptions - Processes command line options to determine if/how
452: the errors of the computed solution are to be viewed.
454: Collective
456: Input Parameter:
457: . eps - the eigensolver context
459: Level: developer
461: .seealso: EPSErrorView()
462: @*/
463: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
464: {
465: PetscViewer viewer;
466: PetscBool flg;
467: static PetscBool incall = PETSC_FALSE;
468: PetscViewerFormat format;
470: PetscFunctionBegin;
471: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
472: incall = PETSC_TRUE;
473: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg));
474: if (flg) {
475: PetscCall(PetscViewerPushFormat(viewer,format));
476: PetscCall(EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer));
477: PetscCall(PetscViewerPopFormat(viewer));
478: PetscCall(PetscViewerDestroy(&viewer));
479: }
480: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg));
481: if (flg) {
482: PetscCall(PetscViewerPushFormat(viewer,format));
483: PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer));
484: PetscCall(PetscViewerPopFormat(viewer));
485: PetscCall(PetscViewerDestroy(&viewer));
486: }
487: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg));
488: if (flg) {
489: PetscCall(PetscViewerPushFormat(viewer,format));
490: PetscCall(EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer));
491: PetscCall(PetscViewerPopFormat(viewer));
492: PetscCall(PetscViewerDestroy(&viewer));
493: }
494: incall = PETSC_FALSE;
495: PetscFunctionReturn(PETSC_SUCCESS);
496: }
498: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
499: {
500: PetscDraw draw;
501: PetscDrawSP drawsp;
502: PetscReal re,im;
503: PetscScalar kr,ki;
504: PetscInt i,nconv;
506: PetscFunctionBegin;
507: if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
508: PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
509: PetscCall(PetscDrawSetTitle(draw,"Computed Eigenvalues"));
510: PetscCall(PetscDrawSPCreate(draw,1,&drawsp));
511: PetscCall(EPS_GetActualConverged(eps,&nconv));
512: for (i=0;i<nconv;i++) {
513: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
514: #if defined(PETSC_USE_COMPLEX)
515: re = PetscRealPart(kr);
516: im = PetscImaginaryPart(kr);
517: #else
518: re = kr;
519: im = ki;
520: #endif
521: PetscCall(PetscDrawSPAddPoint(drawsp,&re,&im));
522: }
523: PetscCall(PetscDrawSPDraw(drawsp,PETSC_TRUE));
524: PetscCall(PetscDrawSPSave(drawsp));
525: PetscCall(PetscDrawSPDestroy(&drawsp));
526: PetscFunctionReturn(PETSC_SUCCESS);
527: }
529: static PetscErrorCode EPSValuesView_BINARY(EPS eps,PetscViewer viewer)
530: {
531: PetscScalar kr,ki;
532: #if defined(PETSC_HAVE_COMPLEX)
533: PetscInt i,nconv;
534: PetscComplex *ev;
535: #endif
537: PetscFunctionBegin;
538: #if defined(PETSC_HAVE_COMPLEX)
539: PetscCall(EPS_GetActualConverged(eps,&nconv));
540: PetscCall(PetscMalloc1(nconv,&ev));
541: for (i=0;i<nconv;i++) {
542: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
543: #if defined(PETSC_USE_COMPLEX)
544: ev[i] = kr;
545: #else
546: ev[i] = PetscCMPLX(kr,ki);
547: #endif
548: }
549: PetscCall(PetscViewerBinaryWrite(viewer,ev,nconv,PETSC_COMPLEX));
550: PetscCall(PetscFree(ev));
551: #endif
552: PetscFunctionReturn(PETSC_SUCCESS);
553: }
555: #if defined(PETSC_HAVE_HDF5)
556: static PetscErrorCode EPSValuesView_HDF5(EPS eps,PetscViewer viewer)
557: {
558: PetscInt i,n,N,nconv;
559: PetscScalar eig;
560: PetscMPIInt rank;
561: Vec v;
562: char vname[30];
563: const char *ename;
565: PetscFunctionBegin;
566: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank));
567: PetscCall(EPS_GetActualConverged(eps,&nconv));
568: N = nconv;
569: n = rank? 0: N;
570: /* create a vector containing the eigenvalues */
571: PetscCall(VecCreateMPI(PetscObjectComm((PetscObject)eps),n,N,&v));
572: PetscCall(PetscObjectGetName((PetscObject)eps,&ename));
573: PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename));
574: PetscCall(PetscObjectSetName((PetscObject)v,vname));
575: if (!rank) {
576: for (i=0;i<nconv;i++) {
577: PetscCall(EPSGetEigenvalue(eps,i,&eig,NULL));
578: PetscCall(VecSetValue(v,i,eig,INSERT_VALUES));
579: }
580: }
581: PetscCall(VecAssemblyBegin(v));
582: PetscCall(VecAssemblyEnd(v));
583: PetscCall(VecView(v,viewer));
584: #if !defined(PETSC_USE_COMPLEX)
585: /* in real scalars write the imaginary part as a separate vector */
586: PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename));
587: PetscCall(PetscObjectSetName((PetscObject)v,vname));
588: if (!rank) {
589: for (i=0;i<nconv;i++) {
590: PetscCall(EPSGetEigenvalue(eps,i,NULL,&eig));
591: PetscCall(VecSetValue(v,i,eig,INSERT_VALUES));
592: }
593: }
594: PetscCall(VecAssemblyBegin(v));
595: PetscCall(VecAssemblyEnd(v));
596: PetscCall(VecView(v,viewer));
597: #endif
598: PetscCall(VecDestroy(&v));
599: PetscFunctionReturn(PETSC_SUCCESS);
600: }
601: #endif
603: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
604: {
605: PetscInt i,nconv;
606: PetscScalar kr,ki;
608: PetscFunctionBegin;
609: PetscCall(PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n"));
610: PetscCall(EPS_GetActualConverged(eps,&nconv));
611: for (i=0;i<nconv;i++) {
612: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
613: PetscCall(PetscViewerASCIIPrintf(viewer," "));
614: PetscCall(SlepcPrintEigenvalueASCII(viewer,kr,ki));
615: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
616: }
617: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
618: PetscFunctionReturn(PETSC_SUCCESS);
619: }
621: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
622: {
623: PetscInt i,nconv;
624: PetscReal re,im;
625: PetscScalar kr,ki;
626: const char *name;
628: PetscFunctionBegin;
629: PetscCall(PetscObjectGetName((PetscObject)eps,&name));
630: PetscCall(PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name));
631: PetscCall(EPS_GetActualConverged(eps,&nconv));
632: for (i=0;i<nconv;i++) {
633: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
634: #if defined(PETSC_USE_COMPLEX)
635: re = PetscRealPart(kr);
636: im = PetscImaginaryPart(kr);
637: #else
638: re = kr;
639: im = ki;
640: #endif
641: if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im));
642: else PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re));
643: }
644: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
645: PetscFunctionReturn(PETSC_SUCCESS);
646: }
648: /*@
649: EPSValuesView - Displays the computed eigenvalues in a viewer.
651: Collective
653: Input Parameters:
654: + eps - the eigensolver context
655: - viewer - the viewer
657: Options Database Key:
658: . -eps_view_values - print computed eigenvalues
660: Level: intermediate
662: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
663: @*/
664: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
665: {
666: PetscBool isascii,isdraw,isbinary;
667: PetscViewerFormat format;
668: #if defined(PETSC_HAVE_HDF5)
669: PetscBool ishdf5;
670: #endif
672: PetscFunctionBegin;
674: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
676: PetscCheckSameComm(eps,1,viewer,2);
677: EPSCheckSolved(eps,1);
678: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
679: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
680: #if defined(PETSC_HAVE_HDF5)
681: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
682: #endif
683: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
684: if (isdraw) PetscCall(EPSValuesView_DRAW(eps,viewer));
685: else if (isbinary) PetscCall(EPSValuesView_BINARY(eps,viewer));
686: #if defined(PETSC_HAVE_HDF5)
687: else if (ishdf5) PetscCall(EPSValuesView_HDF5(eps,viewer));
688: #endif
689: else if (isascii) {
690: PetscCall(PetscViewerGetFormat(viewer,&format));
691: switch (format) {
692: case PETSC_VIEWER_DEFAULT:
693: case PETSC_VIEWER_ASCII_INFO:
694: case PETSC_VIEWER_ASCII_INFO_DETAIL:
695: PetscCall(EPSValuesView_ASCII(eps,viewer));
696: break;
697: case PETSC_VIEWER_ASCII_MATLAB:
698: PetscCall(EPSValuesView_MATLAB(eps,viewer));
699: break;
700: default:
701: PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
702: }
703: }
704: PetscFunctionReturn(PETSC_SUCCESS);
705: }
707: /*@
708: EPSValuesViewFromOptions - Processes command line options to determine if/how
709: the computed eigenvalues are to be viewed.
711: Collective
713: Input Parameters:
714: . eps - the eigensolver context
716: Level: developer
718: .seealso: EPSValuesView()
719: @*/
720: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
721: {
722: PetscViewer viewer;
723: PetscBool flg;
724: static PetscBool incall = PETSC_FALSE;
725: PetscViewerFormat format;
727: PetscFunctionBegin;
728: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
729: incall = PETSC_TRUE;
730: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg));
731: if (flg) {
732: PetscCall(PetscViewerPushFormat(viewer,format));
733: PetscCall(EPSValuesView(eps,viewer));
734: PetscCall(PetscViewerPopFormat(viewer));
735: PetscCall(PetscViewerDestroy(&viewer));
736: }
737: incall = PETSC_FALSE;
738: PetscFunctionReturn(PETSC_SUCCESS);
739: }
741: /*@
742: EPSVectorsView - Outputs computed eigenvectors to a viewer.
744: Collective
746: Input Parameters:
747: + eps - the eigensolver context
748: - viewer - the viewer
750: Options Database Key:
751: . -eps_view_vectors - output eigenvectors.
753: Notes:
754: If PETSc was configured with real scalars, complex conjugate eigenvectors
755: will be viewed as two separate real vectors, one containing the real part
756: and another one containing the imaginary part.
758: If left eigenvectors were computed with a two-sided eigensolver, the right
759: and left eigenvectors are interleaved, that is, the vectors are output in
760: the following order X0, Y0, X1, Y1, X2, Y2, ...
762: Level: intermediate
764: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
765: @*/
766: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
767: {
768: PetscInt i,nconv;
769: Vec xr,xi=NULL;
771: PetscFunctionBegin;
773: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
775: PetscCheckSameComm(eps,1,viewer,2);
776: EPSCheckSolved(eps,1);
777: PetscCall(EPS_GetActualConverged(eps,&nconv));
778: if (nconv) {
779: PetscCall(BVCreateVec(eps->V,&xr));
780: #if !defined(PETSC_USE_COMPLEX)
781: PetscCall(BVCreateVec(eps->V,&xi));
782: #endif
783: for (i=0;i<nconv;i++) {
784: PetscCall(EPSGetEigenvector(eps,i,xr,xi));
785: PetscCall(SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)eps));
786: if (eps->twosided || eps->problem_type==EPS_BSE) {
787: PetscCall(EPSGetLeftEigenvector(eps,i,xr,xi));
788: PetscCall(SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)eps));
789: }
790: }
791: PetscCall(VecDestroy(&xr));
792: #if !defined(PETSC_USE_COMPLEX)
793: PetscCall(VecDestroy(&xi));
794: #endif
795: }
796: PetscFunctionReturn(PETSC_SUCCESS);
797: }
799: /*@
800: EPSVectorsViewFromOptions - Processes command line options to determine if/how
801: the computed eigenvectors are to be viewed.
803: Collective
805: Input Parameter:
806: . eps - the eigensolver context
808: Level: developer
810: .seealso: EPSVectorsView()
811: @*/
812: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
813: {
814: PetscViewer viewer;
815: PetscBool flg = PETSC_FALSE;
816: static PetscBool incall = PETSC_FALSE;
817: PetscViewerFormat format;
819: PetscFunctionBegin;
820: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
821: incall = PETSC_TRUE;
822: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg));
823: if (flg) {
824: PetscCall(PetscViewerPushFormat(viewer,format));
825: PetscCall(EPSVectorsView(eps,viewer));
826: PetscCall(PetscViewerPopFormat(viewer));
827: PetscCall(PetscViewerDestroy(&viewer));
828: }
829: incall = PETSC_FALSE;
830: PetscFunctionReturn(PETSC_SUCCESS);
831: }