Actual source code: nepview.c
slepc-3.22.1 2024-10-28
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
10: /*
11: NEP routines related to various viewers
12: */
14: #include <slepc/private/nepimpl.h>
15: #include <slepc/private/bvimpl.h>
16: #include <petscdraw.h>
18: /*@
19: NEPView - Prints the NEP data structure.
21: Collective
23: Input Parameters:
24: + nep - the nonlinear eigenproblem solver context
25: - viewer - optional visualization context
27: Options Database Key:
28: . -nep_view - Calls NEPView() at end of NEPSolve()
30: Note:
31: The available visualization contexts include
32: + PETSC_VIEWER_STDOUT_SELF - standard output (default)
33: - PETSC_VIEWER_STDOUT_WORLD - synchronized standard
34: output where only the first processor opens
35: the file. All other processors send their
36: data to the first processor to print.
38: The user can open an alternative visualization context with
39: PetscViewerASCIIOpen() - output to a specified file.
41: Level: beginner
43: .seealso: FNView()
44: @*/
45: PetscErrorCode NEPView(NEP nep,PetscViewer viewer)
46: {
47: const char *type=NULL;
48: char str[50];
49: PetscInt i;
50: PetscBool isascii,istrivial;
51: PetscViewer sviewer;
52: MPI_Comm child;
54: PetscFunctionBegin;
56: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
58: PetscCheckSameComm(nep,1,viewer,2);
60: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
61: if (isascii) {
62: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)nep,viewer));
63: PetscCall(PetscViewerASCIIPushTab(viewer));
64: PetscTryTypeMethod(nep,view,viewer);
65: PetscCall(PetscViewerASCIIPopTab(viewer));
66: if (nep->problem_type) {
67: switch (nep->problem_type) {
68: case NEP_GENERAL: type = "general nonlinear eigenvalue problem"; break;
69: case NEP_RATIONAL: type = "rational eigenvalue problem"; break;
70: }
71: } else type = "not yet set";
72: PetscCall(PetscViewerASCIIPrintf(viewer," problem type: %s\n",type));
73: if (nep->fui) {
74: switch (nep->fui) {
75: case NEP_USER_INTERFACE_CALLBACK:
76: PetscCall(PetscViewerASCIIPrintf(viewer," nonlinear operator from user callbacks\n"));
77: break;
78: case NEP_USER_INTERFACE_SPLIT:
79: PetscCall(PetscViewerASCIIPrintf(viewer," nonlinear operator in split form\n"));
80: PetscCall(PetscViewerASCIIPrintf(viewer," number of terms: %" PetscInt_FMT "\n",nep->nt));
81: PetscCall(PetscViewerASCIIPrintf(viewer," nonzero pattern of the matrices: %s\n",MatStructures[nep->mstr]));
82: break;
83: }
84: } else PetscCall(PetscViewerASCIIPrintf(viewer," nonlinear operator not specified yet\n"));
85: PetscCall(PetscViewerASCIIPrintf(viewer," selected portion of the spectrum: "));
86: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
87: PetscCall(SlepcSNPrintfScalar(str,sizeof(str),nep->target,PETSC_FALSE));
88: if (!nep->which) PetscCall(PetscViewerASCIIPrintf(viewer,"not yet set\n"));
89: else switch (nep->which) {
90: case NEP_WHICH_USER:
91: PetscCall(PetscViewerASCIIPrintf(viewer,"user defined\n"));
92: break;
93: case NEP_TARGET_MAGNITUDE:
94: PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str));
95: break;
96: case NEP_TARGET_REAL:
97: PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str));
98: break;
99: case NEP_TARGET_IMAGINARY:
100: PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str));
101: break;
102: case NEP_LARGEST_MAGNITUDE:
103: PetscCall(PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n"));
104: break;
105: case NEP_SMALLEST_MAGNITUDE:
106: PetscCall(PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n"));
107: break;
108: case NEP_LARGEST_REAL:
109: PetscCall(PetscViewerASCIIPrintf(viewer,"largest real parts\n"));
110: break;
111: case NEP_SMALLEST_REAL:
112: PetscCall(PetscViewerASCIIPrintf(viewer,"smallest real parts\n"));
113: break;
114: case NEP_LARGEST_IMAGINARY:
115: PetscCall(PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n"));
116: break;
117: case NEP_SMALLEST_IMAGINARY:
118: PetscCall(PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n"));
119: break;
120: case NEP_ALL:
121: PetscCall(PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n"));
122: break;
123: }
124: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
125: if (nep->twosided) PetscCall(PetscViewerASCIIPrintf(viewer," using two-sided variant (for left eigenvectors)\n"));
126: PetscCall(PetscViewerASCIIPrintf(viewer," number of eigenvalues (nev): %" PetscInt_FMT "\n",nep->nev));
127: PetscCall(PetscViewerASCIIPrintf(viewer," number of column vectors (ncv): %" PetscInt_FMT "\n",nep->ncv));
128: PetscCall(PetscViewerASCIIPrintf(viewer," maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",nep->mpd));
129: PetscCall(PetscViewerASCIIPrintf(viewer," maximum number of iterations: %" PetscInt_FMT "\n",nep->max_it));
130: PetscCall(PetscViewerASCIIPrintf(viewer," tolerance: %g\n",(double)nep->tol));
131: PetscCall(PetscViewerASCIIPrintf(viewer," convergence test: "));
132: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
133: switch (nep->conv) {
134: case NEP_CONV_ABS:
135: PetscCall(PetscViewerASCIIPrintf(viewer,"absolute\n"));break;
136: case NEP_CONV_REL:
137: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n"));break;
138: case NEP_CONV_NORM:
139: PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the matrix norms\n"));
140: if (nep->nrma) {
141: PetscCall(PetscViewerASCIIPrintf(viewer," computed matrix norms: %g",(double)nep->nrma[0]));
142: for (i=1;i<nep->nt;i++) PetscCall(PetscViewerASCIIPrintf(viewer,", %g",(double)nep->nrma[i]));
143: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
144: }
145: break;
146: case NEP_CONV_USER:
147: PetscCall(PetscViewerASCIIPrintf(viewer,"user-defined\n"));break;
148: }
149: PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
150: if (nep->refine) {
151: PetscCall(PetscViewerASCIIPrintf(viewer," iterative refinement: %s, with %s scheme\n",NEPRefineTypes[nep->refine],NEPRefineSchemes[nep->scheme]));
152: PetscCall(PetscViewerASCIIPrintf(viewer," refinement stopping criterion: tol=%g, its=%" PetscInt_FMT "\n",(double)nep->rtol,nep->rits));
153: if (nep->npart>1) PetscCall(PetscViewerASCIIPrintf(viewer," splitting communicator in %" PetscInt_FMT " partitions for refinement\n",nep->npart));
154: }
155: if (nep->nini) PetscCall(PetscViewerASCIIPrintf(viewer," dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(nep->nini)));
156: } else PetscTryTypeMethod(nep,view,viewer);
157: PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
158: if (!nep->V) PetscCall(NEPGetBV(nep,&nep->V));
159: PetscCall(BVView(nep->V,viewer));
160: if (!nep->rg) PetscCall(NEPGetRG(nep,&nep->rg));
161: PetscCall(RGIsTrivial(nep->rg,&istrivial));
162: if (!istrivial) PetscCall(RGView(nep->rg,viewer));
163: if (nep->useds) {
164: if (!nep->ds) PetscCall(NEPGetDS(nep,&nep->ds));
165: PetscCall(DSView(nep->ds,viewer));
166: }
167: PetscCall(PetscViewerPopFormat(viewer));
168: if (nep->refine!=NEP_REFINE_NONE) {
169: PetscCall(PetscViewerASCIIPushTab(viewer));
170: if (nep->npart>1) {
171: if (nep->refinesubc->color==0) {
172: PetscCall(PetscSubcommGetChild(nep->refinesubc,&child));
173: PetscCall(PetscViewerASCIIGetStdout(child,&sviewer));
174: PetscCall(KSPView(nep->refineksp,sviewer));
175: }
176: } else PetscCall(KSPView(nep->refineksp,viewer));
177: PetscCall(PetscViewerASCIIPopTab(viewer));
178: }
179: PetscFunctionReturn(PETSC_SUCCESS);
180: }
182: /*@
183: NEPViewFromOptions - View from options
185: Collective
187: Input Parameters:
188: + nep - the nonlinear eigensolver context
189: . obj - optional object
190: - name - command line option
192: Level: intermediate
194: .seealso: NEPView(), NEPCreate()
195: @*/
196: PetscErrorCode NEPViewFromOptions(NEP nep,PetscObject obj,const char name[])
197: {
198: PetscFunctionBegin;
200: PetscCall(PetscObjectViewFromOptions((PetscObject)nep,obj,name));
201: PetscFunctionReturn(PETSC_SUCCESS);
202: }
204: /*@
205: NEPConvergedReasonView - Displays the reason a NEP solve converged or diverged.
207: Collective
209: Input Parameters:
210: + nep - the nonlinear eigensolver context
211: - viewer - the viewer to display the reason
213: Options Database Keys:
214: . -nep_converged_reason - print reason for convergence, and number of iterations
216: Note:
217: To change the format of the output call PetscViewerPushFormat(viewer,format) before
218: this call. Use PETSC_VIEWER_DEFAULT for the default, use PETSC_VIEWER_FAILED to only
219: display a reason if it fails. The latter can be set in the command line with
220: -nep_converged_reason ::failed
222: Level: intermediate
224: .seealso: NEPSetConvergenceTest(), NEPSetTolerances(), NEPGetIterationNumber(), NEPConvergedReasonViewFromOptions()
225: @*/
226: PetscErrorCode NEPConvergedReasonView(NEP nep,PetscViewer viewer)
227: {
228: PetscBool isAscii;
229: PetscViewerFormat format;
231: PetscFunctionBegin;
232: if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)nep));
233: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
234: if (isAscii) {
235: PetscCall(PetscViewerGetFormat(viewer,&format));
236: PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)nep)->tablevel));
237: if (nep->reason > 0 && format != PETSC_VIEWER_FAILED) PetscCall(PetscViewerASCIIPrintf(viewer,"%s Nonlinear eigensolve converged (%" PetscInt_FMT " eigenpair%s) due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)nep)->prefix?((PetscObject)nep)->prefix:"",nep->nconv,(nep->nconv>1)?"s":"",NEPConvergedReasons[nep->reason],nep->its));
238: else if (nep->reason <= 0) PetscCall(PetscViewerASCIIPrintf(viewer,"%s Nonlinear eigensolve did not converge due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)nep)->prefix?((PetscObject)nep)->prefix:"",NEPConvergedReasons[nep->reason],nep->its));
239: PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)nep)->tablevel));
240: }
241: PetscFunctionReturn(PETSC_SUCCESS);
242: }
244: /*@
245: NEPConvergedReasonViewFromOptions - Processes command line options to determine if/how
246: the NEP converged reason is to be viewed.
248: Collective
250: Input Parameter:
251: . nep - the nonlinear eigensolver context
253: Level: developer
255: .seealso: NEPConvergedReasonView()
256: @*/
257: PetscErrorCode NEPConvergedReasonViewFromOptions(NEP nep)
258: {
259: PetscViewer viewer;
260: PetscBool flg;
261: static PetscBool incall = PETSC_FALSE;
262: PetscViewerFormat format;
264: PetscFunctionBegin;
265: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
266: incall = PETSC_TRUE;
267: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_converged_reason",&viewer,&format,&flg));
268: if (flg) {
269: PetscCall(PetscViewerPushFormat(viewer,format));
270: PetscCall(NEPConvergedReasonView(nep,viewer));
271: PetscCall(PetscViewerPopFormat(viewer));
272: PetscCall(PetscViewerDestroy(&viewer));
273: }
274: incall = PETSC_FALSE;
275: PetscFunctionReturn(PETSC_SUCCESS);
276: }
278: static PetscErrorCode NEPErrorView_ASCII(NEP nep,NEPErrorType etype,PetscViewer viewer)
279: {
280: PetscReal error;
281: PetscInt i,j,k,nvals;
283: PetscFunctionBegin;
284: nvals = (nep->which==NEP_ALL)? nep->nconv: nep->nev;
285: if (nep->which!=NEP_ALL && nep->nconv<nvals) {
286: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " eigenvalues converged\n\n",nep->nev));
287: PetscFunctionReturn(PETSC_SUCCESS);
288: }
289: if (nep->which==NEP_ALL && !nvals) {
290: PetscCall(PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n"));
291: PetscFunctionReturn(PETSC_SUCCESS);
292: }
293: for (i=0;i<nvals;i++) {
294: PetscCall(NEPComputeError(nep,i,etype,&error));
295: if (error>=5.0*nep->tol) {
296: PetscCall(PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",nvals));
297: PetscFunctionReturn(PETSC_SUCCESS);
298: }
299: }
300: if (nep->which==NEP_ALL) PetscCall(PetscViewerASCIIPrintf(viewer," Found %" PetscInt_FMT " eigenvalues, all of them computed up to the required tolerance:",nvals));
301: else PetscCall(PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:"));
302: for (i=0;i<=(nvals-1)/8;i++) {
303: PetscCall(PetscViewerASCIIPrintf(viewer,"\n "));
304: for (j=0;j<PetscMin(8,nvals-8*i);j++) {
305: k = nep->perm[8*i+j];
306: PetscCall(SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]));
307: if (8*i+j+1<nvals) PetscCall(PetscViewerASCIIPrintf(viewer,", "));
308: }
309: }
310: PetscCall(PetscViewerASCIIPrintf(viewer,"\n\n"));
311: PetscFunctionReturn(PETSC_SUCCESS);
312: }
314: static PetscErrorCode NEPErrorView_DETAIL(NEP nep,NEPErrorType etype,PetscViewer viewer)
315: {
316: PetscReal error,re,im;
317: PetscScalar kr,ki;
318: PetscInt i;
319: char ex[30],sep[]=" ---------------------- --------------------\n";
321: PetscFunctionBegin;
322: if (!nep->nconv) PetscFunctionReturn(PETSC_SUCCESS);
323: switch (etype) {
324: case NEP_ERROR_ABSOLUTE:
325: PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||T(k)x||"));
326: break;
327: case NEP_ERROR_RELATIVE:
328: PetscCall(PetscSNPrintf(ex,sizeof(ex)," ||T(k)x||/||kx||"));
329: break;
330: case NEP_ERROR_BACKWARD:
331: PetscCall(PetscSNPrintf(ex,sizeof(ex)," eta(x,k)"));
332: break;
333: }
334: PetscCall(PetscViewerASCIIPrintf(viewer,"%s k %s\n%s",sep,ex,sep));
335: for (i=0;i<nep->nconv;i++) {
336: PetscCall(NEPGetEigenpair(nep,i,&kr,&ki,NULL,NULL));
337: PetscCall(NEPComputeError(nep,i,etype,&error));
338: #if defined(PETSC_USE_COMPLEX)
339: re = PetscRealPart(kr);
340: im = PetscImaginaryPart(kr);
341: #else
342: re = kr;
343: im = ki;
344: #endif
345: if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer," % 9f%+9fi %12g\n",(double)re,(double)im,(double)error));
346: else PetscCall(PetscViewerASCIIPrintf(viewer," % 12f %12g\n",(double)re,(double)error));
347: }
348: PetscCall(PetscViewerASCIIPrintf(viewer,"%s",sep));
349: PetscFunctionReturn(PETSC_SUCCESS);
350: }
352: static PetscErrorCode NEPErrorView_MATLAB(NEP nep,NEPErrorType etype,PetscViewer viewer)
353: {
354: PetscReal error;
355: PetscInt i;
356: const char *name;
358: PetscFunctionBegin;
359: PetscCall(PetscObjectGetName((PetscObject)nep,&name));
360: PetscCall(PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name));
361: for (i=0;i<nep->nconv;i++) {
362: PetscCall(NEPComputeError(nep,i,etype,&error));
363: PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error));
364: }
365: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
366: PetscFunctionReturn(PETSC_SUCCESS);
367: }
369: /*@
370: NEPErrorView - Displays the errors associated with the computed solution
371: (as well as the eigenvalues).
373: Collective
375: Input Parameters:
376: + nep - the nonlinear eigensolver context
377: . etype - error type
378: - viewer - optional visualization context
380: Options Database Keys:
381: + -nep_error_absolute - print absolute errors of each eigenpair
382: . -nep_error_relative - print relative errors of each eigenpair
383: - -nep_error_backward - print backward errors of each eigenpair
385: Notes:
386: By default, this function checks the error of all eigenpairs and prints
387: the eigenvalues if all of them are below the requested tolerance.
388: If the viewer has format=PETSC_VIEWER_ASCII_INFO_DETAIL then a table with
389: eigenvalues and corresponding errors is printed.
391: Level: intermediate
393: .seealso: NEPSolve(), NEPValuesView(), NEPVectorsView()
394: @*/
395: PetscErrorCode NEPErrorView(NEP nep,NEPErrorType etype,PetscViewer viewer)
396: {
397: PetscBool isascii;
398: PetscViewerFormat format;
400: PetscFunctionBegin;
402: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
404: PetscCheckSameComm(nep,1,viewer,3);
405: NEPCheckSolved(nep,1);
406: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
407: if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);
409: PetscCall(PetscViewerGetFormat(viewer,&format));
410: switch (format) {
411: case PETSC_VIEWER_DEFAULT:
412: case PETSC_VIEWER_ASCII_INFO:
413: PetscCall(NEPErrorView_ASCII(nep,etype,viewer));
414: break;
415: case PETSC_VIEWER_ASCII_INFO_DETAIL:
416: PetscCall(NEPErrorView_DETAIL(nep,etype,viewer));
417: break;
418: case PETSC_VIEWER_ASCII_MATLAB:
419: PetscCall(NEPErrorView_MATLAB(nep,etype,viewer));
420: break;
421: default:
422: PetscCall(PetscInfo(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
423: }
424: PetscFunctionReturn(PETSC_SUCCESS);
425: }
427: /*@
428: NEPErrorViewFromOptions - Processes command line options to determine if/how
429: the errors of the computed solution are to be viewed.
431: Collective
433: Input Parameter:
434: . nep - the nonlinear eigensolver context
436: Level: developer
438: .seealso: NEPErrorView()
439: @*/
440: PetscErrorCode NEPErrorViewFromOptions(NEP nep)
441: {
442: PetscViewer viewer;
443: PetscBool flg;
444: static PetscBool incall = PETSC_FALSE;
445: PetscViewerFormat format;
447: PetscFunctionBegin;
448: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
449: incall = PETSC_TRUE;
450: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_absolute",&viewer,&format,&flg));
451: if (flg) {
452: PetscCall(PetscViewerPushFormat(viewer,format));
453: PetscCall(NEPErrorView(nep,NEP_ERROR_ABSOLUTE,viewer));
454: PetscCall(PetscViewerPopFormat(viewer));
455: PetscCall(PetscViewerDestroy(&viewer));
456: }
457: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_relative",&viewer,&format,&flg));
458: if (flg) {
459: PetscCall(PetscViewerPushFormat(viewer,format));
460: PetscCall(NEPErrorView(nep,NEP_ERROR_RELATIVE,viewer));
461: PetscCall(PetscViewerPopFormat(viewer));
462: PetscCall(PetscViewerDestroy(&viewer));
463: }
464: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_error_backward",&viewer,&format,&flg));
465: if (flg) {
466: PetscCall(PetscViewerPushFormat(viewer,format));
467: PetscCall(NEPErrorView(nep,NEP_ERROR_BACKWARD,viewer));
468: PetscCall(PetscViewerPopFormat(viewer));
469: PetscCall(PetscViewerDestroy(&viewer));
470: }
471: incall = PETSC_FALSE;
472: PetscFunctionReturn(PETSC_SUCCESS);
473: }
475: static PetscErrorCode NEPValuesView_DRAW(NEP nep,PetscViewer viewer)
476: {
477: PetscDraw draw;
478: PetscDrawSP drawsp;
479: PetscReal re,im;
480: PetscInt i,k;
482: PetscFunctionBegin;
483: if (!nep->nconv) PetscFunctionReturn(PETSC_SUCCESS);
484: PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
485: PetscCall(PetscDrawSetTitle(draw,"Computed Eigenvalues"));
486: PetscCall(PetscDrawSPCreate(draw,1,&drawsp));
487: for (i=0;i<nep->nconv;i++) {
488: k = nep->perm[i];
489: #if defined(PETSC_USE_COMPLEX)
490: re = PetscRealPart(nep->eigr[k]);
491: im = PetscImaginaryPart(nep->eigr[k]);
492: #else
493: re = nep->eigr[k];
494: im = nep->eigi[k];
495: #endif
496: PetscCall(PetscDrawSPAddPoint(drawsp,&re,&im));
497: }
498: PetscCall(PetscDrawSPDraw(drawsp,PETSC_TRUE));
499: PetscCall(PetscDrawSPSave(drawsp));
500: PetscCall(PetscDrawSPDestroy(&drawsp));
501: PetscFunctionReturn(PETSC_SUCCESS);
502: }
504: static PetscErrorCode NEPValuesView_BINARY(NEP nep,PetscViewer viewer)
505: {
506: #if defined(PETSC_HAVE_COMPLEX)
507: PetscInt i,k;
508: PetscComplex *ev;
509: #endif
511: PetscFunctionBegin;
512: #if defined(PETSC_HAVE_COMPLEX)
513: PetscCall(PetscMalloc1(nep->nconv,&ev));
514: for (i=0;i<nep->nconv;i++) {
515: k = nep->perm[i];
516: #if defined(PETSC_USE_COMPLEX)
517: ev[i] = nep->eigr[k];
518: #else
519: ev[i] = PetscCMPLX(nep->eigr[k],nep->eigi[k]);
520: #endif
521: }
522: PetscCall(PetscViewerBinaryWrite(viewer,ev,nep->nconv,PETSC_COMPLEX));
523: PetscCall(PetscFree(ev));
524: #endif
525: PetscFunctionReturn(PETSC_SUCCESS);
526: }
528: #if defined(PETSC_HAVE_HDF5)
529: static PetscErrorCode NEPValuesView_HDF5(NEP nep,PetscViewer viewer)
530: {
531: PetscInt i,k,n,N;
532: PetscMPIInt rank;
533: Vec v;
534: char vname[30];
535: const char *ename;
537: PetscFunctionBegin;
538: PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)nep),&rank));
539: N = nep->nconv;
540: n = rank? 0: N;
541: /* create a vector containing the eigenvalues */
542: PetscCall(VecCreateMPI(PetscObjectComm((PetscObject)nep),n,N,&v));
543: PetscCall(PetscObjectGetName((PetscObject)nep,&ename));
544: PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename));
545: PetscCall(PetscObjectSetName((PetscObject)v,vname));
546: if (!rank) {
547: for (i=0;i<nep->nconv;i++) {
548: k = nep->perm[i];
549: PetscCall(VecSetValue(v,i,nep->eigr[k],INSERT_VALUES));
550: }
551: }
552: PetscCall(VecAssemblyBegin(v));
553: PetscCall(VecAssemblyEnd(v));
554: PetscCall(VecView(v,viewer));
555: #if !defined(PETSC_USE_COMPLEX)
556: /* in real scalars write the imaginary part as a separate vector */
557: PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename));
558: PetscCall(PetscObjectSetName((PetscObject)v,vname));
559: if (!rank) {
560: for (i=0;i<nep->nconv;i++) {
561: k = nep->perm[i];
562: PetscCall(VecSetValue(v,i,nep->eigi[k],INSERT_VALUES));
563: }
564: }
565: PetscCall(VecAssemblyBegin(v));
566: PetscCall(VecAssemblyEnd(v));
567: PetscCall(VecView(v,viewer));
568: #endif
569: PetscCall(VecDestroy(&v));
570: PetscFunctionReturn(PETSC_SUCCESS);
571: }
572: #endif
574: static PetscErrorCode NEPValuesView_ASCII(NEP nep,PetscViewer viewer)
575: {
576: PetscInt i,k;
578: PetscFunctionBegin;
579: PetscCall(PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n"));
580: for (i=0;i<nep->nconv;i++) {
581: k = nep->perm[i];
582: PetscCall(PetscViewerASCIIPrintf(viewer," "));
583: PetscCall(SlepcPrintEigenvalueASCII(viewer,nep->eigr[k],nep->eigi[k]));
584: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
585: }
586: PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
587: PetscFunctionReturn(PETSC_SUCCESS);
588: }
590: static PetscErrorCode NEPValuesView_MATLAB(NEP nep,PetscViewer viewer)
591: {
592: PetscInt i,k;
593: PetscReal re,im;
594: const char *name;
596: PetscFunctionBegin;
597: PetscCall(PetscObjectGetName((PetscObject)nep,&name));
598: PetscCall(PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name));
599: for (i=0;i<nep->nconv;i++) {
600: k = nep->perm[i];
601: #if defined(PETSC_USE_COMPLEX)
602: re = PetscRealPart(nep->eigr[k]);
603: im = PetscImaginaryPart(nep->eigr[k]);
604: #else
605: re = nep->eigr[k];
606: im = nep->eigi[k];
607: #endif
608: if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im));
609: else PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re));
610: }
611: PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
612: PetscFunctionReturn(PETSC_SUCCESS);
613: }
615: /*@
616: NEPValuesView - Displays the computed eigenvalues in a viewer.
618: Collective
620: Input Parameters:
621: + nep - the nonlinear eigensolver context
622: - viewer - the viewer
624: Options Database Key:
625: . -nep_view_values - print computed eigenvalues
627: Level: intermediate
629: .seealso: NEPSolve(), NEPVectorsView(), NEPErrorView()
630: @*/
631: PetscErrorCode NEPValuesView(NEP nep,PetscViewer viewer)
632: {
633: PetscBool isascii,isdraw,isbinary;
634: PetscViewerFormat format;
635: #if defined(PETSC_HAVE_HDF5)
636: PetscBool ishdf5;
637: #endif
639: PetscFunctionBegin;
641: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
643: PetscCheckSameComm(nep,1,viewer,2);
644: NEPCheckSolved(nep,1);
645: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
646: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
647: #if defined(PETSC_HAVE_HDF5)
648: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
649: #endif
650: PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
651: if (isdraw) PetscCall(NEPValuesView_DRAW(nep,viewer));
652: else if (isbinary) PetscCall(NEPValuesView_BINARY(nep,viewer));
653: #if defined(PETSC_HAVE_HDF5)
654: else if (ishdf5) PetscCall(NEPValuesView_HDF5(nep,viewer));
655: #endif
656: else if (isascii) {
657: PetscCall(PetscViewerGetFormat(viewer,&format));
658: switch (format) {
659: case PETSC_VIEWER_DEFAULT:
660: case PETSC_VIEWER_ASCII_INFO:
661: case PETSC_VIEWER_ASCII_INFO_DETAIL:
662: PetscCall(NEPValuesView_ASCII(nep,viewer));
663: break;
664: case PETSC_VIEWER_ASCII_MATLAB:
665: PetscCall(NEPValuesView_MATLAB(nep,viewer));
666: break;
667: default:
668: PetscCall(PetscInfo(nep,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
669: }
670: }
671: PetscFunctionReturn(PETSC_SUCCESS);
672: }
674: /*@
675: NEPValuesViewFromOptions - Processes command line options to determine if/how
676: the computed eigenvalues are to be viewed.
678: Collective
680: Input Parameter:
681: . nep - the nonlinear eigensolver context
683: Level: developer
685: .seealso: NEPValuesView()
686: @*/
687: PetscErrorCode NEPValuesViewFromOptions(NEP nep)
688: {
689: PetscViewer viewer;
690: PetscBool flg;
691: static PetscBool incall = PETSC_FALSE;
692: PetscViewerFormat format;
694: PetscFunctionBegin;
695: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
696: incall = PETSC_TRUE;
697: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_values",&viewer,&format,&flg));
698: if (flg) {
699: PetscCall(PetscViewerPushFormat(viewer,format));
700: PetscCall(NEPValuesView(nep,viewer));
701: PetscCall(PetscViewerPopFormat(viewer));
702: PetscCall(PetscViewerDestroy(&viewer));
703: }
704: incall = PETSC_FALSE;
705: PetscFunctionReturn(PETSC_SUCCESS);
706: }
708: /*@
709: NEPVectorsView - Outputs computed eigenvectors to a viewer.
711: Collective
713: Input Parameters:
714: + nep - the nonlinear eigensolver context
715: - viewer - the viewer
717: Options Database Key:
718: . -nep_view_vectors - output eigenvectors.
720: Notes:
721: If PETSc was configured with real scalars, complex conjugate eigenvectors
722: will be viewed as two separate real vectors, one containing the real part
723: and another one containing the imaginary part.
725: If left eigenvectors were computed with a two-sided eigensolver, the right
726: and left eigenvectors are interleaved, that is, the vectors are output in
727: the following order X0, Y0, X1, Y1, X2, Y2, ...
729: Level: intermediate
731: .seealso: NEPSolve(), NEPValuesView(), NEPErrorView()
732: @*/
733: PetscErrorCode NEPVectorsView(NEP nep,PetscViewer viewer)
734: {
735: PetscInt i,k;
736: Vec xr,xi=NULL;
738: PetscFunctionBegin;
740: if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)nep),&viewer));
742: PetscCheckSameComm(nep,1,viewer,2);
743: NEPCheckSolved(nep,1);
744: if (nep->nconv) {
745: PetscCall(NEPComputeVectors(nep));
746: PetscCall(BVCreateVec(nep->V,&xr));
747: #if !defined(PETSC_USE_COMPLEX)
748: PetscCall(BVCreateVec(nep->V,&xi));
749: #endif
750: for (i=0;i<nep->nconv;i++) {
751: k = nep->perm[i];
752: PetscCall(BV_GetEigenvector(nep->V,k,nep->eigi[k],xr,xi));
753: PetscCall(SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)nep));
754: if (nep->twosided) {
755: PetscCall(BV_GetEigenvector(nep->W,k,nep->eigi[k],xr,xi));
756: PetscCall(SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)nep));
757: }
758: }
759: PetscCall(VecDestroy(&xr));
760: #if !defined(PETSC_USE_COMPLEX)
761: PetscCall(VecDestroy(&xi));
762: #endif
763: }
764: PetscFunctionReturn(PETSC_SUCCESS);
765: }
767: /*@
768: NEPVectorsViewFromOptions - Processes command line options to determine if/how
769: the computed eigenvectors are to be viewed.
771: Collective
773: Input Parameter:
774: . nep - the nonlinear eigensolver context
776: Level: developer
778: .seealso: NEPVectorsView()
779: @*/
780: PetscErrorCode NEPVectorsViewFromOptions(NEP nep)
781: {
782: PetscViewer viewer;
783: PetscBool flg = PETSC_FALSE;
784: static PetscBool incall = PETSC_FALSE;
785: PetscViewerFormat format;
787: PetscFunctionBegin;
788: if (incall) PetscFunctionReturn(PETSC_SUCCESS);
789: incall = PETSC_TRUE;
790: PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)nep),((PetscObject)nep)->options,((PetscObject)nep)->prefix,"-nep_view_vectors",&viewer,&format,&flg));
791: if (flg) {
792: PetscCall(PetscViewerPushFormat(viewer,format));
793: PetscCall(NEPVectorsView(nep,viewer));
794: PetscCall(PetscViewerPopFormat(viewer));
795: PetscCall(PetscViewerDestroy(&viewer));
796: }
797: incall = PETSC_FALSE;
798: PetscFunctionReturn(PETSC_SUCCESS);
799: }