Actual source code: epsview.c

slepc-3.10.0 2018-09-18
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2018, 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:    EPS routines related to various viewers
 12: */

 14: #include <slepc/private/epsimpl.h>      /*I "slepceps.h" I*/
 15: #include <petscdraw.h>

 17: /*@C
 18:    EPSView - Prints the EPS data structure.

 20:    Collective on EPS

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

 26:    Options Database Key:
 27: .  -eps_view -  Calls EPSView() at end of EPSSolve()

 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: STView(), PetscViewerASCIIOpen()
 43: @*/
 44: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
 45: {
 47:   const char     *type=NULL,*extr=NULL,*bal=NULL;
 48:   char           str[50];
 49:   PetscBool      isascii,isexternal,istrivial;
 50:   PetscInt       tabs;

 54:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));

 58: #if defined(PETSC_USE_COMPLEX)
 59: #define HERM "hermitian"
 60: #else
 61: #define HERM "symmetric"
 62: #endif
 63:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 64:   if (isascii) {
 65:     PetscViewerASCIIGetTab(viewer,&tabs);
 66:     PetscViewerASCIISetTab(viewer,((PetscObject)eps)->tablevel);
 67:     PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer);
 68:     if (eps->ops->view) {
 69:       PetscViewerASCIIPushTab(viewer);
 70:       (*eps->ops->view)(eps,viewer);
 71:       PetscViewerASCIIPopTab(viewer);
 72:     }
 73:     if (eps->problem_type) {
 74:       switch (eps->problem_type) {
 75:         case EPS_HEP:    type = HERM " eigenvalue problem"; break;
 76:         case EPS_GHEP:   type = "generalized " HERM " eigenvalue problem"; break;
 77:         case EPS_NHEP:   type = "non-" HERM " eigenvalue problem"; break;
 78:         case EPS_GNHEP:  type = "generalized non-" HERM " eigenvalue problem"; break;
 79:         case EPS_PGNHEP: type = "generalized non-" HERM " eigenvalue problem with " HERM " positive definite B"; break;
 80:         case EPS_GHIEP:  type = "generalized " HERM "-indefinite eigenvalue problem"; break;
 81:       }
 82:     } else type = "not yet set";
 83:     PetscViewerASCIIPrintf(viewer,"  problem type: %s\n",type);
 84:     if (eps->extraction) {
 85:       switch (eps->extraction) {
 86:         case EPS_RITZ:              extr = "Rayleigh-Ritz"; break;
 87:         case EPS_HARMONIC:          extr = "harmonic Ritz"; break;
 88:         case EPS_HARMONIC_RELATIVE: extr = "relative harmonic Ritz"; break;
 89:         case EPS_HARMONIC_RIGHT:    extr = "right harmonic Ritz"; break;
 90:         case EPS_HARMONIC_LARGEST:  extr = "largest harmonic Ritz"; break;
 91:         case EPS_REFINED:           extr = "refined Ritz"; break;
 92:         case EPS_REFINED_HARMONIC:  extr = "refined harmonic Ritz"; break;
 93:       }
 94:       PetscViewerASCIIPrintf(viewer,"  extraction type: %s\n",extr);
 95:     }
 96:     if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
 97:       switch (eps->balance) {
 98:         case EPS_BALANCE_NONE:    break;
 99:         case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
100:         case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
101:         case EPS_BALANCE_USER:    bal = "user-defined matrix"; break;
102:       }
103:       PetscViewerASCIIPrintf(viewer,"  balancing enabled: %s",bal);
104:       PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
105:       if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) {
106:         PetscViewerASCIIPrintf(viewer,", with its=%D",eps->balance_its);
107:       }
108:       if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) {
109:         PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff);
110:       }
111:       PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
112:       PetscViewerASCIIPrintf(viewer,"\n");
113:     }
114:     PetscViewerASCIIPrintf(viewer,"  selected portion of the spectrum: ");
115:     SlepcSNPrintfScalar(str,50,eps->target,PETSC_FALSE);
116:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
117:     if (!eps->which) {
118:       PetscViewerASCIIPrintf(viewer,"not yet set\n");
119:     } else switch (eps->which) {
120:       case EPS_WHICH_USER:
121:         PetscViewerASCIIPrintf(viewer,"user defined\n");
122:         break;
123:       case EPS_TARGET_MAGNITUDE:
124:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
125:         break;
126:       case EPS_TARGET_REAL:
127:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
128:         break;
129:       case EPS_TARGET_IMAGINARY:
130:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
131:         break;
132:       case EPS_LARGEST_MAGNITUDE:
133:         PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
134:         break;
135:       case EPS_SMALLEST_MAGNITUDE:
136:         PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
137:         break;
138:       case EPS_LARGEST_REAL:
139:         PetscViewerASCIIPrintf(viewer,"largest real parts\n");
140:         break;
141:       case EPS_SMALLEST_REAL:
142:         PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
143:         break;
144:       case EPS_LARGEST_IMAGINARY:
145:         PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
146:         break;
147:       case EPS_SMALLEST_IMAGINARY:
148:         PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
149:         break;
150:       case EPS_ALL:
151:         if (eps->inta || eps->intb) {
152:           PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb);
153:         } else {
154:           PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
155:         }
156:         break;
157:     }
158:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
159:     if (eps->purify) {
160:       PetscViewerASCIIPrintf(viewer,"  postprocessing eigenvectors with purification\n");
161:     }
162:     if (eps->trueres) {
163:       PetscViewerASCIIPrintf(viewer,"  computing true residuals explicitly\n");
164:     }
165:     if (eps->trackall) {
166:       PetscViewerASCIIPrintf(viewer,"  computing all residuals (for tracking convergence)\n");
167:     }
168:     PetscViewerASCIIPrintf(viewer,"  number of eigenvalues (nev): %D\n",eps->nev);
169:     PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %D\n",eps->ncv);
170:     PetscViewerASCIIPrintf(viewer,"  maximum dimension of projected problem (mpd): %D\n",eps->mpd);
171:     PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %D\n",eps->max_it);
172:     PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)eps->tol);
173:     PetscViewerASCIIPrintf(viewer,"  convergence test: ");
174:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
175:     switch (eps->conv) {
176:     case EPS_CONV_ABS:
177:       PetscViewerASCIIPrintf(viewer,"absolute\n");break;
178:     case EPS_CONV_REL:
179:       PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
180:     case EPS_CONV_NORM:
181:       PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");
182:       PetscViewerASCIIPrintf(viewer,"  computed matrix norms: norm(A)=%g",(double)eps->nrma);
183:       if (eps->isgeneralized) {
184:         PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb);
185:       }
186:       PetscViewerASCIIPrintf(viewer,"\n");
187:       break;
188:     case EPS_CONV_USER:
189:       PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
190:     }
191:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
192:     if (eps->nini) {
193:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial space: %D\n",PetscAbs(eps->nini));
194:     }
195:     if (eps->nds) {
196:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided deflation space: %D\n",PetscAbs(eps->nds));
197:     }
198:     PetscViewerASCIISetTab(viewer,tabs);
199:   } else {
200:     if (eps->ops->view) {
201:       (*eps->ops->view)(eps,viewer);
202:     }
203:   }
204:   PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSBLZPACK,EPSFEAST,EPSPRIMME,EPSTRLAN,"");
205:   if (!isexternal) {
206:     PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
207:     if (!eps->V) { EPSGetBV(eps,&eps->V); }
208:     BVView(eps->V,viewer);
209:     if (eps->rg) {
210:       RGIsTrivial(eps->rg,&istrivial);
211:       if (!istrivial) { RGView(eps->rg,viewer); }
212:     }
213:     if (eps->useds) {
214:       if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
215:       DSView(eps->ds,viewer);
216:     }
217:     PetscViewerPopFormat(viewer);
218:   }
219:   if (!eps->st) { EPSGetST(eps,&eps->st); }
220:   STView(eps->st,viewer);
221:   return(0);
222: }

224: /*@C
225:    EPSReasonView - Displays the reason an EPS solve converged or diverged.

227:    Collective on EPS

229:    Parameter:
230: +  eps - the eigensolver context
231: -  viewer - the viewer to display the reason

233:    Options Database Keys:
234: .  -eps_converged_reason - print reason for convergence, and number of iterations

236:    Level: intermediate

238: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber()
239: @*/
240: PetscErrorCode EPSReasonView(EPS eps,PetscViewer viewer)
241: {
243:   PetscBool      isAscii;

246:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
247:   if (isAscii) {
248:     PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
249:     if (eps->reason > 0) {
250:       PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve converged (%D eigenpair%s) due to %s; iterations %D\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",eps->nconv,(eps->nconv>1)?"s":"",EPSConvergedReasons[eps->reason],eps->its);
251:     } else {
252:       PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve did not converge due to %s; iterations %D\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",EPSConvergedReasons[eps->reason],eps->its);
253:     }
254:     PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
255:   }
256:   return(0);
257: }

259: /*@
260:    EPSReasonViewFromOptions - Processes command line options to determine if/how
261:    the EPS converged reason is to be viewed.

263:    Collective on EPS

265:    Input Parameters:
266: .  eps - the eigensolver context

268:    Level: developer
269: @*/
270: PetscErrorCode EPSReasonViewFromOptions(EPS eps)
271: {
272:   PetscErrorCode    ierr;
273:   PetscViewer       viewer;
274:   PetscBool         flg;
275:   static PetscBool  incall = PETSC_FALSE;
276:   PetscViewerFormat format;

279:   if (incall) return(0);
280:   incall = PETSC_TRUE;
281:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg);
282:   if (flg) {
283:     PetscViewerPushFormat(viewer,format);
284:     EPSReasonView(eps,viewer);
285:     PetscViewerPopFormat(viewer);
286:     PetscViewerDestroy(&viewer);
287:   }
288:   incall = PETSC_FALSE;
289:   return(0);
290: }

292: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
293: {
294:   PetscReal      error;
295:   PetscInt       i,j,k,nvals;

299:   nvals = (eps->which==EPS_ALL)? eps->nconv: eps->nev;
300:   if (eps->which!=EPS_ALL && eps->nconv<nvals) {
301:     PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
302:     return(0);
303:   }
304:   if (eps->which==EPS_ALL && !nvals) {
305:     PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
306:     return(0);
307:   }
308:   for (i=0;i<nvals;i++) {
309:     EPSComputeError(eps,i,etype,&error);
310:     if (error>=5.0*eps->tol) {
311:       PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
312:       return(0);
313:     }
314:   }
315:   if (eps->which==EPS_ALL) {
316:     PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
317:   } else {
318:     PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
319:   }
320:   for (i=0;i<=(nvals-1)/8;i++) {
321:     PetscViewerASCIIPrintf(viewer,"\n     ");
322:     for (j=0;j<PetscMin(8,nvals-8*i);j++) {
323:       k = eps->perm[8*i+j];
324:       SlepcPrintEigenvalueASCII(eps->eigr[k],eps->eigi[k]);
325:       if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
326:     }
327:   }
328:   PetscViewerASCIIPrintf(viewer,"\n\n");
329:   return(0);
330: }

332: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
333: {
335:   PetscReal      error,re,im;
336:   PetscScalar    kr,ki;
337:   PetscInt       i;
338: #define EXLEN 30
339:   char           ex[EXLEN],sep[]=" ---------------------- --------------------\n";

342:   if (!eps->nconv) return(0);
343:   switch (etype) {
344:     case EPS_ERROR_ABSOLUTE:
345:       PetscSNPrintf(ex,EXLEN,"   ||Ax-k%sx||",eps->isgeneralized?"B":"");
346:       break;
347:     case EPS_ERROR_RELATIVE:
348:       PetscSNPrintf(ex,EXLEN,"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":"");
349:       break;
350:     case EPS_ERROR_BACKWARD:
351:       PetscSNPrintf(ex,EXLEN,"    eta(x,k)");
352:       break;
353:   }
354:   PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep);
355:   for (i=0;i<eps->nconv;i++) {
356:     EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
357:     EPSComputeError(eps,i,etype,&error);
358: #if defined(PETSC_USE_COMPLEX)
359:     re = PetscRealPart(kr);
360:     im = PetscImaginaryPart(kr);
361: #else
362:     re = kr;
363:     im = ki;
364: #endif
365:     if (im!=0.0) {
366:       PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
367:     } else {
368:       PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
369:     }
370:   }
371:   PetscViewerASCIIPrintf(viewer,sep);
372:   return(0);
373: }

375: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
376: {
378:   PetscReal      error;
379:   PetscInt       i;
380:   const char     *name;

383:   PetscObjectGetName((PetscObject)eps,&name);
384:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
385:   for (i=0;i<eps->nconv;i++) {
386:     EPSComputeError(eps,i,etype,&error);
387:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
388:   }
389:   PetscViewerASCIIPrintf(viewer,"];\n");
390:   return(0);
391: }

393: /*@C
394:    EPSErrorView - Displays the errors associated with the computed solution
395:    (as well as the eigenvalues).

397:    Collective on EPS

399:    Input Parameters:
400: +  eps    - the eigensolver context
401: .  etype  - error type
402: -  viewer - optional visualization context

404:    Options Database Key:
405: +  -eps_error_absolute - print absolute errors of each eigenpair
406: .  -eps_error_relative - print relative errors of each eigenpair
407: -  -eps_error_backward - print backward errors of each eigenpair

409:    Notes:
410:    By default, this function checks the error of all eigenpairs and prints
411:    the eigenvalues if all of them are below the requested tolerance.
412:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
413:    eigenvalues and corresponding errors is printed.

415:    Level: intermediate

417: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
418: @*/
419: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
420: {
421:   PetscBool         isascii;
422:   PetscViewerFormat format;
423:   PetscErrorCode    ierr;

427:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
430:   EPSCheckSolved(eps,1);
431:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
432:   if (!isascii) return(0);

434:   PetscViewerGetFormat(viewer,&format);
435:   switch (format) {
436:     case PETSC_VIEWER_DEFAULT:
437:     case PETSC_VIEWER_ASCII_INFO:
438:       EPSErrorView_ASCII(eps,etype,viewer);
439:       break;
440:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
441:       EPSErrorView_DETAIL(eps,etype,viewer);
442:       break;
443:     case PETSC_VIEWER_ASCII_MATLAB:
444:       EPSErrorView_MATLAB(eps,etype,viewer);
445:       break;
446:     default:
447:       PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
448:   }
449:   return(0);
450: }

452: /*@
453:    EPSErrorViewFromOptions - Processes command line options to determine if/how
454:    the errors of the computed solution are to be viewed.

456:    Collective on EPS

458:    Input Parameters:
459: .  eps - the eigensolver context

461:    Level: developer
462: @*/
463: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
464: {
465:   PetscErrorCode    ierr;
466:   PetscViewer       viewer;
467:   PetscBool         flg;
468:   static PetscBool  incall = PETSC_FALSE;
469:   PetscViewerFormat format;

472:   if (incall) return(0);
473:   incall = PETSC_TRUE;
474:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg);
475:   if (flg) {
476:     PetscViewerPushFormat(viewer,format);
477:     EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer);
478:     PetscViewerPopFormat(viewer);
479:     PetscViewerDestroy(&viewer);
480:   }
481:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg);
482:   if (flg) {
483:     PetscViewerPushFormat(viewer,format);
484:     EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer);
485:     PetscViewerPopFormat(viewer);
486:     PetscViewerDestroy(&viewer);
487:   }
488:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg);
489:   if (flg) {
490:     PetscViewerPushFormat(viewer,format);
491:     EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer);
492:     PetscViewerPopFormat(viewer);
493:     PetscViewerDestroy(&viewer);
494:   }
495:   incall = PETSC_FALSE;
496:   return(0);
497: }

499: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
500: {
502:   PetscDraw      draw;
503:   PetscDrawSP    drawsp;
504:   PetscReal      re,im;
505:   PetscInt       i,k;

508:   if (!eps->nconv) return(0);
509:   PetscViewerDrawGetDraw(viewer,0,&draw);
510:   PetscDrawSetTitle(draw,"Computed Eigenvalues");
511:   PetscDrawSPCreate(draw,1,&drawsp);
512:   for (i=0;i<eps->nconv;i++) {
513:     k = eps->perm[i];
514: #if defined(PETSC_USE_COMPLEX)
515:     re = PetscRealPart(eps->eigr[k]);
516:     im = PetscImaginaryPart(eps->eigr[k]);
517: #else
518:     re = eps->eigr[k];
519:     im = eps->eigi[k];
520: #endif
521:     PetscDrawSPAddPoint(drawsp,&re,&im);
522:   }
523:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
524:   PetscDrawSPSave(drawsp);
525:   PetscDrawSPDestroy(&drawsp);
526:   return(0);
527: }

529: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
530: {
531:   PetscInt       i,k;

535:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
536:   for (i=0;i<eps->nconv;i++) {
537:     k = eps->perm[i];
538:     PetscViewerASCIIPrintf(viewer,"   ");
539:     SlepcPrintEigenvalueASCII(eps->eigr[k],eps->eigi[k]);
540:     PetscViewerASCIIPrintf(viewer,"\n");
541:   }
542:   PetscViewerASCIIPrintf(viewer,"\n");
543:   return(0);
544: }

546: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
547: {
549:   PetscInt       i,k;
550:   PetscReal      re,im;
551:   const char     *name;

554:   PetscObjectGetName((PetscObject)eps,&name);
555:   PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
556:   for (i=0;i<eps->nconv;i++) {
557:     k = eps->perm[i];
558: #if defined(PETSC_USE_COMPLEX)
559:     re = PetscRealPart(eps->eigr[k]);
560:     im = PetscImaginaryPart(eps->eigr[k]);
561: #else
562:     re = eps->eigr[k];
563:     im = eps->eigi[k];
564: #endif
565:     if (im!=0.0) {
566:       PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
567:     } else {
568:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
569:     }
570:   }
571:   PetscViewerASCIIPrintf(viewer,"];\n");
572:   return(0);
573: }

575: /*@C
576:    EPSValuesView - Displays the computed eigenvalues in a viewer.

578:    Collective on EPS

580:    Input Parameters:
581: +  eps    - the eigensolver context
582: -  viewer - the viewer

584:    Options Database Key:
585: .  -eps_view_values - print computed eigenvalues

587:    Level: intermediate

589: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
590: @*/
591: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
592: {
593:   PetscBool         isascii,isdraw;
594:   PetscViewerFormat format;
595:   PetscErrorCode    ierr;

599:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
602:   EPSCheckSolved(eps,1);
603:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
604:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
605:   if (isdraw) {
606:     EPSValuesView_DRAW(eps,viewer);
607:   } else if (isascii) {
608:     PetscViewerGetFormat(viewer,&format);
609:     switch (format) {
610:       case PETSC_VIEWER_DEFAULT:
611:       case PETSC_VIEWER_ASCII_INFO:
612:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
613:         EPSValuesView_ASCII(eps,viewer);
614:         break;
615:       case PETSC_VIEWER_ASCII_MATLAB:
616:         EPSValuesView_MATLAB(eps,viewer);
617:         break;
618:       default:
619:         PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
620:     }
621:   }
622:   return(0);
623: }

625: /*@
626:    EPSValuesViewFromOptions - Processes command line options to determine if/how
627:    the computed eigenvalues are to be viewed.

629:    Collective on EPS

631:    Input Parameters:
632: .  eps - the eigensolver context

634:    Level: developer
635: @*/
636: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
637: {
638:   PetscErrorCode    ierr;
639:   PetscViewer       viewer;
640:   PetscBool         flg;
641:   static PetscBool  incall = PETSC_FALSE;
642:   PetscViewerFormat format;

645:   if (incall) return(0);
646:   incall = PETSC_TRUE;
647:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg);
648:   if (flg) {
649:     PetscViewerPushFormat(viewer,format);
650:     EPSValuesView(eps,viewer);
651:     PetscViewerPopFormat(viewer);
652:     PetscViewerDestroy(&viewer);
653:   }
654:   incall = PETSC_FALSE;
655:   return(0);
656: }

658: /*@C
659:    EPSVectorsView - Outputs computed eigenvectors to a viewer.

661:    Collective on EPS

663:    Parameter:
664: +  eps    - the eigensolver context
665: -  viewer - the viewer

667:    Options Database Keys:
668: .  -eps_view_vectors - output eigenvectors.

670:    Note:
671:    If PETSc was configured with real scalars, complex conjugate eigenvectors
672:    will be viewed as two separate real vectors, one containing the real part
673:    and another one containing the imaginary part.

675:    Level: intermediate

677: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
678: @*/
679: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
680: {
682:   PetscInt       i,k;
683:   Vec            x;
684: #define NMLEN 30
685:   char           vname[NMLEN];
686:   const char     *ename;

690:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
693:   EPSCheckSolved(eps,1);
694:   if (eps->nconv) {
695:     PetscObjectGetName((PetscObject)eps,&ename);
696:     EPSComputeVectors(eps);
697:     for (i=0;i<eps->nconv;i++) {
698:       k = eps->perm[i];
699:       PetscSNPrintf(vname,NMLEN,"V%d_%s",(int)i,ename);
700:       BVGetColumn(eps->V,k,&x);
701:       PetscObjectSetName((PetscObject)x,vname);
702:       VecView(x,viewer);
703:       BVRestoreColumn(eps->V,k,&x);
704:     }
705:   }
706:   return(0);
707: }

709: /*@
710:    EPSVectorsViewFromOptions - Processes command line options to determine if/how
711:    the computed eigenvectors are to be viewed.

713:    Collective on EPS

715:    Input Parameters:
716: .  eps - the eigensolver context

718:    Level: developer
719: @*/
720: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
721: {
722:   PetscErrorCode    ierr;
723:   PetscViewer       viewer;
724:   PetscBool         flg = PETSC_FALSE;
725:   static PetscBool  incall = PETSC_FALSE;
726:   PetscViewerFormat format;

729:   if (incall) return(0);
730:   incall = PETSC_TRUE;
731:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg);
732:   if (flg) {
733:     PetscViewerPushFormat(viewer,format);
734:     EPSVectorsView(eps,viewer);
735:     PetscViewerPopFormat(viewer);
736:     PetscViewerDestroy(&viewer);
737:   }
738:   incall = PETSC_FALSE;
739:   return(0);
740: }