Actual source code: pepview.c

slepc-3.11.2 2019-07-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2019, 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) {
 56:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
 57:   }

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

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

216:    Collective on PEP

218:    Parameter:
219: +  pep - the eigensolver context
220: -  viewer - the viewer to display the reason

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

225:    Level: intermediate

227: .seealso: PEPSetConvergenceTest(), PEPSetTolerances(), PEPGetIterationNumber()
228: @*/
229: PetscErrorCode PEPReasonView(PEP pep,PetscViewer viewer)
230: {
232:   PetscBool      isAscii;

235:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
236:   if (isAscii) {
237:     PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
238:     if (pep->reason > 0) {
239:       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);
240:     } else {
241:       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);
242:     }
243:     PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
244:   }
245:   return(0);
246: }

248: /*@
249:    PEPReasonViewFromOptions - Processes command line options to determine if/how
250:    the PEP converged reason is to be viewed.

252:    Collective on PEP

254:    Input Parameters:
255: .  pep - the eigensolver context

257:    Level: developer
258: @*/
259: PetscErrorCode PEPReasonViewFromOptions(PEP pep)
260: {
261:   PetscErrorCode    ierr;
262:   PetscViewer       viewer;
263:   PetscBool         flg;
264:   static PetscBool  incall = PETSC_FALSE;
265:   PetscViewerFormat format;

268:   if (incall) return(0);
269:   incall = PETSC_TRUE;
270:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_converged_reason",&viewer,&format,&flg);
271:   if (flg) {
272:     PetscViewerPushFormat(viewer,format);
273:     PEPReasonView(pep,viewer);
274:     PetscViewerPopFormat(viewer);
275:     PetscViewerDestroy(&viewer);
276:   }
277:   incall = PETSC_FALSE;
278:   return(0);
279: }

281: static PetscErrorCode PEPErrorView_ASCII(PEP pep,PEPErrorType etype,PetscViewer viewer)
282: {
283:   PetscReal      error;
284:   PetscInt       i,j,k;

288:   if (pep->nconv<pep->nev) {
289:     PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",pep->nev);
290:     return(0);
291:   }
292:   for (i=0;i<pep->nev;i++) {
293:     PEPComputeError(pep,i,etype,&error);
294:     if (error>=5.0*pep->tol) {
295:       PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",pep->nev);
296:       return(0);
297:     }
298:   }
299:   PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
300:   for (i=0;i<=(pep->nev-1)/8;i++) {
301:     PetscViewerASCIIPrintf(viewer,"\n     ");
302:     for (j=0;j<PetscMin(8,pep->nev-8*i);j++) {
303:       k = pep->perm[8*i+j];
304:       SlepcPrintEigenvalueASCII(pep->eigr[k],pep->eigi[k]);
305:       if (8*i+j+1<pep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
306:     }
307:   }
308:   PetscViewerASCIIPrintf(viewer,"\n\n");
309:   return(0);
310: }

312: static PetscErrorCode PEPErrorView_DETAIL(PEP pep,PEPErrorType etype,PetscViewer viewer)
313: {
315:   PetscReal      error,re,im;
316:   PetscScalar    kr,ki;
317:   PetscInt       i;
318: #define EXLEN 30
319:   char           ex[EXLEN],sep[]=" ---------------------- --------------------\n";

322:   if (!pep->nconv) return(0);
323:   switch (etype) {
324:     case PEP_ERROR_ABSOLUTE:
325:       PetscSNPrintf(ex,EXLEN,"   ||P(k)x||");
326:       break;
327:     case PEP_ERROR_RELATIVE:
328:       PetscSNPrintf(ex,EXLEN,"||P(k)x||/||kx||");
329:       break;
330:     case PEP_ERROR_BACKWARD:
331:       PetscSNPrintf(ex,EXLEN,"    eta(x,k)");
332:       break;
333:   }
334:   PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep);
335:   for (i=0;i<pep->nconv;i++) {
336:     PEPGetEigenpair(pep,i,&kr,&ki,NULL,NULL);
337:     PEPComputeError(pep,i,etype,&error);
338: #if defined(PETSC_USE_COMPLEX)
339:     re = PetscRealPart(kr);
340:     im = PetscImaginaryPart(kr);
341: #else
342:     re = kr;
343:     im = ki;
344: #endif
345:     if (im!=0.0) {
346:       PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
347:     } else {
348:       PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
349:     }
350:   }
351:   PetscViewerASCIIPrintf(viewer,sep);
352:   return(0);
353: }

355: static PetscErrorCode PEPErrorView_MATLAB(PEP pep,PEPErrorType etype,PetscViewer viewer)
356: {
358:   PetscReal      error;
359:   PetscInt       i;
360:   const char     *name;

363:   PetscObjectGetName((PetscObject)pep,&name);
364:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
365:   for (i=0;i<pep->nconv;i++) {
366:     PEPComputeError(pep,i,etype,&error);
367:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
368:   }
369:   PetscViewerASCIIPrintf(viewer,"];\n");
370:   return(0);
371: }

373: /*@C
374:    PEPErrorView - Displays the errors associated with the computed solution
375:    (as well as the eigenvalues).

377:    Collective on PEP

379:    Input Parameters:
380: +  pep    - the eigensolver context
381: .  etype  - error type
382: -  viewer - optional visualization context

384:    Options Database Key:
385: +  -pep_error_absolute - print absolute errors of each eigenpair
386: .  -pep_error_relative - print relative errors of each eigenpair
387: -  -pep_error_backward - print backward errors of each eigenpair

389:    Notes:
390:    By default, this function checks the error of all eigenpairs and prints
391:    the eigenvalues if all of them are below the requested tolerance.
392:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
393:    eigenvalues and corresponding errors is printed.

395:    Level: intermediate

397: .seealso: PEPSolve(), PEPValuesView(), PEPVectorsView()
398: @*/
399: PetscErrorCode PEPErrorView(PEP pep,PEPErrorType etype,PetscViewer viewer)
400: {
401:   PetscBool         isascii;
402:   PetscViewerFormat format;
403:   PetscErrorCode    ierr;

407:   if (!viewer) {
408:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
409:   }
412:   PEPCheckSolved(pep,1);
413:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
414:   if (!isascii) return(0);

416:   PetscViewerGetFormat(viewer,&format);
417:   switch (format) {
418:     case PETSC_VIEWER_DEFAULT:
419:     case PETSC_VIEWER_ASCII_INFO:
420:       PEPErrorView_ASCII(pep,etype,viewer);
421:       break;
422:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
423:       PEPErrorView_DETAIL(pep,etype,viewer);
424:       break;
425:     case PETSC_VIEWER_ASCII_MATLAB:
426:       PEPErrorView_MATLAB(pep,etype,viewer);
427:       break;
428:     default:
429:       PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
430:   }
431:   return(0);
432: }

434: /*@
435:    PEPErrorViewFromOptions - Processes command line options to determine if/how
436:    the errors of the computed solution are to be viewed.

438:    Collective on PEP

440:    Input Parameters:
441: .  pep - the eigensolver context

443:    Level: developer
444: @*/
445: PetscErrorCode PEPErrorViewFromOptions(PEP pep)
446: {
447:   PetscErrorCode    ierr;
448:   PetscViewer       viewer;
449:   PetscBool         flg;
450:   static PetscBool  incall = PETSC_FALSE;
451:   PetscViewerFormat format;

454:   if (incall) return(0);
455:   incall = PETSC_TRUE;
456:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_error_absolute",&viewer,&format,&flg);
457:   if (flg) {
458:     PetscViewerPushFormat(viewer,format);
459:     PEPErrorView(pep,PEP_ERROR_ABSOLUTE,viewer);
460:     PetscViewerPopFormat(viewer);
461:     PetscViewerDestroy(&viewer);
462:   }
463:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_error_relative",&viewer,&format,&flg);
464:   if (flg) {
465:     PetscViewerPushFormat(viewer,format);
466:     PEPErrorView(pep,PEP_ERROR_RELATIVE,viewer);
467:     PetscViewerPopFormat(viewer);
468:     PetscViewerDestroy(&viewer);
469:   }
470:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_error_backward",&viewer,&format,&flg);
471:   if (flg) {
472:     PetscViewerPushFormat(viewer,format);
473:     PEPErrorView(pep,PEP_ERROR_BACKWARD,viewer);
474:     PetscViewerPopFormat(viewer);
475:     PetscViewerDestroy(&viewer);
476:   }
477:   incall = PETSC_FALSE;
478:   return(0);
479: }

481: static PetscErrorCode PEPValuesView_DRAW(PEP pep,PetscViewer viewer)
482: {
484:   PetscDraw      draw;
485:   PetscDrawSP    drawsp;
486:   PetscReal      re,im;
487:   PetscInt       i,k;

490:   if (!pep->nconv) return(0);
491:   PetscViewerDrawGetDraw(viewer,0,&draw);
492:   PetscDrawSetTitle(draw,"Computed Eigenvalues");
493:   PetscDrawSPCreate(draw,1,&drawsp);
494:   for (i=0;i<pep->nconv;i++) {
495:     k = pep->perm[i];
496: #if defined(PETSC_USE_COMPLEX)
497:     re = PetscRealPart(pep->eigr[k]);
498:     im = PetscImaginaryPart(pep->eigr[k]);
499: #else
500:     re = pep->eigr[k];
501:     im = pep->eigi[k];
502: #endif
503:     PetscDrawSPAddPoint(drawsp,&re,&im);
504:   }
505:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
506:   PetscDrawSPSave(drawsp);
507:   PetscDrawSPDestroy(&drawsp);
508:   return(0);
509: }

511: static PetscErrorCode PEPValuesView_ASCII(PEP pep,PetscViewer viewer)
512: {
513:   PetscInt       i,k;

517:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
518:   for (i=0;i<pep->nconv;i++) {
519:     k = pep->perm[i];
520:     PetscViewerASCIIPrintf(viewer,"   ");
521:     SlepcPrintEigenvalueASCII(pep->eigr[k],pep->eigi[k]);
522:     PetscViewerASCIIPrintf(viewer,"\n");
523:   }
524:   PetscViewerASCIIPrintf(viewer,"\n");
525:   return(0);
526: }

528: static PetscErrorCode PEPValuesView_MATLAB(PEP pep,PetscViewer viewer)
529: {
531:   PetscInt       i,k;
532:   PetscReal      re,im;
533:   const char     *name;

536:   PetscObjectGetName((PetscObject)pep,&name);
537:   PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
538:   for (i=0;i<pep->nconv;i++) {
539:     k = pep->perm[i];
540: #if defined(PETSC_USE_COMPLEX)
541:     re = PetscRealPart(pep->eigr[k]);
542:     im = PetscImaginaryPart(pep->eigr[k]);
543: #else
544:     re = pep->eigr[k];
545:     im = pep->eigi[k];
546: #endif
547:     if (im!=0.0) {
548:       PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
549:     } else {
550:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
551:     }
552:   }
553:   PetscViewerASCIIPrintf(viewer,"];\n");
554:   return(0);
555: }

557: /*@C
558:    PEPValuesView - Displays the computed eigenvalues in a viewer.

560:    Collective on PEP

562:    Input Parameters:
563: +  pep    - the eigensolver context
564: -  viewer - the viewer

566:    Options Database Key:
567: .  -pep_view_values - print computed eigenvalues

569:    Level: intermediate

571: .seealso: PEPSolve(), PEPVectorsView(), PEPErrorView()
572: @*/
573: PetscErrorCode PEPValuesView(PEP pep,PetscViewer viewer)
574: {
575:   PetscBool         isascii,isdraw;
576:   PetscViewerFormat format;
577:   PetscErrorCode    ierr;

581:   if (!viewer) {
582:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
583:   }
586:   PEPCheckSolved(pep,1);
587:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
588:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
589:   if (isdraw) {
590:     PEPValuesView_DRAW(pep,viewer);
591:   } else if (isascii) {
592:     PetscViewerGetFormat(viewer,&format);
593:     switch (format) {
594:       case PETSC_VIEWER_DEFAULT:
595:       case PETSC_VIEWER_ASCII_INFO:
596:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
597:         PEPValuesView_ASCII(pep,viewer);
598:         break;
599:       case PETSC_VIEWER_ASCII_MATLAB:
600:         PEPValuesView_MATLAB(pep,viewer);
601:         break;
602:       default:
603:         PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
604:     }
605:   }
606:   return(0);
607: }

609: /*@
610:    PEPValuesViewFromOptions - Processes command line options to determine if/how
611:    the computed eigenvalues are to be viewed.

613:    Collective on PEP

615:    Input Parameters:
616: .  pep - the eigensolver context

618:    Level: developer
619: @*/
620: PetscErrorCode PEPValuesViewFromOptions(PEP pep)
621: {
622:   PetscErrorCode    ierr;
623:   PetscViewer       viewer;
624:   PetscBool         flg;
625:   static PetscBool  incall = PETSC_FALSE;
626:   PetscViewerFormat format;

629:   if (incall) return(0);
630:   incall = PETSC_TRUE;
631:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_view_values",&viewer,&format,&flg);
632:   if (flg) {
633:     PetscViewerPushFormat(viewer,format);
634:     PEPValuesView(pep,viewer);
635:     PetscViewerPopFormat(viewer);
636:     PetscViewerDestroy(&viewer);
637:   }
638:   incall = PETSC_FALSE;
639:   return(0);
640: }

642: /*@C
643:    PEPVectorsView - Outputs computed eigenvectors to a viewer.

645:    Collective on PEP

647:    Parameter:
648: +  pep    - the eigensolver context
649: -  viewer - the viewer

651:    Options Database Keys:
652: .  -pep_view_vectors - output eigenvectors.

654:    Note:
655:    If PETSc was configured with real scalars, complex conjugate eigenvectors
656:    will be viewed as two separate real vectors, one containing the real part
657:    and another one containing the imaginary part.

659:    Level: intermediate

661: .seealso: PEPSolve(), PEPValuesView(), PEPErrorView()
662: @*/
663: PetscErrorCode PEPVectorsView(PEP pep,PetscViewer viewer)
664: {
666:   PetscInt       i,k;
667:   Vec            x;
668: #define NMLEN 30
669:   char           vname[NMLEN];
670:   const char     *ename;

674:   if (!viewer) {
675:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
676:   }
679:   PEPCheckSolved(pep,1);
680:   if (pep->nconv) {
681:     PetscObjectGetName((PetscObject)pep,&ename);
682:     PEPComputeVectors(pep);
683:     for (i=0;i<pep->nconv;i++) {
684:       k = pep->perm[i];
685:       PetscSNPrintf(vname,NMLEN,"V%d_%s",(int)i,ename);
686:       BVGetColumn(pep->V,k,&x);
687:       PetscObjectSetName((PetscObject)x,vname);
688:       VecView(x,viewer);
689:       BVRestoreColumn(pep->V,k,&x);
690:     }
691:   }
692:   return(0);
693: }

695: /*@
696:    PEPVectorsViewFromOptions - Processes command line options to determine if/how
697:    the computed eigenvectors are to be viewed.

699:    Collective on PEP

701:    Input Parameters:
702: .  pep - the eigensolver context

704:    Level: developer
705: @*/
706: PetscErrorCode PEPVectorsViewFromOptions(PEP pep)
707: {
708:   PetscErrorCode    ierr;
709:   PetscViewer       viewer;
710:   PetscBool         flg = PETSC_FALSE;
711:   static PetscBool  incall = PETSC_FALSE;
712:   PetscViewerFormat format;

715:   if (incall) return(0);
716:   incall = PETSC_TRUE;
717:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_view_vectors",&viewer,&format,&flg);
718:   if (flg) {
719:     PetscViewerPushFormat(viewer,format);
720:     PEPVectorsView(pep,viewer);
721:     PetscViewerPopFormat(viewer);
722:     PetscViewerDestroy(&viewer);
723:   }
724:   incall = PETSC_FALSE;
725:   return(0);
726: }