SALSA Analysis Modules
iprs.c
Go to the documentation of this file.
00001 /*! \file iprs.c \ingroup categories
00002   \brief Quantities used in the UKY `Intelligent Preconditioner
00003   Recommendation System'
00004 
00005   \section iprs Intelligent Preconditioner Recommendation System
00006 
00007   The University of Kentucky `Intelligent Preconditioner Recommendation
00008   System' (Xu[2003]) uses a number of structural properties of the matrix
00009   that go beyond the properties in our own \ref structure module.
00010 
00011   \section Usage
00012 
00013   Activate this module with
00014 
00015   PetscErrorCode RegisterIprsModules();
00016 
00017   Compute these elements with
00018 
00019   ComputeQuantity("iprs",<element>,A,(AnalysisItem*)&res,&reslen,&flg);
00020 
00021   Available elements are:
00022   - "nnzup" : number of nonzeros in upper triangle
00023   - "nnzlow" : number of nonzeros in lower triangle
00024   - "nnzdia" : number of nonzeros on the diagonal
00025   - "nnz" : total number of nonzeros
00026   - "avgnnzprow" : average nonzeros per row
00027   - "avgdistfromdiag" : average distance of nonzeros to the diagonal
00028   - "relsymm" : relative symmetry, ratio of structural symmetric elements
00029     to total number of nonzeros
00030   - "upband" : bandwidth in the upper triangle
00031   - "loband" : bandwidth in the lower triangle
00032   - "n-nonzero-diags" : number of diagonals that have any nonzero element
00033   - "avg-diag-dist" : average distance of nonzero diagonal to main diagonal
00034   - "sigma-diag-dist" : standard deviation of diag-dist
00035 
00036   \section References
00037 \verbatim
00038 @techreport{Zhang:interim,
00039 author = {Shu{T}ing Xu and Eun-Joo Lee and Jun Zhang},
00040 title = {An Interim Analysis Report on Preconditioners and Matrices},
00041 institution = {University of Kentucky, Lexington; Department of 
00042     Computer Science},
00043 number = {388-03},
00044 year = {2003}
00045 }
00046 \endverbatim
00047 */
00048 #include <stdlib.h>
00049 #include "anamod.h"
00050 #include "anamatrix.h"
00051 #include "anamodsalsamodules.h"
00052 #include "petscerror.h"
00053 #include "petscmat.h"
00054 
00055 #undef __FUNCT__
00056 #define __FUNCT__ "computennz"
00057 /*!
00058   This is the computational routine for all IPRS modules.
00059   Lots of pointer chasing.
00060 */
00061 static PetscErrorCode computennz(Mat A)
00062 {
00063   MPI_Comm comm;
00064   int first,last,M,row,id,max_rowlen,*cols_tmp,
00065     nnzl,gnnzl,nnzu,gnnzu,nnzd,gnnzd, nnz,gnnz,avg, dist,gdist,
00066     p,gp, q,gq;
00067   IS colis,colsis=NULL;
00068   PetscBool has; PetscErrorCode ierr;
00069   PetscFunctionBegin;
00070 
00071   ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr);
00072   ierr = MatGetSize(A,&M,PETSC_NULL); CHKERRQ(ierr);
00073   ierr = MatGetOwnershipRange(A,&first,&last); CHKERRQ(ierr);
00074   ierr = HasQuantity
00075     (A,"structure","max-nnzeros-per-row",&has); CHKERRQ(ierr);
00076   if (has) {
00077     AnalysisItem r;
00078     ierr = MatrixComputeQuantity
00079       (A,"structure","max-nnzeros-per-row",&r,PETSC_NULL,PETSC_NULL); CHKERRQ(ierr);
00080     max_rowlen = r.i;
00081   } else {
00082     ierr = MatGetSize(A,PETSC_NULL,&max_rowlen); CHKERRQ(ierr);
00083   }
00084   ierr = PetscMalloc(max_rowlen*sizeof(int),&cols_tmp); CHKERRQ(ierr);
00085   dist = nnz = nnzl = nnzd = nnzu = p = q = 0;
00086   for (row=first; row<last; row++) {
00087     int ncols,icol; const int *cols;
00088     ierr = MatGetRow(A,row,&ncols,&cols,PETSC_NULL); CHKERRQ(ierr);
00089     if (ncols>max_rowlen) SETERRQ(MPI_COMM_SELF,1,"Max rowlen exceeded");
00090     ierr = PetscMemcpy(cols_tmp,cols,ncols*sizeof(int)); CHKERRQ(ierr);
00091     nnz += ncols;
00092     /* update left and right bandwidth */
00093     if (cols[0]<row && row-cols[0]>p) p = row-cols[0];
00094     if (cols[ncols-1]>row && cols[ncols-1]-row>q) q = cols[ncols-1]-row;
00095     /* compute summed distance from diagonal; count elements in D,L,U */
00096     for (icol=0; icol<ncols; icol++) {
00097       if (cols[icol]<row) {
00098         dist += row-cols[icol]; nnzl++;
00099       } else if (cols[icol]>row) {
00100         dist += cols[icol]-row; nnzu++;
00101       } else nnzd++;
00102     }
00103     /* analyze diagonals */
00104     for (icol=0; icol<ncols; icol++) cols_tmp[icol]-=row;
00105       CHKMEMQ;
00106     if (!colsis) {
00107       ierr = ISCreateGeneral
00108         (PETSC_COMM_SELF,ncols,cols_tmp,PETSC_COPY_VALUES,
00109          &colsis); CHKERRQ(ierr);
00110     } else {
00111       IS tmpis;
00112       ierr = ISCreateGeneral
00113         (PETSC_COMM_SELF,ncols,cols_tmp,PETSC_COPY_VALUES,
00114          &colis); CHKERRQ(ierr);
00115       CHKMEMQ;
00116       ierr = ISSort(colis); CHKERRQ(ierr);
00117       CHKMEMQ;
00118       ierr = ISSum(colsis,colis,&tmpis); CHKERRQ(ierr);
00119       CHKMEMQ;
00120       ierr = ISDestroy(&colsis); CHKERRQ(ierr);
00121       CHKMEMQ;
00122       ierr = ISDestroy(&colis); CHKERRQ(ierr);
00123       colsis = tmpis;
00124     }
00125     for (icol=0; icol<ncols; icol++) cols_tmp[icol]+=row;
00126     ierr = MatRestoreRow(A,row,&ncols,&cols,PETSC_NULL); CHKERRQ(ierr);
00127   }
00128   ierr = PetscFree(cols_tmp); CHKERRQ(ierr);
00129   ierr = MPI_Allreduce((void*)&nnzl,(void*)&gnnzl,1,MPI_INT,MPI_SUM,comm);
00130   ierr = MPI_Allreduce((void*)&nnzu,(void*)&gnnzu,1,MPI_INT,MPI_SUM,comm);
00131   ierr = MPI_Allreduce((void*)&nnzd,(void*)&gnnzd,1,MPI_INT,MPI_SUM,comm);
00132   ierr = MPI_Allreduce((void*)&nnz,(void*)&gnnz,1,MPI_INT,MPI_SUM,comm);
00133   if (gnnz!=(gnnzl+gnnzd+gnnzu))
00134     SETERRQ1(MPI_COMM_WORLD,1,"We missed %d elements",gnnz-(gnnzl+gnnzd+gnnzu));
00135   ierr = MPI_Allreduce((void*)&dist,(void*)&gdist,1,MPI_INT,MPI_SUM,comm);
00136   ierr = MPI_Allreduce((void*)&p,(void*)&gp,1,MPI_INT,MPI_MAX,comm);
00137   ierr = MPI_Allreduce((void*)&q,(void*)&gq,1,MPI_INT,MPI_MAX,comm);
00138 
00139   ierr = GetDataID("iprs","nnzlow",&id,&has); CHKERRQ(ierr);
00140   HASTOEXIST(has);
00141   ierr = PetscObjectComposedDataSetInt
00142     ((PetscObject)A,id,gnnzl); CHKERRQ(ierr);
00143   ierr = GetDataID("iprs","nnzup",&id,&has); CHKERRQ(ierr);
00144   HASTOEXIST(has);
00145   ierr = PetscObjectComposedDataSetInt
00146     ((PetscObject)A,id,gnnzu); CHKERRQ(ierr);
00147   ierr = GetDataID("iprs","nnzdia",&id,&has); CHKERRQ(ierr);
00148   HASTOEXIST(has);
00149   ierr = PetscObjectComposedDataSetInt
00150     ((PetscObject)A,id,gnnzd); CHKERRQ(ierr);
00151   ierr = GetDataID("iprs","nnz",&id,&has); CHKERRQ(ierr);
00152   HASTOEXIST(has);
00153   ierr = PetscObjectComposedDataSetInt
00154     ((PetscObject)A,id,gnnz); CHKERRQ(ierr);
00155 
00156   avg = nnz/M;
00157   ierr = GetDataID("iprs","avgnnzprow",&id,&has); CHKERRQ(ierr);
00158   HASTOEXIST(has);
00159   ierr = PetscObjectComposedDataSetInt
00160     ((PetscObject)A,id,avg); CHKERRQ(ierr);
00161   dist = gdist/M;
00162   ierr = GetDataID("iprs","avgdistfromdiag",&id,&has); CHKERRQ(ierr);
00163   HASTOEXIST(has);
00164   ierr = PetscObjectComposedDataSetInt
00165     ((PetscObject)A,id,dist); CHKERRQ(ierr);
00166 
00167   ierr = GetDataID("iprs","loband",&id,&has); CHKERRQ(ierr);
00168   HASTOEXIST(has);
00169   ierr = PetscObjectComposedDataSetInt((PetscObject)A,id,gp); CHKERRQ(ierr);
00170   ierr = GetDataID("iprs","upband",&id,&has); CHKERRQ(ierr);
00171   HASTOEXIST(has);
00172   ierr = PetscObjectComposedDataSetInt((PetscObject)A,id,gq); CHKERRQ(ierr);
00173 
00174   {
00175     int nd,i; const PetscInt *d; PetscReal avgdist,sigmadist;
00176     ierr = ISGetLocalSize(colsis,&nd); CHKERRQ(ierr);
00177     ierr = GetDataID("iprs","n-nonzero-diags",&id,&has); CHKERRQ(ierr);
00178     HASTOEXIST(has);
00179     ierr = PetscObjectComposedDataSetInt((PetscObject)A,id,nd); CHKERRQ(ierr);
00180     ierr = ISGetIndices(colsis,&d); CHKERRQ(ierr);
00181 
00182     avgdist = 0.;
00183     for (i=0; i<nd; i++) avgdist += PetscAbsScalar(d[i]);
00184     avgdist /= nd;
00185     ierr = GetDataID("iprs","avg-diag-dist",&id,&has); CHKERRQ(ierr);
00186     HASTOEXIST(has);
00187     ierr = PetscObjectComposedDataSetReal
00188       ((PetscObject)A,id,avgdist); CHKERRQ(ierr);
00189 
00190     sigmadist = 0.;
00191     for (i=0; i<nd; i++) sigmadist += (d[i]-avgdist)*(d[i]-avgdist);
00192     sigmadist /= nd;
00193     ierr = GetDataID("iprs","sigma-diag-dist",&id,&has); CHKERRQ(ierr);
00194     HASTOEXIST(has);
00195     ierr = PetscObjectComposedDataSetReal
00196       ((PetscObject)A,id,sigmadist); CHKERRQ(ierr);
00197 
00198     ierr = ISRestoreIndices(colsis,&d); CHKERRQ(ierr);
00199     ierr = ISDestroy(&colsis); CHKERRQ(ierr);
00200   }
00201   PetscFunctionReturn(0);
00202 }
00203 
00204 #undef __FUNCT__
00205 #define __FUNCT__ "NnzUp"
00206 /*!
00207   Compute the number of nonzeroes in the upper triangle of the matrix
00208 */
00209 static PetscErrorCode NnzUp
00210 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00211 {
00212   Mat A = (Mat)prob;
00213   int v=0; PetscBool has; int id; PetscErrorCode ierr;
00214   PetscFunctionBegin;
00215   ierr = GetDataID("iprs","nnzup",&id,&has); CHKERRQ(ierr);
00216   HASTOEXIST(has);
00217   ierr = PetscObjectComposedDataGetInt
00218     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00219   if (lv) *lv = 1;
00220   if (*flg) rv->i = v;
00221   else {
00222     ierr = computennz(A); CHKERRQ(ierr);
00223     ierr = PetscObjectComposedDataGetInt
00224       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00225     if (*flg) rv->i = v;
00226   }
00227   
00228   PetscFunctionReturn(0);
00229 }
00230 
00231 #undef __FUNCT__
00232 #define __FUNCT__ "NnzLow"
00233 static PetscErrorCode NnzLow
00234 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00235 {
00236   Mat A = (Mat)prob;
00237   int v=0; PetscBool has; int id; PetscErrorCode ierr;
00238   PetscFunctionBegin;
00239   ierr = GetDataID("iprs","nnzlow",&id,&has); CHKERRQ(ierr);
00240   HASTOEXIST(has);
00241   ierr = PetscObjectComposedDataGetInt
00242     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00243   if (lv) *lv = 1;
00244   if (*flg) rv->i = v;
00245   else {
00246     ierr = computennz(A); CHKERRQ(ierr);
00247     ierr = PetscObjectComposedDataGetInt
00248       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00249     if (*flg) rv->i = v;
00250   }
00251   
00252   PetscFunctionReturn(0);
00253 }
00254 
00255 #undef __FUNCT__
00256 #define __FUNCT__ "NnzDia"
00257 static PetscErrorCode NnzDia
00258 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00259 {
00260   Mat A = (Mat)prob;
00261   int v=0; PetscBool has; int id; PetscErrorCode ierr;
00262   PetscFunctionBegin;
00263   ierr = GetDataID("iprs","nnzdia",&id,&has); CHKERRQ(ierr);
00264   HASTOEXIST(has);
00265   ierr = PetscObjectComposedDataGetInt
00266     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00267   if (lv) *lv = 1;
00268   if (*flg) rv->i = v;
00269   else {
00270     ierr = computennz(A); CHKERRQ(ierr);
00271     ierr = PetscObjectComposedDataGetInt
00272       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00273     if (*flg) rv->i = v;
00274   }
00275   
00276   PetscFunctionReturn(0);
00277 }
00278 
00279 #undef __FUNCT__
00280 #define __FUNCT__ "Nnz"
00281 static PetscErrorCode Nnz
00282 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00283 {
00284   Mat A = (Mat)prob;
00285   int v=0; PetscBool has; int id; PetscErrorCode ierr;
00286   PetscFunctionBegin;
00287   ierr = GetDataID("iprs","nnz",&id,&has); CHKERRQ(ierr);
00288   HASTOEXIST(has);
00289   ierr = PetscObjectComposedDataGetInt
00290     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00291   if (lv) *lv = 1;
00292   if (*flg) rv->i = v;
00293   else {
00294     ierr = computennz(A); CHKERRQ(ierr);
00295     ierr = PetscObjectComposedDataGetInt
00296       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00297     if (*flg) rv->i = v;
00298   }
00299   
00300   PetscFunctionReturn(0);
00301 }
00302 
00303 #undef __FUNCT__
00304 #define __FUNCT__ "RelSymm"
00305 static PetscErrorCode RelSymm
00306 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00307 {
00308   Mat A = (Mat)prob;
00309   PetscReal v=0.; PetscBool has; int id,ql; PetscErrorCode ierr;
00310   PetscFunctionBegin;
00311   ierr = GetDataID("iprs","relsymm",&id,&has); CHKERRQ(ierr);
00312   HASTOEXIST(has);
00313   ierr = PetscObjectComposedDataGetReal
00314     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00315   if (lv) *lv = 1;
00316   if (*flg) rv->r = v;
00317   else {
00318     AnalysisItem q; int nnz,nnu;
00319     ierr = ComputeQuantity
00320       (prob,"structure","nnzeros",&q,&ql,flg); CHKERRQ(ierr);
00321     if (*flg) {
00322       nnz = q.i;
00323       ierr = ComputeQuantity
00324         (prob,"structure","n-struct-unsymm",&q,&ql,flg); CHKERRQ(ierr);
00325       if (*flg) {
00326         nnu = q.i;
00327         v = (nnz-nnu)/(1.*nnz);
00328         ierr = PetscObjectComposedDataSetReal
00329           ((PetscObject)A,id,v); CHKERRQ(ierr);
00330         rv->r = v;
00331       }
00332     }
00333   }
00334   PetscFunctionReturn(0);
00335 }
00336 
00337 #undef __FUNCT__
00338 #define __FUNCT__ "LoBand"
00339 static PetscErrorCode LoBand
00340 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00341 {
00342   Mat A = (Mat)prob;
00343   int v=0; PetscBool has; int id,ql; PetscErrorCode ierr;
00344   PetscFunctionBegin;
00345   ierr = GetDataID("iprs","loband",&id,&has); CHKERRQ(ierr);
00346   HASTOEXIST(has);
00347   ierr = PetscObjectComposedDataGetInt
00348     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00349   if (lv) *lv = 1;
00350   if (*flg) rv->i = v;
00351   else {
00352     AnalysisItem q;
00353     ierr = HasComputeModule("structure","left-bandwidth",flg); CHKERRQ(ierr);
00354     if (*flg) {
00355       ierr = ComputeQuantity
00356         (prob,"structure","left-bandwidth",&q,&ql,flg); CHKERRQ(ierr);
00357       if (*flg) {
00358         v = q.i;
00359         ierr = PetscObjectComposedDataSetInt
00360           ((PetscObject)A,id,v); CHKERRQ(ierr);
00361         rv->i = v;
00362       }
00363     } else {
00364       ierr = computennz(A); CHKERRQ(ierr);
00365       ierr = PetscObjectComposedDataGetInt
00366         ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00367       if (*flg) rv->i = v;
00368     }
00369   }
00370   
00371   PetscFunctionReturn(0);
00372 }
00373 
00374 #undef __FUNCT__
00375 #define __FUNCT__ "UpBand"
00376 static PetscErrorCode UpBand
00377 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00378 {
00379   Mat A = (Mat)prob;
00380   int v=0; PetscBool has; int id,ql; PetscErrorCode ierr;
00381   PetscFunctionBegin;
00382   ierr = GetDataID("iprs","upband",&id,&has); CHKERRQ(ierr);
00383   HASTOEXIST(has);
00384   ierr = PetscObjectComposedDataGetInt
00385     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00386   if (lv) *lv = 1;
00387   if (*flg) rv->i = v;
00388   else {
00389     ierr = HasComputeModule("structure","right-bandwidth",flg); CHKERRQ(ierr);
00390     if (*flg) {
00391       AnalysisItem q;
00392       ierr = ComputeQuantity
00393         (prob,"structure","right-bandwidth",&q,&ql,flg); CHKERRQ(ierr);
00394       if (*flg) {
00395         v = q.i;
00396         ierr = PetscObjectComposedDataSetInt
00397           ((PetscObject)A,id,v); CHKERRQ(ierr);
00398         rv->i = v;
00399       }
00400     } else {
00401       ierr = computennz(A); CHKERRQ(ierr);
00402       ierr = PetscObjectComposedDataGetInt
00403         ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00404       if (*flg) rv->i = v;
00405     }
00406   }
00407   
00408   PetscFunctionReturn(0);
00409 }
00410 
00411 #undef __FUNCT__
00412 #define __FUNCT__ "AvgNnzpRow"
00413 static PetscErrorCode AvgNnzpRow
00414 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00415 {
00416   Mat A = (Mat)prob;
00417   int v=0; PetscBool has; int id; PetscErrorCode ierr;
00418   PetscFunctionBegin;
00419   ierr = GetDataID("iprs","avgnnzprow",&id,&has); CHKERRQ(ierr);
00420   HASTOEXIST(has);
00421   ierr = PetscObjectComposedDataGetInt
00422     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00423   if (lv) *lv = 1;
00424   if (*flg) rv->i = v;
00425   else {
00426     ierr = computennz(A); CHKERRQ(ierr);
00427     ierr = PetscObjectComposedDataGetInt
00428       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00429     if (*flg) rv->i = v;
00430   }
00431   
00432   PetscFunctionReturn(0);
00433 }
00434 
00435 #if 0
00436 #undef __FUNCT__
00437 #define __FUNCT__ "AvgDistFromDiag"
00438 static PetscErrorCode AvgDistFromDiag
00439 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00440 {
00441   Mat A = (Mat)prob;
00442   int v=0; PetscBool has; int id; PetscErrorCode ierr;
00443   PetscFunctionBegin;
00444   ierr = GetDataID("iprs","avgdistfromdiag",&id,&has); CHKERRQ(ierr);
00445   HASTOEXIST(has);
00446   ierr = PetscObjectComposedDataGetInt
00447     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00448   if (lv) *lv = 1;
00449   if (*flg) rv->i = v;
00450   else {
00451     ierr = computennz(A); CHKERRQ(ierr);
00452     ierr = PetscObjectComposedDataGetInt
00453       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00454     if (*flg) rv->i = v;
00455   }
00456   
00457   PetscFunctionReturn(0);
00458 }
00459 #endif
00460 
00461 #undef __FUNCT__
00462 #define __FUNCT__ "NDiags"
00463 static PetscErrorCode NDiags
00464 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00465 {
00466   Mat A = (Mat)prob;
00467   int v=0; PetscBool has; int id; PetscErrorCode ierr;
00468   PetscFunctionBegin;
00469   ierr = GetDataID("iprs","n-nonzero-diags",&id,&has); CHKERRQ(ierr);
00470   HASTOEXIST(has);
00471   ierr = PetscObjectComposedDataGetInt
00472     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00473   if (lv) *lv = 1;
00474   if (*flg) rv->i = v;
00475   else {
00476     ierr = computennz(A); CHKERRQ(ierr);
00477     ierr = PetscObjectComposedDataGetInt
00478       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00479     if (*flg) rv->i = v;
00480   }
00481   
00482   PetscFunctionReturn(0);
00483 }
00484 
00485 #undef __FUNCT__
00486 #define __FUNCT__ "AvgDiagDist"
00487 static PetscErrorCode AvgDiagDist
00488 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00489 {
00490   Mat A = (Mat)prob;
00491   PetscReal v=0.; PetscBool has; int id; PetscErrorCode ierr;
00492   PetscFunctionBegin;
00493   ierr = GetDataID("iprs","avg-diag-dist",&id,&has); CHKERRQ(ierr);
00494   HASTOEXIST(has);
00495   ierr = PetscObjectComposedDataGetReal
00496     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00497   if (lv) *lv = 1;
00498   if (*flg) rv->r = v;
00499   else {
00500     ierr = computennz(A); CHKERRQ(ierr);
00501     ierr = PetscObjectComposedDataGetReal
00502       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00503     if (*flg) rv->r = v;
00504   }
00505   
00506   PetscFunctionReturn(0);
00507 }
00508 
00509 #undef __FUNCT__
00510 #define __FUNCT__ "SigmaDiagDist"
00511 static PetscErrorCode SigmaDiagDist
00512 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscBool *flg)
00513 {
00514   Mat A = (Mat)prob;
00515   PetscReal v=0.; PetscBool has; int id; PetscErrorCode ierr;
00516   PetscFunctionBegin;
00517   ierr = GetDataID("iprs","sigma-diag-dist",&id,&has); CHKERRQ(ierr);
00518   HASTOEXIST(has);
00519   ierr = PetscObjectComposedDataGetReal
00520     ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00521   if (lv) *lv = 1;
00522   if (*flg) rv->r = v;
00523   else {
00524     ierr = computennz(A); CHKERRQ(ierr);
00525     ierr = PetscObjectComposedDataGetReal
00526       ((PetscObject)A,id,v,*flg); CHKERRQ(ierr);
00527     if (*flg) rv->r = v;
00528   }
00529   
00530   PetscFunctionReturn(0);
00531 }
00532 
00533 #undef __FUNCT__
00534 #define __FUNCT__ "RegisterIprsModules"
00535 PetscErrorCode RegisterIprsModules(void)
00536 {
00537   PetscErrorCode ierr;
00538   PetscFunctionBegin;
00539 
00540   ierr = RegisterModule
00541     ("iprs","nnzup",ANALYSISINTEGER,&NnzUp); CHKERRQ(ierr);
00542   ierr = RegisterModule
00543     ("iprs","nnzlow",ANALYSISINTEGER,&NnzLow); CHKERRQ(ierr);
00544   ierr = RegisterModule
00545     ("iprs","nnzdia",ANALYSISINTEGER,&NnzDia); CHKERRQ(ierr);
00546 
00547   ierr = RegisterModule
00548     ("iprs","nnz",ANALYSISINTEGER,&Nnz); CHKERRQ(ierr);
00549   ierr = RegisterModule
00550     ("iprs","avgnnzprow",ANALYSISINTEGER,&AvgNnzpRow); CHKERRQ(ierr);
00551   ierr = RegisterModule
00552     ("iprs","avgdistfromdiag",ANALYSISINTEGER,&AvgNnzpRow); CHKERRQ(ierr);
00553 
00554   ierr = RegisterModule
00555     ("iprs","relsymm",ANALYSISDOUBLE,&RelSymm); CHKERRQ(ierr);
00556 
00557   ierr = RegisterModule
00558     ("iprs","upband",ANALYSISINTEGER,&UpBand); CHKERRQ(ierr);
00559   ierr = RegisterModule
00560     ("iprs","loband",ANALYSISINTEGER,&LoBand); CHKERRQ(ierr);
00561 
00562   ierr = RegisterModule
00563     ("iprs","n-nonzero-diags",ANALYSISINTEGER,&NDiags); CHKERRQ(ierr);
00564   ierr = RegisterModule
00565     ("iprs","avg-diag-dist",ANALYSISDOUBLE,&AvgDiagDist); CHKERRQ(ierr);
00566   ierr = RegisterModule
00567     ("iprs","sigma-diag-dist",ANALYSISDOUBLE,&SigmaDiagDist); CHKERRQ(ierr);
00568 
00569   PetscFunctionReturn(0);
00570 }
00571 
00572 #undef __FUNCT__
00573 #define __FUNCT__ "DeRegisterIprsModules"
00574 PetscErrorCode DeRegisterIprsModules(void)
00575 {
00576   PetscFunctionBegin;
00577   PetscFunctionReturn(0);
00578 }