Actual source code: nepview.c
slepc-3.6.1 2015-09-03
1: /*
2: The NEP routines related to various viewers.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2015, Universitat Politecnica de Valencia, Spain
8: This file is part of SLEPc.
10: SLEPc is free software: you can redistribute it and/or modify it under the
11: terms of version 3 of the GNU Lesser General Public License as published by
12: the Free Software Foundation.
14: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
15: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
17: more details.
19: You should have received a copy of the GNU Lesser General Public License
20: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
21: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
22: */
24: #include <slepc/private/nepimpl.h> /*I "slepcnep.h" I*/
25: #include <petscdraw.h>
29: /*@C
30: NEPView - Prints the NEP data structure.
32: Collective on NEP
34: Input Parameters:
35: + nep - the nonlinear eigenproblem solver context
36: - viewer - optional visualization context
38: Options Database Key:
39: . -nep_view - Calls NEPView() at end of NEPSolve()
41: Note:
42: The available visualization contexts include
43: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
44: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
45: output where only the first processor opens
46: the file. All other processors send their
47: data to the first processor to print.
49: The user can open an alternative visualization context with
50: PetscViewerASCIIOpen() - output to a specified file.
52: Level: beginner
54: .seealso: PetscViewerASCIIOpen()
55: @*/
56: PetscErrorCode NEPView(NEP nep,PetscViewer viewer)
57: {
59: char str[50];
60: PetscBool isascii,isslp,istrivial,nods;
64: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
68: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
69: if (isascii) {
70: PetscObjectPrintClassNamePrefixType((PetscObject)nep,viewer);
71: if (nep->ops->view) {
72: PetscViewerASCIIPushTab(viewer);
73: (*nep->ops->view)(nep,viewer);
74: PetscViewerASCIIPopTab(viewer);
75: }
76: if (nep->split) {
77: PetscViewerASCIIPrintf(viewer," nonlinear operator in split form\n");
78: } else {
79: PetscViewerASCIIPrintf(viewer," nonlinear operator from user callbacks\n");
80: }
81: PetscViewerASCIIPrintf(viewer," iterative refinement: %s\n",NEPRefineTypes[nep->refine]);
82: if (nep->refine) {
83: PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%D\n",(double)nep->reftol,nep->rits);
84: }
85: if (nep->npart>1) {
86: PetscViewerASCIIPrintf(viewer," splitting communicator in %D partitions for refinement\n",nep->npart);
87: }
88: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
89: SlepcSNPrintfScalar(str,50,nep->target,PETSC_FALSE);
90: if (!nep->which) {
91: PetscViewerASCIIPrintf(viewer,"not yet set\n");
92: } else switch (nep->which) {
93: case NEP_TARGET_MAGNITUDE:
94: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
95: break;
96: case NEP_TARGET_REAL:
97: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
98: break;
99: case NEP_TARGET_IMAGINARY:
100: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
101: break;
102: case NEP_LARGEST_MAGNITUDE:
103: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
104: break;
105: case NEP_SMALLEST_MAGNITUDE:
106: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
107: break;
108: case NEP_LARGEST_REAL:
109: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
110: break;
111: case NEP_SMALLEST_REAL:
112: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
113: break;
114: case NEP_LARGEST_IMAGINARY:
115: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
116: break;
117: case NEP_SMALLEST_IMAGINARY:
118: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
119: break;
120: default: SETERRQ(PetscObjectComm((PetscObject)nep),1,"Wrong value of nep->which");
121: }
122: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",nep->nev);
123: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",nep->ncv);
124: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",nep->mpd);
125: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",nep->max_it);
126: PetscViewerASCIIPrintf(viewer," maximum number of function evaluations: %D\n",nep->max_funcs);
127: PetscViewerASCIIPrintf(viewer," tolerances: relative=%g, absolute=%g, solution=%g\n",(double)nep->rtol,(double)nep->abstol,(double)nep->stol);
128: if (nep->lag) {
129: PetscViewerASCIIPrintf(viewer," updating the preconditioner every %D iterations\n",nep->lag);
130: }
131: if (nep->cctol) {
132: PetscViewerASCIIPrintf(viewer," using a constant tolerance for the linear solver\n");
133: }
134: if (nep->nini) {
135: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(nep->nini));
136: }
137: } else {
138: if (nep->ops->view) {
139: (*nep->ops->view)(nep,viewer);
140: }
141: }
142: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
143: if (!nep->V) { NEPGetBV(nep,&nep->V); }
144: BVView(nep->V,viewer);
145: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
146: RGIsTrivial(nep->rg,&istrivial);
147: if (!istrivial) { RGView(nep->rg,viewer); }
148: PetscObjectTypeCompareAny((PetscObject)nep,&nods,NEPRII,NEPSLP,NEPINTERPOL,"");
149: if (!nods) {
150: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
151: DSView(nep->ds,viewer);
152: }
153: PetscViewerPopFormat(viewer);
154: PetscObjectTypeCompare((PetscObject)nep,NEPSLP,&isslp);
155: if (!isslp) {
156: if (!nep->ksp) { NEPGetKSP(nep,&nep->ksp); }
157: KSPView(nep->ksp,viewer);
158: }
159: return(0);
160: }
164: /*@C
165: NEPReasonView - Displays the reason a NEP solve converged or diverged.
167: Collective on NEP
169: Parameter:
170: + nep - the nonlinear eigensolver context
171: - viewer - the viewer to display the reason
173: Options Database Keys:
174: . -nep_converged_reason - print reason for convergence, and number of iterations
176: Level: intermediate
178: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber()
179: @*/
180: PetscErrorCode NEPReasonView(NEP nep,PetscViewer viewer)
181: {
183: PetscBool isAscii;
186: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
187: if (isAscii) {
188: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
189: if (nep->reason > 0) {
190: 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);
191: } else {
192: 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);
193: }
194: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
195: }
196: return(0);
197: }
201: /*@
202: NEPReasonViewFromOptions - Processes command line options to determine if/how
203: the NEP converged reason is to be viewed.
205: Collective on NEP
207: Input Parameters:
208: . nep - the nonlinear eigensolver context
210: Level: developer
211: @*/
212: PetscErrorCode NEPReasonViewFromOptions(NEP nep)
213: {
214: PetscErrorCode ierr;
215: PetscViewer viewer;
216: PetscBool flg;
217: static PetscBool incall = PETSC_FALSE;
218: PetscViewerFormat format;
221: if (incall) return(0);
222: incall = PETSC_TRUE;
223: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg);
224: if (flg) {
225: PetscViewerPushFormat(viewer,format);
226: NEPReasonView(nep,viewer);
227: PetscViewerPopFormat(viewer);
228: PetscViewerDestroy(&viewer);
229: }
230: incall = PETSC_FALSE;
231: return(0);
232: }
236: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
237: {
238: PetscBool errok;
239: PetscReal error,re,im;
240: PetscScalar kr,ki;
241: PetscInt i,j;
245: if (nep->nconv<nep->nev) {
246: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",nep->nev);
247: return(0);
248: }
249: errok = PETSC_TRUE;
250: for (i=0;i<nep->nev;i++) {
251: NEPComputeError(nep,i,etype,&error);
252: errok = (errok && error<5.0*nep->rtol)? PETSC_TRUE: PETSC_FALSE;
253: }
254: if (!errok) {
255: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nep->nev);
256: return(0);
257: }
258: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
259: for (i=0;i<=(nep->nev-1)/8;i++) {
260: PetscViewerASCIIPrintf(viewer,"\n ");
261: for (j=0;j<PetscMin(8,nep->nev-8*i);j++) {
262: NEPGetEigenpair(nep,8*i+j,&kr,&ki,NULL,NULL);
263: #if defined(PETSC_USE_COMPLEX)
264: re = PetscRealPart(kr);
265: im = PetscImaginaryPart(kr);
266: #else
267: re = kr;
268: im = ki;
269: #endif
270: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
271: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
272: if (im!=0.0) {
273: PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
274: } else {
275: PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
276: }
277: if (8*i+j+1<nep->nev) { PetscViewerASCIIPrintf(viewer,", "); }
278: }
279: }
280: PetscViewerASCIIPrintf(viewer,"\n\n");
281: return(0);
282: }
286: static PetscErrorCode NEPErrorView_DETAIL(NEP nep,NEPErrorType etype,PetscViewer viewer)
287: {
289: PetscReal error,re,im;
290: PetscScalar kr,ki;
291: PetscInt i;
292: #define EXLEN 30
293: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
296: if (!nep->nconv) return(0);
297: switch (etype) {
298: case NEP_ERROR_ABSOLUTE:
299: PetscSNPrintf(ex,EXLEN," ||T(k)x||");
300: break;
301: case NEP_ERROR_RELATIVE:
302: PetscSNPrintf(ex,EXLEN," ||T(k)x||/||kx||");
303: break;
304: }
305: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
306: for (i=0;i<nep->nconv;i++) {
307: NEPGetEigenpair(nep,i,&kr,&ki,NULL,NULL);
308: NEPComputeError(nep,i,etype,&error);
309: #if defined(PETSC_USE_COMPLEX)
310: re = PetscRealPart(kr);
311: im = PetscImaginaryPart(kr);
312: #else
313: re = kr;
314: im = ki;
315: #endif
316: if (im!=0.0) {
317: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
318: } else {
319: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
320: }
321: }
322: PetscViewerASCIIPrintf(viewer,sep);
323: return(0);
324: }
328: static PetscErrorCode NEPErrorView_MATLAB(NEP nep,NEPErrorType etype,PetscViewer viewer)
329: {
331: PetscReal error;
332: PetscInt i;
333: const char *name;
336: PetscObjectGetName((PetscObject)nep,&name);
337: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
338: for (i=0;i<nep->nconv;i++) {
339: NEPComputeError(nep,i,etype,&error);
340: PetscViewerASCIIPrintf(viewer,"%18.16e\n",error);
341: }
342: PetscViewerASCIIPrintf(viewer,"];\n");
343: return(0);
344: }
348: /*@C
349: NEPErrorView - Displays the errors associated with the computed solution
350: (as well as the eigenvalues).
352: Collective on NEP
354: Input Parameters:
355: + nep - the nonlinear eigensolver context
356: . etype - error type
357: - viewer - optional visualization context
359: Options Database Key:
360: + -nep_error_absolute - print absolute errors of each eigenpair
361: - -nep_error_relative - print relative errors of each eigenpair
363: Notes:
364: By default, this function checks the error of all eigenpairs and prints
365: the eigenvalues if all of them are below the requested tolerance.
366: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
367: eigenvalues and corresponding errors is printed.
369: Level: intermediate
371: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
372: @*/
373: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
374: {
375: PetscBool isascii;
376: PetscViewerFormat format;
377: PetscErrorCode ierr;
381: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
384: NEPCheckSolved(nep,1);
385: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
386: if (!isascii) return(0);
388: PetscViewerGetFormat(viewer,&format);
389: switch (format) {
390: case PETSC_VIEWER_DEFAULT:
391: case PETSC_VIEWER_ASCII_INFO:
392: NEPErrorView_ASCII(nep,etype,viewer);
393: break;
394: case PETSC_VIEWER_ASCII_INFO_DETAIL:
395: NEPErrorView_DETAIL(nep,etype,viewer);
396: break;
397: case PETSC_VIEWER_ASCII_MATLAB:
398: NEPErrorView_MATLAB(nep,etype,viewer);
399: break;
400: default:
401: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
402: }
403: return(0);
404: }
408: /*@
409: NEPErrorViewFromOptions - Processes command line options to determine if/how
410: the errors of the computed solution are to be viewed.
412: Collective on NEP
414: Input Parameters:
415: . nep - the nonlinear eigensolver context
417: Level: developer
418: @*/
419: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
420: {
421: PetscErrorCode ierr;
422: PetscViewer viewer;
423: PetscBool flg;
424: static PetscBool incall = PETSC_FALSE;
425: PetscViewerFormat format;
428: if (incall) return(0);
429: incall = PETSC_TRUE;
430: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg);
431: if (flg) {
432: PetscViewerPushFormat(viewer,format);
433: NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer);
434: PetscViewerPopFormat(viewer);
435: PetscViewerDestroy(&viewer);
436: }
437: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg);
438: if (flg) {
439: PetscViewerPushFormat(viewer,format);
440: NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer);
441: PetscViewerPopFormat(viewer);
442: PetscViewerDestroy(&viewer);
443: }
444: incall = PETSC_FALSE;
445: return(0);
446: }
450: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
451: {
453: PetscDraw draw;
454: PetscDrawSP drawsp;
455: PetscReal re,im;
456: PetscInt i,k;
459: if (!nep->nconv) return(0);
460: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Computed Eigenvalues",PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
461: PetscViewerDrawGetDraw(viewer,0,&draw);
462: PetscDrawSPCreate(draw,1,&drawsp);
463: for (i=0;i<nep->nconv;i++) {
464: k = nep->perm[i];
465: #if defined(PETSC_USE_COMPLEX)
466: re = PetscRealPart(nep->eigr[k]);
467: im = PetscImaginaryPart(nep->eigr[k]);
468: #else
469: re = nep->eigr[k];
470: im = nep->eigi[k];
471: #endif
472: PetscDrawSPAddPoint(drawsp,&re,&im);
473: }
474: PetscDrawSPDraw(drawsp,PETSC_TRUE);
475: PetscDrawSPDestroy(&drawsp);
476: PetscViewerDestroy(&viewer);
477: return(0);
478: }
482: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
483: {
484: PetscReal re,im;
485: PetscInt i,k;
489: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
490: for (i=0;i<nep->nconv;i++) {
491: k = nep->perm[i];
492: #if defined(PETSC_USE_COMPLEX)
493: re = PetscRealPart(nep->eigr[k]);
494: im = PetscImaginaryPart(nep->eigr[k]);
495: #else
496: re = nep->eigr[k];
497: im = nep->eigi[k];
498: #endif
499: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
500: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
501: if (im!=0.0) {
502: PetscViewerASCIIPrintf(viewer," %.5f%+.5fi\n",(double)re,(double)im);
503: } else {
504: PetscViewerASCIIPrintf(viewer," %.5f\n",(double)re);
505: }
506: }
507: PetscViewerASCIIPrintf(viewer,"\n");
508: return(0);
509: }
513: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
514: {
516: PetscInt i,k;
517: PetscReal re,im;
518: const char *name;
521: PetscObjectGetName((PetscObject)nep,&name);
522: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
523: for (i=0;i<nep->nconv;i++) {
524: k = nep->perm[i];
525: #if defined(PETSC_USE_COMPLEX)
526: re = PetscRealPart(nep->eigr[k]);
527: im = PetscImaginaryPart(nep->eigr[k]);
528: #else
529: re = nep->eigr[k];
530: im = nep->eigi[k];
531: #endif
532: if (im!=0.0) {
533: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
534: } else {
535: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
536: }
537: }
538: PetscViewerASCIIPrintf(viewer,"];\n");
539: return(0);
540: }
544: /*@C
545: NEPValuesView - Displays the computed eigenvalues in a viewer.
547: Collective on NEP
549: Input Parameters:
550: + nep - the nonlinear eigensolver context
551: - viewer - the viewer
553: Options Database Key:
554: . -nep_view_values - print computed eigenvalues
556: Level: intermediate
558: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
559: @*/
560: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
561: {
562: PetscBool isascii,isdraw;
563: PetscViewerFormat format;
564: PetscErrorCode ierr;
568: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
571: NEPCheckSolved(nep,1);
572: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
573: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
574: if (isdraw) {
575: NEPValuesView_DRAW(nep,viewer);
576: } else if (isascii) {
577: PetscViewerGetFormat(viewer,&format);
578: switch (format) {
579: case PETSC_VIEWER_DEFAULT:
580: case PETSC_VIEWER_ASCII_INFO:
581: case PETSC_VIEWER_ASCII_INFO_DETAIL:
582: NEPValuesView_ASCII(nep,viewer);
583: break;
584: case PETSC_VIEWER_ASCII_MATLAB:
585: NEPValuesView_MATLAB(nep,viewer);
586: break;
587: default:
588: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
589: }
590: }
591: return(0);
592: }
596: /*@
597: NEPValuesViewFromOptions - Processes command line options to determine if/how
598: the computed eigenvalues are to be viewed.
600: Collective on NEP
602: Input Parameters:
603: . nep - the nonlinear eigensolver context
605: Level: developer
606: @*/
607: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
608: {
609: PetscErrorCode ierr;
610: PetscViewer viewer;
611: PetscBool flg;
612: static PetscBool incall = PETSC_FALSE;
613: PetscViewerFormat format;
616: if (incall) return(0);
617: incall = PETSC_TRUE;
618: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg);
619: if (flg) {
620: PetscViewerPushFormat(viewer,format);
621: NEPValuesView(nep,viewer);
622: PetscViewerPopFormat(viewer);
623: PetscViewerDestroy(&viewer);
624: }
625: incall = PETSC_FALSE;
626: return(0);
627: }
631: /*@C
632: NEPVectorsView - Outputs computed eigenvectors to a viewer.
634: Collective on NEP
636: Parameter:
637: + nep - the nonlinear eigensolver context
638: - viewer - the viewer
640: Options Database Keys:
641: . -nep_view_vectors - output eigenvectors.
643: Note:
644: If PETSc was configured with real scalars, complex conjugate eigenvectors
645: will be viewed as two separate real vectors, one containing the real part
646: and another one containing the imaginary part.
648: Level: intermediate
650: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
651: @*/
652: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
653: {
655: PetscInt i,k;
656: Vec x;
657: #define NMLEN 30
658: char vname[NMLEN];
659: const char *ename;
663: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
666: NEPCheckSolved(nep,1);
667: if (nep->nconv) {
668: PetscObjectGetName((PetscObject)nep,&ename);
669: NEPComputeVectors(nep);
670: for (i=0;i<nep->nconv;i++) {
671: k = nep->perm[i];
672: PetscSNPrintf(vname,NMLEN,"V%d_%s",i,ename);
673: BVGetColumn(nep->V,k,&x);
674: PetscObjectSetName((PetscObject)x,vname);
675: VecView(x,viewer);
676: BVRestoreColumn(nep->V,k,&x);
677: }
678: }
679: return(0);
680: }
684: /*@
685: NEPVectorsViewFromOptions - Processes command line options to determine if/how
686: the computed eigenvectors are to be viewed.
688: Collective on NEP
690: Input Parameters:
691: . nep - the nonlinear eigensolver context
693: Level: developer
694: @*/
695: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
696: {
697: PetscErrorCode ierr;
698: PetscViewer viewer;
699: PetscBool flg = PETSC_FALSE;
700: static PetscBool incall = PETSC_FALSE;
701: PetscViewerFormat format;
704: if (incall) return(0);
705: incall = PETSC_TRUE;
706: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg);
707: if (flg) {
708: PetscViewerPushFormat(viewer,format);
709: NEPVectorsView(nep,viewer);
710: PetscViewerPopFormat(viewer);
711: PetscViewerDestroy(&viewer);
712: }
713: incall = PETSC_FALSE;
714: return(0);
715: }