SALSA Analysis Modules
Defines | Functions
simple.c File Reference

Norm-like properties. More...

#include <stdlib.h>
#include "anamod.h"
#include "anamodsalsamodules.h"
#include "petscerror.h"
#include "petscmat.h"
#include "src/mat/impls/aij/seq/aij.h"
Include dependency graph for simple.c:

Go to the source code of this file.

Defines

#define RC_NORM(i, j, v, w)
#define ASSERT3(v, s, i, j)   if (!(v)) SETERRQ3(1,"Violation of <%s> @ (%d,%d): %d",s,i,j);
#define ASSERT(v, s, i, j, t)   if (!(v)) SETERRQ4(1,"Violation of <%s> @ (%d,%d): %d",s,i,j,t);

Functions

static PetscErrorCode computetrace (Mat A)
static PetscErrorCode Trace (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode TraceAbs (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode norm1 (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode normF (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode normInf (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode compute_dd (Mat A, PetscTruth *flg)
static PetscErrorCode DiagonalDominance (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode MatSymmPartNormInf_seq (Mat A, PetscReal *sn_r, PetscReal *an_r, PetscReal *srn_r, PetscReal *arn_r, int *nun_r, PetscTruth *done)
static PetscErrorCode MatSymmPartNormInf (Mat A, PetscTruth *done)
static PetscErrorCode SymmetrySNorm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode SymmetryANorm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode SymmetryFSNorm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode SymmetryFANorm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode NUnstruct (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
PetscErrorCode RegisterSimpleModules ()
PetscErrorCode DeRegisterSimpleModules (void)

Detailed Description

Norm-like properties.

Simple (normlike) quantities

The "simple" module contains simple, normlike quantities, In general, these are properties that can be calculated in time proportional to the number of nonzeros of the matrix. Note that this excludes the 2-norm.

Usage

Activate this module with

PetscErrorCode RegisterSimpleModules();

Compute these elements with

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

Available elements are:

Definition in file simple.c.


Define Documentation

#define ASSERT (   v,
  s,
  i,
  j,
 
)    if (!(v)) SETERRQ4(1,"Violation of <%s> @ (%d,%d): %d",s,i,j,t);

Referenced by MatSymmPartNormInf_seq().

#define ASSERT3 (   v,
  s,
  i,
 
)    if (!(v)) SETERRQ3(1,"Violation of <%s> @ (%d,%d): %d",s,i,j);

Referenced by MatSymmPartNormInf_seq().

#define RC_NORM (   i,
  j,
  v,
 
)
Value:
{                                                                       \
    PetscReal sum = PetscAbsScalar(v+w)/2, dif = PetscAbsScalar(v-w)/2; \
    ASSERT3(!(i==j && v!=w),"same elements on diagonal",i,j);            \
    sronorm[i] += sum; if (i!=j) sronorm[j] += sum;                     \
    aronorm[i] += dif; if (i!=j) aronorm[j] += dif;                     \
    if (i==j) snorm += sum*sum;                                         \
    else {snorm += 2*sum*sum; anorm += 2*dif*dif;}                      \
  }

Referenced by MatSymmPartNormInf_seq().


Function Documentation

static PetscErrorCode compute_dd ( Mat  A,
PetscTruth *  flg 
) [static]

Compute the diagonal dominance of a matrix: minimum value of diagonal element (not absolute!) minus off-diagonal elements absolute.

Definition at line 246 of file simple.c.

References GetDataID(), HASTOEXIST, and id.

Referenced by DiagonalDominance().

{
  PetscReal r,gr;
  MPI_Comm comm;
  int first,last,m,row; PetscErrorCode ierr;

  PetscFunctionBegin;
  ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr);
  ierr = MatGetOwnershipRange(A,&first,&last); CHKERRQ(ierr);
  m = last-first;
  for (row=first; row<last; row++) {
    int col,ncols; const int *cols; 
    const PetscScalar *vals; PetscReal rr=0.;
    ierr = MatGetRow(A,row,&ncols,&cols,&vals); CHKERRQ(ierr);
    for (col=0; col<ncols; col++) {
      PetscReal v = vals[col], vr = PetscRealPart(v), va = PetscAbsScalar(v);
      if (cols[col]==row) rr +=vr; else rr -= va;
    }
    if (row==first) r = rr; else r = PetscMin(r,rr);
    ierr = MatRestoreRow(A,row,&ncols,&cols,&vals); CHKERRQ(ierr);
  }
  MPI_Allreduce((void*)&r,(void*)&gr,1,MPIU_REAL,MPI_MIN,comm);
  {
    int id; PetscTruth has;
    ierr = GetDataID("simple","diagonal-dominance",&id,&has); CHKERRQ(ierr);
    HASTOEXIST(has);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,gr); CHKERRQ(ierr);
  }
  *flg = PETSC_TRUE;
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode computetrace ( Mat  A) [static]

This is the computational routine for trace calculation. It computes both the trace and the absolute trace, that is, using absolute values of matrix elements.

A Petsc vector is created and destroyed.

Definition at line 63 of file simple.c.

References GetDataID(), HASTOEXIST, and id.

Referenced by Trace(), and TraceAbs().

{
  MPI_Comm comm; Vec D; PetscScalar *d;
  PetscTruth has; PetscScalar vsum,gvsum; PetscReal asum,gasum;
  int local_size,i,id;
  PetscErrorCode ierr;

  PetscFunctionBegin;

  ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr);
  ierr = MatGetLocalSize(A,&local_size,PETSC_NULL); CHKERRQ(ierr);
  ierr = VecCreateMPI(comm,local_size,PETSC_DECIDE,&D); CHKERRQ(ierr);
  ierr = MatGetDiagonal(A,D); CHKERRQ(ierr);
  ierr = VecGetArray(D,&d); CHKERRQ(ierr);
  vsum = 0; asum = 0;
  for (i=0; i<local_size; i++) {
    vsum += d[i]; asum += PetscAbsScalar(d[i]);
  }
  ierr = VecRestoreArray(D,&d); CHKERRQ(ierr);
  ierr = VecDestroy(D); CHKERRQ(ierr);
  ierr = MPI_Allreduce
    ((void*)&vsum,(void*)&gvsum,1,MPIU_SCALAR,MPI_SUM,comm); CHKERRQ(ierr);
  ierr = MPI_Allreduce
    ((void*)&asum,(void*)&gasum,1,MPI_DOUBLE,MPI_SUM,comm); CHKERRQ(ierr);

  ierr = GetDataID("simple","trace",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataSetScalar
    ((PetscObject)A,id,gvsum); CHKERRQ(ierr);
  ierr = GetDataID("simple","trace-abs",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataSetReal
    ((PetscObject)A,id,gasum); CHKERRQ(ierr);

  PetscFunctionReturn(0);
}

Here is the call graph for this function:

PetscErrorCode DeRegisterSimpleModules ( void  )

Definition at line 692 of file simple.c.

Referenced by AnaModDeregisterSalsaModules().

{
  PetscFunctionBegin;
  PetscFunctionReturn(0);
}
static PetscErrorCode DiagonalDominance ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Definition at line 280 of file simple.c.

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

Referenced by RegisterSimpleModules().

{
  Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;

  ierr = GetDataID("simple","diagonal-dominance",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = compute_dd(A,flg); CHKERRQ(ierr);
    if (*flg) {
      ierr = PetscObjectComposedDataGetReal
        ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
      if (*flg) rv->r = v;
    }
  }
  
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode MatSymmPartNormInf ( Mat  A,
PetscTruth *  done 
) [static]

This is the computational routine for all quantities relating to symmetric and anti-symmetric parts of the matrix.

This computation can only be done in single-processor mode. If the commandline option (see Commandline options) "-anamod-force sequential" is specified, the matrix is gathered on processor zero to compute this quantity.

See MatSymmPartNormInf_seq() for the actual computation.

Definition at line 443 of file simple.c.

References AnaModGetSequentialMatrix(), GetDataID(), HASTOEXIST, id, and MatSymmPartNormInf_seq().

Referenced by NUnstruct(), SymmetryANorm(), SymmetryFANorm(), SymmetryFSNorm(), and SymmetrySNorm().

{
  MPI_Comm comm; Mat Awork; int id,nun;
  PetscReal snorm,anorm,sronorm,aronorm;
  PetscTruth has,alloc_work,do_work,any_work; PetscErrorCode ierr;

  PetscFunctionBegin;

  ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr);
  ierr = AnaModGetSequentialMatrix
    (A,&Awork,&alloc_work,&do_work,&any_work); CHKERRQ(ierr);
  if (!any_work) {
    *done = PETSC_FALSE; PetscFunctionReturn(0);
  }
  if (do_work) {
    ierr = MatSymmPartNormInf_seq
      (Awork,&snorm,&anorm,&sronorm,&aronorm,&nun,done); CHKERRQ(ierr);
  }
  if (alloc_work) {
    ierr = MatDestroy(Awork); CHKERRQ(ierr);
  }

  MPI_Bcast((void*)done,1,MPI_INT,0,comm);
  MPI_Bcast((void*)&snorm,1,MPIU_REAL,0,comm);
  MPI_Bcast((void*)&anorm,1,MPIU_REAL,0,comm);
  MPI_Bcast((void*)&sronorm,1,MPIU_REAL,0,comm);
  MPI_Bcast((void*)&aronorm,1,MPIU_REAL,0,comm);
  MPI_Bcast((void*)&nun,1,MPI_INT,0,comm);

  /* infinity norms of symmetric / anti-symmetric part */
  ierr = GetDataID("simple","symmetry-snorm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataSetReal
    ((PetscObject)A,id,sronorm); CHKERRQ(ierr);

  ierr = GetDataID("simple","symmetry-anorm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataSetReal
    ((PetscObject)A,id,aronorm); CHKERRQ(ierr);

  /* frobenius norms of symmetric / anti-symmetric part */
  ierr = GetDataID("simple","symmetry-fsnorm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataSetReal
    ((PetscObject)A,id,snorm); CHKERRQ(ierr);

  ierr = GetDataID("simple","symmetry-fanorm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataSetReal
    ((PetscObject)A,id,anorm); CHKERRQ(ierr);

  ierr = GetDataID("structure","n-struct-unsymm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataSetInt
    ((PetscObject)A,id,nun); CHKERRQ(ierr);

  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode MatSymmPartNormInf_seq ( Mat  A,
PetscReal *  sn_r,
PetscReal *  an_r,
PetscReal *  srn_r,
PetscReal *  arn_r,
int *  nun_r,
PetscTruth *  done 
) [static]

This is where the real work of MatSymmPartNormInf() is done

A few temporary arrays of the size of a vector are allocated and freed. Otherwise this is just a whole bunch of pointer chasing.

Definition at line 316 of file simple.c.

References ASSERT, ASSERT3, and RC_NORM.

Referenced by MatSymmPartNormInf().

{
  Mat_SeqAIJ *aij;  
  int *adx,*bdx,*aii,*bii,*aptr,*bptr,ma,na,i,nun; PetscScalar *va,*vb;
  PetscReal snorm,anorm,*sronorm,*aronorm;
  PetscErrorCode ierr;
  PetscFunctionBegin;

  aij = (Mat_SeqAIJ *) A->data;
  ierr = MatGetLocalSize(A,&ma,&na); CHKERRQ(ierr);

  /* two arrays of pointers into the matrix */
  aii = bii = aij->i;
  adx = bdx = aij->j;
  va = vb = aij->a;
  ierr = PetscMalloc(ma*sizeof(int),&aptr); CHKERRQ(ierr);
  ierr = PetscMalloc(ma*sizeof(int),&bptr); CHKERRQ(ierr);
  for (i=0; i<ma; i++) aptr[i] = bptr[i] = aii[i];

  /* temp storage for the inf norm */
  ierr = PetscMalloc(ma*sizeof(PetscReal),&sronorm); CHKERRQ(ierr);
  ierr = PetscMalloc(ma*sizeof(PetscReal),&aronorm); CHKERRQ(ierr);
  ierr = PetscMemzero(sronorm,ma*sizeof(PetscReal)); CHKERRQ(ierr);
  ierr = PetscMemzero(aronorm,ma*sizeof(PetscReal)); CHKERRQ(ierr);

  /* update row & frobenius norms with a(i,j)=v and a(j,i)=w */
#if defined(ANAMODDEBUG)
#define RC_NORM(i,j,v,w) \
  printf("s/a norm calculation: (i,j)=(%d,%d), (vij,vji)=(%e,%e)\n",i,j,v,w); \
  {                                                                      \
    PetscReal sum = PetscAbsScalar(v+w)/2, dif = PetscAbsScalar(v-w)/2;  \
    ASSERT3(!(i==j && v!=w),"same elements on diagonal",i,j);            \
    sronorm[i] += sum; if (i!=j) sronorm[j] += sum;                      \
    aronorm[i] += dif; if (i!=j) aronorm[j] += dif;                      \
    if (i==j) snorm += sum*sum;                                          \
    else {snorm += 2*sum*sum; anorm += 2*dif*dif;}                       \
  }
#else
#define RC_NORM(i,j,v,w) \
  {                                                                     \
    PetscReal sum = PetscAbsScalar(v+w)/2, dif = PetscAbsScalar(v-w)/2; \
    ASSERT3(!(i==j && v!=w),"same elements on diagonal",i,j);            \
    sronorm[i] += sum; if (i!=j) sronorm[j] += sum;                     \
    aronorm[i] += dif; if (i!=j) aronorm[j] += dif;                     \
    if (i==j) snorm += sum*sum;                                         \
    else {snorm += 2*sum*sum; anorm += 2*dif*dif;}                      \
  }
#endif
#define ASSERT3(v,s,i,j) \
    if (!(v)) SETERRQ3(1,"Violation of <%s> @ (%d,%d): %d",s,i,j);
#define ASSERT(v,s,i,j,t) \
    if (!(v)) SETERRQ4(1,"Violation of <%s> @ (%d,%d): %d",s,i,j,t);

  snorm = anorm = 0; nun = 0;
  for (i=0; i<ma; i++) {
    int idc,idr; PetscScalar vc,vr;
    /* get the a-pointer into the upper triangle */
    while (aptr[i]<aii[i+1] && adx[aptr[i]]<i) aptr[i]++;
    ASSERT3(aptr[i]==aii[i+1]||adx[aptr[i]]>=i,"A in upper0",i,adx[aptr[i]]);
    /* incorporate L elements that are not structurally symmetric */
    while (bptr[i]<bii[i+1] && bdx[bptr[i]]<i) {
      vr = PetscAbsScalar(vb[bptr[i]]); nun++;
      RC_NORM(i,bdx[bptr[i]],vr,0.);
      bptr[i]++;
    }
    ASSERT3(bptr[i]==bii[i+1] || bdx[bptr[i]]>=i,"B part done",i,bdx[bptr[i]]);
    while (aptr[i]<aii[i+1]) {
      /* in row i we are at element j = idc @ adx( aptr[i] );
         vc = A( i, j ) */
      idc = adx[aptr[i]]; vc = va[aptr[i]];
      ASSERT(idc>=i,"A in upper",i,idc,i-idc);
      /* in row j=idc look at elements i' = idr @ bdx(bptr[idc]); vr = B(j,i')
         first catch up with structurally unsymmetric elements */
      idr = bdx[bptr[idc]]; vr = vb[bptr[idc]];
      while (idr<i) {
        RC_NORM(idc,idr,vr,0.);
        bptr[idc]++; idr = bdx[bptr[idc]]; vr = vb[bptr[idc]]; nun++;
      }     
      /* now see if (i,j) == (j,i') */
      if (idr==i) {
        /* structurally symmetric; process these elements of L and U */
        RC_NORM(i,idc,vc,vr);
        bptr[idc]++;
      } else {
        /* element (j,i') in L corresponds to a later row in U,
           so we have a structurally unsymmetric element in U */
        nun++; RC_NORM(i,idc,vc,0.);
      }
      aptr[i]++;
    }
    ASSERT(aptr[i]==aii[i+1],"Aptr has traversed row i",i+1,adx[aptr[i]],aii[i+1]);
  }

  ierr = PetscFree(aptr); CHKERRQ(ierr);
  ierr = PetscFree(bptr); CHKERRQ(ierr);

  snorm = sqrt(snorm); anorm = sqrt(anorm);

  for (i=1; i<ma; i++) {
    if (sronorm[i]>sronorm[0]) sronorm[0] = sronorm[i];
    if (aronorm[i]>aronorm[0]) aronorm[0] = aronorm[i];
  }

  *sn_r = snorm; *an_r = anorm; *srn_r = sronorm[0]; *arn_r = aronorm[0];
  *nun_r = nun;
  *done = PETSC_TRUE;

  ierr = PetscFree(sronorm); CHKERRQ(ierr);
  ierr = PetscFree(aronorm); CHKERRQ(ierr);

  PetscFunctionReturn(0);
}
static PetscErrorCode norm1 ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute the 1-norm of a matrix.

Definition at line 159 of file simple.c.

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

Referenced by RegisterSimpleModules().

{
  Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;

  ierr = GetDataID("simple","norm1",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = MatNorm(A,NORM_1,&v); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,v); CHKERRQ(ierr);
    *flg = PETSC_TRUE;
    rv->r = v;
  }
  
  *flg = PETSC_TRUE;
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode normF ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute the Frobenius norm of a matrix.

Definition at line 186 of file simple.c.

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

Referenced by RegisterSimpleModules().

{
  Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;

  ierr = GetDataID("simple","normF",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = MatNorm(A,NORM_FROBENIUS,&v); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,v); CHKERRQ(ierr);
    *flg = PETSC_TRUE;
    rv->r = v;
  }
  
  *flg = PETSC_TRUE;
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode normInf ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute the infinity-norm of a matrix.

Definition at line 213 of file simple.c.

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

Referenced by RegisterSimpleModules().

{
  Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;

  ierr = GetDataID("simple","normInf",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = MatNorm(A,NORM_INFINITY,&v); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataSetReal
      ((PetscObject)A,id,v); CHKERRQ(ierr);
    *flg = PETSC_TRUE;
    rv->r = v;
  }
  
  *flg = PETSC_TRUE;
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode NUnstruct ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute and store the number of structurally unsymmetric elements of a matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.

The "n-struct-unsymm" feature is declared to be in category "structure", rather than simple.

Definition at line 628 of file simple.c.

References GetDataID(), HASTOEXIST, AnalysisItem::i, id, and MatSymmPartNormInf().

Referenced by RegisterSimpleModules().

{
  Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;
  ierr = GetDataID("structure","n-struct-unsymm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetInt
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->i = v;
  else {
    ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr);
    if (*flg) {
      ierr = PetscObjectComposedDataGetInt
        ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
      if (*flg) rv->i = v;
    }
  }
  
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

PetscErrorCode RegisterSimpleModules ( void  )

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

Definition at line 655 of file simple.c.

References ANALYSISDOUBLE, ANALYSISINTEGER, DiagonalDominance(), norm1(), normF(), normInf(), NUnstruct(), RegisterModule(), SymmetryANorm(), SymmetryFANorm(), SymmetryFSNorm(), SymmetrySNorm(), Trace(), and TraceAbs().

Referenced by AnaModRegisterSalsaModules(), and AnaModRegisterStandardModules().

{
  PetscErrorCode ierr;
  PetscFunctionBegin;

  ierr = RegisterModule
    ("simple","trace",ANALYSISDOUBLE,&Trace); CHKERRQ(ierr);
  ierr = RegisterModule
    ("simple","trace-abs",ANALYSISDOUBLE,&TraceAbs); CHKERRQ(ierr);

  ierr = RegisterModule
    ("simple","norm1",ANALYSISDOUBLE,&norm1); CHKERRQ(ierr);
  ierr = RegisterModule
    ("simple","normInf",ANALYSISDOUBLE,&normInf); CHKERRQ(ierr);
  ierr = RegisterModule
    ("simple","normF",ANALYSISDOUBLE,&normF); CHKERRQ(ierr);

  ierr = RegisterModule
    ("simple","diagonal-dominance",ANALYSISDOUBLE,&DiagonalDominance); CHKERRQ(ierr);

  ierr = RegisterModule
    ("simple","symmetry-snorm",ANALYSISDOUBLE,&SymmetrySNorm); CHKERRQ(ierr);
  ierr = RegisterModule
    ("simple","symmetry-anorm",ANALYSISDOUBLE,&SymmetryANorm); CHKERRQ(ierr);
  ierr = RegisterModule
    ("simple","symmetry-fsnorm",ANALYSISDOUBLE,&SymmetryFSNorm); CHKERRQ(ierr);
  ierr = RegisterModule
    ("simple","symmetry-fanorm",ANALYSISDOUBLE,&SymmetryFANorm); CHKERRQ(ierr);

  ierr = RegisterModule
    ("structure","n-struct-unsymm",ANALYSISINTEGER,&NUnstruct); CHKERRQ(ierr);

  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode SymmetryANorm ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute and store the infinity norm of the antisymmetric part of the matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.

Definition at line 538 of file simple.c.

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

Referenced by RegisterSimpleModules().

{
    Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;
  ierr = GetDataID("simple","symmetry-anorm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr);
    if (*flg) {
      ierr = PetscObjectComposedDataGetReal
        ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
      if (*flg) rv->r = v;
    }
  }
  
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode SymmetryFANorm ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute and store the Frobenius norm of the antisymmetric part of the matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.

Definition at line 596 of file simple.c.

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

Referenced by RegisterSimpleModules().

{
    Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;
  ierr = GetDataID("simple","symmetry-fanorm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr);
    if (*flg) {
      ierr = PetscObjectComposedDataGetReal
        ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
      if (*flg) rv->r = v;
    }
  }
  
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode SymmetryFSNorm ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute and store the Frobenius norm of the symmetric part of the matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.

Definition at line 567 of file simple.c.

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

Referenced by RegisterSimpleModules().

{
    Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;
  ierr = GetDataID("simple","symmetry-fsnorm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr);
    if (*flg) {
      ierr = PetscObjectComposedDataGetReal
        ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
      if (*flg) rv->r = v;
    }
  }
  
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode SymmetrySNorm ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute and store the infinity norm of the symmetric part of the matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.

Definition at line 509 of file simple.c.

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

Referenced by RegisterSimpleModules().

{
  Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;
  ierr = GetDataID("simple","symmetry-snorm",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr);
    if (*flg) {
      ierr = PetscObjectComposedDataGetReal
        ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
      if (*flg) rv->r = v;
    }
  }
  
  PetscFunctionReturn(0);
}

Here is the call graph for this function:

static PetscErrorCode Trace ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute and store the trace of a matrix; the computation is done in computetrace().

Definition at line 106 of file simple.c.

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

Referenced by RegisterSimpleModules().

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

Here is the call graph for this function:

static PetscErrorCode TraceAbs ( AnaModNumericalProblem  prob,
AnalysisItem rv,
int *  lv,
PetscTruth *  flg 
) [static]

Compute and store the absolute trace of a matrix; the computation is done in computetrace().

Definition at line 135 of file simple.c.

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

Referenced by RegisterSimpleModules().

{
  Mat A = (Mat)prob;
  PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr;
  PetscFunctionBegin;
  ierr = GetDataID("simple","trace-abs",&id,&has); CHKERRQ(ierr);
  HASTOEXIST(has);
  ierr = PetscObjectComposedDataGetReal
    ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
  if (*flg) rv->r = v;
  else {
    ierr = computetrace(A); CHKERRQ(ierr);
    ierr = PetscObjectComposedDataGetReal
      ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
    if (*flg) rv->r = v;
  }
  
  PetscFunctionReturn(0);
}

Here is the call graph for this function: