Actual source code: nepview.c
slepc-3.16.3 2022-04-11
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2021, 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 <slepc/private/bvimpl.h>
16: #include <petscdraw.h>
18: /*@C
19: NEPView - Prints the NEP data structure.
21: Collective on nep
23: Input Parameters:
24: + nep - the nonlinear eigenproblem solver context
25: - viewer - optional visualization context
27: Options Database Key:
28: . -nep_view - Calls NEPView() at end of NEPSolve()
30: Note:
31: The available visualization contexts include
32: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
33: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
34: output where only the first processor opens
35: the file. All other processors send their
36: data to the first processor to print.
38: The user can open an alternative visualization context with
39: PetscViewerASCIIOpen() - output to a specified file.
41: Level: beginner
42: @*/
43: PetscErrorCode NEPView(NEP nep,PetscViewer viewer)
44: {
46: const char *type=NULL;
47: char str[50];
48: PetscInt i;
49: PetscBool isascii,istrivial;
50: PetscViewer sviewer;
54: if (!viewer) {
55: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
56: }
60: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
61: if (isascii) {
62: PetscObjectPrintClassNamePrefixType((PetscObject)nep,viewer);
63: if (nep->ops->view) {
64: PetscViewerASCIIPushTab(viewer);
65: (*nep->ops->view)(nep,viewer);
66: PetscViewerASCIIPopTab(viewer);
67: }
68: if (nep->problem_type) {
69: switch (nep->problem_type) {
70: case NEP_GENERAL: type = "general nonlinear eigenvalue problem"; break;
71: case NEP_RATIONAL: type = "rational eigenvalue problem"; break;
72: }
73: } else type = "not yet set";
74: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
75: if (nep->fui) {
76: switch (nep->fui) {
77: case NEP_USER_INTERFACE_CALLBACK:
78: PetscViewerASCIIPrintf(viewer," nonlinear operator from user callbacks\n");
79: break;
80: case NEP_USER_INTERFACE_SPLIT:
81: PetscViewerASCIIPrintf(viewer," nonlinear operator in split form, with %D terms\n",nep->nt);
82: break;
83: }
84: } else {
85: PetscViewerASCIIPrintf(viewer," nonlinear operator not specified yet\n");
86: }
87: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
88: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
89: SlepcSNPrintfScalar(str,sizeof(str),nep->target,PETSC_FALSE);
90: if (!nep->which) {
91: PetscViewerASCIIPrintf(viewer,"not yet set\n");
92: } else switch (nep->which) {
93: case NEP_WHICH_USER:
94: PetscViewerASCIIPrintf(viewer,"user defined\n");
95: break;
96: case NEP_TARGET_MAGNITUDE:
97: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
98: break;
99: case NEP_TARGET_REAL:
100: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
101: break;
102: case NEP_TARGET_IMAGINARY:
103: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
104: break;
105: case NEP_LARGEST_MAGNITUDE:
106: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
107: break;
108: case NEP_SMALLEST_MAGNITUDE:
109: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
110: break;
111: case NEP_LARGEST_REAL:
112: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
113: break;
114: case NEP_SMALLEST_REAL:
115: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
116: break;
117: case NEP_LARGEST_IMAGINARY:
118: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
119: break;
120: case NEP_SMALLEST_IMAGINARY:
121: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
122: break;
123: case NEP_ALL:
124: PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
125: break;
126: }
127: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
128: if (nep->twosided) {
129: PetscViewerASCIIPrintf(viewer," using two-sided variant (for left eigenvectors)\n");
130: }
131: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",nep->nev);
132: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",nep->ncv);
133: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",nep->mpd);
134: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",nep->max_it);
135: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)nep->tol);
136: PetscViewerASCIIPrintf(viewer," convergence test: ");
137: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
138: switch (nep->conv) {
139: case NEP_CONV_ABS:
140: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
141: case NEP_CONV_REL:
142: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
143: case NEP_CONV_NORM:
144: PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
145: if (nep->nrma) {
146: PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)nep->nrma[0]);
147: for (i=1;i<nep->nt;i++) {
148: PetscViewerASCIIPrintf(viewer,", %g",(double)nep->nrma[i]);
149: }
150: PetscViewerASCIIPrintf(viewer,"\n");
151: }
152: break;
153: case NEP_CONV_USER:
154: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
155: }
156: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
157: if (nep->refine) {
158: PetscViewerASCIIPrintf(viewer," iterative refinement: %s, with %s scheme\n",NEPRefineTypes[nep->refine],NEPRefineSchemes[nep->scheme]);
159: PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%D\n",(double)nep->rtol,nep->rits);
160: if (nep->npart>1) {
161: PetscViewerASCIIPrintf(viewer," splitting communicator in %D partitions for refinement\n",nep->npart);
162: }
163: }
164: if (nep->nini) {
165: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(nep->nini));
166: }
167: } else {
168: if (nep->ops->view) {
169: (*nep->ops->view)(nep,viewer);
170: }
171: }
172: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
173: if (!nep->V) { NEPGetBV(nep,&nep->V); }
174: BVView(nep->V,viewer);
175: if (!nep->rg) { NEPGetRG(nep,&nep->rg); }
176: RGIsTrivial(nep->rg,&istrivial);
177: if (!istrivial) { RGView(nep->rg,viewer); }
178: if (nep->useds) {
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: PetscViewerASCIIPushTab(viewer);
185: if (nep->npart>1) {
186: if (nep->refinesubc->color==0) {
187: PetscViewerASCIIGetStdout(PetscSubcommChild(nep->refinesubc),&sviewer);
188: KSPView(nep->refineksp,sviewer);
189: }
190: } else {
191: KSPView(nep->refineksp,viewer);
192: }
193: PetscViewerASCIIPopTab(viewer);
194: }
195: return(0);
196: }
198: /*@C
199: NEPViewFromOptions - View from options
201: Collective on NEP
203: Input Parameters:
204: + nep - the nonlinear eigensolver context
205: . obj - optional object
206: - name - command line option
208: Level: intermediate
210: .seealso: NEPView(), NEPCreate()
211: @*/
212: PetscErrorCode NEPViewFromOptions(NEP nep,PetscObject obj,const char name[])
213: {
218: PetscObjectViewFromOptions((PetscObject)nep,obj,name);
219: return(0);
220: }
222: /*@C
223: NEPConvergedReasonView - Displays the reason a NEP solve converged or diverged.
225: Collective on nep
227: Input Parameters:
228: + nep - the nonlinear eigensolver context
229: - viewer - the viewer to display the reason
231: Options Database Keys:
232: . -nep_converged_reason - print reason for convergence, and number of iterations
234: Note:
235: To change the format of the output call PetscViewerPushFormat(viewer,format) before
236: this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
237: display a reason if it fails. The latter can be set in the command line with
238: -nep_converged_reason ::failed
240: Level: intermediate
242: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber(), NEPConvergedReasonViewFromOptions(
243: @*/
244: PetscErrorCode NEPConvergedReasonView(NEP nep,PetscViewer viewer)
245: {
246: PetscErrorCode ierr;
247: PetscBool isAscii;
248: PetscViewerFormat format;
251: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
252: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
253: if (isAscii) {
254: PetscViewerGetFormat(viewer,&format);
255: PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel);
256: if (nep->reason > 0 && format != PETSC_VIEWER_FAILED) {
257: 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);
258: } else if (nep->reason <= 0) {
259: 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);
260: }
261: PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel);
262: }
263: return(0);
264: }
266: /*@
267: NEPConvergedReasonViewFromOptions - Processes command line options to determine if/how
268: the NEP converged reason is to be viewed.
270: Collective on nep
272: Input Parameter:
273: . nep - the nonlinear eigensolver context
275: Level: developer
277: .seealso: NEPConvergedReasonView()
278: @*/
279: PetscErrorCode NEPConvergedReasonViewFromOptions(NEP nep)
280: {
281: PetscErrorCode ierr;
282: PetscViewer viewer;
283: PetscBool flg;
284: static PetscBool incall = PETSC_FALSE;
285: PetscViewerFormat format;
288: if (incall) return(0);
289: incall = PETSC_TRUE;
290: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg);
291: if (flg) {
292: PetscViewerPushFormat(viewer,format);
293: NEPConvergedReasonView(nep,viewer);
294: PetscViewerPopFormat(viewer);
295: PetscViewerDestroy(&viewer);
296: }
297: incall = PETSC_FALSE;
298: return(0);
299: }
301: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
302: {
303: PetscReal error;
304: PetscInt i,j,k,nvals;
308: nvals = (nep->which==NEP_ALL)? nep->nconv: nep->nev;
309: if (nep->which!=NEP_ALL && nep->nconv<nvals) {
310: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",nep->nev);
311: return(0);
312: }
313: if (nep->which==NEP_ALL && !nvals) {
314: PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
315: return(0);
316: }
317: for (i=0;i<nvals;i++) {
318: NEPComputeError(nep,i,etype,&error);
319: if (error>=5.0*nep->tol) {
320: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
321: return(0);
322: }
323: }
324: if (nep->which==NEP_ALL) {
325: PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
326: } else {
327: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
328: }
329: for (i=0;i<=(nvals-1)/8;i++) {
330: PetscViewerASCIIPrintf(viewer,"\n ");
331: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
332: k = nep->perm[8*i+j];
333: SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]);
334: if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
335: }
336: }
337: PetscViewerASCIIPrintf(viewer,"\n\n");
338: return(0);
339: }
341: static PetscErrorCode NEPErrorView_DETAIL(NEP nep,NEPErrorType etype,PetscViewer viewer)
342: {
344: PetscReal error,re,im;
345: PetscScalar kr,ki;
346: PetscInt i;
347: char ex[30],sep[]=" ---------------------- --------------------\n";
350: if (!nep->nconv) return(0);
351: switch (etype) {
352: case NEP_ERROR_ABSOLUTE:
353: PetscSNPrintf(ex,sizeof(ex)," ||T(k)x||");
354: break;
355: case NEP_ERROR_RELATIVE:
356: PetscSNPrintf(ex,sizeof(ex)," ||T(k)x||/||kx||");
357: break;
358: case NEP_ERROR_BACKWARD:
359: PetscSNPrintf(ex,sizeof(ex)," eta(x,k)");
360: break;
361: }
362: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
363: for (i=0;i<nep->nconv;i++) {
364: NEPGetEigenpair(nep,i,&kr,&ki,NULL,NULL);
365: NEPComputeError(nep,i,etype,&error);
366: #if defined(PETSC_USE_COMPLEX)
367: re = PetscRealPart(kr);
368: im = PetscImaginaryPart(kr);
369: #else
370: re = kr;
371: im = ki;
372: #endif
373: if (im!=0.0) {
374: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
375: } else {
376: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
377: }
378: }
379: PetscViewerASCIIPrintf(viewer,"%s",sep);
380: return(0);
381: }
383: static PetscErrorCode NEPErrorView_MATLAB(NEP nep,NEPErrorType etype,PetscViewer viewer)
384: {
386: PetscReal error;
387: PetscInt i;
388: const char *name;
391: PetscObjectGetName((PetscObject)nep,&name);
392: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
393: for (i=0;i<nep->nconv;i++) {
394: NEPComputeError(nep,i,etype,&error);
395: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
396: }
397: PetscViewerASCIIPrintf(viewer,"];\n");
398: return(0);
399: }
401: /*@C
402: NEPErrorView - Displays the errors associated with the computed solution
403: (as well as the eigenvalues).
405: Collective on nep
407: Input Parameters:
408: + nep - the nonlinear eigensolver context
409: . etype - error type
410: - viewer - optional visualization context
412: Options Database Key:
413: + -nep_error_absolute - print absolute errors of each eigenpair
414: . -nep_error_relative - print relative errors of each eigenpair
415: - -nep_error_backward - print backward errors of each eigenpair
417: Notes:
418: By default, this function checks the error of all eigenpairs and prints
419: the eigenvalues if all of them are below the requested tolerance.
420: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
421: eigenvalues and corresponding errors is printed.
423: Level: intermediate
425: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
426: @*/
427: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
428: {
429: PetscBool isascii;
430: PetscViewerFormat format;
431: PetscErrorCode ierr;
435: if (!viewer) {
436: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
437: }
440: NEPCheckSolved(nep,1);
441: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
442: if (!isascii) return(0);
444: PetscViewerGetFormat(viewer,&format);
445: switch (format) {
446: case PETSC_VIEWER_DEFAULT:
447: case PETSC_VIEWER_ASCII_INFO:
448: NEPErrorView_ASCII(nep,etype,viewer);
449: break;
450: case PETSC_VIEWER_ASCII_INFO_DETAIL:
451: NEPErrorView_DETAIL(nep,etype,viewer);
452: break;
453: case PETSC_VIEWER_ASCII_MATLAB:
454: NEPErrorView_MATLAB(nep,etype,viewer);
455: break;
456: default:
457: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
458: }
459: return(0);
460: }
462: /*@
463: NEPErrorViewFromOptions - Processes command line options to determine if/how
464: the errors of the computed solution are to be viewed.
466: Collective on nep
468: Input Parameter:
469: . nep - the nonlinear eigensolver context
471: Level: developer
472: @*/
473: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
474: {
475: PetscErrorCode ierr;
476: PetscViewer viewer;
477: PetscBool flg;
478: static PetscBool incall = PETSC_FALSE;
479: PetscViewerFormat format;
482: if (incall) return(0);
483: incall = PETSC_TRUE;
484: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg);
485: if (flg) {
486: PetscViewerPushFormat(viewer,format);
487: NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer);
488: PetscViewerPopFormat(viewer);
489: PetscViewerDestroy(&viewer);
490: }
491: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg);
492: if (flg) {
493: PetscViewerPushFormat(viewer,format);
494: NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer);
495: PetscViewerPopFormat(viewer);
496: PetscViewerDestroy(&viewer);
497: }
498: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_backward",&viewer,&format,&flg);
499: if (flg) {
500: PetscViewerPushFormat(viewer,format);
501: NEPErrorView(nep,NEP_ERROR_BACKWARD,viewer);
502: PetscViewerPopFormat(viewer);
503: PetscViewerDestroy(&viewer);
504: }
505: incall = PETSC_FALSE;
506: return(0);
507: }
509: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
510: {
512: PetscDraw draw;
513: PetscDrawSP drawsp;
514: PetscReal re,im;
515: PetscInt i,k;
518: if (!nep->nconv) return(0);
519: PetscViewerDrawGetDraw(viewer,0,&draw);
520: PetscDrawSetTitle(draw,"Computed Eigenvalues");
521: PetscDrawSPCreate(draw,1,&drawsp);
522: for (i=0;i<nep->nconv;i++) {
523: k = nep->perm[i];
524: #if defined(PETSC_USE_COMPLEX)
525: re = PetscRealPart(nep->eigr[k]);
526: im = PetscImaginaryPart(nep->eigr[k]);
527: #else
528: re = nep->eigr[k];
529: im = nep->eigi[k];
530: #endif
531: PetscDrawSPAddPoint(drawsp,&re,&im);
532: }
533: PetscDrawSPDraw(drawsp,PETSC_TRUE);
534: PetscDrawSPSave(drawsp);
535: PetscDrawSPDestroy(&drawsp);
536: return(0);
537: }
539: static PetscErrorCode NEPValuesView_BINARY(NEP nep,PetscViewer viewer)
540: {
541: #if defined(PETSC_HAVE_COMPLEX)
542: PetscInt i,k;
543: PetscComplex *ev;
545: #endif
548: #if defined(PETSC_HAVE_COMPLEX)
549: PetscMalloc1(nep->nconv,&ev);
550: for (i=0;i<nep->nconv;i++) {
551: k = nep->perm[i];
552: #if defined(PETSC_USE_COMPLEX)
553: ev[i] = nep->eigr[k];
554: #else
555: ev[i] = PetscCMPLX(nep->eigr[k],nep->eigi[k]);
556: #endif
557: }
558: PetscViewerBinaryWrite(viewer,ev,nep->nconv,PETSC_COMPLEX);
559: PetscFree(ev);
560: #endif
561: return(0);
562: }
564: #if defined(PETSC_HAVE_HDF5)
565: static PetscErrorCode NEPValuesView_HDF5(NEP nep,PetscViewer viewer)
566: {
568: PetscInt i,k,n,N;
569: PetscMPIInt rank;
570: Vec v;
571: char vname[30];
572: const char *ename;
575: MPI_Comm_rank(PetscObjectComm((PetscObject)nep),&rank);
576: N = nep->nconv;
577: n = rank? 0: N;
578: /* create a vector containing the eigenvalues */
579: VecCreateMPI(PetscObjectComm((PetscObject)nep),n,N,&v);
580: PetscObjectGetName((PetscObject)nep,&ename);
581: PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename);
582: PetscObjectSetName((PetscObject)v,vname);
583: if (!rank) {
584: for (i=0;i<nep->nconv;i++) {
585: k = nep->perm[i];
586: VecSetValue(v,i,nep->eigr[k],INSERT_VALUES);
587: }
588: }
589: VecAssemblyBegin(v);
590: VecAssemblyEnd(v);
591: VecView(v,viewer);
592: #if !defined(PETSC_USE_COMPLEX)
593: /* in real scalars write the imaginary part as a separate vector */
594: PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename);
595: PetscObjectSetName((PetscObject)v,vname);
596: if (!rank) {
597: for (i=0;i<nep->nconv;i++) {
598: k = nep->perm[i];
599: VecSetValue(v,i,nep->eigi[k],INSERT_VALUES);
600: }
601: }
602: VecAssemblyBegin(v);
603: VecAssemblyEnd(v);
604: VecView(v,viewer);
605: #endif
606: VecDestroy(&v);
607: return(0);
608: }
609: #endif
611: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
612: {
613: PetscInt i,k;
617: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
618: for (i=0;i<nep->nconv;i++) {
619: k = nep->perm[i];
620: PetscViewerASCIIPrintf(viewer," ");
621: SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]);
622: PetscViewerASCIIPrintf(viewer,"\n");
623: }
624: PetscViewerASCIIPrintf(viewer,"\n");
625: return(0);
626: }
628: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
629: {
631: PetscInt i,k;
632: PetscReal re,im;
633: const char *name;
636: PetscObjectGetName((PetscObject)nep,&name);
637: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
638: for (i=0;i<nep->nconv;i++) {
639: k = nep->perm[i];
640: #if defined(PETSC_USE_COMPLEX)
641: re = PetscRealPart(nep->eigr[k]);
642: im = PetscImaginaryPart(nep->eigr[k]);
643: #else
644: re = nep->eigr[k];
645: im = nep->eigi[k];
646: #endif
647: if (im!=0.0) {
648: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
649: } else {
650: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
651: }
652: }
653: PetscViewerASCIIPrintf(viewer,"];\n");
654: return(0);
655: }
657: /*@C
658: NEPValuesView - Displays the computed eigenvalues in a viewer.
660: Collective on nep
662: Input Parameters:
663: + nep - the nonlinear eigensolver context
664: - viewer - the viewer
666: Options Database Key:
667: . -nep_view_values - print computed eigenvalues
669: Level: intermediate
671: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
672: @*/
673: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
674: {
675: PetscBool isascii,isdraw,isbinary;
676: PetscViewerFormat format;
677: PetscErrorCode ierr;
678: #if defined(PETSC_HAVE_HDF5)
679: PetscBool ishdf5;
680: #endif
684: if (!viewer) {
685: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
686: }
689: NEPCheckSolved(nep,1);
690: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
691: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
692: #if defined(PETSC_HAVE_HDF5)
693: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
694: #endif
695: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
696: if (isdraw) {
697: NEPValuesView_DRAW(nep,viewer);
698: } else if (isbinary) {
699: NEPValuesView_BINARY(nep,viewer);
700: #if defined(PETSC_HAVE_HDF5)
701: } else if (ishdf5) {
702: NEPValuesView_HDF5(nep,viewer);
703: #endif
704: } else if (isascii) {
705: PetscViewerGetFormat(viewer,&format);
706: switch (format) {
707: case PETSC_VIEWER_DEFAULT:
708: case PETSC_VIEWER_ASCII_INFO:
709: case PETSC_VIEWER_ASCII_INFO_DETAIL:
710: NEPValuesView_ASCII(nep,viewer);
711: break;
712: case PETSC_VIEWER_ASCII_MATLAB:
713: NEPValuesView_MATLAB(nep,viewer);
714: break;
715: default:
716: PetscInfo1(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
717: }
718: }
719: return(0);
720: }
722: /*@
723: NEPValuesViewFromOptions - Processes command line options to determine if/how
724: the computed eigenvalues are to be viewed.
726: Collective on nep
728: Input Parameter:
729: . nep - the nonlinear eigensolver context
731: Level: developer
732: @*/
733: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
734: {
735: PetscErrorCode ierr;
736: PetscViewer viewer;
737: PetscBool flg;
738: static PetscBool incall = PETSC_FALSE;
739: PetscViewerFormat format;
742: if (incall) return(0);
743: incall = PETSC_TRUE;
744: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg);
745: if (flg) {
746: PetscViewerPushFormat(viewer,format);
747: NEPValuesView(nep,viewer);
748: PetscViewerPopFormat(viewer);
749: PetscViewerDestroy(&viewer);
750: }
751: incall = PETSC_FALSE;
752: return(0);
753: }
755: /*@C
756: NEPVectorsView - Outputs computed eigenvectors to a viewer.
758: Collective on nep
760: Input Parameters:
761: + nep - the nonlinear eigensolver context
762: - viewer - the viewer
764: Options Database Keys:
765: . -nep_view_vectors - output eigenvectors.
767: Notes:
768: If PETSc was configured with real scalars, complex conjugate eigenvectors
769: will be viewed as two separate real vectors, one containing the real part
770: and another one containing the imaginary part.
772: If left eigenvectors were computed with a two-sided eigensolver, the right
773: and left eigenvectors are interleaved, that is, the vectors are output in
774: the following order: X0, Y0, X1, Y1, X2, Y2, ...
776: Level: intermediate
778: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
779: @*/
780: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
781: {
783: PetscInt i,k;
784: Vec xr,xi=NULL;
788: if (!viewer) {
789: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer);
790: }
793: NEPCheckSolved(nep,1);
794: if (nep->nconv) {
795: NEPComputeVectors(nep);
796: BVCreateVec(nep->V,&xr);
797: #if !defined(PETSC_USE_COMPLEX)
798: BVCreateVec(nep->V,&xi);
799: #endif
800: for (i=0;i<nep->nconv;i++) {
801: k = nep->perm[i];
802: BV_GetEigenvector(nep->V,k,nep->eigi[k],xr,xi);
803: SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)nep);
804: if (nep->twosided) {
805: BV_GetEigenvector(nep->W,k,nep->eigi[k],xr,xi);
806: SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)nep);
807: }
808: }
809: VecDestroy(&xr);
810: #if !defined(PETSC_USE_COMPLEX)
811: VecDestroy(&xi);
812: #endif
813: }
814: return(0);
815: }
817: /*@
818: NEPVectorsViewFromOptions - Processes command line options to determine if/how
819: the computed eigenvectors are to be viewed.
821: Collective on nep
823: Input Parameter:
824: . nep - the nonlinear eigensolver context
826: Level: developer
827: @*/
828: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
829: {
830: PetscErrorCode ierr;
831: PetscViewer viewer;
832: PetscBool flg = PETSC_FALSE;
833: static PetscBool incall = PETSC_FALSE;
834: PetscViewerFormat format;
837: if (incall) return(0);
838: incall = PETSC_TRUE;
839: PetscOptionsGetViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg);
840: if (flg) {
841: PetscViewerPushFormat(viewer,format);
842: NEPVectorsView(nep,viewer);
843: PetscViewerPopFormat(viewer);
844: PetscViewerDestroy(&viewer);
845: }
846: incall = PETSC_FALSE;
847: return(0);
848: }