Actual source code: nepview.c

slepc-main 2024-11-15
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: /*@
 19:    NEPView - Prints the NEP data structure.

 21:    Collective

 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;

 54:   PetscFunctionBegin;
 56:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
 58:   PetscCheckSameComm(nep,1,viewer,2);

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

182: /*@
183:    NEPViewFromOptions - View from options

185:    Collective

187:    Input Parameters:
188: +  nep  - the nonlinear eigensolver context
189: .  obj  - optional object
190: -  name - command line option

192:    Level: intermediate

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

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

207:    Collective

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

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

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

222:    Level: intermediate

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

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

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

248:    Collective

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

253:    Level: developer

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

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

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

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

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

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

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

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

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

373:    Collective

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

380:    Options Database Keys:
381: +  -nep_error_absolute - print absolute errors of each eigenpair
382: .  -nep_error_relative - print relative errors of each eigenpair
383: -  -nep_error_backward - print backward errors of each eigenpair

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

391:    Level: intermediate

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

400:   PetscFunctionBegin;
402:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
404:   PetscCheckSameComm(nep,1,viewer,3);
405:   NEPCheckSolved(nep,1);
406:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
407:   if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);

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

427: /*@
428:    NEPErrorViewFromOptions - Processes command line options to determine if/how
429:    the errors of the computed solution are to be viewed.

431:    Collective

433:    Input Parameter:
434: .  nep - the nonlinear eigensolver context

436:    Level: developer

438: .seealso: NEPErrorView()
439: @*/
440: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
441: {
442:   PetscViewer       viewer;
443:   PetscBool         flg;
444:   static PetscBool  incall = PETSC_FALSE;
445:   PetscViewerFormat format;

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

475: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
476: {
477:   PetscDraw      draw;
478:   PetscDrawSP    drawsp;
479:   PetscReal      re,im;
480:   PetscInt       i,k;

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

504: static PetscErrorCode NEPValuesView_BINARY(NEP nep,PetscViewer viewer)
505: {
506: #if defined(PETSC_HAVE_COMPLEX)
507:   PetscInt       i,k;
508:   PetscComplex   *ev;
509: #endif

511:   PetscFunctionBegin;
512: #if defined(PETSC_HAVE_COMPLEX)
513:   PetscCall(PetscMalloc1(nep->nconv,&ev));
514:   for (i=0;i<nep->nconv;i++) {
515:     k = nep->perm[i];
516: #if defined(PETSC_USE_COMPLEX)
517:     ev[i] = nep->eigr[k];
518: #else
519:     ev[i] = PetscCMPLX(nep->eigr[k],nep->eigi[k]);
520: #endif
521:   }
522:   PetscCall(PetscViewerBinaryWrite(viewer,ev,nep->nconv,PETSC_COMPLEX));
523:   PetscCall(PetscFree(ev));
524: #endif
525:   PetscFunctionReturn(PETSC_SUCCESS);
526: }

528: #if defined(PETSC_HAVE_HDF5)
529: static PetscErrorCode NEPValuesView_HDF5(NEP nep,PetscViewer viewer)
530: {
531:   PetscInt       i,k,n,N;
532:   PetscMPIInt    rank;
533:   Vec            v;
534:   char           vname[30];
535:   const char     *ename;

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

574: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
575: {
576:   PetscInt       i,k;

578:   PetscFunctionBegin;
579:   PetscCall(PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n"));
580:   for (i=0;i<nep->nconv;i++) {
581:     k = nep->perm[i];
582:     PetscCall(PetscViewerASCIIPrintf(viewer,"   "));
583:     PetscCall(SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]));
584:     PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
585:   }
586:   PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
587:   PetscFunctionReturn(PETSC_SUCCESS);
588: }

590: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
591: {
592:   PetscInt       i,k;
593:   PetscReal      re,im;
594:   const char     *name;

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

615: /*@
616:    NEPValuesView - Displays the computed eigenvalues in a viewer.

618:    Collective

620:    Input Parameters:
621: +  nep    - the nonlinear eigensolver context
622: -  viewer - the viewer

624:    Options Database Key:
625: .  -nep_view_values - print computed eigenvalues

627:    Level: intermediate

629: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
630: @*/
631: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
632: {
633:   PetscBool         isascii,isdraw,isbinary;
634:   PetscViewerFormat format;
635: #if defined(PETSC_HAVE_HDF5)
636:   PetscBool         ishdf5;
637: #endif

639:   PetscFunctionBegin;
641:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
643:   PetscCheckSameComm(nep,1,viewer,2);
644:   NEPCheckSolved(nep,1);
645:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
646:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
647: #if defined(PETSC_HAVE_HDF5)
648:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
649: #endif
650:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
651:   if (isdraw) PetscCall(NEPValuesView_DRAW(nep,viewer));
652:   else if (isbinary) PetscCall(NEPValuesView_BINARY(nep,viewer));
653: #if defined(PETSC_HAVE_HDF5)
654:   else if (ishdf5) PetscCall(NEPValuesView_HDF5(nep,viewer));
655: #endif
656:   else if (isascii) {
657:     PetscCall(PetscViewerGetFormat(viewer,&format));
658:     switch (format) {
659:       case PETSC_VIEWER_DEFAULT:
660:       case PETSC_VIEWER_ASCII_INFO:
661:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
662:         PetscCall(NEPValuesView_ASCII(nep,viewer));
663:         break;
664:       case PETSC_VIEWER_ASCII_MATLAB:
665:         PetscCall(NEPValuesView_MATLAB(nep,viewer));
666:         break;
667:       default:
668:         PetscCall(PetscInfo(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
669:     }
670:   }
671:   PetscFunctionReturn(PETSC_SUCCESS);
672: }

674: /*@
675:    NEPValuesViewFromOptions - Processes command line options to determine if/how
676:    the computed eigenvalues are to be viewed.

678:    Collective

680:    Input Parameter:
681: .  nep - the nonlinear eigensolver context

683:    Level: developer

685: .seealso: NEPValuesView()
686: @*/
687: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
688: {
689:   PetscViewer       viewer;
690:   PetscBool         flg;
691:   static PetscBool  incall = PETSC_FALSE;
692:   PetscViewerFormat format;

694:   PetscFunctionBegin;
695:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
696:   incall = PETSC_TRUE;
697:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg));
698:   if (flg) {
699:     PetscCall(PetscViewerPushFormat(viewer,format));
700:     PetscCall(NEPValuesView(nep,viewer));
701:     PetscCall(PetscViewerPopFormat(viewer));
702:     PetscCall(PetscViewerDestroy(&viewer));
703:   }
704:   incall = PETSC_FALSE;
705:   PetscFunctionReturn(PETSC_SUCCESS);
706: }

708: /*@
709:    NEPVectorsView - Outputs computed eigenvectors to a viewer.

711:    Collective

713:    Input Parameters:
714: +  nep    - the nonlinear eigensolver context
715: -  viewer - the viewer

717:    Options Database Key:
718: .  -nep_view_vectors - output eigenvectors.

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

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

729:    Level: intermediate

731: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
732: @*/
733: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
734: {
735:   PetscInt       i,k;
736:   Vec            xr,xi=NULL;

738:   PetscFunctionBegin;
740:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
742:   PetscCheckSameComm(nep,1,viewer,2);
743:   NEPCheckSolved(nep,1);
744:   if (nep->nconv) {
745:     PetscCall(NEPComputeVectors(nep));
746:     PetscCall(BVCreateVec(nep->V,&xr));
747: #if !defined(PETSC_USE_COMPLEX)
748:     PetscCall(BVCreateVec(nep->V,&xi));
749: #endif
750:     for (i=0;i<nep->nconv;i++) {
751:       k = nep->perm[i];
752:       PetscCall(BV_GetEigenvector(nep->V,k,nep->eigi[k],xr,xi));
753:       PetscCall(SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)nep));
754:       if (nep->twosided) {
755:         PetscCall(BV_GetEigenvector(nep->W,k,nep->eigi[k],xr,xi));
756:         PetscCall(SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)nep));
757:       }
758:     }
759:     PetscCall(VecDestroy(&xr));
760: #if !defined(PETSC_USE_COMPLEX)
761:     PetscCall(VecDestroy(&xi));
762: #endif
763:   }
764:   PetscFunctionReturn(PETSC_SUCCESS);
765: }

767: /*@
768:    NEPVectorsViewFromOptions - Processes command line options to determine if/how
769:    the computed eigenvectors are to be viewed.

771:    Collective

773:    Input Parameter:
774: .  nep - the nonlinear eigensolver context

776:    Level: developer

778: .seealso: NEPVectorsView()
779: @*/
780: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
781: {
782:   PetscViewer       viewer;
783:   PetscBool         flg = PETSC_FALSE;
784:   static PetscBool  incall = PETSC_FALSE;
785:   PetscViewerFormat format;

787:   PetscFunctionBegin;
788:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
789:   incall = PETSC_TRUE;
790:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg));
791:   if (flg) {
792:     PetscCall(PetscViewerPushFormat(viewer,format));
793:     PetscCall(NEPVectorsView(nep,viewer));
794:     PetscCall(PetscViewerPopFormat(viewer));
795:     PetscCall(PetscViewerDestroy(&viewer));
796:   }
797:   incall = PETSC_FALSE;
798:   PetscFunctionReturn(PETSC_SUCCESS);
799: }