Actual source code: nepview.c
slepc-3.8.2 2017-12-01
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2017, 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,nods;
51: PetscViewer sviewer;
55: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
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: case NEP_USER_INTERFACE_DERIVATIVES:
83: PetscViewerASCIIPrintf(viewer," nonlinear operator from user callbacks for the successive derivatives\n");
84: break;
85: }
86: } else {
87: PetscViewerASCIIPrintf(viewer," nonlinear operator not specified yet\n");
88: }
89: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
90: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
91: SlepcSNPrintfScalar(str,50,nep->target,PETSC_FALSE);
92: if (!nep->which) {
93: PetscViewerASCIIPrintf(viewer,"not yet set\n");
94: } else switch (nep->which) {
95: case NEP_WHICH_USER:
96: PetscViewerASCIIPrintf(viewer,"user defined\n");
97: break;
98: case NEP_TARGET_MAGNITUDE:
99: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
100: break;
101: case NEP_TARGET_REAL:
102: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
103: break;
104: case NEP_TARGET_IMAGINARY:
105: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
106: break;
107: case NEP_LARGEST_MAGNITUDE:
108: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
109: break;
110: case NEP_SMALLEST_MAGNITUDE:
111: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
112: break;
113: case NEP_LARGEST_REAL:
114: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
115: break;
116: case NEP_SMALLEST_REAL:
117: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
118: break;
119: case NEP_LARGEST_IMAGINARY:
120: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
121: break;
122: case NEP_SMALLEST_IMAGINARY:
123: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
124: break;
125: case NEP_ALL:
126: PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
127: break;
128: }
129: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
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: PetscObjectTypeCompareAny((PetscObject)nep,&nods,NEPRII,NEPSLP,NEPINTERPOL,"");
178: if (!nods) {
179: if (!nep->ds) { NEPGetDS(nep,&nep->ds); }
180: DSView(nep->ds,viewer);
181: }
182: PetscViewerPopFormat(viewer);
183: if (nep->refine!=NEP_REFINE_NONE) {
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: }
193: return(0);
194: }
196: /*@C
197: NEPReasonView - Displays the reason a NEP solve converged or diverged.
199: Collective on NEP
201: Parameter:
202: + nep - the nonlinear eigensolver context
203: - viewer - the viewer to display the reason
205: Options Database Keys:
206: . -nep_converged_reason - print reason for convergence, and number of iterations
208: Level: intermediate
210: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber()
211: @*/
212: PetscErrorCode NEPReasonView(NEP nep,PetscViewer viewer)
213: {
215: PetscBool isAscii;
218: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
219: if (isAscii) {
220: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
221: if (nep->reason > 0) {
222: 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);
223: } else {
224: 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);
225: }
226: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
227: }
228: return(0);
229: }
231: /*@
232: NEPReasonViewFromOptions - Processes command line options to determine if/how
233: the NEP converged reason is to be viewed.
235: Collective on NEP
237: Input Parameters:
238: . nep - the nonlinear eigensolver context
240: Level: developer
241: @*/
242: PetscErrorCode NEPReasonViewFromOptions(NEP nep)
243: {
244: PetscErrorCode ierr;
245: PetscViewer viewer;
246: PetscBool flg;
247: static PetscBool incall = PETSC_FALSE;
248: PetscViewerFormat format;
251: if (incall) return(0);
252: incall = PETSC_TRUE;
253: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg);
254: if (flg) {
255: PetscViewerPushFormat(viewer,format);
256: NEPReasonView(nep,viewer);
257: PetscViewerPopFormat(viewer);
258: PetscViewerDestroy(&viewer);
259: }
260: incall = PETSC_FALSE;
261: return(0);
262: }
264: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
265: {
266: PetscReal error;
267: PetscInt i,j,k,nvals;
271: nvals = (nep->which==NEP_ALL)? nep->nconv: nep->nev;
272: if (nep->which!=NEP_ALL && nep->nconv<nvals) {
273: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",nep->nev);
274: return(0);
275: }
276: if (nep->which==NEP_ALL && !nvals) {
277: PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
278: return(0);
279: }
280: for (i=0;i<nvals;i++) {
281: NEPComputeError(nep,i,etype,&error);
282: if (error>=5.0*nep->tol) {
283: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
284: return(0);
285: }
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: k = nep->perm[8*i+j];
296: SlepcPrintEigenvalueASCII(nep->eigr[k],nep->eigi[k]);
297: if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
298: }
299: }
300: PetscViewerASCIIPrintf(viewer,"\n\n");
301: return(0);
302: }
304: static PetscErrorCode NEPErrorView_DETAIL(NEP nep,NEPErrorType etype,PetscViewer viewer)
305: {
307: PetscReal error,re,im;
308: PetscScalar kr,ki;
309: PetscInt i;
310: #define EXLEN 30
311: char ex[EXLEN],sep[]=" ---------------------- --------------------\n";
314: if (!nep->nconv) return(0);
315: switch (etype) {
316: case NEP_ERROR_ABSOLUTE:
317: PetscSNPrintf(ex,EXLEN," ||T(k)x||");
318: break;
319: case NEP_ERROR_RELATIVE:
320: PetscSNPrintf(ex,EXLEN," ||T(k)x||/||kx||");
321: break;
322: case NEP_ERROR_BACKWARD:
323: PetscSNPrintf(ex,EXLEN," eta(x,k)");
324: break;
325: }
326: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
327: for (i=0;i<nep->nconv;i++) {
328: NEPGetEigenpair(nep,i,&kr,&ki,NULL,NULL);
329: NEPComputeError(nep,i,etype,&error);
330: #if defined(PETSC_USE_COMPLEX)
331: re = PetscRealPart(kr);
332: im = PetscImaginaryPart(kr);
333: #else
334: re = kr;
335: im = ki;
336: #endif
337: if (im!=0.0) {
338: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
339: } else {
340: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
341: }
342: }
343: PetscViewerASCIIPrintf(viewer,sep);
344: return(0);
345: }
347: static PetscErrorCode NEPErrorView_MATLAB(NEP nep,NEPErrorType etype,PetscViewer viewer)
348: {
350: PetscReal error;
351: PetscInt i;
352: const char *name;
355: PetscObjectGetName((PetscObject)nep,&name);
356: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
357: for (i=0;i<nep->nconv;i++) {
358: NEPComputeError(nep,i,etype,&error);
359: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
360: }
361: PetscViewerASCIIPrintf(viewer,"];\n");
362: return(0);
363: }
365: /*@C
366: NEPErrorView - Displays the errors associated with the computed solution
367: (as well as the eigenvalues).
369: Collective on NEP
371: Input Parameters:
372: + nep - the nonlinear eigensolver context
373: . etype - error type
374: - viewer - optional visualization context
376: Options Database Key:
377: + -nep_error_absolute - print absolute errors of each eigenpair
378: . -nep_error_relative - print relative errors of each eigenpair
379: - -nep_error_backward - print backward errors of each eigenpair
381: Notes:
382: By default, this function checks the error of all eigenpairs and prints
383: the eigenvalues if all of them are below the requested tolerance.
384: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
385: eigenvalues and corresponding errors is printed.
387: Level: intermediate
389: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
390: @*/
391: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
392: {
393: PetscBool isascii;
394: PetscViewerFormat format;
395: PetscErrorCode ierr;
399: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
402: NEPCheckSolved(nep,1);
403: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
404: if (!isascii) return(0);
406: PetscViewerGetFormat(viewer,&format);
407: switch (format) {
408: case PETSC_VIEWER_DEFAULT:
409: case PETSC_VIEWER_ASCII_INFO:
410: NEPErrorView_ASCII(nep,etype,viewer);
411: break;
412: case PETSC_VIEWER_ASCII_INFO_DETAIL:
413: NEPErrorView_DETAIL(nep,etype,viewer);
414: break;
415: case PETSC_VIEWER_ASCII_MATLAB:
416: NEPErrorView_MATLAB(nep,etype,viewer);
417: break;
418: default:
419: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
420: }
421: return(0);
422: }
424: /*@
425: NEPErrorViewFromOptions - Processes command line options to determine if/how
426: the errors of the computed solution are to be viewed.
428: Collective on NEP
430: Input Parameters:
431: . nep - the nonlinear eigensolver context
433: Level: developer
434: @*/
435: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
436: {
437: PetscErrorCode ierr;
438: PetscViewer viewer;
439: PetscBool flg;
440: static PetscBool incall = PETSC_FALSE;
441: PetscViewerFormat format;
444: if (incall) return(0);
445: incall = PETSC_TRUE;
446: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg);
447: if (flg) {
448: PetscViewerPushFormat(viewer,format);
449: NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer);
450: PetscViewerPopFormat(viewer);
451: PetscViewerDestroy(&viewer);
452: }
453: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg);
454: if (flg) {
455: PetscViewerPushFormat(viewer,format);
456: NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer);
457: PetscViewerPopFormat(viewer);
458: PetscViewerDestroy(&viewer);
459: }
460: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_error_backward",&viewer,&format,&flg);
461: if (flg) {
462: PetscViewerPushFormat(viewer,format);
463: NEPErrorView(nep,NEP_ERROR_BACKWARD,viewer);
464: PetscViewerPopFormat(viewer);
465: PetscViewerDestroy(&viewer);
466: }
467: incall = PETSC_FALSE;
468: return(0);
469: }
471: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
472: {
474: PetscDraw draw;
475: PetscDrawSP drawsp;
476: PetscReal re,im;
477: PetscInt i,k;
480: if (!nep->nconv) return(0);
481: PetscViewerDrawGetDraw(viewer,0,&draw);
482: PetscDrawSetTitle(draw,"Computed Eigenvalues");
483: PetscDrawSPCreate(draw,1,&drawsp);
484: for (i=0;i<nep->nconv;i++) {
485: k = nep->perm[i];
486: #if defined(PETSC_USE_COMPLEX)
487: re = PetscRealPart(nep->eigr[k]);
488: im = PetscImaginaryPart(nep->eigr[k]);
489: #else
490: re = nep->eigr[k];
491: im = nep->eigi[k];
492: #endif
493: PetscDrawSPAddPoint(drawsp,&re,&im);
494: }
495: PetscDrawSPDraw(drawsp,PETSC_TRUE);
496: PetscDrawSPSave(drawsp);
497: PetscDrawSPDestroy(&drawsp);
498: return(0);
499: }
501: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
502: {
503: PetscInt i,k;
507: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
508: for (i=0;i<nep->nconv;i++) {
509: k = nep->perm[i];
510: PetscViewerASCIIPrintf(viewer," ");
511: SlepcPrintEigenvalueASCII(nep->eigr[k],nep->eigi[k]);
512: PetscViewerASCIIPrintf(viewer,"\n");
513: }
514: PetscViewerASCIIPrintf(viewer,"\n");
515: return(0);
516: }
518: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
519: {
521: PetscInt i,k;
522: PetscReal re,im;
523: const char *name;
526: PetscObjectGetName((PetscObject)nep,&name);
527: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
528: for (i=0;i<nep->nconv;i++) {
529: k = nep->perm[i];
530: #if defined(PETSC_USE_COMPLEX)
531: re = PetscRealPart(nep->eigr[k]);
532: im = PetscImaginaryPart(nep->eigr[k]);
533: #else
534: re = nep->eigr[k];
535: im = nep->eigi[k];
536: #endif
537: if (im!=0.0) {
538: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
539: } else {
540: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
541: }
542: }
543: PetscViewerASCIIPrintf(viewer,"];\n");
544: return(0);
545: }
547: /*@C
548: NEPValuesView - Displays the computed eigenvalues in a viewer.
550: Collective on NEP
552: Input Parameters:
553: + nep - the nonlinear eigensolver context
554: - viewer - the viewer
556: Options Database Key:
557: . -nep_view_values - print computed eigenvalues
559: Level: intermediate
561: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
562: @*/
563: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
564: {
565: PetscBool isascii,isdraw;
566: PetscViewerFormat format;
567: PetscErrorCode ierr;
571: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
574: NEPCheckSolved(nep,1);
575: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
576: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
577: if (isdraw) {
578: NEPValuesView_DRAW(nep,viewer);
579: } else if (isascii) {
580: PetscViewerGetFormat(viewer,&format);
581: switch (format) {
582: case PETSC_VIEWER_DEFAULT:
583: case PETSC_VIEWER_ASCII_INFO:
584: case PETSC_VIEWER_ASCII_INFO_DETAIL:
585: NEPValuesView_ASCII(nep,viewer);
586: break;
587: case PETSC_VIEWER_ASCII_MATLAB:
588: NEPValuesView_MATLAB(nep,viewer);
589: break;
590: default:
591: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
592: }
593: }
594: return(0);
595: }
597: /*@
598: NEPValuesViewFromOptions - Processes command line options to determine if/how
599: the computed eigenvalues are to be viewed.
601: Collective on NEP
603: Input Parameters:
604: . nep - the nonlinear eigensolver context
606: Level: developer
607: @*/
608: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
609: {
610: PetscErrorCode ierr;
611: PetscViewer viewer;
612: PetscBool flg;
613: static PetscBool incall = PETSC_FALSE;
614: PetscViewerFormat format;
617: if (incall) return(0);
618: incall = PETSC_TRUE;
619: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg);
620: if (flg) {
621: PetscViewerPushFormat(viewer,format);
622: NEPValuesView(nep,viewer);
623: PetscViewerPopFormat(viewer);
624: PetscViewerDestroy(&viewer);
625: }
626: incall = PETSC_FALSE;
627: return(0);
628: }
630: /*@C
631: NEPVectorsView - Outputs computed eigenvectors to a viewer.
633: Collective on NEP
635: Parameter:
636: + nep - the nonlinear eigensolver context
637: - viewer - the viewer
639: Options Database Keys:
640: . -nep_view_vectors - output eigenvectors.
642: Note:
643: If PETSc was configured with real scalars, complex conjugate eigenvectors
644: will be viewed as two separate real vectors, one containing the real part
645: and another one containing the imaginary part.
647: Level: intermediate
649: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
650: @*/
651: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
652: {
654: PetscInt i,k;
655: Vec x;
656: #define NMLEN 30
657: char vname[NMLEN];
658: const char *ename;
662: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
665: NEPCheckSolved(nep,1);
666: if (nep->nconv) {
667: PetscObjectGetName((PetscObject)nep,&ename);
668: NEPComputeVectors(nep);
669: for (i=0;i<nep->nconv;i++) {
670: k = nep->perm[i];
671: PetscSNPrintf(vname,NMLEN,"V%d_%s",(int)i,ename);
672: BVGetColumn(nep->V,k,&x);
673: PetscObjectSetName((PetscObject)x,vname);
674: VecView(x,viewer);
675: BVRestoreColumn(nep->V,k,&x);
676: }
677: }
678: return(0);
679: }
681: /*@
682: NEPVectorsViewFromOptions - Processes command line options to determine if/how
683: the computed eigenvectors are to be viewed.
685: Collective on NEP
687: Input Parameters:
688: . nep - the nonlinear eigensolver context
690: Level: developer
691: @*/
692: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
693: {
694: PetscErrorCode ierr;
695: PetscViewer viewer;
696: PetscBool flg = PETSC_FALSE;
697: static PetscBool incall = PETSC_FALSE;
698: PetscViewerFormat format;
701: if (incall) return(0);
702: incall = PETSC_TRUE;
703: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg);
704: if (flg) {
705: PetscViewerPushFormat(viewer,format);
706: NEPVectorsView(nep,viewer);
707: PetscViewerPopFormat(viewer);
708: PetscViewerDestroy(&viewer);
709: }
710: incall = PETSC_FALSE;
711: return(0);
712: }