Actual source code: pepview.c
slepc-3.6.3 2016-03-29
1: /*
2: The PEP routines related to various viewers.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2015, 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/pepimpl.h> /*I "slepcpep.h" I*/
25: #include <petscdraw.h>
29: /*@C
30: PEPView - Prints the PEP data structure.
32: Collective on PEP
34: Input Parameters:
35: + pep - the polynomial eigenproblem solver context
36: - viewer - optional visualization context
38: Options Database Key:
39: . -pep_view - Calls PEPView() at end of PEPSolve()
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: PetscViewerASCIIOpen()
55: @*/
56: PetscErrorCode PEPView(PEP pep,PetscViewer viewer)
57: {
59: const char *type;
60: char str[50];
61: PetscBool isascii,islinear,istrivial;
62: PetscInt i;
63: PetscViewer sviewer;
67: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
71: #if defined(PETSC_USE_COMPLEX)
72: #define HERM "hermitian"
73: #else
74: #define HERM "symmetric"
75: #endif
76: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
77: if (isascii) {
78: PetscObjectPrintClassNamePrefixType((PetscObject)pep,viewer);
79: if (pep->ops->view) {
80: PetscViewerASCIIPushTab(viewer);
81: (*pep->ops->view)(pep,viewer);
82: PetscViewerASCIIPopTab(viewer);
83: }
84: if (pep->problem_type) {
85: switch (pep->problem_type) {
86: case PEP_GENERAL: type = "general polynomial eigenvalue problem"; break;
87: case PEP_HERMITIAN: type = HERM " polynomial eigenvalue problem"; break;
88: case PEP_GYROSCOPIC: type = "gyroscopic polynomial eigenvalue problem"; break;
89: default: SETERRQ(PetscObjectComm((PetscObject)pep),1,"Wrong value of pep->problem_type");
90: }
91: } else type = "not yet set";
92: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
93: PetscViewerASCIIPrintf(viewer," polynomial represented in %s basis\n",PEPBasisTypes[pep->basis]);
94: switch (pep->scale) {
95: case PEP_SCALE_NONE:
96: break;
97: case PEP_SCALE_SCALAR:
98: PetscViewerASCIIPrintf(viewer," scalar balancing enabled, with scaling factor=%g\n",(double)pep->sfactor);
99: break;
100: case PEP_SCALE_DIAGONAL:
101: PetscViewerASCIIPrintf(viewer," diagonal balancing enabled, with its=%D and lambda=%g\n",pep->sits,(double)pep->slambda);
102: break;
103: case PEP_SCALE_BOTH:
104: PetscViewerASCIIPrintf(viewer," scalar & diagonal balancing enabled, with scaling factor=%g, its=%D and lambda=%g\n",(double)pep->sfactor,pep->sits,(double)pep->slambda);
105: break;
106: }
107: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",PEPExtractTypes[pep->extract]);
108: PetscViewerASCIIPrintf(viewer," iterative refinement: %s%s\n",PEPRefineTypes[pep->refine],pep->schur?", with a Schur complement approach":"");
109: if (pep->refine) {
110: PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%D\n",(double)pep->rtol,pep->rits);
111: if (pep->npart>1) {
112: PetscViewerASCIIPrintf(viewer," splitting communicator in %D partitions for refinement\n",pep->npart);
113: }
114: }
115: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
116: SlepcSNPrintfScalar(str,50,pep->target,PETSC_FALSE);
117: if (!pep->which) {
118: PetscViewerASCIIPrintf(viewer,"not yet set\n");
119: } else switch (pep->which) {
120: case PEP_WHICH_USER:
121: PetscViewerASCIIPrintf(viewer,"user defined\n");
122: break;
123: case PEP_TARGET_MAGNITUDE:
124: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
125: break;
126: case PEP_TARGET_REAL:
127: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
128: break;
129: case PEP_TARGET_IMAGINARY:
130: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
131: break;
132: case PEP_LARGEST_MAGNITUDE:
133: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
134: break;
135: case PEP_SMALLEST_MAGNITUDE:
136: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
137: break;
138: case PEP_LARGEST_REAL:
139: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
140: break;
141: case PEP_SMALLEST_REAL:
142: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
143: break;
144: case PEP_LARGEST_IMAGINARY:
145: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
146: break;
147: case PEP_SMALLEST_IMAGINARY:
148: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
149: break;
150: default: SETERRQ(PetscObjectComm((PetscObject)pep),1,"Wrong value of pep->which");
151: }
152: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",pep->nev);
153: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",pep->ncv);
154: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",pep->mpd);
155: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",pep->max_it);
156: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)pep->tol);
157: PetscViewerASCIIPrintf(viewer," convergence test: ");
158: switch (pep->conv) {
159: case PEP_CONV_ABS:
160: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
161: case PEP_CONV_EIG:
162: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
163: case PEP_CONV_LINEAR:
164: PetscViewerASCIIPrintf(viewer,"related to the linearized eigenproblem\n");
165: if (pep->nrma) {
166: PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)pep->nrma[0]);
167: for (i=1;i<pep->nmat;i++) {
168: PetscViewerASCIIPrintf(viewer,", %g",(double)pep->nrma[i]);
169: }
170: PetscViewerASCIIPrintf(viewer,"\n");
171: }
172: break;
173: case PEP_CONV_NORM:
174: PetscViewerASCIIPrintf(viewer,"related to the matrix norms\n");
175: if (pep->nrma) {
176: PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)pep->nrma[0]);
177: for (i=1;i<pep->nmat;i++) {
178: PetscViewerASCIIPrintf(viewer,", %g",(double)pep->nrma[i]);
179: }
180: PetscViewerASCIIPrintf(viewer,"\n");
181: }
182: break;
183: case PEP_CONV_USER:
184: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
185: }
186: if (pep->nini) {
187: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(pep->nini));
188: }
189: } else {
190: if (pep->ops->view) {
191: (*pep->ops->view)(pep,viewer);
192: }
193: }
194: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
195: if (!pep->V) { PEPGetBV(pep,&pep->V); }
196: BVView(pep->V,viewer);
197: if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
198: RGIsTrivial(pep->rg,&istrivial);
199: if (!istrivial) { RGView(pep->rg,viewer); }
200: PetscObjectTypeCompare((PetscObject)pep,PEPLINEAR,&islinear);
201: if (!islinear) {
202: if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
203: DSView(pep->ds,viewer);
204: }
205: PetscViewerPopFormat(viewer);
206: if (!pep->st) { PEPGetST(pep,&pep->st); }
207: STView(pep->st,viewer);
208: if (pep->refine!=PEP_REFINE_NONE) {
209: if (pep->npart>1) {
210: if (pep->refinesubc->color==0) {
211: PetscViewerASCIIGetStdout(PetscSubcommChild(pep->refinesubc),&sviewer);
212: KSPView(pep->refineksp,sviewer);
213: }
214: } else {
215: KSPView(pep->refineksp,viewer);
216: }
217: }
218: return(0);
219: }
223: /*@C
224: PEPReasonView - Displays the reason a PEP solve converged or diverged.
226: Collective on PEP
228: Parameter:
229: + pep - the eigensolver context
230: - viewer - the viewer to display the reason
232: Options Database Keys:
233: . -pep_converged_reason - print reason for convergence, and number of iterations
235: Level: intermediate
237: .seealso: PEPSetConvergenceTest(), PEPSetTolerances(), PEPGetIterationNumber()
238: @*/
239: PetscErrorCode PEPReasonView(PEP pep,PetscViewer viewer)
240: {
242: PetscBool isAscii;
245: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
246: if (isAscii) {
247: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
248: if (pep->reason > 0) {
249: 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);
250: } else {
251: 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);
252: }
253: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
254: }
255: return(0);
256: }
260: /*@
261: PEPReasonViewFromOptions - Processes command line options to determine if/how
262: the PEP converged reason is to be viewed.
264: Collective on PEP
266: Input Parameters:
267: . pep - the eigensolver context
269: Level: developer
270: @*/
271: PetscErrorCode PEPReasonViewFromOptions(PEP pep)
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)pep),((PetscObject)pep)->prefix,"-pep_converged_reason",&viewer,&format,&flg);
283: if (flg) {
284: PetscViewerPushFormat(viewer,format);
285: PEPReasonView(pep,viewer);
286: PetscViewerPopFormat(viewer);
287: PetscViewerDestroy(&viewer);
288: }
289: incall = PETSC_FALSE;
290: return(0);
291: }
295: static PetscErrorCode PEPErrorView_ASCII(PEP pep,PEPErrorType etype,PetscViewer viewer)
296: {
297: PetscBool errok;
298: PetscReal error,re,im;
299: PetscScalar kr,ki;
300: PetscInt i,j;
304: if (pep->nconv<pep->nev) {
305: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",pep->nev);
306: return(0);
307: }
308: errok = PETSC_TRUE;
309: for (i=0;i<pep->nev;i++) {
310: PEPComputeError(pep,i,etype,&error);
311: errok = (errok && error<5.0*pep->tol)? PETSC_TRUE: PETSC_FALSE;
312: }
313: if (!errok) {
314: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",pep->nev);
315: return(0);
316: }
317: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
318: for (i=0;i<=(pep->nev-1)/8;i++) {
319: PetscViewerASCIIPrintf(viewer,"\n ");
320: for (j=0;j<PetscMin(8,pep->nev-8*i);j++) {
321: PEPGetEigenpair(pep,8*i+j,&kr,&ki,NULL,NULL);
322: #if defined(PETSC_USE_COMPLEX)
323: re = PetscRealPart(kr);
324: im = PetscImaginaryPart(kr);
325: #else
326: re = kr;
327: im = ki;
328: #endif
329: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
330: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
331: if (im!=0.0) {
332: PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
333: } else {
334: PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
335: }
336: if (8*i+j+1<pep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
337: }
338: }
339: PetscViewerASCIIPrintf(viewer,"\n\n");
340: return(0);
341: }
345: static PetscErrorCode PEPErrorView_DETAIL(PEP pep,PEPErrorType etype,PetscViewer viewer)
346: {
348: PetscReal error,re,im;
349: PetscScalar kr,ki;
350: PetscInt i;
351: #define EXLEN 30
352: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
355: if (!pep->nconv) return(0);
356: switch (etype) {
357: case PEP_ERROR_ABSOLUTE:
358: PetscSNPrintf(ex,EXLEN," ||P(k)x||");
359: break;
360: case PEP_ERROR_RELATIVE:
361: PetscSNPrintf(ex,EXLEN,"||P(k)x||/||kx||");
362: break;
363: case PEP_ERROR_BACKWARD:
364: PetscSNPrintf(ex,EXLEN," eta(x,k)");
365: break;
366: }
367: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
368: for (i=0;i<pep->nconv;i++) {
369: PEPGetEigenpair(pep,i,&kr,&ki,NULL,NULL);
370: PEPComputeError(pep,i,etype,&error);
371: #if defined(PETSC_USE_COMPLEX)
372: re = PetscRealPart(kr);
373: im = PetscImaginaryPart(kr);
374: #else
375: re = kr;
376: im = ki;
377: #endif
378: if (im!=0.0) {
379: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
380: } else {
381: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
382: }
383: }
384: PetscViewerASCIIPrintf(viewer,sep);
385: return(0);
386: }
390: static PetscErrorCode PEPErrorView_MATLAB(PEP pep,PEPErrorType etype,PetscViewer viewer)
391: {
393: PetscReal error;
394: PetscInt i;
395: const char *name;
398: PetscObjectGetName((PetscObject)pep,&name);
399: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
400: for (i=0;i<pep->nconv;i++) {
401: PEPComputeError(pep,i,etype,&error);
402: PetscViewerASCIIPrintf(viewer,"%18.16e\n",error);
403: }
404: PetscViewerASCIIPrintf(viewer,"];\n");
405: return(0);
406: }
410: /*@C
411: PEPErrorView - Displays the errors associated with the computed solution
412: (as well as the eigenvalues).
414: Collective on PEP
416: Input Parameters:
417: + pep - the eigensolver context
418: . etype - error type
419: - viewer - optional visualization context
421: Options Database Key:
422: + -pep_error_absolute - print absolute errors of each eigenpair
423: . -pep_error_relative - print relative errors of each eigenpair
424: - -pep_error_backward - print backward errors of each eigenpair
426: Notes:
427: By default, this function checks the error of all eigenpairs and prints
428: the eigenvalues if all of them are below the requested tolerance.
429: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
430: eigenvalues and corresponding errors is printed.
432: Level: intermediate
434: .seealso: PEPSolve(), PEPValuesView(), PEPVectorsView()
435: @*/
436: PetscErrorCode PEPErrorView(PEP pep,PEPErrorType etype,PetscViewer viewer)
437: {
438: PetscBool isascii;
439: PetscViewerFormat format;
440: PetscErrorCode ierr;
444: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
447: PEPCheckSolved(pep,1);
448: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
449: if (!isascii) return(0);
451: PetscViewerGetFormat(viewer,&format);
452: switch (format) {
453: case PETSC_VIEWER_DEFAULT:
454: case PETSC_VIEWER_ASCII_INFO:
455: PEPErrorView_ASCII(pep,etype,viewer);
456: break;
457: case PETSC_VIEWER_ASCII_INFO_DETAIL:
458: PEPErrorView_DETAIL(pep,etype,viewer);
459: break;
460: case PETSC_VIEWER_ASCII_MATLAB:
461: PEPErrorView_MATLAB(pep,etype,viewer);
462: break;
463: default:
464: PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
465: }
466: return(0);
467: }
471: /*@
472: PEPErrorViewFromOptions - Processes command line options to determine if/how
473: the errors of the computed solution are to be viewed.
475: Collective on PEP
477: Input Parameters:
478: . pep - the eigensolver context
480: Level: developer
481: @*/
482: PetscErrorCode PEPErrorViewFromOptions(PEP pep)
483: {
484: PetscErrorCode ierr;
485: PetscViewer viewer;
486: PetscBool flg;
487: static PetscBool incall = PETSC_FALSE;
488: PetscViewerFormat format;
491: if (incall) return(0);
492: incall = PETSC_TRUE;
493: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_error_absolute",&viewer,&format,&flg);
494: if (flg) {
495: PetscViewerPushFormat(viewer,format);
496: PEPErrorView(pep,PEP_ERROR_ABSOLUTE,viewer);
497: PetscViewerPopFormat(viewer);
498: PetscViewerDestroy(&viewer);
499: }
500: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_error_relative",&viewer,&format,&flg);
501: if (flg) {
502: PetscViewerPushFormat(viewer,format);
503: PEPErrorView(pep,PEP_ERROR_RELATIVE,viewer);
504: PetscViewerPopFormat(viewer);
505: PetscViewerDestroy(&viewer);
506: }
507: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_error_backward",&viewer,&format,&flg);
508: if (flg) {
509: PetscViewerPushFormat(viewer,format);
510: PEPErrorView(pep,PEP_ERROR_BACKWARD,viewer);
511: PetscViewerPopFormat(viewer);
512: PetscViewerDestroy(&viewer);
513: }
514: incall = PETSC_FALSE;
515: return(0);
516: }
520: static PetscErrorCode PEPValuesView_DRAW(PEP pep,PetscViewer viewer)
521: {
523: PetscDraw draw;
524: PetscDrawSP drawsp;
525: PetscReal re,im;
526: PetscInt i,k;
529: if (!pep->nconv) return(0);
530: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Computed Eigenvalues",PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
531: PetscViewerDrawGetDraw(viewer,0,&draw);
532: PetscDrawSPCreate(draw,1,&drawsp);
533: for (i=0;i<pep->nconv;i++) {
534: k = pep->perm[i];
535: #if defined(PETSC_USE_COMPLEX)
536: re = PetscRealPart(pep->eigr[k]);
537: im = PetscImaginaryPart(pep->eigr[k]);
538: #else
539: re = pep->eigr[k];
540: im = pep->eigi[k];
541: #endif
542: PetscDrawSPAddPoint(drawsp,&re,&im);
543: }
544: PetscDrawSPDraw(drawsp,PETSC_TRUE);
545: PetscDrawSPDestroy(&drawsp);
546: PetscViewerDestroy(&viewer);
547: return(0);
548: }
552: static PetscErrorCode PEPValuesView_ASCII(PEP pep,PetscViewer viewer)
553: {
554: PetscReal re,im;
555: PetscInt i,k;
559: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
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 (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
570: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
571: if (im!=0.0) {
572: PetscViewerASCIIPrintf(viewer," %.5f%+.5fi\n",(double)re,(double)im);
573: } else {
574: PetscViewerASCIIPrintf(viewer," %.5f\n",(double)re);
575: }
576: }
577: PetscViewerASCIIPrintf(viewer,"\n");
578: return(0);
579: }
583: static PetscErrorCode PEPValuesView_MATLAB(PEP pep,PetscViewer viewer)
584: {
586: PetscInt i,k;
587: PetscReal re,im;
588: const char *name;
591: PetscObjectGetName((PetscObject)pep,&name);
592: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
593: for (i=0;i<pep->nconv;i++) {
594: k = pep->perm[i];
595: #if defined(PETSC_USE_COMPLEX)
596: re = PetscRealPart(pep->eigr[k]);
597: im = PetscImaginaryPart(pep->eigr[k]);
598: #else
599: re = pep->eigr[k];
600: im = pep->eigi[k];
601: #endif
602: if (im!=0.0) {
603: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
604: } else {
605: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
606: }
607: }
608: PetscViewerASCIIPrintf(viewer,"];\n");
609: return(0);
610: }
614: /*@C
615: PEPValuesView - Displays the computed eigenvalues in a viewer.
617: Collective on PEP
619: Input Parameters:
620: + pep - the eigensolver context
621: - viewer - the viewer
623: Options Database Key:
624: . -pep_view_values - print computed eigenvalues
626: Level: intermediate
628: .seealso: PEPSolve(), PEPVectorsView(), PEPErrorView()
629: @*/
630: PetscErrorCode PEPValuesView(PEP pep,PetscViewer viewer)
631: {
632: PetscBool isascii,isdraw;
633: PetscViewerFormat format;
634: PetscErrorCode ierr;
638: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
641: PEPCheckSolved(pep,1);
642: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
643: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
644: if (isdraw) {
645: PEPValuesView_DRAW(pep,viewer);
646: } else if (isascii) {
647: PetscViewerGetFormat(viewer,&format);
648: switch (format) {
649: case PETSC_VIEWER_DEFAULT:
650: case PETSC_VIEWER_ASCII_INFO:
651: case PETSC_VIEWER_ASCII_INFO_DETAIL:
652: PEPValuesView_ASCII(pep,viewer);
653: break;
654: case PETSC_VIEWER_ASCII_MATLAB:
655: PEPValuesView_MATLAB(pep,viewer);
656: break;
657: default:
658: PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
659: }
660: }
661: return(0);
662: }
666: /*@
667: PEPValuesViewFromOptions - Processes command line options to determine if/how
668: the computed eigenvalues are to be viewed.
670: Collective on PEP
672: Input Parameters:
673: . pep - the eigensolver context
675: Level: developer
676: @*/
677: PetscErrorCode PEPValuesViewFromOptions(PEP pep)
678: {
679: PetscErrorCode ierr;
680: PetscViewer viewer;
681: PetscBool flg;
682: static PetscBool incall = PETSC_FALSE;
683: PetscViewerFormat format;
686: if (incall) return(0);
687: incall = PETSC_TRUE;
688: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_view_values",&viewer,&format,&flg);
689: if (flg) {
690: PetscViewerPushFormat(viewer,format);
691: PEPValuesView(pep,viewer);
692: PetscViewerPopFormat(viewer);
693: PetscViewerDestroy(&viewer);
694: }
695: incall = PETSC_FALSE;
696: return(0);
697: }
701: /*@C
702: PEPVectorsView - Outputs computed eigenvectors to a viewer.
704: Collective on PEP
706: Parameter:
707: + pep - the eigensolver context
708: - viewer - the viewer
710: Options Database Keys:
711: . -pep_view_vectors - output eigenvectors.
713: Note:
714: If PETSc was configured with real scalars, complex conjugate eigenvectors
715: will be viewed as two separate real vectors, one containing the real part
716: and another one containing the imaginary part.
718: Level: intermediate
720: .seealso: PEPSolve(), PEPValuesView(), PEPErrorView()
721: @*/
722: PetscErrorCode PEPVectorsView(PEP pep,PetscViewer viewer)
723: {
725: PetscInt i,k;
726: Vec x;
727: #define NMLEN 30
728: char vname[NMLEN];
729: const char *ename;
733: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
736: PEPCheckSolved(pep,1);
737: if (pep->nconv) {
738: PetscObjectGetName((PetscObject)pep,&ename);
739: PEPComputeVectors(pep);
740: for (i=0;i<pep->nconv;i++) {
741: k = pep->perm[i];
742: PetscSNPrintf(vname,NMLEN,"V%d_%s",i,ename);
743: BVGetColumn(pep->V,k,&x);
744: PetscObjectSetName((PetscObject)x,vname);
745: VecView(x,viewer);
746: BVRestoreColumn(pep->V,k,&x);
747: }
748: }
749: return(0);
750: }
754: /*@
755: PEPVectorsViewFromOptions - Processes command line options to determine if/how
756: the computed eigenvectors are to be viewed.
758: Collective on PEP
760: Input Parameters:
761: . pep - the eigensolver context
763: Level: developer
764: @*/
765: PetscErrorCode PEPVectorsViewFromOptions(PEP pep)
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)pep),((PetscObject)pep)->prefix,"-pep_view_vectors",&viewer,&format,&flg);
777: if (flg) {
778: PetscViewerPushFormat(viewer,format);
779: PEPVectorsView(pep,viewer);
780: PetscViewerPopFormat(viewer);
781: PetscViewerDestroy(&viewer);
782: }
783: incall = PETSC_FALSE;
784: return(0);
785: }