Actual source code: nepview.c

slepc-3.11.2 2019-07-30
Report Typos and Errors
  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-2019, 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>      /*I "slepcnep.h" I*/
 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

 42: .seealso: PetscViewerASCIIOpen()
 43: @*/
 44: PetscErrorCode NEPView(NEP nep,PetscViewer viewer)
 45: {
 47:   const char     *type=NULL;
 48:   char           str[50];
 49:   PetscInt       i;
 50:   PetscBool      isascii,istrivial;
 51:   PetscViewer    sviewer;

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

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

202: /*@C
203:    NEPReasonView - Displays the reason a NEP solve converged or diverged.

205:    Collective on NEP

207:    Parameter:
208: +  nep - the nonlinear eigensolver context
209: -  viewer - the viewer to display the reason

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

214:    Level: intermediate

216: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber()
217: @*/
218: PetscErrorCode NEPReasonView(NEP nep,PetscViewer viewer)
219: {
221:   PetscBool      isAscii;

224:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
225:   if (isAscii) {
226:     PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
227:     if (nep->reason > 0) {
228:       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);
229:     } else {
230:       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);
231:     }
232:     PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
233:   }
234:   return(0);
235: }

237: /*@
238:    NEPReasonViewFromOptions - Processes command line options to determine if/how
239:    the NEP converged reason is to be viewed.

241:    Collective on NEP

243:    Input Parameters:
244: .  nep - the nonlinear eigensolver context

246:    Level: developer
247: @*/
248: PetscErrorCode NEPReasonViewFromOptions(NEP nep)
249: {
250:   PetscErrorCode    ierr;
251:   PetscViewer       viewer;
252:   PetscBool         flg;
253:   static PetscBool  incall = PETSC_FALSE;
254:   PetscViewerFormat format;

257:   if (incall) return(0);
258:   incall = PETSC_TRUE;
259:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg);
260:   if (flg) {
261:     PetscViewerPushFormat(viewer,format);
262:     NEPReasonView(nep,viewer);
263:     PetscViewerPopFormat(viewer);
264:     PetscViewerDestroy(&viewer);
265:   }
266:   incall = PETSC_FALSE;
267:   return(0);
268: }

270: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
271: {
272:   PetscReal      error;
273:   PetscInt       i,j,k,nvals;

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

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

320:   if (!nep->nconv) return(0);
321:   switch (etype) {
322:     case NEP_ERROR_ABSOLUTE:
323:       PetscSNPrintf(ex,EXLEN,"    ||T(k)x||");
324:       break;
325:     case NEP_ERROR_RELATIVE:
326:       PetscSNPrintf(ex,EXLEN," ||T(k)x||/||kx||");
327:       break;
328:     case NEP_ERROR_BACKWARD:
329:       PetscSNPrintf(ex,EXLEN,"    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) {
344:       PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error);
345:     } else {
346:       PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error);
347:     }
348:   }
349:   PetscViewerASCIIPrintf(viewer,sep);
350:   return(0);
351: }

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

361:   PetscObjectGetName((PetscObject)nep,&name);
362:   PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
363:   for (i=0;i<nep->nconv;i++) {
364:     NEPComputeError(nep,i,etype,&error);
365:     PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
366:   }
367:   PetscViewerASCIIPrintf(viewer,"];\n");
368:   return(0);
369: }

371: /*@C
372:    NEPErrorView - Displays the errors associated with the computed solution
373:    (as well as the eigenvalues).

375:    Collective on NEP

377:    Input Parameters:
378: +  nep    - the nonlinear eigensolver context
379: .  etype  - error type
380: -  viewer - optional visualization context

382:    Options Database Key:
383: +  -nep_error_absolute - print absolute errors of each eigenpair
384: .  -nep_error_relative - print relative errors of each eigenpair
385: -  -nep_error_backward - print backward errors of each eigenpair

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

393:    Level: intermediate

395: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
396: @*/
397: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
398: {
399:   PetscBool         isascii;
400:   PetscViewerFormat format;
401:   PetscErrorCode    ierr;

405:   if (!viewer) {
406:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
407:   }
410:   NEPCheckSolved(nep,1);
411:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
412:   if (!isascii) return(0);

414:   PetscViewerGetFormat(viewer,&format);
415:   switch (format) {
416:     case PETSC_VIEWER_DEFAULT:
417:     case PETSC_VIEWER_ASCII_INFO:
418:       NEPErrorView_ASCII(nep,etype,viewer);
419:       break;
420:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
421:       NEPErrorView_DETAIL(nep,etype,viewer);
422:       break;
423:     case PETSC_VIEWER_ASCII_MATLAB:
424:       NEPErrorView_MATLAB(nep,etype,viewer);
425:       break;
426:     default:
427:       PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
428:   }
429:   return(0);
430: }

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

436:    Collective on NEP

438:    Input Parameters:
439: .  nep - the nonlinear eigensolver context

441:    Level: developer
442: @*/
443: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
444: {
445:   PetscErrorCode    ierr;
446:   PetscViewer       viewer;
447:   PetscBool         flg;
448:   static PetscBool  incall = PETSC_FALSE;
449:   PetscViewerFormat format;

452:   if (incall) return(0);
453:   incall = PETSC_TRUE;
454:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg);
455:   if (flg) {
456:     PetscViewerPushFormat(viewer,format);
457:     NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer);
458:     PetscViewerPopFormat(viewer);
459:     PetscViewerDestroy(&viewer);
460:   }
461:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg);
462:   if (flg) {
463:     PetscViewerPushFormat(viewer,format);
464:     NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer);
465:     PetscViewerPopFormat(viewer);
466:     PetscViewerDestroy(&viewer);
467:   }
468:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_backward",&viewer,&format,&flg);
469:   if (flg) {
470:     PetscViewerPushFormat(viewer,format);
471:     NEPErrorView(nep,NEP_ERROR_BACKWARD,viewer);
472:     PetscViewerPopFormat(viewer);
473:     PetscViewerDestroy(&viewer);
474:   }
475:   incall = PETSC_FALSE;
476:   return(0);
477: }

479: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
480: {
482:   PetscDraw      draw;
483:   PetscDrawSP    drawsp;
484:   PetscReal      re,im;
485:   PetscInt       i,k;

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

509: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
510: {
511:   PetscInt       i,k;

515:   PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
516:   for (i=0;i<nep->nconv;i++) {
517:     k = nep->perm[i];
518:     PetscViewerASCIIPrintf(viewer,"   ");
519:     SlepcPrintEigenvalueASCII(nep->eigr[k],nep->eigi[k]);
520:     PetscViewerASCIIPrintf(viewer,"\n");
521:   }
522:   PetscViewerASCIIPrintf(viewer,"\n");
523:   return(0);
524: }

526: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
527: {
529:   PetscInt       i,k;
530:   PetscReal      re,im;
531:   const char     *name;

534:   PetscObjectGetName((PetscObject)nep,&name);
535:   PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
536:   for (i=0;i<nep->nconv;i++) {
537:     k = nep->perm[i];
538: #if defined(PETSC_USE_COMPLEX)
539:     re = PetscRealPart(nep->eigr[k]);
540:     im = PetscImaginaryPart(nep->eigr[k]);
541: #else
542:     re = nep->eigr[k];
543:     im = nep->eigi[k];
544: #endif
545:     if (im!=0.0) {
546:       PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
547:     } else {
548:       PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
549:     }
550:   }
551:   PetscViewerASCIIPrintf(viewer,"];\n");
552:   return(0);
553: }

555: /*@C
556:    NEPValuesView - Displays the computed eigenvalues in a viewer.

558:    Collective on NEP

560:    Input Parameters:
561: +  nep    - the nonlinear eigensolver context
562: -  viewer - the viewer

564:    Options Database Key:
565: .  -nep_view_values - print computed eigenvalues

567:    Level: intermediate

569: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
570: @*/
571: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
572: {
573:   PetscBool         isascii,isdraw;
574:   PetscViewerFormat format;
575:   PetscErrorCode    ierr;

579:   if (!viewer) {
580:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
581:   }
584:   NEPCheckSolved(nep,1);
585:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
586:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
587:   if (isdraw) {
588:     NEPValuesView_DRAW(nep,viewer);
589:   } else if (isascii) {
590:     PetscViewerGetFormat(viewer,&format);
591:     switch (format) {
592:       case PETSC_VIEWER_DEFAULT:
593:       case PETSC_VIEWER_ASCII_INFO:
594:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
595:         NEPValuesView_ASCII(nep,viewer);
596:         break;
597:       case PETSC_VIEWER_ASCII_MATLAB:
598:         NEPValuesView_MATLAB(nep,viewer);
599:         break;
600:       default:
601:         PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
602:     }
603:   }
604:   return(0);
605: }

607: /*@
608:    NEPValuesViewFromOptions - Processes command line options to determine if/how
609:    the computed eigenvalues are to be viewed.

611:    Collective on NEP

613:    Input Parameters:
614: .  nep - the nonlinear eigensolver context

616:    Level: developer
617: @*/
618: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
619: {
620:   PetscErrorCode    ierr;
621:   PetscViewer       viewer;
622:   PetscBool         flg;
623:   static PetscBool  incall = PETSC_FALSE;
624:   PetscViewerFormat format;

627:   if (incall) return(0);
628:   incall = PETSC_TRUE;
629:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg);
630:   if (flg) {
631:     PetscViewerPushFormat(viewer,format);
632:     NEPValuesView(nep,viewer);
633:     PetscViewerPopFormat(viewer);
634:     PetscViewerDestroy(&viewer);
635:   }
636:   incall = PETSC_FALSE;
637:   return(0);
638: }

640: /*@C
641:    NEPVectorsView - Outputs computed eigenvectors to a viewer.

643:    Collective on NEP

645:    Parameter:
646: +  nep    - the nonlinear eigensolver context
647: -  viewer - the viewer

649:    Options Database Keys:
650: .  -nep_view_vectors - output eigenvectors.

652:    Notes:
653:    If PETSc was configured with real scalars, complex conjugate eigenvectors
654:    will be viewed as two separate real vectors, one containing the real part
655:    and another one containing the imaginary part.

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

661:    Level: intermediate

663: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
664: @*/
665: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
666: {
668:   PetscInt       i,k;
669:   Vec            x;
670: #define NMLEN 30
671:   char           vname[NMLEN];
672:   const char     *ename;

676:   if (!viewer) {
677:     PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
678:   }
681:   NEPCheckSolved(nep,1);
682:   if (nep->nconv) {
683:     PetscObjectGetName((PetscObject)nep,&ename);
684:     NEPComputeVectors(nep);
685:     for (i=0;i<nep->nconv;i++) {
686:       k = nep->perm[i];
687:       PetscSNPrintf(vname,NMLEN,"X%d_%s",(int)i,ename);
688:       BVGetColumn(nep->V,k,&x);
689:       PetscObjectSetName((PetscObject)x,vname);
690:       VecView(x,viewer);
691:       BVRestoreColumn(nep->V,k,&x);
692:       if (nep->twosided) {
693:         PetscSNPrintf(vname,NMLEN,"Y%d_%s",(int)i,ename);
694:         BVGetColumn(nep->W,k,&x);
695:         PetscObjectSetName((PetscObject)x,vname);
696:         VecView(x,viewer);
697:         BVRestoreColumn(nep->W,k,&x);
698:       }
699:     }
700:   }
701:   return(0);
702: }

704: /*@
705:    NEPVectorsViewFromOptions - Processes command line options to determine if/how
706:    the computed eigenvectors are to be viewed.

708:    Collective on NEP

710:    Input Parameters:
711: .  nep - the nonlinear eigensolver context

713:    Level: developer
714: @*/
715: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
716: {
717:   PetscErrorCode    ierr;
718:   PetscViewer       viewer;
719:   PetscBool         flg = PETSC_FALSE;
720:   static PetscBool  incall = PETSC_FALSE;
721:   PetscViewerFormat format;

724:   if (incall) return(0);
725:   incall = PETSC_TRUE;
726:   PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg);
727:   if (flg) {
728:     PetscViewerPushFormat(viewer,format);
729:     NEPVectorsView(nep,viewer);
730:     PetscViewerPopFormat(viewer);
731:     PetscViewerDestroy(&viewer);
732:   }
733:   incall = PETSC_FALSE;
734:   return(0);
735: }