LCOV - code coverage report
Current view: top level - pep/interface - pepview.c (source / functions) Hit Total Coverage
Test: SLEPc Lines: 321 389 82.5 %
Date: 2024-04-19 00:31:36 Functions: 17 17 100.0 %
Legend: Lines: hit not hit

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

Generated by: LCOV version 1.14