Actual source code: nepview.c

slepc-3.14.1 2020-12-08
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2020, 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:    NEP routines related to various viewers
 12: */

 14: #include <slepc/private/nepimpl.h>
 15: #include <petscdraw.h>

 17: /*@C
 18:    NEPView - Prints the NEP data structure.

 20:    Collective on nep

 22:    Input Parameters:
 23: +  nep - the nonlinear eigenproblem solver context
 24: -  viewer - optional visualization context

 26:    Options Database Key:
 27: .  -nep_view -  Calls NEPView() at end of NEPSolve()

 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
 41: @*/
 42: PetscErrorCode NEPView(NEP nep,PetscViewer viewer)
 43: {
 45:   const char     *type=NULL;
 46:   char           str[50];
 47:   PetscInt       i;
 48:   PetscBool      isascii,istrivial;
 49:   PetscViewer    sviewer;

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

 59:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
 60:   if (isascii) {
 61:     PetscObjectPrintClassNamePrefixType((PetscObject)nep,viewer);
 62:     if (nep->ops->view) {
 63:       PetscViewerASCIIPushTab(viewer);
 64:       (*nep->ops->view)(nep,viewer);
 65:       PetscViewerASCIIPopTab(viewer);
 66:     }
 67:     if (nep->problem_type) {
 68:       switch (nep->problem_type) {
 69:         case NEP_GENERAL:  type = "general nonlinear eigenvalue problem"; break;
 70:         case NEP_RATIONAL: type = "rational eigenvalue problem"; break;
 71:       }
 72:     } else type = "not yet set";
 73:     PetscViewerASCIIPrintf(viewer,"  problem type: %s\n",type);
 74:     if (nep->fui) {
 75:       switch (nep->fui) {
 76:       case NEP_USER_INTERFACE_CALLBACK:
 77:         PetscViewerASCIIPrintf(viewer,"  nonlinear operator from user callbacks\n");
 78:         break;
 79:       case NEP_USER_INTERFACE_SPLIT:
 80:         PetscViewerASCIIPrintf(viewer,"  nonlinear operator in split form, with %D terms\n",nep->nt);
 81:         break;
 82:       }
 83:     } else {
 84:       PetscViewerASCIIPrintf(viewer,"  nonlinear operator not specified yet\n");
 85:     }
 86:     PetscViewerASCIIPrintf(viewer,"  selected portion of the spectrum: ");
 87:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
 88:     SlepcSNPrintfScalar(str,sizeof(str),nep->target,PETSC_FALSE);
 89:     if (!nep->which) {
 90:       PetscViewerASCIIPrintf(viewer,"not yet set\n");
 91:     } else switch (nep->which) {
 92:       case NEP_WHICH_USER:
 93:         PetscViewerASCIIPrintf(viewer,"user defined\n");
 94:         break;
 95:       case NEP_TARGET_MAGNITUDE:
 96:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
 97:         break;
 98:       case NEP_TARGET_REAL:
 99:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
100:         break;
101:       case NEP_TARGET_IMAGINARY:
102:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
103:         break;
104:       case NEP_LARGEST_MAGNITUDE:
105:         PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
106:         break;
107:       case NEP_SMALLEST_MAGNITUDE:
108:         PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
109:         break;
110:       case NEP_LARGEST_REAL:
111:         PetscViewerASCIIPrintf(viewer,"largest real parts\n");
112:         break;
113:       case NEP_SMALLEST_REAL:
114:         PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
115:         break;
116:       case NEP_LARGEST_IMAGINARY:
117:         PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
118:         break;
119:       case NEP_SMALLEST_IMAGINARY:
120:         PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
121:         break;
122:       case NEP_ALL:
123:         PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
124:         break;
125:     }
126:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
127:     if (nep->twosided) {
128:       PetscViewerASCIIPrintf(viewer,"  using two-sided variant (for left eigenvectors)\n");
129:     }
130:     PetscViewerASCIIPrintf(viewer,"  number of eigenvalues (nev): %D\n",nep->nev);
131:     PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %D\n",nep->ncv);
132:     PetscViewerASCIIPrintf(viewer,"  maximum dimension of projected problem (mpd): %D\n",nep->mpd);
133:     PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %D\n",nep->max_it);
134:     PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)nep->tol);
135:     PetscViewerASCIIPrintf(viewer,"  convergence test: ");
136:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
137:     switch (nep->conv) {
138:     case NEP_CONV_ABS:
139:       PetscViewerASCIIPrintf(viewer,"absolute\n");break;
140:     case NEP_CONV_REL:
141:       PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
142:     case NEP_CONV_NORM:
143:       PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
144:       if (nep->nrma) {
145:         PetscViewerASCIIPrintf(viewer,"  computed matrix norms: %g",(double)nep->nrma[0]);
146:         for (i=1;i<nep->nt;i++) {
147:           PetscViewerASCIIPrintf(viewer,", %g",(double)nep->nrma[i]);
148:         }
149:         PetscViewerASCIIPrintf(viewer,"\n");
150:       }
151:       break;
152:     case NEP_CONV_USER:
153:       PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
154:     }
155:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
156:     if (nep->refine) {
157:       PetscViewerASCIIPrintf(viewer,"  iterative refinement: %s, with %s scheme\n",NEPRefineTypes[nep->refine],NEPRefineSchemes[nep->scheme]);
158:       PetscViewerASCIIPrintf(viewer,"  refinement stopping criterion: tol=%g, its=%D\n",(double)nep->rtol,nep->rits);
159:       if (nep->npart>1) {
160:         PetscViewerASCIIPrintf(viewer,"  splitting communicator in %D partitions for refinement\n",nep->npart);
161:       }
162:     }
163:     if (nep->nini) {
164:       PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial space: %D\n",PetscAbs(nep->nini));
165:     }
166:   } else {
167:     if (nep->ops->view) {
168:       (*nep->ops->view)(nep,viewer);
169:     }
170:   }
171:   PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
172:   if (!nep->V) { NEPGetBV(nep,&nep->V); }
173:   BVView(nep->V,viewer);
174:   if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
175:   RGIsTrivial(nep->rg,&istrivial);
176:   if (!istrivial) { RGView(nep->rg,viewer); }
177:   if (nep->useds) {
178:     if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
179:     DSView(nep->ds,viewer);
180:   }
181:   PetscViewerPopFormat(viewer);
182:   if (nep->refine!=NEP_REFINE_NONE) {
183:     PetscViewerASCIIPushTab(viewer);
184:     if (nep->npart>1) {
185:       if (nep->refinesubc->color==0) {
186:         PetscViewerASCIIGetStdout(PetscSubcommChild(nep->refinesubc),&sviewer);
187:         KSPView(nep->refineksp,sviewer);
188:       }
189:     } else {
190:       KSPView(nep->refineksp,viewer);
191:     }
192:     PetscViewerASCIIPopTab(viewer);
193:   }
194:   return(0);
195: }

197: /*@C
198:    NEPViewFromOptions - View from options

200:    Collective on NEP

202:    Input Parameters:
203: +  nep  - the nonlinear eigensolver context
204: .  obj  - optional object
205: -  name - command line option

207:    Level: intermediate

209: .seealso: NEPView(), NEPCreate()
210: @*/
211: PetscErrorCode NEPViewFromOptions(NEP nep,PetscObject obj,const char name[])
212: {

217:   PetscObjectViewFromOptions((PetscObject)nep,obj,name);
218:   return(0);
219: }

221: /*@C
222:    NEPConvergedReasonView - Displays the reason a NEP solve converged or diverged.

224:    Collective on nep

226:    Input Parameters:
227: +  nep - the nonlinear eigensolver context
228: -  viewer - the viewer to display the reason

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

233:    Note:
234:    To change the format of the output call PetscViewerPushFormat(viewer,format) before
235:    this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
236:    display a reason if it fails. The latter can be set in the command line with
237:    -nep_converged_reason ::failed

239:    Level: intermediate

241: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber(), NEPConvergedReasonViewFromOptions(
242: @*/
243: PetscErrorCode NEPConvergedReasonView(NEP nep,PetscViewer viewer)
244: {
245:   PetscErrorCode    ierr;
246:   PetscBool         isAscii;
247:   PetscViewerFormat format;

250:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
251:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
252:   if (isAscii) {
253:     PetscViewerGetFormat(viewer,&format);
254:     PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
255:     if (nep->reason > 0 && format != PETSC_VIEWER_FAILED) {
256:       PetscViewerASCIIPrintf(viewer,"%s Nonlinear eigensolve converged (%D eigenpair%s) due to %s; iterations %D\n",((PetscObject)nep)->prefix?((PetscObject)nep)->prefix:"",nep->nconv,(nep->nconv>1)?"s":"",NEPConvergedReasons[nep->reason],nep->its);
257:     } else if (nep->reason <= 0) {
258:       PetscViewerASCIIPrintf(viewer,"%s Nonlinear eigensolve did not converge due to %s; iterations %D\n",((PetscObject)nep)->prefix?((PetscObject)nep)->prefix:"",NEPConvergedReasons[nep->reason],nep->its);
259:     }
260:     PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
261:   }
262:   return(0);
263: }

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

269:    Collective on nep

271:    Input Parameter:
272: .  nep - the nonlinear eigensolver context

274:    Level: developer

276: .seealso: NEPConvergedReasonView()
277: @*/
278: PetscErrorCode NEPConvergedReasonViewFromOptions(NEP nep)
279: {
280:   PetscErrorCode    ierr;
281:   PetscViewer       viewer;
282:   PetscBool         flg;
283:   static PetscBool  incall = PETSC_FALSE;
284:   PetscViewerFormat format;

287:   if (incall) return(0);
288:   incall = PETSC_TRUE;
289:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg);
290:   if (flg) {
291:     PetscViewerPushFormat(viewer,format);
292:     NEPConvergedReasonView(nep,viewer);
293:     PetscViewerPopFormat(viewer);
294:     PetscViewerDestroy(&viewer);
295:   }
296:   incall = PETSC_FALSE;
297:   return(0);
298: }

300: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
301: {
302:   PetscReal      error;
303:   PetscInt       i,j,k,nvals;

307:   nvals = (nep->which==NEP_ALL)? nep->nconv: nep->nev;
308:   if (nep->which!=NEP_ALL && nep->nconv<nvals) {
309:     PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",nep->nev);
310:     return(0);
311:   }
312:   if (nep->which==NEP_ALL && !nvals) {
313:     PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
314:     return(0);
315:   }
316:   for (i=0;i<nvals;i++) {
317:     NEPComputeError(nep,i,etype,&error);
318:     if (error>=5.0*nep->tol) {
319:       PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
320:       return(0);
321:     }
322:   }
323:   if (nep->which==NEP_ALL) {
324:     PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
325:   } else {
326:     PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
327:   }
328:   for (i=0;i<=(nvals-1)/8;i++) {
329:     PetscViewerASCIIPrintf(viewer,"\n     ");
330:     for (j=0;j<PetscMin(8,nvals-8*i);j++) {
331:       k = nep->perm[8*i+j];
332:       SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]);
333:       if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
334:     }
335:   }
336:   PetscViewerASCIIPrintf(viewer,"\n\n");
337:   return(0);
338: }

340: static PetscErrorCode NEPErrorView_DETAIL(NEP nep,NEPErrorType etype,PetscViewer viewer)
341: {
343:   PetscReal      error,re,im;
344:   PetscScalar    kr,ki;
345:   PetscInt       i;
346:   char           ex[30],sep[]=" ---------------------- --------------------\n";

349:   if (!nep->nconv) return(0);
350:   switch (etype) {
351:     case NEP_ERROR_ABSOLUTE:
352:       PetscSNPrintf(ex,sizeof(ex),"    ||T(k)x||");
353:       break;
354:     case NEP_ERROR_RELATIVE:
355:       PetscSNPrintf(ex,sizeof(ex)," ||T(k)x||/||kx||");
356:       break;
357:     case NEP_ERROR_BACKWARD:
358:       PetscSNPrintf(ex,sizeof(ex),"    eta(x,k)");
359:       break;
360:   }
361:   PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep);
362:   for (i=0;i<nep->nconv;i++) {
363:     NEPGetEigenpair(nep,i,&kr,&ki,NULL,NULL);
364:     NEPComputeError(nep,i,etype,&error);
365: #if defined(PETSC_USE_COMPLEX)
366:     re = PetscRealPart(kr);
367:     im = PetscImaginaryPart(kr);
368: #else
369:     re = kr;
370:     im = ki;
371: #endif
372:     if (im!=0.0) {
373:       PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
374:     } else {
375:       PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
376:     }
377:   }
378:   PetscViewerASCIIPrintf(viewer,sep);
379:   return(0);
380: }

382: static PetscErrorCode NEPErrorView_MATLAB(NEP nep,NEPErrorType etype,PetscViewer viewer)
383: {
385:   PetscReal      error;
386:   PetscInt       i;
387:   const char     *name;

390:   PetscObjectGetName((PetscObject)nep,&name);
391:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
392:   for (i=0;i<nep->nconv;i++) {
393:     NEPComputeError(nep,i,etype,&error);
394:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
395:   }
396:   PetscViewerASCIIPrintf(viewer,"];\n");
397:   return(0);
398: }

400: /*@C
401:    NEPErrorView - Displays the errors associated with the computed solution
402:    (as well as the eigenvalues).

404:    Collective on nep

406:    Input Parameters:
407: +  nep    - the nonlinear eigensolver context
408: .  etype  - error type
409: -  viewer - optional visualization context

411:    Options Database Key:
412: +  -nep_error_absolute - print absolute errors of each eigenpair
413: .  -nep_error_relative - print relative errors of each eigenpair
414: -  -nep_error_backward - print backward errors of each eigenpair

416:    Notes:
417:    By default, this function checks the error of all eigenpairs and prints
418:    the eigenvalues if all of them are below the requested tolerance.
419:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
420:    eigenvalues and corresponding errors is printed.

422:    Level: intermediate

424: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
425: @*/
426: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
427: {
428:   PetscBool         isascii;
429:   PetscViewerFormat format;
430:   PetscErrorCode    ierr;

434:   if (!viewer) {
435:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
436:   }
439:   NEPCheckSolved(nep,1);
440:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
441:   if (!isascii) return(0);

443:   PetscViewerGetFormat(viewer,&format);
444:   switch (format) {
445:     case PETSC_VIEWER_DEFAULT:
446:     case PETSC_VIEWER_ASCII_INFO:
447:       NEPErrorView_ASCII(nep,etype,viewer);
448:       break;
449:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
450:       NEPErrorView_DETAIL(nep,etype,viewer);
451:       break;
452:     case PETSC_VIEWER_ASCII_MATLAB:
453:       NEPErrorView_MATLAB(nep,etype,viewer);
454:       break;
455:     default:
456:       PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
457:   }
458:   return(0);
459: }

461: /*@
462:    NEPErrorViewFromOptions - Processes command line options to determine if/how
463:    the errors of the computed solution are to be viewed.

465:    Collective on nep

467:    Input Parameter:
468: .  nep - the nonlinear eigensolver context

470:    Level: developer
471: @*/
472: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
473: {
474:   PetscErrorCode    ierr;
475:   PetscViewer       viewer;
476:   PetscBool         flg;
477:   static PetscBool  incall = PETSC_FALSE;
478:   PetscViewerFormat format;

481:   if (incall) return(0);
482:   incall = PETSC_TRUE;
483:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg);
484:   if (flg) {
485:     PetscViewerPushFormat(viewer,format);
486:     NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer);
487:     PetscViewerPopFormat(viewer);
488:     PetscViewerDestroy(&viewer);
489:   }
490:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg);
491:   if (flg) {
492:     PetscViewerPushFormat(viewer,format);
493:     NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer);
494:     PetscViewerPopFormat(viewer);
495:     PetscViewerDestroy(&viewer);
496:   }
497:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_backward",&viewer,&format,&flg);
498:   if (flg) {
499:     PetscViewerPushFormat(viewer,format);
500:     NEPErrorView(nep,NEP_ERROR_BACKWARD,viewer);
501:     PetscViewerPopFormat(viewer);
502:     PetscViewerDestroy(&viewer);
503:   }
504:   incall = PETSC_FALSE;
505:   return(0);
506: }

508: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
509: {
511:   PetscDraw      draw;
512:   PetscDrawSP    drawsp;
513:   PetscReal      re,im;
514:   PetscInt       i,k;

517:   if (!nep->nconv) return(0);
518:   PetscViewerDrawGetDraw(viewer,0,&draw);
519:   PetscDrawSetTitle(draw,"Computed Eigenvalues");
520:   PetscDrawSPCreate(draw,1,&drawsp);
521:   for (i=0;i<nep->nconv;i++) {
522:     k = nep->perm[i];
523: #if defined(PETSC_USE_COMPLEX)
524:     re = PetscRealPart(nep->eigr[k]);
525:     im = PetscImaginaryPart(nep->eigr[k]);
526: #else
527:     re = nep->eigr[k];
528:     im = nep->eigi[k];
529: #endif
530:     PetscDrawSPAddPoint(drawsp,&re,&im);
531:   }
532:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
533:   PetscDrawSPSave(drawsp);
534:   PetscDrawSPDestroy(&drawsp);
535:   return(0);
536: }

538: static PetscErrorCode NEPValuesView_BINARY(NEP nep,PetscViewer viewer)
539: {
540: #if defined(PETSC_HAVE_COMPLEX)
541:   PetscInt       i,k;
542:   PetscComplex   *ev;
544: #endif

547: #if defined(PETSC_HAVE_COMPLEX)
548:   PetscMalloc1(nep->nconv,&ev);
549:   for (i=0;i<nep->nconv;i++) {
550:     k = nep->perm[i];
551: #if defined(PETSC_USE_COMPLEX)
552:     ev[i] = nep->eigr[k];
553: #else
554:     ev[i] = PetscCMPLX(nep->eigr[k],nep->eigi[k]);
555: #endif
556:   }
557:   PetscViewerBinaryWrite(viewer,ev,nep->nconv,PETSC_COMPLEX);
558:   PetscFree(ev);
559: #endif
560:   return(0);
561: }

563: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
564: {
565:   PetscInt       i,k;

569:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
570:   for (i=0;i<nep->nconv;i++) {
571:     k = nep->perm[i];
572:     PetscViewerASCIIPrintf(viewer,"   ");
573:     SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]);
574:     PetscViewerASCIIPrintf(viewer,"\n");
575:   }
576:   PetscViewerASCIIPrintf(viewer,"\n");
577:   return(0);
578: }

580: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
581: {
583:   PetscInt       i,k;
584:   PetscReal      re,im;
585:   const char     *name;

588:   PetscObjectGetName((PetscObject)nep,&name);
589:   PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
590:   for (i=0;i<nep->nconv;i++) {
591:     k = nep->perm[i];
592: #if defined(PETSC_USE_COMPLEX)
593:     re = PetscRealPart(nep->eigr[k]);
594:     im = PetscImaginaryPart(nep->eigr[k]);
595: #else
596:     re = nep->eigr[k];
597:     im = nep->eigi[k];
598: #endif
599:     if (im!=0.0) {
600:       PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
601:     } else {
602:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
603:     }
604:   }
605:   PetscViewerASCIIPrintf(viewer,"];\n");
606:   return(0);
607: }

609: /*@C
610:    NEPValuesView - Displays the computed eigenvalues in a viewer.

612:    Collective on nep

614:    Input Parameters:
615: +  nep    - the nonlinear eigensolver context
616: -  viewer - the viewer

618:    Options Database Key:
619: .  -nep_view_values - print computed eigenvalues

621:    Level: intermediate

623: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
624: @*/
625: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
626: {
627:   PetscBool         isascii,isdraw,isbinary;
628:   PetscViewerFormat format;
629:   PetscErrorCode    ierr;

633:   if (!viewer) {
634:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
635:   }
638:   NEPCheckSolved(nep,1);
639:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
640:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
641:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
642:   if (isdraw) {
643:     NEPValuesView_DRAW(nep,viewer);
644:   } else if (isbinary) {
645:     NEPValuesView_BINARY(nep,viewer);
646:   } else if (isascii) {
647:     PetscViewerGetFormat(viewer,&format);
648:     switch (format) {
649:       case PETSC_VIEWER_DEFAULT:
650:       case PETSC_VIEWER_ASCII_INFO:
651:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
652:         NEPValuesView_ASCII(nep,viewer);
653:         break;
654:       case PETSC_VIEWER_ASCII_MATLAB:
655:         NEPValuesView_MATLAB(nep,viewer);
656:         break;
657:       default:
658:         PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
659:     }
660:   }
661:   return(0);
662: }

664: /*@
665:    NEPValuesViewFromOptions - Processes command line options to determine if/how
666:    the computed eigenvalues are to be viewed.

668:    Collective on nep

670:    Input Parameter:
671: .  nep - the nonlinear eigensolver context

673:    Level: developer
674: @*/
675: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
676: {
677:   PetscErrorCode    ierr;
678:   PetscViewer       viewer;
679:   PetscBool         flg;
680:   static PetscBool  incall = PETSC_FALSE;
681:   PetscViewerFormat format;

684:   if (incall) return(0);
685:   incall = PETSC_TRUE;
686:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg);
687:   if (flg) {
688:     PetscViewerPushFormat(viewer,format);
689:     NEPValuesView(nep,viewer);
690:     PetscViewerPopFormat(viewer);
691:     PetscViewerDestroy(&viewer);
692:   }
693:   incall = PETSC_FALSE;
694:   return(0);
695: }

697: /*@C
698:    NEPVectorsView - Outputs computed eigenvectors to a viewer.

700:    Collective on nep

702:    Input Parameters:
703: +  nep    - the nonlinear eigensolver context
704: -  viewer - the viewer

706:    Options Database Keys:
707: .  -nep_view_vectors - output eigenvectors.

709:    Notes:
710:    If PETSc was configured with real scalars, complex conjugate eigenvectors
711:    will be viewed as two separate real vectors, one containing the real part
712:    and another one containing the imaginary part.

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

718:    Level: intermediate

720: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
721: @*/
722: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
723: {
725:   PetscInt       i,k;
726:   Vec            x;
727:   char           vname[30];
728:   const char     *ename;

732:   if (!viewer) {
733:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
734:   }
737:   NEPCheckSolved(nep,1);
738:   if (nep->nconv) {
739:     PetscObjectGetName((PetscObject)nep,&ename);
740:     NEPComputeVectors(nep);
741:     for (i=0;i<nep->nconv;i++) {
742:       k = nep->perm[i];
743:       PetscSNPrintf(vname,sizeof(vname),"X%d_%s",(int)i,ename);
744:       BVGetColumn(nep->V,k,&x);
745:       PetscObjectSetName((PetscObject)x,vname);
746:       VecView(x,viewer);
747:       BVRestoreColumn(nep->V,k,&x);
748:       if (nep->twosided) {
749:         PetscSNPrintf(vname,sizeof(vname),"Y%d_%s",(int)i,ename);
750:         BVGetColumn(nep->W,k,&x);
751:         PetscObjectSetName((PetscObject)x,vname);
752:         VecView(x,viewer);
753:         BVRestoreColumn(nep->W,k,&x);
754:       }
755:     }
756:   }
757:   return(0);
758: }

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

764:    Collective on nep

766:    Input Parameter:
767: .  nep - the nonlinear eigensolver context

769:    Level: developer
770: @*/
771: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
772: {
773:   PetscErrorCode    ierr;
774:   PetscViewer       viewer;
775:   PetscBool         flg = PETSC_FALSE;
776:   static PetscBool  incall = PETSC_FALSE;
777:   PetscViewerFormat format;

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