Actual source code: pepview.c

slepc-3.14.0 2020-09-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2020, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    The PEP routines related to various viewers
 12: */

 14: #include <slepc/private/pepimpl.h>
 15: #include <petscdraw.h>

 17: /*@C
 18:    PEPView - Prints the PEP data structure.

 20:    Collective on pep

 22:    Input Parameters:
 23: +  pep - the polynomial eigenproblem solver context
 24: -  viewer - optional visualization context

 26:    Options Database Key:
 27: .  -pep_view -  Calls PEPView() at end of PEPSolve()

 29:    Note:
 30:    The available visualization contexts include
 31: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
 32: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
 33:          output where only the first processor opens
 34:          the file.  All other processors send their
 35:          data to the first processor to print.

 37:    The user can open an alternative visualization context with
 38:    PetscViewerASCIIOpen() - output to a specified file.

 40:    Level: beginner
 41: @*/
 42: PetscErrorCode PEPView(PEP pep,PetscViewer viewer)
 43: {
 45:   const char     *type=NULL;
 46:   char           str[50];
 47:   PetscBool      isascii,islinear,istrivial;
 48:   PetscInt       i;
 49:   PetscViewer    sviewer;

 53:   if (!viewer) {
 54:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
 55:   }

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

206: /*@C
207:    PEPViewFromOptions - View from options

209:    Collective on PEP

211:    Input Parameters:
212: +  pep  - the eigensolver context
213: .  obj  - optional object
214: -  name - command line option

216:    Level: intermediate

218: .seealso: PEPView(), PEPCreate()
219: @*/
220: PetscErrorCode PEPViewFromOptions(PEP pep,PetscObject obj,const char name[])
221: {

226:   PetscObjectViewFromOptions((PetscObject)pep,obj,name);
227:   return(0);
228: }

230: /*@C
231:    PEPConvergedReasonView - Displays the reason a PEP solve converged or diverged.

233:    Collective on pep

235:    Input Parameters:
236: +  pep - the eigensolver context
237: -  viewer - the viewer to display the reason

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

242:    Note:
243:    To change the format of the output call PetscViewerPushFormat(viewer,format) before
244:    this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
245:    display a reason if it fails. The latter can be set in the command line with
246:    -pep_converged_reason ::failed

248:    Level: intermediate

250: .seealso: PEPSetConvergenceTest(), PEPSetTolerances(), PEPGetIterationNumber(), PEPConvergedReasonViewFromOptions()
251: @*/
252: PetscErrorCode PEPConvergedReasonView(PEP pep,PetscViewer viewer)
253: {
254:   PetscErrorCode    ierr;
255:   PetscBool         isAscii;
256:   PetscViewerFormat format;

259:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
260:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
261:   if (isAscii) {
262:     PetscViewerGetFormat(viewer,&format);
263:     PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
264:     if (pep->reason > 0 && format != PETSC_VIEWER_FAILED) {
265:       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);
266:     } else if (pep->reason <= 0) {
267:       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);
268:     }
269:     PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
270:   }
271:   return(0);
272: }

274: /*@
275:    PEPConvergedReasonViewFromOptions - Processes command line options to determine if/how
276:    the PEP converged reason is to be viewed.

278:    Collective on pep

280:    Input Parameter:
281: .  pep - the eigensolver context

283:    Level: developer

285: .seealso: PEPConvergedReasonView()
286: @*/
287: PetscErrorCode PEPConvergedReasonViewFromOptions(PEP pep)
288: {
289:   PetscErrorCode    ierr;
290:   PetscViewer       viewer;
291:   PetscBool         flg;
292:   static PetscBool  incall = PETSC_FALSE;
293:   PetscViewerFormat format;

296:   if (incall) return(0);
297:   incall = PETSC_TRUE;
298:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_converged_reason",&viewer,&format,&flg);
299:   if (flg) {
300:     PetscViewerPushFormat(viewer,format);
301:     PEPConvergedReasonView(pep,viewer);
302:     PetscViewerPopFormat(viewer);
303:     PetscViewerDestroy(&viewer);
304:   }
305:   incall = PETSC_FALSE;
306:   return(0);
307: }

309: static PetscErrorCode PEPErrorView_ASCII(PEP pep,PEPErrorType etype,PetscViewer viewer)
310: {
311:   PetscReal      error;
312:   PetscInt       i,j,k;

316:   if (pep->nconv<pep->nev) {
317:     PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",pep->nev);
318:     return(0);
319:   }
320:   for (i=0;i<pep->nev;i++) {
321:     PEPComputeError(pep,i,etype,&error);
322:     if (error>=5.0*pep->tol) {
323:       PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",pep->nev);
324:       return(0);
325:     }
326:   }
327:   PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
328:   for (i=0;i<=(pep->nev-1)/8;i++) {
329:     PetscViewerASCIIPrintf(viewer,"\n     ");
330:     for (j=0;j<PetscMin(8,pep->nev-8*i);j++) {
331:       k = pep->perm[8*i+j];
332:       SlepcPrintEigenvalueASCII(viewer,pep->eigr[k],pep->eigi[k]);
333:       if (8*i+j+1<pep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
334:     }
335:   }
336:   PetscViewerASCIIPrintf(viewer,"\n\n");
337:   return(0);
338: }

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

349:   if (!pep->nconv) return(0);
350:   switch (etype) {
351:     case PEP_ERROR_ABSOLUTE:
352:       PetscSNPrintf(ex,sizeof(ex),"   ||P(k)x||");
353:       break;
354:     case PEP_ERROR_RELATIVE:
355:       PetscSNPrintf(ex,sizeof(ex),"||P(k)x||/||kx||");
356:       break;
357:     case PEP_ERROR_BACKWARD:
358:       PetscSNPrintf(ex,sizeof(ex),"    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: }

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

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

400: /*@C
401:    PEPErrorView - Displays the errors associated with the computed solution
402:    (as well as the eigenvalues).

404:    Collective on pep

406:    Input Parameters:
407: +  pep    - the eigensolver context
408: .  etype  - error type
409: -  viewer - optional visualization context

411:    Options Database Key:
412: +  -pep_error_absolute - print absolute errors of each eigenpair
413: .  -pep_error_relative - print relative errors of each eigenpair
414: -  -pep_error_backward - print backward errors of each eigenpair

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

422:    Level: intermediate

424: .seealso: PEPSolve(), PEPValuesView(), PEPVectorsView()
425: @*/
426: PetscErrorCode PEPErrorView(PEP pep,PEPErrorType etype,PetscViewer viewer)
427: {
428:   PetscBool         isascii;
429:   PetscViewerFormat format;
430:   PetscErrorCode    ierr;

434:   if (!viewer) {
435:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
436:   }
439:   PEPCheckSolved(pep,1);
440:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
441:   if (!isascii) return(0);

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

461: /*@
462:    PEPErrorViewFromOptions - Processes command line options to determine if/how
463:    the errors of the computed solution are to be viewed.

465:    Collective on pep

467:    Input Parameter:
468: .  pep - the eigensolver context

470:    Level: developer
471: @*/
472: PetscErrorCode PEPErrorViewFromOptions(PEP pep)
473: {
474:   PetscErrorCode    ierr;
475:   PetscViewer       viewer;
476:   PetscBool         flg;
477:   static PetscBool  incall = PETSC_FALSE;
478:   PetscViewerFormat format;

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

508: static PetscErrorCode PEPValuesView_DRAW(PEP pep,PetscViewer viewer)
509: {
511:   PetscDraw      draw;
512:   PetscDrawSP    drawsp;
513:   PetscReal      re,im;
514:   PetscInt       i,k;

517:   if (!pep->nconv) return(0);
518:   PetscViewerDrawGetDraw(viewer,0,&draw);
519:   PetscDrawSetTitle(draw,"Computed Eigenvalues");
520:   PetscDrawSPCreate(draw,1,&drawsp);
521:   for (i=0;i<pep->nconv;i++) {
522:     k = pep->perm[i];
523: #if defined(PETSC_USE_COMPLEX)
524:     re = PetscRealPart(pep->eigr[k]);
525:     im = PetscImaginaryPart(pep->eigr[k]);
526: #else
527:     re = pep->eigr[k];
528:     im = pep->eigi[k];
529: #endif
530:     PetscDrawSPAddPoint(drawsp,&re,&im);
531:   }
532:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
533:   PetscDrawSPSave(drawsp);
534:   PetscDrawSPDestroy(&drawsp);
535:   return(0);
536: }

538: static PetscErrorCode PEPValuesView_BINARY(PEP pep,PetscViewer viewer)
539: {
540: #if defined(PETSC_HAVE_COMPLEX)
541:   PetscInt       i,k;
542:   PetscComplex   *ev;
544: #endif

547: #if defined(PETSC_HAVE_COMPLEX)
548:   PetscMalloc1(pep->nconv,&ev);
549:   for (i=0;i<pep->nconv;i++) {
550:     k = pep->perm[i];
551: #if defined(PETSC_USE_COMPLEX)
552:     ev[i] = pep->eigr[k];
553: #else
554:     ev[i] = PetscCMPLX(pep->eigr[k],pep->eigi[k]);
555: #endif
556:   }
557:   PetscViewerBinaryWrite(viewer,ev,pep->nconv,PETSC_COMPLEX);
558:   PetscFree(ev);
559: #endif
560:   return(0);
561: }

563: static PetscErrorCode PEPValuesView_ASCII(PEP pep,PetscViewer viewer)
564: {
565:   PetscInt       i,k;

569:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
570:   for (i=0;i<pep->nconv;i++) {
571:     k = pep->perm[i];
572:     PetscViewerASCIIPrintf(viewer,"   ");
573:     SlepcPrintEigenvalueASCII(viewer,pep->eigr[k],pep->eigi[k]);
574:     PetscViewerASCIIPrintf(viewer,"\n");
575:   }
576:   PetscViewerASCIIPrintf(viewer,"\n");
577:   return(0);
578: }

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

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

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,isbinary;
628:   PetscViewerFormat format;
629:   PetscErrorCode    ierr;

633:   if (!viewer) {
634:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
635:   }
638:   PEPCheckSolved(pep,1);
639:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
640:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
641:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
642:   if (isdraw) {
643:     PEPValuesView_DRAW(pep,viewer);
644:   } else if (isbinary) {
645:     PEPValuesView_BINARY(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: }

664: /*@
665:    PEPValuesViewFromOptions - Processes command line options to determine if/how
666:    the computed eigenvalues are to be viewed.

668:    Collective on pep

670:    Input Parameter:
671: .  pep - the eigensolver context

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

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

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

700:    Collective on pep

702:    Input Parameters:
703: +  pep    - the eigensolver context
704: -  viewer - the viewer

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

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

714:    Level: intermediate

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

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

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 Parameter:
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)->options,((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: }