Actual source code: pepview.c
slepc-3.16.2 2022-02-01
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: The PEP routines related to various viewers
12: */
14: #include <slepc/private/pepimpl.h>
15: #include <slepc/private/bvimpl.h>
16: #include <petscdraw.h>
18: /*@C
19: PEPView - Prints the PEP data structure.
21: Collective on pep
23: Input Parameters:
24: + pep - the polynomial eigenproblem solver context
25: - viewer - optional visualization context
27: Options Database Key:
28: . -pep_view - Calls PEPView() at end of PEPSolve()
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 PEPView(PEP pep,PetscViewer viewer)
44: {
46: const char *type=NULL;
47: char str[50];
48: PetscBool isascii,islinear,istrivial;
49: PetscInt i;
50: PetscViewer sviewer;
54: if (!viewer) {
55: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
56: }
60: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
61: if (isascii) {
62: PetscObjectPrintClassNamePrefixType((PetscObject)pep,viewer);
63: if (pep->ops->view) {
64: PetscViewerASCIIPushTab(viewer);
65: (*pep->ops->view)(pep,viewer);
66: PetscViewerASCIIPopTab(viewer);
67: }
68: if (pep->problem_type) {
69: switch (pep->problem_type) {
70: case PEP_GENERAL: type = "general polynomial eigenvalue problem"; break;
71: case PEP_HERMITIAN: type = SLEPC_STRING_HERMITIAN " polynomial eigenvalue problem"; break;
72: case PEP_HYPERBOLIC: type = "hyperbolic polynomial eigenvalue problem"; break;
73: case PEP_GYROSCOPIC: type = "gyroscopic polynomial eigenvalue problem"; break;
74: }
75: } else type = "not yet set";
76: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
77: PetscViewerASCIIPrintf(viewer," polynomial represented in %s basis\n",PEPBasisTypes[pep->basis]);
78: switch (pep->scale) {
79: case PEP_SCALE_NONE:
80: break;
81: case PEP_SCALE_SCALAR:
82: PetscViewerASCIIPrintf(viewer," parameter scaling enabled, with scaling factor=%g\n",(double)pep->sfactor);
83: break;
84: case PEP_SCALE_DIAGONAL:
85: PetscViewerASCIIPrintf(viewer," diagonal balancing enabled, with its=%D and lambda=%g\n",pep->sits,(double)pep->slambda);
86: break;
87: case PEP_SCALE_BOTH:
88: PetscViewerASCIIPrintf(viewer," parameter scaling & diagonal balancing enabled, with scaling factor=%g, its=%D and lambda=%g\n",(double)pep->sfactor,pep->sits,(double)pep->slambda);
89: break;
90: }
91: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
92: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
93: SlepcSNPrintfScalar(str,sizeof(str),pep->target,PETSC_FALSE);
94: if (!pep->which) {
95: PetscViewerASCIIPrintf(viewer,"not yet set\n");
96: } else switch (pep->which) {
97: case PEP_WHICH_USER:
98: PetscViewerASCIIPrintf(viewer,"user defined\n");
99: break;
100: case PEP_TARGET_MAGNITUDE:
101: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
102: break;
103: case PEP_TARGET_REAL:
104: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
105: break;
106: case PEP_TARGET_IMAGINARY:
107: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
108: break;
109: case PEP_LARGEST_MAGNITUDE:
110: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
111: break;
112: case PEP_SMALLEST_MAGNITUDE:
113: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
114: break;
115: case PEP_LARGEST_REAL:
116: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
117: break;
118: case PEP_SMALLEST_REAL:
119: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
120: break;
121: case PEP_LARGEST_IMAGINARY:
122: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
123: break;
124: case PEP_SMALLEST_IMAGINARY:
125: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
126: break;
127: case PEP_ALL:
128: if (pep->inta || pep->intb) {
129: PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)pep->inta,(double)pep->intb);
130: } else {
131: PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
132: }
133: break;
134: }
135: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
136: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %D\n",pep->nev);
137: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %D\n",pep->ncv);
138: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %D\n",pep->mpd);
139: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %D\n",pep->max_it);
140: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)pep->tol);
141: PetscViewerASCIIPrintf(viewer," convergence test: ");
142: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
143: switch (pep->conv) {
144: case PEP_CONV_ABS:
145: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
146: case PEP_CONV_REL:
147: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
148: case PEP_CONV_NORM:
149: PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n");
150: if (pep->nrma) {
151: PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)pep->nrma[0]);
152: for (i=1;i<pep->nmat;i++) {
153: PetscViewerASCIIPrintf(viewer,", %g",(double)pep->nrma[i]);
154: }
155: PetscViewerASCIIPrintf(viewer,"\n");
156: }
157: break;
158: case PEP_CONV_USER:
159: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
160: }
161: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
162: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",PEPExtractTypes[pep->extract]);
163: if (pep->refine) {
164: PetscViewerASCIIPrintf(viewer," iterative refinement: %s, with %s scheme\n",PEPRefineTypes[pep->refine],PEPRefineSchemes[pep->scheme]);
165: PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%D\n",(double)pep->rtol,pep->rits);
166: if (pep->npart>1) {
167: PetscViewerASCIIPrintf(viewer," splitting communicator in %D partitions for refinement\n",pep->npart);
168: }
169: }
170: if (pep->nini) {
171: PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %D\n",PetscAbs(pep->nini));
172: }
173: } else {
174: if (pep->ops->view) {
175: (*pep->ops->view)(pep,viewer);
176: }
177: }
178: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
179: if (!pep->V) { PEPGetBV(pep,&pep->V); }
180: BVView(pep->V,viewer);
181: if (!pep->rg) { PEPGetRG(pep,&pep->rg); }
182: RGIsTrivial(pep->rg,&istrivial);
183: if (!istrivial) { RGView(pep->rg,viewer); }
184: PetscObjectTypeCompare((PetscObject)pep,PEPLINEAR,&islinear);
185: if (!islinear) {
186: if (!pep->ds) { PEPGetDS(pep,&pep->ds); }
187: DSView(pep->ds,viewer);
188: }
189: PetscViewerPopFormat(viewer);
190: if (!pep->st) { PEPGetST(pep,&pep->st); }
191: STView(pep->st,viewer);
192: if (pep->refine!=PEP_REFINE_NONE) {
193: PetscViewerASCIIPushTab(viewer);
194: if (pep->npart>1) {
195: if (pep->refinesubc->color==0) {
196: PetscViewerASCIIGetStdout(PetscSubcommChild(pep->refinesubc),&sviewer);
197: KSPView(pep->refineksp,sviewer);
198: }
199: } else {
200: KSPView(pep->refineksp,viewer);
201: }
202: PetscViewerASCIIPopTab(viewer);
203: }
204: return(0);
205: }
207: /*@C
208: PEPViewFromOptions - View from options
210: Collective on PEP
212: Input Parameters:
213: + pep - the eigensolver context
214: . obj - optional object
215: - name - command line option
217: Level: intermediate
219: .seealso: PEPView(), PEPCreate()
220: @*/
221: PetscErrorCode PEPViewFromOptions(PEP pep,PetscObject obj,const char name[])
222: {
227: PetscObjectViewFromOptions((PetscObject)pep,obj,name);
228: return(0);
229: }
231: /*@C
232: PEPConvergedReasonView - Displays the reason a PEP solve converged or diverged.
234: Collective on pep
236: Input Parameters:
237: + pep - the eigensolver context
238: - viewer - the viewer to display the reason
240: Options Database Keys:
241: . -pep_converged_reason - print reason for convergence, and number of iterations
243: Note:
244: To change the format of the output call PetscViewerPushFormat(viewer,format) before
245: this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
246: display a reason if it fails. The latter can be set in the command line with
247: -pep_converged_reason ::failed
249: Level: intermediate
251: .seealso: PEPSetConvergenceTest(), PEPSetTolerances(), PEPGetIterationNumber(), PEPConvergedReasonViewFromOptions()
252: @*/
253: PetscErrorCode PEPConvergedReasonView(PEP pep,PetscViewer viewer)
254: {
255: PetscErrorCode ierr;
256: PetscBool isAscii;
257: PetscViewerFormat format;
260: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)pep));
261: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
262: if (isAscii) {
263: PetscViewerGetFormat(viewer,&format);
264: PetscViewerASCIIAddTab(viewer,((PetscObject)pep)->tablevel);
265: if (pep->reason > 0 && format != PETSC_VIEWER_FAILED) {
266: PetscViewerASCIIPrintf(viewer,"%s Polynomial eigensolve converged (%D eigenpair%s) due to %s; iterations %D\n",((PetscObject)pep)->prefix?((PetscObject)pep)->prefix:"",pep->nconv,(pep->nconv>1)?"s":"",PEPConvergedReasons[pep->reason],pep->its);
267: } else if (pep->reason <= 0) {
268: PetscViewerASCIIPrintf(viewer,"%s Polynomial eigensolve did not converge due to %s; iterations %D\n",((PetscObject)pep)->prefix?((PetscObject)pep)->prefix:"",PEPConvergedReasons[pep->reason],pep->its);
269: }
270: PetscViewerASCIISubtractTab(viewer,((PetscObject)pep)->tablevel);
271: }
272: return(0);
273: }
275: /*@
276: PEPConvergedReasonViewFromOptions - Processes command line options to determine if/how
277: the PEP converged reason is to be viewed.
279: Collective on pep
281: Input Parameter:
282: . pep - the eigensolver context
284: Level: developer
286: .seealso: PEPConvergedReasonView()
287: @*/
288: PetscErrorCode PEPConvergedReasonViewFromOptions(PEP pep)
289: {
290: PetscErrorCode ierr;
291: PetscViewer viewer;
292: PetscBool flg;
293: static PetscBool incall = PETSC_FALSE;
294: PetscViewerFormat format;
297: if (incall) return(0);
298: incall = PETSC_TRUE;
299: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_converged_reason",&viewer,&format,&flg);
300: if (flg) {
301: PetscViewerPushFormat(viewer,format);
302: PEPConvergedReasonView(pep,viewer);
303: PetscViewerPopFormat(viewer);
304: PetscViewerDestroy(&viewer);
305: }
306: incall = PETSC_FALSE;
307: return(0);
308: }
310: static PetscErrorCode PEPErrorView_ASCII(PEP pep,PEPErrorType etype,PetscViewer viewer)
311: {
312: PetscReal error;
313: PetscInt i,j,k,nvals;
317: nvals = (pep->which==PEP_ALL)? pep->nconv: pep->nev;
318: if (pep->which!=PEP_ALL && pep->nconv<pep->nev) {
319: PetscViewerASCIIPrintf(viewer," Problem: less than %D eigenvalues converged\n\n",pep->nev);
320: return(0);
321: }
322: if (pep->which==PEP_ALL && !nvals) {
323: PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
324: return(0);
325: }
326: for (i=0;i<nvals;i++) {
327: PEPComputeError(pep,i,etype,&error);
328: if (error>=5.0*pep->tol) {
329: PetscViewerASCIIPrintf(viewer," Problem: some of the first %D relative errors are higher than the tolerance\n\n",nvals);
330: return(0);
331: }
332: }
333: if (pep->which==PEP_ALL) {
334: PetscViewerASCIIPrintf(viewer," Found %D eigenvalues, all of them computed up to the required tolerance:",nvals);
335: } else {
336: PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
337: }
338: for (i=0;i<=(nvals-1)/8;i++) {
339: PetscViewerASCIIPrintf(viewer,"\n ");
340: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
341: k = pep->perm[8*i+j];
342: SlepcPrintEigenvalueASCII(viewer,pep->eigr[k],pep->eigi[k]);
343: if (8*i+j+1<nvals) { PetscViewerASCIIPrintf(viewer,", "); }
344: }
345: }
346: PetscViewerASCIIPrintf(viewer,"\n\n");
347: return(0);
348: }
350: static PetscErrorCode PEPErrorView_DETAIL(PEP pep,PEPErrorType etype,PetscViewer viewer)
351: {
353: PetscReal error,re,im;
354: PetscScalar kr,ki;
355: PetscInt i;
356: char ex[30],sep[]=" ---------------------- --------------------\n";
359: if (!pep->nconv) return(0);
360: switch (etype) {
361: case PEP_ERROR_ABSOLUTE:
362: PetscSNPrintf(ex,sizeof(ex)," ||P(k)x||");
363: break;
364: case PEP_ERROR_RELATIVE:
365: PetscSNPrintf(ex,sizeof(ex),"||P(k)x||/||kx||");
366: break;
367: case PEP_ERROR_BACKWARD:
368: PetscSNPrintf(ex,sizeof(ex)," eta(x,k)");
369: break;
370: }
371: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
372: for (i=0;i<pep->nconv;i++) {
373: PEPGetEigenpair(pep,i,&kr,&ki,NULL,NULL);
374: PEPComputeError(pep,i,etype,&error);
375: #if defined(PETSC_USE_COMPLEX)
376: re = PetscRealPart(kr);
377: im = PetscImaginaryPart(kr);
378: #else
379: re = kr;
380: im = ki;
381: #endif
382: if (im!=0.0) {
383: PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
384: } else {
385: PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
386: }
387: }
388: PetscViewerASCIIPrintf(viewer,"%s",sep);
389: return(0);
390: }
392: static PetscErrorCode PEPErrorView_MATLAB(PEP pep,PEPErrorType etype,PetscViewer viewer)
393: {
395: PetscReal error;
396: PetscInt i;
397: const char *name;
400: PetscObjectGetName((PetscObject)pep,&name);
401: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
402: for (i=0;i<pep->nconv;i++) {
403: PEPComputeError(pep,i,etype,&error);
404: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
405: }
406: PetscViewerASCIIPrintf(viewer,"];\n");
407: return(0);
408: }
410: /*@C
411: PEPErrorView - Displays the errors associated with the computed solution
412: (as well as the eigenvalues).
414: Collective on pep
416: Input Parameters:
417: + pep - the eigensolver context
418: . etype - error type
419: - viewer - optional visualization context
421: Options Database Key:
422: + -pep_error_absolute - print absolute errors of each eigenpair
423: . -pep_error_relative - print relative errors of each eigenpair
424: - -pep_error_backward - print backward errors of each eigenpair
426: Notes:
427: By default, this function checks the error of all eigenpairs and prints
428: the eigenvalues if all of them are below the requested tolerance.
429: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
430: eigenvalues and corresponding errors is printed.
432: Level: intermediate
434: .seealso: PEPSolve(), PEPValuesView(), PEPVectorsView()
435: @*/
436: PetscErrorCode PEPErrorView(PEP pep,PEPErrorType etype,PetscViewer viewer)
437: {
438: PetscBool isascii;
439: PetscViewerFormat format;
440: PetscErrorCode ierr;
444: if (!viewer) {
445: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
446: }
449: PEPCheckSolved(pep,1);
450: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
451: if (!isascii) return(0);
453: PetscViewerGetFormat(viewer,&format);
454: switch (format) {
455: case PETSC_VIEWER_DEFAULT:
456: case PETSC_VIEWER_ASCII_INFO:
457: PEPErrorView_ASCII(pep,etype,viewer);
458: break;
459: case PETSC_VIEWER_ASCII_INFO_DETAIL:
460: PEPErrorView_DETAIL(pep,etype,viewer);
461: break;
462: case PETSC_VIEWER_ASCII_MATLAB:
463: PEPErrorView_MATLAB(pep,etype,viewer);
464: break;
465: default:
466: PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
467: }
468: return(0);
469: }
471: /*@
472: PEPErrorViewFromOptions - Processes command line options to determine if/how
473: the errors of the computed solution are to be viewed.
475: Collective on pep
477: Input Parameter:
478: . pep - the eigensolver context
480: Level: developer
481: @*/
482: PetscErrorCode PEPErrorViewFromOptions(PEP pep)
483: {
484: PetscErrorCode ierr;
485: PetscViewer viewer;
486: PetscBool flg;
487: static PetscBool incall = PETSC_FALSE;
488: PetscViewerFormat format;
491: if (incall) return(0);
492: incall = PETSC_TRUE;
493: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_error_absolute",&viewer,&format,&flg);
494: if (flg) {
495: PetscViewerPushFormat(viewer,format);
496: PEPErrorView(pep,PEP_ERROR_ABSOLUTE,viewer);
497: PetscViewerPopFormat(viewer);
498: PetscViewerDestroy(&viewer);
499: }
500: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_error_relative",&viewer,&format,&flg);
501: if (flg) {
502: PetscViewerPushFormat(viewer,format);
503: PEPErrorView(pep,PEP_ERROR_RELATIVE,viewer);
504: PetscViewerPopFormat(viewer);
505: PetscViewerDestroy(&viewer);
506: }
507: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_error_backward",&viewer,&format,&flg);
508: if (flg) {
509: PetscViewerPushFormat(viewer,format);
510: PEPErrorView(pep,PEP_ERROR_BACKWARD,viewer);
511: PetscViewerPopFormat(viewer);
512: PetscViewerDestroy(&viewer);
513: }
514: incall = PETSC_FALSE;
515: return(0);
516: }
518: static PetscErrorCode PEPValuesView_DRAW(PEP pep,PetscViewer viewer)
519: {
521: PetscDraw draw;
522: PetscDrawSP drawsp;
523: PetscReal re,im;
524: PetscInt i,k;
527: if (!pep->nconv) return(0);
528: PetscViewerDrawGetDraw(viewer,0,&draw);
529: PetscDrawSetTitle(draw,"Computed Eigenvalues");
530: PetscDrawSPCreate(draw,1,&drawsp);
531: for (i=0;i<pep->nconv;i++) {
532: k = pep->perm[i];
533: #if defined(PETSC_USE_COMPLEX)
534: re = PetscRealPart(pep->eigr[k]);
535: im = PetscImaginaryPart(pep->eigr[k]);
536: #else
537: re = pep->eigr[k];
538: im = pep->eigi[k];
539: #endif
540: PetscDrawSPAddPoint(drawsp,&re,&im);
541: }
542: PetscDrawSPDraw(drawsp,PETSC_TRUE);
543: PetscDrawSPSave(drawsp);
544: PetscDrawSPDestroy(&drawsp);
545: return(0);
546: }
548: static PetscErrorCode PEPValuesView_BINARY(PEP pep,PetscViewer viewer)
549: {
550: #if defined(PETSC_HAVE_COMPLEX)
551: PetscInt i,k;
552: PetscComplex *ev;
554: #endif
557: #if defined(PETSC_HAVE_COMPLEX)
558: PetscMalloc1(pep->nconv,&ev);
559: for (i=0;i<pep->nconv;i++) {
560: k = pep->perm[i];
561: #if defined(PETSC_USE_COMPLEX)
562: ev[i] = pep->eigr[k];
563: #else
564: ev[i] = PetscCMPLX(pep->eigr[k],pep->eigi[k]);
565: #endif
566: }
567: PetscViewerBinaryWrite(viewer,ev,pep->nconv,PETSC_COMPLEX);
568: PetscFree(ev);
569: #endif
570: return(0);
571: }
573: #if defined(PETSC_HAVE_HDF5)
574: static PetscErrorCode PEPValuesView_HDF5(PEP pep,PetscViewer viewer)
575: {
577: PetscInt i,k,n,N;
578: PetscMPIInt rank;
579: Vec v;
580: char vname[30];
581: const char *ename;
584: MPI_Comm_rank(PetscObjectComm((PetscObject)pep),&rank);
585: N = pep->nconv;
586: n = rank? 0: N;
587: /* create a vector containing the eigenvalues */
588: VecCreateMPI(PetscObjectComm((PetscObject)pep),n,N,&v);
589: PetscObjectGetName((PetscObject)pep,&ename);
590: PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename);
591: PetscObjectSetName((PetscObject)v,vname);
592: if (!rank) {
593: for (i=0;i<pep->nconv;i++) {
594: k = pep->perm[i];
595: VecSetValue(v,i,pep->eigr[k],INSERT_VALUES);
596: }
597: }
598: VecAssemblyBegin(v);
599: VecAssemblyEnd(v);
600: VecView(v,viewer);
601: #if !defined(PETSC_USE_COMPLEX)
602: /* in real scalars write the imaginary part as a separate vector */
603: PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename);
604: PetscObjectSetName((PetscObject)v,vname);
605: if (!rank) {
606: for (i=0;i<pep->nconv;i++) {
607: k = pep->perm[i];
608: VecSetValue(v,i,pep->eigi[k],INSERT_VALUES);
609: }
610: }
611: VecAssemblyBegin(v);
612: VecAssemblyEnd(v);
613: VecView(v,viewer);
614: #endif
615: VecDestroy(&v);
616: return(0);
617: }
618: #endif
620: static PetscErrorCode PEPValuesView_ASCII(PEP pep,PetscViewer viewer)
621: {
622: PetscInt i,k;
626: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
627: for (i=0;i<pep->nconv;i++) {
628: k = pep->perm[i];
629: PetscViewerASCIIPrintf(viewer," ");
630: SlepcPrintEigenvalueASCII(viewer,pep->eigr[k],pep->eigi[k]);
631: PetscViewerASCIIPrintf(viewer,"\n");
632: }
633: PetscViewerASCIIPrintf(viewer,"\n");
634: return(0);
635: }
637: static PetscErrorCode PEPValuesView_MATLAB(PEP pep,PetscViewer viewer)
638: {
640: PetscInt i,k;
641: PetscReal re,im;
642: const char *name;
645: PetscObjectGetName((PetscObject)pep,&name);
646: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
647: for (i=0;i<pep->nconv;i++) {
648: k = pep->perm[i];
649: #if defined(PETSC_USE_COMPLEX)
650: re = PetscRealPart(pep->eigr[k]);
651: im = PetscImaginaryPart(pep->eigr[k]);
652: #else
653: re = pep->eigr[k];
654: im = pep->eigi[k];
655: #endif
656: if (im!=0.0) {
657: PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
658: } else {
659: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
660: }
661: }
662: PetscViewerASCIIPrintf(viewer,"];\n");
663: return(0);
664: }
666: /*@C
667: PEPValuesView - Displays the computed eigenvalues in a viewer.
669: Collective on pep
671: Input Parameters:
672: + pep - the eigensolver context
673: - viewer - the viewer
675: Options Database Key:
676: . -pep_view_values - print computed eigenvalues
678: Level: intermediate
680: .seealso: PEPSolve(), PEPVectorsView(), PEPErrorView()
681: @*/
682: PetscErrorCode PEPValuesView(PEP pep,PetscViewer viewer)
683: {
684: PetscBool isascii,isdraw,isbinary;
685: PetscViewerFormat format;
686: PetscErrorCode ierr;
687: #if defined(PETSC_HAVE_HDF5)
688: PetscBool ishdf5;
689: #endif
693: if (!viewer) {
694: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
695: }
698: PEPCheckSolved(pep,1);
699: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
700: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
701: #if defined(PETSC_HAVE_HDF5)
702: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
703: #endif
704: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
705: if (isdraw) {
706: PEPValuesView_DRAW(pep,viewer);
707: } else if (isbinary) {
708: PEPValuesView_BINARY(pep,viewer);
709: #if defined(PETSC_HAVE_HDF5)
710: } else if (ishdf5) {
711: PEPValuesView_HDF5(pep,viewer);
712: #endif
713: } else if (isascii) {
714: PetscViewerGetFormat(viewer,&format);
715: switch (format) {
716: case PETSC_VIEWER_DEFAULT:
717: case PETSC_VIEWER_ASCII_INFO:
718: case PETSC_VIEWER_ASCII_INFO_DETAIL:
719: PEPValuesView_ASCII(pep,viewer);
720: break;
721: case PETSC_VIEWER_ASCII_MATLAB:
722: PEPValuesView_MATLAB(pep,viewer);
723: break;
724: default:
725: PetscInfo1(pep,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
726: }
727: }
728: return(0);
729: }
731: /*@
732: PEPValuesViewFromOptions - Processes command line options to determine if/how
733: the computed eigenvalues are to be viewed.
735: Collective on pep
737: Input Parameter:
738: . pep - the eigensolver context
740: Level: developer
741: @*/
742: PetscErrorCode PEPValuesViewFromOptions(PEP pep)
743: {
744: PetscErrorCode ierr;
745: PetscViewer viewer;
746: PetscBool flg;
747: static PetscBool incall = PETSC_FALSE;
748: PetscViewerFormat format;
751: if (incall) return(0);
752: incall = PETSC_TRUE;
753: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_view_values",&viewer,&format,&flg);
754: if (flg) {
755: PetscViewerPushFormat(viewer,format);
756: PEPValuesView(pep,viewer);
757: PetscViewerPopFormat(viewer);
758: PetscViewerDestroy(&viewer);
759: }
760: incall = PETSC_FALSE;
761: return(0);
762: }
764: /*@C
765: PEPVectorsView - Outputs computed eigenvectors to a viewer.
767: Collective on pep
769: Input Parameters:
770: + pep - the eigensolver context
771: - viewer - the viewer
773: Options Database Keys:
774: . -pep_view_vectors - output eigenvectors.
776: Note:
777: If PETSc was configured with real scalars, complex conjugate eigenvectors
778: will be viewed as two separate real vectors, one containing the real part
779: and another one containing the imaginary part.
781: Level: intermediate
783: .seealso: PEPSolve(), PEPValuesView(), PEPErrorView()
784: @*/
785: PetscErrorCode PEPVectorsView(PEP pep,PetscViewer viewer)
786: {
788: PetscInt i,k;
789: Vec xr,xi=NULL;
793: if (!viewer) {
794: PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)pep),&viewer);
795: }
798: PEPCheckSolved(pep,1);
799: if (pep->nconv) {
800: PEPComputeVectors(pep);
801: BVCreateVec(pep->V,&xr);
802: #if !defined(PETSC_USE_COMPLEX)
803: BVCreateVec(pep->V,&xi);
804: #endif
805: for (i=0;i<pep->nconv;i++) {
806: k = pep->perm[i];
807: BV_GetEigenvector(pep->V,k,pep->eigi[k],xr,xi);
808: SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)pep);
809: }
810: VecDestroy(&xr);
811: #if !defined(PETSC_USE_COMPLEX)
812: VecDestroy(&xi);
813: #endif
814: }
815: return(0);
816: }
818: /*@
819: PEPVectorsViewFromOptions - Processes command line options to determine if/how
820: the computed eigenvectors are to be viewed.
822: Collective on pep
824: Input Parameter:
825: . pep - the eigensolver context
827: Level: developer
828: @*/
829: PetscErrorCode PEPVectorsViewFromOptions(PEP pep)
830: {
831: PetscErrorCode ierr;
832: PetscViewer viewer;
833: PetscBool flg = PETSC_FALSE;
834: static PetscBool incall = PETSC_FALSE;
835: PetscViewerFormat format;
838: if (incall) return(0);
839: incall = PETSC_TRUE;
840: PetscOptionsGetViewer(PetscObjectComm((PetscObject)pep),((PetscObject)pep)->options,((PetscObject)pep)->prefix,"-pep_view_vectors",&viewer,&format,&flg);
841: if (flg) {
842: PetscViewerPushFormat(viewer,format);
843: PEPVectorsView(pep,viewer);
844: PetscViewerPopFormat(viewer);
845: PetscViewerDestroy(&viewer);
846: }
847: incall = PETSC_FALSE;
848: return(0);
849: }