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