Actual source code: epsview.c

slepc-3.8.2 2017-12-01
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2017, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    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;

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

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

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

223:    Collective on EPS

225:    Parameter:
226: +  eps - the eigensolver context
227: -  viewer - the viewer to display the reason

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

232:    Level: intermediate

234: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber()
235: @*/
236: PetscErrorCode EPSReasonView(EPS eps,PetscViewer viewer)
237: {
239:   PetscBool      isAscii;

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

255: /*@
256:    EPSReasonViewFromOptions - Processes command line options to determine if/how
257:    the EPS converged reason is to be viewed.

259:    Collective on EPS

261:    Input Parameters:
262: .  eps - the eigensolver context

264:    Level: developer
265: @*/
266: PetscErrorCode EPSReasonViewFromOptions(EPS eps)
267: {
268:   PetscErrorCode    ierr;
269:   PetscViewer       viewer;
270:   PetscBool         flg;
271:   static PetscBool  incall = PETSC_FALSE;
272:   PetscViewerFormat format;

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

288: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
289: {
290:   PetscReal      error;
291:   PetscInt       i,j,k,nvals;

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

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

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

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

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

389: /*@C
390:    EPSErrorView - Displays the errors associated with the computed solution
391:    (as well as the eigenvalues).

393:    Collective on EPS

395:    Input Parameters:
396: +  eps    - the eigensolver context
397: .  etype  - error type
398: -  viewer - optional visualization context

400:    Options Database Key:
401: +  -eps_error_absolute - print absolute errors of each eigenpair
402: .  -eps_error_relative - print relative errors of each eigenpair
403: -  -eps_error_backward - print backward errors of each eigenpair

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

411:    Level: intermediate

413: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
414: @*/
415: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
416: {
417:   PetscBool         isascii;
418:   PetscViewerFormat format;
419:   PetscErrorCode    ierr;

423:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
426:   EPSCheckSolved(eps,1);
427:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
428:   if (!isascii) return(0);

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

448: /*@
449:    EPSErrorViewFromOptions - Processes command line options to determine if/how
450:    the errors of the computed solution are to be viewed.

452:    Collective on EPS

454:    Input Parameters:
455: .  eps - the eigensolver context

457:    Level: developer
458: @*/
459: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
460: {
461:   PetscErrorCode    ierr;
462:   PetscViewer       viewer;
463:   PetscBool         flg;
464:   static PetscBool  incall = PETSC_FALSE;
465:   PetscViewerFormat format;

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

495: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
496: {
498:   PetscDraw      draw;
499:   PetscDrawSP    drawsp;
500:   PetscReal      re,im;
501:   PetscInt       i,k;

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

525: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
526: {
527:   PetscInt       i,k;

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

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

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

571: /*@C
572:    EPSValuesView - Displays the computed eigenvalues in a viewer.

574:    Collective on EPS

576:    Input Parameters:
577: +  eps    - the eigensolver context
578: -  viewer - the viewer

580:    Options Database Key:
581: .  -eps_view_values - print computed eigenvalues

583:    Level: intermediate

585: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
586: @*/
587: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
588: {
589:   PetscBool         isascii,isdraw;
590:   PetscViewerFormat format;
591:   PetscErrorCode    ierr;

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

621: /*@
622:    EPSValuesViewFromOptions - Processes command line options to determine if/how
623:    the computed eigenvalues are to be viewed.

625:    Collective on EPS

627:    Input Parameters:
628: .  eps - the eigensolver context

630:    Level: developer
631: @*/
632: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
633: {
634:   PetscErrorCode    ierr;
635:   PetscViewer       viewer;
636:   PetscBool         flg;
637:   static PetscBool  incall = PETSC_FALSE;
638:   PetscViewerFormat format;

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

654: /*@C
655:    EPSVectorsView - Outputs computed eigenvectors to a viewer.

657:    Collective on EPS

659:    Parameter:
660: +  eps    - the eigensolver context
661: -  viewer - the viewer

663:    Options Database Keys:
664: .  -eps_view_vectors - output eigenvectors.

666:    Note:
667:    If PETSc was configured with real scalars, complex conjugate eigenvectors
668:    will be viewed as two separate real vectors, one containing the real part
669:    and another one containing the imaginary part.

671:    Level: intermediate

673: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
674: @*/
675: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
676: {
678:   PetscInt       i,k;
679:   Vec            x;
680: #define NMLEN 30
681:   char           vname[NMLEN];
682:   const char     *ename;

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

705: /*@
706:    EPSVectorsViewFromOptions - Processes command line options to determine if/how
707:    the computed eigenvectors are to be viewed.

709:    Collective on EPS

711:    Input Parameters:
712: .  eps - the eigensolver context

714:    Level: developer
715: @*/
716: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
717: {
718:   PetscErrorCode    ierr;
719:   PetscViewer       viewer;
720:   PetscBool         flg = PETSC_FALSE;
721:   static PetscBool  incall = PETSC_FALSE;
722:   PetscViewerFormat format;

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