Actual source code: nepview.c
slepc-3.12.1 2019-11-08
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2019, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: NEP routines related to various viewers
12: */
14: #include <slepc/private/nepimpl.h> /*I "slepcnep.h" I*/
15: #include <petscdraw.h>
17: /*@C
18: NEPView - Prints the NEP data structure.
20: Collective on nep
22: Input Parameters:
23: + nep - the nonlinear eigenproblem solver context
24: - viewer - optional visualization context
26: Options Database Key:
27: . -nep_view - Calls NEPView() at end of NEPSolve()
29: Note:
30: The available visualization contexts include
31: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
32: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
33: output where only the first processor opens
34: the file. All other processors send their
35: data to the first processor to print.
37: The user can open an alternative visualization context with
38: PetscViewerASCIIOpen() - output to a specified file.
40: Level: beginner
42: .seealso: PetscViewerASCIIOpen()
43: @*/
44: PetscErrorCode NEPView(NEP nep,PetscViewer viewer)
45: {
47: const char *type=NULL;
48: char str[50];
49: PetscInt i;
50: PetscBool isascii,istrivial;
51: PetscViewer sviewer;
55: if (!viewer) {
56: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
57: }
61: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
62: if (isascii) {
63: PetscObjectPrintClassNamePrefixType((PetscObject)nep,viewer);
64: if (nep->ops->view) {
65: PetscViewerASCIIPushTab(viewer);
66: (*nep->ops->view)(nep,viewer);
67: PetscViewerASCIIPopTab(viewer);
68: }
69: if (nep->problem_type) {
70: switch (nep->problem_type) {
71: case NEP_GENERAL: type = "general nonlinear eigenvalue problem"; break;
72: case NEP_RATIONAL: type = "rational eigenvalue problem"; break;
73: }
74: } else type = "not yet set";
75: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
76: if (nep->fui) {
77: switch (nep->fui) {
78: case NEP_USER_INTERFACE_CALLBACK:
79: PetscViewerASCIIPrintf(viewer," nonlinear operator from user callbacks\n");
80: break;
81: case NEP_USER_INTERFACE_SPLIT:
82: PetscViewerASCIIPrintf(viewer," nonlinear operator in split form, with %D terms\n",nep->nt);
83: break;
84: }
85: } else {
86: PetscViewerASCIIPrintf(viewer," nonlinear operator not specified yet\n");
87: }
88: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
89: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
90: SlepcSNPrintfScalar(str,50,nep->target,PETSC_FALSE);
91: if (!nep->which) {
92: PetscViewerASCIIPrintf(viewer,"not yet set\n");
93: } else switch (nep->which) {
94: case NEP_WHICH_USER:
95: PetscViewerASCIIPrintf(viewer,"user defined\n");
96: break;
97: case NEP_TARGET_MAGNITUDE:
98: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
99: break;
100: case NEP_TARGET_REAL:
101: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
102: break;
103: case NEP_TARGET_IMAGINARY:
104: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
105: break;
106: case NEP_LARGEST_MAGNITUDE:
107: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
108: break;
109: case NEP_SMALLEST_MAGNITUDE:
110: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
111: break;
112: case NEP_LARGEST_REAL:
113: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
114: break;
115: case NEP_SMALLEST_REAL:
116: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
117: break;
118: case NEP_LARGEST_IMAGINARY:
119: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
120: break;
121: case NEP_SMALLEST_IMAGINARY:
122: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
123: break;
124: case NEP_ALL:
125: PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
126: break;
127: }
128: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
129: if (nep->twosided) {
130: PetscViewerASCIIPrintf(viewer," using two-sided variant (for left eigenvectors)\n");
131: }
132: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",nep->nev);
133: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",nep->ncv);
134: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",nep->mpd);
135: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",nep->max_it);
136: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)nep->tol);
137: PetscViewerASCIIPrintf(viewer," convergence test: ");
138: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
139: switch (nep->conv) {
140: case NEP_CONV_ABS:
141: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
142: case NEP_CONV_REL:
143: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
144: case NEP_CONV_NORM:
145: PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
146: if (nep->nrma) {
147: PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)nep->nrma[0]);
148: for (i=1;i<nep->nt;i++) {
149: PetscViewerASCIIPrintf(viewer,", %g",(double)nep->nrma[i]);
150: }
151: PetscViewerASCIIPrintf(viewer,"\n");
152: }
153: break;
154: case NEP_CONV_USER:
155: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
156: }
157: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
158: if (nep->refine) {
159: PetscViewerASCIIPrintf(viewer," iterative refinement: %s, with %s scheme\n",NEPRefineTypes[nep->refine],NEPRefineSchemes[nep->scheme]);
160: PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%D\n",(double)nep->rtol,nep->rits);
161: if (nep->npart>1) {
162: PetscViewerASCIIPrintf(viewer," splitting communicator in %D partitions for refinement\n",nep->npart);
163: }
164: }
165: if (nep->nini) {
166: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(nep->nini));
167: }
168: } else {
169: if (nep->ops->view) {
170: (*nep->ops->view)(nep,viewer);
171: }
172: }
173: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
174: if (!nep->V) { NEPGetBV(nep,&nep->V); }
175: BVView(nep->V,viewer);
176: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
177: RGIsTrivial(nep->rg,&istrivial);
178: if (!istrivial) { RGView(nep->rg,viewer); }
179: if (nep->useds) {
180: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
181: DSView(nep->ds,viewer);
182: }
183: PetscViewerPopFormat(viewer);
184: if (nep->refine!=NEP_REFINE_NONE) {
185: PetscViewerASCIIPushTab(viewer);
186: if (nep->npart>1) {
187: if (nep->refinesubc->color==0) {
188: PetscViewerASCIIGetStdout(PetscSubcommChild(nep->refinesubc),&sviewer);
189: KSPView(nep->refineksp,sviewer);
190: }
191: } else {
192: KSPView(nep->refineksp,viewer);
193: }
194: PetscViewerASCIIPopTab(viewer);
195: }
196: return(0);
197: }
199: /*@C
200: NEPReasonView - Displays the reason a NEP solve converged or diverged.
202: Collective on nep
204: Parameter:
205: + nep - the nonlinear eigensolver context
206: - viewer - the viewer to display the reason
208: Options Database Keys:
209: . -nep_converged_reason - print reason for convergence, and number of iterations
211: Level: intermediate
213: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber()
214: @*/
215: PetscErrorCode NEPReasonView(NEP nep,PetscViewer viewer)
216: {
218: PetscBool isAscii;
221: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
222: if (isAscii) {
223: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
224: if (nep->reason > 0) {
225: PetscViewerASCIIPrintf(viewer,"%s Nonlinear eigensolve converged (%D eigenpair%s) due to %s; iterations %D\n",((PetscObject)nep)->prefix?((PetscObject)nep)->prefix:"",nep->nconv,(nep->nconv>1)?"s":"",NEPConvergedReasons[nep->reason],nep->its);
226: } else {
227: PetscViewerASCIIPrintf(viewer,"%s Nonlinear eigensolve did not converge due to %s; iterations %D\n",((PetscObject)nep)->prefix?((PetscObject)nep)->prefix:"",NEPConvergedReasons[nep->reason],nep->its);
228: }
229: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
230: }
231: return(0);
232: }
234: /*@
235: NEPReasonViewFromOptions - Processes command line options to determine if/how
236: the NEP converged reason is to be viewed.
238: Collective on nep
240: Input Parameters:
241: . nep - the nonlinear eigensolver context
243: Level: developer
244: @*/
245: PetscErrorCode NEPReasonViewFromOptions(NEP nep)
246: {
247: PetscErrorCode ierr;
248: PetscViewer viewer;
249: PetscBool flg;
250: static PetscBool incall = PETSC_FALSE;
251: PetscViewerFormat format;
254: if (incall) return(0);
255: incall = PETSC_TRUE;
256: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg);
257: if (flg) {
258: PetscViewerPushFormat(viewer,format);
259: NEPReasonView(nep,viewer);
260: PetscViewerPopFormat(viewer);
261: PetscViewerDestroy(&viewer);
262: }
263: incall = PETSC_FALSE;
264: return(0);
265: }
267: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
268: {
269: PetscReal error;
270: PetscInt i,j,k,nvals;
274: nvals = (nep->which==NEP_ALL)? nep->nconv: nep->nev;
275: if (nep->which!=NEP_ALL && nep->nconv<nvals) {
276: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",nep->nev);
277: return(0);
278: }
279: if (nep->which==NEP_ALL && !nvals) {
280: PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
281: return(0);
282: }
283: for (i=0;i<nvals;i++) {
284: NEPComputeError(nep,i,etype,&error);
285: if (error>=5.0*nep->tol) {
286: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
287: return(0);
288: }
289: }
290: if (nep->which==NEP_ALL) {
291: PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
292: } else {
293: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
294: }
295: for (i=0;i<=(nvals-1)/8;i++) {
296: PetscViewerASCIIPrintf(viewer,"\n ");
297: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
298: k = nep->perm[8*i+j];
299: SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]);
300: if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
301: }
302: }
303: PetscViewerASCIIPrintf(viewer,"\n\n");
304: return(0);
305: }
307: static PetscErrorCode NEPErrorView_DETAIL(NEP nep,NEPErrorType etype,PetscViewer viewer)
308: {
310: PetscReal error,re,im;
311: PetscScalar kr,ki;
312: PetscInt i;
313: #define EXLEN 30
314: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
317: if (!nep->nconv) return(0);
318: switch (etype) {
319: case NEP_ERROR_ABSOLUTE:
320: PetscSNPrintf(ex,EXLEN," ||T(k)x||");
321: break;
322: case NEP_ERROR_RELATIVE:
323: PetscSNPrintf(ex,EXLEN," ||T(k)x||/||kx||");
324: break;
325: case NEP_ERROR_BACKWARD:
326: PetscSNPrintf(ex,EXLEN," eta(x,k)");
327: break;
328: }
329: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
330: for (i=0;i<nep->nconv;i++) {
331: NEPGetEigenpair(nep,i,&kr,&ki,NULL,NULL);
332: NEPComputeError(nep,i,etype,&error);
333: #if defined(PETSC_USE_COMPLEX)
334: re = PetscRealPart(kr);
335: im = PetscImaginaryPart(kr);
336: #else
337: re = kr;
338: im = ki;
339: #endif
340: if (im!=0.0) {
341: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
342: } else {
343: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
344: }
345: }
346: PetscViewerASCIIPrintf(viewer,sep);
347: return(0);
348: }
350: static PetscErrorCode NEPErrorView_MATLAB(NEP nep,NEPErrorType etype,PetscViewer viewer)
351: {
353: PetscReal error;
354: PetscInt i;
355: const char *name;
358: PetscObjectGetName((PetscObject)nep,&name);
359: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
360: for (i=0;i<nep->nconv;i++) {
361: NEPComputeError(nep,i,etype,&error);
362: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
363: }
364: PetscViewerASCIIPrintf(viewer,"];\n");
365: return(0);
366: }
368: /*@C
369: NEPErrorView - Displays the errors associated with the computed solution
370: (as well as the eigenvalues).
372: Collective on nep
374: Input Parameters:
375: + nep - the nonlinear eigensolver context
376: . etype - error type
377: - viewer - optional visualization context
379: Options Database Key:
380: + -nep_error_absolute - print absolute errors of each eigenpair
381: . -nep_error_relative - print relative errors of each eigenpair
382: - -nep_error_backward - print backward errors of each eigenpair
384: Notes:
385: By default, this function checks the error of all eigenpairs and prints
386: the eigenvalues if all of them are below the requested tolerance.
387: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
388: eigenvalues and corresponding errors is printed.
390: Level: intermediate
392: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
393: @*/
394: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
395: {
396: PetscBool isascii;
397: PetscViewerFormat format;
398: PetscErrorCode ierr;
402: if (!viewer) {
403: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
404: }
407: NEPCheckSolved(nep,1);
408: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
409: if (!isascii) return(0);
411: PetscViewerGetFormat(viewer,&format);
412: switch (format) {
413: case PETSC_VIEWER_DEFAULT:
414: case PETSC_VIEWER_ASCII_INFO:
415: NEPErrorView_ASCII(nep,etype,viewer);
416: break;
417: case PETSC_VIEWER_ASCII_INFO_DETAIL:
418: NEPErrorView_DETAIL(nep,etype,viewer);
419: break;
420: case PETSC_VIEWER_ASCII_MATLAB:
421: NEPErrorView_MATLAB(nep,etype,viewer);
422: break;
423: default:
424: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
425: }
426: return(0);
427: }
429: /*@
430: NEPErrorViewFromOptions - Processes command line options to determine if/how
431: the errors of the computed solution are to be viewed.
433: Collective on nep
435: Input Parameters:
436: . nep - the nonlinear eigensolver context
438: Level: developer
439: @*/
440: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
441: {
442: PetscErrorCode ierr;
443: PetscViewer viewer;
444: PetscBool flg;
445: static PetscBool incall = PETSC_FALSE;
446: PetscViewerFormat format;
449: if (incall) return(0);
450: incall = PETSC_TRUE;
451: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg);
452: if (flg) {
453: PetscViewerPushFormat(viewer,format);
454: NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer);
455: PetscViewerPopFormat(viewer);
456: PetscViewerDestroy(&viewer);
457: }
458: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg);
459: if (flg) {
460: PetscViewerPushFormat(viewer,format);
461: NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer);
462: PetscViewerPopFormat(viewer);
463: PetscViewerDestroy(&viewer);
464: }
465: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_backward",&viewer,&format,&flg);
466: if (flg) {
467: PetscViewerPushFormat(viewer,format);
468: NEPErrorView(nep,NEP_ERROR_BACKWARD,viewer);
469: PetscViewerPopFormat(viewer);
470: PetscViewerDestroy(&viewer);
471: }
472: incall = PETSC_FALSE;
473: return(0);
474: }
476: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
477: {
479: PetscDraw draw;
480: PetscDrawSP drawsp;
481: PetscReal re,im;
482: PetscInt i,k;
485: if (!nep->nconv) return(0);
486: PetscViewerDrawGetDraw(viewer,0,&draw);
487: PetscDrawSetTitle(draw,"Computed Eigenvalues");
488: PetscDrawSPCreate(draw,1,&drawsp);
489: for (i=0;i<nep->nconv;i++) {
490: k = nep->perm[i];
491: #if defined(PETSC_USE_COMPLEX)
492: re = PetscRealPart(nep->eigr[k]);
493: im = PetscImaginaryPart(nep->eigr[k]);
494: #else
495: re = nep->eigr[k];
496: im = nep->eigi[k];
497: #endif
498: PetscDrawSPAddPoint(drawsp,&re,&im);
499: }
500: PetscDrawSPDraw(drawsp,PETSC_TRUE);
501: PetscDrawSPSave(drawsp);
502: PetscDrawSPDestroy(&drawsp);
503: return(0);
504: }
506: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
507: {
508: PetscInt i,k;
512: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
513: for (i=0;i<nep->nconv;i++) {
514: k = nep->perm[i];
515: PetscViewerASCIIPrintf(viewer," ");
516: SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]);
517: PetscViewerASCIIPrintf(viewer,"\n");
518: }
519: PetscViewerASCIIPrintf(viewer,"\n");
520: return(0);
521: }
523: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
524: {
526: PetscInt i,k;
527: PetscReal re,im;
528: const char *name;
531: PetscObjectGetName((PetscObject)nep,&name);
532: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
533: for (i=0;i<nep->nconv;i++) {
534: k = nep->perm[i];
535: #if defined(PETSC_USE_COMPLEX)
536: re = PetscRealPart(nep->eigr[k]);
537: im = PetscImaginaryPart(nep->eigr[k]);
538: #else
539: re = nep->eigr[k];
540: im = nep->eigi[k];
541: #endif
542: if (im!=0.0) {
543: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
544: } else {
545: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
546: }
547: }
548: PetscViewerASCIIPrintf(viewer,"];\n");
549: return(0);
550: }
552: /*@C
553: NEPValuesView - Displays the computed eigenvalues in a viewer.
555: Collective on nep
557: Input Parameters:
558: + nep - the nonlinear eigensolver context
559: - viewer - the viewer
561: Options Database Key:
562: . -nep_view_values - print computed eigenvalues
564: Level: intermediate
566: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
567: @*/
568: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
569: {
570: PetscBool isascii,isdraw;
571: PetscViewerFormat format;
572: PetscErrorCode ierr;
576: if (!viewer) {
577: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
578: }
581: NEPCheckSolved(nep,1);
582: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
583: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
584: if (isdraw) {
585: NEPValuesView_DRAW(nep,viewer);
586: } else if (isascii) {
587: PetscViewerGetFormat(viewer,&format);
588: switch (format) {
589: case PETSC_VIEWER_DEFAULT:
590: case PETSC_VIEWER_ASCII_INFO:
591: case PETSC_VIEWER_ASCII_INFO_DETAIL:
592: NEPValuesView_ASCII(nep,viewer);
593: break;
594: case PETSC_VIEWER_ASCII_MATLAB:
595: NEPValuesView_MATLAB(nep,viewer);
596: break;
597: default:
598: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
599: }
600: }
601: return(0);
602: }
604: /*@
605: NEPValuesViewFromOptions - Processes command line options to determine if/how
606: the computed eigenvalues are to be viewed.
608: Collective on nep
610: Input Parameters:
611: . nep - the nonlinear eigensolver context
613: Level: developer
614: @*/
615: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
616: {
617: PetscErrorCode ierr;
618: PetscViewer viewer;
619: PetscBool flg;
620: static PetscBool incall = PETSC_FALSE;
621: PetscViewerFormat format;
624: if (incall) return(0);
625: incall = PETSC_TRUE;
626: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg);
627: if (flg) {
628: PetscViewerPushFormat(viewer,format);
629: NEPValuesView(nep,viewer);
630: PetscViewerPopFormat(viewer);
631: PetscViewerDestroy(&viewer);
632: }
633: incall = PETSC_FALSE;
634: return(0);
635: }
637: /*@C
638: NEPVectorsView - Outputs computed eigenvectors to a viewer.
640: Collective on nep
642: Parameter:
643: + nep - the nonlinear eigensolver context
644: - viewer - the viewer
646: Options Database Keys:
647: . -nep_view_vectors - output eigenvectors.
649: Notes:
650: If PETSc was configured with real scalars, complex conjugate eigenvectors
651: will be viewed as two separate real vectors, one containing the real part
652: and another one containing the imaginary part.
654: If left eigenvectors were computed with a two-sided eigensolver, the right
655: and left eigenvectors are interleaved, that is, the vectors are output in
656: the following order: X0, Y0, X1, Y1, X2, Y2, ...
658: Level: intermediate
660: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
661: @*/
662: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
663: {
665: PetscInt i,k;
666: Vec x;
667: #define NMLEN 30
668: char vname[NMLEN];
669: const char *ename;
673: if (!viewer) {
674: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
675: }
678: NEPCheckSolved(nep,1);
679: if (nep->nconv) {
680: PetscObjectGetName((PetscObject)nep,&ename);
681: NEPComputeVectors(nep);
682: for (i=0;i<nep->nconv;i++) {
683: k = nep->perm[i];
684: PetscSNPrintf(vname,NMLEN,"X%d_%s",(int)i,ename);
685: BVGetColumn(nep->V,k,&x);
686: PetscObjectSetName((PetscObject)x,vname);
687: VecView(x,viewer);
688: BVRestoreColumn(nep->V,k,&x);
689: if (nep->twosided) {
690: PetscSNPrintf(vname,NMLEN,"Y%d_%s",(int)i,ename);
691: BVGetColumn(nep->W,k,&x);
692: PetscObjectSetName((PetscObject)x,vname);
693: VecView(x,viewer);
694: BVRestoreColumn(nep->W,k,&x);
695: }
696: }
697: }
698: return(0);
699: }
701: /*@
702: NEPVectorsViewFromOptions - Processes command line options to determine if/how
703: the computed eigenvectors are to be viewed.
705: Collective on nep
707: Input Parameters:
708: . nep - the nonlinear eigensolver context
710: Level: developer
711: @*/
712: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
713: {
714: PetscErrorCode ierr;
715: PetscViewer viewer;
716: PetscBool flg = PETSC_FALSE;
717: static PetscBool incall = PETSC_FALSE;
718: PetscViewerFormat format;
721: if (incall) return(0);
722: incall = PETSC_TRUE;
723: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg);
724: if (flg) {
725: PetscViewerPushFormat(viewer,format);
726: NEPVectorsView(nep,viewer);
727: PetscViewerPopFormat(viewer);
728: PetscViewerDestroy(&viewer);
729: }
730: incall = PETSC_FALSE;
731: return(0);
732: }