Actual source code: epsview.c

slepc-3.15.2 2021-09-20
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2021, 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>
 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()
 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) {
 54:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
 55:   }

 59:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 60:   if (isascii) {
 61:     PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer);
 62:     if (eps->ops->view) {
 63:       PetscViewerASCIIPushTab(viewer);
 64:       (*eps->ops->view)(eps,viewer);
 65:       PetscViewerASCIIPopTab(viewer);
 66:     }
 67:     if (eps->problem_type) {
 68:       switch (eps->problem_type) {
 69:         case EPS_HEP:    type = SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
 70:         case EPS_GHEP:   type = "generalized " SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
 71:         case EPS_NHEP:   type = "non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
 72:         case EPS_GNHEP:  type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
 73:         case EPS_PGNHEP: type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem with " SLEPC_STRING_HERMITIAN " positive definite B"; break;
 74:         case EPS_GHIEP:  type = "generalized " SLEPC_STRING_HERMITIAN "-indefinite eigenvalue problem"; break;
 75:       }
 76:     } else type = "not yet set";
 77:     PetscViewerASCIIPrintf(viewer,"  problem type: %s\n",type);
 78:     if (eps->extraction) {
 79:       switch (eps->extraction) {
 80:         case EPS_RITZ:              extr = "Rayleigh-Ritz"; break;
 81:         case EPS_HARMONIC:          extr = "harmonic Ritz"; break;
 82:         case EPS_HARMONIC_RELATIVE: extr = "relative harmonic Ritz"; break;
 83:         case EPS_HARMONIC_RIGHT:    extr = "right harmonic Ritz"; break;
 84:         case EPS_HARMONIC_LARGEST:  extr = "largest harmonic Ritz"; break;
 85:         case EPS_REFINED:           extr = "refined Ritz"; break;
 86:         case EPS_REFINED_HARMONIC:  extr = "refined harmonic Ritz"; break;
 87:       }
 88:       PetscViewerASCIIPrintf(viewer,"  extraction type: %s\n",extr);
 89:     }
 90:     if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
 91:       switch (eps->balance) {
 92:         case EPS_BALANCE_NONE:    break;
 93:         case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
 94:         case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
 95:         case EPS_BALANCE_USER:    bal = "user-defined matrix"; break;
 96:       }
 97:       PetscViewerASCIIPrintf(viewer,"  balancing enabled: %s",bal);
 98:       PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
 99:       if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) {
100:         PetscViewerASCIIPrintf(viewer,", with its=%D",eps->balance_its);
101:       }
102:       if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) {
103:         PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff);
104:       }
105:       PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
106:       PetscViewerASCIIPrintf(viewer,"\n");
107:     }
108:     PetscViewerASCIIPrintf(viewer,"  selected portion of the spectrum: ");
109:     SlepcSNPrintfScalar(str,sizeof(str),eps->target,PETSC_FALSE);
110:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
111:     if (!eps->which) {
112:       PetscViewerASCIIPrintf(viewer,"not yet set\n");
113:     } else switch (eps->which) {
114:       case EPS_WHICH_USER:
115:         PetscViewerASCIIPrintf(viewer,"user defined\n");
116:         break;
117:       case EPS_TARGET_MAGNITUDE:
118:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
119:         break;
120:       case EPS_TARGET_REAL:
121:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
122:         break;
123:       case EPS_TARGET_IMAGINARY:
124:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
125:         break;
126:       case EPS_LARGEST_MAGNITUDE:
127:         PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
128:         break;
129:       case EPS_SMALLEST_MAGNITUDE:
130:         PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
131:         break;
132:       case EPS_LARGEST_REAL:
133:         PetscViewerASCIIPrintf(viewer,"largest real parts\n");
134:         break;
135:       case EPS_SMALLEST_REAL:
136:         PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
137:         break;
138:       case EPS_LARGEST_IMAGINARY:
139:         PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
140:         break;
141:       case EPS_SMALLEST_IMAGINARY:
142:         PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
143:         break;
144:       case EPS_ALL:
145:         if (eps->inta || eps->intb) {
146:           PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb);
147:         } else {
148:           PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
149:         }
150:         break;
151:     }
152:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
153:     if (eps->twosided && eps->problem_type!=EPS_HEP && eps->problem_type!=EPS_GHEP) {
154:       PetscViewerASCIIPrintf(viewer,"  using two-sided variant (for left eigenvectors)\n");
155:     }
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->ninil) {
193:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided left initial space: %D\n",PetscAbs(eps->ninil));
194:     }
195:     if (eps->nds) {
196:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided deflation space: %D\n",PetscAbs(eps->nds));
197:     }
198:   } else {
199:     if (eps->ops->view) {
200:       (*eps->ops->view)(eps,viewer);
201:     }
202:   }
203:   PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSBLZPACK,EPSELEMENTAL,EPSFEAST,EPSPRIMME,EPSSCALAPACK,EPSELPA,EPSEVSL,EPSTRLAN,"");
204:   if (!isexternal) {
205:     PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
206:     if (!eps->V) { EPSGetBV(eps,&eps->V); }
207:     BVView(eps->V,viewer);
208:     if (eps->rg) {
209:       RGIsTrivial(eps->rg,&istrivial);
210:       if (!istrivial) { RGView(eps->rg,viewer); }
211:     }
212:     if (eps->useds) {
213:       if (!eps->ds) { EPSGetDS(eps,&eps->ds); }
214:       DSView(eps->ds,viewer);
215:     }
216:     PetscViewerPopFormat(viewer);
217:   }
218:   if (!eps->st) { EPSGetST(eps,&eps->st); }
219:   STView(eps->st,viewer);
220:   return(0);
221: }

223: /*@C
224:    EPSViewFromOptions - View from options

226:    Collective on EPS

228:    Input Parameters:
229: +  eps  - the eigensolver context
230: .  obj  - optional object
231: -  name - command line option

233:    Level: intermediate

235: .seealso: EPSView(), EPSCreate()
236: @*/
237: PetscErrorCode EPSViewFromOptions(EPS eps,PetscObject obj,const char name[])
238: {

243:   PetscObjectViewFromOptions((PetscObject)eps,obj,name);
244:   return(0);
245: }

247: /*@C
248:    EPSConvergedReasonView - Displays the reason an EPS solve converged or diverged.

250:    Collective on eps

252:    Input Parameters:
253: +  eps - the eigensolver context
254: -  viewer - the viewer to display the reason

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

259:    Note:
260:    To change the format of the output call PetscViewerPushFormat(viewer,format) before
261:    this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
262:    display a reason if it fails. The latter can be set in the command line with
263:    -eps_converged_reason ::failed

265:    Level: intermediate

267: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber(), EPSConvergedReasonViewFromOptions()
268: @*/
269: PetscErrorCode EPSConvergedReasonView(EPS eps,PetscViewer viewer)
270: {
271:   PetscErrorCode    ierr;
272:   PetscBool         isAscii;
273:   PetscViewerFormat format;

276:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
277:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
278:   if (isAscii) {
279:     PetscViewerGetFormat(viewer,&format);
280:     PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
281:     if (eps->reason > 0 && format != PETSC_VIEWER_FAILED) {
282:       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);
283:     } else if (eps->reason <= 0) {
284:       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);
285:     }
286:     PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
287:   }
288:   return(0);
289: }

291: /*@
292:    EPSConvergedReasonViewFromOptions - Processes command line options to determine if/how
293:    the EPS converged reason is to be viewed.

295:    Collective on eps

297:    Input Parameter:
298: .  eps - the eigensolver context

300:    Level: developer

302: .seealso: EPSConvergedReasonView()
303: @*/
304: PetscErrorCode EPSConvergedReasonViewFromOptions(EPS eps)
305: {
306:   PetscErrorCode    ierr;
307:   PetscViewer       viewer;
308:   PetscBool         flg;
309:   static PetscBool  incall = PETSC_FALSE;
310:   PetscViewerFormat format;

313:   if (incall) return(0);
314:   incall = PETSC_TRUE;
315:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg);
316:   if (flg) {
317:     PetscViewerPushFormat(viewer,format);
318:     EPSConvergedReasonView(eps,viewer);
319:     PetscViewerPopFormat(viewer);
320:     PetscViewerDestroy(&viewer);
321:   }
322:   incall = PETSC_FALSE;
323:   return(0);
324: }

326: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
327: {
328:   PetscReal      error;
329:   PetscInt       i,j,k,nvals;

333:   nvals = (eps->which==EPS_ALL)? eps->nconv: eps->nev;
334:   if (eps->which!=EPS_ALL && eps->nconv<nvals) {
335:     PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
336:     return(0);
337:   }
338:   if (eps->which==EPS_ALL && !nvals) {
339:     PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
340:     return(0);
341:   }
342:   for (i=0;i<nvals;i++) {
343:     EPSComputeError(eps,i,etype,&error);
344:     if (error>=5.0*eps->tol) {
345:       PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
346:       return(0);
347:     }
348:   }
349:   if (eps->which==EPS_ALL) {
350:     PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
351:   } else {
352:     PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
353:   }
354:   for (i=0;i<=(nvals-1)/8;i++) {
355:     PetscViewerASCIIPrintf(viewer,"\n     ");
356:     for (j=0;j<PetscMin(8,nvals-8*i);j++) {
357:       k = eps->perm[8*i+j];
358:       SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]);
359:       if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
360:     }
361:   }
362:   PetscViewerASCIIPrintf(viewer,"\n\n");
363:   return(0);
364: }

366: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
367: {
369:   PetscReal      error,re,im;
370:   PetscScalar    kr,ki;
371:   PetscInt       i;
372:   char           ex[30],sep[]=" ---------------------- --------------------\n";

375:   if (!eps->nconv) return(0);
376:   switch (etype) {
377:     case EPS_ERROR_ABSOLUTE:
378:       PetscSNPrintf(ex,sizeof(ex),"   ||Ax-k%sx||",eps->isgeneralized?"B":"");
379:       break;
380:     case EPS_ERROR_RELATIVE:
381:       PetscSNPrintf(ex,sizeof(ex),"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":"");
382:       break;
383:     case EPS_ERROR_BACKWARD:
384:       PetscSNPrintf(ex,sizeof(ex),"    eta(x,k)");
385:       break;
386:   }
387:   PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep);
388:   for (i=0;i<eps->nconv;i++) {
389:     EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
390:     EPSComputeError(eps,i,etype,&error);
391: #if defined(PETSC_USE_COMPLEX)
392:     re = PetscRealPart(kr);
393:     im = PetscImaginaryPart(kr);
394: #else
395:     re = kr;
396:     im = ki;
397: #endif
398:     if (im!=0.0) {
399:       PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
400:     } else {
401:       PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
402:     }
403:   }
404:   PetscViewerASCIIPrintf(viewer,sep);
405:   return(0);
406: }

408: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
409: {
411:   PetscReal      error;
412:   PetscInt       i;
413:   const char     *name;

416:   PetscObjectGetName((PetscObject)eps,&name);
417:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
418:   for (i=0;i<eps->nconv;i++) {
419:     EPSComputeError(eps,i,etype,&error);
420:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
421:   }
422:   PetscViewerASCIIPrintf(viewer,"];\n");
423:   return(0);
424: }

426: /*@C
427:    EPSErrorView - Displays the errors associated with the computed solution
428:    (as well as the eigenvalues).

430:    Collective on eps

432:    Input Parameters:
433: +  eps    - the eigensolver context
434: .  etype  - error type
435: -  viewer - optional visualization context

437:    Options Database Key:
438: +  -eps_error_absolute - print absolute errors of each eigenpair
439: .  -eps_error_relative - print relative errors of each eigenpair
440: -  -eps_error_backward - print backward errors of each eigenpair

442:    Notes:
443:    By default, this function checks the error of all eigenpairs and prints
444:    the eigenvalues if all of them are below the requested tolerance.
445:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
446:    eigenvalues and corresponding errors is printed.

448:    Level: intermediate

450: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
451: @*/
452: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
453: {
454:   PetscBool         isascii;
455:   PetscViewerFormat format;
456:   PetscErrorCode    ierr;

460:   if (!viewer) {
461:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
462:   }
465:   EPSCheckSolved(eps,1);
466:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
467:   if (!isascii) return(0);

469:   PetscViewerGetFormat(viewer,&format);
470:   switch (format) {
471:     case PETSC_VIEWER_DEFAULT:
472:     case PETSC_VIEWER_ASCII_INFO:
473:       EPSErrorView_ASCII(eps,etype,viewer);
474:       break;
475:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
476:       EPSErrorView_DETAIL(eps,etype,viewer);
477:       break;
478:     case PETSC_VIEWER_ASCII_MATLAB:
479:       EPSErrorView_MATLAB(eps,etype,viewer);
480:       break;
481:     default:
482:       PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
483:   }
484:   return(0);
485: }

487: /*@
488:    EPSErrorViewFromOptions - Processes command line options to determine if/how
489:    the errors of the computed solution are to be viewed.

491:    Collective on eps

493:    Input Parameter:
494: .  eps - the eigensolver context

496:    Level: developer
497: @*/
498: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
499: {
500:   PetscErrorCode    ierr;
501:   PetscViewer       viewer;
502:   PetscBool         flg;
503:   static PetscBool  incall = PETSC_FALSE;
504:   PetscViewerFormat format;

507:   if (incall) return(0);
508:   incall = PETSC_TRUE;
509:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg);
510:   if (flg) {
511:     PetscViewerPushFormat(viewer,format);
512:     EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer);
513:     PetscViewerPopFormat(viewer);
514:     PetscViewerDestroy(&viewer);
515:   }
516:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg);
517:   if (flg) {
518:     PetscViewerPushFormat(viewer,format);
519:     EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer);
520:     PetscViewerPopFormat(viewer);
521:     PetscViewerDestroy(&viewer);
522:   }
523:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg);
524:   if (flg) {
525:     PetscViewerPushFormat(viewer,format);
526:     EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer);
527:     PetscViewerPopFormat(viewer);
528:     PetscViewerDestroy(&viewer);
529:   }
530:   incall = PETSC_FALSE;
531:   return(0);
532: }

534: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
535: {
537:   PetscDraw      draw;
538:   PetscDrawSP    drawsp;
539:   PetscReal      re,im;
540:   PetscInt       i,k;

543:   if (!eps->nconv) return(0);
544:   PetscViewerDrawGetDraw(viewer,0,&draw);
545:   PetscDrawSetTitle(draw,"Computed Eigenvalues");
546:   PetscDrawSPCreate(draw,1,&drawsp);
547:   for (i=0;i<eps->nconv;i++) {
548:     k = eps->perm[i];
549: #if defined(PETSC_USE_COMPLEX)
550:     re = PetscRealPart(eps->eigr[k]);
551:     im = PetscImaginaryPart(eps->eigr[k]);
552: #else
553:     re = eps->eigr[k];
554:     im = eps->eigi[k];
555: #endif
556:     PetscDrawSPAddPoint(drawsp,&re,&im);
557:   }
558:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
559:   PetscDrawSPSave(drawsp);
560:   PetscDrawSPDestroy(&drawsp);
561:   return(0);
562: }

564: static PetscErrorCode EPSValuesView_BINARY(EPS eps,PetscViewer viewer)
565: {
566: #if defined(PETSC_HAVE_COMPLEX)
567:   PetscInt       i,k;
568:   PetscComplex   *ev;
570: #endif

573: #if defined(PETSC_HAVE_COMPLEX)
574:   PetscMalloc1(eps->nconv,&ev);
575:   for (i=0;i<eps->nconv;i++) {
576:     k = eps->perm[i];
577: #if defined(PETSC_USE_COMPLEX)
578:     ev[i] = eps->eigr[k];
579: #else
580:     ev[i] = PetscCMPLX(eps->eigr[k],eps->eigi[k]);
581: #endif
582:   }
583:   PetscViewerBinaryWrite(viewer,ev,eps->nconv,PETSC_COMPLEX);
584:   PetscFree(ev);
585: #endif
586:   return(0);
587: }

589: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
590: {
591:   PetscInt       i,k;

595:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
596:   for (i=0;i<eps->nconv;i++) {
597:     k = eps->perm[i];
598:     PetscViewerASCIIPrintf(viewer,"   ");
599:     SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]);
600:     PetscViewerASCIIPrintf(viewer,"\n");
601:   }
602:   PetscViewerASCIIPrintf(viewer,"\n");
603:   return(0);
604: }

606: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
607: {
609:   PetscInt       i,k;
610:   PetscReal      re,im;
611:   const char     *name;

614:   PetscObjectGetName((PetscObject)eps,&name);
615:   PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
616:   for (i=0;i<eps->nconv;i++) {
617:     k = eps->perm[i];
618: #if defined(PETSC_USE_COMPLEX)
619:     re = PetscRealPart(eps->eigr[k]);
620:     im = PetscImaginaryPart(eps->eigr[k]);
621: #else
622:     re = eps->eigr[k];
623:     im = eps->eigi[k];
624: #endif
625:     if (im!=0.0) {
626:       PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
627:     } else {
628:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
629:     }
630:   }
631:   PetscViewerASCIIPrintf(viewer,"];\n");
632:   return(0);
633: }

635: /*@C
636:    EPSValuesView - Displays the computed eigenvalues in a viewer.

638:    Collective on eps

640:    Input Parameters:
641: +  eps    - the eigensolver context
642: -  viewer - the viewer

644:    Options Database Key:
645: .  -eps_view_values - print computed eigenvalues

647:    Level: intermediate

649: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
650: @*/
651: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
652: {
653:   PetscBool         isascii,isdraw,isbinary;
654:   PetscViewerFormat format;
655:   PetscErrorCode    ierr;

659:   if (!viewer) {
660:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
661:   }
664:   EPSCheckSolved(eps,1);
665:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
666:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
667:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
668:   if (isdraw) {
669:     EPSValuesView_DRAW(eps,viewer);
670:   } else if (isbinary) {
671:     EPSValuesView_BINARY(eps,viewer);
672:   } else if (isascii) {
673:     PetscViewerGetFormat(viewer,&format);
674:     switch (format) {
675:       case PETSC_VIEWER_DEFAULT:
676:       case PETSC_VIEWER_ASCII_INFO:
677:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
678:         EPSValuesView_ASCII(eps,viewer);
679:         break;
680:       case PETSC_VIEWER_ASCII_MATLAB:
681:         EPSValuesView_MATLAB(eps,viewer);
682:         break;
683:       default:
684:         PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
685:     }
686:   }
687:   return(0);
688: }

690: /*@
691:    EPSValuesViewFromOptions - Processes command line options to determine if/how
692:    the computed eigenvalues are to be viewed.

694:    Collective on eps

696:    Input Parameters:
697: .  eps - the eigensolver context

699:    Level: developer
700: @*/
701: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
702: {
703:   PetscErrorCode    ierr;
704:   PetscViewer       viewer;
705:   PetscBool         flg;
706:   static PetscBool  incall = PETSC_FALSE;
707:   PetscViewerFormat format;

710:   if (incall) return(0);
711:   incall = PETSC_TRUE;
712:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg);
713:   if (flg) {
714:     PetscViewerPushFormat(viewer,format);
715:     EPSValuesView(eps,viewer);
716:     PetscViewerPopFormat(viewer);
717:     PetscViewerDestroy(&viewer);
718:   }
719:   incall = PETSC_FALSE;
720:   return(0);
721: }

723: /*@C
724:    EPSVectorsView - Outputs computed eigenvectors to a viewer.

726:    Collective on eps

728:    Parameter:
729: +  eps    - the eigensolver context
730: -  viewer - the viewer

732:    Options Database Keys:
733: .  -eps_view_vectors - output eigenvectors.

735:    Notes:
736:    If PETSc was configured with real scalars, complex conjugate eigenvectors
737:    will be viewed as two separate real vectors, one containing the real part
738:    and another one containing the imaginary part.

740:    If left eigenvectors were computed with a two-sided eigensolver, the right
741:    and left eigenvectors are interleaved, that is, the vectors are output in
742:    the following order: X0, Y0, X1, Y1, X2, Y2, ...

744:    Level: intermediate

746: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
747: @*/
748: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
749: {
751:   PetscInt       i,k;
752:   Vec            x;
753:   char           vname[30];
754:   const char     *ename;

758:   if (!viewer) {
759:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
760:   }
763:   EPSCheckSolved(eps,1);
764:   if (eps->nconv) {
765:     PetscObjectGetName((PetscObject)eps,&ename);
766:     EPSComputeVectors(eps);
767:     for (i=0;i<eps->nconv;i++) {
768:       k = eps->perm[i];
769:       PetscSNPrintf(vname,sizeof(vname),"X%d_%s",(int)i,ename);
770:       BVGetColumn(eps->V,k,&x);
771:       PetscObjectSetName((PetscObject)x,vname);
772:       VecView(x,viewer);
773:       BVRestoreColumn(eps->V,k,&x);
774:       if (eps->twosided) {
775:         PetscSNPrintf(vname,sizeof(vname),"Y%d_%s",(int)i,ename);
776:         BVGetColumn(eps->W,k,&x);
777:         PetscObjectSetName((PetscObject)x,vname);
778:         VecView(x,viewer);
779:         BVRestoreColumn(eps->W,k,&x);
780:       }
781:     }
782:   }
783:   return(0);
784: }

786: /*@
787:    EPSVectorsViewFromOptions - Processes command line options to determine if/how
788:    the computed eigenvectors are to be viewed.

790:    Collective on eps

792:    Input Parameter:
793: .  eps - the eigensolver context

795:    Level: developer
796: @*/
797: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
798: {
799:   PetscErrorCode    ierr;
800:   PetscViewer       viewer;
801:   PetscBool         flg = PETSC_FALSE;
802:   static PetscBool  incall = PETSC_FALSE;
803:   PetscViewerFormat format;

806:   if (incall) return(0);
807:   incall = PETSC_TRUE;
808:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg);
809:   if (flg) {
810:     PetscViewerPushFormat(viewer,format);
811:     EPSVectorsView(eps,viewer);
812:     PetscViewerPopFormat(viewer);
813:     PetscViewerDestroy(&viewer);
814:   }
815:   incall = PETSC_FALSE;
816:   return(0);
817: }