SALSA Analysis Modules
Defines | Functions
lapack.c File Reference

Dense routines from Lapack: eigenvalue and schur stuff. More...

#include <stdlib.h>
#include "anamod.h"
#include "anamodsalsamodules.h"
#include "petscerror.h"
#include "petscmat.h"
Include dependency graph for lapack.c:

Go to the source code of this file.

Defines

#define ABS(x)   (x>0?x:-x)
#define EIGENVALUE(re, im)

Functions

void LAPACK_DGEESX (const char *jobvs, const char *sort, char *sel, const char *sense, const int *n, double *a, const int *lda, int *sdim, double *wr, double *wi, double *vs, const int *ldvs, double *rconde, double *rcondv, double *work, const int *lwork, int *iwork, const int *liwork, int *bwork, int *info)
static PetscErrorCode eigenvaluecomp (Mat A)
static PetscErrorCode Departure (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg)
static PetscErrorCode MaxEVbyMagRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg)
static PetscErrorCode MaxEVbyMagIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg)
static PetscErrorCode MinEVbyMagRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg)
static PetscErrorCode MinEVbyMagIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg)
static PetscErrorCode MaxEVbyRealRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg)
static PetscErrorCode MaxEVbyRealIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg)
static PetscErrorCode MaxEVbyImRe (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg)
static PetscErrorCode MaxEVbyImIm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg)
PetscErrorCode RegisterLapackModules ()

Detailed Description

Dense routines from Lapack: eigenvalue and schur stuff.

Dense calculations from Lapack

Usage

Activate this module with

PetscErrorCode RegisterLapackModules();

Compute these elements with

ComputeQuantity("lapack",<element>,A,(AnalysisItem*)&res,&flg);

Available elements are:

Definition in file lapack.c.


Define Documentation

#define ABS (   x)    (x>0?x:-x)

Referenced by eigenvaluecomp().

#define EIGENVALUE (   re,
  im 
)
Value:
{ \
          PetscReal mag = sqrt(re*re+im*im); \
          if (!maxbyreset) {maxbyreset = 1; maxbyrere = re; maxbyreim = im; } \
          else if (re>maxbyrere) { maxbyrere = re; maxbyreim = im; }    \
          if (!maxbyimset) {maxbyimset = 1; maxbyimre = re; maxbyimim = im; } \
          else if (im>maxbyimre) { maxbyimre = re; maxbyimim = im; }    \
          if (!maxbymagset) {maxbymagset = 1; maxmag = mag;             \
            maxbymagre = re; maxbymagim = im; }                         \
          else if (mag>maxmag) { maxbymagre = re; maxbymagim = im; }    \
          if (!minbymagset) {minbymagset = 1; minmag = mag;             \
            minbymagre = re; minbymagim = im; }                         \
          else if (mag<minmag) { minbymagre = re; minbymagim = im; }    \
      }

Referenced by eigenvaluecomp().


Function Documentation

static PetscErrorCode Departure ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscBool *  flg 
) [static]

Definition at line 196 of file lapack.c.

References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.

Referenced by RegisterLapackModules().

{
  Mat A = (Mat)prob;
  PetscReal v=0; PetscBool has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;
  ierr = GetDataID("lapack","departure",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetScalar
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) {
    rv->r = v;
  } else {
    ierr = eigenvaluecomp(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetScalar
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) {
      rv->r = v;
    } else SETERRQ(MPI_COMM_WORLD,1,"Should have computed departure by now");
  }
  
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode eigenvaluecomp ( Mat  A) [static]

This is the computational routine for lapack eigenvalue calculation.

Definition at line 41 of file lapack.c.

References ABS, EIGENVALUE, GetDataID(), HASTOEXIST, id, and LAPACK_DGEESX().

Referenced by Departure(), MaxEVbyImIm(), MaxEVbyImRe(), MaxEVbyMagIm(), MaxEVbyMagRe(), MaxEVbyRealIm(), MaxEVbyRealRe(), MinEVbyMagIm(), and MinEVbyMagRe().

{
  MPI_Comm comm; Mat Adense; 
  int maxbyreset=0,maxbyimset=0,maxbymagset=0,minbymagset=0;
  PetscReal dep=0,
    maxbymagre,maxbymagim, minbymagre,minbymagim,
    maxbyrere,maxbyreim, maxbyimre,maxbyimim, 
    maxmag,minmag;
  PetscErrorCode ierr;

  PetscFunctionBegin;

  ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr);
  ierr = MatConvert(A,MATSEQDENSE,MAT_INITIAL_MATRIX,&Adense); CHKERRQ(ierr);

  {
    char JOBVS = 'N', SORT = 'N', SENSE = 'N', SELECT='X';
    int N,LDA,SDIM=0,LDVS=1,LWORK,LIWORK,*IWORK,INFO;
    PetscScalar *AA,*VS=NULL,*RCONDE=NULL,*RCONDV=NULL;
    PetscReal *WR,*WI,*WORK,*BWORK=NULL;

    ierr = MatGetSize(A,&N,PETSC_NULL); CHKERRQ(ierr);
    LDA = N;
    ierr = MatGetArray(Adense,&AA); CHKERRQ(ierr);
    ierr = PetscMalloc(N*sizeof(PetscReal),&WR); CHKERRQ(ierr);
    ierr = PetscMalloc(N*sizeof(PetscReal),&WI); CHKERRQ(ierr);
    LWORK = 3*N+1;
    ierr = PetscMalloc(LWORK*sizeof(PetscReal),&WORK); CHKERRQ(ierr);
    LIWORK = 10;
    ierr = PetscMalloc(LIWORK*sizeof(int),&IWORK); CHKERRQ(ierr);
    CHKMEMQ;
    LAPACK_DGEESX(&JOBVS, &SORT, &SELECT, &SENSE, &N, AA, &LDA, &SDIM,
                  WR, WI, VS, &LDVS, RCONDE, RCONDV, WORK, &LWORK,
                  IWORK, &LIWORK, (int*)BWORK, &INFO);
    if (INFO<0) {
      INFO = -INFO;
      SETERRQ1(MPI_COMM_WORLD,1,"DGEESX Problem with parameter %d",INFO);
    } else if (INFO>0) {
      SETERRQ1(MPI_COMM_WORLD,1,"DGEESX error %d",INFO);
    }
    CHKMEMQ;
    {
      int i,j,rind=0; PetscReal f=0;
      for (i=0; i<N; i++) {
#define ABS(x) (x>0?x:-x)
#define EIGENVALUE(re,im) { \
          PetscReal mag = sqrt(re*re+im*im); \
          if (!maxbyreset) {maxbyreset = 1; maxbyrere = re; maxbyreim = im; } \
          else if (re>maxbyrere) { maxbyrere = re; maxbyreim = im; }    \
          if (!maxbyimset) {maxbyimset = 1; maxbyimre = re; maxbyimim = im; } \
          else if (im>maxbyimre) { maxbyimre = re; maxbyimim = im; }    \
          if (!maxbymagset) {maxbymagset = 1; maxmag = mag;             \
            maxbymagre = re; maxbymagim = im; }                         \
          else if (mag>maxmag) { maxbymagre = re; maxbymagim = im; }    \
          if (!minbymagset) {minbymagset = 1; minmag = mag;             \
            minbymagre = re; minbymagim = im; }                         \
          else if (mag<minmag) { minbymagre = re; minbymagim = im; }    \
      }
        PetscReal dd = AA[rind],db = ABS(AA[rind+1]),dt;
        dt = 1.e-14 * ABS(dd);
        /*printf("row %d, diagonal %e, below %e, test %e\n",i,dd,db,dt);*/
        if (db<dt) {
          /* case of 1x1 eigenvalue block */
          int ind=rind+N;
          for (j=i+1; j<N; j++) {
            PetscReal aa = AA[ind];
            f += aa*aa; ind += N;
          }
          EIGENVALUE(dd,0.);
        } else {
          /* case of 2x2 eigenvalue block */
          int ind;
          PetscReal dd = AA[rind],do1=AA[rind+1],do2=AA[rind+N],dn=AA[rind+N+1];
          if (ABS((dd-dn)/dd)>1.e-13)
            SETERRQ4(PETSC_COMM_SELF,1,"Strange 2x2 block %e,%e,%e,%e\n",dd,do1,do2,dn);
          EIGENVALUE(dd,sqrt(do1*do2));
          ind = rind+2*N;
          for (j=i+2; j<N; j++) {
            PetscReal aa = AA[ind];
            f += aa*aa; ind += N;
          }
          i++; rind += N+1; ind = rind+N;
          for (j=i+1; j<N; j++) {
            PetscReal aa = AA[ind];
            f += aa*aa; ind += N;
          }
        }
        rind += N+1;
      }
      dep = sqrt(f);
    }
    ierr = PetscFree(WR); CHKERRQ(ierr);
    ierr = PetscFree(WI); CHKERRQ(ierr);
  }

  ierr = MatDestroy(&Adense); CHKERRQ(ierr);

  {
    int id; PetscBool has;
    ierr = GetDataID("lapack","departure",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetScalar
      ((PetscObject)A,id,dep); CHKERRQ(ierr);

    ierr = GetDataID
      ("lapack","lambda-max-by-magnitude-re",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,maxbymagre); CHKERRQ(ierr);
    ierr = GetDataID
      ("lapack","lambda-max-by-magnitude-im",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,maxbymagim); CHKERRQ(ierr);

    ierr = GetDataID
      ("lapack","lambda-min-by-magnitude-re",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,minbymagre); CHKERRQ(ierr);
    ierr = GetDataID
      ("lapack","lambda-min-by-magnitude-im",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,minbymagim); CHKERRQ(ierr);

    ierr = GetDataID
      ("lapack","lambda-max-by-real-part-re",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,maxbyrere); CHKERRQ(ierr);
    ierr = GetDataID
      ("lapack","lambda-max-by-real-part-im",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,maxbyreim); CHKERRQ(ierr);

    ierr = GetDataID
      ("lapack","lambda-max-by-im-part-re",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,maxbyimre); CHKERRQ(ierr);
    ierr = GetDataID
      ("lapack","lambda-max-by-im-part-im",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,maxbyimim); CHKERRQ(ierr);
  }

  PetscFunctionReturn(0);
}

Here is the call graph for this function:

void LAPACK_DGEESX ( const char *  jobvs,
const char *  sort,
char *  sel,
const char *  sense,
const int *  n,
double *  a,
const int *  lda,
int *  sdim,
double *  wr,
double *  wi,
double *  vs,
const int *  ldvs,
double *  rconde,
double *  rcondv,
double *  work,
const int *  lwork,
int *  iwork,
const int *  liwork,
int *  bwork,
int *  info 
)

Referenced by eigenvaluecomp().

static PetscErrorCode MaxEVbyImIm ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscBool *  flg 
) [static]

Definition at line 397 of file lapack.c.

References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.

Referenced by RegisterLapackModules().

{
  Mat A = (Mat)prob;
  PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0;
  PetscFunctionBegin;
  ierr = GetDataID
    ("lapack","lambda-max-by-im-part-im",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = eigenvaluecomp(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetScalar
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) {
      rv->r = v;
    } else SETERRQ(MPI_COMM_WORLD,1,"Should have computed lambda-max-by-im-part-im by now");
  }
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode MaxEVbyImRe ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscBool *  flg 
) [static]

Definition at line 372 of file lapack.c.

References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.

Referenced by RegisterLapackModules().

{
  Mat A = (Mat)prob;
  PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0;
  PetscFunctionBegin;
  ierr = GetDataID
    ("lapack","lambda-max-by-im-part-re",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = eigenvaluecomp(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetScalar
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) {
      rv->r = v;
    } else SETERRQ(MPI_COMM_WORLD,1,"Should have computed lambda-max-by-im-part-re by now");
  }
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode MaxEVbyMagIm ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscBool *  flg 
) [static]

Definition at line 247 of file lapack.c.

References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.

Referenced by RegisterLapackModules().

{
  Mat A = (Mat)prob;
  PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0;
  PetscFunctionBegin;
  ierr = GetDataID
    ("lapack","lambda-max-by-magnitude-im",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = eigenvaluecomp(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetScalar
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) {
      rv->r = v;
    } else SETERRQ(MPI_COMM_WORLD,1,"Should have computed lambda-max-by-magnitude-im by now");
  }
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode MaxEVbyMagRe ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscBool *  flg 
) [static]

Definition at line 222 of file lapack.c.

References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.

Referenced by RegisterLapackModules().

{
  Mat A = (Mat)prob;
  PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0;
  PetscFunctionBegin;
  ierr = GetDataID
    ("lapack","lambda-max-by-magnitude-re",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = eigenvaluecomp(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetScalar
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) {
      rv->r = v;
    } else SETERRQ(MPI_COMM_WORLD,1,"Should have computed lambda-max-by-magnitude-re by now");
  }
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode MaxEVbyRealIm ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscBool *  flg 
) [static]

Definition at line 347 of file lapack.c.

References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.

Referenced by RegisterLapackModules().

{
  Mat A = (Mat)prob;
  PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0;
  PetscFunctionBegin;
  ierr = GetDataID
    ("lapack","lambda-max-by-real-part-im",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = eigenvaluecomp(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetScalar
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) {
      rv->r = v;
    } else SETERRQ(MPI_COMM_WORLD,1,"Should have computed lambda-max-by-real-part-im by now");
  }
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode MaxEVbyRealRe ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscBool *  flg 
) [static]

Definition at line 322 of file lapack.c.

References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.

Referenced by RegisterLapackModules().

{
  Mat A = (Mat)prob;
  PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0;
  PetscFunctionBegin;
  ierr = GetDataID
    ("lapack","lambda-max-by-real-part-re",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = eigenvaluecomp(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetScalar
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) {
      rv->r = v;
    } else SETERRQ(MPI_COMM_WORLD,1,"Should have computed lambda-max-by-real-part-re by now");
  }
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode MinEVbyMagIm ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscBool *  flg 
) [static]

Definition at line 297 of file lapack.c.

References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.

Referenced by RegisterLapackModules().

{
  Mat A = (Mat)prob;
  PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0;
  PetscFunctionBegin;
  ierr = GetDataID
    ("lapack","lambda-min-by-magnitude-im",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = eigenvaluecomp(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetScalar
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) {
      rv->r = v;
    } else SETERRQ(MPI_COMM_WORLD,1,"Should have computed lambda-min-by-magnitude-im by now");
  }
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode MinEVbyMagRe ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscBool *  flg 
) [static]

Definition at line 272 of file lapack.c.

References eigenvaluecomp(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.

Referenced by RegisterLapackModules().

{
  Mat A = (Mat)prob;
  PetscBool has; int id; PetscErrorCode ierr; PetscReal v = 0;
  PetscFunctionBegin;
  ierr = GetDataID
    ("lapack","lambda-min-by-magnitude-re",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = eigenvaluecomp(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetScalar
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) {
      rv->r = v;
    } else SETERRQ(MPI_COMM_WORLD,1,"Should have computed lambda-min-by-magnitude-re by now");
  }
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

PetscErrorCode RegisterLapackModules ( )

Declare norm-like modules that can be performed with lapack calculations.

Definition at line 424 of file lapack.c.

References ANALYSISDOUBLE, Departure(), MaxEVbyImIm(), MaxEVbyImRe(), MaxEVbyMagIm(), MaxEVbyMagRe(), MaxEVbyRealIm(), MaxEVbyRealRe(), MinEVbyMagIm(), MinEVbyMagRe(), and RegisterModule().

Referenced by AnaModRegisterSalsaModules(), and AnaModRegisterStandardModules().

{
  PetscErrorCode ierr;
  PetscFunctionBegin;

  ierr = RegisterModule
    ("lapack","departure",ANALYSISDOUBLE,&Departure); CHKERRQ(ierr);

  ierr = RegisterModule
    ("lapack","lambda-max-by-magnitude-re",
     ANALYSISDOUBLE,&MaxEVbyMagRe); CHKERRQ(ierr);
  ierr = RegisterModule
    ("lapack","lambda-max-by-magnitude-im",
     ANALYSISDOUBLE,&MaxEVbyMagIm); CHKERRQ(ierr);
  ierr = RegisterModule
    ("lapack","lambda-min-by-magnitude-re",
     ANALYSISDOUBLE,&MinEVbyMagRe); CHKERRQ(ierr);
  ierr = RegisterModule
    ("lapack","lambda-min-by-magnitude-im",
     ANALYSISDOUBLE,&MinEVbyMagIm); CHKERRQ(ierr);
  ierr = RegisterModule
    ("lapack","lambda-max-by-real-part-re",
     ANALYSISDOUBLE,&MaxEVbyRealRe); CHKERRQ(ierr);
  ierr = RegisterModule
    ("lapack","lambda-max-by-real-part-im",
     ANALYSISDOUBLE,&MaxEVbyRealIm); CHKERRQ(ierr);
  ierr = RegisterModule
    ("lapack","lambda-max-by-im-part-re",
     ANALYSISDOUBLE,&MaxEVbyImRe); CHKERRQ(ierr);
  ierr = RegisterModule
    ("lapack","lambda-max-by-im-part-im",
     ANALYSISDOUBLE,&MaxEVbyImIm); CHKERRQ(ierr);

  PetscFunctionReturn(0);
}

Here is the call graph for this function: