Actual source code: epsview.c
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: /*@
19: EPSView - Prints the `EPS` data structure.
21: Collective
23: Input Parameters:
24: + eps - the linear eigensolver context
25: - viewer - optional visualization context
27: Options Database Key:
28: . -eps_view - calls `EPSView()` at end of `EPSSolve()`
30: Notes:
31: The available visualization contexts include
32: + `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
33: - `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard output where only the
34: first process opens the file; all other processes send their data to the
35: first one to print
37: The user can open an alternative visualization context with `PetscViewerASCIIOpen()`
38: to output to a specified file.
40: Level: beginner
42: .seealso: [](ch:eps), `EPSCreate()`, `EPSViewFromOptions()`, `STView()`
43: @*/
44: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
45: {
46: const char *type=NULL,*extr=NULL,*bal=NULL;
47: char str[50];
48: PetscBool isascii,isexternal,istrivial,isstruct=PETSC_FALSE,reduced=PETSC_FALSE,flg;
49: Mat A;
51: PetscFunctionBegin;
53: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
55: PetscCheckSameComm(eps,1,viewer,2);
57: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
58: if (isascii) {
59: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer));
60: PetscCall(PetscViewerASCIIPushTab(viewer));
61: PetscTryTypeMethod(eps,view,viewer);
62: PetscCall(PetscViewerASCIIPopTab(viewer));
63: if (eps->problem_type) {
64: switch (eps->problem_type) {
65: case EPS_HEP: type = SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
66: case EPS_GHEP: type = "generalized " SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
67: case EPS_NHEP: type = "non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
68: case EPS_GNHEP: type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
69: case EPS_PGNHEP: type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem with " SLEPC_STRING_HERMITIAN " positive definite B"; break;
70: case EPS_GHIEP: type = "generalized " SLEPC_STRING_HERMITIAN "-indefinite eigenvalue problem"; break;
71: case EPS_BSE: type = "structured Bethe-Salpeter eigenvalue problem"; break;
72: case EPS_HAMILT: type = "structured Hamiltonian eigenvalue problem"; break;
73: case EPS_LREP: type = "structured Linear Response eigenvalue problem"; break;
74: }
75: } else type = "not yet set";
76: PetscCall(PetscViewerASCIIPrintf(viewer," problem type: %s\n",type));
77: PetscCall(EPSGetOperators(eps,&A,NULL));
78: if (A) PetscCall(SlepcCheckMatStruct(A,0,&isstruct));
79: if (isstruct) {
80: PetscCall(SlepcCheckMatStruct(A,SLEPC_MAT_STRUCT_BSE,&flg));
81: if (flg) PetscCall(PetscViewerASCIIPrintf(viewer," matrix A has a Bethe-Salpeter structure\n"));
82: PetscCall(SlepcCheckMatStruct(A,SLEPC_MAT_STRUCT_HAMILT,&flg));
83: if (flg) PetscCall(PetscViewerASCIIPrintf(viewer," matrix A has a Hamiltonian structure\n"));
84: PetscCall(SlepcCheckMatStruct(A,SLEPC_MAT_STRUCT_LREP,&flg));
85: if (flg) {
86: PetscCall(SlepcCheckMatLREPReduced(A,&reduced));
87: PetscCall(PetscViewerASCIIPrintf(viewer," matrix A has a%s Linear Response structure\n",reduced?" reduced":""));
88: }
89: }
90: if (eps->extraction) {
91: switch (eps->extraction) {
92: case EPS_RITZ: extr = "Rayleigh-Ritz"; break;
93: case EPS_HARMONIC: extr = "harmonic Ritz"; break;
94: case EPS_HARMONIC_RELATIVE: extr = "relative harmonic Ritz"; break;
95: case EPS_HARMONIC_RIGHT: extr = "right harmonic Ritz"; break;
96: case EPS_HARMONIC_LARGEST: extr = "largest harmonic Ritz"; break;
97: case EPS_REFINED: extr = "refined Ritz"; break;
98: case EPS_REFINED_HARMONIC: extr = "refined harmonic Ritz"; break;
99: }
100: PetscCall(PetscViewerASCIIPrintf(viewer," extraction type: %s\n",extr));
101: }
102: if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
103: switch (eps->balance) {
104: case EPS_BALANCE_NONE: break;
105: case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
106: case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
107: case EPS_BALANCE_USER: bal = "user-defined matrix"; break;
108: }
109: PetscCall(PetscViewerASCIIPrintf(viewer," balancing enabled: %s",bal));
110: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
111: if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) PetscCall(PetscViewerASCIIPrintf(viewer,", with its=%" PetscInt_FMT,eps->balance_its));
112: if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff));
113: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
114: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
115: }
116: PetscCall(PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: "));
117: PetscCall(SlepcSNPrintfScalar(str,sizeof(str),eps->target,PETSC_FALSE));
118: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
119: if (!eps->which) PetscCall(PetscViewerASCIIPrintf(viewer,"not yet set\n"));
120: else switch (eps->which) {
121: case EPS_WHICH_USER:
122: PetscCall(PetscViewerASCIIPrintf(viewer,"user defined\n"));
123: break;
124: case EPS_TARGET_MAGNITUDE:
125: PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str));
126: break;
127: case EPS_TARGET_REAL:
128: PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str));
129: break;
130: case EPS_TARGET_IMAGINARY:
131: PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str));
132: break;
133: case EPS_LARGEST_MAGNITUDE:
134: PetscCall(PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n"));
135: break;
136: case EPS_SMALLEST_MAGNITUDE:
137: PetscCall(PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n"));
138: break;
139: case EPS_LARGEST_REAL:
140: PetscCall(PetscViewerASCIIPrintf(viewer,"largest real parts\n"));
141: break;
142: case EPS_SMALLEST_REAL:
143: PetscCall(PetscViewerASCIIPrintf(viewer,"smallest real parts\n"));
144: break;
145: case EPS_LARGEST_IMAGINARY:
146: PetscCall(PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n"));
147: break;
148: case EPS_SMALLEST_IMAGINARY:
149: PetscCall(PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n"));
150: break;
151: case EPS_ALL:
152: if (eps->inta || eps->intb) PetscCall(PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb));
153: else PetscCall(PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n"));
154: break;
155: }
156: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
157: if (eps->twosided && eps->problem_type!=EPS_HEP && eps->problem_type!=EPS_GHEP) PetscCall(PetscViewerASCIIPrintf(viewer," using two-sided variant (for left eigenvectors)\n"));
158: if (eps->purify) PetscCall(PetscViewerASCIIPrintf(viewer," postprocessing eigenvectors with purification\n"));
159: if (eps->trueres) PetscCall(PetscViewerASCIIPrintf(viewer," computing true residuals explicitly\n"));
160: if (eps->trackall) PetscCall(PetscViewerASCIIPrintf(viewer," computing all residuals (for tracking convergence)\n"));
161: if (eps->stop==EPS_STOP_THRESHOLD) PetscCall(PetscViewerASCIIPrintf(viewer," computing eigenvalues %s the threshold: %g%s\n",(eps->which==EPS_SMALLEST_MAGNITUDE||eps->which==EPS_SMALLEST_REAL)?"below":"above",(double)eps->thres,eps->threlative?" (relative)":""));
162: if (eps->nev) PetscCall(PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %" PetscInt_FMT "\n",eps->nev));
163: PetscCall(PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",eps->ncv));
164: PetscCall(PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",eps->mpd));
165: PetscCall(PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",eps->max_it));
166: PetscCall(PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)eps->tol));
167: PetscCall(PetscViewerASCIIPrintf(viewer," convergence test: "));
168: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
169: switch (eps->conv) {
170: case EPS_CONV_ABS:
171: PetscCall(PetscViewerASCIIPrintf(viewer,"absolute\n"));break;
172: case EPS_CONV_REL:
173: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n"));break;
174: case EPS_CONV_NORM:
175: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n"));
176: PetscCall(PetscViewerASCIIPrintf(viewer," computed matrix norms: norm(A)=%g",(double)eps->nrma));
177: if (eps->isgeneralized) PetscCall(PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb));
178: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
179: break;
180: case EPS_CONV_USER:
181: PetscCall(PetscViewerASCIIPrintf(viewer,"user-defined\n"));break;
182: }
183: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
184: if (eps->nini) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(eps->nini)));
185: if (eps->ninil) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided left initial space: %" PetscInt_FMT "\n",PetscAbs(eps->ninil)));
186: if (eps->nds) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided deflation space: %" PetscInt_FMT "\n",PetscAbs(eps->nds)));
187: } else PetscTryTypeMethod(eps,view,viewer);
188: PetscCall(PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSELEMENTAL,EPSFEAST,EPSPRIMME,EPSSCALAPACK,EPSELPA,EPSEVSL,""));
189: if (!isexternal) {
190: PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
191: if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
192: PetscCall(BVView(eps->V,viewer));
193: if (eps->rg) {
194: PetscCall(RGIsTrivial(eps->rg,&istrivial));
195: if (!istrivial) PetscCall(RGView(eps->rg,viewer));
196: }
197: if (eps->useds) {
198: if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
199: PetscCall(DSView(eps->ds,viewer));
200: }
201: PetscCall(PetscViewerPopFormat(viewer));
202: }
203: if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
204: PetscCall(STView(eps->st,viewer));
205: PetscFunctionReturn(PETSC_SUCCESS);
206: }
208: /*@
209: EPSViewFromOptions - View (print) an `EPS` object based on values in the options database.
211: Collective
213: Input Parameters:
214: + eps - the linear eigensolver context
215: . obj - optional object that provides the options prefix used to query the options database
216: - name - command line option
218: Level: intermediate
220: .seealso: [](ch:eps), `EPSView()`, `EPSCreate()`, `PetscObjectViewFromOptions()`
221: @*/
222: PetscErrorCode EPSViewFromOptions(EPS eps,PetscObject obj,const char name[])
223: {
224: PetscFunctionBegin;
226: PetscCall(PetscObjectViewFromOptions((PetscObject)eps,obj,name));
227: PetscFunctionReturn(PETSC_SUCCESS);
228: }
230: /*@
231: EPSConvergedReasonView - Displays the reason an `EPS` solve converged or diverged.
233: Collective
235: Input Parameters:
236: + eps - the linear eigensolver context
237: - viewer - the viewer to display the reason
239: Options Database Key:
240: . -eps_converged_reason - print reason for convergence/divergence, and number of iterations
242: Notes:
243: Use `EPSConvergedReasonViewFromOptions()` to display the reason based on values
244: in the options database.
246: To change the format of the output call `PetscViewerPushFormat()` before this
247: call. Use `PETSC_VIEWER_DEFAULT` for the default, or `PETSC_VIEWER_FAILED` to only
248: display a reason if it fails. The latter can be set in the command line with
249: `-eps_converged_reason ::failed`.
251: Level: intermediate
253: .seealso: [](ch:eps), `EPSSetConvergenceTest()`, `EPSSetTolerances()`, `EPSGetIterationNumber()`, `EPSConvergedReasonViewFromOptions()`
254: @*/
255: PetscErrorCode EPSConvergedReasonView(EPS eps,PetscViewer viewer)
256: {
257: PetscBool isAscii;
258: PetscViewerFormat format;
259: PetscInt nconv;
261: PetscFunctionBegin;
262: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
263: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
264: if (isAscii) {
265: PetscCall(PetscViewerGetFormat(viewer,&format));
266: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
267: PetscCall(EPS_GetActualConverged(eps,&nconv));
268: if (eps->reason > 0 && format != PETSC_VIEWER_FAILED) PetscCall(PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve converged (%" PetscInt_FMT " eigenpair%s) due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",nconv,(nconv>1)?"s":"",EPSConvergedReasons[eps->reason],eps->its));
269: else if (eps->reason <= 0) PetscCall(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));
270: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
271: }
272: PetscFunctionReturn(PETSC_SUCCESS);
273: }
275: /*@
276: EPSConvergedReasonViewFromOptions - Processes command line options to determine if/how
277: the `EPS` converged reason is to be viewed.
279: Collective
281: Input Parameter:
282: . eps - the linear eigensolver context
284: Level: intermediate
286: .seealso: [](ch:eps), `EPSConvergedReasonView()`
287: @*/
288: PetscErrorCode EPSConvergedReasonViewFromOptions(EPS eps)
289: {
290: PetscViewer viewer;
291: PetscBool flg;
292: static PetscBool incall = PETSC_FALSE;
293: PetscViewerFormat format;
295: PetscFunctionBegin;
296: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
297: incall = PETSC_TRUE;
298: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg));
299: if (flg) {
300: PetscCall(PetscViewerPushFormat(viewer,format));
301: PetscCall(EPSConvergedReasonView(eps,viewer));
302: PetscCall(PetscViewerPopFormat(viewer));
303: PetscCall(PetscViewerDestroy(&viewer));
304: }
305: incall = PETSC_FALSE;
306: PetscFunctionReturn(PETSC_SUCCESS);
307: }
309: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
310: {
311: PetscReal error;
312: PetscScalar kr,ki;
313: PetscInt i,j,nvals,nconv;
314: PetscBool requestednev=PETSC_TRUE; /* user has specified nev */
316: PetscFunctionBegin;
317: PetscCall(EPS_GetActualConverged(eps,&nconv));
318: if (eps->which==EPS_ALL) {
319: requestednev = PETSC_FALSE;
320: if (eps->stop==EPS_STOP_THRESHOLD && eps->nev>0) requestednev = PETSC_TRUE;
321: } else if (eps->stop==EPS_STOP_THRESHOLD && eps->nev==0) requestednev = PETSC_FALSE;
322: nvals = requestednev? eps->nev: nconv;
324: if (requestednev && nconv<nvals) {
325: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " eigenvalues converged\n\n",eps->nev));
326: PetscFunctionReturn(PETSC_SUCCESS);
327: }
328: if (!requestednev && !nvals) {
329: PetscCall(PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n"));
330: PetscFunctionReturn(PETSC_SUCCESS);
331: }
332: for (i=0;i<nvals;i++) {
333: PetscCall(EPSComputeError(eps,i,etype,&error));
334: if (error>=5.0*eps->tol) {
335: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",nvals));
336: PetscFunctionReturn(PETSC_SUCCESS);
337: }
338: }
339: if (!requestednev) PetscCall(PetscViewerASCIIPrintf(viewer," Found %" PetscInt_FMT " eigenvalues, all of them computed up to the required tolerance:",nvals));
340: else PetscCall(PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:"));
341: for (i=0;i<=(nvals-1)/8;i++) {
342: PetscCall(PetscViewerASCIIPrintf(viewer,"\n "));
343: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
344: PetscCall(EPSGetEigenvalue(eps,8*i+j,&kr,&ki));
345: PetscCall(SlepcPrintEigenvalueASCII(viewer,kr,ki));
346: if (8*i+j+1<nvals) PetscCall(PetscViewerASCIIPrintf(viewer,", "));
347: }
348: }
349: PetscCall(PetscViewerASCIIPrintf(viewer,"\n\n"));
350: PetscFunctionReturn(PETSC_SUCCESS);
351: }
353: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
354: {
355: PetscReal error,re,im;
356: PetscScalar kr,ki;
357: PetscInt i,nconv;
358: char ex[30],sep[]=" ---------------------- --------------------\n";
360: PetscFunctionBegin;
361: if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
362: switch (etype) {
363: case EPS_ERROR_ABSOLUTE:
364: PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||Ax-k%sx||",eps->isgeneralized?"B":""));
365: break;
366: case EPS_ERROR_RELATIVE:
367: PetscCall(PetscSNPrintf(ex,sizeof(ex),"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":""));
368: break;
369: case EPS_ERROR_BACKWARD:
370: PetscCall(PetscSNPrintf(ex,sizeof(ex)," eta(x,k)"));
371: break;
372: }
373: PetscCall(PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep));
374: PetscCall(EPS_GetActualConverged(eps,&nconv));
375: for (i=0;i<nconv;i++) {
376: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
377: PetscCall(EPSComputeError(eps,i,etype,&error));
378: #if defined(PETSC_USE_COMPLEX)
379: re = PetscRealPart(kr);
380: im = PetscImaginaryPart(kr);
381: #else
382: re = kr;
383: im = ki;
384: #endif
385: if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error));
386: else PetscCall(PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error));
387: }
388: PetscCall(PetscViewerASCIIPrintf(viewer,"%s",sep));
389: PetscFunctionReturn(PETSC_SUCCESS);
390: }
392: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
393: {
394: PetscReal error;
395: PetscInt i,nconv;
396: const char *name;
398: PetscFunctionBegin;
399: PetscCall(PetscObjectGetName((PetscObject)eps,&name));
400: PetscCall(PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name));
401: PetscCall(EPS_GetActualConverged(eps,&nconv));
402: for (i=0;i<nconv;i++) {
403: PetscCall(EPSComputeError(eps,i,etype,&error));
404: PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error));
405: }
406: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
407: PetscFunctionReturn(PETSC_SUCCESS);
408: }
410: /*@
411: EPSErrorView - Displays the errors associated with the computed solution
412: (as well as the eigenvalues).
414: Collective
416: Input Parameters:
417: + eps - the linear eigensolver context
418: . etype - error type
419: - viewer - optional visualization context
421: Options Database Keys:
422: + -eps_error_absolute - print absolute errors of each eigenpair
423: . -eps_error_relative - print relative errors of each eigenpair
424: - -eps_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: All the command-line options listed above admit an optional argument
433: specifying the viewer type and options. For instance, use
434: `-eps_error_relative :myerr.m:ascii_matlab` to save the errors in a file
435: that can be executed in Matlab.
437: Level: intermediate
439: .seealso: [](ch:eps), `EPSSolve()`, `EPSValuesView()`, `EPSVectorsView()`
440: @*/
441: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
442: {
443: PetscBool isascii;
444: PetscViewerFormat format;
446: PetscFunctionBegin;
448: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
450: PetscCheckSameComm(eps,1,viewer,3);
451: EPSCheckSolved(eps,1);
452: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
453: if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);
455: PetscCall(PetscViewerGetFormat(viewer,&format));
456: switch (format) {
457: case PETSC_VIEWER_DEFAULT:
458: case PETSC_VIEWER_ASCII_INFO:
459: PetscCall(EPSErrorView_ASCII(eps,etype,viewer));
460: break;
461: case PETSC_VIEWER_ASCII_INFO_DETAIL:
462: PetscCall(EPSErrorView_DETAIL(eps,etype,viewer));
463: break;
464: case PETSC_VIEWER_ASCII_MATLAB:
465: PetscCall(EPSErrorView_MATLAB(eps,etype,viewer));
466: break;
467: default:
468: PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
469: }
470: PetscFunctionReturn(PETSC_SUCCESS);
471: }
473: /*@
474: EPSErrorViewFromOptions - Processes command line options to determine if/how
475: the errors of the computed solution are to be viewed.
477: Collective
479: Input Parameter:
480: . eps - the linear eigensolver context
482: Level: developer
484: .seealso: [](ch:eps), `EPSErrorView()`
485: @*/
486: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
487: {
488: PetscViewer viewer;
489: PetscBool flg;
490: static PetscBool incall = PETSC_FALSE;
491: PetscViewerFormat format;
493: PetscFunctionBegin;
494: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
495: incall = PETSC_TRUE;
496: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg));
497: if (flg) {
498: PetscCall(PetscViewerPushFormat(viewer,format));
499: PetscCall(EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer));
500: PetscCall(PetscViewerPopFormat(viewer));
501: PetscCall(PetscViewerDestroy(&viewer));
502: }
503: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg));
504: if (flg) {
505: PetscCall(PetscViewerPushFormat(viewer,format));
506: PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer));
507: PetscCall(PetscViewerPopFormat(viewer));
508: PetscCall(PetscViewerDestroy(&viewer));
509: }
510: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg));
511: if (flg) {
512: PetscCall(PetscViewerPushFormat(viewer,format));
513: PetscCall(EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer));
514: PetscCall(PetscViewerPopFormat(viewer));
515: PetscCall(PetscViewerDestroy(&viewer));
516: }
517: incall = PETSC_FALSE;
518: PetscFunctionReturn(PETSC_SUCCESS);
519: }
521: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
522: {
523: PetscDraw draw;
524: PetscDrawSP drawsp;
525: PetscReal re,im;
526: PetscScalar kr,ki;
527: PetscInt i,nconv;
529: PetscFunctionBegin;
530: if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
531: PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
532: PetscCall(PetscDrawSetTitle(draw,"Computed Eigenvalues"));
533: PetscCall(PetscDrawSPCreate(draw,1,&drawsp));
534: PetscCall(EPS_GetActualConverged(eps,&nconv));
535: for (i=0;i<nconv;i++) {
536: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
537: #if defined(PETSC_USE_COMPLEX)
538: re = PetscRealPart(kr);
539: im = PetscImaginaryPart(kr);
540: #else
541: re = kr;
542: im = ki;
543: #endif
544: PetscCall(PetscDrawSPAddPoint(drawsp,&re,&im));
545: }
546: PetscCall(PetscDrawSPDraw(drawsp,PETSC_TRUE));
547: PetscCall(PetscDrawSPSave(drawsp));
548: PetscCall(PetscDrawSPDestroy(&drawsp));
549: PetscFunctionReturn(PETSC_SUCCESS);
550: }
552: static PetscErrorCode EPSValuesView_BINARY(EPS eps,PetscViewer viewer)
553: {
554: PetscScalar kr,ki;
555: #if defined(PETSC_HAVE_COMPLEX)
556: PetscInt i,nconv;
557: PetscComplex *ev;
558: #endif
560: PetscFunctionBegin;
561: #if defined(PETSC_HAVE_COMPLEX)
562: PetscCall(EPS_GetActualConverged(eps,&nconv));
563: PetscCall(PetscMalloc1(nconv,&ev));
564: for (i=0;i<nconv;i++) {
565: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
566: #if defined(PETSC_USE_COMPLEX)
567: ev[i] = kr;
568: #else
569: ev[i] = PetscCMPLX(kr,ki);
570: #endif
571: }
572: PetscCall(PetscViewerBinaryWrite(viewer,ev,nconv,PETSC_COMPLEX));
573: PetscCall(PetscFree(ev));
574: #endif
575: PetscFunctionReturn(PETSC_SUCCESS);
576: }
578: #if defined(PETSC_HAVE_HDF5)
579: static PetscErrorCode EPSValuesView_HDF5(EPS eps,PetscViewer viewer)
580: {
581: PetscInt i,n,N,nconv;
582: PetscScalar eig;
583: PetscMPIInt rank;
584: Vec v;
585: char vname[30];
586: const char *ename;
588: PetscFunctionBegin;
589: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank));
590: PetscCall(EPS_GetActualConverged(eps,&nconv));
591: N = nconv;
592: n = rank? 0: N;
593: /* create a vector containing the eigenvalues */
594: PetscCall(VecCreateMPI(PetscObjectComm((PetscObject)eps),n,N,&v));
595: PetscCall(PetscObjectGetName((PetscObject)eps,&ename));
596: PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename));
597: PetscCall(PetscObjectSetName((PetscObject)v,vname));
598: if (!rank) {
599: for (i=0;i<nconv;i++) {
600: PetscCall(EPSGetEigenvalue(eps,i,&eig,NULL));
601: PetscCall(VecSetValue(v,i,eig,INSERT_VALUES));
602: }
603: }
604: PetscCall(VecAssemblyBegin(v));
605: PetscCall(VecAssemblyEnd(v));
606: PetscCall(VecView(v,viewer));
607: #if !defined(PETSC_USE_COMPLEX)
608: /* in real scalars write the imaginary part as a separate vector */
609: PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename));
610: PetscCall(PetscObjectSetName((PetscObject)v,vname));
611: if (!rank) {
612: for (i=0;i<nconv;i++) {
613: PetscCall(EPSGetEigenvalue(eps,i,NULL,&eig));
614: PetscCall(VecSetValue(v,i,eig,INSERT_VALUES));
615: }
616: }
617: PetscCall(VecAssemblyBegin(v));
618: PetscCall(VecAssemblyEnd(v));
619: PetscCall(VecView(v,viewer));
620: #endif
621: PetscCall(VecDestroy(&v));
622: PetscFunctionReturn(PETSC_SUCCESS);
623: }
624: #endif
626: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
627: {
628: PetscInt i,nconv;
629: PetscScalar kr,ki;
631: PetscFunctionBegin;
632: PetscCall(PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n"));
633: PetscCall(EPS_GetActualConverged(eps,&nconv));
634: for (i=0;i<nconv;i++) {
635: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
636: PetscCall(PetscViewerASCIIPrintf(viewer," "));
637: PetscCall(SlepcPrintEigenvalueASCII(viewer,kr,ki));
638: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
639: }
640: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
641: PetscFunctionReturn(PETSC_SUCCESS);
642: }
644: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
645: {
646: PetscInt i,nconv;
647: PetscReal re,im;
648: PetscScalar kr,ki;
649: const char *name;
651: PetscFunctionBegin;
652: PetscCall(PetscObjectGetName((PetscObject)eps,&name));
653: PetscCall(PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name));
654: PetscCall(EPS_GetActualConverged(eps,&nconv));
655: for (i=0;i<nconv;i++) {
656: PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
657: #if defined(PETSC_USE_COMPLEX)
658: re = PetscRealPart(kr);
659: im = PetscImaginaryPart(kr);
660: #else
661: re = kr;
662: im = ki;
663: #endif
664: if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im));
665: else PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re));
666: }
667: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
668: PetscFunctionReturn(PETSC_SUCCESS);
669: }
671: /*@
672: EPSValuesView - Displays the computed eigenvalues in a viewer.
674: Collective
676: Input Parameters:
677: + eps - the linear eigensolver context
678: - viewer - the viewer
680: Options Database Key:
681: . -eps_view_values - print computed eigenvalues
683: Note:
684: The command-line option listed above admits an optional argument
685: specifying the viewer type and options. For instance, use
686: `-eps_view_values :evals.m:ascii_matlab` to save the values in a file
687: that can be executed in Matlab.
688: See `PetscObjectViewFromOptions()` for more details.
690: Level: intermediate
692: .seealso: [](ch:eps), `EPSSolve()`, `EPSVectorsView()`, `EPSErrorView()`
693: @*/
694: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
695: {
696: PetscBool isascii,isdraw,isbinary;
697: PetscViewerFormat format;
698: #if defined(PETSC_HAVE_HDF5)
699: PetscBool ishdf5;
700: #endif
702: PetscFunctionBegin;
704: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
706: PetscCheckSameComm(eps,1,viewer,2);
707: EPSCheckSolved(eps,1);
708: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
709: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
710: #if defined(PETSC_HAVE_HDF5)
711: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
712: #endif
713: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
714: if (isdraw) PetscCall(EPSValuesView_DRAW(eps,viewer));
715: else if (isbinary) PetscCall(EPSValuesView_BINARY(eps,viewer));
716: #if defined(PETSC_HAVE_HDF5)
717: else if (ishdf5) PetscCall(EPSValuesView_HDF5(eps,viewer));
718: #endif
719: else if (isascii) {
720: PetscCall(PetscViewerGetFormat(viewer,&format));
721: switch (format) {
722: case PETSC_VIEWER_DEFAULT:
723: case PETSC_VIEWER_ASCII_INFO:
724: case PETSC_VIEWER_ASCII_INFO_DETAIL:
725: PetscCall(EPSValuesView_ASCII(eps,viewer));
726: break;
727: case PETSC_VIEWER_ASCII_MATLAB:
728: PetscCall(EPSValuesView_MATLAB(eps,viewer));
729: break;
730: default:
731: PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
732: }
733: }
734: PetscFunctionReturn(PETSC_SUCCESS);
735: }
737: /*@
738: EPSValuesViewFromOptions - Processes command line options to determine if/how
739: the computed eigenvalues are to be viewed.
741: Collective
743: Input Parameter:
744: . eps - the linear eigensolver context
746: Level: developer
748: .seealso: [](ch:eps), `EPSValuesView()`
749: @*/
750: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
751: {
752: PetscViewer viewer;
753: PetscBool flg;
754: static PetscBool incall = PETSC_FALSE;
755: PetscViewerFormat format;
757: PetscFunctionBegin;
758: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
759: incall = PETSC_TRUE;
760: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg));
761: if (flg) {
762: PetscCall(PetscViewerPushFormat(viewer,format));
763: PetscCall(EPSValuesView(eps,viewer));
764: PetscCall(PetscViewerPopFormat(viewer));
765: PetscCall(PetscViewerDestroy(&viewer));
766: }
767: incall = PETSC_FALSE;
768: PetscFunctionReturn(PETSC_SUCCESS);
769: }
771: /*@
772: EPSVectorsView - Outputs computed eigenvectors to a viewer.
774: Collective
776: Input Parameters:
777: + eps - the linear eigensolver context
778: - viewer - the viewer
780: Options Database Key:
781: . -eps_view_vectors - output eigenvectors
783: Notes:
784: If PETSc was configured with real scalars, complex conjugate eigenvectors
785: will be viewed as two separate real vectors, one containing the real part
786: and another one containing the imaginary part.
788: If left eigenvectors were computed with a two-sided eigensolver, the right
789: and left eigenvectors are interleaved, that is, the vectors are output in
790: the following order\: `X0, Y0, X1, Y1, X2, Y2, ...`
792: The command-line option listed above admits an optional argument
793: specifying the viewer type and options. For instance, use
794: `-eps_view_vectors binary:evecs.bin` to save the vectors in a binary file.
795: See `PetscObjectViewFromOptions()` for more details.
797: Level: intermediate
799: .seealso: [](ch:eps), `EPSSolve()`, `EPSValuesView()`, `EPSErrorView()`
800: @*/
801: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
802: {
803: PetscInt i,nconv;
804: Vec xr,xi=NULL;
806: PetscFunctionBegin;
808: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
810: PetscCheckSameComm(eps,1,viewer,2);
811: EPSCheckSolved(eps,1);
812: PetscCall(EPS_GetActualConverged(eps,&nconv));
813: if (nconv) {
814: PetscCall(BVCreateVec(eps->V,&xr));
815: #if !defined(PETSC_USE_COMPLEX)
816: PetscCall(BVCreateVec(eps->V,&xi));
817: #endif
818: for (i=0;i<nconv;i++) {
819: PetscCall(EPSGetEigenvector(eps,i,xr,xi));
820: PetscCall(SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)eps));
821: if (eps->twosided || eps->problem_type==EPS_BSE) {
822: PetscCall(EPSGetLeftEigenvector(eps,i,xr,xi));
823: PetscCall(SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)eps));
824: }
825: }
826: PetscCall(VecDestroy(&xr));
827: #if !defined(PETSC_USE_COMPLEX)
828: PetscCall(VecDestroy(&xi));
829: #endif
830: }
831: PetscFunctionReturn(PETSC_SUCCESS);
832: }
834: /*@
835: EPSVectorsViewFromOptions - Processes command line options to determine if/how
836: the computed eigenvectors are to be viewed.
838: Collective
840: Input Parameter:
841: . eps - the linear eigensolver context
843: Level: developer
845: .seealso: [](ch:eps), `EPSVectorsView()`
846: @*/
847: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
848: {
849: PetscViewer viewer;
850: PetscBool flg = PETSC_FALSE;
851: static PetscBool incall = PETSC_FALSE;
852: PetscViewerFormat format;
854: PetscFunctionBegin;
855: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
856: incall = PETSC_TRUE;
857: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg));
858: if (flg) {
859: PetscCall(PetscViewerPushFormat(viewer,format));
860: PetscCall(EPSVectorsView(eps,viewer));
861: PetscCall(PetscViewerPopFormat(viewer));
862: PetscCall(PetscViewerDestroy(&viewer));
863: }
864: incall = PETSC_FALSE;
865: PetscFunctionReturn(PETSC_SUCCESS);
866: }