Actual source code: nepview.c

slepc-3.17.0 2022-03-31
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:    NEP routines related to various viewers
 12: */

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

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

 21:    Collective on nep

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

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

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

 55:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);

 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\n");
 81:         PetscViewerASCIIPrintf(viewer,"    number of terms: %" PetscInt_FMT "\n",nep->nt);
 82:         PetscViewerASCIIPrintf(viewer,"    nonzero pattern of the matrices: %s\n",MatStructures[nep->mstr]);
 83:         break;
 84:       }
 85:     } else PetscViewerASCIIPrintf(viewer,"  nonlinear operator not specified yet\n");
 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) PetscViewerASCIIPrintf(viewer,"not yet set\n");
 90:     else switch (nep->which) {
 91:       case NEP_WHICH_USER:
 92:         PetscViewerASCIIPrintf(viewer,"user defined\n");
 93:         break;
 94:       case NEP_TARGET_MAGNITUDE:
 95:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
 96:         break;
 97:       case NEP_TARGET_REAL:
 98:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
 99:         break;
100:       case NEP_TARGET_IMAGINARY:
101:         PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
102:         break;
103:       case NEP_LARGEST_MAGNITUDE:
104:         PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
105:         break;
106:       case NEP_SMALLEST_MAGNITUDE:
107:         PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
108:         break;
109:       case NEP_LARGEST_REAL:
110:         PetscViewerASCIIPrintf(viewer,"largest real parts\n");
111:         break;
112:       case NEP_SMALLEST_REAL:
113:         PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
114:         break;
115:       case NEP_LARGEST_IMAGINARY:
116:         PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
117:         break;
118:       case NEP_SMALLEST_IMAGINARY:
119:         PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
120:         break;
121:       case NEP_ALL:
122:         PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
123:         break;
124:     }
125:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
126:     if (nep->twosided) PetscViewerASCIIPrintf(viewer,"  using two-sided variant (for left eigenvectors)\n");
127:     PetscViewerASCIIPrintf(viewer,"  number of eigenvalues (nev): %" PetscInt_FMT "\n",nep->nev);
128:     PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %" PetscInt_FMT "\n",nep->ncv);
129:     PetscViewerASCIIPrintf(viewer,"  maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",nep->mpd);
130:     PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %" PetscInt_FMT "\n",nep->max_it);
131:     PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)nep->tol);
132:     PetscViewerASCIIPrintf(viewer,"  convergence test: ");
133:     PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
134:     switch (nep->conv) {
135:     case NEP_CONV_ABS:
136:       PetscViewerASCIIPrintf(viewer,"absolute\n");break;
137:     case NEP_CONV_REL:
138:       PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
139:     case NEP_CONV_NORM:
140:       PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
141:       if (nep->nrma) {
142:         PetscViewerASCIIPrintf(viewer,"  computed matrix norms: %g",(double)nep->nrma[0]);
143:         for (i=1;i<nep->nt;i++) PetscViewerASCIIPrintf(viewer,", %g",(double)nep->nrma[i]);
144:         PetscViewerASCIIPrintf(viewer,"\n");
145:       }
146:       break;
147:     case NEP_CONV_USER:
148:       PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
149:     }
150:     PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
151:     if (nep->refine) {
152:       PetscViewerASCIIPrintf(viewer,"  iterative refinement: %s, with %s scheme\n",NEPRefineTypes[nep->refine],NEPRefineSchemes[nep->scheme]);
153:       PetscViewerASCIIPrintf(viewer,"  refinement stopping criterion: tol=%g, its=%" PetscInt_FMT "\n",(double)nep->rtol,nep->rits);
154:       if (nep->npart>1) PetscViewerASCIIPrintf(viewer,"  splitting communicator in %" PetscInt_FMT " partitions for refinement\n",nep->npart);
155:     }
156:     if (nep->nini) PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(nep->nini));
157:   } else {
158:     if (nep->ops->view) (*nep->ops->view)(nep,viewer);
159:   }
160:   PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
161:   if (!nep->V) NEPGetBV(nep,&nep->V);
162:   BVView(nep->V,viewer);
163:   if (!nep->rg) NEPGetRG(nep,&nep->rg);
164:   RGIsTrivial(nep->rg,&istrivial);
165:   if (!istrivial) RGView(nep->rg,viewer);
166:   if (nep->useds) {
167:     if (!nep->ds) NEPGetDS(nep,&nep->ds);
168:     DSView(nep->ds,viewer);
169:   }
170:   PetscViewerPopFormat(viewer);
171:   if (nep->refine!=NEP_REFINE_NONE) {
172:     PetscViewerASCIIPushTab(viewer);
173:     if (nep->npart>1) {
174:       if (nep->refinesubc->color==0) {
175:         PetscSubcommGetChild(nep->refinesubc,&child);
176:         PetscViewerASCIIGetStdout(child,&sviewer);
177:         KSPView(nep->refineksp,sviewer);
178:       }
179:     } else KSPView(nep->refineksp,viewer);
180:     PetscViewerASCIIPopTab(viewer);
181:   }
182:   PetscFunctionReturn(0);
183: }

185: /*@C
186:    NEPViewFromOptions - View from options

188:    Collective on NEP

190:    Input Parameters:
191: +  nep  - the nonlinear eigensolver context
192: .  obj  - optional object
193: -  name - command line option

195:    Level: intermediate

197: .seealso: NEPView(), NEPCreate()
198: @*/
199: PetscErrorCode NEPViewFromOptions(NEP nep,PetscObject obj,const char name[])
200: {
202:   PetscObjectViewFromOptions((PetscObject)nep,obj,name);
203:   PetscFunctionReturn(0);
204: }

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

209:    Collective on nep

211:    Input Parameters:
212: +  nep - the nonlinear eigensolver context
213: -  viewer - the viewer to display the reason

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

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

224:    Level: intermediate

226: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber(), NEPConvergedReasonViewFromOptions(
227: @*/
228: PetscErrorCode NEPConvergedReasonView(NEP nep,PetscViewer viewer)
229: {
230:   PetscBool         isAscii;
231:   PetscViewerFormat format;

233:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
234:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
235:   if (isAscii) {
236:     PetscViewerGetFormat(viewer,&format);
237:     PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
238:     if (nep->reason > 0 && format != PETSC_VIEWER_FAILED) PetscViewerASCIIPrintf(viewer,"%s Nonlinear eigensolve converged (%" PetscInt_FMT " eigenpair%s) due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)nep)->prefix?((PetscObject)nep)->prefix:"",nep->nconv,(nep->nconv>1)?"s":"",NEPConvergedReasons[nep->reason],nep->its);
239:     else if (nep->reason <= 0) PetscViewerASCIIPrintf(viewer,"%s Nonlinear eigensolve did not converge due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)nep)->prefix?((PetscObject)nep)->prefix:"",NEPConvergedReasons[nep->reason],nep->its);
240:     PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
241:   }
242:   PetscFunctionReturn(0);
243: }

245: /*@
246:    NEPConvergedReasonViewFromOptions - Processes command line options to determine if/how
247:    the NEP converged reason is to be viewed.

249:    Collective on nep

251:    Input Parameter:
252: .  nep - the nonlinear eigensolver context

254:    Level: developer

256: .seealso: NEPConvergedReasonView()
257: @*/
258: PetscErrorCode NEPConvergedReasonViewFromOptions(NEP nep)
259: {
260:   PetscViewer       viewer;
261:   PetscBool         flg;
262:   static PetscBool  incall = PETSC_FALSE;
263:   PetscViewerFormat format;

265:   if (incall) PetscFunctionReturn(0);
266:   incall = PETSC_TRUE;
267:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg);
268:   if (flg) {
269:     PetscViewerPushFormat(viewer,format);
270:     NEPConvergedReasonView(nep,viewer);
271:     PetscViewerPopFormat(viewer);
272:     PetscViewerDestroy(&viewer);
273:   }
274:   incall = PETSC_FALSE;
275:   PetscFunctionReturn(0);
276: }

278: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
279: {
280:   PetscReal      error;
281:   PetscInt       i,j,k,nvals;

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

313: static PetscErrorCode NEPErrorView_DETAIL(NEP nep,NEPErrorType etype,PetscViewer viewer)
314: {
315:   PetscReal      error,re,im;
316:   PetscScalar    kr,ki;
317:   PetscInt       i;
318:   char           ex[30],sep[]=" ---------------------- --------------------\n";

320:   if (!nep->nconv) PetscFunctionReturn(0);
321:   switch (etype) {
322:     case NEP_ERROR_ABSOLUTE:
323:       PetscSNPrintf(ex,sizeof(ex),"    ||T(k)x||");
324:       break;
325:     case NEP_ERROR_RELATIVE:
326:       PetscSNPrintf(ex,sizeof(ex)," ||T(k)x||/||kx||");
327:       break;
328:     case NEP_ERROR_BACKWARD:
329:       PetscSNPrintf(ex,sizeof(ex),"    eta(x,k)");
330:       break;
331:   }
332:   PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep);
333:   for (i=0;i<nep->nconv;i++) {
334:     NEPGetEigenpair(nep,i,&kr,&ki,NULL,NULL);
335:     NEPComputeError(nep,i,etype,&error);
336: #if defined(PETSC_USE_COMPLEX)
337:     re = PetscRealPart(kr);
338:     im = PetscImaginaryPart(kr);
339: #else
340:     re = kr;
341:     im = ki;
342: #endif
343:     if (im!=0.0) PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
344:     else PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
345:   }
346:   PetscViewerASCIIPrintf(viewer,"%s",sep);
347:   PetscFunctionReturn(0);
348: }

350: static PetscErrorCode NEPErrorView_MATLAB(NEP nep,NEPErrorType etype,PetscViewer viewer)
351: {
352:   PetscReal      error;
353:   PetscInt       i;
354:   const char     *name;

356:   PetscObjectGetName((PetscObject)nep,&name);
357:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
358:   for (i=0;i<nep->nconv;i++) {
359:     NEPComputeError(nep,i,etype,&error);
360:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
361:   }
362:   PetscViewerASCIIPrintf(viewer,"];\n");
363:   PetscFunctionReturn(0);
364: }

366: /*@C
367:    NEPErrorView - Displays the errors associated with the computed solution
368:    (as well as the eigenvalues).

370:    Collective on nep

372:    Input Parameters:
373: +  nep    - the nonlinear eigensolver context
374: .  etype  - error type
375: -  viewer - optional visualization context

377:    Options Database Keys:
378: +  -nep_error_absolute - print absolute errors of each eigenpair
379: .  -nep_error_relative - print relative errors of each eigenpair
380: -  -nep_error_backward - print backward errors of each eigenpair

382:    Notes:
383:    By default, this function checks the error of all eigenpairs and prints
384:    the eigenvalues if all of them are below the requested tolerance.
385:    If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
386:    eigenvalues and corresponding errors is printed.

388:    Level: intermediate

390: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
391: @*/
392: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
393: {
394:   PetscBool         isascii;
395:   PetscViewerFormat format;

398:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
401:   NEPCheckSolved(nep,1);
402:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
403:   if (!isascii) PetscFunctionReturn(0);

405:   PetscViewerGetFormat(viewer,&format);
406:   switch (format) {
407:     case PETSC_VIEWER_DEFAULT:
408:     case PETSC_VIEWER_ASCII_INFO:
409:       NEPErrorView_ASCII(nep,etype,viewer);
410:       break;
411:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
412:       NEPErrorView_DETAIL(nep,etype,viewer);
413:       break;
414:     case PETSC_VIEWER_ASCII_MATLAB:
415:       NEPErrorView_MATLAB(nep,etype,viewer);
416:       break;
417:     default:
418:       PetscInfo(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
419:   }
420:   PetscFunctionReturn(0);
421: }

423: /*@
424:    NEPErrorViewFromOptions - Processes command line options to determine if/how
425:    the errors of the computed solution are to be viewed.

427:    Collective on nep

429:    Input Parameter:
430: .  nep - the nonlinear eigensolver context

432:    Level: developer

434: .seealso: NEPErrorView()
435: @*/
436: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
437: {
438:   PetscViewer       viewer;
439:   PetscBool         flg;
440:   static PetscBool  incall = PETSC_FALSE;
441:   PetscViewerFormat format;

443:   if (incall) PetscFunctionReturn(0);
444:   incall = PETSC_TRUE;
445:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg);
446:   if (flg) {
447:     PetscViewerPushFormat(viewer,format);
448:     NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer);
449:     PetscViewerPopFormat(viewer);
450:     PetscViewerDestroy(&viewer);
451:   }
452:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg);
453:   if (flg) {
454:     PetscViewerPushFormat(viewer,format);
455:     NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer);
456:     PetscViewerPopFormat(viewer);
457:     PetscViewerDestroy(&viewer);
458:   }
459:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_backward",&viewer,&format,&flg);
460:   if (flg) {
461:     PetscViewerPushFormat(viewer,format);
462:     NEPErrorView(nep,NEP_ERROR_BACKWARD,viewer);
463:     PetscViewerPopFormat(viewer);
464:     PetscViewerDestroy(&viewer);
465:   }
466:   incall = PETSC_FALSE;
467:   PetscFunctionReturn(0);
468: }

470: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
471: {
472:   PetscDraw      draw;
473:   PetscDrawSP    drawsp;
474:   PetscReal      re,im;
475:   PetscInt       i,k;

477:   if (!nep->nconv) PetscFunctionReturn(0);
478:   PetscViewerDrawGetDraw(viewer,0,&draw);
479:   PetscDrawSetTitle(draw,"Computed Eigenvalues");
480:   PetscDrawSPCreate(draw,1,&drawsp);
481:   for (i=0;i<nep->nconv;i++) {
482:     k = nep->perm[i];
483: #if defined(PETSC_USE_COMPLEX)
484:     re = PetscRealPart(nep->eigr[k]);
485:     im = PetscImaginaryPart(nep->eigr[k]);
486: #else
487:     re = nep->eigr[k];
488:     im = nep->eigi[k];
489: #endif
490:     PetscDrawSPAddPoint(drawsp,&re,&im);
491:   }
492:   PetscDrawSPDraw(drawsp,PETSC_TRUE);
493:   PetscDrawSPSave(drawsp);
494:   PetscDrawSPDestroy(&drawsp);
495:   PetscFunctionReturn(0);
496: }

498: static PetscErrorCode NEPValuesView_BINARY(NEP nep,PetscViewer viewer)
499: {
500: #if defined(PETSC_HAVE_COMPLEX)
501:   PetscInt       i,k;
502:   PetscComplex   *ev;
503: #endif

505: #if defined(PETSC_HAVE_COMPLEX)
506:   PetscMalloc1(nep->nconv,&ev);
507:   for (i=0;i<nep->nconv;i++) {
508:     k = nep->perm[i];
509: #if defined(PETSC_USE_COMPLEX)
510:     ev[i] = nep->eigr[k];
511: #else
512:     ev[i] = PetscCMPLX(nep->eigr[k],nep->eigi[k]);
513: #endif
514:   }
515:   PetscViewerBinaryWrite(viewer,ev,nep->nconv,PETSC_COMPLEX);
516:   PetscFree(ev);
517: #endif
518:   PetscFunctionReturn(0);
519: }

521: #if defined(PETSC_HAVE_HDF5)
522: static PetscErrorCode NEPValuesView_HDF5(NEP nep,PetscViewer viewer)
523: {
524:   PetscInt       i,k,n,N;
525:   PetscMPIInt    rank;
526:   Vec            v;
527:   char           vname[30];
528:   const char     *ename;

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

566: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
567: {
568:   PetscInt       i,k;

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

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

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

605: /*@C
606:    NEPValuesView - Displays the computed eigenvalues in a viewer.

608:    Collective on nep

610:    Input Parameters:
611: +  nep    - the nonlinear eigensolver context
612: -  viewer - the viewer

614:    Options Database Key:
615: .  -nep_view_values - print computed eigenvalues

617:    Level: intermediate

619: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
620: @*/
621: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
622: {
623:   PetscBool         isascii,isdraw,isbinary;
624:   PetscViewerFormat format;
625: #if defined(PETSC_HAVE_HDF5)
626:   PetscBool         ishdf5;
627: #endif

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

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

667:    Collective on nep

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

672:    Level: developer

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

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

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

699:    Collective on nep

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

705:    Options Database Key:
706: .  -nep_view_vectors - output eigenvectors.

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

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

717:    Level: intermediate

719: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
720: @*/
721: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
722: {
723:   PetscInt       i,k;
724:   Vec            xr,xi=NULL;

727:   if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
730:   NEPCheckSolved(nep,1);
731:   if (nep->nconv) {
732:     NEPComputeVectors(nep);
733:     BVCreateVec(nep->V,&xr);
734: #if !defined(PETSC_USE_COMPLEX)
735:     BVCreateVec(nep->V,&xi);
736: #endif
737:     for (i=0;i<nep->nconv;i++) {
738:       k = nep->perm[i];
739:       BV_GetEigenvector(nep->V,k,nep->eigi[k],xr,xi);
740:       SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)nep);
741:       if (nep->twosided) {
742:         BV_GetEigenvector(nep->W,k,nep->eigi[k],xr,xi);
743:         SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)nep);
744:       }
745:     }
746:     VecDestroy(&xr);
747: #if !defined(PETSC_USE_COMPLEX)
748:     VecDestroy(&xi);
749: #endif
750:   }
751:   PetscFunctionReturn(0);
752: }

754: /*@
755:    NEPVectorsViewFromOptions - Processes command line options to determine if/how
756:    the computed eigenvectors are to be viewed.

758:    Collective on nep

760:    Input Parameter:
761: .  nep - the nonlinear eigensolver context

763:    Level: developer

765: .seealso: NEPVectorsView()
766: @*/
767: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
768: {
769:   PetscViewer       viewer;
770:   PetscBool         flg = PETSC_FALSE;
771:   static PetscBool  incall = PETSC_FALSE;
772:   PetscViewerFormat format;

774:   if (incall) PetscFunctionReturn(0);
775:   incall = PETSC_TRUE;
776:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg);
777:   if (flg) {
778:     PetscViewerPushFormat(viewer,format);
779:     NEPVectorsView(nep,viewer);
780:     PetscViewerPopFormat(viewer);
781:     PetscViewerDestroy(&viewer);
782:   }
783:   incall = PETSC_FALSE;
784:   PetscFunctionReturn(0);
785: }