Actual source code: epsview.c

slepc-3.6.1 2015-09-03
Report Typos and Errors
  1: /*
  2:    The EPS routines related to various viewers.

  4:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  5:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  6:    Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain

  8:    This file is part of SLEPc.

 10:    SLEPc is free software: you can redistribute it and/or modify it under  the
 11:    terms of version 3 of the GNU Lesser General Public License as published by
 12:    the Free Software Foundation.

 14:    SLEPc  is  distributed in the hope that it will be useful, but WITHOUT  ANY
 15:    WARRANTY;  without even the implied warranty of MERCHANTABILITY or  FITNESS
 16:    FOR  A  PARTICULAR PURPOSE. See the GNU Lesser General Public  License  for
 17:    more details.

 19:    You  should have received a copy of the GNU Lesser General  Public  License
 20:    along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
 21:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
 22: */

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

 29: /*@C
 30:    EPSView - Prints the EPS data structure.

 32:    Collective on EPS

 34:    Input Parameters:
 35: +  eps - the eigenproblem solver context
 36: -  viewer - optional visualization context

 38:    Options Database Key:
 39: .  -eps_view -  Calls EPSView() at end of EPSSolve()

 41:    Note:
 42:    The available visualization contexts include
 43: +     PETSC_VIEWER_STDOUT_SELF - standard output (default)
 44: -     PETSC_VIEWER_STDOUT_WORLD - synchronized standard
 45:          output where only the first processor opens
 46:          the file.  All other processors send their
 47:          data to the first processor to print.

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

 52:    Level: beginner

 54: .seealso: STView(), PetscViewerASCIIOpen()
 55: @*/
 56: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
 57: {
 59:   const char     *type,*extr,*bal;
 60:   char           str[50];
 61:   PetscBool      isascii,ispower,isexternal,istrivial;

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

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

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

231:    Collective on EPS

233:    Parameter:
234: +  eps - the eigensolver context
235: -  viewer - the viewer to display the reason

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

240:    Level: intermediate

242: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber()
243: @*/
244: PetscErrorCode EPSReasonView(EPS eps,PetscViewer viewer)
245: {
247:   PetscBool      isAscii;

250:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
251:   if (isAscii) {
252:     PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
253:     if (eps->reason > 0) {
254:       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);
255:     } else {
256:       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);
257:     }
258:     PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
259:   }
260:   return(0);
261: }

265: /*@
266:    EPSReasonViewFromOptions - Processes command line options to determine if/how
267:    the EPS converged reason is to be viewed. 

269:    Collective on EPS

271:    Input Parameters:
272: .  eps - the eigensolver context

274:    Level: developer
275: @*/
276: PetscErrorCode EPSReasonViewFromOptions(EPS eps)
277: {
278:   PetscErrorCode    ierr;
279:   PetscViewer       viewer;
280:   PetscBool         flg;
281:   static PetscBool  incall = PETSC_FALSE;
282:   PetscViewerFormat format;

285:   if (incall) return(0);
286:   incall = PETSC_TRUE;
287:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg);
288:   if (flg) {
289:     PetscViewerPushFormat(viewer,format);
290:     EPSReasonView(eps,viewer);
291:     PetscViewerPopFormat(viewer);
292:     PetscViewerDestroy(&viewer);
293:   }
294:   incall = PETSC_FALSE;
295:   return(0);
296: }

300: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
301: {
302:   PetscBool      errok;
303:   PetscReal      error,re,im;
304:   PetscScalar    kr,ki;
305:   PetscInt       i,j;

309:   if (eps->nconv<eps->nev) {
310:     PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",eps->nev);
311:     return(0);
312:   }
313:   errok = PETSC_TRUE;
314:   for (i=0;i<eps->nev;i++) {
315:     EPSComputeError(eps,i,etype,&error);
316:     errok = (errok && error<5.0*eps->tol)? PETSC_TRUE: PETSC_FALSE;
317:   }
318:   if (!errok) {
319:     PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",eps->nev);
320:     return(0);
321:   }
322:   PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
323:   for (i=0;i<=(eps->nev-1)/8;i++) {
324:     PetscViewerASCIIPrintf(viewer,"\n     ");
325:     for (j=0;j<PetscMin(8,eps->nev-8*i);j++) {
326:       EPSGetEigenpair(eps,8*i+j,&kr,&ki,NULL,NULL);
327: #if defined(PETSC_USE_COMPLEX)
328:       re = PetscRealPart(kr);
329:       im = PetscImaginaryPart(kr);
330: #else
331:       re = kr;
332:       im = ki;
333: #endif
334:       if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
335:       if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
336:       if (im!=0.0) {
337:         PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
338:       } else {
339:         PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
340:       }
341:       if (8*i+j+1<eps->nev) { PetscViewerASCIIPrintf(viewer,", "); }
342:     }
343:   }
344:   PetscViewerASCIIPrintf(viewer,"\n\n");
345:   return(0);
346: }

350: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
351: {
353:   PetscReal      error,re,im;
354:   PetscScalar    kr,ki;
355:   PetscInt       i;
356: #define EXLEN 30
357:   char           ex[EXLEN],sep[]=" ---------------------- --------------------\n";

360:   if (!eps->nconv) return(0);
361:   switch (etype) {
362:     case EPS_ERROR_ABSOLUTE:
363:       PetscSNPrintf(ex,EXLEN,"   ||Ax-k%sx||",eps->isgeneralized?"B":"");
364:       break;
365:     case EPS_ERROR_RELATIVE:
366:       PetscSNPrintf(ex,EXLEN,"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":"");
367:       break;
368:     case EPS_ERROR_BACKWARD:
369:       PetscSNPrintf(ex,EXLEN,"    eta(x,k)");
370:       break;
371:   }
372:   PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep);
373:   for (i=0;i<eps->nconv;i++) {
374:     EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
375:     EPSComputeError(eps,i,etype,&error);
376: #if defined(PETSC_USE_COMPLEX)
377:     re = PetscRealPart(kr);
378:     im = PetscImaginaryPart(kr);
379: #else
380:     re = kr;
381:     im = ki;
382: #endif
383:     if (im!=0.0) {
384:       PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
385:     } else {
386:       PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
387:     }
388:   }
389:   PetscViewerASCIIPrintf(viewer,sep);
390:   return(0);
391: }

395: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
396: {
398:   PetscReal      error;
399:   PetscInt       i;
400:   const char     *name;

403:   PetscObjectGetName((PetscObject)eps,&name);
404:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
405:   for (i=0;i<eps->nconv;i++) {
406:     EPSComputeError(eps,i,etype,&error);
407:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",error);
408:   }
409:   PetscViewerASCIIPrintf(viewer,"];\n");
410:   return(0);
411: }

415: /*@C
416:    EPSErrorView - Displays the errors associated with the computed solution
417:    (as well as the eigenvalues).

419:    Collective on EPS

421:    Input Parameters:
422: +  eps    - the eigensolver context
423: .  etype  - error type
424: -  viewer - optional visualization context

426:    Options Database Key:
427: +  -eps_error_absolute - print absolute errors of each eigenpair
428: .  -eps_error_relative - print relative errors of each eigenpair
429: -  -eps_error_backward - print backward errors of each eigenpair

431:    Notes:
432:    By default, this function checks the error of all eigenpairs and prints
433:    the eigenvalues if all of them are below the requested tolerance.
434:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
435:    eigenvalues and corresponding errors is printed.

437:    Level: intermediate

439: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
440: @*/
441: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
442: {
443:   PetscBool         isascii;
444:   PetscViewerFormat format;
445:   PetscErrorCode    ierr;

449:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
452:   EPSCheckSolved(eps,1);
453:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
454:   if (!isascii) return(0);

456:   PetscViewerGetFormat(viewer,&format);
457:   switch (format) {
458:     case PETSC_VIEWER_DEFAULT:
459:     case PETSC_VIEWER_ASCII_INFO:
460:       EPSErrorView_ASCII(eps,etype,viewer);
461:       break;
462:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
463:       EPSErrorView_DETAIL(eps,etype,viewer);
464:       break;
465:     case PETSC_VIEWER_ASCII_MATLAB:
466:       EPSErrorView_MATLAB(eps,etype,viewer);
467:       break;
468:     default:
469:       PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
470:   }
471:   return(0);
472: }

476: /*@
477:    EPSErrorViewFromOptions - Processes command line options to determine if/how
478:    the errors of the computed solution are to be viewed. 

480:    Collective on EPS

482:    Input Parameters:
483: .  eps - the eigensolver context

485:    Level: developer
486: @*/
487: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
488: {
489:   PetscErrorCode    ierr;
490:   PetscViewer       viewer;
491:   PetscBool         flg;
492:   static PetscBool  incall = PETSC_FALSE;
493:   PetscViewerFormat format;

496:   if (incall) return(0);
497:   incall = PETSC_TRUE;
498:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg);
499:   if (flg) {
500:     PetscViewerPushFormat(viewer,format);
501:     EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer);
502:     PetscViewerPopFormat(viewer);
503:     PetscViewerDestroy(&viewer);
504:   }
505:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg);
506:   if (flg) {
507:     PetscViewerPushFormat(viewer,format);
508:     EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer);
509:     PetscViewerPopFormat(viewer);
510:     PetscViewerDestroy(&viewer);
511:   }
512:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg);
513:   if (flg) {
514:     PetscViewerPushFormat(viewer,format);
515:     EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer);
516:     PetscViewerPopFormat(viewer);
517:     PetscViewerDestroy(&viewer);
518:   }
519:   incall = PETSC_FALSE;
520:   return(0);
521: }

525: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
526: {
528:   PetscDraw      draw;
529:   PetscDrawSP    drawsp;
530:   PetscReal      re,im;
531:   PetscInt       i,k;

534:   if (!eps->nconv) return(0);
535:   PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Computed Eigenvalues",PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
536:   PetscViewerDrawGetDraw(viewer,0,&draw);
537:   PetscDrawSPCreate(draw,1,&drawsp);
538:   for (i=0;i<eps->nconv;i++) {
539:     k = eps->perm[i];
540: #if defined(PETSC_USE_COMPLEX)
541:     re = PetscRealPart(eps->eigr[k]);
542:     im = PetscImaginaryPart(eps->eigr[k]);
543: #else
544:     re = eps->eigr[k];
545:     im = eps->eigi[k];
546: #endif
547:     PetscDrawSPAddPoint(drawsp,&re,&im);
548:   }
549:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
550:   PetscDrawSPDestroy(&drawsp);
551:   PetscViewerDestroy(&viewer);
552:   return(0);
553: }

557: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
558: {
559:   PetscReal      re,im;
560:   PetscInt       i,k;

564:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
565:   for (i=0;i<eps->nconv;i++) {
566:     k = eps->perm[i];
567: #if defined(PETSC_USE_COMPLEX)
568:     re = PetscRealPart(eps->eigr[k]);
569:     im = PetscImaginaryPart(eps->eigr[k]);
570: #else
571:     re = eps->eigr[k];
572:     im = eps->eigi[k];
573: #endif
574:     if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
575:     if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
576:     if (im!=0.0) {
577:       PetscViewerASCIIPrintf(viewer,"   %.5f%+.5fi\n",(double)re,(double)im);
578:     } else {
579:       PetscViewerASCIIPrintf(viewer,"   %.5f\n",(double)re);
580:     }
581:   }
582:   PetscViewerASCIIPrintf(viewer,"\n");
583:   return(0);
584: }

588: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
589: {
591:   PetscInt       i,k;
592:   PetscReal      re,im;
593:   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) {
608:       PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
609:     } else {
610:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
611:     }
612:   }
613:   PetscViewerASCIIPrintf(viewer,"];\n");
614:   return(0);
615: }

619: /*@C
620:    EPSValuesView - Displays the computed eigenvalues in a viewer.

622:    Collective on EPS

624:    Input Parameters:
625: +  eps    - the eigensolver context
626: -  viewer - the viewer

628:    Options Database Key:
629: .  -eps_view_values - print computed eigenvalues

631:    Level: intermediate

633: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
634: @*/
635: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
636: {
637:   PetscBool         isascii,isdraw;
638:   PetscViewerFormat format;
639:   PetscErrorCode    ierr;

643:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
646:   EPSCheckSolved(eps,1);
647:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
648:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
649:   if (isdraw) {
650:     EPSValuesView_DRAW(eps,viewer);
651:   } else if (isascii) {
652:     PetscViewerGetFormat(viewer,&format);
653:     switch (format) {
654:       case PETSC_VIEWER_DEFAULT:
655:       case PETSC_VIEWER_ASCII_INFO:
656:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
657:         EPSValuesView_ASCII(eps,viewer);
658:         break;
659:       case PETSC_VIEWER_ASCII_MATLAB:
660:         EPSValuesView_MATLAB(eps,viewer);
661:         break;
662:       default:
663:         PetscInfo1(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
664:     }
665:   }
666:   return(0);
667: }

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

675:    Collective on EPS

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

680:    Level: developer
681: @*/
682: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
683: {
684:   PetscErrorCode    ierr;
685:   PetscViewer       viewer;
686:   PetscBool         flg;
687:   static PetscBool  incall = PETSC_FALSE;
688:   PetscViewerFormat format;

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

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

709:    Collective on EPS

711:    Parameter:
712: +  eps    - the eigensolver context
713: -  viewer - the viewer

715:    Options Database Keys:
716: .  -eps_view_vectors - output eigenvectors.

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

723:    Level: intermediate

725: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
726: @*/
727: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
728: {
730:   PetscInt       i,k;
731:   Vec            x;
732: #define NMLEN 30
733:   char           vname[NMLEN];
734:   const char     *ename;

738:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
741:   EPSCheckSolved(eps,1);
742:   if (eps->nconv) {
743:     PetscObjectGetName((PetscObject)eps,&ename);
744:     EPSComputeVectors(eps);
745:     for (i=0;i<eps->nconv;i++) {
746:       k = eps->perm[i];
747:       PetscSNPrintf(vname,NMLEN,"V%d_%s",i,ename);
748:       BVGetColumn(eps->V,k,&x);
749:       PetscObjectSetName((PetscObject)x,vname);
750:       VecView(x,viewer);
751:       BVRestoreColumn(eps->V,k,&x);
752:     }
753:   }
754:   return(0);
755: }

759: /*@
760:    EPSVectorsViewFromOptions - Processes command line options to determine if/how
761:    the computed eigenvectors are to be viewed. 

763:    Collective on EPS

765:    Input Parameters:
766: .  eps - the eigensolver context

768:    Level: developer
769: @*/
770: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
771: {
772:   PetscErrorCode    ierr;
773:   PetscViewer       viewer;
774:   PetscBool         flg = PETSC_FALSE;
775:   static PetscBool  incall = PETSC_FALSE;
776:   PetscViewerFormat format;

779:   if (incall) return(0);
780:   incall = PETSC_TRUE;
781:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg);
782:   if (flg) {
783:     PetscViewerPushFormat(viewer,format);
784:     EPSVectorsView(eps,viewer);
785:     PetscViewerPopFormat(viewer);
786:     PetscViewerDestroy(&viewer);
787:   }
788:   incall = PETSC_FALSE;
789:   return(0);
790: }