Actual source code: epsview.c

slepc-3.16.1 2021-11-17
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 <slepc/private/bvimpl.h>
 16: #include <petscdraw.h>

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

 21:    Collective on eps

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

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

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

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

 41:    Level: beginner

 43: .seealso: STView()
 44: @*/
 45: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
 46: {
 48:   const char     *type=NULL,*extr=NULL,*bal=NULL;
 49:   char           str[50];
 50:   PetscBool      isascii,isexternal,istrivial;

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

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

227:    Collective on EPS

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

234:    Level: intermediate

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

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

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

251:    Collective on eps

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

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

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

266:    Level: intermediate

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

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

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

296:    Collective on eps

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

301:    Level: developer

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

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

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

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

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

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

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

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

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

431:    Collective on eps

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

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

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

449:    Level: intermediate

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

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

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

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

492:    Collective on eps

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

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

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

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

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

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

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

590: #if defined(PETSC_HAVE_HDF5)
591: static PetscErrorCode EPSValuesView_HDF5(EPS eps,PetscViewer viewer)
592: {
594:   PetscInt       i,k,n,N;
595:   PetscMPIInt    rank;
596:   Vec            v;
597:   char           vname[30];
598:   const char     *ename;

601:   MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank);
602:   N = eps->nconv;
603:   n = rank? 0: N;
604:   /* create a vector containing the eigenvalues */
605:   VecCreateMPI(PetscObjectComm((PetscObject)eps),n,N,&v);
606:   PetscObjectGetName((PetscObject)eps,&ename);
607:   PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename);
608:   PetscObjectSetName((PetscObject)v,vname);
609:   if (!rank) {
610:     for (i=0;i<eps->nconv;i++) {
611:       k = eps->perm[i];
612:       VecSetValue(v,i,eps->eigr[k],INSERT_VALUES);
613:     }
614:   }
615:   VecAssemblyBegin(v);
616:   VecAssemblyEnd(v);
617:   VecView(v,viewer);
618: #if !defined(PETSC_USE_COMPLEX)
619:   /* in real scalars write the imaginary part as a separate vector */
620:   PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename);
621:   PetscObjectSetName((PetscObject)v,vname);
622:   if (!rank) {
623:     for (i=0;i<eps->nconv;i++) {
624:       k = eps->perm[i];
625:       VecSetValue(v,i,eps->eigi[k],INSERT_VALUES);
626:     }
627:   }
628:   VecAssemblyBegin(v);
629:   VecAssemblyEnd(v);
630:   VecView(v,viewer);
631: #endif
632:   VecDestroy(&v);
633:   return(0);
634: }
635: #endif

637: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
638: {
639:   PetscInt       i,k;

643:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
644:   for (i=0;i<eps->nconv;i++) {
645:     k = eps->perm[i];
646:     PetscViewerASCIIPrintf(viewer,"   ");
647:     SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]);
648:     PetscViewerASCIIPrintf(viewer,"\n");
649:   }
650:   PetscViewerASCIIPrintf(viewer,"\n");
651:   return(0);
652: }

654: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
655: {
657:   PetscInt       i,k;
658:   PetscReal      re,im;
659:   const char     *name;

662:   PetscObjectGetName((PetscObject)eps,&name);
663:   PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
664:   for (i=0;i<eps->nconv;i++) {
665:     k = eps->perm[i];
666: #if defined(PETSC_USE_COMPLEX)
667:     re = PetscRealPart(eps->eigr[k]);
668:     im = PetscImaginaryPart(eps->eigr[k]);
669: #else
670:     re = eps->eigr[k];
671:     im = eps->eigi[k];
672: #endif
673:     if (im!=0.0) {
674:       PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
675:     } else {
676:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
677:     }
678:   }
679:   PetscViewerASCIIPrintf(viewer,"];\n");
680:   return(0);
681: }

683: /*@C
684:    EPSValuesView - Displays the computed eigenvalues in a viewer.

686:    Collective on eps

688:    Input Parameters:
689: +  eps    - the eigensolver context
690: -  viewer - the viewer

692:    Options Database Key:
693: .  -eps_view_values - print computed eigenvalues

695:    Level: intermediate

697: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
698: @*/
699: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
700: {
701:   PetscBool         isascii,isdraw,isbinary;
702:   PetscViewerFormat format;
703:   PetscErrorCode    ierr;
704: #if defined(PETSC_HAVE_HDF5)
705:   PetscBool         ishdf5;
706: #endif

710:   if (!viewer) {
711:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
712:   }
715:   EPSCheckSolved(eps,1);
716:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
717:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
718: #if defined(PETSC_HAVE_HDF5)
719:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
720: #endif
721:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
722:   if (isdraw) {
723:     EPSValuesView_DRAW(eps,viewer);
724:   } else if (isbinary) {
725:     EPSValuesView_BINARY(eps,viewer);
726: #if defined(PETSC_HAVE_HDF5)
727:   } else if (ishdf5) {
728:     EPSValuesView_HDF5(eps,viewer);
729: #endif
730:   } else if (isascii) {
731:     PetscViewerGetFormat(viewer,&format);
732:     switch (format) {
733:       case PETSC_VIEWER_DEFAULT:
734:       case PETSC_VIEWER_ASCII_INFO:
735:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
736:         EPSValuesView_ASCII(eps,viewer);
737:         break;
738:       case PETSC_VIEWER_ASCII_MATLAB:
739:         EPSValuesView_MATLAB(eps,viewer);
740:         break;
741:       default:
742:         PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
743:     }
744:   }
745:   return(0);
746: }

748: /*@
749:    EPSValuesViewFromOptions - Processes command line options to determine if/how
750:    the computed eigenvalues are to be viewed.

752:    Collective on eps

754:    Input Parameters:
755: .  eps - the eigensolver context

757:    Level: developer
758: @*/
759: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
760: {
761:   PetscErrorCode    ierr;
762:   PetscViewer       viewer;
763:   PetscBool         flg;
764:   static PetscBool  incall = PETSC_FALSE;
765:   PetscViewerFormat format;

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

781: /*@C
782:    EPSVectorsView - Outputs computed eigenvectors to a viewer.

784:    Collective on eps

786:    Parameter:
787: +  eps    - the eigensolver context
788: -  viewer - the viewer

790:    Options Database Keys:
791: .  -eps_view_vectors - output eigenvectors.

793:    Notes:
794:    If PETSc was configured with real scalars, complex conjugate eigenvectors
795:    will be viewed as two separate real vectors, one containing the real part
796:    and another one containing the imaginary part.

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

802:    Level: intermediate

804: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
805: @*/
806: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
807: {
809:   PetscInt       i,k;
810:   Vec            xr,xi=NULL;

814:   if (!viewer) {
815:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
816:   }
819:   EPSCheckSolved(eps,1);
820:   if (eps->nconv) {
821:     EPSComputeVectors(eps);
822:     BVCreateVec(eps->V,&xr);
823: #if !defined(PETSC_USE_COMPLEX)
824:     BVCreateVec(eps->V,&xi);
825: #endif
826:     for (i=0;i<eps->nconv;i++) {
827:       k = eps->perm[i];
828:       BV_GetEigenvector(eps->V,k,eps->eigi[k],xr,xi);
829:       SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)eps);
830:       if (eps->twosided) {
831:         BV_GetEigenvector(eps->W,k,eps->eigi[k],xr,xi);
832:         SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)eps);
833:       }
834:     }
835:     VecDestroy(&xr);
836: #if !defined(PETSC_USE_COMPLEX)
837:     VecDestroy(&xi);
838: #endif
839:   }
840:   return(0);
841: }

843: /*@
844:    EPSVectorsViewFromOptions - Processes command line options to determine if/how
845:    the computed eigenvectors are to be viewed.

847:    Collective on eps

849:    Input Parameter:
850: .  eps - the eigensolver context

852:    Level: developer
853: @*/
854: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
855: {
856:   PetscErrorCode    ierr;
857:   PetscViewer       viewer;
858:   PetscBool         flg = PETSC_FALSE;
859:   static PetscBool  incall = PETSC_FALSE;
860:   PetscViewerFormat format;

863:   if (incall) return(0);
864:   incall = PETSC_TRUE;
865:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg);
866:   if (flg) {
867:     PetscViewerPushFormat(viewer,format);
868:     EPSVectorsView(eps,viewer);
869:     PetscViewerPopFormat(viewer);
870:     PetscViewerDestroy(&viewer);
871:   }
872:   incall = PETSC_FALSE;
873:   return(0);
874: }