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: }