Actual source code: epsview.c

  1: /*
  2:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  3:    SLEPc - Scalable Library for Eigenvalue Problem Computations
  4:    Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain

  6:    This file is part of SLEPc.
  7:    SLEPc is distributed under a 2-clause BSD license (see LICENSE).
  8:    - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
  9: */
 10: /*
 11:    EPS routines related to various viewers
 12: */

 14: #include <slepc/private/epsimpl.h>
 15: #include <slepc/private/bvimpl.h>
 16: #include <petscdraw.h>

 18: /*@
 19:    EPSView - Prints the `EPS` data structure.

 21:    Collective

 23:    Input Parameters:
 24: +  eps - the linear eigensolver context
 25: -  viewer - optional visualization context

 27:    Options Database Key:
 28: .  -eps_view - calls `EPSView()` at end of `EPSSolve()`

 30:    Notes:
 31:    The available visualization contexts include
 32: +     `PETSC_VIEWER_STDOUT_SELF` - standard output (default)
 33: -     `PETSC_VIEWER_STDOUT_WORLD` - synchronized standard output where only the
 34:          first process opens the file; all other processes send their data to the
 35:          first one to print

 37:    The user can open an alternative visualization context with `PetscViewerASCIIOpen()`
 38:    to output to a specified file.

 40:    Level: beginner

 42: .seealso: [](ch:eps), `EPSCreate()`, `EPSViewFromOptions()`, `STView()`
 43: @*/
 44: PetscErrorCode EPSView(EPS eps,PetscViewer viewer)
 45: {
 46:   const char     *type=NULL,*extr=NULL,*bal=NULL;
 47:   char           str[50];
 48:   PetscBool      isascii,isexternal,istrivial,isstruct=PETSC_FALSE,reduced=PETSC_FALSE,flg;
 49:   Mat            A;

 51:   PetscFunctionBegin;
 53:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
 55:   PetscCheckSameComm(eps,1,viewer,2);

 57:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
 58:   if (isascii) {
 59:     PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)eps,viewer));
 60:     PetscCall(PetscViewerASCIIPushTab(viewer));
 61:     PetscTryTypeMethod(eps,view,viewer);
 62:     PetscCall(PetscViewerASCIIPopTab(viewer));
 63:     if (eps->problem_type) {
 64:       switch (eps->problem_type) {
 65:         case EPS_HEP:    type = SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
 66:         case EPS_GHEP:   type = "generalized " SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
 67:         case EPS_NHEP:   type = "non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
 68:         case EPS_GNHEP:  type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem"; break;
 69:         case EPS_PGNHEP: type = "generalized non-" SLEPC_STRING_HERMITIAN " eigenvalue problem with " SLEPC_STRING_HERMITIAN " positive definite B"; break;
 70:         case EPS_GHIEP:  type = "generalized " SLEPC_STRING_HERMITIAN "-indefinite eigenvalue problem"; break;
 71:         case EPS_BSE:    type = "structured Bethe-Salpeter eigenvalue problem"; break;
 72:         case EPS_HAMILT: type = "structured Hamiltonian eigenvalue problem"; break;
 73:         case EPS_LREP:   type = "structured Linear Response eigenvalue problem"; break;
 74:       }
 75:     } else type = "not yet set";
 76:     PetscCall(PetscViewerASCIIPrintf(viewer,"  problem type: %s\n",type));
 77:     PetscCall(EPSGetOperators(eps,&A,NULL));
 78:     if (A) PetscCall(SlepcCheckMatStruct(A,0,&isstruct));
 79:     if (isstruct) {
 80:       PetscCall(SlepcCheckMatStruct(A,SLEPC_MAT_STRUCT_BSE,&flg));
 81:       if (flg) PetscCall(PetscViewerASCIIPrintf(viewer,"  matrix A has a Bethe-Salpeter structure\n"));
 82:       PetscCall(SlepcCheckMatStruct(A,SLEPC_MAT_STRUCT_HAMILT,&flg));
 83:       if (flg) PetscCall(PetscViewerASCIIPrintf(viewer,"  matrix A has a Hamiltonian structure\n"));
 84:       PetscCall(SlepcCheckMatStruct(A,SLEPC_MAT_STRUCT_LREP,&flg));
 85:       if (flg) {
 86:         PetscCall(SlepcCheckMatLREPReduced(A,&reduced));
 87:         PetscCall(PetscViewerASCIIPrintf(viewer,"  matrix A has a%s Linear Response structure\n",reduced?" reduced":""));
 88:       }
 89:     }
 90:     if (eps->extraction) {
 91:       switch (eps->extraction) {
 92:         case EPS_RITZ:              extr = "Rayleigh-Ritz"; break;
 93:         case EPS_HARMONIC:          extr = "harmonic Ritz"; break;
 94:         case EPS_HARMONIC_RELATIVE: extr = "relative harmonic Ritz"; break;
 95:         case EPS_HARMONIC_RIGHT:    extr = "right harmonic Ritz"; break;
 96:         case EPS_HARMONIC_LARGEST:  extr = "largest harmonic Ritz"; break;
 97:         case EPS_REFINED:           extr = "refined Ritz"; break;
 98:         case EPS_REFINED_HARMONIC:  extr = "refined harmonic Ritz"; break;
 99:       }
100:       PetscCall(PetscViewerASCIIPrintf(viewer,"  extraction type: %s\n",extr));
101:     }
102:     if (!eps->ishermitian && eps->balance!=EPS_BALANCE_NONE) {
103:       switch (eps->balance) {
104:         case EPS_BALANCE_NONE:    break;
105:         case EPS_BALANCE_ONESIDE: bal = "one-sided Krylov"; break;
106:         case EPS_BALANCE_TWOSIDE: bal = "two-sided Krylov"; break;
107:         case EPS_BALANCE_USER:    bal = "user-defined matrix"; break;
108:       }
109:       PetscCall(PetscViewerASCIIPrintf(viewer,"  balancing enabled: %s",bal));
110:       PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
111:       if (eps->balance==EPS_BALANCE_ONESIDE || eps->balance==EPS_BALANCE_TWOSIDE) PetscCall(PetscViewerASCIIPrintf(viewer,", with its=%" PetscInt_FMT,eps->balance_its));
112:       if (eps->balance==EPS_BALANCE_TWOSIDE && eps->balance_cutoff!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer," and cutoff=%g",(double)eps->balance_cutoff));
113:       PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
114:       PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
115:     }
116:     PetscCall(PetscViewerASCIIPrintf(viewer,"  selected portion of the spectrum: "));
117:     PetscCall(SlepcSNPrintfScalar(str,sizeof(str),eps->target,PETSC_FALSE));
118:     PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
119:     if (!eps->which) PetscCall(PetscViewerASCIIPrintf(viewer,"not yet set\n"));
120:     else switch (eps->which) {
121:       case EPS_WHICH_USER:
122:         PetscCall(PetscViewerASCIIPrintf(viewer,"user defined\n"));
123:         break;
124:       case EPS_TARGET_MAGNITUDE:
125:         PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (in magnitude)\n",str));
126:         break;
127:       case EPS_TARGET_REAL:
128:         PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the real axis)\n",str));
129:         break;
130:       case EPS_TARGET_IMAGINARY:
131:         PetscCall(PetscViewerASCIIPrintf(viewer,"closest to target: %s (along the imaginary axis)\n",str));
132:         break;
133:       case EPS_LARGEST_MAGNITUDE:
134:         PetscCall(PetscViewerASCIIPrintf(viewer,"largest eigenvalues in magnitude\n"));
135:         break;
136:       case EPS_SMALLEST_MAGNITUDE:
137:         PetscCall(PetscViewerASCIIPrintf(viewer,"smallest eigenvalues in magnitude\n"));
138:         break;
139:       case EPS_LARGEST_REAL:
140:         PetscCall(PetscViewerASCIIPrintf(viewer,"largest real parts\n"));
141:         break;
142:       case EPS_SMALLEST_REAL:
143:         PetscCall(PetscViewerASCIIPrintf(viewer,"smallest real parts\n"));
144:         break;
145:       case EPS_LARGEST_IMAGINARY:
146:         PetscCall(PetscViewerASCIIPrintf(viewer,"largest imaginary parts\n"));
147:         break;
148:       case EPS_SMALLEST_IMAGINARY:
149:         PetscCall(PetscViewerASCIIPrintf(viewer,"smallest imaginary parts\n"));
150:         break;
151:       case EPS_ALL:
152:         if (eps->inta || eps->intb) PetscCall(PetscViewerASCIIPrintf(viewer,"all eigenvalues in interval [%g,%g]\n",(double)eps->inta,(double)eps->intb));
153:         else PetscCall(PetscViewerASCIIPrintf(viewer,"all eigenvalues in the region\n"));
154:         break;
155:     }
156:     PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
157:     if (eps->twosided && eps->problem_type!=EPS_HEP && eps->problem_type!=EPS_GHEP) PetscCall(PetscViewerASCIIPrintf(viewer,"  using two-sided variant (for left eigenvectors)\n"));
158:     if (eps->purify) PetscCall(PetscViewerASCIIPrintf(viewer,"  postprocessing eigenvectors with purification\n"));
159:     if (eps->trueres) PetscCall(PetscViewerASCIIPrintf(viewer,"  computing true residuals explicitly\n"));
160:     if (eps->trackall) PetscCall(PetscViewerASCIIPrintf(viewer,"  computing all residuals (for tracking convergence)\n"));
161:     if (eps->stop==EPS_STOP_THRESHOLD) PetscCall(PetscViewerASCIIPrintf(viewer,"  computing eigenvalues %s the threshold: %g%s\n",(eps->which==EPS_SMALLEST_MAGNITUDE||eps->which==EPS_SMALLEST_REAL)?"below":"above",(double)eps->thres,eps->threlative?" (relative)":""));
162:     if (eps->nev) PetscCall(PetscViewerASCIIPrintf(viewer,"  number of eigenvalues (nev): %" PetscInt_FMT "\n",eps->nev));
163:     PetscCall(PetscViewerASCIIPrintf(viewer,"  number of column vectors (ncv): %" PetscInt_FMT "\n",eps->ncv));
164:     PetscCall(PetscViewerASCIIPrintf(viewer,"  maximum dimension of projected problem (mpd): %" PetscInt_FMT "\n",eps->mpd));
165:     PetscCall(PetscViewerASCIIPrintf(viewer,"  maximum number of iterations: %" PetscInt_FMT "\n",eps->max_it));
166:     PetscCall(PetscViewerASCIIPrintf(viewer,"  tolerance: %g\n",(double)eps->tol));
167:     PetscCall(PetscViewerASCIIPrintf(viewer,"  convergence test: "));
168:     PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_FALSE));
169:     switch (eps->conv) {
170:     case EPS_CONV_ABS:
171:       PetscCall(PetscViewerASCIIPrintf(viewer,"absolute\n"));break;
172:     case EPS_CONV_REL:
173:       PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue\n"));break;
174:     case EPS_CONV_NORM:
175:       PetscCall(PetscViewerASCIIPrintf(viewer,"relative to the eigenvalue and matrix norms\n"));
176:       PetscCall(PetscViewerASCIIPrintf(viewer,"  computed matrix norms: norm(A)=%g",(double)eps->nrma));
177:       if (eps->isgeneralized) PetscCall(PetscViewerASCIIPrintf(viewer,", norm(B)=%g",(double)eps->nrmb));
178:       PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
179:       break;
180:     case EPS_CONV_USER:
181:       PetscCall(PetscViewerASCIIPrintf(viewer,"user-defined\n"));break;
182:     }
183:     PetscCall(PetscViewerASCIIUseTabs(viewer,PETSC_TRUE));
184:     if (eps->nini) PetscCall(PetscViewerASCIIPrintf(viewer,"  dimension of user-provided initial space: %" PetscInt_FMT "\n",PetscAbs(eps->nini)));
185:     if (eps->ninil) PetscCall(PetscViewerASCIIPrintf(viewer,"  dimension of user-provided left initial space: %" PetscInt_FMT "\n",PetscAbs(eps->ninil)));
186:     if (eps->nds) PetscCall(PetscViewerASCIIPrintf(viewer,"  dimension of user-provided deflation space: %" PetscInt_FMT "\n",PetscAbs(eps->nds)));
187:   } else PetscTryTypeMethod(eps,view,viewer);
188:   PetscCall(PetscObjectTypeCompareAny((PetscObject)eps,&isexternal,EPSARPACK,EPSBLOPEX,EPSELEMENTAL,EPSFEAST,EPSPRIMME,EPSSCALAPACK,EPSELPA,EPSEVSL,""));
189:   if (!isexternal) {
190:     PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO));
191:     if (!eps->V) PetscCall(EPSGetBV(eps,&eps->V));
192:     PetscCall(BVView(eps->V,viewer));
193:     if (eps->rg) {
194:       PetscCall(RGIsTrivial(eps->rg,&istrivial));
195:       if (!istrivial) PetscCall(RGView(eps->rg,viewer));
196:     }
197:     if (eps->useds) {
198:       if (!eps->ds) PetscCall(EPSGetDS(eps,&eps->ds));
199:       PetscCall(DSView(eps->ds,viewer));
200:     }
201:     PetscCall(PetscViewerPopFormat(viewer));
202:   }
203:   if (!eps->st) PetscCall(EPSGetST(eps,&eps->st));
204:   PetscCall(STView(eps->st,viewer));
205:   PetscFunctionReturn(PETSC_SUCCESS);
206: }

208: /*@
209:    EPSViewFromOptions - View (print) an `EPS` object based on values in the options database.

211:    Collective

213:    Input Parameters:
214: +  eps  - the linear eigensolver context
215: .  obj  - optional object that provides the options prefix used to query the options database
216: -  name - command line option

218:    Level: intermediate

220: .seealso: [](ch:eps), `EPSView()`, `EPSCreate()`, `PetscObjectViewFromOptions()`
221: @*/
222: PetscErrorCode EPSViewFromOptions(EPS eps,PetscObject obj,const char name[])
223: {
224:   PetscFunctionBegin;
226:   PetscCall(PetscObjectViewFromOptions((PetscObject)eps,obj,name));
227:   PetscFunctionReturn(PETSC_SUCCESS);
228: }

230: /*@
231:    EPSConvergedReasonView - Displays the reason an `EPS` solve converged or diverged.

233:    Collective

235:    Input Parameters:
236: +  eps - the linear eigensolver context
237: -  viewer - the viewer to display the reason

239:    Options Database Key:
240: .  -eps_converged_reason - print reason for convergence/divergence, and number of iterations

242:    Notes:
243:    Use `EPSConvergedReasonViewFromOptions()` to display the reason based on values
244:    in the options database.

246:    To change the format of the output call `PetscViewerPushFormat()` before this
247:    call. Use `PETSC_VIEWER_DEFAULT` for the default, or `PETSC_VIEWER_FAILED` to only
248:    display a reason if it fails. The latter can be set in the command line with
249:    `-eps_converged_reason ::failed`.

251:    Level: intermediate

253: .seealso: [](ch:eps), `EPSSetConvergenceTest()`, `EPSSetTolerances()`, `EPSGetIterationNumber()`, `EPSConvergedReasonViewFromOptions()`
254: @*/
255: PetscErrorCode EPSConvergedReasonView(EPS eps,PetscViewer viewer)
256: {
257:   PetscBool         isAscii;
258:   PetscViewerFormat format;
259:   PetscInt          nconv;

261:   PetscFunctionBegin;
262:   if (!viewer) viewer = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)eps));
263:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isAscii));
264:   if (isAscii) {
265:     PetscCall(PetscViewerGetFormat(viewer,&format));
266:     PetscCall(PetscViewerASCIIAddTab(viewer,((PetscObject)eps)->tablevel));
267:     PetscCall(EPS_GetActualConverged(eps,&nconv));
268:     if (eps->reason > 0 && format != PETSC_VIEWER_FAILED) PetscCall(PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve converged (%" PetscInt_FMT " eigenpair%s) due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",nconv,(nconv>1)?"s":"",EPSConvergedReasons[eps->reason],eps->its));
269:     else if (eps->reason <= 0) PetscCall(PetscViewerASCIIPrintf(viewer,"%s Linear eigensolve did not converge due to %s; iterations %" PetscInt_FMT "\n",((PetscObject)eps)->prefix?((PetscObject)eps)->prefix:"",EPSConvergedReasons[eps->reason],eps->its));
270:     PetscCall(PetscViewerASCIISubtractTab(viewer,((PetscObject)eps)->tablevel));
271:   }
272:   PetscFunctionReturn(PETSC_SUCCESS);
273: }

275: /*@
276:    EPSConvergedReasonViewFromOptions - Processes command line options to determine if/how
277:    the `EPS` converged reason is to be viewed.

279:    Collective

281:    Input Parameter:
282: .  eps - the linear eigensolver context

284:    Level: intermediate

286: .seealso: [](ch:eps), `EPSConvergedReasonView()`
287: @*/
288: PetscErrorCode EPSConvergedReasonViewFromOptions(EPS eps)
289: {
290:   PetscViewer       viewer;
291:   PetscBool         flg;
292:   static PetscBool  incall = PETSC_FALSE;
293:   PetscViewerFormat format;

295:   PetscFunctionBegin;
296:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
297:   incall = PETSC_TRUE;
298:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_converged_reason",&viewer,&format,&flg));
299:   if (flg) {
300:     PetscCall(PetscViewerPushFormat(viewer,format));
301:     PetscCall(EPSConvergedReasonView(eps,viewer));
302:     PetscCall(PetscViewerPopFormat(viewer));
303:     PetscCall(PetscViewerDestroy(&viewer));
304:   }
305:   incall = PETSC_FALSE;
306:   PetscFunctionReturn(PETSC_SUCCESS);
307: }

309: static PetscErrorCode EPSErrorView_ASCII(EPS eps,EPSErrorType etype,PetscViewer viewer)
310: {
311:   PetscReal      error;
312:   PetscScalar    kr,ki;
313:   PetscInt       i,j,nvals,nconv;
314:   PetscBool      requestednev=PETSC_TRUE; /* user has specified nev */

316:   PetscFunctionBegin;
317:   PetscCall(EPS_GetActualConverged(eps,&nconv));
318:   if (eps->which==EPS_ALL) {
319:     requestednev = PETSC_FALSE;
320:     if (eps->stop==EPS_STOP_THRESHOLD && eps->nev>0) requestednev = PETSC_TRUE;
321:   } else if (eps->stop==EPS_STOP_THRESHOLD && eps->nev==0) requestednev = PETSC_FALSE;
322:   nvals = requestednev? eps->nev: nconv;

324:   if (requestednev && nconv<nvals) {
325:     PetscCall(PetscViewerASCIIPrintf(viewer," Problem: less than %" PetscInt_FMT " eigenvalues converged\n\n",eps->nev));
326:     PetscFunctionReturn(PETSC_SUCCESS);
327:   }
328:   if (!requestednev && !nvals) {
329:     PetscCall(PetscViewerASCIIPrintf(viewer," No eigenvalues have been found\n\n"));
330:     PetscFunctionReturn(PETSC_SUCCESS);
331:   }
332:   for (i=0;i<nvals;i++) {
333:     PetscCall(EPSComputeError(eps,i,etype,&error));
334:     if (error>=5.0*eps->tol) {
335:       PetscCall(PetscViewerASCIIPrintf(viewer," Problem: some of the first %" PetscInt_FMT " relative errors are higher than the tolerance\n\n",nvals));
336:       PetscFunctionReturn(PETSC_SUCCESS);
337:     }
338:   }
339:   if (!requestednev) PetscCall(PetscViewerASCIIPrintf(viewer," Found %" PetscInt_FMT " eigenvalues, all of them computed up to the required tolerance:",nvals));
340:   else PetscCall(PetscViewerASCIIPrintf(viewer," All requested eigenvalues computed up to the required tolerance:"));
341:   for (i=0;i<=(nvals-1)/8;i++) {
342:     PetscCall(PetscViewerASCIIPrintf(viewer,"\n     "));
343:     for (j=0;j<PetscMin(8,nvals-8*i);j++) {
344:       PetscCall(EPSGetEigenvalue(eps,8*i+j,&kr,&ki));
345:       PetscCall(SlepcPrintEigenvalueASCII(viewer,kr,ki));
346:       if (8*i+j+1<nvals) PetscCall(PetscViewerASCIIPrintf(viewer,", "));
347:     }
348:   }
349:   PetscCall(PetscViewerASCIIPrintf(viewer,"\n\n"));
350:   PetscFunctionReturn(PETSC_SUCCESS);
351: }

353: static PetscErrorCode EPSErrorView_DETAIL(EPS eps,EPSErrorType etype,PetscViewer viewer)
354: {
355:   PetscReal      error,re,im;
356:   PetscScalar    kr,ki;
357:   PetscInt       i,nconv;
358:   char           ex[30],sep[]=" ---------------------- --------------------\n";

360:   PetscFunctionBegin;
361:   if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
362:   switch (etype) {
363:     case EPS_ERROR_ABSOLUTE:
364:       PetscCall(PetscSNPrintf(ex,sizeof(ex),"   ||Ax-k%sx||",eps->isgeneralized?"B":""));
365:       break;
366:     case EPS_ERROR_RELATIVE:
367:       PetscCall(PetscSNPrintf(ex,sizeof(ex),"||Ax-k%sx||/||kx||",eps->isgeneralized?"B":""));
368:       break;
369:     case EPS_ERROR_BACKWARD:
370:       PetscCall(PetscSNPrintf(ex,sizeof(ex),"    eta(x,k)"));
371:       break;
372:   }
373:   PetscCall(PetscViewerASCIIPrintf(viewer,"%s            k             %s\n%s",sep,ex,sep));
374:   PetscCall(EPS_GetActualConverged(eps,&nconv));
375:   for (i=0;i<nconv;i++) {
376:     PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
377:     PetscCall(EPSComputeError(eps,i,etype,&error));
378: #if defined(PETSC_USE_COMPLEX)
379:     re = PetscRealPart(kr);
380:     im = PetscImaginaryPart(kr);
381: #else
382:     re = kr;
383:     im = ki;
384: #endif
385:     if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"  % 9f%+9fi      %12g\n",(double)re,(double)im,(double)error));
386:     else PetscCall(PetscViewerASCIIPrintf(viewer,"    % 12f           %12g\n",(double)re,(double)error));
387:   }
388:   PetscCall(PetscViewerASCIIPrintf(viewer,"%s",sep));
389:   PetscFunctionReturn(PETSC_SUCCESS);
390: }

392: static PetscErrorCode EPSErrorView_MATLAB(EPS eps,EPSErrorType etype,PetscViewer viewer)
393: {
394:   PetscReal      error;
395:   PetscInt       i,nconv;
396:   const char     *name;

398:   PetscFunctionBegin;
399:   PetscCall(PetscObjectGetName((PetscObject)eps,&name));
400:   PetscCall(PetscViewerASCIIPrintf(viewer,"Error_%s = [\n",name));
401:   PetscCall(EPS_GetActualConverged(eps,&nconv));
402:   for (i=0;i<nconv;i++) {
403:     PetscCall(EPSComputeError(eps,i,etype,&error));
404:     PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)error));
405:   }
406:   PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
407:   PetscFunctionReturn(PETSC_SUCCESS);
408: }

410: /*@
411:    EPSErrorView - Displays the errors associated with the computed solution
412:    (as well as the eigenvalues).

414:    Collective

416:    Input Parameters:
417: +  eps    - the linear eigensolver context
418: .  etype  - error type
419: -  viewer - optional visualization context

421:    Options Database Keys:
422: +  -eps_error_absolute - print absolute errors of each eigenpair
423: .  -eps_error_relative - print relative errors of each eigenpair
424: -  -eps_error_backward - print backward errors of each eigenpair

426:    Notes:
427:    By default, this function checks the error of all eigenpairs and prints
428:    the eigenvalues if all of them are below the requested tolerance.
429:    If the viewer has format `PETSC_VIEWER_ASCII_INFO_DETAIL` then a table with
430:    eigenvalues and corresponding errors is printed.

432:    All the command-line options listed above admit an optional argument
433:    specifying the viewer type and options. For instance, use
434:    `-eps_error_relative :myerr.m:ascii_matlab` to save the errors in a file
435:    that can be executed in Matlab.

437:    Level: intermediate

439: .seealso: [](ch:eps), `EPSSolve()`, `EPSValuesView()`, `EPSVectorsView()`
440: @*/
441: PetscErrorCode EPSErrorView(EPS eps,EPSErrorType etype,PetscViewer viewer)
442: {
443:   PetscBool         isascii;
444:   PetscViewerFormat format;

446:   PetscFunctionBegin;
448:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
450:   PetscCheckSameComm(eps,1,viewer,3);
451:   EPSCheckSolved(eps,1);
452:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
453:   if (!isascii) PetscFunctionReturn(PETSC_SUCCESS);

455:   PetscCall(PetscViewerGetFormat(viewer,&format));
456:   switch (format) {
457:     case PETSC_VIEWER_DEFAULT:
458:     case PETSC_VIEWER_ASCII_INFO:
459:       PetscCall(EPSErrorView_ASCII(eps,etype,viewer));
460:       break;
461:     case PETSC_VIEWER_ASCII_INFO_DETAIL:
462:       PetscCall(EPSErrorView_DETAIL(eps,etype,viewer));
463:       break;
464:     case PETSC_VIEWER_ASCII_MATLAB:
465:       PetscCall(EPSErrorView_MATLAB(eps,etype,viewer));
466:       break;
467:     default:
468:       PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
469:   }
470:   PetscFunctionReturn(PETSC_SUCCESS);
471: }

473: /*@
474:    EPSErrorViewFromOptions - Processes command line options to determine if/how
475:    the errors of the computed solution are to be viewed.

477:    Collective

479:    Input Parameter:
480: .  eps - the linear eigensolver context

482:    Level: developer

484: .seealso: [](ch:eps), `EPSErrorView()`
485: @*/
486: PetscErrorCode EPSErrorViewFromOptions(EPS eps)
487: {
488:   PetscViewer       viewer;
489:   PetscBool         flg;
490:   static PetscBool  incall = PETSC_FALSE;
491:   PetscViewerFormat format;

493:   PetscFunctionBegin;
494:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
495:   incall = PETSC_TRUE;
496:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_absolute",&viewer,&format,&flg));
497:   if (flg) {
498:     PetscCall(PetscViewerPushFormat(viewer,format));
499:     PetscCall(EPSErrorView(eps,EPS_ERROR_ABSOLUTE,viewer));
500:     PetscCall(PetscViewerPopFormat(viewer));
501:     PetscCall(PetscViewerDestroy(&viewer));
502:   }
503:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_relative",&viewer,&format,&flg));
504:   if (flg) {
505:     PetscCall(PetscViewerPushFormat(viewer,format));
506:     PetscCall(EPSErrorView(eps,EPS_ERROR_RELATIVE,viewer));
507:     PetscCall(PetscViewerPopFormat(viewer));
508:     PetscCall(PetscViewerDestroy(&viewer));
509:   }
510:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_error_backward",&viewer,&format,&flg));
511:   if (flg) {
512:     PetscCall(PetscViewerPushFormat(viewer,format));
513:     PetscCall(EPSErrorView(eps,EPS_ERROR_BACKWARD,viewer));
514:     PetscCall(PetscViewerPopFormat(viewer));
515:     PetscCall(PetscViewerDestroy(&viewer));
516:   }
517:   incall = PETSC_FALSE;
518:   PetscFunctionReturn(PETSC_SUCCESS);
519: }

521: static PetscErrorCode EPSValuesView_DRAW(EPS eps,PetscViewer viewer)
522: {
523:   PetscDraw      draw;
524:   PetscDrawSP    drawsp;
525:   PetscReal      re,im;
526:   PetscScalar    kr,ki;
527:   PetscInt       i,nconv;

529:   PetscFunctionBegin;
530:   if (!eps->nconv) PetscFunctionReturn(PETSC_SUCCESS);
531:   PetscCall(PetscViewerDrawGetDraw(viewer,0,&draw));
532:   PetscCall(PetscDrawSetTitle(draw,"Computed Eigenvalues"));
533:   PetscCall(PetscDrawSPCreate(draw,1,&drawsp));
534:   PetscCall(EPS_GetActualConverged(eps,&nconv));
535:   for (i=0;i<nconv;i++) {
536:     PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
537: #if defined(PETSC_USE_COMPLEX)
538:     re = PetscRealPart(kr);
539:     im = PetscImaginaryPart(kr);
540: #else
541:     re = kr;
542:     im = ki;
543: #endif
544:     PetscCall(PetscDrawSPAddPoint(drawsp,&re,&im));
545:   }
546:   PetscCall(PetscDrawSPDraw(drawsp,PETSC_TRUE));
547:   PetscCall(PetscDrawSPSave(drawsp));
548:   PetscCall(PetscDrawSPDestroy(&drawsp));
549:   PetscFunctionReturn(PETSC_SUCCESS);
550: }

552: static PetscErrorCode EPSValuesView_BINARY(EPS eps,PetscViewer viewer)
553: {
554:   PetscScalar    kr,ki;
555: #if defined(PETSC_HAVE_COMPLEX)
556:   PetscInt       i,nconv;
557:   PetscComplex   *ev;
558: #endif

560:   PetscFunctionBegin;
561: #if defined(PETSC_HAVE_COMPLEX)
562:   PetscCall(EPS_GetActualConverged(eps,&nconv));
563:   PetscCall(PetscMalloc1(nconv,&ev));
564:   for (i=0;i<nconv;i++) {
565:     PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
566: #if defined(PETSC_USE_COMPLEX)
567:     ev[i] = kr;
568: #else
569:     ev[i] = PetscCMPLX(kr,ki);
570: #endif
571:   }
572:   PetscCall(PetscViewerBinaryWrite(viewer,ev,nconv,PETSC_COMPLEX));
573:   PetscCall(PetscFree(ev));
574: #endif
575:   PetscFunctionReturn(PETSC_SUCCESS);
576: }

578: #if defined(PETSC_HAVE_HDF5)
579: static PetscErrorCode EPSValuesView_HDF5(EPS eps,PetscViewer viewer)
580: {
581:   PetscInt       i,n,N,nconv;
582:   PetscScalar    eig;
583:   PetscMPIInt    rank;
584:   Vec            v;
585:   char           vname[30];
586:   const char     *ename;

588:   PetscFunctionBegin;
589:   PetscCallMPI(MPI_Comm_rank(PetscObjectComm((PetscObject)eps),&rank));
590:   PetscCall(EPS_GetActualConverged(eps,&nconv));
591:   N = nconv;
592:   n = rank? 0: N;
593:   /* create a vector containing the eigenvalues */
594:   PetscCall(VecCreateMPI(PetscObjectComm((PetscObject)eps),n,N,&v));
595:   PetscCall(PetscObjectGetName((PetscObject)eps,&ename));
596:   PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigr_%s",ename));
597:   PetscCall(PetscObjectSetName((PetscObject)v,vname));
598:   if (!rank) {
599:     for (i=0;i<nconv;i++) {
600:       PetscCall(EPSGetEigenvalue(eps,i,&eig,NULL));
601:       PetscCall(VecSetValue(v,i,eig,INSERT_VALUES));
602:     }
603:   }
604:   PetscCall(VecAssemblyBegin(v));
605:   PetscCall(VecAssemblyEnd(v));
606:   PetscCall(VecView(v,viewer));
607: #if !defined(PETSC_USE_COMPLEX)
608:   /* in real scalars write the imaginary part as a separate vector */
609:   PetscCall(PetscSNPrintf(vname,sizeof(vname),"eigi_%s",ename));
610:   PetscCall(PetscObjectSetName((PetscObject)v,vname));
611:   if (!rank) {
612:     for (i=0;i<nconv;i++) {
613:       PetscCall(EPSGetEigenvalue(eps,i,NULL,&eig));
614:       PetscCall(VecSetValue(v,i,eig,INSERT_VALUES));
615:     }
616:   }
617:   PetscCall(VecAssemblyBegin(v));
618:   PetscCall(VecAssemblyEnd(v));
619:   PetscCall(VecView(v,viewer));
620: #endif
621:   PetscCall(VecDestroy(&v));
622:   PetscFunctionReturn(PETSC_SUCCESS);
623: }
624: #endif

626: static PetscErrorCode EPSValuesView_ASCII(EPS eps,PetscViewer viewer)
627: {
628:   PetscInt       i,nconv;
629:   PetscScalar    kr,ki;

631:   PetscFunctionBegin;
632:   PetscCall(PetscViewerASCIIPrintf(viewer,"Eigenvalues = \n"));
633:   PetscCall(EPS_GetActualConverged(eps,&nconv));
634:   for (i=0;i<nconv;i++) {
635:     PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
636:     PetscCall(PetscViewerASCIIPrintf(viewer,"   "));
637:     PetscCall(SlepcPrintEigenvalueASCII(viewer,kr,ki));
638:     PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
639:   }
640:   PetscCall(PetscViewerASCIIPrintf(viewer,"\n"));
641:   PetscFunctionReturn(PETSC_SUCCESS);
642: }

644: static PetscErrorCode EPSValuesView_MATLAB(EPS eps,PetscViewer viewer)
645: {
646:   PetscInt       i,nconv;
647:   PetscReal      re,im;
648:   PetscScalar    kr,ki;
649:   const char     *name;

651:   PetscFunctionBegin;
652:   PetscCall(PetscObjectGetName((PetscObject)eps,&name));
653:   PetscCall(PetscViewerASCIIPrintf(viewer,"Lambda_%s = [\n",name));
654:   PetscCall(EPS_GetActualConverged(eps,&nconv));
655:   for (i=0;i<nconv;i++) {
656:     PetscCall(EPSGetEigenvalue(eps,i,&kr,&ki));
657: #if defined(PETSC_USE_COMPLEX)
658:     re = PetscRealPart(kr);
659:     im = PetscImaginaryPart(kr);
660: #else
661:     re = kr;
662:     im = ki;
663: #endif
664:     if (im!=0.0) PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e%+18.16ei\n",(double)re,(double)im));
665:     else PetscCall(PetscViewerASCIIPrintf(viewer,"%18.16e\n",(double)re));
666:   }
667:   PetscCall(PetscViewerASCIIPrintf(viewer,"];\n"));
668:   PetscFunctionReturn(PETSC_SUCCESS);
669: }

671: /*@
672:    EPSValuesView - Displays the computed eigenvalues in a viewer.

674:    Collective

676:    Input Parameters:
677: +  eps    - the linear eigensolver context
678: -  viewer - the viewer

680:    Options Database Key:
681: .  -eps_view_values - print computed eigenvalues

683:    Note:
684:    The command-line option listed above admits an optional argument
685:    specifying the viewer type and options. For instance, use
686:    `-eps_view_values :evals.m:ascii_matlab` to save the values in a file
687:    that can be executed in Matlab.
688:    See `PetscObjectViewFromOptions()` for more details.

690:    Level: intermediate

692: .seealso: [](ch:eps), `EPSSolve()`, `EPSVectorsView()`, `EPSErrorView()`
693: @*/
694: PetscErrorCode EPSValuesView(EPS eps,PetscViewer viewer)
695: {
696:   PetscBool         isascii,isdraw,isbinary;
697:   PetscViewerFormat format;
698: #if defined(PETSC_HAVE_HDF5)
699:   PetscBool         ishdf5;
700: #endif

702:   PetscFunctionBegin;
704:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
706:   PetscCheckSameComm(eps,1,viewer,2);
707:   EPSCheckSolved(eps,1);
708:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERDRAW,&isdraw));
709:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary));
710: #if defined(PETSC_HAVE_HDF5)
711:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERHDF5,&ishdf5));
712: #endif
713:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii));
714:   if (isdraw) PetscCall(EPSValuesView_DRAW(eps,viewer));
715:   else if (isbinary) PetscCall(EPSValuesView_BINARY(eps,viewer));
716: #if defined(PETSC_HAVE_HDF5)
717:   else if (ishdf5) PetscCall(EPSValuesView_HDF5(eps,viewer));
718: #endif
719:   else if (isascii) {
720:     PetscCall(PetscViewerGetFormat(viewer,&format));
721:     switch (format) {
722:       case PETSC_VIEWER_DEFAULT:
723:       case PETSC_VIEWER_ASCII_INFO:
724:       case PETSC_VIEWER_ASCII_INFO_DETAIL:
725:         PetscCall(EPSValuesView_ASCII(eps,viewer));
726:         break;
727:       case PETSC_VIEWER_ASCII_MATLAB:
728:         PetscCall(EPSValuesView_MATLAB(eps,viewer));
729:         break;
730:       default:
731:         PetscCall(PetscInfo(eps,"Unsupported viewer format %s\n",PetscViewerFormats[format]));
732:     }
733:   }
734:   PetscFunctionReturn(PETSC_SUCCESS);
735: }

737: /*@
738:    EPSValuesViewFromOptions - Processes command line options to determine if/how
739:    the computed eigenvalues are to be viewed.

741:    Collective

743:    Input Parameter:
744: .  eps - the linear eigensolver context

746:    Level: developer

748: .seealso: [](ch:eps), `EPSValuesView()`
749: @*/
750: PetscErrorCode EPSValuesViewFromOptions(EPS eps)
751: {
752:   PetscViewer       viewer;
753:   PetscBool         flg;
754:   static PetscBool  incall = PETSC_FALSE;
755:   PetscViewerFormat format;

757:   PetscFunctionBegin;
758:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
759:   incall = PETSC_TRUE;
760:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_values",&viewer,&format,&flg));
761:   if (flg) {
762:     PetscCall(PetscViewerPushFormat(viewer,format));
763:     PetscCall(EPSValuesView(eps,viewer));
764:     PetscCall(PetscViewerPopFormat(viewer));
765:     PetscCall(PetscViewerDestroy(&viewer));
766:   }
767:   incall = PETSC_FALSE;
768:   PetscFunctionReturn(PETSC_SUCCESS);
769: }

771: /*@
772:    EPSVectorsView - Outputs computed eigenvectors to a viewer.

774:    Collective

776:    Input Parameters:
777: +  eps    - the linear eigensolver context
778: -  viewer - the viewer

780:    Options Database Key:
781: .  -eps_view_vectors - output eigenvectors

783:    Notes:
784:    If PETSc was configured with real scalars, complex conjugate eigenvectors
785:    will be viewed as two separate real vectors, one containing the real part
786:    and another one containing the imaginary part.

788:    If left eigenvectors were computed with a two-sided eigensolver, the right
789:    and left eigenvectors are interleaved, that is, the vectors are output in
790:    the following order\: `X0, Y0, X1, Y1, X2, Y2, ...`

792:    The command-line option listed above admits an optional argument
793:    specifying the viewer type and options. For instance, use
794:    `-eps_view_vectors binary:evecs.bin` to save the vectors in a binary file.
795:    See `PetscObjectViewFromOptions()` for more details.

797:    Level: intermediate

799: .seealso: [](ch:eps), `EPSSolve()`, `EPSValuesView()`, `EPSErrorView()`
800: @*/
801: PetscErrorCode EPSVectorsView(EPS eps,PetscViewer viewer)
802: {
803:   PetscInt       i,nconv;
804:   Vec            xr,xi=NULL;

806:   PetscFunctionBegin;
808:   if (!viewer) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)eps),&viewer));
810:   PetscCheckSameComm(eps,1,viewer,2);
811:   EPSCheckSolved(eps,1);
812:   PetscCall(EPS_GetActualConverged(eps,&nconv));
813:   if (nconv) {
814:     PetscCall(BVCreateVec(eps->V,&xr));
815: #if !defined(PETSC_USE_COMPLEX)
816:     PetscCall(BVCreateVec(eps->V,&xi));
817: #endif
818:     for (i=0;i<nconv;i++) {
819:       PetscCall(EPSGetEigenvector(eps,i,xr,xi));
820:       PetscCall(SlepcViewEigenvector(viewer,xr,xi,"X",i,(PetscObject)eps));
821:       if (eps->twosided || eps->problem_type==EPS_BSE) {
822:         PetscCall(EPSGetLeftEigenvector(eps,i,xr,xi));
823:         PetscCall(SlepcViewEigenvector(viewer,xr,xi,"Y",i,(PetscObject)eps));
824:       }
825:     }
826:     PetscCall(VecDestroy(&xr));
827: #if !defined(PETSC_USE_COMPLEX)
828:     PetscCall(VecDestroy(&xi));
829: #endif
830:   }
831:   PetscFunctionReturn(PETSC_SUCCESS);
832: }

834: /*@
835:    EPSVectorsViewFromOptions - Processes command line options to determine if/how
836:    the computed eigenvectors are to be viewed.

838:    Collective

840:    Input Parameter:
841: .  eps - the linear eigensolver context

843:    Level: developer

845: .seealso: [](ch:eps), `EPSVectorsView()`
846: @*/
847: PetscErrorCode EPSVectorsViewFromOptions(EPS eps)
848: {
849:   PetscViewer       viewer;
850:   PetscBool         flg = PETSC_FALSE;
851:   static PetscBool  incall = PETSC_FALSE;
852:   PetscViewerFormat format;

854:   PetscFunctionBegin;
855:   if (incall) PetscFunctionReturn(PETSC_SUCCESS);
856:   incall = PETSC_TRUE;
857:   PetscCall(PetscOptionsCreateViewer(PetscObjectComm((PetscObject)eps),((PetscObject)eps)->options,((PetscObject)eps)->prefix,"-eps_view_vectors",&viewer,&format,&flg));
858:   if (flg) {
859:     PetscCall(PetscViewerPushFormat(viewer,format));
860:     PetscCall(EPSVectorsView(eps,viewer));
861:     PetscCall(PetscViewerPopFormat(viewer));
862:     PetscCall(PetscViewerDestroy(&viewer));
863:   }
864:   incall = PETSC_FALSE;
865:   PetscFunctionReturn(PETSC_SUCCESS);
866: }