Actual source code: epsview.c
slepc-3.18.2 2023-01-26
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, 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: EPS routines related to various viewers
12: */
14: #include <slepc/private/epsimpl.h>
15: #include <slepc/private/bvimpl.h>
16: #include <petscdraw.h>
18: /*@C
19: EPSView - Prints the EPS data structure.
21: Collective on eps
23: Input Parameters:
24: + eps - the eigenproblem solver context
25: - viewer - optional visualization context
27: Options Database Key:
28: . -eps_view - Calls EPSView() at end of EPSSolve()
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
43: .seealso: STView()
44: @*/
45: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
46: {
47: const char *type=NULL,*extr=NULL,*bal=NULL;
48: char str[50];
49: PetscBool isascii,isexternal,istrivial;
52: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
56: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
57: if (isascii) {
58: PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer);
59: PetscViewerASCIIPushTab(viewer);
60: PetscTryTypeMethod(eps,view,viewer);
61: PetscViewerASCIIPopTab(viewer);
62: if (eps->problem_type) {
63: switch (eps->problem_type) {
64: case EPS_HEP: type = SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
65: case EPS_GHEP: type = "generalized " SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
66: case EPS_NHEP: type = "non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
67: case EPS_GNHEP: type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
68: case EPS_PGNHEP: type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem with " SLEPC_STRING_HERMITIAN " positive definite B"; break;
69: case EPS_GHIEP: type = "generalized " SLEPC_STRING_HERMITIAN "-indefinite eigenvalue problem"; break;
70: }
71: } else type = "not yet set";
72: PetscViewerASCIIPrintf(viewer," problem type: %s\n",type);
73: if (eps->extraction) {
74: switch (eps->extraction) {
75: case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
76: case EPS_HARMONIC: extr = "harmonic Ritz"; break;
77: case EPS_HARMONIC_RELATIVE: extr = "relative harmonic Ritz"; break;
78: case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
79: case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
80: case EPS_REFINED: extr = "refined Ritz"; break;
81: case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
82: }
83: PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr);
84: }
85: if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
86: switch (eps->balance) {
87: case EPS_BALANCE_NONE: break;
88: case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
89: case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
90: case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
91: }
92: PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal);
93: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
94: if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) PetscViewerASCIIPrintf(viewer,", with its=%" PetscInt_FMT,eps->balance_its);
95: if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff);
96: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
97: PetscViewerASCIIPrintf(viewer,"\n");
98: }
99: PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: ");
100: SlepcSNPrintfScalar(str,sizeof(str),eps->target,PETSC_FALSE);
101: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
102: if (!eps->which) PetscViewerASCIIPrintf(viewer,"not yet set\n");
103: else switch (eps->which) {
104: case EPS_WHICH_USER:
105: PetscViewerASCIIPrintf(viewer,"user defined\n");
106: break;
107: case EPS_TARGET_MAGNITUDE:
108: PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str);
109: break;
110: case EPS_TARGET_REAL:
111: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str);
112: break;
113: case EPS_TARGET_IMAGINARY:
114: PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str);
115: break;
116: case EPS_LARGEST_MAGNITUDE:
117: PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n");
118: break;
119: case EPS_SMALLEST_MAGNITUDE:
120: PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n");
121: break;
122: case EPS_LARGEST_REAL:
123: PetscViewerASCIIPrintf(viewer,"largest real parts\n");
124: break;
125: case EPS_SMALLEST_REAL:
126: PetscViewerASCIIPrintf(viewer,"smallest real parts\n");
127: break;
128: case EPS_LARGEST_IMAGINARY:
129: PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n");
130: break;
131: case EPS_SMALLEST_IMAGINARY:
132: PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n");
133: break;
134: case EPS_ALL:
135: if (eps->inta || eps->intb) PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb);
136: else PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n");
137: break;
138: }
139: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
140: if (eps->twosided && eps->problem_type!=EPS_HEP && eps->problem_type!=EPS_GHEP) PetscViewerASCIIPrintf(viewer," using two-sided variant (for left eigenvectors)\n");
141: if (eps->purify) PetscViewerASCIIPrintf(viewer," postprocessing eigenvectors with purification\n");
142: if (eps->trueres) PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n");
143: if (eps->trackall) PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n");
144: PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %" PetscInt_FMT "\n",eps->nev);
145: PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",eps->ncv);
146: PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",eps->mpd);
147: PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",eps->max_it);
148: PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)eps->tol);
149: PetscViewerASCIIPrintf(viewer," convergence test: ");
150: PetscViewerASCIIUseTabs(viewer,PETSC_FALSE);
151: switch (eps->conv) {
152: case EPS_CONV_ABS:
153: PetscViewerASCIIPrintf(viewer,"absolute\n");break;
154: case EPS_CONV_REL:
155: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n");break;
156: case EPS_CONV_NORM:
157: PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n");
158: PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)eps->nrma);
159: if (eps->isgeneralized) PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb);
160: PetscViewerASCIIPrintf(viewer,"\n");
161: break;
162: case EPS_CONV_USER:
163: PetscViewerASCIIPrintf(viewer,"user-defined\n");break;
164: }
165: PetscViewerASCIIUseTabs(viewer,PETSC_TRUE);
166: if (eps->nini) PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(eps->nini));
167: if (eps->ninil) PetscViewerASCIIPrintf(viewer," dimension of user-provided left initial space: %" PetscInt_FMT "\n",PetscAbs(eps->ninil));
168: if (eps->nds) PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %" PetscInt_FMT "\n",PetscAbs(eps->nds));
169: } else PetscTryTypeMethod(eps,view,viewer);
170: PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSELEMENTAL,EPSFEAST,EPSPRIMME,EPSSCALAPACK,EPSELPA,EPSEVSL,EPSTRLAN,"");
171: if (!isexternal) {
172: PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO);
173: if (!eps->V) EPSGetBV(eps,&eps->V);
174: BVView(eps->V,viewer);
175: if (eps->rg) {
176: RGIsTrivial(eps->rg,&istrivial);
177: if (!istrivial) RGView(eps->rg,viewer);
178: }
179: if (eps->useds) {
180: if (!eps->ds) EPSGetDS(eps,&eps->ds);
181: DSView(eps->ds,viewer);
182: }
183: PetscViewerPopFormat(viewer);
184: }
185: if (!eps->st) EPSGetST(eps,&eps->st);
186: STView(eps->st,viewer);
187: return 0;
188: }
190: /*@C
191: EPSViewFromOptions - View from options
193: Collective on EPS
195: Input Parameters:
196: + eps - the eigensolver context
197: . obj - optional object
198: - name - command line option
200: Level: intermediate
202: .seealso: EPSView(), EPSCreate()
203: @*/
204: PetscErrorCode EPSViewFromOptions(EPS eps,PetscObject obj,const char name[])
205: {
207: PetscObjectViewFromOptions((PetscObject)eps,obj,name);
208: return 0;
209: }
211: /*@C
212: EPSConvergedReasonView - Displays the reason an EPS solve converged or diverged.
214: Collective on eps
216: Input Parameters:
217: + eps - the eigensolver context
218: - viewer - the viewer to display the reason
220: Options Database Keys:
221: . -eps_converged_reason - print reason for convergence, and number of iterations
223: Note:
224: To change the format of the output call PetscViewerPushFormat(viewer,format) before
225: this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
226: display a reason if it fails. The latter can be set in the command line with
227: -eps_converged_reason ::failed
229: Level: intermediate
231: .seealso: EPSSetConvergenceTest(), EPSSetTolerances(), EPSGetIterationNumber(), EPSConvergedReasonViewFromOptions()
232: @*/
233: PetscErrorCode EPSConvergedReasonView(EPS eps,PetscViewer viewer)
234: {
235: PetscBool isAscii;
236: PetscViewerFormat format;
238: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
239: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii);
240: if (isAscii) {
241: PetscViewerGetFormat(viewer,&format);
242: PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel);
243: if (eps->reason > 0 && format != PETSC_VIEWER_FAILED) PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve converged (%" PetscInt_FMT " eigenpair%s) due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",eps->nconv,(eps->nconv>1)?"s":"",EPSConvergedReasons[eps->reason],eps->its);
244: else if (eps->reason <= 0) PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve did not converge due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",EPSConvergedReasons[eps->reason],eps->its);
245: PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel);
246: }
247: return 0;
248: }
250: /*@
251: EPSConvergedReasonViewFromOptions - Processes command line options to determine if/how
252: the EPS converged reason is to be viewed.
254: Collective on eps
256: Input Parameter:
257: . eps - the eigensolver context
259: Level: developer
261: .seealso: EPSConvergedReasonView()
262: @*/
263: PetscErrorCode EPSConvergedReasonViewFromOptions(EPS eps)
264: {
265: PetscViewer viewer;
266: PetscBool flg;
267: static PetscBool incall = PETSC_FALSE;
268: PetscViewerFormat format;
270: if (incall) return 0;
271: incall = PETSC_TRUE;
272: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg);
273: if (flg) {
274: PetscViewerPushFormat(viewer,format);
275: EPSConvergedReasonView(eps,viewer);
276: PetscViewerPopFormat(viewer);
277: PetscViewerDestroy(&viewer);
278: }
279: incall = PETSC_FALSE;
280: return 0;
281: }
283: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
284: {
285: PetscReal error;
286: PetscInt i,j,k,nvals;
288: nvals = (eps->which==EPS_ALL)? eps->nconv: eps->nev;
289: if (eps->which!=EPS_ALL && eps->nconv<nvals) {
290: PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " eigenvalues converged\n\n",eps->nev);
291: return 0;
292: }
293: if (eps->which==EPS_ALL && !nvals) {
294: PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n");
295: return 0;
296: }
297: for (i=0;i<nvals;i++) {
298: EPSComputeError(eps,i,etype,&error);
299: if (error>=5.0*eps->tol) {
300: PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",nvals);
301: return 0;
302: }
303: }
304: if (eps->which==EPS_ALL) PetscViewerASCIIPrintf(viewer," Found %" PetscInt_FMT " eigenvalues, all of them computed up to the required tolerance:",nvals);
305: else PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:");
306: for (i=0;i<=(nvals-1)/8;i++) {
307: PetscViewerASCIIPrintf(viewer,"\n ");
308: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
309: k = eps->perm[8*i+j];
310: SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]);
311: if (8*i+j+1<nvals) PetscViewerASCIIPrintf(viewer,", ");
312: }
313: }
314: PetscViewerASCIIPrintf(viewer,"\n\n");
315: return 0;
316: }
318: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
319: {
320: PetscReal error,re,im;
321: PetscScalar kr,ki;
322: PetscInt i;
323: char ex[30],sep[]=" ---------------------- --------------------\n";
325: if (!eps->nconv) return 0;
326: switch (etype) {
327: case EPS_ERROR_ABSOLUTE:
328: PetscSNPrintf(ex,sizeof(ex)," ||Ax-k%sx||",eps->isgeneralized?"B":"");
329: break;
330: case EPS_ERROR_RELATIVE:
331: PetscSNPrintf(ex,sizeof(ex),"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":"");
332: break;
333: case EPS_ERROR_BACKWARD:
334: PetscSNPrintf(ex,sizeof(ex)," eta(x,k)");
335: break;
336: }
337: PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep);
338: for (i=0;i<eps->nconv;i++) {
339: EPSGetEigenpair(eps,i,&kr,&ki,NULL,NULL);
340: EPSComputeError(eps,i,etype,&error);
341: #if defined(PETSC_USE_COMPLEX)
342: re = PetscRealPart(kr);
343: im = PetscImaginaryPart(kr);
344: #else
345: re = kr;
346: im = ki;
347: #endif
348: if (im!=0.0) PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error);
349: else PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error);
350: }
351: PetscViewerASCIIPrintf(viewer,"%s",sep);
352: return 0;
353: }
355: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
356: {
357: PetscReal error;
358: PetscInt i;
359: const char *name;
361: PetscObjectGetName((PetscObject)eps,&name);
362: PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name);
363: for (i=0;i<eps->nconv;i++) {
364: EPSComputeError(eps,i,etype,&error);
365: PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error);
366: }
367: PetscViewerASCIIPrintf(viewer,"];\n");
368: return 0;
369: }
371: /*@C
372: EPSErrorView - Displays the errors associated with the computed solution
373: (as well as the eigenvalues).
375: Collective on eps
377: Input Parameters:
378: + eps - the eigensolver context
379: . etype - error type
380: - viewer - optional visualization context
382: Options Database Keys:
383: + -eps_error_absolute - print absolute errors of each eigenpair
384: . -eps_error_relative - print relative errors of each eigenpair
385: - -eps_error_backward - print backward errors of each eigenpair
387: Notes:
388: By default, this function checks the error of all eigenpairs and prints
389: the eigenvalues if all of them are below the requested tolerance.
390: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
391: eigenvalues and corresponding errors is printed.
393: Level: intermediate
395: .seealso: EPSSolve(), EPSValuesView(), EPSVectorsView()
396: @*/
397: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
398: {
399: PetscBool isascii;
400: PetscViewerFormat format;
403: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
406: EPSCheckSolved(eps,1);
407: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
408: if (!isascii) return 0;
410: PetscViewerGetFormat(viewer,&format);
411: switch (format) {
412: case PETSC_VIEWER_DEFAULT:
413: case PETSC_VIEWER_ASCII_INFO:
414: EPSErrorView_ASCII(eps,etype,viewer);
415: break;
416: case PETSC_VIEWER_ASCII_INFO_DETAIL:
417: EPSErrorView_DETAIL(eps,etype,viewer);
418: break;
419: case PETSC_VIEWER_ASCII_MATLAB:
420: EPSErrorView_MATLAB(eps,etype,viewer);
421: break;
422: default:
423: PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
424: }
425: return 0;
426: }
428: /*@
429: EPSErrorViewFromOptions - Processes command line options to determine if/how
430: the errors of the computed solution are to be viewed.
432: Collective on eps
434: Input Parameter:
435: . eps - the eigensolver context
437: Level: developer
439: .seealso: EPSErrorView()
440: @*/
441: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
442: {
443: PetscViewer viewer;
444: PetscBool flg;
445: static PetscBool incall = PETSC_FALSE;
446: PetscViewerFormat format;
448: if (incall) return 0;
449: incall = PETSC_TRUE;
450: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg);
451: if (flg) {
452: PetscViewerPushFormat(viewer,format);
453: EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer);
454: PetscViewerPopFormat(viewer);
455: PetscViewerDestroy(&viewer);
456: }
457: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg);
458: if (flg) {
459: PetscViewerPushFormat(viewer,format);
460: EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer);
461: PetscViewerPopFormat(viewer);
462: PetscViewerDestroy(&viewer);
463: }
464: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg);
465: if (flg) {
466: PetscViewerPushFormat(viewer,format);
467: EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer);
468: PetscViewerPopFormat(viewer);
469: PetscViewerDestroy(&viewer);
470: }
471: incall = PETSC_FALSE;
472: return 0;
473: }
475: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
476: {
477: PetscDraw draw;
478: PetscDrawSP drawsp;
479: PetscReal re,im;
480: PetscInt i,k;
482: if (!eps->nconv) return 0;
483: PetscViewerDrawGetDraw(viewer,0,&draw);
484: PetscDrawSetTitle(draw,"Computed Eigenvalues");
485: PetscDrawSPCreate(draw,1,&drawsp);
486: for (i=0;i<eps->nconv;i++) {
487: k = eps->perm[i];
488: #if defined(PETSC_USE_COMPLEX)
489: re = PetscRealPart(eps->eigr[k]);
490: im = PetscImaginaryPart(eps->eigr[k]);
491: #else
492: re = eps->eigr[k];
493: im = eps->eigi[k];
494: #endif
495: PetscDrawSPAddPoint(drawsp,&re,&im);
496: }
497: PetscDrawSPDraw(drawsp,PETSC_TRUE);
498: PetscDrawSPSave(drawsp);
499: PetscDrawSPDestroy(&drawsp);
500: return 0;
501: }
503: static PetscErrorCode EPSValuesView_BINARY(EPS eps,PetscViewer viewer)
504: {
505: #if defined(PETSC_HAVE_COMPLEX)
506: PetscInt i,k;
507: PetscComplex *ev;
508: #endif
510: #if defined(PETSC_HAVE_COMPLEX)
511: PetscMalloc1(eps->nconv,&ev);
512: for (i=0;i<eps->nconv;i++) {
513: k = eps->perm[i];
514: #if defined(PETSC_USE_COMPLEX)
515: ev[i] = eps->eigr[k];
516: #else
517: ev[i] = PetscCMPLX(eps->eigr[k],eps->eigi[k]);
518: #endif
519: }
520: PetscViewerBinaryWrite(viewer,ev,eps->nconv,PETSC_COMPLEX);
521: PetscFree(ev);
522: #endif
523: return 0;
524: }
526: #if defined(PETSC_HAVE_HDF5)
527: static PetscErrorCode EPSValuesView_HDF5(EPS eps,PetscViewer viewer)
528: {
529: PetscInt i,k,n,N;
530: PetscMPIInt rank;
531: Vec v;
532: char vname[30];
533: const char *ename;
535: MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank);
536: N = eps->nconv;
537: n = rank? 0: N;
538: /* create a vector containing the eigenvalues */
539: VecCreateMPI(PetscObjectComm((PetscObject)eps),n,N,&v);
540: PetscObjectGetName((PetscObject)eps,&ename);
541: PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename);
542: PetscObjectSetName((PetscObject)v,vname);
543: if (!rank) {
544: for (i=0;i<eps->nconv;i++) {
545: k = eps->perm[i];
546: VecSetValue(v,i,eps->eigr[k],INSERT_VALUES);
547: }
548: }
549: VecAssemblyBegin(v);
550: VecAssemblyEnd(v);
551: VecView(v,viewer);
552: #if !defined(PETSC_USE_COMPLEX)
553: /* in real scalars write the imaginary part as a separate vector */
554: PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename);
555: PetscObjectSetName((PetscObject)v,vname);
556: if (!rank) {
557: for (i=0;i<eps->nconv;i++) {
558: k = eps->perm[i];
559: VecSetValue(v,i,eps->eigi[k],INSERT_VALUES);
560: }
561: }
562: VecAssemblyBegin(v);
563: VecAssemblyEnd(v);
564: VecView(v,viewer);
565: #endif
566: VecDestroy(&v);
567: return 0;
568: }
569: #endif
571: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
572: {
573: PetscInt i,k;
575: PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n");
576: for (i=0;i<eps->nconv;i++) {
577: k = eps->perm[i];
578: PetscViewerASCIIPrintf(viewer," ");
579: SlepcPrintEigenvalueASCII(viewer,eps->eigr[k],eps->eigi[k]);
580: PetscViewerASCIIPrintf(viewer,"\n");
581: }
582: PetscViewerASCIIPrintf(viewer,"\n");
583: return 0;
584: }
586: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
587: {
588: PetscInt i,k;
589: PetscReal re,im;
590: const char *name;
592: PetscObjectGetName((PetscObject)eps,&name);
593: PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name);
594: for (i=0;i<eps->nconv;i++) {
595: k = eps->perm[i];
596: #if defined(PETSC_USE_COMPLEX)
597: re = PetscRealPart(eps->eigr[k]);
598: im = PetscImaginaryPart(eps->eigr[k]);
599: #else
600: re = eps->eigr[k];
601: im = eps->eigi[k];
602: #endif
603: if (im!=0.0) PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im);
604: else PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re);
605: }
606: PetscViewerASCIIPrintf(viewer,"];\n");
607: return 0;
608: }
610: /*@C
611: EPSValuesView - Displays the computed eigenvalues in a viewer.
613: Collective on eps
615: Input Parameters:
616: + eps - the eigensolver context
617: - viewer - the viewer
619: Options Database Key:
620: . -eps_view_values - print computed eigenvalues
622: Level: intermediate
624: .seealso: EPSSolve(), EPSVectorsView(), EPSErrorView()
625: @*/
626: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
627: {
628: PetscBool isascii,isdraw,isbinary;
629: PetscViewerFormat format;
630: #if defined(PETSC_HAVE_HDF5)
631: PetscBool ishdf5;
632: #endif
635: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
638: EPSCheckSolved(eps,1);
639: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw);
640: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
641: #if defined(PETSC_HAVE_HDF5)
642: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5);
643: #endif
644: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
645: if (isdraw) EPSValuesView_DRAW(eps,viewer);
646: else if (isbinary) EPSValuesView_BINARY(eps,viewer);
647: #if defined(PETSC_HAVE_HDF5)
648: else if (ishdf5) EPSValuesView_HDF5(eps,viewer);
649: #endif
650: else if (isascii) {
651: PetscViewerGetFormat(viewer,&format);
652: switch (format) {
653: case PETSC_VIEWER_DEFAULT:
654: case PETSC_VIEWER_ASCII_INFO:
655: case PETSC_VIEWER_ASCII_INFO_DETAIL:
656: EPSValuesView_ASCII(eps,viewer);
657: break;
658: case PETSC_VIEWER_ASCII_MATLAB:
659: EPSValuesView_MATLAB(eps,viewer);
660: break;
661: default:
662: PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]);
663: }
664: }
665: return 0;
666: }
668: /*@
669: EPSValuesViewFromOptions - Processes command line options to determine if/how
670: the computed eigenvalues are to be viewed.
672: Collective on eps
674: Input Parameters:
675: . eps - the eigensolver context
677: Level: developer
679: .seealso: EPSValuesView()
680: @*/
681: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
682: {
683: PetscViewer viewer;
684: PetscBool flg;
685: static PetscBool incall = PETSC_FALSE;
686: PetscViewerFormat format;
688: if (incall) return 0;
689: incall = PETSC_TRUE;
690: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg);
691: if (flg) {
692: PetscViewerPushFormat(viewer,format);
693: EPSValuesView(eps,viewer);
694: PetscViewerPopFormat(viewer);
695: PetscViewerDestroy(&viewer);
696: }
697: incall = PETSC_FALSE;
698: return 0;
699: }
701: /*@C
702: EPSVectorsView - Outputs computed eigenvectors to a viewer.
704: Collective on eps
706: Input Parameters:
707: + eps - the eigensolver context
708: - viewer - the viewer
710: Options Database Key:
711: . -eps_view_vectors - output eigenvectors.
713: Notes:
714: If PETSc was configured with real scalars, complex conjugate eigenvectors
715: will be viewed as two separate real vectors, one containing the real part
716: and another one containing the imaginary part.
718: If left eigenvectors were computed with a two-sided eigensolver, the right
719: and left eigenvectors are interleaved, that is, the vectors are output in
720: the following order X0, Y0, X1, Y1, X2, Y2, ...
722: Level: intermediate
724: .seealso: EPSSolve(), EPSValuesView(), EPSErrorView()
725: @*/
726: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
727: {
728: PetscInt i,k;
729: Vec xr,xi=NULL;
732: if (!viewer) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer);
735: EPSCheckSolved(eps,1);
736: if (eps->nconv) {
737: EPSComputeVectors(eps);
738: BVCreateVec(eps->V,&xr);
739: #if !defined(PETSC_USE_COMPLEX)
740: BVCreateVec(eps->V,&xi);
741: #endif
742: for (i=0;i<eps->nconv;i++) {
743: k = eps->perm[i];
744: BV_GetEigenvector(eps->V,k,eps->eigi[k],xr,xi);
745: SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)eps);
746: if (eps->twosided) {
747: BV_GetEigenvector(eps->W,k,eps->eigi[k],xr,xi);
748: SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)eps);
749: }
750: }
751: VecDestroy(&xr);
752: #if !defined(PETSC_USE_COMPLEX)
753: VecDestroy(&xi);
754: #endif
755: }
756: return 0;
757: }
759: /*@
760: EPSVectorsViewFromOptions - Processes command line options to determine if/how
761: the computed eigenvectors are to be viewed.
763: Collective on eps
765: Input Parameter:
766: . eps - the eigensolver context
768: Level: developer
770: .seealso: EPSVectorsView()
771: @*/
772: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
773: {
774: PetscViewer viewer;
775: PetscBool flg = PETSC_FALSE;
776: static PetscBool incall = PETSC_FALSE;
777: PetscViewerFormat format;
779: if (incall) return 0;
780: incall = PETSC_TRUE;
781: PetscOptionsGetViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg);
782: if (flg) {
783: PetscViewerPushFormat(viewer,format);
784: EPSVectorsView(eps,viewer);
785: PetscViewerPopFormat(viewer);
786: PetscViewerDestroy(&viewer);
787: }
788: incall = PETSC_FALSE;
789: return 0;
790: }