Actual source code: pepview.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: The PEP routines related to various viewers
12: */
14: #include <slepc/private/pepimpl.h> /*I "slepcpep.h" I*/
15: #include <petscdraw.h>
17: /*@C
18: PEPView - Prints the PEP data structure.
20: Collective on pep
22: Input Parameters:
23: + pep - the polynomial eigenproblem solver context
24: - viewer - optional visualization context
26: Options Database Key:
27: . -pep_view - Calls PEPView() at end of PEPSolve()
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
41: @*/
42: PetscErrorCode PEPView(PEP pep,PetscViewer viewer)
43: {
45: const char *type=NULL;
46: char str[50];
47: PetscBool isascii,islinear,istrivial;
48: PetscInt i;
49: PetscViewer sviewer;
53: if (!viewer) {
54: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&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)pep,viewer);
67: if (pep->ops->view) {
68: PetscViewerASCIIPushTab(viewer);
69: (*pep->ops->view)(pep,viewer);
70: PetscViewerASCIIPopTab(viewer);
71: }
72: if (pep->problem_type) {
73: switch (pep->problem_type) {
74: case PEP_GENERAL: type = "general polynomial eigenvalue problem"; break;
75: case PEP_HERMITIAN: type = HERM " polynomial eigenvalue problem"; break;
76: case PEP_HYPERBOLIC: type = "hyperbolic polynomial eigenvalue problem"; break;
77: case PEP_GYROSCOPIC: type = "gyroscopic polynomial eigenvalue problem"; break;
78: }
79: } else type = "not yet set";
80: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
81: PetscViewerASCIIPrintf(viewer," polynomial represented in %s basis\n",PEPBasisTypes[pep->basis]);
82: switch (pep->scale) {
83: case PEP_SCALE_NONE:
84: break;
85: case PEP_SCALE_SCALAR:
86: PetscViewerASCIIPrintf(viewer," parameter scaling enabled, with scaling factor=%g\n",(double)pep->sfactor);
87: break;
88: case PEP_SCALE_DIAGONAL:
89: PetscViewerASCIIPrintf(viewer," diagonal balancing enabled, with its=%D and lambda=%g\n",pep->sits,(double)pep->slambda);
90: break;
91: case PEP_SCALE_BOTH:
92: PetscViewerASCIIPrintf(viewer," parameter scaling & diagonal balancing enabled, with scaling factor=%g, its=%D and lambda=%g\n",(double)pep->sfactor,pep->sits,(double)pep->slambda);
93: break;
94: }
95: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
96: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
97: SlepcSNPrintfScalar(str,50,pep->target,PETSC_FALSE);
98: if (!pep->which) {
99: PetscViewerASCIIPrintf(viewer,"not yet set\n");
100: } else switch (pep->which) {
101: case PEP_WHICH_USER:
102: PetscViewerASCIIPrintf(viewer,"user defined\n");
103: break;
104: case PEP_TARGET_MAGNITUDE:
105: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
106: break;
107: case PEP_TARGET_REAL:
108: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
109: break;
110: case PEP_TARGET_IMAGINARY:
111: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
112: break;
113: case PEP_LARGEST_MAGNITUDE:
114: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
115: break;
116: case PEP_SMALLEST_MAGNITUDE:
117: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
118: break;
119: case PEP_LARGEST_REAL:
120: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
121: break;
122: case PEP_SMALLEST_REAL:
123: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
124: break;
125: case PEP_LARGEST_IMAGINARY:
126: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
127: break;
128: case PEP_SMALLEST_IMAGINARY:
129: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
130: break;
131: case PEP_ALL:
132: if (pep->inta || pep->intb) {
133: PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)pep->inta,(double)pep->intb);
134: } else {
135: PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
136: }
137: break;
138: }
139: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
140: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",pep->nev);
141: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",pep->ncv);
142: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",pep->mpd);
143: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",pep->max_it);
144: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)pep->tol);
145: PetscViewerASCIIPrintf(viewer," convergence test: ");
146: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
147: switch (pep->conv) {
148: case PEP_CONV_ABS:
149: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
150: case PEP_CONV_REL:
151: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
152: case PEP_CONV_NORM:
153: PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
154: if (pep->nrma) {
155: PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)pep->nrma[0]);
156: for (i=1;i<pep->nmat;i++) {
157: PetscViewerASCIIPrintf(viewer,", %g",(double)pep->nrma[i]);
158: }
159: PetscViewerASCIIPrintf(viewer,"\n");
160: }
161: break;
162: case PEP_CONV_USER:
163: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
164: }
165: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
166: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",PEPExtractTypes[pep->extract]);
167: if (pep->refine) {
168: PetscViewerASCIIPrintf(viewer," iterative refinement: %s, with %s scheme\n",PEPRefineTypes[pep->refine],PEPRefineSchemes[pep->scheme]);
169: PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%D\n",(double)pep->rtol,pep->rits);
170: if (pep->npart>1) {
171: PetscViewerASCIIPrintf(viewer," splitting communicator in %D partitions for refinement\n",pep->npart);
172: }
173: }
174: if (pep->nini) {
175: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(pep->nini));
176: }
177: } else {
178: if (pep->ops->view) {
179: (*pep->ops->view)(pep,viewer);
180: }
181: }
182: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
183: if (!pep->V) { PEPGetBV(pep,&pep->V); }
184: BVView(pep->V,viewer);
185: if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
186: RGIsTrivial(pep->rg,&istrivial);
187: if (!istrivial) { RGView(pep->rg,viewer); }
188: PetscObjectTypeCompare((PetscObject)pep,PEPLINEAR,&islinear);
189: if (!islinear) {
190: if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
191: DSView(pep->ds,viewer);
192: }
193: PetscViewerPopFormat(viewer);
194: if (!pep->st) { PEPGetST(pep,&pep->st); }
195: STView(pep->st,viewer);
196: if (pep->refine!=PEP_REFINE_NONE) {
197: PetscViewerASCIIPushTab(viewer);
198: if (pep->npart>1) {
199: if (pep->refinesubc->color==0) {
200: PetscViewerASCIIGetStdout(PetscSubcommChild(pep->refinesubc),&sviewer);
201: KSPView(pep->refineksp,sviewer);
202: }
203: } else {
204: KSPView(pep->refineksp,viewer);
205: }
206: PetscViewerASCIIPopTab(viewer);
207: }
208: return(0);
209: }
211: /*@C
212: PEPViewFromOptions - View from options
214: Collective on PEP
216: Input Parameters:
217: + pep - the eigensolver context
218: . obj - optional object
219: - name - command line option
221: Level: intermediate
223: .seealso: PEPView(), PEPCreate()
224: @*/
225: PetscErrorCode PEPViewFromOptions(PEP pep,PetscObject obj,const char name[])
226: {
231: PetscObjectViewFromOptions((PetscObject)pep,obj,name);
232: return(0);
233: }
235: /*@C
236: PEPReasonView - Displays the reason a PEP solve converged or diverged.
238: Collective on pep
240: Input Parameters:
241: + pep - the eigensolver context
242: - viewer - the viewer to display the reason
244: Options Database Keys:
245: . -pep_converged_reason - print reason for convergence, and number of iterations
247: Level: intermediate
249: .seealso: PEPSetConvergenceTest(), PEPSetTolerances(), PEPGetIterationNumber()
250: @*/
251: PetscErrorCode PEPReasonView(PEP pep,PetscViewer viewer)
252: {
254: PetscBool isAscii;
257: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
258: if (isAscii) {
259: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
260: if (pep->reason > 0) {
261: PetscViewerASCIIPrintf(viewer,"%s Polynomial eigensolve converged (%D eigenpair%s) due to %s; iterations %D\n",((PetscObject)pep)->prefix?((PetscObject)pep)->prefix:"",pep->nconv,(pep->nconv>1)?"s":"",PEPConvergedReasons[pep->reason],pep->its);
262: } else {
263: PetscViewerASCIIPrintf(viewer,"%s Polynomial eigensolve did not converge due to %s; iterations %D\n",((PetscObject)pep)->prefix?((PetscObject)pep)->prefix:"",PEPConvergedReasons[pep->reason],pep->its);
264: }
265: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
266: }
267: return(0);
268: }
270: /*@
271: PEPReasonViewFromOptions - Processes command line options to determine if/how
272: the PEP converged reason is to be viewed.
274: Collective on pep
276: Input Parameter:
277: . pep - the eigensolver context
279: Level: developer
280: @*/
281: PetscErrorCode PEPReasonViewFromOptions(PEP pep)
282: {
283: PetscErrorCode ierr;
284: PetscViewer viewer;
285: PetscBool flg;
286: static PetscBool incall = PETSC_FALSE;
287: PetscViewerFormat format;
290: if (incall) return(0);
291: incall = PETSC_TRUE;
292: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_converged_reason",&viewer,&format,&flg);
293: if (flg) {
294: PetscViewerPushFormat(viewer,format);
295: PEPReasonView(pep,viewer);
296: PetscViewerPopFormat(viewer);
297: PetscViewerDestroy(&viewer);
298: }
299: incall = PETSC_FALSE;
300: return(0);
301: }
303: static PetscErrorCode PEPErrorView_ASCII(PEP pep,PEPErrorType etype,PetscViewer viewer)
304: {
305: PetscReal error;
306: PetscInt i,j,k;
310: if (pep->nconv<pep->nev) {
311: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",pep->nev);
312: return(0);
313: }
314: for (i=0;i<pep->nev;i++) {
315: PEPComputeError(pep,i,etype,&error);
316: if (error>=5.0*pep->tol) {
317: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",pep->nev);
318: return(0);
319: }
320: }
321: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
322: for (i=0;i<=(pep->nev-1)/8;i++) {
323: PetscViewerASCIIPrintf(viewer,"\n ");
324: for (j=0;j<PetscMin(8,pep->nev-8*i);j++) {
325: k = pep->perm[8*i+j];
326: SlepcPrintEigenvalueASCII(viewer,pep->eigr[k],pep->eigi[k]);
327: if (8*i+j+1<pep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
328: }
329: }
330: PetscViewerASCIIPrintf(viewer,"\n\n");
331: return(0);
332: }
334: static PetscErrorCode PEPErrorView_DETAIL(PEP pep,PEPErrorType etype,PetscViewer viewer)
335: {
337: PetscReal error,re,im;
338: PetscScalar kr,ki;
339: PetscInt i;
340: #define EXLEN 30
341: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
344: if (!pep->nconv) return(0);
345: switch (etype) {
346: case PEP_ERROR_ABSOLUTE:
347: PetscSNPrintf(ex,EXLEN," ||P(k)x||");
348: break;
349: case PEP_ERROR_RELATIVE:
350: PetscSNPrintf(ex,EXLEN,"||P(k)x||/||kx||");
351: break;
352: case PEP_ERROR_BACKWARD:
353: PetscSNPrintf(ex,EXLEN," eta(x,k)");
354: break;
355: }
356: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
357: for (i=0;i<pep->nconv;i++) {
358: PEPGetEigenpair(pep,i,&kr,&ki,NULL,NULL);
359: PEPComputeError(pep,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) {
368: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
369: } else {
370: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
371: }
372: }
373: PetscViewerASCIIPrintf(viewer,sep);
374: return(0);
375: }
377: static PetscErrorCode PEPErrorView_MATLAB(PEP pep,PEPErrorType etype,PetscViewer viewer)
378: {
380: PetscReal error;
381: PetscInt i;
382: const char *name;
385: PetscObjectGetName((PetscObject)pep,&name);
386: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
387: for (i=0;i<pep->nconv;i++) {
388: PEPComputeError(pep,i,etype,&error);
389: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
390: }
391: PetscViewerASCIIPrintf(viewer,"];\n");
392: return(0);
393: }
395: /*@C
396: PEPErrorView - Displays the errors associated with the computed solution
397: (as well as the eigenvalues).
399: Collective on pep
401: Input Parameters:
402: + pep - the eigensolver context
403: . etype - error type
404: - viewer - optional visualization context
406: Options Database Key:
407: + -pep_error_absolute - print absolute errors of each eigenpair
408: . -pep_error_relative - print relative errors of each eigenpair
409: - -pep_error_backward - print backward errors of each eigenpair
411: Notes:
412: By default, this function checks the error of all eigenpairs and prints
413: the eigenvalues if all of them are below the requested tolerance.
414: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
415: eigenvalues and corresponding errors is printed.
417: Level: intermediate
419: .seealso: PEPSolve(), PEPValuesView(), PEPVectorsView()
420: @*/
421: PetscErrorCode PEPErrorView(PEP pep,PEPErrorType etype,PetscViewer viewer)
422: {
423: PetscBool isascii;
424: PetscViewerFormat format;
425: PetscErrorCode ierr;
429: if (!viewer) {
430: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
431: }
434: PEPCheckSolved(pep,1);
435: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
436: if (!isascii) return(0);
438: PetscViewerGetFormat(viewer,&format);
439: switch (format) {
440: case PETSC_VIEWER_DEFAULT:
441: case PETSC_VIEWER_ASCII_INFO:
442: PEPErrorView_ASCII(pep,etype,viewer);
443: break;
444: case PETSC_VIEWER_ASCII_INFO_DETAIL:
445: PEPErrorView_DETAIL(pep,etype,viewer);
446: break;
447: case PETSC_VIEWER_ASCII_MATLAB:
448: PEPErrorView_MATLAB(pep,etype,viewer);
449: break;
450: default:
451: PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
452: }
453: return(0);
454: }
456: /*@
457: PEPErrorViewFromOptions - Processes command line options to determine if/how
458: the errors of the computed solution are to be viewed.
460: Collective on pep
462: Input Parameter:
463: . pep - the eigensolver context
465: Level: developer
466: @*/
467: PetscErrorCode PEPErrorViewFromOptions(PEP pep)
468: {
469: PetscErrorCode ierr;
470: PetscViewer viewer;
471: PetscBool flg;
472: static PetscBool incall = PETSC_FALSE;
473: PetscViewerFormat format;
476: if (incall) return(0);
477: incall = PETSC_TRUE;
478: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_error_absolute",&viewer,&format,&flg);
479: if (flg) {
480: PetscViewerPushFormat(viewer,format);
481: PEPErrorView(pep,PEP_ERROR_ABSOLUTE,viewer);
482: PetscViewerPopFormat(viewer);
483: PetscViewerDestroy(&viewer);
484: }
485: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_error_relative",&viewer,&format,&flg);
486: if (flg) {
487: PetscViewerPushFormat(viewer,format);
488: PEPErrorView(pep,PEP_ERROR_RELATIVE,viewer);
489: PetscViewerPopFormat(viewer);
490: PetscViewerDestroy(&viewer);
491: }
492: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_error_backward",&viewer,&format,&flg);
493: if (flg) {
494: PetscViewerPushFormat(viewer,format);
495: PEPErrorView(pep,PEP_ERROR_BACKWARD,viewer);
496: PetscViewerPopFormat(viewer);
497: PetscViewerDestroy(&viewer);
498: }
499: incall = PETSC_FALSE;
500: return(0);
501: }
503: static PetscErrorCode PEPValuesView_DRAW(PEP pep,PetscViewer viewer)
504: {
506: PetscDraw draw;
507: PetscDrawSP drawsp;
508: PetscReal re,im;
509: PetscInt i,k;
512: if (!pep->nconv) return(0);
513: PetscViewerDrawGetDraw(viewer,0,&draw);
514: PetscDrawSetTitle(draw,"Computed Eigenvalues");
515: PetscDrawSPCreate(draw,1,&drawsp);
516: for (i=0;i<pep->nconv;i++) {
517: k = pep->perm[i];
518: #if defined(PETSC_USE_COMPLEX)
519: re = PetscRealPart(pep->eigr[k]);
520: im = PetscImaginaryPart(pep->eigr[k]);
521: #else
522: re = pep->eigr[k];
523: im = pep->eigi[k];
524: #endif
525: PetscDrawSPAddPoint(drawsp,&re,&im);
526: }
527: PetscDrawSPDraw(drawsp,PETSC_TRUE);
528: PetscDrawSPSave(drawsp);
529: PetscDrawSPDestroy(&drawsp);
530: return(0);
531: }
533: static PetscErrorCode PEPValuesView_ASCII(PEP pep,PetscViewer viewer)
534: {
535: PetscInt i,k;
539: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
540: for (i=0;i<pep->nconv;i++) {
541: k = pep->perm[i];
542: PetscViewerASCIIPrintf(viewer," ");
543: SlepcPrintEigenvalueASCII(viewer,pep->eigr[k],pep->eigi[k]);
544: PetscViewerASCIIPrintf(viewer,"\n");
545: }
546: PetscViewerASCIIPrintf(viewer,"\n");
547: return(0);
548: }
550: static PetscErrorCode PEPValuesView_MATLAB(PEP pep,PetscViewer viewer)
551: {
553: PetscInt i,k;
554: PetscReal re,im;
555: const char *name;
558: PetscObjectGetName((PetscObject)pep,&name);
559: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
560: for (i=0;i<pep->nconv;i++) {
561: k = pep->perm[i];
562: #if defined(PETSC_USE_COMPLEX)
563: re = PetscRealPart(pep->eigr[k]);
564: im = PetscImaginaryPart(pep->eigr[k]);
565: #else
566: re = pep->eigr[k];
567: im = pep->eigi[k];
568: #endif
569: if (im!=0.0) {
570: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
571: } else {
572: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
573: }
574: }
575: PetscViewerASCIIPrintf(viewer,"];\n");
576: return(0);
577: }
579: /*@C
580: PEPValuesView - Displays the computed eigenvalues in a viewer.
582: Collective on pep
584: Input Parameters:
585: + pep - the eigensolver context
586: - viewer - the viewer
588: Options Database Key:
589: . -pep_view_values - print computed eigenvalues
591: Level: intermediate
593: .seealso: PEPSolve(), PEPVectorsView(), PEPErrorView()
594: @*/
595: PetscErrorCode PEPValuesView(PEP pep,PetscViewer viewer)
596: {
597: PetscBool isascii,isdraw;
598: PetscViewerFormat format;
599: PetscErrorCode ierr;
603: if (!viewer) {
604: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
605: }
608: PEPCheckSolved(pep,1);
609: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
610: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
611: if (isdraw) {
612: PEPValuesView_DRAW(pep,viewer);
613: } else if (isascii) {
614: PetscViewerGetFormat(viewer,&format);
615: switch (format) {
616: case PETSC_VIEWER_DEFAULT:
617: case PETSC_VIEWER_ASCII_INFO:
618: case PETSC_VIEWER_ASCII_INFO_DETAIL:
619: PEPValuesView_ASCII(pep,viewer);
620: break;
621: case PETSC_VIEWER_ASCII_MATLAB:
622: PEPValuesView_MATLAB(pep,viewer);
623: break;
624: default:
625: PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
626: }
627: }
628: return(0);
629: }
631: /*@
632: PEPValuesViewFromOptions - Processes command line options to determine if/how
633: the computed eigenvalues are to be viewed.
635: Collective on pep
637: Input Parameter:
638: . pep - the eigensolver context
640: Level: developer
641: @*/
642: PetscErrorCode PEPValuesViewFromOptions(PEP pep)
643: {
644: PetscErrorCode ierr;
645: PetscViewer viewer;
646: PetscBool flg;
647: static PetscBool incall = PETSC_FALSE;
648: PetscViewerFormat format;
651: if (incall) return(0);
652: incall = PETSC_TRUE;
653: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_view_values",&viewer,&format,&flg);
654: if (flg) {
655: PetscViewerPushFormat(viewer,format);
656: PEPValuesView(pep,viewer);
657: PetscViewerPopFormat(viewer);
658: PetscViewerDestroy(&viewer);
659: }
660: incall = PETSC_FALSE;
661: return(0);
662: }
664: /*@C
665: PEPVectorsView - Outputs computed eigenvectors to a viewer.
667: Collective on pep
669: Input Parameters:
670: + pep - the eigensolver context
671: - viewer - the viewer
673: Options Database Keys:
674: . -pep_view_vectors - output eigenvectors.
676: Note:
677: If PETSc was configured with real scalars, complex conjugate eigenvectors
678: will be viewed as two separate real vectors, one containing the real part
679: and another one containing the imaginary part.
681: Level: intermediate
683: .seealso: PEPSolve(), PEPValuesView(), PEPErrorView()
684: @*/
685: PetscErrorCode PEPVectorsView(PEP pep,PetscViewer viewer)
686: {
688: PetscInt i,k;
689: Vec x;
690: #define NMLEN 30
691: char vname[NMLEN];
692: const char *ename;
696: if (!viewer) {
697: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
698: }
701: PEPCheckSolved(pep,1);
702: if (pep->nconv) {
703: PetscObjectGetName((PetscObject)pep,&ename);
704: PEPComputeVectors(pep);
705: for (i=0;i<pep->nconv;i++) {
706: k = pep->perm[i];
707: PetscSNPrintf(vname,NMLEN,"V%d_%s",(int)i,ename);
708: BVGetColumn(pep->V,k,&x);
709: PetscObjectSetName((PetscObject)x,vname);
710: VecView(x,viewer);
711: BVRestoreColumn(pep->V,k,&x);
712: }
713: }
714: return(0);
715: }
717: /*@
718: PEPVectorsViewFromOptions - Processes command line options to determine if/how
719: the computed eigenvectors are to be viewed.
721: Collective on pep
723: Input Parameter:
724: . pep - the eigensolver context
726: Level: developer
727: @*/
728: PetscErrorCode PEPVectorsViewFromOptions(PEP pep)
729: {
730: PetscErrorCode ierr;
731: PetscViewer viewer;
732: PetscBool flg = PETSC_FALSE;
733: static PetscBool incall = PETSC_FALSE;
734: PetscViewerFormat format;
737: if (incall) return(0);
738: incall = PETSC_TRUE;
739: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_view_vectors",&viewer,&format,&flg);
740: if (flg) {
741: PetscViewerPushFormat(viewer,format);
742: PEPVectorsView(pep,viewer);
743: PetscViewerPopFormat(viewer);
744: PetscViewerDestroy(&viewer);
745: }
746: incall = PETSC_FALSE;
747: return(0);
748: }