Actual source code: epsview.c
slepc-3.11.2 2019-07-30
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, 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> /*I "slepceps.h" I*/
15: #include <petscdraw.h>
17: /*@C
18: EPSView - Prints the EPS data structure.
20: Collective on EPS
22: Input Parameters:
23: + eps - the eigenproblem solver context
24: - viewer - optional visualization context
26: Options Database Key:
27: . -eps_view - Calls EPSView() at end of EPSSolve()
29: Note:
30: The available visualization contexts include
31: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
32: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
33: output where only the first processor opens
34: the file. All other processors send their
35: data to the first processor to print.
37: The user can open an alternative visualization context with
38: PetscViewerASCIIOpen() - output to a specified file.
40: Level: beginner
42: .seealso: STView(), PetscViewerASCIIOpen()
43: @*/
44: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
45: {
47: const char *type=NULL,*extr=NULL,*bal=NULL;
48: char str[50];
49: PetscBool isascii,isexternal,istrivial;
53: if (!viewer) {
54: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
55: }
59: #if defined(PETSC_USE_COMPLEX)
60: #define HERM "hermitian"
61: #else
62: #define HERM "symmetric"
63: #endif
64: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
65: if (isascii) {
66: PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer);
67: if (eps->ops->view) {
68: PetscViewerASCIIPushTab(viewer);
69: (*eps->ops->view)(eps,viewer);
70: PetscViewerASCIIPopTab(viewer);
71: }
72: if (eps->problem_type) {
73: switch (eps->problem_type) {
74: case EPS_HEP: type = HERM " eigenvalue problem"; break;
75: case EPS_GHEP: type = "generalized " HERM " eigenvalue problem"; break;
76: case EPS_NHEP: type = "non-" HERM " eigenvalue problem"; break;
77: case EPS_GNHEP: type = "generalized non-" HERM " eigenvalue problem"; break;
78: case EPS_PGNHEP: type = "generalized non-" HERM " eigenvalue problem with " HERM " positive definite B"; break;
79: case EPS_GHIEP: type = "generalized " HERM "-indefinite eigenvalue problem"; break;
80: }
81: } else type = "not yet set";
82: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
83: if (eps->extraction) {
84: switch (eps->extraction) {
85: case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
86: case EPS_HARMONIC: extr = "harmonic Ritz"; break;
87: case EPS_HARMONIC_RELATIVE: extr = "relative harmonic Ritz"; break;
88: case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
89: case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
90: case EPS_REFINED: extr = "refined Ritz"; break;
91: case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
92: }
93: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr);
94: }
95: if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
96: switch (eps->balance) {
97: case EPS_BALANCE_NONE: break;
98: case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
99: case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
100: case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
101: }
102: PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal);
103: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
104: if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) {
105: PetscViewerASCIIPrintf(viewer,", with its=%D",eps->balance_its);
106: }
107: if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) {
108: PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff);
109: }
110: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
111: PetscViewerASCIIPrintf(viewer,"\n");
112: }
113: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
114: SlepcSNPrintfScalar(str,50,eps->target,PETSC_FALSE);
115: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
116: if (!eps->which) {
117: PetscViewerASCIIPrintf(viewer,"not yet set\n");
118: } else switch (eps->which) {
119: case EPS_WHICH_USER:
120: PetscViewerASCIIPrintf(viewer,"user defined\n");
121: break;
122: case EPS_TARGET_MAGNITUDE:
123: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
124: break;
125: case EPS_TARGET_REAL:
126: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
127: break;
128: case EPS_TARGET_IMAGINARY:
129: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
130: break;
131: case EPS_LARGEST_MAGNITUDE:
132: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
133: break;
134: case EPS_SMALLEST_MAGNITUDE:
135: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
136: break;
137: case EPS_LARGEST_REAL:
138: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
139: break;
140: case EPS_SMALLEST_REAL:
141: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
142: break;
143: case EPS_LARGEST_IMAGINARY:
144: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
145: break;
146: case EPS_SMALLEST_IMAGINARY:
147: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
148: break;
149: case EPS_ALL:
150: if (eps->inta || eps->intb) {
151: PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb);
152: } else {
153: PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
154: }
155: break;
156: }
157: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
158: if (eps->twosided && eps->problem_type!=EPS_HEP && eps->problem_type!=EPS_GHEP) {
159: PetscViewerASCIIPrintf(viewer," using two-sided variant (for left eigenvectors)\n");
160: }
161: if (eps->purify) {
162: PetscViewerASCIIPrintf(viewer," postprocessing eigenvectors with purification\n");
163: }
164: if (eps->trueres) {
165: PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n");
166: }
167: if (eps->trackall) {
168: PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n");
169: }
170: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",eps->nev);
171: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",eps->ncv);
172: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",eps->mpd);
173: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",eps->max_it);
174: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)eps->tol);
175: PetscViewerASCIIPrintf(viewer," convergence test: ");
176: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
177: switch (eps->conv) {
178: case EPS_CONV_ABS:
179: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
180: case EPS_CONV_REL:
181: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
182: case EPS_CONV_NORM:
183: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");
184: PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)eps->nrma);
185: if (eps->isgeneralized) {
186: PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb);
187: }
188: PetscViewerASCIIPrintf(viewer,"\n");
189: break;
190: case EPS_CONV_USER:
191: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
192: }
193: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
194: if (eps->nini) {
195: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(eps->nini));
196: }
197: if (eps->nds) {
198: PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %D\n",PetscAbs(eps->nds));
199: }
200: } else {
201: if (eps->ops->view) {
202: (*eps->ops->view)(eps,viewer);
203: }
204: }
205: PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSBLZPACK,EPSFEAST,EPSPRIMME,EPSTRLAN,"");
206: if (!isexternal) {
207: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
208: if (!eps->V) { EPSGetBV(eps,&eps->V); }
209: BVView(eps->V,viewer);
210: if (eps->rg) {
211: RGIsTrivial(eps->rg,&istrivial);
212: if (!istrivial) { RGView(eps->rg,viewer); }
213: }
214: if (eps->useds) {
215: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
216: DSView(eps->ds,viewer);
217: }
218: PetscViewerPopFormat(viewer);
219: }
220: if (!eps->st) { EPSGetST(eps,&eps->st); }
221: STView(eps->st,viewer);
222: return(0);
223: }
225: /*@C
226: EPSReasonView - Displays the reason an EPS solve converged or diverged.
228: Collective on EPS
230: Parameter:
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: Level: intermediate
239: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber()
240: @*/
241: PetscErrorCode EPSReasonView(EPS eps,PetscViewer viewer)
242: {
244: PetscBool isAscii;
247: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
248: if (isAscii) {
249: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
250: if (eps->reason > 0) {
251: PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve converged (%D eigenpair%s) due to %s; iterations %D\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",eps->nconv,(eps->nconv>1)?"s":"",EPSConvergedReasons[eps->reason],eps->its);
252: } else {
253: PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve did not converge due to %s; iterations %D\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",EPSConvergedReasons[eps->reason],eps->its);
254: }
255: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
256: }
257: return(0);
258: }
260: /*@
261: EPSReasonViewFromOptions - Processes command line options to determine if/how
262: the EPS converged reason is to be viewed.
264: Collective on EPS
266: Input Parameters:
267: . eps - the eigensolver context
269: Level: developer
270: @*/
271: PetscErrorCode EPSReasonViewFromOptions(EPS eps)
272: {
273: PetscErrorCode ierr;
274: PetscViewer viewer;
275: PetscBool flg;
276: static PetscBool incall = PETSC_FALSE;
277: PetscViewerFormat format;
280: if (incall) return(0);
281: incall = PETSC_TRUE;
282: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg);
283: if (flg) {
284: PetscViewerPushFormat(viewer,format);
285: EPSReasonView(eps,viewer);
286: PetscViewerPopFormat(viewer);
287: PetscViewerDestroy(&viewer);
288: }
289: incall = PETSC_FALSE;
290: return(0);
291: }
293: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
294: {
295: PetscReal error;
296: PetscInt i,j,k,nvals;
300: nvals = (eps->which==EPS_ALL)? eps->nconv: eps->nev;
301: if (eps->which!=EPS_ALL && eps->nconv<nvals) {
302: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
303: return(0);
304: }
305: if (eps->which==EPS_ALL && !nvals) {
306: PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
307: return(0);
308: }
309: for (i=0;i<nvals;i++) {
310: EPSComputeError(eps,i,etype,&error);
311: if (error>=5.0*eps->tol) {
312: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
313: return(0);
314: }
315: }
316: if (eps->which==EPS_ALL) {
317: PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
318: } else {
319: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
320: }
321: for (i=0;i<=(nvals-1)/8;i++) {
322: PetscViewerASCIIPrintf(viewer,"\n ");
323: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
324: k = eps->perm[8*i+j];
325: SlepcPrintEigenvalueASCII(eps->eigr[k],eps->eigi[k]);
326: if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
327: }
328: }
329: PetscViewerASCIIPrintf(viewer,"\n\n");
330: return(0);
331: }
333: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
334: {
336: PetscReal error,re,im;
337: PetscScalar kr,ki;
338: PetscInt i;
339: #define EXLEN 30
340: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
343: if (!eps->nconv) return(0);
344: switch (etype) {
345: case EPS_ERROR_ABSOLUTE:
346: PetscSNPrintf(ex,EXLEN," ||Ax-k%sx||",eps->isgeneralized?"B":"");
347: break;
348: case EPS_ERROR_RELATIVE:
349: PetscSNPrintf(ex,EXLEN,"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":"");
350: break;
351: case EPS_ERROR_BACKWARD:
352: PetscSNPrintf(ex,EXLEN," eta(x,k)");
353: break;
354: }
355: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
356: for (i=0;i<eps->nconv;i++) {
357: EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
358: EPSComputeError(eps,i,etype,&error);
359: #if defined(PETSC_USE_COMPLEX)
360: re = PetscRealPart(kr);
361: im = PetscImaginaryPart(kr);
362: #else
363: re = kr;
364: im = ki;
365: #endif
366: if (im!=0.0) {
367: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
368: } else {
369: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
370: }
371: }
372: PetscViewerASCIIPrintf(viewer,sep);
373: return(0);
374: }
376: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
377: {
379: PetscReal error;
380: PetscInt i;
381: const char *name;
384: PetscObjectGetName((PetscObject)eps,&name);
385: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
386: for (i=0;i<eps->nconv;i++) {
387: EPSComputeError(eps,i,etype,&error);
388: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
389: }
390: PetscViewerASCIIPrintf(viewer,"];\n");
391: return(0);
392: }
394: /*@C
395: EPSErrorView - Displays the errors associated with the computed solution
396: (as well as the eigenvalues).
398: Collective on EPS
400: Input Parameters:
401: + eps - the eigensolver context
402: . etype - error type
403: - viewer - optional visualization context
405: Options Database Key:
406: + -eps_error_absolute - print absolute errors of each eigenpair
407: . -eps_error_relative - print relative errors of each eigenpair
408: - -eps_error_backward - print backward errors of each eigenpair
410: Notes:
411: By default, this function checks the error of all eigenpairs and prints
412: the eigenvalues if all of them are below the requested tolerance.
413: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
414: eigenvalues and corresponding errors is printed.
416: Level: intermediate
418: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
419: @*/
420: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
421: {
422: PetscBool isascii;
423: PetscViewerFormat format;
424: PetscErrorCode ierr;
428: if (!viewer) {
429: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
430: }
433: EPSCheckSolved(eps,1);
434: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
435: if (!isascii) return(0);
437: PetscViewerGetFormat(viewer,&format);
438: switch (format) {
439: case PETSC_VIEWER_DEFAULT:
440: case PETSC_VIEWER_ASCII_INFO:
441: EPSErrorView_ASCII(eps,etype,viewer);
442: break;
443: case PETSC_VIEWER_ASCII_INFO_DETAIL:
444: EPSErrorView_DETAIL(eps,etype,viewer);
445: break;
446: case PETSC_VIEWER_ASCII_MATLAB:
447: EPSErrorView_MATLAB(eps,etype,viewer);
448: break;
449: default:
450: PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
451: }
452: return(0);
453: }
455: /*@
456: EPSErrorViewFromOptions - Processes command line options to determine if/how
457: the errors of the computed solution are to be viewed.
459: Collective on EPS
461: Input Parameters:
462: . eps - the eigensolver context
464: Level: developer
465: @*/
466: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
467: {
468: PetscErrorCode ierr;
469: PetscViewer viewer;
470: PetscBool flg;
471: static PetscBool incall = PETSC_FALSE;
472: PetscViewerFormat format;
475: if (incall) return(0);
476: incall = PETSC_TRUE;
477: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg);
478: if (flg) {
479: PetscViewerPushFormat(viewer,format);
480: EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer);
481: PetscViewerPopFormat(viewer);
482: PetscViewerDestroy(&viewer);
483: }
484: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg);
485: if (flg) {
486: PetscViewerPushFormat(viewer,format);
487: EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer);
488: PetscViewerPopFormat(viewer);
489: PetscViewerDestroy(&viewer);
490: }
491: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg);
492: if (flg) {
493: PetscViewerPushFormat(viewer,format);
494: EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer);
495: PetscViewerPopFormat(viewer);
496: PetscViewerDestroy(&viewer);
497: }
498: incall = PETSC_FALSE;
499: return(0);
500: }
502: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
503: {
505: PetscDraw draw;
506: PetscDrawSP drawsp;
507: PetscReal re,im;
508: PetscInt i,k;
511: if (!eps->nconv) return(0);
512: PetscViewerDrawGetDraw(viewer,0,&draw);
513: PetscDrawSetTitle(draw,"Computed Eigenvalues");
514: PetscDrawSPCreate(draw,1,&drawsp);
515: for (i=0;i<eps->nconv;i++) {
516: k = eps->perm[i];
517: #if defined(PETSC_USE_COMPLEX)
518: re = PetscRealPart(eps->eigr[k]);
519: im = PetscImaginaryPart(eps->eigr[k]);
520: #else
521: re = eps->eigr[k];
522: im = eps->eigi[k];
523: #endif
524: PetscDrawSPAddPoint(drawsp,&re,&im);
525: }
526: PetscDrawSPDraw(drawsp,PETSC_TRUE);
527: PetscDrawSPSave(drawsp);
528: PetscDrawSPDestroy(&drawsp);
529: return(0);
530: }
532: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
533: {
534: PetscInt i,k;
538: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
539: for (i=0;i<eps->nconv;i++) {
540: k = eps->perm[i];
541: PetscViewerASCIIPrintf(viewer," ");
542: SlepcPrintEigenvalueASCII(eps->eigr[k],eps->eigi[k]);
543: PetscViewerASCIIPrintf(viewer,"\n");
544: }
545: PetscViewerASCIIPrintf(viewer,"\n");
546: return(0);
547: }
549: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
550: {
552: PetscInt i,k;
553: PetscReal re,im;
554: const char *name;
557: PetscObjectGetName((PetscObject)eps,&name);
558: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
559: for (i=0;i<eps->nconv;i++) {
560: k = eps->perm[i];
561: #if defined(PETSC_USE_COMPLEX)
562: re = PetscRealPart(eps->eigr[k]);
563: im = PetscImaginaryPart(eps->eigr[k]);
564: #else
565: re = eps->eigr[k];
566: im = eps->eigi[k];
567: #endif
568: if (im!=0.0) {
569: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
570: } else {
571: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
572: }
573: }
574: PetscViewerASCIIPrintf(viewer,"];\n");
575: return(0);
576: }
578: /*@C
579: EPSValuesView - Displays the computed eigenvalues in a viewer.
581: Collective on EPS
583: Input Parameters:
584: + eps - the eigensolver context
585: - viewer - the viewer
587: Options Database Key:
588: . -eps_view_values - print computed eigenvalues
590: Level: intermediate
592: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
593: @*/
594: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
595: {
596: PetscBool isascii,isdraw;
597: PetscViewerFormat format;
598: PetscErrorCode ierr;
602: if (!viewer) {
603: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
604: }
607: EPSCheckSolved(eps,1);
608: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
609: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
610: if (isdraw) {
611: EPSValuesView_DRAW(eps,viewer);
612: } else if (isascii) {
613: PetscViewerGetFormat(viewer,&format);
614: switch (format) {
615: case PETSC_VIEWER_DEFAULT:
616: case PETSC_VIEWER_ASCII_INFO:
617: case PETSC_VIEWER_ASCII_INFO_DETAIL:
618: EPSValuesView_ASCII(eps,viewer);
619: break;
620: case PETSC_VIEWER_ASCII_MATLAB:
621: EPSValuesView_MATLAB(eps,viewer);
622: break;
623: default:
624: PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
625: }
626: }
627: return(0);
628: }
630: /*@
631: EPSValuesViewFromOptions - Processes command line options to determine if/how
632: the computed eigenvalues are to be viewed.
634: Collective on EPS
636: Input Parameters:
637: . eps - the eigensolver context
639: Level: developer
640: @*/
641: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
642: {
643: PetscErrorCode ierr;
644: PetscViewer viewer;
645: PetscBool flg;
646: static PetscBool incall = PETSC_FALSE;
647: PetscViewerFormat format;
650: if (incall) return(0);
651: incall = PETSC_TRUE;
652: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg);
653: if (flg) {
654: PetscViewerPushFormat(viewer,format);
655: EPSValuesView(eps,viewer);
656: PetscViewerPopFormat(viewer);
657: PetscViewerDestroy(&viewer);
658: }
659: incall = PETSC_FALSE;
660: return(0);
661: }
663: /*@C
664: EPSVectorsView - Outputs computed eigenvectors to a viewer.
666: Collective on EPS
668: Parameter:
669: + eps - the eigensolver context
670: - viewer - the viewer
672: Options Database Keys:
673: . -eps_view_vectors - output eigenvectors.
675: Notes:
676: If PETSc was configured with real scalars, complex conjugate eigenvectors
677: will be viewed as two separate real vectors, one containing the real part
678: and another one containing the imaginary part.
680: If left eigenvectors were computed with a two-sided eigensolver, the right
681: and left eigenvectors are interleaved, that is, the vectors are output in
682: the following order: X0, Y0, X1, Y1, X2, Y2, ...
684: Level: intermediate
686: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
687: @*/
688: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
689: {
691: PetscInt i,k;
692: Vec x;
693: #define NMLEN 30
694: char vname[NMLEN];
695: const char *ename;
699: if (!viewer) {
700: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
701: }
704: EPSCheckSolved(eps,1);
705: if (eps->nconv) {
706: PetscObjectGetName((PetscObject)eps,&ename);
707: EPSComputeVectors(eps);
708: for (i=0;i<eps->nconv;i++) {
709: k = eps->perm[i];
710: PetscSNPrintf(vname,NMLEN,"X%d_%s",(int)i,ename);
711: BVGetColumn(eps->V,k,&x);
712: PetscObjectSetName((PetscObject)x,vname);
713: VecView(x,viewer);
714: BVRestoreColumn(eps->V,k,&x);
715: if (eps->twosided) {
716: PetscSNPrintf(vname,NMLEN,"Y%d_%s",(int)i,ename);
717: BVGetColumn(eps->W,k,&x);
718: PetscObjectSetName((PetscObject)x,vname);
719: VecView(x,viewer);
720: BVRestoreColumn(eps->W,k,&x);
721: }
722: }
723: }
724: return(0);
725: }
727: /*@
728: EPSVectorsViewFromOptions - Processes command line options to determine if/how
729: the computed eigenvectors are to be viewed.
731: Collective on EPS
733: Input Parameters:
734: . eps - the eigensolver context
736: Level: developer
737: @*/
738: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
739: {
740: PetscErrorCode ierr;
741: PetscViewer viewer;
742: PetscBool flg = PETSC_FALSE;
743: static PetscBool incall = PETSC_FALSE;
744: PetscViewerFormat format;
747: if (incall) return(0);
748: incall = PETSC_TRUE;
749: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg);
750: if (flg) {
751: PetscViewerPushFormat(viewer,format);
752: EPSVectorsView(eps,viewer);
753: PetscViewerPopFormat(viewer);
754: PetscViewerDestroy(&viewer);
755: }
756: incall = PETSC_FALSE;
757: return(0);
758: }