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