Actual source code: nepview.c
slepc-3.13.3 2020-06-14
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, 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 <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
41: @*/
42: PetscErrorCode NEPView(NEP nep,PetscViewer viewer)
43: {
45: const char *type=NULL;
46: char str[50];
47: PetscInt i;
48: PetscBool isascii,istrivial;
49: PetscViewer sviewer;
53: if (!viewer) {
54: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
55: }
59: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
60: if (isascii) {
61: PetscObjectPrintClassNamePrefixType((PetscObject)nep,viewer);
62: if (nep->ops->view) {
63: PetscViewerASCIIPushTab(viewer);
64: (*nep->ops->view)(nep,viewer);
65: PetscViewerASCIIPopTab(viewer);
66: }
67: if (nep->problem_type) {
68: switch (nep->problem_type) {
69: case NEP_GENERAL: type = "general nonlinear eigenvalue problem"; break;
70: case NEP_RATIONAL: type = "rational eigenvalue problem"; break;
71: }
72: } else type = "not yet set";
73: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
74: if (nep->fui) {
75: switch (nep->fui) {
76: case NEP_USER_INTERFACE_CALLBACK:
77: PetscViewerASCIIPrintf(viewer," nonlinear operator from user callbacks\n");
78: break;
79: case NEP_USER_INTERFACE_SPLIT:
80: PetscViewerASCIIPrintf(viewer," nonlinear operator in split form, with %D terms\n",nep->nt);
81: break;
82: }
83: } else {
84: PetscViewerASCIIPrintf(viewer," nonlinear operator not specified yet\n");
85: }
86: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
87: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
88: SlepcSNPrintfScalar(str,50,nep->target,PETSC_FALSE);
89: if (!nep->which) {
90: PetscViewerASCIIPrintf(viewer,"not yet set\n");
91: } else switch (nep->which) {
92: case NEP_WHICH_USER:
93: PetscViewerASCIIPrintf(viewer,"user defined\n");
94: break;
95: case NEP_TARGET_MAGNITUDE:
96: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
97: break;
98: case NEP_TARGET_REAL:
99: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
100: break;
101: case NEP_TARGET_IMAGINARY:
102: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
103: break;
104: case NEP_LARGEST_MAGNITUDE:
105: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
106: break;
107: case NEP_SMALLEST_MAGNITUDE:
108: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
109: break;
110: case NEP_LARGEST_REAL:
111: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
112: break;
113: case NEP_SMALLEST_REAL:
114: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
115: break;
116: case NEP_LARGEST_IMAGINARY:
117: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
118: break;
119: case NEP_SMALLEST_IMAGINARY:
120: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
121: break;
122: case NEP_ALL:
123: PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
124: break;
125: }
126: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
127: if (nep->twosided) {
128: PetscViewerASCIIPrintf(viewer," using two-sided variant (for left eigenvectors)\n");
129: }
130: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",nep->nev);
131: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",nep->ncv);
132: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",nep->mpd);
133: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",nep->max_it);
134: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)nep->tol);
135: PetscViewerASCIIPrintf(viewer," convergence test: ");
136: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
137: switch (nep->conv) {
138: case NEP_CONV_ABS:
139: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
140: case NEP_CONV_REL:
141: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
142: case NEP_CONV_NORM:
143: PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
144: if (nep->nrma) {
145: PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)nep->nrma[0]);
146: for (i=1;i<nep->nt;i++) {
147: PetscViewerASCIIPrintf(viewer,", %g",(double)nep->nrma[i]);
148: }
149: PetscViewerASCIIPrintf(viewer,"\n");
150: }
151: break;
152: case NEP_CONV_USER:
153: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
154: }
155: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
156: if (nep->refine) {
157: PetscViewerASCIIPrintf(viewer," iterative refinement: %s, with %s scheme\n",NEPRefineTypes[nep->refine],NEPRefineSchemes[nep->scheme]);
158: PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%D\n",(double)nep->rtol,nep->rits);
159: if (nep->npart>1) {
160: PetscViewerASCIIPrintf(viewer," splitting communicator in %D partitions for refinement\n",nep->npart);
161: }
162: }
163: if (nep->nini) {
164: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(nep->nini));
165: }
166: } else {
167: if (nep->ops->view) {
168: (*nep->ops->view)(nep,viewer);
169: }
170: }
171: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
172: if (!nep->V) { NEPGetBV(nep,&nep->V); }
173: BVView(nep->V,viewer);
174: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
175: RGIsTrivial(nep->rg,&istrivial);
176: if (!istrivial) { RGView(nep->rg,viewer); }
177: if (nep->useds) {
178: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
179: DSView(nep->ds,viewer);
180: }
181: PetscViewerPopFormat(viewer);
182: if (nep->refine!=NEP_REFINE_NONE) {
183: PetscViewerASCIIPushTab(viewer);
184: if (nep->npart>1) {
185: if (nep->refinesubc->color==0) {
186: PetscViewerASCIIGetStdout(PetscSubcommChild(nep->refinesubc),&sviewer);
187: KSPView(nep->refineksp,sviewer);
188: }
189: } else {
190: KSPView(nep->refineksp,viewer);
191: }
192: PetscViewerASCIIPopTab(viewer);
193: }
194: return(0);
195: }
197: /*@C
198: NEPViewFromOptions - View from options
200: Collective on NEP
202: Input Parameters:
203: + nep - the nonlinear eigensolver context
204: . obj - optional object
205: - name - command line option
207: Level: intermediate
209: .seealso: NEPView(), NEPCreate()
210: @*/
211: PetscErrorCode NEPViewFromOptions(NEP nep,PetscObject obj,const char name[])
212: {
217: PetscObjectViewFromOptions((PetscObject)nep,obj,name);
218: return(0);
219: }
221: /*@C
222: NEPReasonView - Displays the reason a NEP solve converged or diverged.
224: Collective on nep
226: Input Parameters:
227: + nep - the nonlinear eigensolver context
228: - viewer - the viewer to display the reason
230: Options Database Keys:
231: . -nep_converged_reason - print reason for convergence, and number of iterations
233: Level: intermediate
235: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber()
236: @*/
237: PetscErrorCode NEPReasonView(NEP nep,PetscViewer viewer)
238: {
240: PetscBool isAscii;
243: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
244: if (isAscii) {
245: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
246: if (nep->reason > 0) {
247: 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);
248: } else {
249: 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);
250: }
251: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
252: }
253: return(0);
254: }
256: /*@
257: NEPReasonViewFromOptions - Processes command line options to determine if/how
258: the NEP converged reason is to be viewed.
260: Collective on nep
262: Input Parameter:
263: . nep - the nonlinear eigensolver context
265: Level: developer
266: @*/
267: PetscErrorCode NEPReasonViewFromOptions(NEP nep)
268: {
269: PetscErrorCode ierr;
270: PetscViewer viewer;
271: PetscBool flg;
272: static PetscBool incall = PETSC_FALSE;
273: PetscViewerFormat format;
276: if (incall) return(0);
277: incall = PETSC_TRUE;
278: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg);
279: if (flg) {
280: PetscViewerPushFormat(viewer,format);
281: NEPReasonView(nep,viewer);
282: PetscViewerPopFormat(viewer);
283: PetscViewerDestroy(&viewer);
284: }
285: incall = PETSC_FALSE;
286: return(0);
287: }
289: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
290: {
291: PetscReal error;
292: PetscInt i,j,k,nvals;
296: nvals = (nep->which==NEP_ALL)? nep->nconv: nep->nev;
297: if (nep->which!=NEP_ALL && nep->nconv<nvals) {
298: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",nep->nev);
299: return(0);
300: }
301: if (nep->which==NEP_ALL && !nvals) {
302: PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
303: return(0);
304: }
305: for (i=0;i<nvals;i++) {
306: NEPComputeError(nep,i,etype,&error);
307: if (error>=5.0*nep->tol) {
308: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
309: return(0);
310: }
311: }
312: if (nep->which==NEP_ALL) {
313: PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
314: } else {
315: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
316: }
317: for (i=0;i<=(nvals-1)/8;i++) {
318: PetscViewerASCIIPrintf(viewer,"\n ");
319: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
320: k = nep->perm[8*i+j];
321: SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]);
322: if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
323: }
324: }
325: PetscViewerASCIIPrintf(viewer,"\n\n");
326: return(0);
327: }
329: static PetscErrorCode NEPErrorView_DETAIL(NEP nep,NEPErrorType etype,PetscViewer viewer)
330: {
332: PetscReal error,re,im;
333: PetscScalar kr,ki;
334: PetscInt i;
335: #define EXLEN 30
336: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
339: if (!nep->nconv) return(0);
340: switch (etype) {
341: case NEP_ERROR_ABSOLUTE:
342: PetscSNPrintf(ex,EXLEN," ||T(k)x||");
343: break;
344: case NEP_ERROR_RELATIVE:
345: PetscSNPrintf(ex,EXLEN," ||T(k)x||/||kx||");
346: break;
347: case NEP_ERROR_BACKWARD:
348: PetscSNPrintf(ex,EXLEN," eta(x,k)");
349: break;
350: }
351: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
352: for (i=0;i<nep->nconv;i++) {
353: NEPGetEigenpair(nep,i,&kr,&ki,NULL,NULL);
354: NEPComputeError(nep,i,etype,&error);
355: #if defined(PETSC_USE_COMPLEX)
356: re = PetscRealPart(kr);
357: im = PetscImaginaryPart(kr);
358: #else
359: re = kr;
360: im = ki;
361: #endif
362: if (im!=0.0) {
363: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
364: } else {
365: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
366: }
367: }
368: PetscViewerASCIIPrintf(viewer,sep);
369: return(0);
370: }
372: static PetscErrorCode NEPErrorView_MATLAB(NEP nep,NEPErrorType etype,PetscViewer viewer)
373: {
375: PetscReal error;
376: PetscInt i;
377: const char *name;
380: PetscObjectGetName((PetscObject)nep,&name);
381: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
382: for (i=0;i<nep->nconv;i++) {
383: NEPComputeError(nep,i,etype,&error);
384: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
385: }
386: PetscViewerASCIIPrintf(viewer,"];\n");
387: return(0);
388: }
390: /*@C
391: NEPErrorView - Displays the errors associated with the computed solution
392: (as well as the eigenvalues).
394: Collective on nep
396: Input Parameters:
397: + nep - the nonlinear eigensolver context
398: . etype - error type
399: - viewer - optional visualization context
401: Options Database Key:
402: + -nep_error_absolute - print absolute errors of each eigenpair
403: . -nep_error_relative - print relative errors of each eigenpair
404: - -nep_error_backward - print backward errors of each eigenpair
406: Notes:
407: By default, this function checks the error of all eigenpairs and prints
408: the eigenvalues if all of them are below the requested tolerance.
409: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
410: eigenvalues and corresponding errors is printed.
412: Level: intermediate
414: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
415: @*/
416: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
417: {
418: PetscBool isascii;
419: PetscViewerFormat format;
420: PetscErrorCode ierr;
424: if (!viewer) {
425: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
426: }
429: NEPCheckSolved(nep,1);
430: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
431: if (!isascii) return(0);
433: PetscViewerGetFormat(viewer,&format);
434: switch (format) {
435: case PETSC_VIEWER_DEFAULT:
436: case PETSC_VIEWER_ASCII_INFO:
437: NEPErrorView_ASCII(nep,etype,viewer);
438: break;
439: case PETSC_VIEWER_ASCII_INFO_DETAIL:
440: NEPErrorView_DETAIL(nep,etype,viewer);
441: break;
442: case PETSC_VIEWER_ASCII_MATLAB:
443: NEPErrorView_MATLAB(nep,etype,viewer);
444: break;
445: default:
446: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
447: }
448: return(0);
449: }
451: /*@
452: NEPErrorViewFromOptions - Processes command line options to determine if/how
453: the errors of the computed solution are to be viewed.
455: Collective on nep
457: Input Parameter:
458: . nep - the nonlinear eigensolver context
460: Level: developer
461: @*/
462: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
463: {
464: PetscErrorCode ierr;
465: PetscViewer viewer;
466: PetscBool flg;
467: static PetscBool incall = PETSC_FALSE;
468: PetscViewerFormat format;
471: if (incall) return(0);
472: incall = PETSC_TRUE;
473: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg);
474: if (flg) {
475: PetscViewerPushFormat(viewer,format);
476: NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer);
477: PetscViewerPopFormat(viewer);
478: PetscViewerDestroy(&viewer);
479: }
480: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg);
481: if (flg) {
482: PetscViewerPushFormat(viewer,format);
483: NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer);
484: PetscViewerPopFormat(viewer);
485: PetscViewerDestroy(&viewer);
486: }
487: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_backward",&viewer,&format,&flg);
488: if (flg) {
489: PetscViewerPushFormat(viewer,format);
490: NEPErrorView(nep,NEP_ERROR_BACKWARD,viewer);
491: PetscViewerPopFormat(viewer);
492: PetscViewerDestroy(&viewer);
493: }
494: incall = PETSC_FALSE;
495: return(0);
496: }
498: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
499: {
501: PetscDraw draw;
502: PetscDrawSP drawsp;
503: PetscReal re,im;
504: PetscInt i,k;
507: if (!nep->nconv) return(0);
508: PetscViewerDrawGetDraw(viewer,0,&draw);
509: PetscDrawSetTitle(draw,"Computed Eigenvalues");
510: PetscDrawSPCreate(draw,1,&drawsp);
511: for (i=0;i<nep->nconv;i++) {
512: k = nep->perm[i];
513: #if defined(PETSC_USE_COMPLEX)
514: re = PetscRealPart(nep->eigr[k]);
515: im = PetscImaginaryPart(nep->eigr[k]);
516: #else
517: re = nep->eigr[k];
518: im = nep->eigi[k];
519: #endif
520: PetscDrawSPAddPoint(drawsp,&re,&im);
521: }
522: PetscDrawSPDraw(drawsp,PETSC_TRUE);
523: PetscDrawSPSave(drawsp);
524: PetscDrawSPDestroy(&drawsp);
525: return(0);
526: }
528: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
529: {
530: PetscInt i,k;
534: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
535: for (i=0;i<nep->nconv;i++) {
536: k = nep->perm[i];
537: PetscViewerASCIIPrintf(viewer," ");
538: SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]);
539: PetscViewerASCIIPrintf(viewer,"\n");
540: }
541: PetscViewerASCIIPrintf(viewer,"\n");
542: return(0);
543: }
545: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
546: {
548: PetscInt i,k;
549: PetscReal re,im;
550: const char *name;
553: PetscObjectGetName((PetscObject)nep,&name);
554: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
555: for (i=0;i<nep->nconv;i++) {
556: k = nep->perm[i];
557: #if defined(PETSC_USE_COMPLEX)
558: re = PetscRealPart(nep->eigr[k]);
559: im = PetscImaginaryPart(nep->eigr[k]);
560: #else
561: re = nep->eigr[k];
562: im = nep->eigi[k];
563: #endif
564: if (im!=0.0) {
565: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
566: } else {
567: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
568: }
569: }
570: PetscViewerASCIIPrintf(viewer,"];\n");
571: return(0);
572: }
574: /*@C
575: NEPValuesView - Displays the computed eigenvalues in a viewer.
577: Collective on nep
579: Input Parameters:
580: + nep - the nonlinear eigensolver context
581: - viewer - the viewer
583: Options Database Key:
584: . -nep_view_values - print computed eigenvalues
586: Level: intermediate
588: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
589: @*/
590: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
591: {
592: PetscBool isascii,isdraw;
593: PetscViewerFormat format;
594: PetscErrorCode ierr;
598: if (!viewer) {
599: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
600: }
603: NEPCheckSolved(nep,1);
604: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
605: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
606: if (isdraw) {
607: NEPValuesView_DRAW(nep,viewer);
608: } else if (isascii) {
609: PetscViewerGetFormat(viewer,&format);
610: switch (format) {
611: case PETSC_VIEWER_DEFAULT:
612: case PETSC_VIEWER_ASCII_INFO:
613: case PETSC_VIEWER_ASCII_INFO_DETAIL:
614: NEPValuesView_ASCII(nep,viewer);
615: break;
616: case PETSC_VIEWER_ASCII_MATLAB:
617: NEPValuesView_MATLAB(nep,viewer);
618: break;
619: default:
620: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
621: }
622: }
623: return(0);
624: }
626: /*@
627: NEPValuesViewFromOptions - Processes command line options to determine if/how
628: the computed eigenvalues are to be viewed.
630: Collective on nep
632: Input Parameter:
633: . nep - the nonlinear eigensolver context
635: Level: developer
636: @*/
637: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
638: {
639: PetscErrorCode ierr;
640: PetscViewer viewer;
641: PetscBool flg;
642: static PetscBool incall = PETSC_FALSE;
643: PetscViewerFormat format;
646: if (incall) return(0);
647: incall = PETSC_TRUE;
648: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg);
649: if (flg) {
650: PetscViewerPushFormat(viewer,format);
651: NEPValuesView(nep,viewer);
652: PetscViewerPopFormat(viewer);
653: PetscViewerDestroy(&viewer);
654: }
655: incall = PETSC_FALSE;
656: return(0);
657: }
659: /*@C
660: NEPVectorsView - Outputs computed eigenvectors to a viewer.
662: Collective on nep
664: Input Parameters:
665: + nep - the nonlinear eigensolver context
666: - viewer - the viewer
668: Options Database Keys:
669: . -nep_view_vectors - output eigenvectors.
671: Notes:
672: If PETSc was configured with real scalars, complex conjugate eigenvectors
673: will be viewed as two separate real vectors, one containing the real part
674: and another one containing the imaginary part.
676: If left eigenvectors were computed with a two-sided eigensolver, the right
677: and left eigenvectors are interleaved, that is, the vectors are output in
678: the following order: X0, Y0, X1, Y1, X2, Y2, ...
680: Level: intermediate
682: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
683: @*/
684: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
685: {
687: PetscInt i,k;
688: Vec x;
689: #define NMLEN 30
690: char vname[NMLEN];
691: const char *ename;
695: if (!viewer) {
696: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
697: }
700: NEPCheckSolved(nep,1);
701: if (nep->nconv) {
702: PetscObjectGetName((PetscObject)nep,&ename);
703: NEPComputeVectors(nep);
704: for (i=0;i<nep->nconv;i++) {
705: k = nep->perm[i];
706: PetscSNPrintf(vname,NMLEN,"X%d_%s",(int)i,ename);
707: BVGetColumn(nep->V,k,&x);
708: PetscObjectSetName((PetscObject)x,vname);
709: VecView(x,viewer);
710: BVRestoreColumn(nep->V,k,&x);
711: if (nep->twosided) {
712: PetscSNPrintf(vname,NMLEN,"Y%d_%s",(int)i,ename);
713: BVGetColumn(nep->W,k,&x);
714: PetscObjectSetName((PetscObject)x,vname);
715: VecView(x,viewer);
716: BVRestoreColumn(nep->W,k,&x);
717: }
718: }
719: }
720: return(0);
721: }
723: /*@
724: NEPVectorsViewFromOptions - Processes command line options to determine if/how
725: the computed eigenvectors are to be viewed.
727: Collective on nep
729: Input Parameter:
730: . nep - the nonlinear eigensolver context
732: Level: developer
733: @*/
734: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
735: {
736: PetscErrorCode ierr;
737: PetscViewer viewer;
738: PetscBool flg = PETSC_FALSE;
739: static PetscBool incall = PETSC_FALSE;
740: PetscViewerFormat format;
743: if (incall) return(0);
744: incall = PETSC_TRUE;
745: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg);
746: if (flg) {
747: PetscViewerPushFormat(viewer,format);
748: NEPVectorsView(nep,viewer);
749: PetscViewerPopFormat(viewer);
750: PetscViewerDestroy(&viewer);
751: }
752: incall = PETSC_FALSE;
753: return(0);
754: }