Actual source code: epsview.c

slepc-3.17.2 2022-08-09
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, 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: {
 47:   const char     *type=NULL,*extr=NULL,*bal=NULL;
 48:   char           str[50];
 49:   PetscBool      isascii,isexternal,istrivial;

 52:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);

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

194: /*@C
195:    EPSViewFromOptions - View from options

197:    Collective on EPS

199:    Input Parameters:
200: +  eps  - the eigensolver context
201: .  obj  - optional object
202: -  name - command line option

204:    Level: intermediate

206: .seealso: EPSView(), EPSCreate()
207: @*/
208: PetscErrorCode EPSViewFromOptions(EPS eps,PetscObject obj,const char name[])
209: {
211:   PetscObjectViewFromOptions((PetscObject)eps,obj,name);
212:   PetscFunctionReturn(0);
213: }

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

218:    Collective on eps

220:    Input Parameters:
221: +  eps - the eigensolver context
222: -  viewer - the viewer to display the reason

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

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

233:    Level: intermediate

235: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber(), EPSConvergedReasonViewFromOptions()
236: @*/
237: PetscErrorCode EPSConvergedReasonView(EPS eps,PetscViewer viewer)
238: {
239:   PetscBool         isAscii;
240:   PetscViewerFormat format;

242:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
243:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
244:   if (isAscii) {
245:     PetscViewerGetFormat(viewer,&format);
246:     PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
247:     if (eps->reason > 0 && format != PETSC_VIEWER_FAILED) PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve converged (%" PetscInt_FMT " eigenpair%s) due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",eps->nconv,(eps->nconv>1)?"s":"",EPSConvergedReasons[eps->reason],eps->its);
248:     else if (eps->reason <= 0) PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve did not converge due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",EPSConvergedReasons[eps->reason],eps->its);
249:     PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
250:   }
251:   PetscFunctionReturn(0);
252: }

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

258:    Collective on eps

260:    Input Parameter:
261: .  eps - the eigensolver context

263:    Level: developer

265: .seealso: EPSConvergedReasonView()
266: @*/
267: PetscErrorCode EPSConvergedReasonViewFromOptions(EPS eps)
268: {
269:   PetscViewer       viewer;
270:   PetscBool         flg;
271:   static PetscBool  incall = PETSC_FALSE;
272:   PetscViewerFormat format;

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

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

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

322: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
323: {
324:   PetscReal      error,re,im;
325:   PetscScalar    kr,ki;
326:   PetscInt       i;
327:   char           ex[30],sep[]=" ---------------------- --------------------\n";

329:   if (!eps->nconv) PetscFunctionReturn(0);
330:   switch (etype) {
331:     case EPS_ERROR_ABSOLUTE:
332:       PetscSNPrintf(ex,sizeof(ex),"   ||Ax-k%sx||",eps->isgeneralized?"B":"");
333:       break;
334:     case EPS_ERROR_RELATIVE:
335:       PetscSNPrintf(ex,sizeof(ex),"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":"");
336:       break;
337:     case EPS_ERROR_BACKWARD:
338:       PetscSNPrintf(ex,sizeof(ex),"    eta(x,k)");
339:       break;
340:   }
341:   PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep);
342:   for (i=0;i<eps->nconv;i++) {
343:     EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
344:     EPSComputeError(eps,i,etype,&error);
345: #if defined(PETSC_USE_COMPLEX)
346:     re = PetscRealPart(kr);
347:     im = PetscImaginaryPart(kr);
348: #else
349:     re = kr;
350:     im = ki;
351: #endif
352:     if (im!=0.0) PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
353:     else PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
354:   }
355:   PetscViewerASCIIPrintf(viewer,"%s",sep);
356:   PetscFunctionReturn(0);
357: }

359: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
360: {
361:   PetscReal      error;
362:   PetscInt       i;
363:   const char     *name;

365:   PetscObjectGetName((PetscObject)eps,&name);
366:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
367:   for (i=0;i<eps->nconv;i++) {
368:     EPSComputeError(eps,i,etype,&error);
369:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
370:   }
371:   PetscViewerASCIIPrintf(viewer,"];\n");
372:   PetscFunctionReturn(0);
373: }

375: /*@C
376:    EPSErrorView - Displays the errors associated with the computed solution
377:    (as well as the eigenvalues).

379:    Collective on eps

381:    Input Parameters:
382: +  eps    - the eigensolver context
383: .  etype  - error type
384: -  viewer - optional visualization context

386:    Options Database Keys:
387: +  -eps_error_absolute - print absolute errors of each eigenpair
388: .  -eps_error_relative - print relative errors of each eigenpair
389: -  -eps_error_backward - print backward errors of each eigenpair

391:    Notes:
392:    By default, this function checks the error of all eigenpairs and prints
393:    the eigenvalues if all of them are below the requested tolerance.
394:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
395:    eigenvalues and corresponding errors is printed.

397:    Level: intermediate

399: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
400: @*/
401: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
402: {
403:   PetscBool         isascii;
404:   PetscViewerFormat format;

407:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
410:   EPSCheckSolved(eps,1);
411:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
412:   if (!isascii) PetscFunctionReturn(0);

414:   PetscViewerGetFormat(viewer,&format);
415:   switch (format) {
416:     case PETSC_VIEWER_DEFAULT:
417:     case PETSC_VIEWER_ASCII_INFO:
418:       EPSErrorView_ASCII(eps,etype,viewer);
419:       break;
420:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
421:       EPSErrorView_DETAIL(eps,etype,viewer);
422:       break;
423:     case PETSC_VIEWER_ASCII_MATLAB:
424:       EPSErrorView_MATLAB(eps,etype,viewer);
425:       break;
426:     default:
427:       PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
428:   }
429:   PetscFunctionReturn(0);
430: }

432: /*@
433:    EPSErrorViewFromOptions - Processes command line options to determine if/how
434:    the errors of the computed solution are to be viewed.

436:    Collective on eps

438:    Input Parameter:
439: .  eps - the eigensolver context

441:    Level: developer

443: .seealso: EPSErrorView()
444: @*/
445: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
446: {
447:   PetscViewer       viewer;
448:   PetscBool         flg;
449:   static PetscBool  incall = PETSC_FALSE;
450:   PetscViewerFormat format;

452:   if (incall) PetscFunctionReturn(0);
453:   incall = PETSC_TRUE;
454:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg);
455:   if (flg) {
456:     PetscViewerPushFormat(viewer,format);
457:     EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer);
458:     PetscViewerPopFormat(viewer);
459:     PetscViewerDestroy(&viewer);
460:   }
461:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg);
462:   if (flg) {
463:     PetscViewerPushFormat(viewer,format);
464:     EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer);
465:     PetscViewerPopFormat(viewer);
466:     PetscViewerDestroy(&viewer);
467:   }
468:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg);
469:   if (flg) {
470:     PetscViewerPushFormat(viewer,format);
471:     EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer);
472:     PetscViewerPopFormat(viewer);
473:     PetscViewerDestroy(&viewer);
474:   }
475:   incall = PETSC_FALSE;
476:   PetscFunctionReturn(0);
477: }

479: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
480: {
481:   PetscDraw      draw;
482:   PetscDrawSP    drawsp;
483:   PetscReal      re,im;
484:   PetscInt       i,k;

486:   if (!eps->nconv) PetscFunctionReturn(0);
487:   PetscViewerDrawGetDraw(viewer,0,&draw);
488:   PetscDrawSetTitle(draw,"Computed Eigenvalues");
489:   PetscDrawSPCreate(draw,1,&drawsp);
490:   for (i=0;i<eps->nconv;i++) {
491:     k = eps->perm[i];
492: #if defined(PETSC_USE_COMPLEX)
493:     re = PetscRealPart(eps->eigr[k]);
494:     im = PetscImaginaryPart(eps->eigr[k]);
495: #else
496:     re = eps->eigr[k];
497:     im = eps->eigi[k];
498: #endif
499:     PetscDrawSPAddPoint(drawsp,&re,&im);
500:   }
501:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
502:   PetscDrawSPSave(drawsp);
503:   PetscDrawSPDestroy(&drawsp);
504:   PetscFunctionReturn(0);
505: }

507: static PetscErrorCode EPSValuesView_BINARY(EPS eps,PetscViewer viewer)
508: {
509: #if defined(PETSC_HAVE_COMPLEX)
510:   PetscInt       i,k;
511:   PetscComplex   *ev;
512: #endif

514: #if defined(PETSC_HAVE_COMPLEX)
515:   PetscMalloc1(eps->nconv,&ev);
516:   for (i=0;i<eps->nconv;i++) {
517:     k = eps->perm[i];
518: #if defined(PETSC_USE_COMPLEX)
519:     ev[i] = eps->eigr[k];
520: #else
521:     ev[i] = PetscCMPLX(eps->eigr[k],eps->eigi[k]);
522: #endif
523:   }
524:   PetscViewerBinaryWrite(viewer,ev,eps->nconv,PETSC_COMPLEX);
525:   PetscFree(ev);
526: #endif
527:   PetscFunctionReturn(0);
528: }

530: #if defined(PETSC_HAVE_HDF5)
531: static PetscErrorCode EPSValuesView_HDF5(EPS eps,PetscViewer viewer)
532: {
533:   PetscInt       i,k,n,N;
534:   PetscMPIInt    rank;
535:   Vec            v;
536:   char           vname[30];
537:   const char     *ename;

539:   MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank);
540:   N = eps->nconv;
541:   n = rank? 0: N;
542:   /* create a vector containing the eigenvalues */
543:   VecCreateMPI(PetscObjectComm((PetscObject)eps),n,N,&v);
544:   PetscObjectGetName((PetscObject)eps,&ename);
545:   PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename);
546:   PetscObjectSetName((PetscObject)v,vname);
547:   if (!rank) {
548:     for (i=0;i<eps->nconv;i++) {
549:       k = eps->perm[i];
550:       VecSetValue(v,i,eps->eigr[k],INSERT_VALUES);
551:     }
552:   }
553:   VecAssemblyBegin(v);
554:   VecAssemblyEnd(v);
555:   VecView(v,viewer);
556: #if !defined(PETSC_USE_COMPLEX)
557:   /* in real scalars write the imaginary part as a separate vector */
558:   PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename);
559:   PetscObjectSetName((PetscObject)v,vname);
560:   if (!rank) {
561:     for (i=0;i<eps->nconv;i++) {
562:       k = eps->perm[i];
563:       VecSetValue(v,i,eps->eigi[k],INSERT_VALUES);
564:     }
565:   }
566:   VecAssemblyBegin(v);
567:   VecAssemblyEnd(v);
568:   VecView(v,viewer);
569: #endif
570:   VecDestroy(&v);
571:   PetscFunctionReturn(0);
572: }
573: #endif

575: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
576: {
577:   PetscInt       i,k;

579:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
580:   for (i=0;i<eps->nconv;i++) {
581:     k = eps->perm[i];
582:     PetscViewerASCIIPrintf(viewer,"   ");
583:     SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]);
584:     PetscViewerASCIIPrintf(viewer,"\n");
585:   }
586:   PetscViewerASCIIPrintf(viewer,"\n");
587:   PetscFunctionReturn(0);
588: }

590: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
591: {
592:   PetscInt       i,k;
593:   PetscReal      re,im;
594:   const char     *name;

596:   PetscObjectGetName((PetscObject)eps,&name);
597:   PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
598:   for (i=0;i<eps->nconv;i++) {
599:     k = eps->perm[i];
600: #if defined(PETSC_USE_COMPLEX)
601:     re = PetscRealPart(eps->eigr[k]);
602:     im = PetscImaginaryPart(eps->eigr[k]);
603: #else
604:     re = eps->eigr[k];
605:     im = eps->eigi[k];
606: #endif
607:     if (im!=0.0) PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
608:     else PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
609:   }
610:   PetscViewerASCIIPrintf(viewer,"];\n");
611:   PetscFunctionReturn(0);
612: }

614: /*@C
615:    EPSValuesView - Displays the computed eigenvalues in a viewer.

617:    Collective on eps

619:    Input Parameters:
620: +  eps    - the eigensolver context
621: -  viewer - the viewer

623:    Options Database Key:
624: .  -eps_view_values - print computed eigenvalues

626:    Level: intermediate

628: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
629: @*/
630: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
631: {
632:   PetscBool         isascii,isdraw,isbinary;
633:   PetscViewerFormat format;
634: #if defined(PETSC_HAVE_HDF5)
635:   PetscBool         ishdf5;
636: #endif

639:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
642:   EPSCheckSolved(eps,1);
643:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
644:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
645: #if defined(PETSC_HAVE_HDF5)
646:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
647: #endif
648:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
649:   if (isdraw) EPSValuesView_DRAW(eps,viewer);
650:   else if (isbinary) EPSValuesView_BINARY(eps,viewer);
651: #if defined(PETSC_HAVE_HDF5)
652:   else if (ishdf5) EPSValuesView_HDF5(eps,viewer);
653: #endif
654:   else if (isascii) {
655:     PetscViewerGetFormat(viewer,&format);
656:     switch (format) {
657:       case PETSC_VIEWER_DEFAULT:
658:       case PETSC_VIEWER_ASCII_INFO:
659:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
660:         EPSValuesView_ASCII(eps,viewer);
661:         break;
662:       case PETSC_VIEWER_ASCII_MATLAB:
663:         EPSValuesView_MATLAB(eps,viewer);
664:         break;
665:       default:
666:         PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
667:     }
668:   }
669:   PetscFunctionReturn(0);
670: }

672: /*@
673:    EPSValuesViewFromOptions - Processes command line options to determine if/how
674:    the computed eigenvalues are to be viewed.

676:    Collective on eps

678:    Input Parameters:
679: .  eps - the eigensolver context

681:    Level: developer

683: .seealso: EPSValuesView()
684: @*/
685: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
686: {
687:   PetscViewer       viewer;
688:   PetscBool         flg;
689:   static PetscBool  incall = PETSC_FALSE;
690:   PetscViewerFormat format;

692:   if (incall) PetscFunctionReturn(0);
693:   incall = PETSC_TRUE;
694:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg);
695:   if (flg) {
696:     PetscViewerPushFormat(viewer,format);
697:     EPSValuesView(eps,viewer);
698:     PetscViewerPopFormat(viewer);
699:     PetscViewerDestroy(&viewer);
700:   }
701:   incall = PETSC_FALSE;
702:   PetscFunctionReturn(0);
703: }

705: /*@C
706:    EPSVectorsView - Outputs computed eigenvectors to a viewer.

708:    Collective on eps

710:    Input Parameters:
711: +  eps    - the eigensolver context
712: -  viewer - the viewer

714:    Options Database Key:
715: .  -eps_view_vectors - output eigenvectors.

717:    Notes:
718:    If PETSc was configured with real scalars, complex conjugate eigenvectors
719:    will be viewed as two separate real vectors, one containing the real part
720:    and another one containing the imaginary part.

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

726:    Level: intermediate

728: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
729: @*/
730: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
731: {
732:   PetscInt       i,k;
733:   Vec            xr,xi=NULL;

736:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
739:   EPSCheckSolved(eps,1);
740:   if (eps->nconv) {
741:     EPSComputeVectors(eps);
742:     BVCreateVec(eps->V,&xr);
743: #if !defined(PETSC_USE_COMPLEX)
744:     BVCreateVec(eps->V,&xi);
745: #endif
746:     for (i=0;i<eps->nconv;i++) {
747:       k = eps->perm[i];
748:       BV_GetEigenvector(eps->V,k,eps->eigi[k],xr,xi);
749:       SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)eps);
750:       if (eps->twosided) {
751:         BV_GetEigenvector(eps->W,k,eps->eigi[k],xr,xi);
752:         SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)eps);
753:       }
754:     }
755:     VecDestroy(&xr);
756: #if !defined(PETSC_USE_COMPLEX)
757:     VecDestroy(&xi);
758: #endif
759:   }
760:   PetscFunctionReturn(0);
761: }

763: /*@
764:    EPSVectorsViewFromOptions - Processes command line options to determine if/how
765:    the computed eigenvectors are to be viewed.

767:    Collective on eps

769:    Input Parameter:
770: .  eps - the eigensolver context

772:    Level: developer

774: .seealso: EPSVectorsView()
775: @*/
776: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
777: {
778:   PetscViewer       viewer;
779:   PetscBool         flg = PETSC_FALSE;
780:   static PetscBool  incall = PETSC_FALSE;
781:   PetscViewerFormat format;

783:   if (incall) PetscFunctionReturn(0);
784:   incall = PETSC_TRUE;
785:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg);
786:   if (flg) {
787:     PetscViewerPushFormat(viewer,format);
788:     EPSVectorsView(eps,viewer);
789:     PetscViewerPopFormat(viewer);
790:     PetscViewerDestroy(&viewer);
791:   }
792:   incall = PETSC_FALSE;
793:   PetscFunctionReturn(0);
794: }