Actual source code: nepview.c

  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: /*@
 19:    NEPView - Prints the `NEP` data structure.

 21:    Collective

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

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

 30:    Notes:
 31:    The available visualization contexts include
 32: +     `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
 33: -     `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard output where only the
 34:          first process opens the file; all other processes send their data to the
 35:          first one to print

 37:    The user can open an alternative visualization context with `PetscViewerASCIIOpen()`
 38:    to output to a specified file.

 40:    Level: beginner

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

 53:   PetscFunctionBegin;
 55:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
 57:   PetscCheckSameComm(nep,1,viewer,2);

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

181: /*@
182:    NEPViewFromOptions - View (print) a `NEP` object based on values in the options database.

184:    Collective

186:    Input Parameters:
187: +  nep  - the nonlinear eigensolver context
188: .  obj  - optional object that provides the options prefix used to query the options database
189: -  name - command line option

191:    Level: intermediate

193: .seealso: [](ch:nep), `NEPView()`, `NEPCreate()`, `PetscObjectViewFromOptions()`
194: @*/
195: PetscErrorCode NEPViewFromOptions(NEP nep,PetscObject obj,const char name[])
196: {
197:   PetscFunctionBegin;
199:   PetscCall(PetscObjectViewFromOptions((PetscObject)nep,obj,name));
200:   PetscFunctionReturn(PETSC_SUCCESS);
201: }

203: /*@
204:    NEPConvergedReasonView - Displays the reason a `NEP` solve converged or diverged.

206:    Collective

208:    Input Parameters:
209: +  nep - the nonlinear eigensolver context
210: -  viewer - the viewer to display the reason

212:    Options Database Key:
213: .  -nep_converged_reason - print reason for convergence/divergence, and number of iterations

215:    Notes:
216:    Use `NEPConvergedReasonViewFromOptions()` to display the reason based on values
217:    in the options database.

219:    To change the format of the output call `PetscViewerPushFormat()` before this
220:    call. Use `PETSC_VIEWER_DEFAULT` for the default, or `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: [](ch:nep), `NEPSetConvergenceTest()`, `NEPSetTolerances()`, `NEPGetIterationNumber()`, `NEPConvergedReasonViewFromOptions()`
227: @*/
228: PetscErrorCode NEPConvergedReasonView(NEP nep,PetscViewer viewer)
229: {
230:   PetscBool         isAscii;
231:   PetscViewerFormat format;

233:   PetscFunctionBegin;
234:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
235:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
236:   if (isAscii) {
237:     PetscCall(PetscViewerGetFormat(viewer,&format));
238:     PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
239:     if (nep->reason > 0 && format != PETSC_VIEWER_FAILED) PetscCall(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));
240:     else if (nep->reason <= 0) PetscCall(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));
241:     PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
242:   }
243:   PetscFunctionReturn(PETSC_SUCCESS);
244: }

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

250:    Collective

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

255:    Level: intermediate

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

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

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

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

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

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

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

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

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

375:    Collective

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

382:    Options Database Keys:
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:    All the command-line options listed above admit an optional argument
394:    specifying the viewer type and options. For instance, use
395:    `-nep_error_relative :myerr.m:ascii_matlab` to save the errors in a file
396:    that can be executed in Matlab.

398:    Level: intermediate

400: .seealso: [](ch:nep), `NEPSolve()`, `NEPValuesView()`, `NEPVectorsView()`
401: @*/
402: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
403: {
404:   PetscBool         isascii;
405:   PetscViewerFormat format;

407:   PetscFunctionBegin;
409:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
411:   PetscCheckSameComm(nep,1,viewer,3);
412:   NEPCheckSolved(nep,1);
413:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
414:   if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);

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

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

438:    Collective

440:    Input Parameter:
441: .  nep - the nonlinear eigensolver context

443:    Level: developer

445: .seealso: [](ch:nep), `NEPErrorView()`
446: @*/
447: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
448: {
449:   PetscViewer       viewer;
450:   PetscBool         flg;
451:   static PetscBool  incall = PETSC_FALSE;
452:   PetscViewerFormat format;

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

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

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

511: static PetscErrorCode NEPValuesView_BINARY(NEP nep,PetscViewer viewer)
512: {
513: #if defined(PETSC_HAVE_COMPLEX)
514:   PetscInt       i,k;
515:   PetscComplex   *ev;
516: #endif

518:   PetscFunctionBegin;
519: #if defined(PETSC_HAVE_COMPLEX)
520:   PetscCall(PetscMalloc1(nep->nconv,&ev));
521:   for (i=0;i<nep->nconv;i++) {
522:     k = nep->perm[i];
523: #if defined(PETSC_USE_COMPLEX)
524:     ev[i] = nep->eigr[k];
525: #else
526:     ev[i] = PetscCMPLX(nep->eigr[k],nep->eigi[k]);
527: #endif
528:   }
529:   PetscCall(PetscViewerBinaryWrite(viewer,ev,nep->nconv,PETSC_COMPLEX));
530:   PetscCall(PetscFree(ev));
531: #endif
532:   PetscFunctionReturn(PETSC_SUCCESS);
533: }

535: #if defined(PETSC_HAVE_HDF5)
536: static PetscErrorCode NEPValuesView_HDF5(NEP nep,PetscViewer viewer)
537: {
538:   PetscInt       i,k,n,N;
539:   PetscMPIInt    rank;
540:   Vec            v;
541:   char           vname[30];
542:   const char     *ename;

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

581: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
582: {
583:   PetscInt       i,k;

585:   PetscFunctionBegin;
586:   PetscCall(PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n"));
587:   for (i=0;i<nep->nconv;i++) {
588:     k = nep->perm[i];
589:     PetscCall(PetscViewerASCIIPrintf(viewer,"   "));
590:     PetscCall(SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]));
591:     PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
592:   }
593:   PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
594:   PetscFunctionReturn(PETSC_SUCCESS);
595: }

597: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
598: {
599:   PetscInt       i,k;
600:   PetscReal      re,im;
601:   const char     *name;

603:   PetscFunctionBegin;
604:   PetscCall(PetscObjectGetName((PetscObject)nep,&name));
605:   PetscCall(PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name));
606:   for (i=0;i<nep->nconv;i++) {
607:     k = nep->perm[i];
608: #if defined(PETSC_USE_COMPLEX)
609:     re = PetscRealPart(nep->eigr[k]);
610:     im = PetscImaginaryPart(nep->eigr[k]);
611: #else
612:     re = nep->eigr[k];
613:     im = nep->eigi[k];
614: #endif
615:     if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im));
616:     else PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re));
617:   }
618:   PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
619:   PetscFunctionReturn(PETSC_SUCCESS);
620: }

622: /*@
623:    NEPValuesView - Displays the computed eigenvalues in a viewer.

625:    Collective

627:    Input Parameters:
628: +  nep    - the nonlinear eigensolver context
629: -  viewer - the viewer

631:    Options Database Key:
632: .  -nep_view_values - print computed eigenvalues

634:    Note:
635:    The command-line option listed above admits an optional argument
636:    specifying the viewer type and options. For instance, use
637:    `-nep_view_values :evals.m:ascii_matlab` to save the values in a file
638:    that can be executed in Matlab.

640:    Level: intermediate

642: .seealso: [](ch:nep), `NEPSolve()`, `NEPVectorsView()`, `NEPErrorView()`
643: @*/
644: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
645: {
646:   PetscBool         isascii,isdraw,isbinary;
647:   PetscViewerFormat format;
648: #if defined(PETSC_HAVE_HDF5)
649:   PetscBool         ishdf5;
650: #endif

652:   PetscFunctionBegin;
654:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
656:   PetscCheckSameComm(nep,1,viewer,2);
657:   NEPCheckSolved(nep,1);
658:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
659:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
660: #if defined(PETSC_HAVE_HDF5)
661:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
662: #endif
663:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
664:   if (isdraw) PetscCall(NEPValuesView_DRAW(nep,viewer));
665:   else if (isbinary) PetscCall(NEPValuesView_BINARY(nep,viewer));
666: #if defined(PETSC_HAVE_HDF5)
667:   else if (ishdf5) PetscCall(NEPValuesView_HDF5(nep,viewer));
668: #endif
669:   else if (isascii) {
670:     PetscCall(PetscViewerGetFormat(viewer,&format));
671:     switch (format) {
672:       case PETSC_VIEWER_DEFAULT:
673:       case PETSC_VIEWER_ASCII_INFO:
674:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
675:         PetscCall(NEPValuesView_ASCII(nep,viewer));
676:         break;
677:       case PETSC_VIEWER_ASCII_MATLAB:
678:         PetscCall(NEPValuesView_MATLAB(nep,viewer));
679:         break;
680:       default:
681:         PetscCall(PetscInfo(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
682:     }
683:   }
684:   PetscFunctionReturn(PETSC_SUCCESS);
685: }

687: /*@
688:    NEPValuesViewFromOptions - Processes command line options to determine if/how
689:    the computed eigenvalues are to be viewed.

691:    Collective

693:    Input Parameter:
694: .  nep - the nonlinear eigensolver context

696:    Level: developer

698: .seealso: [](ch:nep), `NEPValuesView()`
699: @*/
700: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
701: {
702:   PetscViewer       viewer;
703:   PetscBool         flg;
704:   static PetscBool  incall = PETSC_FALSE;
705:   PetscViewerFormat format;

707:   PetscFunctionBegin;
708:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
709:   incall = PETSC_TRUE;
710:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg));
711:   if (flg) {
712:     PetscCall(PetscViewerPushFormat(viewer,format));
713:     PetscCall(NEPValuesView(nep,viewer));
714:     PetscCall(PetscViewerPopFormat(viewer));
715:     PetscCall(PetscViewerDestroy(&viewer));
716:   }
717:   incall = PETSC_FALSE;
718:   PetscFunctionReturn(PETSC_SUCCESS);
719: }

721: /*@
722:    NEPVectorsView - Outputs computed eigenvectors to a viewer.

724:    Collective

726:    Input Parameters:
727: +  nep    - the nonlinear eigensolver context
728: -  viewer - the viewer

730:    Options Database Key:
731: .  -nep_view_vectors - output eigenvectors

733:    Notes:
734:    If PETSc was configured with real scalars, complex conjugate eigenvectors
735:    will be viewed as two separate real vectors, one containing the real part
736:    and another one containing the imaginary part.

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

742:    The command-line option listed above admits an optional argument
743:    specifying the viewer type and options. For instance, use
744:    `-nep_view_vectors binary:evecs.bin` to save the vectors in a binary file.

746:    Level: intermediate

748: .seealso: [](ch:nep), `NEPSolve()`, `NEPValuesView()`, `NEPErrorView()`
749: @*/
750: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
751: {
752:   PetscInt       i,k;
753:   Vec            xr,xi=NULL;

755:   PetscFunctionBegin;
757:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
759:   PetscCheckSameComm(nep,1,viewer,2);
760:   NEPCheckSolved(nep,1);
761:   if (nep->nconv) {
762:     PetscCall(NEPComputeVectors(nep));
763:     PetscCall(BVCreateVec(nep->V,&xr));
764: #if !defined(PETSC_USE_COMPLEX)
765:     PetscCall(BVCreateVec(nep->V,&xi));
766: #endif
767:     for (i=0;i<nep->nconv;i++) {
768:       k = nep->perm[i];
769:       PetscCall(BV_GetEigenvector(nep->V,k,nep->eigi[k],xr,xi));
770:       PetscCall(SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)nep));
771:       if (nep->twosided) {
772:         PetscCall(BV_GetEigenvector(nep->W,k,nep->eigi[k],xr,xi));
773:         PetscCall(SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)nep));
774:       }
775:     }
776:     PetscCall(VecDestroy(&xr));
777: #if !defined(PETSC_USE_COMPLEX)
778:     PetscCall(VecDestroy(&xi));
779: #endif
780:   }
781:   PetscFunctionReturn(PETSC_SUCCESS);
782: }

784: /*@
785:    NEPVectorsViewFromOptions - Processes command line options to determine if/how
786:    the computed eigenvectors are to be viewed.

788:    Collective

790:    Input Parameter:
791: .  nep - the nonlinear eigensolver context

793:    Level: developer

795: .seealso: [](ch:nep), `NEPVectorsView()`
796: @*/
797: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
798: {
799:   PetscViewer       viewer;
800:   PetscBool         flg = PETSC_FALSE;
801:   static PetscBool  incall = PETSC_FALSE;
802:   PetscViewerFormat format;

804:   PetscFunctionBegin;
805:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
806:   incall = PETSC_TRUE;
807:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg));
808:   if (flg) {
809:     PetscCall(PetscViewerPushFormat(viewer,format));
810:     PetscCall(NEPVectorsView(nep,viewer));
811:     PetscCall(PetscViewerPopFormat(viewer));
812:     PetscCall(PetscViewerDestroy(&viewer));
813:   }
814:   incall = PETSC_FALSE;
815:   PetscFunctionReturn(PETSC_SUCCESS);
816: }