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: }