Actual source code: nepview.c
slepc-3.7.4 2017-05-17
1: /*
2: The NEP routines related to various viewers.
4: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
5: SLEPc - Scalable Library for Eigenvalue Problem Computations
6: Copyright (c) 2002-2016, 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: PetscInt i;
61: PetscBool isascii,istrivial,nods;
65: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
69: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
70: if (isascii) {
71: PetscObjectPrintClassNamePrefixType((PetscObject)nep,viewer);
72: if (nep->ops->view) {
73: PetscViewerASCIIPushTab(viewer);
74: (*nep->ops->view)(nep,viewer);
75: PetscViewerASCIIPopTab(viewer);
76: }
77: if (nep->fui) {
78: switch (nep->fui) {
79: case NEP_USER_INTERFACE_CALLBACK:
80: PetscViewerASCIIPrintf(viewer," nonlinear operator from user callbacks\n");
81: break;
82: case NEP_USER_INTERFACE_SPLIT:
83: PetscViewerASCIIPrintf(viewer," nonlinear operator in split form\n");
84: break;
85: case NEP_USER_INTERFACE_DERIVATIVES:
86: PetscViewerASCIIPrintf(viewer," nonlinear operator from user callbacks for the successive derivatives\n");
87: break;
88: }
89: } else {
90: PetscViewerASCIIPrintf(viewer," nonlinear operator not specified yet\n");
91: }
92: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
93: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
94: SlepcSNPrintfScalar(str,50,nep->target,PETSC_FALSE);
95: if (!nep->which) {
96: PetscViewerASCIIPrintf(viewer,"not yet set\n");
97: } else switch (nep->which) {
98: case NEP_WHICH_USER:
99: PetscViewerASCIIPrintf(viewer,"user defined\n");
100: break;
101: case NEP_TARGET_MAGNITUDE:
102: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
103: break;
104: case NEP_TARGET_REAL:
105: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
106: break;
107: case NEP_TARGET_IMAGINARY:
108: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
109: break;
110: case NEP_LARGEST_MAGNITUDE:
111: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
112: break;
113: case NEP_SMALLEST_MAGNITUDE:
114: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
115: break;
116: case NEP_LARGEST_REAL:
117: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
118: break;
119: case NEP_SMALLEST_REAL:
120: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
121: break;
122: case NEP_LARGEST_IMAGINARY:
123: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
124: break;
125: case NEP_SMALLEST_IMAGINARY:
126: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
127: break;
128: case NEP_ALL:
129: PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
130: break;
131: default: SETERRQ(PetscObjectComm((PetscObject)nep),1,"Wrong value of nep->which");
132: }
133: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
134: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",nep->nev);
135: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",nep->ncv);
136: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",nep->mpd);
137: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",nep->max_it);
138: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)nep->tol);
139: PetscViewerASCIIPrintf(viewer," convergence test: ");
140: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
141: switch (nep->conv) {
142: case NEP_CONV_ABS:
143: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
144: case NEP_CONV_REL:
145: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
146: case NEP_CONV_NORM:
147: PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
148: if (nep->nrma) {
149: PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)nep->nrma[0]);
150: for (i=1;i<nep->nt;i++) {
151: PetscViewerASCIIPrintf(viewer,", %g",(double)nep->nrma[i]);
152: }
153: PetscViewerASCIIPrintf(viewer,"\n");
154: }
155: break;
156: case NEP_CONV_USER:
157: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
158: }
159: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
160: if (nep->refine) {
161: PetscViewerASCIIPrintf(viewer," iterative refinement: %s, with %s scheme\n",NEPRefineTypes[nep->refine],NEPRefineSchemes[nep->scheme]);
162: PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%D\n",(double)nep->rtol,nep->rits);
163: if (nep->npart>1) {
164: PetscViewerASCIIPrintf(viewer," splitting communicator in %D partitions for refinement\n",nep->npart);
165: }
166: }
167: if (nep->nini) {
168: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(nep->nini));
169: }
170: } else {
171: if (nep->ops->view) {
172: (*nep->ops->view)(nep,viewer);
173: }
174: }
175: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
176: if (!nep->V) { NEPGetBV(nep,&nep->V); }
177: BVView(nep->V,viewer);
178: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
179: RGIsTrivial(nep->rg,&istrivial);
180: if (!istrivial) { RGView(nep->rg,viewer); }
181: PetscObjectTypeCompareAny((PetscObject)nep,&nods,NEPRII,NEPSLP,NEPINTERPOL,"");
182: if (!nods) {
183: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
184: DSView(nep->ds,viewer);
185: }
186: PetscViewerPopFormat(viewer);
187: return(0);
188: }
192: /*@C
193: NEPReasonView - Displays the reason a NEP solve converged or diverged.
195: Collective on NEP
197: Parameter:
198: + nep - the nonlinear eigensolver context
199: - viewer - the viewer to display the reason
201: Options Database Keys:
202: . -nep_converged_reason - print reason for convergence, and number of iterations
204: Level: intermediate
206: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber()
207: @*/
208: PetscErrorCode NEPReasonView(NEP nep,PetscViewer viewer)
209: {
211: PetscBool isAscii;
214: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
215: if (isAscii) {
216: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
217: if (nep->reason > 0) {
218: 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);
219: } else {
220: 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);
221: }
222: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
223: }
224: return(0);
225: }
229: /*@
230: NEPReasonViewFromOptions - Processes command line options to determine if/how
231: the NEP converged reason is to be viewed.
233: Collective on NEP
235: Input Parameters:
236: . nep - the nonlinear eigensolver context
238: Level: developer
239: @*/
240: PetscErrorCode NEPReasonViewFromOptions(NEP nep)
241: {
242: PetscErrorCode ierr;
243: PetscViewer viewer;
244: PetscBool flg;
245: static PetscBool incall = PETSC_FALSE;
246: PetscViewerFormat format;
249: if (incall) return(0);
250: incall = PETSC_TRUE;
251: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg);
252: if (flg) {
253: PetscViewerPushFormat(viewer,format);
254: NEPReasonView(nep,viewer);
255: PetscViewerPopFormat(viewer);
256: PetscViewerDestroy(&viewer);
257: }
258: incall = PETSC_FALSE;
259: return(0);
260: }
264: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
265: {
266: PetscBool errok;
267: PetscReal error,re,im;
268: PetscScalar kr,ki;
269: PetscInt i,j,nvals;
273: if (nep->which!=NEP_ALL && nep->nconv<nep->nev) {
274: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",nep->nev);
275: return(0);
276: }
277: errok = PETSC_TRUE;
278: nvals = (nep->which==NEP_ALL)? nep->nconv: nep->nev;
279: for (i=0;i<nvals;i++) {
280: NEPComputeError(nep,i,etype,&error);
281: errok = (errok && error<5.0*nep->tol)? PETSC_TRUE: PETSC_FALSE;
282: }
283: if (!errok) {
284: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
285: return(0);
286: }
287: if (nep->which==NEP_ALL) {
288: PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
289: } else {
290: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
291: }
292: for (i=0;i<=(nvals-1)/8;i++) {
293: PetscViewerASCIIPrintf(viewer,"\n ");
294: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
295: NEPGetEigenpair(nep,8*i+j,&kr,&ki,NULL,NULL);
296: #if defined(PETSC_USE_COMPLEX)
297: re = PetscRealPart(kr);
298: im = PetscImaginaryPart(kr);
299: #else
300: re = kr;
301: im = ki;
302: #endif
303: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
304: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
305: if (im!=0.0) {
306: PetscViewerASCIIPrintf(viewer,"%.5f%+.5fi",(double)re,(double)im);
307: } else {
308: PetscViewerASCIIPrintf(viewer,"%.5f",(double)re);
309: }
310: if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
311: }
312: }
313: PetscViewerASCIIPrintf(viewer,"\n\n");
314: return(0);
315: }
319: static PetscErrorCode NEPErrorView_DETAIL(NEP nep,NEPErrorType etype,PetscViewer viewer)
320: {
322: PetscReal error,re,im;
323: PetscScalar kr,ki;
324: PetscInt i;
325: #define EXLEN 30
326: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
329: if (!nep->nconv) return(0);
330: switch (etype) {
331: case NEP_ERROR_ABSOLUTE:
332: PetscSNPrintf(ex,EXLEN," ||T(k)x||");
333: break;
334: case NEP_ERROR_RELATIVE:
335: PetscSNPrintf(ex,EXLEN," ||T(k)x||/||kx||");
336: break;
337: case NEP_ERROR_BACKWARD:
338: PetscSNPrintf(ex,EXLEN," eta(x,k)");
339: break;
340: }
341: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
342: for (i=0;i<nep->nconv;i++) {
343: NEPGetEigenpair(nep,i,&kr,&ki,NULL,NULL);
344: NEPComputeError(nep,i,etype,&error);
345: #if defined(PETSC_USE_COMPLEX)
346: re = PetscRealPart(kr);
347: im = PetscImaginaryPart(kr);
348: #else
349: re = kr;
350: im = ki;
351: #endif
352: if (im!=0.0) {
353: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
354: } else {
355: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
356: }
357: }
358: PetscViewerASCIIPrintf(viewer,sep);
359: return(0);
360: }
364: static PetscErrorCode NEPErrorView_MATLAB(NEP nep,NEPErrorType etype,PetscViewer viewer)
365: {
367: PetscReal error;
368: PetscInt i;
369: const char *name;
372: PetscObjectGetName((PetscObject)nep,&name);
373: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
374: for (i=0;i<nep->nconv;i++) {
375: NEPComputeError(nep,i,etype,&error);
376: PetscViewerASCIIPrintf(viewer,"%18.16e\n",error);
377: }
378: PetscViewerASCIIPrintf(viewer,"];\n");
379: return(0);
380: }
384: /*@C
385: NEPErrorView - Displays the errors associated with the computed solution
386: (as well as the eigenvalues).
388: Collective on NEP
390: Input Parameters:
391: + nep - the nonlinear eigensolver context
392: . etype - error type
393: - viewer - optional visualization context
395: Options Database Key:
396: + -nep_error_absolute - print absolute errors of each eigenpair
397: . -nep_error_relative - print relative errors of each eigenpair
398: - -nep_error_backward - print backward errors of each eigenpair
400: Notes:
401: By default, this function checks the error of all eigenpairs and prints
402: the eigenvalues if all of them are below the requested tolerance.
403: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
404: eigenvalues and corresponding errors is printed.
406: Level: intermediate
408: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
409: @*/
410: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
411: {
412: PetscBool isascii;
413: PetscViewerFormat format;
414: PetscErrorCode ierr;
418: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
421: NEPCheckSolved(nep,1);
422: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
423: if (!isascii) return(0);
425: PetscViewerGetFormat(viewer,&format);
426: switch (format) {
427: case PETSC_VIEWER_DEFAULT:
428: case PETSC_VIEWER_ASCII_INFO:
429: NEPErrorView_ASCII(nep,etype,viewer);
430: break;
431: case PETSC_VIEWER_ASCII_INFO_DETAIL:
432: NEPErrorView_DETAIL(nep,etype,viewer);
433: break;
434: case PETSC_VIEWER_ASCII_MATLAB:
435: NEPErrorView_MATLAB(nep,etype,viewer);
436: break;
437: default:
438: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
439: }
440: return(0);
441: }
445: /*@
446: NEPErrorViewFromOptions - Processes command line options to determine if/how
447: the errors of the computed solution are to be viewed.
449: Collective on NEP
451: Input Parameters:
452: . nep - the nonlinear eigensolver context
454: Level: developer
455: @*/
456: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
457: {
458: PetscErrorCode ierr;
459: PetscViewer viewer;
460: PetscBool flg;
461: static PetscBool incall = PETSC_FALSE;
462: PetscViewerFormat format;
465: if (incall) return(0);
466: incall = PETSC_TRUE;
467: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg);
468: if (flg) {
469: PetscViewerPushFormat(viewer,format);
470: NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer);
471: PetscViewerPopFormat(viewer);
472: PetscViewerDestroy(&viewer);
473: }
474: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg);
475: if (flg) {
476: PetscViewerPushFormat(viewer,format);
477: NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer);
478: PetscViewerPopFormat(viewer);
479: PetscViewerDestroy(&viewer);
480: }
481: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_error_backward",&viewer,&format,&flg);
482: if (flg) {
483: PetscViewerPushFormat(viewer,format);
484: NEPErrorView(nep,NEP_ERROR_BACKWARD,viewer);
485: PetscViewerPopFormat(viewer);
486: PetscViewerDestroy(&viewer);
487: }
488: incall = PETSC_FALSE;
489: return(0);
490: }
494: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
495: {
497: PetscDraw draw;
498: PetscDrawSP drawsp;
499: PetscReal re,im;
500: PetscInt i,k;
503: if (!nep->nconv) return(0);
504: PetscViewerDrawOpen(PETSC_COMM_SELF,0,"Computed Eigenvalues",PETSC_DECIDE,PETSC_DECIDE,300,300,&viewer);
505: PetscViewerDrawGetDraw(viewer,0,&draw);
506: PetscDrawSPCreate(draw,1,&drawsp);
507: for (i=0;i<nep->nconv;i++) {
508: k = nep->perm[i];
509: #if defined(PETSC_USE_COMPLEX)
510: re = PetscRealPart(nep->eigr[k]);
511: im = PetscImaginaryPart(nep->eigr[k]);
512: #else
513: re = nep->eigr[k];
514: im = nep->eigi[k];
515: #endif
516: PetscDrawSPAddPoint(drawsp,&re,&im);
517: }
518: PetscDrawSPDraw(drawsp,PETSC_TRUE);
519: PetscDrawSPSave(drawsp);
520: PetscDrawSPDestroy(&drawsp);
521: PetscViewerDestroy(&viewer);
522: return(0);
523: }
527: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
528: {
529: PetscReal re,im;
530: PetscInt i,k;
534: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
535: for (i=0;i<nep->nconv;i++) {
536: k = nep->perm[i];
537: #if defined(PETSC_USE_COMPLEX)
538: re = PetscRealPart(nep->eigr[k]);
539: im = PetscImaginaryPart(nep->eigr[k]);
540: #else
541: re = nep->eigr[k];
542: im = nep->eigi[k];
543: #endif
544: if (PetscAbs(re)/PetscAbs(im)<PETSC_SMALL) re = 0.0;
545: if (PetscAbs(im)/PetscAbs(re)<PETSC_SMALL) im = 0.0;
546: if (im!=0.0) {
547: PetscViewerASCIIPrintf(viewer," %.5f%+.5fi\n",(double)re,(double)im);
548: } else {
549: PetscViewerASCIIPrintf(viewer," %.5f\n",(double)re);
550: }
551: }
552: PetscViewerASCIIPrintf(viewer,"\n");
553: return(0);
554: }
558: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
559: {
561: PetscInt i,k;
562: PetscReal re,im;
563: const char *name;
566: PetscObjectGetName((PetscObject)nep,&name);
567: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
568: for (i=0;i<nep->nconv;i++) {
569: k = nep->perm[i];
570: #if defined(PETSC_USE_COMPLEX)
571: re = PetscRealPart(nep->eigr[k]);
572: im = PetscImaginaryPart(nep->eigr[k]);
573: #else
574: re = nep->eigr[k];
575: im = nep->eigi[k];
576: #endif
577: if (im!=0.0) {
578: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
579: } else {
580: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
581: }
582: }
583: PetscViewerASCIIPrintf(viewer,"];\n");
584: return(0);
585: }
589: /*@C
590: NEPValuesView - Displays the computed eigenvalues in a viewer.
592: Collective on NEP
594: Input Parameters:
595: + nep - the nonlinear eigensolver context
596: - viewer - the viewer
598: Options Database Key:
599: . -nep_view_values - print computed eigenvalues
601: Level: intermediate
603: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
604: @*/
605: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
606: {
607: PetscBool isascii,isdraw;
608: PetscViewerFormat format;
609: PetscErrorCode ierr;
613: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
616: NEPCheckSolved(nep,1);
617: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
618: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
619: if (isdraw) {
620: NEPValuesView_DRAW(nep,viewer);
621: } else if (isascii) {
622: PetscViewerGetFormat(viewer,&format);
623: switch (format) {
624: case PETSC_VIEWER_DEFAULT:
625: case PETSC_VIEWER_ASCII_INFO:
626: case PETSC_VIEWER_ASCII_INFO_DETAIL:
627: NEPValuesView_ASCII(nep,viewer);
628: break;
629: case PETSC_VIEWER_ASCII_MATLAB:
630: NEPValuesView_MATLAB(nep,viewer);
631: break;
632: default:
633: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
634: }
635: }
636: return(0);
637: }
641: /*@
642: NEPValuesViewFromOptions - Processes command line options to determine if/how
643: the computed eigenvalues are to be viewed.
645: Collective on NEP
647: Input Parameters:
648: . nep - the nonlinear eigensolver context
650: Level: developer
651: @*/
652: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
653: {
654: PetscErrorCode ierr;
655: PetscViewer viewer;
656: PetscBool flg;
657: static PetscBool incall = PETSC_FALSE;
658: PetscViewerFormat format;
661: if (incall) return(0);
662: incall = PETSC_TRUE;
663: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg);
664: if (flg) {
665: PetscViewerPushFormat(viewer,format);
666: NEPValuesView(nep,viewer);
667: PetscViewerPopFormat(viewer);
668: PetscViewerDestroy(&viewer);
669: }
670: incall = PETSC_FALSE;
671: return(0);
672: }
676: /*@C
677: NEPVectorsView - Outputs computed eigenvectors to a viewer.
679: Collective on NEP
681: Parameter:
682: + nep - the nonlinear eigensolver context
683: - viewer - the viewer
685: Options Database Keys:
686: . -nep_view_vectors - output eigenvectors.
688: Note:
689: If PETSc was configured with real scalars, complex conjugate eigenvectors
690: will be viewed as two separate real vectors, one containing the real part
691: and another one containing the imaginary part.
693: Level: intermediate
695: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
696: @*/
697: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
698: {
700: PetscInt i,k;
701: Vec x;
702: #define NMLEN 30
703: char vname[NMLEN];
704: const char *ename;
708: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
711: NEPCheckSolved(nep,1);
712: if (nep->nconv) {
713: PetscObjectGetName((PetscObject)nep,&ename);
714: NEPComputeVectors(nep);
715: for (i=0;i<nep->nconv;i++) {
716: k = nep->perm[i];
717: PetscSNPrintf(vname,NMLEN,"V%d_%s",(int)i,ename);
718: BVGetColumn(nep->V,k,&x);
719: PetscObjectSetName((PetscObject)x,vname);
720: VecView(x,viewer);
721: BVRestoreColumn(nep->V,k,&x);
722: }
723: }
724: return(0);
725: }
729: /*@
730: NEPVectorsViewFromOptions - Processes command line options to determine if/how
731: the computed eigenvectors are to be viewed.
733: Collective on NEP
735: Input Parameters:
736: . nep - the nonlinear eigensolver context
738: Level: developer
739: @*/
740: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
741: {
742: PetscErrorCode ierr;
743: PetscViewer viewer;
744: PetscBool flg = PETSC_FALSE;
745: static PetscBool incall = PETSC_FALSE;
746: PetscViewerFormat format;
749: if (incall) return(0);
750: incall = PETSC_TRUE;
751: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg);
752: if (flg) {
753: PetscViewerPushFormat(viewer,format);
754: NEPVectorsView(nep,viewer);
755: PetscViewerPopFormat(viewer);
756: PetscViewerDestroy(&viewer);
757: }
758: incall = PETSC_FALSE;
759: return(0);
760: }