Actual source code: epsview.c
slepc-3.7.4 2017-05-17
1: /*
2: The EPS routines related to various viewers.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include <slepc/private/epsimpl.h> /*I "slepceps.h" I*/
25: #include <petscdraw.h>
29: /*@C
30: EPSView - Prints the EPS data structure.
32: Collective on EPS
34: Input Parameters:
35: + eps - the eigenproblem solver context
36: - viewer - optional visualization context
38: Options Database Key:
39: . -eps_view - Calls EPSView() at end of EPSSolve()
41: Note:
42: The available visualization contexts include
43: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
44: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
45: output where only the first processor opens
46: the file. All other processors send their
47: data to the first processor to print.
49: The user can open an alternative visualization context with
50: PetscViewerASCIIOpen() - output to a specified file.
52: Level: beginner
54: .seealso: STView(), PetscViewerASCIIOpen()
55: @*/
56: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
57: {
59: const char *type,*extr,*bal;
60: char str[50];
61: PetscBool isascii,ispower,isexternal,istrivial;
65: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
69: #if defined(PETSC_USE_COMPLEX)
70: #define HERM "hermitian"
71: #else
72: #define HERM "symmetric"
73: #endif
74: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
75: if (isascii) {
76: PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer);
77: if (eps->ops->view) {
78: PetscViewerASCIIPushTab(viewer);
79: (*eps->ops->view)(eps,viewer);
80: PetscViewerASCIIPopTab(viewer);
81: }
82: if (eps->problem_type) {
83: switch (eps->problem_type) {
84: case EPS_HEP: type = HERM " eigenvalue problem"; break;
85: case EPS_GHEP: type = "generalized " HERM " eigenvalue problem"; break;
86: case EPS_NHEP: type = "non-" HERM " eigenvalue problem"; break;
87: case EPS_GNHEP: type = "generalized non-" HERM " eigenvalue problem"; break;
88: case EPS_PGNHEP: type = "generalized non-" HERM " eigenvalue problem with " HERM " positive definite B"; break;
89: case EPS_GHIEP: type = "generalized " HERM "-indefinite eigenvalue problem"; break;
90: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->problem_type");
91: }
92: } else type = "not yet set";
93: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
94: if (eps->extraction) {
95: switch (eps->extraction) {
96: case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
97: case EPS_HARMONIC: extr = "harmonic Ritz"; break;
98: case EPS_HARMONIC_RELATIVE:extr = "relative harmonic Ritz"; break;
99: case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
100: case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
101: case EPS_REFINED: extr = "refined Ritz"; break;
102: case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
103: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->extraction");
104: }
105: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr);
106: }
107: if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
108: switch (eps->balance) {
109: case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
110: case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
111: case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
112: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->balance");
113: }
114: PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal);
115: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
116: if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) {
117: PetscViewerASCIIPrintf(viewer,", with its=%D",eps->balance_its);
118: }
119: if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) {
120: PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff);
121: }
122: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
123: PetscViewerASCIIPrintf(viewer,"\n");
124: }
125: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
126: SlepcSNPrintfScalar(str,50,eps->target,PETSC_FALSE);
127: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
128: if (!eps->which) {
129: PetscViewerASCIIPrintf(viewer,"not yet set\n");
130: } else switch (eps->which) {
131: case EPS_WHICH_USER:
132: PetscViewerASCIIPrintf(viewer,"user defined\n");
133: break;
134: case EPS_TARGET_MAGNITUDE:
135: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
136: break;
137: case EPS_TARGET_REAL:
138: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
139: break;
140: case EPS_TARGET_IMAGINARY:
141: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
142: break;
143: case EPS_LARGEST_MAGNITUDE:
144: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
145: break;
146: case EPS_SMALLEST_MAGNITUDE:
147: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
148: break;
149: case EPS_LARGEST_REAL:
150: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
151: break;
152: case EPS_SMALLEST_REAL:
153: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
154: break;
155: case EPS_LARGEST_IMAGINARY:
156: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
157: break;
158: case EPS_SMALLEST_IMAGINARY:
159: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
160: break;
161: case EPS_ALL:
162: if (eps->inta || eps->intb) {
163: PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb);
164: } else {
165: PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
166: }
167: break;
168: default: SETERRQ(PetscObjectComm((PetscObject)eps),1,"Wrong value of eps->which");
169: }
170: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
171: if (eps->isgeneralized && eps->ishermitian && eps->purify) {
172: PetscViewerASCIIPrintf(viewer," postprocessing eigenvectors with purification\n");
173: }
174: if (eps->trueres) {
175: PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n");
176: }
177: if (eps->trackall) {
178: PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n");
179: }
180: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",eps->nev);
181: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",eps->ncv);
182: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",eps->mpd);
183: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",eps->max_it);
184: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)eps->tol);
185: PetscViewerASCIIPrintf(viewer," convergence test: ");
186: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
187: switch (eps->conv) {
188: case EPS_CONV_ABS:
189: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
190: case EPS_CONV_REL:
191: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
192: case EPS_CONV_NORM:
193: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");
194: PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)eps->nrma);
195: if (eps->isgeneralized) {
196: PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb);
197: }
198: PetscViewerASCIIPrintf(viewer,"\n");
199: break;
200: case EPS_CONV_USER:
201: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
202: }
203: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
204: if (eps->nini) {
205: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(eps->nini));
206: }
207: if (eps->nds) {
208: PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %D\n",PetscAbs(eps->nds));
209: }
210: } else {
211: if (eps->ops->view) {
212: (*eps->ops->view)(eps,viewer);
213: }
214: }
215: PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLZPACK,EPSTRLAN,EPSBLOPEX,EPSPRIMME,"");
216: if (!isexternal) {
217: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
218: if (!eps->V) { EPSGetBV(eps,&eps->V); }
219: BVView(eps->V,viewer);
220: if (eps->rg) {
221: RGIsTrivial(eps->rg,&istrivial);
222: if (!istrivial) { RGView(eps->rg,viewer); }
223: }
224: PetscObjectTypeCompare((PetscObject)eps,EPSPOWER,&ispower);
225: if (!ispower) {
226: if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
227: DSView(eps->ds,viewer);
228: }
229: PetscViewerPopFormat(viewer);
230: }
231: if (!eps->st) { EPSGetST(eps,&eps->st); }
232: STView(eps->st,viewer);
233: return(0);
234: }
238: /*@C
239: EPSReasonView - Displays the reason an EPS solve converged or diverged.
241: Collective on EPS
243: Parameter:
244: + eps - the eigensolver context
245: - viewer - the viewer to display the reason
247: Options Database Keys:
248: . -eps_converged_reason - print reason for convergence, and number of iterations
250: Level: intermediate
252: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber()
253: @*/
254: PetscErrorCode EPSReasonView(EPS eps,PetscViewer viewer)
255: {
257: PetscBool isAscii;
260: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
261: if (isAscii) {
262: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
263: if (eps->reason > 0) {
264: 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);
265: } else {
266: 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);
267: }
268: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
269: }
270: return(0);
271: }
275: /*@
276: EPSReasonViewFromOptions - Processes command line options to determine if/how
277: the EPS converged reason is to be viewed.
279: Collective on EPS
281: Input Parameters:
282: . eps - the eigensolver context
284: Level: developer
285: @*/
286: PetscErrorCode EPSReasonViewFromOptions(EPS eps)
287: {
288: PetscErrorCode ierr;
289: PetscViewer viewer;
290: PetscBool flg;
291: static PetscBool incall = PETSC_FALSE;
292: PetscViewerFormat format;
295: if (incall) return(0);
296: incall = PETSC_TRUE;
297: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg);
298: if (flg) {
299: PetscViewerPushFormat(viewer,format);
300: EPSReasonView(eps,viewer);
301: PetscViewerPopFormat(viewer);
302: PetscViewerDestroy(&viewer);
303: }
304: incall = PETSC_FALSE;
305: return(0);
306: }
310: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
311: {
312: PetscBool errok;
313: PetscReal error,re,im;
314: PetscScalar kr,ki;
315: PetscInt i,j,nvals;
319: if (eps->which!=EPS_ALL && eps->nconv<eps->nev) {
320: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
321: return(0);
322: }
323: errok = PETSC_TRUE;
324: nvals = (eps->which==EPS_ALL)? eps->nconv: eps->nev;
325: for (i=0;i<nvals;i++) {
326: EPSComputeError(eps,i,etype,&error);
327: errok = (errok && error<5.0*eps->tol)? PETSC_TRUE: PETSC_FALSE;
328: }
329: if (!errok) {
330: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
331: return(0);
332: }
333: if (eps->which==EPS_ALL) {
334: PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
335: } else {
336: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
337: }
338: for (i=0;i<=(nvals-1)/8;i++) {
339: PetscViewerASCIIPrintf(viewer,"\n ");
340: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
341: EPSGetEigenpair(eps,8*i+j,&kr,&ki,NULL,NULL);
342: #if defined(PETSC_USE_COMPLEX)
343: re = PetscRealPart(kr);
344: im = PetscImaginaryPart(kr);
345: #else
346: re = kr;
347: im = ki;
348: #endif
349: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
350: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
351: if (im!=0.0) {
352: PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
353: } else {
354: PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
355: }
356: if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
357: }
358: }
359: PetscViewerASCIIPrintf(viewer,"\n\n");
360: return(0);
361: }
365: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
366: {
368: PetscReal error,re,im;
369: PetscScalar kr,ki;
370: PetscInt i;
371: #define EXLEN 30
372: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
375: if (!eps->nconv) return(0);
376: switch (etype) {
377: case EPS_ERROR_ABSOLUTE:
378: PetscSNPrintf(ex,EXLEN," ||Ax-k%sx||",eps->isgeneralized?"B":"");
379: break;
380: case EPS_ERROR_RELATIVE:
381: PetscSNPrintf(ex,EXLEN,"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":"");
382: break;
383: case EPS_ERROR_BACKWARD:
384: PetscSNPrintf(ex,EXLEN," eta(x,k)");
385: break;
386: }
387: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
388: for (i=0;i<eps->nconv;i++) {
389: EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
390: EPSComputeError(eps,i,etype,&error);
391: #if defined(PETSC_USE_COMPLEX)
392: re = PetscRealPart(kr);
393: im = PetscImaginaryPart(kr);
394: #else
395: re = kr;
396: im = ki;
397: #endif
398: if (im!=0.0) {
399: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
400: } else {
401: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
402: }
403: }
404: PetscViewerASCIIPrintf(viewer,sep);
405: return(0);
406: }
410: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
411: {
413: PetscReal error;
414: PetscInt i;
415: const char *name;
418: PetscObjectGetName((PetscObject)eps,&name);
419: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
420: for (i=0;i<eps->nconv;i++) {
421: EPSComputeError(eps,i,etype,&error);
422: PetscViewerASCIIPrintf(viewer,"%18.16e\n",error);
423: }
424: PetscViewerASCIIPrintf(viewer,"];\n");
425: return(0);
426: }
430: /*@C
431: EPSErrorView - Displays the errors associated with the computed solution
432: (as well as the eigenvalues).
434: Collective on EPS
436: Input Parameters:
437: + eps - the eigensolver context
438: . etype - error type
439: - viewer - optional visualization context
441: Options Database Key:
442: + -eps_error_absolute - print absolute errors of each eigenpair
443: . -eps_error_relative - print relative errors of each eigenpair
444: - -eps_error_backward - print backward errors of each eigenpair
446: Notes:
447: By default, this function checks the error of all eigenpairs and prints
448: the eigenvalues if all of them are below the requested tolerance.
449: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
450: eigenvalues and corresponding errors is printed.
452: Level: intermediate
454: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
455: @*/
456: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
457: {
458: PetscBool isascii;
459: PetscViewerFormat format;
460: PetscErrorCode ierr;
464: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
467: EPSCheckSolved(eps,1);
468: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
469: if (!isascii) return(0);
471: PetscViewerGetFormat(viewer,&format);
472: switch (format) {
473: case PETSC_VIEWER_DEFAULT:
474: case PETSC_VIEWER_ASCII_INFO:
475: EPSErrorView_ASCII(eps,etype,viewer);
476: break;
477: case PETSC_VIEWER_ASCII_INFO_DETAIL:
478: EPSErrorView_DETAIL(eps,etype,viewer);
479: break;
480: case PETSC_VIEWER_ASCII_MATLAB:
481: EPSErrorView_MATLAB(eps,etype,viewer);
482: break;
483: default:
484: PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
485: }
486: return(0);
487: }
491: /*@
492: EPSErrorViewFromOptions - Processes command line options to determine if/how
493: the errors of the computed solution are to be viewed.
495: Collective on EPS
497: Input Parameters:
498: . eps - the eigensolver context
500: Level: developer
501: @*/
502: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
503: {
504: PetscErrorCode ierr;
505: PetscViewer viewer;
506: PetscBool flg;
507: static PetscBool incall = PETSC_FALSE;
508: PetscViewerFormat format;
511: if (incall) return(0);
512: incall = PETSC_TRUE;
513: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg);
514: if (flg) {
515: PetscViewerPushFormat(viewer,format);
516: EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer);
517: PetscViewerPopFormat(viewer);
518: PetscViewerDestroy(&viewer);
519: }
520: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg);
521: if (flg) {
522: PetscViewerPushFormat(viewer,format);
523: EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer);
524: PetscViewerPopFormat(viewer);
525: PetscViewerDestroy(&viewer);
526: }
527: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg);
528: if (flg) {
529: PetscViewerPushFormat(viewer,format);
530: EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer);
531: PetscViewerPopFormat(viewer);
532: PetscViewerDestroy(&viewer);
533: }
534: incall = PETSC_FALSE;
535: return(0);
536: }
540: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
541: {
543: PetscDraw draw;
544: PetscDrawSP drawsp;
545: PetscReal re,im;
546: PetscInt i,k;
549: if (!eps->nconv) return(0);
550: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Computed Eigenvalues",PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
551: PetscViewerDrawGetDraw(viewer,0,&draw);
552: PetscDrawSPCreate(draw,1,&drawsp);
553: for (i=0;i<eps->nconv;i++) {
554: k = eps->perm[i];
555: #if defined(PETSC_USE_COMPLEX)
556: re = PetscRealPart(eps->eigr[k]);
557: im = PetscImaginaryPart(eps->eigr[k]);
558: #else
559: re = eps->eigr[k];
560: im = eps->eigi[k];
561: #endif
562: PetscDrawSPAddPoint(drawsp,&re,&im);
563: }
564: PetscDrawSPDraw(drawsp,PETSC_TRUE);
565: PetscDrawSPSave(drawsp);
566: PetscDrawSPDestroy(&drawsp);
567: PetscViewerDestroy(&viewer);
568: return(0);
569: }
573: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
574: {
575: PetscReal re,im;
576: PetscInt i,k;
580: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
581: for (i=0;i<eps->nconv;i++) {
582: k = eps->perm[i];
583: #if defined(PETSC_USE_COMPLEX)
584: re = PetscRealPart(eps->eigr[k]);
585: im = PetscImaginaryPart(eps->eigr[k]);
586: #else
587: re = eps->eigr[k];
588: im = eps->eigi[k];
589: #endif
590: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
591: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
592: if (im!=0.0) {
593: PetscViewerASCIIPrintf(viewer," %.5f%+.5fi\n",(double)re,(double)im);
594: } else {
595: PetscViewerASCIIPrintf(viewer," %.5f\n",(double)re);
596: }
597: }
598: PetscViewerASCIIPrintf(viewer,"\n");
599: return(0);
600: }
604: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
605: {
607: PetscInt i,k;
608: PetscReal re,im;
609: const char *name;
612: PetscObjectGetName((PetscObject)eps,&name);
613: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
614: for (i=0;i<eps->nconv;i++) {
615: k = eps->perm[i];
616: #if defined(PETSC_USE_COMPLEX)
617: re = PetscRealPart(eps->eigr[k]);
618: im = PetscImaginaryPart(eps->eigr[k]);
619: #else
620: re = eps->eigr[k];
621: im = eps->eigi[k];
622: #endif
623: if (im!=0.0) {
624: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
625: } else {
626: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
627: }
628: }
629: PetscViewerASCIIPrintf(viewer,"];\n");
630: return(0);
631: }
635: /*@C
636: EPSValuesView - Displays the computed eigenvalues in a viewer.
638: Collective on EPS
640: Input Parameters:
641: + eps - the eigensolver context
642: - viewer - the viewer
644: Options Database Key:
645: . -eps_view_values - print computed eigenvalues
647: Level: intermediate
649: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
650: @*/
651: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
652: {
653: PetscBool isascii,isdraw;
654: PetscViewerFormat format;
655: PetscErrorCode ierr;
659: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
662: EPSCheckSolved(eps,1);
663: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
664: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
665: if (isdraw) {
666: EPSValuesView_DRAW(eps,viewer);
667: } else if (isascii) {
668: PetscViewerGetFormat(viewer,&format);
669: switch (format) {
670: case PETSC_VIEWER_DEFAULT:
671: case PETSC_VIEWER_ASCII_INFO:
672: case PETSC_VIEWER_ASCII_INFO_DETAIL:
673: EPSValuesView_ASCII(eps,viewer);
674: break;
675: case PETSC_VIEWER_ASCII_MATLAB:
676: EPSValuesView_MATLAB(eps,viewer);
677: break;
678: default:
679: PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
680: }
681: }
682: return(0);
683: }
687: /*@
688: EPSValuesViewFromOptions - Processes command line options to determine if/how
689: the computed eigenvalues are to be viewed.
691: Collective on EPS
693: Input Parameters:
694: . eps - the eigensolver context
696: Level: developer
697: @*/
698: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
699: {
700: PetscErrorCode ierr;
701: PetscViewer viewer;
702: PetscBool flg;
703: static PetscBool incall = PETSC_FALSE;
704: PetscViewerFormat format;
707: if (incall) return(0);
708: incall = PETSC_TRUE;
709: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg);
710: if (flg) {
711: PetscViewerPushFormat(viewer,format);
712: EPSValuesView(eps,viewer);
713: PetscViewerPopFormat(viewer);
714: PetscViewerDestroy(&viewer);
715: }
716: incall = PETSC_FALSE;
717: return(0);
718: }
722: /*@C
723: EPSVectorsView - Outputs computed eigenvectors to a viewer.
725: Collective on EPS
727: Parameter:
728: + eps - the eigensolver context
729: - viewer - the viewer
731: Options Database Keys:
732: . -eps_view_vectors - output eigenvectors.
734: Note:
735: If PETSc was configured with real scalars, complex conjugate eigenvectors
736: will be viewed as two separate real vectors, one containing the real part
737: and another one containing the imaginary part.
739: Level: intermediate
741: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
742: @*/
743: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
744: {
746: PetscInt i,k;
747: Vec x;
748: #define NMLEN 30
749: char vname[NMLEN];
750: const char *ename;
754: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
757: EPSCheckSolved(eps,1);
758: if (eps->nconv) {
759: PetscObjectGetName((PetscObject)eps,&ename);
760: EPSComputeVectors(eps);
761: for (i=0;i<eps->nconv;i++) {
762: k = eps->perm[i];
763: PetscSNPrintf(vname,NMLEN,"V%d_%s",(int)i,ename);
764: BVGetColumn(eps->V,k,&x);
765: PetscObjectSetName((PetscObject)x,vname);
766: VecView(x,viewer);
767: BVRestoreColumn(eps->V,k,&x);
768: }
769: }
770: return(0);
771: }
775: /*@
776: EPSVectorsViewFromOptions - Processes command line options to determine if/how
777: the computed eigenvectors are to be viewed.
779: Collective on EPS
781: Input Parameters:
782: . eps - the eigensolver context
784: Level: developer
785: @*/
786: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
787: {
788: PetscErrorCode ierr;
789: PetscViewer viewer;
790: PetscBool flg = PETSC_FALSE;
791: static PetscBool incall = PETSC_FALSE;
792: PetscViewerFormat format;
795: if (incall) return(0);
796: incall = PETSC_TRUE;
797: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg);
798: if (flg) {
799: PetscViewerPushFormat(viewer,format);
800: EPSVectorsView(eps,viewer);
801: PetscViewerPopFormat(viewer);
802: PetscViewerDestroy(&viewer);
803: }
804: incall = PETSC_FALSE;
805: return(0);
806: }