SALSA Analysis Modules
|
Various estimates of how `wild' a matrix is. More...
#include <stdlib.h>
#include "anamod.h"
#include "anamodsalsamodules.h"
#include "petscerror.h"
#include "petscmat.h"
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) |
Various estimates of how `wild' a matrix is.
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.
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.
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().
{ Mat A = (Mat)prob; PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("variance","col-variability",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { ierr = ComputeVariability(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } PetscFunctionReturn(0); }
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().
{ int sign; Vec V; PetscScalar *v; MPI_Comm comm; PetscReal average,g_average,variance,g_variance; int local,nn,nz,np,gnn,gnz,gnp,row,id, N; PetscTruth has; PetscErrorCode ierr; PetscFunctionBegin; ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr); ierr = MatGetLocalSize(A,&local,PETSC_NULL); CHKERRQ(ierr); ierr = MatGetSize(A,&N,PETSC_NULL); CHKERRQ(ierr); ierr = VecCreateMPI(comm,local,PETSC_DECIDE,&V); CHKERRQ(ierr); ierr = MatGetDiagonal(A,V); CHKERRQ(ierr); ierr = VecGetArray(V,&v); CHKERRQ(ierr); average = 0; nn = np = nz = 0; for (row=0; row<local; row++) { if (v[row]>0) np++; else if (v[row]<0) nn++; else nz++; average += PetscAbsScalar(v[row]); } MPI_Allreduce((void*)&nn,(void*)&gnn,1,MPI_INT,MPI_SUM,comm); MPI_Allreduce((void*)&nz,(void*)&gnz,1,MPI_INT,MPI_SUM,comm); MPI_Allreduce((void*)&np,(void*)&gnp,1,MPI_INT,MPI_SUM,comm); if (np==N) sign = DIAGONAL_POSITIVE; else if (nn==0) { if (np==0) sign = DIAGONAL_ZERO; else sign = DIAGONAL_NONNEGATIVE;} else if (nn==N) sign = DIAGONAL_NEGATIVE; else if (np==0) sign = DIAGONAL_NONPOSITIVE; else sign = DIAGONAL_INDEFINITE; MPI_Allreduce((void*)&average,(void*)&g_average,1,MPIU_REAL,MPI_SUM,comm); average = g_average / N; variance = 0; for (row=0; row<local; row++) { PetscReal d = PetscAbsScalar(v[row])-average; variance += d*d; } MPI_Allreduce((void*)&variance,(void*)&g_variance,1,MPIU_REAL,MPI_SUM,comm); variance = sqrt(g_variance / N); ierr = GetDataID("variance","diagonal-average",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,average); CHKERRQ(ierr); ierr = GetDataID("variance","diagonal-variance",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,variance); CHKERRQ(ierr); ierr = GetDataID("variance","diagonal-sign",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,sign); CHKERRQ(ierr); ierr = VecRestoreArray(V,&v); CHKERRQ(ierr); ierr = VecDestroy(V); CHKERRQ(ierr); PetscFunctionReturn(0); }
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().
{ MPI_Comm comm; PetscReal rr,cr; PetscReal *rmax,*rmin,*cmax,*cmin,*cmaxx,*cminn,rr_local; int m,n,M,N,first,last,i,id,row; PetscTruth has; PetscErrorCode ierr; PetscFunctionBegin; ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr); ierr = MatGetSize(A,&M,&N); CHKERRQ(ierr); ierr = MatGetLocalSize(A,&m,&n); CHKERRQ(ierr); ierr = MatGetOwnershipRange(A,&first,&last); CHKERRQ(ierr); ierr = PetscMalloc(m*sizeof(PetscReal),&rmax); CHKERRQ(ierr); ierr = PetscMalloc(m*sizeof(PetscReal),&rmin); CHKERRQ(ierr); ierr = PetscMalloc(N*sizeof(PetscReal),&cmax); CHKERRQ(ierr); ierr = PetscMalloc(N*sizeof(PetscReal),&cmin); CHKERRQ(ierr); ierr = PetscMalloc(N*sizeof(PetscReal),&cmaxx); CHKERRQ(ierr); ierr = PetscMalloc(N*sizeof(PetscReal),&cminn); CHKERRQ(ierr); /* init; we take logs, so 1.e+5 is beyond infinity */ for (i=0; i<m; i++) {rmax[i] = 0; rmin[i] = 1.e+5;} for (i=0; i<N; i++) {cmax[i] = 0; cmin[i] = 1.e+5;} for (row=first; row<last; row++) { int irow=row-first,ncols,icol; const int *cols; const PetscScalar *vals; ierr = MatGetRow(A,row,&ncols,&cols,&vals); CHKERRQ(ierr); for (icol=0; icol<ncols; icol++) { int col=cols[icol]; if (vals[icol]) { PetscReal logval = (PetscReal) log10(PetscAbsScalar(vals[icol])); if (logval>rmax[irow]) rmax[irow] = logval; if (logval<rmin[irow]) rmin[irow] = logval; if (logval>cmax[col]) cmax[col] = logval; if (logval<cmin[col]) cmin[col] = logval; } } ierr = MatRestoreRow(A,row,&ncols,&cols,&vals); CHKERRQ(ierr); } /* * now get global row spread */ for (i=0; i<m; i++) rmax[i] -= rmin[i]; for (i=1; i<m; i++) if (rmax[i]>rmax[0]) rmax[0] = rmax[i]; rr_local = rmax[0]; MPI_Allreduce(&rr_local,&rr,1,MPIU_SCALAR,MPI_MAX,comm); ierr = GetDataID("variance","row-variability",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,rr); CHKERRQ(ierr); /* * global column spread: reduce the full length row of values, * then local max is global */ MPI_Allreduce(cmax,cmaxx,N,MPIU_SCALAR,MPI_MAX,comm); MPI_Allreduce(cmin,cminn,N,MPIU_SCALAR,MPI_MIN,comm); for (i=0; i<N; i++) cmaxx[i] -= cminn[i]; for (i=1; i<N; i++) if (cmaxx[i]>cmaxx[0]) cmaxx[0] = cmaxx[i]; cr = cmaxx[0]; ierr = GetDataID("variance","col-variability",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,cr); CHKERRQ(ierr); ierr = PetscFree(rmax); CHKERRQ(ierr); ierr = PetscFree(cmax); CHKERRQ(ierr); ierr = PetscFree(cmaxx); CHKERRQ(ierr); ierr = PetscFree(rmin); CHKERRQ(ierr); ierr = PetscFree(cmin); CHKERRQ(ierr); ierr = PetscFree(cminn); CHKERRQ(ierr); PetscFunctionReturn(0); }
PetscErrorCode DeRegisterVarianceModules | ( | void | ) |
Definition at line 333 of file variance.c.
Referenced by AnaModDeregisterSalsaModules().
{ PetscFunctionBegin; PetscFunctionReturn(0); }
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().
{ Mat A = (Mat)prob; PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("variance","diagonal-average",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { ierr = ComputeDiagonal(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } PetscFunctionReturn(0); }
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().
{ Mat A = (Mat)prob; int v = 0; PetscTruth has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("variance","diagonal-sign",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; else { ierr = ComputeDiagonal(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } PetscFunctionReturn(0); }
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().
{ Mat A = (Mat)prob; PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("variance","diagonal-variance",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { ierr = ComputeDiagonal(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } PetscFunctionReturn(0); }
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().
{ PetscErrorCode ierr; PetscFunctionBegin; ierr = RegisterModule ("variance","row-variability",ANALYSISDOUBLE,&RowVariability); CHKERRQ(ierr); ierr = RegisterModule ("variance","col-variability",ANALYSISDOUBLE,&ColVariability); CHKERRQ(ierr); ierr = RegisterModule ("variance","diagonal-average",ANALYSISDOUBLE,&DiagonalAverage); CHKERRQ(ierr); ierr = RegisterModule ("variance","diagonal-variance",ANALYSISDOUBLE,&DiagonalVariance); CHKERRQ(ierr); ierr = RegisterModule ("variance","diagonal-sign",ANALYSISINTEGER,&DiagonalSign); CHKERRQ(ierr); PetscFunctionReturn(0); }
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().
{ Mat A = (Mat)prob; PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("variance","row-variability",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; else { ierr = ComputeVariability(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } PetscFunctionReturn(0); }