Actual source code: pepview.c

slepc-3.8.0 2017-10-20
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2017, 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

 42: .seealso: PetscViewerASCIIOpen()
 43: @*/
 44: PetscErrorCode PEPView(PEP pep,PetscViewer viewer)
 45: {
 47:   const char     *type=NULL;
 48:   char           str[50];
 49:   PetscBool      isascii,islinear,istrivial;
 50:   PetscInt       i;
 51:   PetscViewer    sviewer;

 55:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));

 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_GYROSCOPIC: type = "gyroscopic polynomial eigenvalue problem"; break;
 77:       }
 78:     } else type = "not yet set";
 79:     PetscViewerASCIIPrintf(viewer,"  problem type: %s\n",type);
 80:     PetscViewerASCIIPrintf(viewer,"  polynomial represented in %s basis\n",PEPBasisTypes[pep->basis]);
 81:     switch (pep->scale) {
 82:       case PEP_SCALE_NONE:
 83:         break;
 84:       case PEP_SCALE_SCALAR:
 85:         PetscViewerASCIIPrintf(viewer,"  parameter scaling enabled, with scaling factor=%g\n",(double)pep->sfactor);
 86:         break;
 87:       case PEP_SCALE_DIAGONAL:
 88:         PetscViewerASCIIPrintf(viewer,"  diagonal balancing enabled, with its=%D and lambda=%g\n",pep->sits,(double)pep->slambda);
 89:         break;
 90:       case PEP_SCALE_BOTH:
 91:         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);
 92:         break;
 93:     }
 94:     PetscViewerASCIIPrintf(viewer,"  selected portion of the spectrum: ");
 95:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
 96:     SlepcSNPrintfScalar(str,50,pep->target,PETSC_FALSE);
 97:     if (!pep->which) {
 98:       PetscViewerASCIIPrintf(viewer,"not yet set\n");
 99:     } else switch (pep->which) {
100:       case PEP_WHICH_USER:
101:         PetscViewerASCIIPrintf(viewer,"user defined\n");
102:         break;
103:       case PEP_TARGET_MAGNITUDE:
104:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
105:         break;
106:       case PEP_TARGET_REAL:
107:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
108:         break;
109:       case PEP_TARGET_IMAGINARY:
110:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
111:         break;
112:       case PEP_LARGEST_MAGNITUDE:
113:         PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
114:         break;
115:       case PEP_SMALLEST_MAGNITUDE:
116:         PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
117:         break;
118:       case PEP_LARGEST_REAL:
119:         PetscViewerASCIIPrintf(viewer,"largest real parts\n");
120:         break;
121:       case PEP_SMALLEST_REAL:
122:         PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
123:         break;
124:       case PEP_LARGEST_IMAGINARY:
125:         PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
126:         break;
127:       case PEP_SMALLEST_IMAGINARY:
128:         PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
129:         break;
130:     }
131:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
132:     PetscViewerASCIIPrintf(viewer,"  number of eigenvalues (nev): %D\n",pep->nev);
133:     PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %D\n",pep->ncv);
134:     PetscViewerASCIIPrintf(viewer,"  maximum dimension of projected problem (mpd): %D\n",pep->mpd);
135:     PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %D\n",pep->max_it);
136:     PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)pep->tol);
137:     PetscViewerASCIIPrintf(viewer,"  convergence test: ");
138:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
139:     switch (pep->conv) {
140:     case PEP_CONV_ABS:
141:       PetscViewerASCIIPrintf(viewer,"absolute\n");break;
142:     case PEP_CONV_REL:
143:       PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
144:     case PEP_CONV_NORM:
145:       PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
146:       if (pep->nrma) {
147:         PetscViewerASCIIPrintf(viewer,"  computed matrix norms: %g",(double)pep->nrma[0]);
148:         for (i=1;i<pep->nmat;i++) {
149:           PetscViewerASCIIPrintf(viewer,", %g",(double)pep->nrma[i]);
150:         }
151:         PetscViewerASCIIPrintf(viewer,"\n");
152:       }
153:       break;
154:     case PEP_CONV_USER:
155:       PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
156:     }
157:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
158:     PetscViewerASCIIPrintf(viewer,"  extraction type: %s\n",PEPExtractTypes[pep->extract]);
159:     if (pep->refine) {
160:       PetscViewerASCIIPrintf(viewer,"  iterative refinement: %s, with %s scheme\n",PEPRefineTypes[pep->refine],PEPRefineSchemes[pep->scheme]);
161:       PetscViewerASCIIPrintf(viewer,"  refinement stopping criterion: tol=%g, its=%D\n",(double)pep->rtol,pep->rits);
162:       if (pep->npart>1) {
163:         PetscViewerASCIIPrintf(viewer,"  splitting communicator in %D partitions for refinement\n",pep->npart);
164:       }
165:     }
166:     if (pep->nini) {
167:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial space: %D\n",PetscAbs(pep->nini));
168:     }
169:   } else {
170:     if (pep->ops->view) {
171:       (*pep->ops->view)(pep,viewer);
172:     }
173:   }
174:   PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
175:   if (!pep->V) { PEPGetBV(pep,&pep->V); }
176:   BVView(pep->V,viewer);
177:   if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
178:   RGIsTrivial(pep->rg,&istrivial);
179:   if (!istrivial) { RGView(pep->rg,viewer); }
180:   PetscObjectTypeCompare((PetscObject)pep,PEPLINEAR,&islinear);
181:   if (!islinear) {
182:     if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
183:     DSView(pep->ds,viewer);
184:   }
185:   PetscViewerPopFormat(viewer);
186:   if (!pep->st) { PEPGetST(pep,&pep->st); }
187:   STView(pep->st,viewer);
188:   if (pep->refine!=PEP_REFINE_NONE) {
189:     if (pep->npart>1) {
190:       if (pep->refinesubc->color==0) {
191:         PetscViewerASCIIGetStdout(PetscSubcommChild(pep->refinesubc),&sviewer);
192:         KSPView(pep->refineksp,sviewer);
193:       }
194:     } else {
195:       KSPView(pep->refineksp,viewer);
196:     }
197:   }
198:   return(0);
199: }

201: /*@C
202:    PEPReasonView - Displays the reason a PEP solve converged or diverged.

204:    Collective on PEP

206:    Parameter:
207: +  pep - the eigensolver context
208: -  viewer - the viewer to display the reason

210:    Options Database Keys:
211: .  -pep_converged_reason - print reason for convergence, and number of iterations

213:    Level: intermediate

215: .seealso: PEPSetConvergenceTest(), PEPSetTolerances(), PEPGetIterationNumber()
216: @*/
217: PetscErrorCode PEPReasonView(PEP pep,PetscViewer viewer)
218: {
220:   PetscBool      isAscii;

223:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
224:   if (isAscii) {
225:     PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
226:     if (pep->reason > 0) {
227:       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);
228:     } else {
229:       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);
230:     }
231:     PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
232:   }
233:   return(0);
234: }

236: /*@
237:    PEPReasonViewFromOptions - Processes command line options to determine if/how
238:    the PEP converged reason is to be viewed.

240:    Collective on PEP

242:    Input Parameters:
243: .  pep - the eigensolver context

245:    Level: developer
246: @*/
247: PetscErrorCode PEPReasonViewFromOptions(PEP pep)
248: {
249:   PetscErrorCode    ierr;
250:   PetscViewer       viewer;
251:   PetscBool         flg;
252:   static PetscBool  incall = PETSC_FALSE;
253:   PetscViewerFormat format;

256:   if (incall) return(0);
257:   incall = PETSC_TRUE;
258:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_converged_reason",&viewer,&format,&flg);
259:   if (flg) {
260:     PetscViewerPushFormat(viewer,format);
261:     PEPReasonView(pep,viewer);
262:     PetscViewerPopFormat(viewer);
263:     PetscViewerDestroy(&viewer);
264:   }
265:   incall = PETSC_FALSE;
266:   return(0);
267: }

269: static PetscErrorCode PEPErrorView_ASCII(PEP pep,PEPErrorType etype,PetscViewer viewer)
270: {
271:   PetscReal      error;
272:   PetscInt       i,j,k;

276:   if (pep->nconv<pep->nev) {
277:     PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",pep->nev);
278:     return(0);
279:   }
280:   for (i=0;i<pep->nev;i++) {
281:     PEPComputeError(pep,i,etype,&error);
282:     if (error>=5.0*pep->tol) {
283:       PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",pep->nev);
284:       return(0);
285:     }
286:   }
287:   PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
288:   for (i=0;i<=(pep->nev-1)/8;i++) {
289:     PetscViewerASCIIPrintf(viewer,"\n     ");
290:     for (j=0;j<PetscMin(8,pep->nev-8*i);j++) {
291:       k = pep->perm[8*i+j];
292:       SlepcPrintEigenvalueASCII(pep->eigr[k],pep->eigi[k]);
293:       if (8*i+j+1<pep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
294:     }
295:   }
296:   PetscViewerASCIIPrintf(viewer,"\n\n");
297:   return(0);
298: }

300: static PetscErrorCode PEPErrorView_DETAIL(PEP pep,PEPErrorType etype,PetscViewer viewer)
301: {
303:   PetscReal      error,re,im;
304:   PetscScalar    kr,ki;
305:   PetscInt       i;
306: #define EXLEN 30
307:   char           ex[EXLEN],sep[]=" ---------------------- --------------------\n";

310:   if (!pep->nconv) return(0);
311:   switch (etype) {
312:     case PEP_ERROR_ABSOLUTE:
313:       PetscSNPrintf(ex,EXLEN,"   ||P(k)x||");
314:       break;
315:     case PEP_ERROR_RELATIVE:
316:       PetscSNPrintf(ex,EXLEN,"||P(k)x||/||kx||");
317:       break;
318:     case PEP_ERROR_BACKWARD:
319:       PetscSNPrintf(ex,EXLEN,"    eta(x,k)");
320:       break;
321:   }
322:   PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep);
323:   for (i=0;i<pep->nconv;i++) {
324:     PEPGetEigenpair(pep,i,&kr,&ki,NULL,NULL);
325:     PEPComputeError(pep,i,etype,&error);
326: #if defined(PETSC_USE_COMPLEX)
327:     re = PetscRealPart(kr);
328:     im = PetscImaginaryPart(kr);
329: #else
330:     re = kr;
331:     im = ki;
332: #endif
333:     if (im!=0.0) {
334:       PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
335:     } else {
336:       PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
337:     }
338:   }
339:   PetscViewerASCIIPrintf(viewer,sep);
340:   return(0);
341: }

343: static PetscErrorCode PEPErrorView_MATLAB(PEP pep,PEPErrorType etype,PetscViewer viewer)
344: {
346:   PetscReal      error;
347:   PetscInt       i;
348:   const char     *name;

351:   PetscObjectGetName((PetscObject)pep,&name);
352:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
353:   for (i=0;i<pep->nconv;i++) {
354:     PEPComputeError(pep,i,etype,&error);
355:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
356:   }
357:   PetscViewerASCIIPrintf(viewer,"];\n");
358:   return(0);
359: }

361: /*@C
362:    PEPErrorView - Displays the errors associated with the computed solution
363:    (as well as the eigenvalues).

365:    Collective on PEP

367:    Input Parameters:
368: +  pep    - the eigensolver context
369: .  etype  - error type
370: -  viewer - optional visualization context

372:    Options Database Key:
373: +  -pep_error_absolute - print absolute errors of each eigenpair
374: .  -pep_error_relative - print relative errors of each eigenpair
375: -  -pep_error_backward - print backward errors of each eigenpair

377:    Notes:
378:    By default, this function checks the error of all eigenpairs and prints
379:    the eigenvalues if all of them are below the requested tolerance.
380:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
381:    eigenvalues and corresponding errors is printed.

383:    Level: intermediate

385: .seealso: PEPSolve(), PEPValuesView(), PEPVectorsView()
386: @*/
387: PetscErrorCode PEPErrorView(PEP pep,PEPErrorType etype,PetscViewer viewer)
388: {
389:   PetscBool         isascii;
390:   PetscViewerFormat format;
391:   PetscErrorCode    ierr;

395:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
398:   PEPCheckSolved(pep,1);
399:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
400:   if (!isascii) return(0);

402:   PetscViewerGetFormat(viewer,&format);
403:   switch (format) {
404:     case PETSC_VIEWER_DEFAULT:
405:     case PETSC_VIEWER_ASCII_INFO:
406:       PEPErrorView_ASCII(pep,etype,viewer);
407:       break;
408:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
409:       PEPErrorView_DETAIL(pep,etype,viewer);
410:       break;
411:     case PETSC_VIEWER_ASCII_MATLAB:
412:       PEPErrorView_MATLAB(pep,etype,viewer);
413:       break;
414:     default:
415:       PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
416:   }
417:   return(0);
418: }

420: /*@
421:    PEPErrorViewFromOptions - Processes command line options to determine if/how
422:    the errors of the computed solution are to be viewed.

424:    Collective on PEP

426:    Input Parameters:
427: .  pep - the eigensolver context

429:    Level: developer
430: @*/
431: PetscErrorCode PEPErrorViewFromOptions(PEP pep)
432: {
433:   PetscErrorCode    ierr;
434:   PetscViewer       viewer;
435:   PetscBool         flg;
436:   static PetscBool  incall = PETSC_FALSE;
437:   PetscViewerFormat format;

440:   if (incall) return(0);
441:   incall = PETSC_TRUE;
442:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_error_absolute",&viewer,&format,&flg);
443:   if (flg) {
444:     PetscViewerPushFormat(viewer,format);
445:     PEPErrorView(pep,PEP_ERROR_ABSOLUTE,viewer);
446:     PetscViewerPopFormat(viewer);
447:     PetscViewerDestroy(&viewer);
448:   }
449:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_error_relative",&viewer,&format,&flg);
450:   if (flg) {
451:     PetscViewerPushFormat(viewer,format);
452:     PEPErrorView(pep,PEP_ERROR_RELATIVE,viewer);
453:     PetscViewerPopFormat(viewer);
454:     PetscViewerDestroy(&viewer);
455:   }
456:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_error_backward",&viewer,&format,&flg);
457:   if (flg) {
458:     PetscViewerPushFormat(viewer,format);
459:     PEPErrorView(pep,PEP_ERROR_BACKWARD,viewer);
460:     PetscViewerPopFormat(viewer);
461:     PetscViewerDestroy(&viewer);
462:   }
463:   incall = PETSC_FALSE;
464:   return(0);
465: }

467: static PetscErrorCode PEPValuesView_DRAW(PEP pep,PetscViewer viewer)
468: {
470:   PetscDraw      draw;
471:   PetscDrawSP    drawsp;
472:   PetscReal      re,im;
473:   PetscInt       i,k;

476:   if (!pep->nconv) return(0);
477:   PetscViewerDrawGetDraw(viewer,0,&draw);
478:   PetscDrawSetTitle(draw,"Computed Eigenvalues");
479:   PetscDrawSPCreate(draw,1,&drawsp);
480:   for (i=0;i<pep->nconv;i++) {
481:     k = pep->perm[i];
482: #if defined(PETSC_USE_COMPLEX)
483:     re = PetscRealPart(pep->eigr[k]);
484:     im = PetscImaginaryPart(pep->eigr[k]);
485: #else
486:     re = pep->eigr[k];
487:     im = pep->eigi[k];
488: #endif
489:     PetscDrawSPAddPoint(drawsp,&re,&im);
490:   }
491:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
492:   PetscDrawSPSave(drawsp);
493:   PetscDrawSPDestroy(&drawsp);
494:   return(0);
495: }

497: static PetscErrorCode PEPValuesView_ASCII(PEP pep,PetscViewer viewer)
498: {
499:   PetscInt       i,k;

503:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
504:   for (i=0;i<pep->nconv;i++) {
505:     k = pep->perm[i];
506:     PetscViewerASCIIPrintf(viewer,"   ");
507:     SlepcPrintEigenvalueASCII(pep->eigr[k],pep->eigi[k]);
508:     PetscViewerASCIIPrintf(viewer,"\n");
509:   }
510:   PetscViewerASCIIPrintf(viewer,"\n");
511:   return(0);
512: }

514: static PetscErrorCode PEPValuesView_MATLAB(PEP pep,PetscViewer viewer)
515: {
517:   PetscInt       i,k;
518:   PetscReal      re,im;
519:   const char     *name;

522:   PetscObjectGetName((PetscObject)pep,&name);
523:   PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
524:   for (i=0;i<pep->nconv;i++) {
525:     k = pep->perm[i];
526: #if defined(PETSC_USE_COMPLEX)
527:     re = PetscRealPart(pep->eigr[k]);
528:     im = PetscImaginaryPart(pep->eigr[k]);
529: #else
530:     re = pep->eigr[k];
531:     im = pep->eigi[k];
532: #endif
533:     if (im!=0.0) {
534:       PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
535:     } else {
536:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
537:     }
538:   }
539:   PetscViewerASCIIPrintf(viewer,"];\n");
540:   return(0);
541: }

543: /*@C
544:    PEPValuesView - Displays the computed eigenvalues in a viewer.

546:    Collective on PEP

548:    Input Parameters:
549: +  pep    - the eigensolver context
550: -  viewer - the viewer

552:    Options Database Key:
553: .  -pep_view_values - print computed eigenvalues

555:    Level: intermediate

557: .seealso: PEPSolve(), PEPVectorsView(), PEPErrorView()
558: @*/
559: PetscErrorCode PEPValuesView(PEP pep,PetscViewer viewer)
560: {
561:   PetscBool         isascii,isdraw;
562:   PetscViewerFormat format;
563:   PetscErrorCode    ierr;

567:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
570:   PEPCheckSolved(pep,1);
571:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
572:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
573:   if (isdraw) {
574:     PEPValuesView_DRAW(pep,viewer);
575:   } else if (isascii) {
576:     PetscViewerGetFormat(viewer,&format);
577:     switch (format) {
578:       case PETSC_VIEWER_DEFAULT:
579:       case PETSC_VIEWER_ASCII_INFO:
580:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
581:         PEPValuesView_ASCII(pep,viewer);
582:         break;
583:       case PETSC_VIEWER_ASCII_MATLAB:
584:         PEPValuesView_MATLAB(pep,viewer);
585:         break;
586:       default:
587:         PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
588:     }
589:   }
590:   return(0);
591: }

593: /*@
594:    PEPValuesViewFromOptions - Processes command line options to determine if/how
595:    the computed eigenvalues are to be viewed.

597:    Collective on PEP

599:    Input Parameters:
600: .  pep - the eigensolver context

602:    Level: developer
603: @*/
604: PetscErrorCode PEPValuesViewFromOptions(PEP pep)
605: {
606:   PetscErrorCode    ierr;
607:   PetscViewer       viewer;
608:   PetscBool         flg;
609:   static PetscBool  incall = PETSC_FALSE;
610:   PetscViewerFormat format;

613:   if (incall) return(0);
614:   incall = PETSC_TRUE;
615:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_view_values",&viewer,&format,&flg);
616:   if (flg) {
617:     PetscViewerPushFormat(viewer,format);
618:     PEPValuesView(pep,viewer);
619:     PetscViewerPopFormat(viewer);
620:     PetscViewerDestroy(&viewer);
621:   }
622:   incall = PETSC_FALSE;
623:   return(0);
624: }

626: /*@C
627:    PEPVectorsView - Outputs computed eigenvectors to a viewer.

629:    Collective on PEP

631:    Parameter:
632: +  pep    - the eigensolver context
633: -  viewer - the viewer

635:    Options Database Keys:
636: .  -pep_view_vectors - output eigenvectors.

638:    Note:
639:    If PETSc was configured with real scalars, complex conjugate eigenvectors
640:    will be viewed as two separate real vectors, one containing the real part
641:    and another one containing the imaginary part.

643:    Level: intermediate

645: .seealso: PEPSolve(), PEPValuesView(), PEPErrorView()
646: @*/
647: PetscErrorCode PEPVectorsView(PEP pep,PetscViewer viewer)
648: {
650:   PetscInt       i,k;
651:   Vec            x;
652: #define NMLEN 30
653:   char           vname[NMLEN];
654:   const char     *ename;

658:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
661:   PEPCheckSolved(pep,1);
662:   if (pep->nconv) {
663:     PetscObjectGetName((PetscObject)pep,&ename);
664:     PEPComputeVectors(pep);
665:     for (i=0;i<pep->nconv;i++) {
666:       k = pep->perm[i];
667:       PetscSNPrintf(vname,NMLEN,"V%d_%s",(int)i,ename);
668:       BVGetColumn(pep->V,k,&x);
669:       PetscObjectSetName((PetscObject)x,vname);
670:       VecView(x,viewer);
671:       BVRestoreColumn(pep->V,k,&x);
672:     }
673:   }
674:   return(0);
675: }

677: /*@
678:    PEPVectorsViewFromOptions - Processes command line options to determine if/how
679:    the computed eigenvectors are to be viewed.

681:    Collective on PEP

683:    Input Parameters:
684: .  pep - the eigensolver context

686:    Level: developer
687: @*/
688: PetscErrorCode PEPVectorsViewFromOptions(PEP pep)
689: {
690:   PetscErrorCode    ierr;
691:   PetscViewer       viewer;
692:   PetscBool         flg = PETSC_FALSE;
693:   static PetscBool  incall = PETSC_FALSE;
694:   PetscViewerFormat format;

697:   if (incall) return(0);
698:   incall = PETSC_TRUE;
699:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_view_vectors",&viewer,&format,&flg);
700:   if (flg) {
701:     PetscViewerPushFormat(viewer,format);
702:     PEPVectorsView(pep,viewer);
703:     PetscViewerPopFormat(viewer);
704:     PetscViewerDestroy(&viewer);
705:   }
706:   incall = PETSC_FALSE;
707:   return(0);
708: }