Actual source code: pepview.c

slepc-3.7.1 2016-05-27
Report Typos and Errors
  1: /*
  2:    The PEP routines related to various viewers.

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2016, Universitat Politecnica de Valencia, Spain

  8:    This file is part of SLEPc.

 10:    SLEPc is free software: you can redistribute it and/or modify it under  the
 11:    terms of version 3 of the GNU Lesser General Public License as published by
 12:    the Free Software Foundation.

 14:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY
 15:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS
 16:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for
 17:    more details.

 19:    You  should have received a copy of the GNU Lesser General  Public  License
 20:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 21:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 22: */

 24: #include <slepc/private/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,"  parameter scaling 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,"  parameter scaling & 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,"  selected portion of the spectrum: ");
108:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
109:     SlepcSNPrintfScalar(str,50,pep->target,PETSC_FALSE);
110:     if (!pep->which) {
111:       PetscViewerASCIIPrintf(viewer,"not yet set\n");
112:     } else switch (pep->which) {
113:       case PEP_WHICH_USER:
114:         PetscViewerASCIIPrintf(viewer,"user defined\n");
115:         break;
116:       case PEP_TARGET_MAGNITUDE:
117:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
118:         break;
119:       case PEP_TARGET_REAL:
120:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
121:         break;
122:       case PEP_TARGET_IMAGINARY:
123:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
124:         break;
125:       case PEP_LARGEST_MAGNITUDE:
126:         PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
127:         break;
128:       case PEP_SMALLEST_MAGNITUDE:
129:         PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
130:         break;
131:       case PEP_LARGEST_REAL:
132:         PetscViewerASCIIPrintf(viewer,"largest real parts\n");
133:         break;
134:       case PEP_SMALLEST_REAL:
135:         PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
136:         break;
137:       case PEP_LARGEST_IMAGINARY:
138:         PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
139:         break;
140:       case PEP_SMALLEST_IMAGINARY:
141:         PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
142:         break;
143:       default: SETERRQ(PetscObjectComm((PetscObject)pep),1,"Wrong value of pep->which");
144:     }
145:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
146:     PetscViewerASCIIPrintf(viewer,"  number of eigenvalues (nev): %D\n",pep->nev);
147:     PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %D\n",pep->ncv);
148:     PetscViewerASCIIPrintf(viewer,"  maximum dimension of projected problem (mpd): %D\n",pep->mpd);
149:     PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %D\n",pep->max_it);
150:     PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)pep->tol);
151:     PetscViewerASCIIPrintf(viewer,"  convergence test: ");
152:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
153:     switch (pep->conv) {
154:     case PEP_CONV_ABS:
155:       PetscViewerASCIIPrintf(viewer,"absolute\n");break;
156:     case PEP_CONV_REL:
157:       PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
158:     case PEP_CONV_NORM:
159:       PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
160:       if (pep->nrma) {
161:         PetscViewerASCIIPrintf(viewer,"  computed matrix norms: %g",(double)pep->nrma[0]);
162:         for (i=1;i<pep->nmat;i++) {
163:           PetscViewerASCIIPrintf(viewer,", %g",(double)pep->nrma[i]);
164:         }
165:         PetscViewerASCIIPrintf(viewer,"\n");
166:       }
167:       break;
168:     case PEP_CONV_USER:
169:       PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
170:     }
171:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
172:     PetscViewerASCIIPrintf(viewer,"  extraction type: %s\n",PEPExtractTypes[pep->extract]);
173:     if (pep->refine) {
174:       PetscViewerASCIIPrintf(viewer,"  iterative refinement: %s, with %s scheme\n",PEPRefineTypes[pep->refine],PEPRefineSchemes[pep->scheme]);
175:       PetscViewerASCIIPrintf(viewer,"  refinement stopping criterion: tol=%g, its=%D\n",(double)pep->rtol,pep->rits);
176:       if (pep->npart>1) {
177:         PetscViewerASCIIPrintf(viewer,"  splitting communicator in %D partitions for refinement\n",pep->npart);
178:       }
179:     }
180:     if (pep->nini) {
181:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial space: %D\n",PetscAbs(pep->nini));
182:     }
183:   } else {
184:     if (pep->ops->view) {
185:       (*pep->ops->view)(pep,viewer);
186:     }
187:   }
188:   PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
189:   if (!pep->V) { PEPGetBV(pep,&pep->V); }
190:   BVView(pep->V,viewer);
191:   if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
192:   RGIsTrivial(pep->rg,&istrivial);
193:   if (!istrivial) { RGView(pep->rg,viewer); }
194:   PetscObjectTypeCompare((PetscObject)pep,PEPLINEAR,&islinear);
195:   if (!islinear) {
196:     if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
197:     DSView(pep->ds,viewer);
198:   }
199:   PetscViewerPopFormat(viewer);
200:   if (!pep->st) { PEPGetST(pep,&pep->st); }
201:   STView(pep->st,viewer);
202:   if (pep->refine!=PEP_REFINE_NONE) {
203:     if (pep->npart>1) {
204:       if (pep->refinesubc->color==0) {
205:         PetscViewerASCIIGetStdout(PetscSubcommChild(pep->refinesubc),&sviewer);
206:         KSPView(pep->refineksp,sviewer);
207:       }
208:     } else {
209:       KSPView(pep->refineksp,viewer);
210:     }
211:   } 
212:   return(0);
213: }

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

220:    Collective on PEP

222:    Parameter:
223: +  pep - the eigensolver context
224: -  viewer - the viewer to display the reason

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

229:    Level: intermediate

231: .seealso: PEPSetConvergenceTest(), PEPSetTolerances(), PEPGetIterationNumber()
232: @*/
233: PetscErrorCode PEPReasonView(PEP pep,PetscViewer viewer)
234: {
236:   PetscBool      isAscii;

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

254: /*@
255:    PEPReasonViewFromOptions - Processes command line options to determine if/how
256:    the PEP converged reason is to be viewed. 

258:    Collective on PEP

260:    Input Parameters:
261: .  pep - the eigensolver context

263:    Level: developer
264: @*/
265: PetscErrorCode PEPReasonViewFromOptions(PEP pep)
266: {
267:   PetscErrorCode    ierr;
268:   PetscViewer       viewer;
269:   PetscBool         flg;
270:   static PetscBool  incall = PETSC_FALSE;
271:   PetscViewerFormat format;

274:   if (incall) return(0);
275:   incall = PETSC_TRUE;
276:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_converged_reason",&viewer,&format,&flg);
277:   if (flg) {
278:     PetscViewerPushFormat(viewer,format);
279:     PEPReasonView(pep,viewer);
280:     PetscViewerPopFormat(viewer);
281:     PetscViewerDestroy(&viewer);
282:   }
283:   incall = PETSC_FALSE;
284:   return(0);
285: }

289: static PetscErrorCode PEPErrorView_ASCII(PEP pep,PEPErrorType etype,PetscViewer viewer)
290: {
291:   PetscBool      errok;
292:   PetscReal      error,re,im;
293:   PetscScalar    kr,ki;
294:   PetscInt       i,j;

298:   if (pep->nconv<pep->nev) {
299:     PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",pep->nev);
300:     return(0);
301:   }
302:   errok = PETSC_TRUE;
303:   for (i=0;i<pep->nev;i++) {
304:     PEPComputeError(pep,i,etype,&error);
305:     errok = (errok && error<5.0*pep->tol)? PETSC_TRUE: PETSC_FALSE;
306:   }
307:   if (!errok) {
308:     PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",pep->nev);
309:     return(0);
310:   }
311:   PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
312:   for (i=0;i<=(pep->nev-1)/8;i++) {
313:     PetscViewerASCIIPrintf(viewer,"\n     ");
314:     for (j=0;j<PetscMin(8,pep->nev-8*i);j++) {
315:       PEPGetEigenpair(pep,8*i+j,&kr,&ki,NULL,NULL);
316: #if defined(PETSC_USE_COMPLEX)
317:       re = PetscRealPart(kr);
318:       im = PetscImaginaryPart(kr);
319: #else
320:       re = kr;
321:       im = ki;
322: #endif
323:       if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
324:       if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
325:       if (im!=0.0) {
326:         PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
327:       } else {
328:         PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
329:       }
330:       if (8*i+j+1<pep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
331:     }
332:   }
333:   PetscViewerASCIIPrintf(viewer,"\n\n");
334:   return(0);
335: }

339: static PetscErrorCode PEPErrorView_DETAIL(PEP pep,PEPErrorType etype,PetscViewer viewer)
340: {
342:   PetscReal      error,re,im;
343:   PetscScalar    kr,ki;
344:   PetscInt       i;
345: #define EXLEN 30
346:   char           ex[EXLEN],sep[]=" ---------------------- --------------------\n";

349:   if (!pep->nconv) return(0);
350:   switch (etype) {
351:     case PEP_ERROR_ABSOLUTE:
352:       PetscSNPrintf(ex,EXLEN,"   ||P(k)x||");
353:       break;
354:     case PEP_ERROR_RELATIVE:
355:       PetscSNPrintf(ex,EXLEN,"||P(k)x||/||kx||");
356:       break;
357:     case PEP_ERROR_BACKWARD:
358:       PetscSNPrintf(ex,EXLEN,"    eta(x,k)");
359:       break;
360:   }
361:   PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep);
362:   for (i=0;i<pep->nconv;i++) {
363:     PEPGetEigenpair(pep,i,&kr,&ki,NULL,NULL);
364:     PEPComputeError(pep,i,etype,&error);
365: #if defined(PETSC_USE_COMPLEX)
366:     re = PetscRealPart(kr);
367:     im = PetscImaginaryPart(kr);
368: #else
369:     re = kr;
370:     im = ki;
371: #endif
372:     if (im!=0.0) {
373:       PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
374:     } else {
375:       PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
376:     }
377:   }
378:   PetscViewerASCIIPrintf(viewer,sep);
379:   return(0);
380: }

384: static PetscErrorCode PEPErrorView_MATLAB(PEP pep,PEPErrorType etype,PetscViewer viewer)
385: {
387:   PetscReal      error;
388:   PetscInt       i;
389:   const char     *name;

392:   PetscObjectGetName((PetscObject)pep,&name);
393:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
394:   for (i=0;i<pep->nconv;i++) {
395:     PEPComputeError(pep,i,etype,&error);
396:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",error);
397:   }
398:   PetscViewerASCIIPrintf(viewer,"];\n");
399:   return(0);
400: }

404: /*@C
405:    PEPErrorView - Displays the errors associated with the computed solution
406:    (as well as the eigenvalues).

408:    Collective on PEP

410:    Input Parameters:
411: +  pep    - the eigensolver context
412: .  etype  - error type
413: -  viewer - optional visualization context

415:    Options Database Key:
416: +  -pep_error_absolute - print absolute errors of each eigenpair
417: .  -pep_error_relative - print relative errors of each eigenpair
418: -  -pep_error_backward - print backward errors of each eigenpair

420:    Notes:
421:    By default, this function checks the error of all eigenpairs and prints
422:    the eigenvalues if all of them are below the requested tolerance.
423:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
424:    eigenvalues and corresponding errors is printed.

426:    Level: intermediate

428: .seealso: PEPSolve(), PEPValuesView(), PEPVectorsView()
429: @*/
430: PetscErrorCode PEPErrorView(PEP pep,PEPErrorType etype,PetscViewer viewer)
431: {
432:   PetscBool         isascii;
433:   PetscViewerFormat format;
434:   PetscErrorCode    ierr;

438:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
441:   PEPCheckSolved(pep,1);
442:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
443:   if (!isascii) return(0);

445:   PetscViewerGetFormat(viewer,&format);
446:   switch (format) {
447:     case PETSC_VIEWER_DEFAULT:
448:     case PETSC_VIEWER_ASCII_INFO:
449:       PEPErrorView_ASCII(pep,etype,viewer);
450:       break;
451:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
452:       PEPErrorView_DETAIL(pep,etype,viewer);
453:       break;
454:     case PETSC_VIEWER_ASCII_MATLAB:
455:       PEPErrorView_MATLAB(pep,etype,viewer);
456:       break;
457:     default:
458:       PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
459:   }
460:   return(0);
461: }

465: /*@
466:    PEPErrorViewFromOptions - Processes command line options to determine if/how
467:    the errors of the computed solution are to be viewed. 

469:    Collective on PEP

471:    Input Parameters:
472: .  pep - the eigensolver context

474:    Level: developer
475: @*/
476: PetscErrorCode PEPErrorViewFromOptions(PEP pep)
477: {
478:   PetscErrorCode    ierr;
479:   PetscViewer       viewer;
480:   PetscBool         flg;
481:   static PetscBool  incall = PETSC_FALSE;
482:   PetscViewerFormat format;

485:   if (incall) return(0);
486:   incall = PETSC_TRUE;
487:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_error_absolute",&viewer,&format,&flg);
488:   if (flg) {
489:     PetscViewerPushFormat(viewer,format);
490:     PEPErrorView(pep,PEP_ERROR_ABSOLUTE,viewer);
491:     PetscViewerPopFormat(viewer);
492:     PetscViewerDestroy(&viewer);
493:   }
494:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_error_relative",&viewer,&format,&flg);
495:   if (flg) {
496:     PetscViewerPushFormat(viewer,format);
497:     PEPErrorView(pep,PEP_ERROR_RELATIVE,viewer);
498:     PetscViewerPopFormat(viewer);
499:     PetscViewerDestroy(&viewer);
500:   }
501:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_error_backward",&viewer,&format,&flg);
502:   if (flg) {
503:     PetscViewerPushFormat(viewer,format);
504:     PEPErrorView(pep,PEP_ERROR_BACKWARD,viewer);
505:     PetscViewerPopFormat(viewer);
506:     PetscViewerDestroy(&viewer);
507:   }
508:   incall = PETSC_FALSE;
509:   return(0);
510: }

514: static PetscErrorCode PEPValuesView_DRAW(PEP pep,PetscViewer viewer)
515: {
517:   PetscDraw      draw;
518:   PetscDrawSP    drawsp;
519:   PetscReal      re,im;
520:   PetscInt       i,k;

523:   if (!pep->nconv) return(0);
524:   PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Computed Eigenvalues",PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
525:   PetscViewerDrawGetDraw(viewer,0,&draw);
526:   PetscDrawSPCreate(draw,1,&drawsp);
527:   for (i=0;i<pep->nconv;i++) {
528:     k = pep->perm[i];
529: #if defined(PETSC_USE_COMPLEX)
530:     re = PetscRealPart(pep->eigr[k]);
531:     im = PetscImaginaryPart(pep->eigr[k]);
532: #else
533:     re = pep->eigr[k];
534:     im = pep->eigi[k];
535: #endif
536:     PetscDrawSPAddPoint(drawsp,&re,&im);
537:   }
538:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
539:   PetscDrawSPSave(drawsp);
540:   PetscDrawSPDestroy(&drawsp);
541:   PetscViewerDestroy(&viewer);
542:   return(0);
543: }

547: static PetscErrorCode PEPValuesView_ASCII(PEP pep,PetscViewer viewer)
548: {
549:   PetscReal      re,im;
550:   PetscInt       i,k;

554:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
555:   for (i=0;i<pep->nconv;i++) {
556:     k = pep->perm[i];
557: #if defined(PETSC_USE_COMPLEX)
558:     re = PetscRealPart(pep->eigr[k]);
559:     im = PetscImaginaryPart(pep->eigr[k]);
560: #else
561:     re = pep->eigr[k];
562:     im = pep->eigi[k];
563: #endif
564:     if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
565:     if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
566:     if (im!=0.0) {
567:       PetscViewerASCIIPrintf(viewer,"   %.5f%+.5fi\n",(double)re,(double)im);
568:     } else {
569:       PetscViewerASCIIPrintf(viewer,"   %.5f\n",(double)re);
570:     }
571:   }
572:   PetscViewerASCIIPrintf(viewer,"\n");
573:   return(0);
574: }

578: static PetscErrorCode PEPValuesView_MATLAB(PEP pep,PetscViewer viewer)
579: {
581:   PetscInt       i,k;
582:   PetscReal      re,im;
583:   const char     *name;

586:   PetscObjectGetName((PetscObject)pep,&name);
587:   PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
588:   for (i=0;i<pep->nconv;i++) {
589:     k = pep->perm[i];
590: #if defined(PETSC_USE_COMPLEX)
591:     re = PetscRealPart(pep->eigr[k]);
592:     im = PetscImaginaryPart(pep->eigr[k]);
593: #else
594:     re = pep->eigr[k];
595:     im = pep->eigi[k];
596: #endif
597:     if (im!=0.0) {
598:       PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
599:     } else {
600:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
601:     }
602:   }
603:   PetscViewerASCIIPrintf(viewer,"];\n");
604:   return(0);
605: }

609: /*@C
610:    PEPValuesView - Displays the computed eigenvalues in a viewer.

612:    Collective on PEP

614:    Input Parameters:
615: +  pep    - the eigensolver context
616: -  viewer - the viewer

618:    Options Database Key:
619: .  -pep_view_values - print computed eigenvalues

621:    Level: intermediate

623: .seealso: PEPSolve(), PEPVectorsView(), PEPErrorView()
624: @*/
625: PetscErrorCode PEPValuesView(PEP pep,PetscViewer viewer)
626: {
627:   PetscBool         isascii,isdraw;
628:   PetscViewerFormat format;
629:   PetscErrorCode    ierr;

633:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
636:   PEPCheckSolved(pep,1);
637:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
638:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
639:   if (isdraw) {
640:     PEPValuesView_DRAW(pep,viewer);
641:   } else if (isascii) {
642:     PetscViewerGetFormat(viewer,&format);
643:     switch (format) {
644:       case PETSC_VIEWER_DEFAULT:
645:       case PETSC_VIEWER_ASCII_INFO:
646:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
647:         PEPValuesView_ASCII(pep,viewer);
648:         break;
649:       case PETSC_VIEWER_ASCII_MATLAB:
650:         PEPValuesView_MATLAB(pep,viewer);
651:         break;
652:       default:
653:         PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
654:     }
655:   }
656:   return(0);
657: }

661: /*@
662:    PEPValuesViewFromOptions - Processes command line options to determine if/how
663:    the computed eigenvalues are to be viewed. 

665:    Collective on PEP

667:    Input Parameters:
668: .  pep - the eigensolver context

670:    Level: developer
671: @*/
672: PetscErrorCode PEPValuesViewFromOptions(PEP pep)
673: {
674:   PetscErrorCode    ierr;
675:   PetscViewer       viewer;
676:   PetscBool         flg;
677:   static PetscBool  incall = PETSC_FALSE;
678:   PetscViewerFormat format;

681:   if (incall) return(0);
682:   incall = PETSC_TRUE;
683:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_view_values",&viewer,&format,&flg);
684:   if (flg) {
685:     PetscViewerPushFormat(viewer,format);
686:     PEPValuesView(pep,viewer);
687:     PetscViewerPopFormat(viewer);
688:     PetscViewerDestroy(&viewer);
689:   }
690:   incall = PETSC_FALSE;
691:   return(0);
692: }

696: /*@C
697:    PEPVectorsView - Outputs computed eigenvectors to a viewer.

699:    Collective on PEP

701:    Parameter:
702: +  pep    - the eigensolver context
703: -  viewer - the viewer

705:    Options Database Keys:
706: .  -pep_view_vectors - output eigenvectors.

708:    Note:
709:    If PETSc was configured with real scalars, complex conjugate eigenvectors
710:    will be viewed as two separate real vectors, one containing the real part
711:    and another one containing the imaginary part.

713:    Level: intermediate

715: .seealso: PEPSolve(), PEPValuesView(), PEPErrorView()
716: @*/
717: PetscErrorCode PEPVectorsView(PEP pep,PetscViewer viewer)
718: {
720:   PetscInt       i,k;
721:   Vec            x;
722: #define NMLEN 30
723:   char           vname[NMLEN];
724:   const char     *ename;

728:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
731:   PEPCheckSolved(pep,1);
732:   if (pep->nconv) {
733:     PetscObjectGetName((PetscObject)pep,&ename);
734:     PEPComputeVectors(pep);
735:     for (i=0;i<pep->nconv;i++) {
736:       k = pep->perm[i];
737:       PetscSNPrintf(vname,NMLEN,"V%d_%s",(int)i,ename);
738:       BVGetColumn(pep->V,k,&x);
739:       PetscObjectSetName((PetscObject)x,vname);
740:       VecView(x,viewer);
741:       BVRestoreColumn(pep->V,k,&x);
742:     }
743:   }
744:   return(0);
745: }

749: /*@
750:    PEPVectorsViewFromOptions - Processes command line options to determine if/how
751:    the computed eigenvectors are to be viewed. 

753:    Collective on PEP

755:    Input Parameters:
756: .  pep - the eigensolver context

758:    Level: developer
759: @*/
760: PetscErrorCode PEPVectorsViewFromOptions(PEP pep)
761: {
762:   PetscErrorCode    ierr;
763:   PetscViewer       viewer;
764:   PetscBool         flg = PETSC_FALSE;
765:   static PetscBool  incall = PETSC_FALSE;
766:   PetscViewerFormat format;

769:   if (incall) return(0);
770:   incall = PETSC_TRUE;
771:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->prefix,"-pep_view_vectors",&viewer,&format,&flg);
772:   if (flg) {
773:     PetscViewerPushFormat(viewer,format);
774:     PEPVectorsView(pep,viewer);
775:     PetscViewerPopFormat(viewer);
776:     PetscViewerDestroy(&viewer);
777:   }
778:   incall = PETSC_FALSE;
779:   return(0);
780: }