variance.c File Reference

Various estimates of how `wild' a matrix is. More...

#include <stdlib.h>
#include "anamod.h"
#include "anamodsalsamodules.h"
#include "petscerror.h"
#include "petscmat.h"

Include dependency graph for variance.c:

Go to the source code of this file.

Functions

static PetscErrorCode ComputeVariability (Mat A)
static PetscErrorCode RowVariability (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode ColVariability (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode ComputeDiagonal (Mat A)
static PetscErrorCode DiagonalAverage (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode DiagonalVariance (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
static PetscErrorCode DiagonalSign (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg)
PetscErrorCode RegisterVarianceModules ()
PetscErrorCode DeRegisterVarianceModules (void)


Detailed Description

Various estimates of how `wild' a matrix is.

Measurements of element variance

This module contains various estimates of how different elements in a matrix are. These measures are completely heuristic, and do not relate to any known mathematical theory.

Usage

Activate this module with

PetscErrorCode RegisterVarianceModules();

Compute these elements with

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

Available elements are:

Definition in file variance.c.


Function Documentation

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

This routine relies on a call to ComputeVariability()

Definition at line 146 of file variance.c.

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

Referenced by RegisterVarianceModules().

00147 {
00148   Mat A = (Mat)prob;
00149   PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr;
00150   PetscFunctionBegin;
00151   ierr = GetDataID("variance","col-variability",&id,&has); CHKERRQ(ierr);
00152   HASTOEXIST(has);
00153   ierr = PetscObjectComposedDataGetReal
00154     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00155   if (*flg) rv->r = v;
00156   else {
00157     ierr = ComputeVariability(A); CHKERRQ(ierr);
00158     ierr = PetscObjectComposedDataGetReal
00159       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00160     if (*flg) rv->r = v;
00161   }
00162   
00163   PetscFunctionReturn(0);
00164 }

Here is the call graph for this function:

static PetscErrorCode ComputeDiagonal ( Mat  A  )  [static]

This is a computational routine; it computes the value of several modules

Definition at line 173 of file variance.c.

References DIAGONAL_INDEFINITE, DIAGONAL_NEGATIVE, DIAGONAL_NONNEGATIVE, DIAGONAL_NONPOSITIVE, DIAGONAL_POSITIVE, DIAGONAL_ZERO, GetDataID(), HASTOEXIST, and id.

Referenced by DiagonalAverage(), DiagonalSign(), and DiagonalVariance().

00174 {
00175   int sign;
00176   Vec V; PetscScalar *v; MPI_Comm comm;
00177   PetscReal average,g_average,variance,g_variance;
00178   int local,nn,nz,np,gnn,gnz,gnp,row,id, N;
00179   PetscTruth has; PetscErrorCode ierr;
00180 
00181   PetscFunctionBegin;
00182   ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr);
00183   ierr = MatGetLocalSize(A,&local,PETSC_NULL); CHKERRQ(ierr);
00184   ierr = MatGetSize(A,&N,PETSC_NULL); CHKERRQ(ierr);
00185   ierr = VecCreateMPI(comm,local,PETSC_DECIDE,&V); CHKERRQ(ierr);
00186   ierr = MatGetDiagonal(A,V); CHKERRQ(ierr);
00187   ierr = VecGetArray(V,&v); CHKERRQ(ierr);
00188 
00189   average = 0; nn = np = nz = 0;
00190   for (row=0; row<local; row++) {
00191     if (v[row]>0) np++; else if (v[row]<0) nn++; else nz++;
00192     average += PetscAbsScalar(v[row]);
00193   }
00194   MPI_Allreduce((void*)&nn,(void*)&gnn,1,MPI_INT,MPI_SUM,comm);
00195   MPI_Allreduce((void*)&nz,(void*)&gnz,1,MPI_INT,MPI_SUM,comm);
00196   MPI_Allreduce((void*)&np,(void*)&gnp,1,MPI_INT,MPI_SUM,comm);
00197   if (np==N) sign = DIAGONAL_POSITIVE;
00198   else if (nn==0) {
00199     if (np==0) sign = DIAGONAL_ZERO; else sign = DIAGONAL_NONNEGATIVE;}
00200   else if (nn==N) sign = DIAGONAL_NEGATIVE;
00201   else if (np==0) sign = DIAGONAL_NONPOSITIVE;
00202   else sign = DIAGONAL_INDEFINITE;
00203   MPI_Allreduce((void*)&average,(void*)&g_average,1,MPIU_REAL,MPI_SUM,comm);
00204   average = g_average / N; variance = 0;
00205   for (row=0; row<local; row++) {
00206     PetscReal d = PetscAbsScalar(v[row])-average;
00207     variance += d*d;
00208   }
00209   MPI_Allreduce((void*)&variance,(void*)&g_variance,1,MPIU_REAL,MPI_SUM,comm);
00210   variance = sqrt(g_variance / N);
00211   ierr = GetDataID("variance","diagonal-average",&id,&has); CHKERRQ(ierr);
00212   HASTOEXIST(has);
00213   ierr = PetscObjectComposedDataSetReal
00214     ((PetscObject)A,id,average); CHKERRQ(ierr);
00215 
00216   ierr = GetDataID("variance","diagonal-variance",&id,&has); CHKERRQ(ierr);
00217   HASTOEXIST(has);
00218   ierr = PetscObjectComposedDataSetReal
00219     ((PetscObject)A,id,variance); CHKERRQ(ierr);
00220 
00221   ierr = GetDataID("variance","diagonal-sign",&id,&has); CHKERRQ(ierr);
00222   HASTOEXIST(has);
00223   ierr = PetscObjectComposedDataSetInt
00224     ((PetscObject)A,id,sign); CHKERRQ(ierr);
00225 
00226   ierr = VecRestoreArray(V,&v); CHKERRQ(ierr);
00227   ierr = VecDestroy(V); CHKERRQ(ierr);
00228   PetscFunctionReturn(0);
00229 }

Here is the call graph for this function:

static PetscErrorCode ComputeVariability ( Mat  A  )  [static]

Variability in rows and columns.

This is a computational routine.

Definition at line 42 of file variance.c.

References GetDataID(), HASTOEXIST, and id.

Referenced by ColVariability(), and RowVariability().

00043 {
00044   MPI_Comm comm;
00045   PetscReal rr,cr;
00046   PetscReal *rmax,*rmin,*cmax,*cmin,*cmaxx,*cminn,rr_local;
00047   int m,n,M,N,first,last,i,id,row; 
00048   PetscTruth has; PetscErrorCode ierr;
00049 
00050   PetscFunctionBegin;
00051   ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr);
00052   ierr = MatGetSize(A,&M,&N); CHKERRQ(ierr);
00053   ierr = MatGetLocalSize(A,&m,&n); CHKERRQ(ierr);
00054   ierr = MatGetOwnershipRange(A,&first,&last); CHKERRQ(ierr);
00055   ierr = PetscMalloc(m*sizeof(PetscReal),&rmax); CHKERRQ(ierr);
00056   ierr = PetscMalloc(m*sizeof(PetscReal),&rmin); CHKERRQ(ierr);
00057   ierr = PetscMalloc(N*sizeof(PetscReal),&cmax); CHKERRQ(ierr);
00058   ierr = PetscMalloc(N*sizeof(PetscReal),&cmin); CHKERRQ(ierr);
00059   ierr = PetscMalloc(N*sizeof(PetscReal),&cmaxx); CHKERRQ(ierr);
00060   ierr = PetscMalloc(N*sizeof(PetscReal),&cminn); CHKERRQ(ierr);
00061 
00062   /* init; we take logs, so 1.e+5 is beyond infinity */
00063   for (i=0; i<m; i++) {rmax[i] = 0; rmin[i] = 1.e+5;}
00064   for (i=0; i<N; i++) {cmax[i] = 0; cmin[i] = 1.e+5;}
00065   for (row=first; row<last; row++) {
00066     int irow=row-first,ncols,icol; const int *cols; const PetscScalar *vals;
00067     ierr = MatGetRow(A,row,&ncols,&cols,&vals); CHKERRQ(ierr);
00068     for (icol=0; icol<ncols; icol++) {
00069       int col=cols[icol];
00070       if (vals[icol]) {
00071         PetscReal logval = (PetscReal) log10(PetscAbsScalar(vals[icol]));
00072         if (logval>rmax[irow]) rmax[irow] = logval;
00073         if (logval<rmin[irow]) rmin[irow] = logval;
00074         if (logval>cmax[col]) cmax[col] = logval;
00075         if (logval<cmin[col]) cmin[col] = logval;
00076       }
00077     }
00078     ierr = MatRestoreRow(A,row,&ncols,&cols,&vals); CHKERRQ(ierr);
00079   }
00080 
00081   /*
00082    * now get global row spread
00083    */
00084   for (i=0; i<m; i++) rmax[i] -= rmin[i];
00085   for (i=1; i<m; i++) if (rmax[i]>rmax[0]) rmax[0] = rmax[i];
00086   rr_local = rmax[0];
00087   MPI_Allreduce(&rr_local,&rr,1,MPIU_SCALAR,MPI_MAX,comm);
00088   ierr = GetDataID("variance","row-variability",&id,&has); CHKERRQ(ierr);
00089   HASTOEXIST(has);
00090   ierr = PetscObjectComposedDataSetReal
00091     ((PetscObject)A,id,rr); CHKERRQ(ierr);
00092 
00093   /*
00094    * global column spread: reduce the full length row of values,
00095    * then local max is global
00096    */
00097   MPI_Allreduce(cmax,cmaxx,N,MPIU_SCALAR,MPI_MAX,comm);
00098   MPI_Allreduce(cmin,cminn,N,MPIU_SCALAR,MPI_MIN,comm);
00099   for (i=0; i<N; i++) cmaxx[i] -= cminn[i];
00100   for (i=1; i<N; i++) if (cmaxx[i]>cmaxx[0]) cmaxx[0] = cmaxx[i];
00101   cr = cmaxx[0];
00102   ierr = GetDataID("variance","col-variability",&id,&has); CHKERRQ(ierr);
00103   HASTOEXIST(has);
00104   ierr = PetscObjectComposedDataSetReal
00105     ((PetscObject)A,id,cr); CHKERRQ(ierr);
00106 
00107   ierr = PetscFree(rmax); CHKERRQ(ierr);
00108   ierr = PetscFree(cmax); CHKERRQ(ierr);
00109   ierr = PetscFree(cmaxx); CHKERRQ(ierr);
00110   ierr = PetscFree(rmin); CHKERRQ(ierr);
00111   ierr = PetscFree(cmin); CHKERRQ(ierr);
00112   ierr = PetscFree(cminn); CHKERRQ(ierr);
00113   PetscFunctionReturn(0);
00114 }

Here is the call graph for this function:

PetscErrorCode DeRegisterVarianceModules ( void   ) 

Definition at line 333 of file variance.c.

Referenced by AnaModDeregisterSalsaModules().

00334 {
00335   PetscFunctionBegin;
00336   PetscFunctionReturn(0);
00337 }

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

This routine relies on a call to ComputeDiagonal()

Definition at line 237 of file variance.c.

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

Referenced by RegisterVarianceModules().

00238 {
00239   Mat A = (Mat)prob;
00240   PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr;
00241   PetscFunctionBegin;
00242   ierr = GetDataID("variance","diagonal-average",&id,&has); CHKERRQ(ierr);
00243   HASTOEXIST(has);
00244   ierr = PetscObjectComposedDataGetReal
00245     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00246   if (*flg) rv->r = v;
00247   else {
00248     ierr = ComputeDiagonal(A); CHKERRQ(ierr);
00249     ierr = PetscObjectComposedDataGetReal
00250       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00251     if (*flg) rv->r = v;
00252   }
00253   
00254   PetscFunctionReturn(0);
00255 }

Here is the call graph for this function:

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

This routine relies on a call to ComputeDiagonal()

Definition at line 289 of file variance.c.

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

Referenced by RegisterVarianceModules().

00290 {
00291   Mat A = (Mat)prob;
00292   int v = 0; PetscTruth has; int id; PetscErrorCode ierr;
00293   PetscFunctionBegin;
00294   ierr = GetDataID("variance","diagonal-sign",&id,&has); CHKERRQ(ierr);
00295   HASTOEXIST(has);
00296   ierr = PetscObjectComposedDataGetInt
00297     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00298   if (*flg) rv->i = v;
00299   else {
00300     ierr = ComputeDiagonal(A); CHKERRQ(ierr);
00301     ierr = PetscObjectComposedDataGetInt
00302       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00303     if (*flg) rv->i = v;
00304   }
00305   
00306   PetscFunctionReturn(0);
00307 }

Here is the call graph for this function:

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

This routine relies on a call to ComputeDiagonal()

Definition at line 263 of file variance.c.

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

Referenced by RegisterVarianceModules().

00264 {
00265   Mat A = (Mat)prob;
00266   PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr;
00267   PetscFunctionBegin;
00268   ierr = GetDataID("variance","diagonal-variance",&id,&has); CHKERRQ(ierr);
00269   HASTOEXIST(has);
00270   ierr = PetscObjectComposedDataGetReal
00271     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00272   if (*flg) rv->r = v;
00273   else {
00274     ierr = ComputeDiagonal(A); CHKERRQ(ierr);
00275     ierr = PetscObjectComposedDataGetReal
00276       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00277     if (*flg) rv->r = v;
00278   }
00279   
00280   PetscFunctionReturn(0);
00281 }

Here is the call graph for this function:

PetscErrorCode RegisterVarianceModules ( void   ) 

Definition at line 311 of file variance.c.

References ANALYSISDOUBLE, ANALYSISINTEGER, ColVariability(), DiagonalAverage(), DiagonalSign(), DiagonalVariance(), RegisterModule(), and RowVariability().

Referenced by AnaModRegisterSalsaModules(), and AnaModRegisterStandardModules().

00312 {
00313   PetscErrorCode ierr;
00314   PetscFunctionBegin;
00315 
00316   ierr = RegisterModule
00317     ("variance","row-variability",ANALYSISDOUBLE,&RowVariability); CHKERRQ(ierr);
00318   ierr = RegisterModule
00319     ("variance","col-variability",ANALYSISDOUBLE,&ColVariability); CHKERRQ(ierr);
00320 
00321   ierr = RegisterModule
00322     ("variance","diagonal-average",ANALYSISDOUBLE,&DiagonalAverage); CHKERRQ(ierr);
00323   ierr = RegisterModule
00324     ("variance","diagonal-variance",ANALYSISDOUBLE,&DiagonalVariance); CHKERRQ(ierr);
00325   ierr = RegisterModule
00326     ("variance","diagonal-sign",ANALYSISINTEGER,&DiagonalSign); CHKERRQ(ierr);
00327 
00328   PetscFunctionReturn(0);
00329 }

Here is the call graph for this function:

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

This routine relies on a call to ComputeVariability()

Definition at line 121 of file variance.c.

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

Referenced by RegisterVarianceModules().

00122 {
00123   Mat A = (Mat)prob;
00124   PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr;
00125   PetscFunctionBegin;
00126   ierr = GetDataID("variance","row-variability",&id,&has); CHKERRQ(ierr);
00127   HASTOEXIST(has);
00128   ierr = PetscObjectComposedDataGetReal
00129     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00130   if (*flg) rv->r = v;
00131   else {
00132     ierr = ComputeVariability(A); CHKERRQ(ierr);
00133     ierr = PetscObjectComposedDataGetReal
00134       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00135     if (*flg) rv->r = v;
00136   }
00137   
00138   PetscFunctionReturn(0);
00139 }

Here is the call graph for this function:


Generated on Sun Oct 4 04:01:25 2009 for SALSA Analysis Modules by  doxygen 1.5.9