Actual source code: epsview.c
slepc-3.13.1 2020-04-12
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, 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()
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->ninil) {
198: PetscViewerASCIIPrintf(viewer," dimension of user-provided left initial space: %D\n",PetscAbs(eps->ninil));
199: }
200: if (eps->nds) {
201: PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %D\n",PetscAbs(eps->nds));
202: }
203: } else {
204: if (eps->ops->view) {
205: (*eps->ops->view)(eps,viewer);
206: }
207: }
208: PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSBLZPACK,EPSPRIMME,EPSTRLAN,"");
209: if (!isexternal) {
210: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
211: if (!eps->V) { EPSGetBV(eps,&eps->V); }
212: BVView(eps->V,viewer);
213: if (eps->rg) {
214: RGIsTrivial(eps->rg,&istrivial);
215: if (!istrivial) { RGView(eps->rg,viewer); }
216: }
217: if (eps->useds) {
218: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
219: DSView(eps->ds,viewer);
220: }
221: PetscViewerPopFormat(viewer);
222: }
223: if (!eps->st) { EPSGetST(eps,&eps->st); }
224: STView(eps->st,viewer);
225: return(0);
226: }
228: /*@C
229: EPSViewFromOptions - View from options
231: Collective on EPS
233: Input Parameters:
234: + eps - the eigensolver context
235: . obj - optional object
236: - name - command line option
238: Level: intermediate
240: .seealso: EPSView(), EPSCreate()
241: @*/
242: PetscErrorCode EPSViewFromOptions(EPS eps,PetscObject obj,const char name[])
243: {
248: PetscObjectViewFromOptions((PetscObject)eps,obj,name);
249: return(0);
250: }
252: /*@C
253: EPSReasonView - Displays the reason an EPS solve converged or diverged.
255: Collective on eps
257: Input Parameters:
258: + eps - the eigensolver context
259: - viewer - the viewer to display the reason
261: Options Database Keys:
262: . -eps_converged_reason - print reason for convergence, and number of iterations
264: Level: intermediate
266: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber()
267: @*/
268: PetscErrorCode EPSReasonView(EPS eps,PetscViewer viewer)
269: {
271: PetscBool isAscii;
274: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
275: if (isAscii) {
276: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
277: if (eps->reason > 0) {
278: 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);
279: } else {
280: 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);
281: }
282: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
283: }
284: return(0);
285: }
287: /*@
288: EPSReasonViewFromOptions - Processes command line options to determine if/how
289: the EPS converged reason is to be viewed.
291: Collective on eps
293: Input Parameter:
294: . eps - the eigensolver context
296: Level: developer
297: @*/
298: PetscErrorCode EPSReasonViewFromOptions(EPS eps)
299: {
300: PetscErrorCode ierr;
301: PetscViewer viewer;
302: PetscBool flg;
303: static PetscBool incall = PETSC_FALSE;
304: PetscViewerFormat format;
307: if (incall) return(0);
308: incall = PETSC_TRUE;
309: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg);
310: if (flg) {
311: PetscViewerPushFormat(viewer,format);
312: EPSReasonView(eps,viewer);
313: PetscViewerPopFormat(viewer);
314: PetscViewerDestroy(&viewer);
315: }
316: incall = PETSC_FALSE;
317: return(0);
318: }
320: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
321: {
322: PetscReal error;
323: PetscInt i,j,k,nvals;
327: nvals = (eps->which==EPS_ALL)? eps->nconv: eps->nev;
328: if (eps->which!=EPS_ALL && eps->nconv<nvals) {
329: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
330: return(0);
331: }
332: if (eps->which==EPS_ALL && !nvals) {
333: PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
334: return(0);
335: }
336: for (i=0;i<nvals;i++) {
337: EPSComputeError(eps,i,etype,&error);
338: if (error>=5.0*eps->tol) {
339: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
340: return(0);
341: }
342: }
343: if (eps->which==EPS_ALL) {
344: PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
345: } else {
346: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
347: }
348: for (i=0;i<=(nvals-1)/8;i++) {
349: PetscViewerASCIIPrintf(viewer,"\n ");
350: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
351: k = eps->perm[8*i+j];
352: SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]);
353: if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
354: }
355: }
356: PetscViewerASCIIPrintf(viewer,"\n\n");
357: return(0);
358: }
360: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
361: {
363: PetscReal error,re,im;
364: PetscScalar kr,ki;
365: PetscInt i;
366: #define EXLEN 30
367: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
370: if (!eps->nconv) return(0);
371: switch (etype) {
372: case EPS_ERROR_ABSOLUTE:
373: PetscSNPrintf(ex,EXLEN," ||Ax-k%sx||",eps->isgeneralized?"B":"");
374: break;
375: case EPS_ERROR_RELATIVE:
376: PetscSNPrintf(ex,EXLEN,"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":"");
377: break;
378: case EPS_ERROR_BACKWARD:
379: PetscSNPrintf(ex,EXLEN," eta(x,k)");
380: break;
381: }
382: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
383: for (i=0;i<eps->nconv;i++) {
384: EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
385: EPSComputeError(eps,i,etype,&error);
386: #if defined(PETSC_USE_COMPLEX)
387: re = PetscRealPart(kr);
388: im = PetscImaginaryPart(kr);
389: #else
390: re = kr;
391: im = ki;
392: #endif
393: if (im!=0.0) {
394: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
395: } else {
396: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
397: }
398: }
399: PetscViewerASCIIPrintf(viewer,sep);
400: return(0);
401: }
403: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
404: {
406: PetscReal error;
407: PetscInt i;
408: const char *name;
411: PetscObjectGetName((PetscObject)eps,&name);
412: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
413: for (i=0;i<eps->nconv;i++) {
414: EPSComputeError(eps,i,etype,&error);
415: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
416: }
417: PetscViewerASCIIPrintf(viewer,"];\n");
418: return(0);
419: }
421: /*@C
422: EPSErrorView - Displays the errors associated with the computed solution
423: (as well as the eigenvalues).
425: Collective on eps
427: Input Parameters:
428: + eps - the eigensolver context
429: . etype - error type
430: - viewer - optional visualization context
432: Options Database Key:
433: + -eps_error_absolute - print absolute errors of each eigenpair
434: . -eps_error_relative - print relative errors of each eigenpair
435: - -eps_error_backward - print backward errors of each eigenpair
437: Notes:
438: By default, this function checks the error of all eigenpairs and prints
439: the eigenvalues if all of them are below the requested tolerance.
440: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
441: eigenvalues and corresponding errors is printed.
443: Level: intermediate
445: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
446: @*/
447: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
448: {
449: PetscBool isascii;
450: PetscViewerFormat format;
451: PetscErrorCode ierr;
455: if (!viewer) {
456: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
457: }
460: EPSCheckSolved(eps,1);
461: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
462: if (!isascii) return(0);
464: PetscViewerGetFormat(viewer,&format);
465: switch (format) {
466: case PETSC_VIEWER_DEFAULT:
467: case PETSC_VIEWER_ASCII_INFO:
468: EPSErrorView_ASCII(eps,etype,viewer);
469: break;
470: case PETSC_VIEWER_ASCII_INFO_DETAIL:
471: EPSErrorView_DETAIL(eps,etype,viewer);
472: break;
473: case PETSC_VIEWER_ASCII_MATLAB:
474: EPSErrorView_MATLAB(eps,etype,viewer);
475: break;
476: default:
477: PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
478: }
479: return(0);
480: }
482: /*@
483: EPSErrorViewFromOptions - Processes command line options to determine if/how
484: the errors of the computed solution are to be viewed.
486: Collective on eps
488: Input Parameter:
489: . eps - the eigensolver context
491: Level: developer
492: @*/
493: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
494: {
495: PetscErrorCode ierr;
496: PetscViewer viewer;
497: PetscBool flg;
498: static PetscBool incall = PETSC_FALSE;
499: PetscViewerFormat format;
502: if (incall) return(0);
503: incall = PETSC_TRUE;
504: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg);
505: if (flg) {
506: PetscViewerPushFormat(viewer,format);
507: EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer);
508: PetscViewerPopFormat(viewer);
509: PetscViewerDestroy(&viewer);
510: }
511: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg);
512: if (flg) {
513: PetscViewerPushFormat(viewer,format);
514: EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer);
515: PetscViewerPopFormat(viewer);
516: PetscViewerDestroy(&viewer);
517: }
518: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg);
519: if (flg) {
520: PetscViewerPushFormat(viewer,format);
521: EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer);
522: PetscViewerPopFormat(viewer);
523: PetscViewerDestroy(&viewer);
524: }
525: incall = PETSC_FALSE;
526: return(0);
527: }
529: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
530: {
532: PetscDraw draw;
533: PetscDrawSP drawsp;
534: PetscReal re,im;
535: PetscInt i,k;
538: if (!eps->nconv) return(0);
539: PetscViewerDrawGetDraw(viewer,0,&draw);
540: PetscDrawSetTitle(draw,"Computed Eigenvalues");
541: PetscDrawSPCreate(draw,1,&drawsp);
542: for (i=0;i<eps->nconv;i++) {
543: k = eps->perm[i];
544: #if defined(PETSC_USE_COMPLEX)
545: re = PetscRealPart(eps->eigr[k]);
546: im = PetscImaginaryPart(eps->eigr[k]);
547: #else
548: re = eps->eigr[k];
549: im = eps->eigi[k];
550: #endif
551: PetscDrawSPAddPoint(drawsp,&re,&im);
552: }
553: PetscDrawSPDraw(drawsp,PETSC_TRUE);
554: PetscDrawSPSave(drawsp);
555: PetscDrawSPDestroy(&drawsp);
556: return(0);
557: }
559: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
560: {
561: PetscInt i,k;
565: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
566: for (i=0;i<eps->nconv;i++) {
567: k = eps->perm[i];
568: PetscViewerASCIIPrintf(viewer," ");
569: SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]);
570: PetscViewerASCIIPrintf(viewer,"\n");
571: }
572: PetscViewerASCIIPrintf(viewer,"\n");
573: return(0);
574: }
576: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
577: {
579: PetscInt i,k;
580: PetscReal re,im;
581: const char *name;
584: PetscObjectGetName((PetscObject)eps,&name);
585: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
586: for (i=0;i<eps->nconv;i++) {
587: k = eps->perm[i];
588: #if defined(PETSC_USE_COMPLEX)
589: re = PetscRealPart(eps->eigr[k]);
590: im = PetscImaginaryPart(eps->eigr[k]);
591: #else
592: re = eps->eigr[k];
593: im = eps->eigi[k];
594: #endif
595: if (im!=0.0) {
596: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
597: } else {
598: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
599: }
600: }
601: PetscViewerASCIIPrintf(viewer,"];\n");
602: return(0);
603: }
605: /*@C
606: EPSValuesView - Displays the computed eigenvalues in a viewer.
608: Collective on eps
610: Input Parameters:
611: + eps - the eigensolver context
612: - viewer - the viewer
614: Options Database Key:
615: . -eps_view_values - print computed eigenvalues
617: Level: intermediate
619: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
620: @*/
621: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
622: {
623: PetscBool isascii,isdraw;
624: PetscViewerFormat format;
625: PetscErrorCode ierr;
629: if (!viewer) {
630: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
631: }
634: EPSCheckSolved(eps,1);
635: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
636: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
637: if (isdraw) {
638: EPSValuesView_DRAW(eps,viewer);
639: } else if (isascii) {
640: PetscViewerGetFormat(viewer,&format);
641: switch (format) {
642: case PETSC_VIEWER_DEFAULT:
643: case PETSC_VIEWER_ASCII_INFO:
644: case PETSC_VIEWER_ASCII_INFO_DETAIL:
645: EPSValuesView_ASCII(eps,viewer);
646: break;
647: case PETSC_VIEWER_ASCII_MATLAB:
648: EPSValuesView_MATLAB(eps,viewer);
649: break;
650: default:
651: PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
652: }
653: }
654: return(0);
655: }
657: /*@
658: EPSValuesViewFromOptions - Processes command line options to determine if/how
659: the computed eigenvalues are to be viewed.
661: Collective on eps
663: Input Parameters:
664: . eps - the eigensolver context
666: Level: developer
667: @*/
668: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
669: {
670: PetscErrorCode ierr;
671: PetscViewer viewer;
672: PetscBool flg;
673: static PetscBool incall = PETSC_FALSE;
674: PetscViewerFormat format;
677: if (incall) return(0);
678: incall = PETSC_TRUE;
679: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg);
680: if (flg) {
681: PetscViewerPushFormat(viewer,format);
682: EPSValuesView(eps,viewer);
683: PetscViewerPopFormat(viewer);
684: PetscViewerDestroy(&viewer);
685: }
686: incall = PETSC_FALSE;
687: return(0);
688: }
690: /*@C
691: EPSVectorsView - Outputs computed eigenvectors to a viewer.
693: Collective on eps
695: Parameter:
696: + eps - the eigensolver context
697: - viewer - the viewer
699: Options Database Keys:
700: . -eps_view_vectors - output eigenvectors.
702: Notes:
703: If PETSc was configured with real scalars, complex conjugate eigenvectors
704: will be viewed as two separate real vectors, one containing the real part
705: and another one containing the imaginary part.
707: If left eigenvectors were computed with a two-sided eigensolver, the right
708: and left eigenvectors are interleaved, that is, the vectors are output in
709: the following order: X0, Y0, X1, Y1, X2, Y2, ...
711: Level: intermediate
713: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
714: @*/
715: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
716: {
718: PetscInt i,k;
719: Vec x;
720: #define NMLEN 30
721: char vname[NMLEN];
722: const char *ename;
726: if (!viewer) {
727: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
728: }
731: EPSCheckSolved(eps,1);
732: if (eps->nconv) {
733: PetscObjectGetName((PetscObject)eps,&ename);
734: EPSComputeVectors(eps);
735: for (i=0;i<eps->nconv;i++) {
736: k = eps->perm[i];
737: PetscSNPrintf(vname,NMLEN,"X%d_%s",(int)i,ename);
738: BVGetColumn(eps->V,k,&x);
739: PetscObjectSetName((PetscObject)x,vname);
740: VecView(x,viewer);
741: BVRestoreColumn(eps->V,k,&x);
742: if (eps->twosided) {
743: PetscSNPrintf(vname,NMLEN,"Y%d_%s",(int)i,ename);
744: BVGetColumn(eps->W,k,&x);
745: PetscObjectSetName((PetscObject)x,vname);
746: VecView(x,viewer);
747: BVRestoreColumn(eps->W,k,&x);
748: }
749: }
750: }
751: return(0);
752: }
754: /*@
755: EPSVectorsViewFromOptions - Processes command line options to determine if/how
756: the computed eigenvectors are to be viewed.
758: Collective on eps
760: Input Parameter:
761: . eps - the eigensolver context
763: Level: developer
764: @*/
765: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
766: {
767: PetscErrorCode ierr;
768: PetscViewer viewer;
769: PetscBool flg = PETSC_FALSE;
770: static PetscBool incall = PETSC_FALSE;
771: PetscViewerFormat format;
774: if (incall) return(0);
775: incall = PETSC_TRUE;
776: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg);
777: if (flg) {
778: PetscViewerPushFormat(viewer,format);
779: EPSVectorsView(eps,viewer);
780: PetscViewerPopFormat(viewer);
781: PetscViewerDestroy(&viewer);
782: }
783: incall = PETSC_FALSE;
784: return(0);
785: }