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