SALSA Analysis Modules
|
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 PetscTruth 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 ierr = PetscMemcpy(cols_tmp,cols,ncols*sizeof(int)); CHKERRQ(ierr); 00090 nnz += ncols; 00091 /* left and right bandwidth */ 00092 if (cols[0]<row && row-cols[0]>p) p = row-cols[0]; 00093 if (cols[ncols-1]>row && cols[ncols-1]-row>q) q = cols[ncols-1]-row; 00094 /* distance from diagonal, count elements in L,U */ 00095 for (icol=0; icol<ncols; icol++) { 00096 if (cols[icol]<row) { 00097 dist += row-cols[icol]; nnzl++; 00098 } else if (cols[icol]>row) { 00099 dist += cols[icol]-row; nnzu++; 00100 } else nnzd++; 00101 } 00102 /* analyze diagonals */ 00103 for (icol=0; icol<ncols; icol++) cols_tmp[icol]-=row; 00104 if (!colsis) { 00105 ierr = ISCreateGeneral 00106 (PETSC_COMM_SELF,ncols,cols_tmp,&colsis); CHKERRQ(ierr); 00107 } else { 00108 IS tmpis; 00109 ierr = ISCreateGeneralWithArray 00110 (PETSC_COMM_SELF,ncols,cols_tmp,&colis); CHKERRQ(ierr); 00111 ierr = ISSum(colsis,colis,&tmpis); CHKERRQ(ierr); 00112 ierr = ISDestroy(colsis); CHKERRQ(ierr); 00113 ierr = ISDestroy(colis); CHKERRQ(ierr); 00114 colsis = tmpis; 00115 } 00116 for (icol=0; icol<ncols; icol++) cols_tmp[icol]+=row; 00117 ierr = MatRestoreRow(A,row,&ncols,&cols,PETSC_NULL); CHKERRQ(ierr); 00118 } 00119 ierr = PetscFree(cols_tmp); CHKERRQ(ierr); 00120 MPI_Allreduce((void*)&nnzl,(void*)&gnnzl,1,MPI_INT,MPI_SUM,comm); 00121 MPI_Allreduce((void*)&nnzu,(void*)&gnnzu,1,MPI_INT,MPI_SUM,comm); 00122 MPI_Allreduce((void*)&nnzd,(void*)&gnnzd,1,MPI_INT,MPI_SUM,comm); 00123 MPI_Allreduce((void*)&nnz,(void*)&gnnz,1,MPI_INT,MPI_SUM,comm); 00124 if (gnnz!=(gnnzl+gnnzd+gnnzu)) 00125 SETERRQ1(1,"We missed %d elements",gnnz-(gnnzl+gnnzd+gnnzu)); 00126 MPI_Allreduce((void*)&dist,(void*)&gdist,1,MPI_INT,MPI_SUM,comm); 00127 MPI_Allreduce((void*)&p,(void*)&gp,1,MPI_INT,MPI_MAX,comm); 00128 MPI_Allreduce((void*)&q,(void*)&gq,1,MPI_INT,MPI_MAX,comm); 00129 00130 ierr = GetDataID("iprs","nnzlow",&id,&has); CHKERRQ(ierr); 00131 HASTOEXIST(has); 00132 ierr = PetscObjectComposedDataSetInt 00133 ((PetscObject)A,id,gnnzl); CHKERRQ(ierr); 00134 ierr = GetDataID("iprs","nnzup",&id,&has); CHKERRQ(ierr); 00135 HASTOEXIST(has); 00136 ierr = PetscObjectComposedDataSetInt 00137 ((PetscObject)A,id,gnnzu); CHKERRQ(ierr); 00138 ierr = GetDataID("iprs","nnzdia",&id,&has); CHKERRQ(ierr); 00139 HASTOEXIST(has); 00140 ierr = PetscObjectComposedDataSetInt 00141 ((PetscObject)A,id,gnnzd); CHKERRQ(ierr); 00142 ierr = GetDataID("iprs","nnz",&id,&has); CHKERRQ(ierr); 00143 HASTOEXIST(has); 00144 ierr = PetscObjectComposedDataSetInt 00145 ((PetscObject)A,id,gnnz); CHKERRQ(ierr); 00146 00147 avg = nnz/M; 00148 ierr = GetDataID("iprs","avgnnzprow",&id,&has); CHKERRQ(ierr); 00149 HASTOEXIST(has); 00150 ierr = PetscObjectComposedDataSetInt 00151 ((PetscObject)A,id,avg); CHKERRQ(ierr); 00152 dist = gdist/M; 00153 ierr = GetDataID("iprs","avgdistfromdiag",&id,&has); CHKERRQ(ierr); 00154 HASTOEXIST(has); 00155 ierr = PetscObjectComposedDataSetInt 00156 ((PetscObject)A,id,dist); CHKERRQ(ierr); 00157 00158 ierr = GetDataID("iprs","loband",&id,&has); CHKERRQ(ierr); 00159 HASTOEXIST(has); 00160 ierr = PetscObjectComposedDataSetInt((PetscObject)A,id,gp); CHKERRQ(ierr); 00161 ierr = GetDataID("iprs","upband",&id,&has); CHKERRQ(ierr); 00162 HASTOEXIST(has); 00163 ierr = PetscObjectComposedDataSetInt((PetscObject)A,id,gq); CHKERRQ(ierr); 00164 00165 { 00166 int nd,i; const PetscInt *d; PetscReal avgdist,sigmadist; 00167 ierr = ISGetLocalSize(colsis,&nd); CHKERRQ(ierr); 00168 ierr = GetDataID("iprs","n-nonzero-diags",&id,&has); CHKERRQ(ierr); 00169 HASTOEXIST(has); 00170 ierr = PetscObjectComposedDataSetInt((PetscObject)A,id,nd); CHKERRQ(ierr); 00171 ierr = ISGetIndices(colsis,&d); CHKERRQ(ierr); 00172 00173 avgdist = 0.; 00174 for (i=0; i<nd; i++) avgdist += PetscAbsScalar(d[i]); 00175 avgdist /= nd; 00176 ierr = GetDataID("iprs","avg-diag-dist",&id,&has); CHKERRQ(ierr); 00177 HASTOEXIST(has); 00178 ierr = PetscObjectComposedDataSetReal 00179 ((PetscObject)A,id,avgdist); CHKERRQ(ierr); 00180 00181 sigmadist = 0.; 00182 for (i=0; i<nd; i++) sigmadist += (d[i]-avgdist)*(d[i]-avgdist); 00183 sigmadist /= nd; 00184 ierr = GetDataID("iprs","sigma-diag-dist",&id,&has); CHKERRQ(ierr); 00185 HASTOEXIST(has); 00186 ierr = PetscObjectComposedDataSetReal 00187 ((PetscObject)A,id,sigmadist); CHKERRQ(ierr); 00188 00189 ierr = ISRestoreIndices(colsis,&d); CHKERRQ(ierr); 00190 ierr = ISDestroy(colsis); CHKERRQ(ierr); 00191 } 00192 PetscFunctionReturn(0); 00193 } 00194 00195 #undef __FUNCT__ 00196 #define __FUNCT__ "NnzUp" 00197 /*! 00198 Compute the number of nonzeroes in the upper triangle of the matrix 00199 */ 00200 static PetscErrorCode NnzUp 00201 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00202 { 00203 Mat A = (Mat)prob; 00204 int v=0; PetscTruth has; int id; PetscErrorCode ierr; 00205 PetscFunctionBegin; 00206 ierr = GetDataID("iprs","nnzup",&id,&has); CHKERRQ(ierr); 00207 HASTOEXIST(has); 00208 ierr = PetscObjectComposedDataGetInt 00209 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00210 if (lv) *lv = 1; 00211 if (*flg) rv->i = v; 00212 else { 00213 ierr = computennz(A); CHKERRQ(ierr); 00214 ierr = PetscObjectComposedDataGetInt 00215 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00216 if (*flg) rv->i = v; 00217 } 00218 00219 PetscFunctionReturn(0); 00220 } 00221 00222 #undef __FUNCT__ 00223 #define __FUNCT__ "NnzLow" 00224 static PetscErrorCode NnzLow 00225 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00226 { 00227 Mat A = (Mat)prob; 00228 int v=0; PetscTruth has; int id; PetscErrorCode ierr; 00229 PetscFunctionBegin; 00230 ierr = GetDataID("iprs","nnzlow",&id,&has); CHKERRQ(ierr); 00231 HASTOEXIST(has); 00232 ierr = PetscObjectComposedDataGetInt 00233 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00234 if (lv) *lv = 1; 00235 if (*flg) rv->i = v; 00236 else { 00237 ierr = computennz(A); CHKERRQ(ierr); 00238 ierr = PetscObjectComposedDataGetInt 00239 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00240 if (*flg) rv->i = v; 00241 } 00242 00243 PetscFunctionReturn(0); 00244 } 00245 00246 #undef __FUNCT__ 00247 #define __FUNCT__ "NnzDia" 00248 static PetscErrorCode NnzDia 00249 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00250 { 00251 Mat A = (Mat)prob; 00252 int v=0; PetscTruth has; int id; PetscErrorCode ierr; 00253 PetscFunctionBegin; 00254 ierr = GetDataID("iprs","nnzdia",&id,&has); CHKERRQ(ierr); 00255 HASTOEXIST(has); 00256 ierr = PetscObjectComposedDataGetInt 00257 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00258 if (lv) *lv = 1; 00259 if (*flg) rv->i = v; 00260 else { 00261 ierr = computennz(A); CHKERRQ(ierr); 00262 ierr = PetscObjectComposedDataGetInt 00263 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00264 if (*flg) rv->i = v; 00265 } 00266 00267 PetscFunctionReturn(0); 00268 } 00269 00270 #undef __FUNCT__ 00271 #define __FUNCT__ "Nnz" 00272 static PetscErrorCode Nnz 00273 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00274 { 00275 Mat A = (Mat)prob; 00276 int v=0; PetscTruth has; int id; PetscErrorCode ierr; 00277 PetscFunctionBegin; 00278 ierr = GetDataID("iprs","nnz",&id,&has); CHKERRQ(ierr); 00279 HASTOEXIST(has); 00280 ierr = PetscObjectComposedDataGetInt 00281 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00282 if (lv) *lv = 1; 00283 if (*flg) rv->i = v; 00284 else { 00285 ierr = computennz(A); CHKERRQ(ierr); 00286 ierr = PetscObjectComposedDataGetInt 00287 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00288 if (*flg) rv->i = v; 00289 } 00290 00291 PetscFunctionReturn(0); 00292 } 00293 00294 #undef __FUNCT__ 00295 #define __FUNCT__ "RelSymm" 00296 static PetscErrorCode RelSymm 00297 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00298 { 00299 Mat A = (Mat)prob; 00300 PetscReal v=0.; PetscTruth has; int id,ql; PetscErrorCode ierr; 00301 PetscFunctionBegin; 00302 ierr = GetDataID("iprs","relsymm",&id,&has); CHKERRQ(ierr); 00303 HASTOEXIST(has); 00304 ierr = PetscObjectComposedDataGetReal 00305 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00306 if (lv) *lv = 1; 00307 if (*flg) rv->r = v; 00308 else { 00309 AnalysisItem q; int nnz,nnu; 00310 ierr = ComputeQuantity 00311 (prob,"structure","nnzeros",&q,&ql,flg); CHKERRQ(ierr); 00312 if (*flg) { 00313 nnz = q.i; 00314 ierr = ComputeQuantity 00315 (prob,"structure","n-struct-unsymm",&q,&ql,flg); CHKERRQ(ierr); 00316 if (*flg) { 00317 nnu = q.i; 00318 v = (nnz-nnu)/(1.*nnz); 00319 ierr = PetscObjectComposedDataSetReal 00320 ((PetscObject)A,id,v); CHKERRQ(ierr); 00321 rv->r = v; 00322 } 00323 } 00324 } 00325 PetscFunctionReturn(0); 00326 } 00327 00328 #undef __FUNCT__ 00329 #define __FUNCT__ "LoBand" 00330 static PetscErrorCode LoBand 00331 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00332 { 00333 Mat A = (Mat)prob; 00334 int v=0; PetscTruth has; int id,ql; PetscErrorCode ierr; 00335 PetscFunctionBegin; 00336 ierr = GetDataID("iprs","loband",&id,&has); CHKERRQ(ierr); 00337 HASTOEXIST(has); 00338 ierr = PetscObjectComposedDataGetInt 00339 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00340 if (lv) *lv = 1; 00341 if (*flg) rv->i = v; 00342 else { 00343 AnalysisItem q; 00344 ierr = HasComputeModule("structure","left-bandwidth",flg); CHKERRQ(ierr); 00345 if (*flg) { 00346 ierr = ComputeQuantity 00347 (prob,"structure","left-bandwidth",&q,&ql,flg); CHKERRQ(ierr); 00348 if (*flg) { 00349 v = q.i; 00350 ierr = PetscObjectComposedDataSetInt 00351 ((PetscObject)A,id,v); CHKERRQ(ierr); 00352 rv->i = v; 00353 } 00354 } else { 00355 ierr = computennz(A); CHKERRQ(ierr); 00356 ierr = PetscObjectComposedDataGetInt 00357 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00358 if (*flg) rv->i = v; 00359 } 00360 } 00361 00362 PetscFunctionReturn(0); 00363 } 00364 00365 #undef __FUNCT__ 00366 #define __FUNCT__ "UpBand" 00367 static PetscErrorCode UpBand 00368 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00369 { 00370 Mat A = (Mat)prob; 00371 int v=0; PetscTruth has; int id,ql; PetscErrorCode ierr; 00372 PetscFunctionBegin; 00373 ierr = GetDataID("iprs","upband",&id,&has); CHKERRQ(ierr); 00374 HASTOEXIST(has); 00375 ierr = PetscObjectComposedDataGetInt 00376 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00377 if (lv) *lv = 1; 00378 if (*flg) rv->i = v; 00379 else { 00380 ierr = HasComputeModule("structure","right-bandwidth",flg); CHKERRQ(ierr); 00381 if (*flg) { 00382 AnalysisItem q; 00383 ierr = ComputeQuantity 00384 (prob,"structure","right-bandwidth",&q,&ql,flg); CHKERRQ(ierr); 00385 if (*flg) { 00386 v = q.i; 00387 ierr = PetscObjectComposedDataSetInt 00388 ((PetscObject)A,id,v); CHKERRQ(ierr); 00389 rv->i = v; 00390 } 00391 } else { 00392 ierr = computennz(A); CHKERRQ(ierr); 00393 ierr = PetscObjectComposedDataGetInt 00394 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00395 if (*flg) rv->i = v; 00396 } 00397 } 00398 00399 PetscFunctionReturn(0); 00400 } 00401 00402 #undef __FUNCT__ 00403 #define __FUNCT__ "AvgNnzpRow" 00404 static PetscErrorCode AvgNnzpRow 00405 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00406 { 00407 Mat A = (Mat)prob; 00408 int v=0; PetscTruth has; int id; PetscErrorCode ierr; 00409 PetscFunctionBegin; 00410 ierr = GetDataID("iprs","avgnnzprow",&id,&has); CHKERRQ(ierr); 00411 HASTOEXIST(has); 00412 ierr = PetscObjectComposedDataGetInt 00413 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00414 if (lv) *lv = 1; 00415 if (*flg) rv->i = v; 00416 else { 00417 ierr = computennz(A); CHKERRQ(ierr); 00418 ierr = PetscObjectComposedDataGetInt 00419 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00420 if (*flg) rv->i = v; 00421 } 00422 00423 PetscFunctionReturn(0); 00424 } 00425 00426 #undef __FUNCT__ 00427 #define __FUNCT__ "AvgDistFromDiag" 00428 static PetscErrorCode AvgDistFromDiag 00429 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00430 { 00431 Mat A = (Mat)prob; 00432 int v=0; PetscTruth has; int id; PetscErrorCode ierr; 00433 PetscFunctionBegin; 00434 ierr = GetDataID("iprs","avgdistfromdiag",&id,&has); CHKERRQ(ierr); 00435 HASTOEXIST(has); 00436 ierr = PetscObjectComposedDataGetInt 00437 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00438 if (lv) *lv = 1; 00439 if (*flg) rv->i = v; 00440 else { 00441 ierr = computennz(A); CHKERRQ(ierr); 00442 ierr = PetscObjectComposedDataGetInt 00443 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00444 if (*flg) rv->i = v; 00445 } 00446 00447 PetscFunctionReturn(0); 00448 } 00449 00450 #undef __FUNCT__ 00451 #define __FUNCT__ "NDiags" 00452 static PetscErrorCode NDiags 00453 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00454 { 00455 Mat A = (Mat)prob; 00456 int v=0; PetscTruth has; int id; PetscErrorCode ierr; 00457 PetscFunctionBegin; 00458 ierr = GetDataID("iprs","n-nonzero-diags",&id,&has); CHKERRQ(ierr); 00459 HASTOEXIST(has); 00460 ierr = PetscObjectComposedDataGetInt 00461 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00462 if (lv) *lv = 1; 00463 if (*flg) rv->i = v; 00464 else { 00465 ierr = computennz(A); CHKERRQ(ierr); 00466 ierr = PetscObjectComposedDataGetInt 00467 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00468 if (*flg) rv->i = v; 00469 } 00470 00471 PetscFunctionReturn(0); 00472 } 00473 00474 #undef __FUNCT__ 00475 #define __FUNCT__ "AvgDiagDist" 00476 static PetscErrorCode AvgDiagDist 00477 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00478 { 00479 Mat A = (Mat)prob; 00480 PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr; 00481 PetscFunctionBegin; 00482 ierr = GetDataID("iprs","avg-diag-dist",&id,&has); CHKERRQ(ierr); 00483 HASTOEXIST(has); 00484 ierr = PetscObjectComposedDataGetReal 00485 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00486 if (lv) *lv = 1; 00487 if (*flg) rv->r = v; 00488 else { 00489 ierr = computennz(A); CHKERRQ(ierr); 00490 ierr = PetscObjectComposedDataGetReal 00491 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00492 if (*flg) rv->r = v; 00493 } 00494 00495 PetscFunctionReturn(0); 00496 } 00497 00498 #undef __FUNCT__ 00499 #define __FUNCT__ "SigmaDiagDist" 00500 static PetscErrorCode SigmaDiagDist 00501 (AnaModNumericalProblem prob,AnalysisItem *rv,int *lv,PetscTruth *flg) 00502 { 00503 Mat A = (Mat)prob; 00504 PetscReal v=0.; PetscTruth has; int id; PetscErrorCode ierr; 00505 PetscFunctionBegin; 00506 ierr = GetDataID("iprs","sigma-diag-dist",&id,&has); CHKERRQ(ierr); 00507 HASTOEXIST(has); 00508 ierr = PetscObjectComposedDataGetReal 00509 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00510 if (lv) *lv = 1; 00511 if (*flg) rv->r = v; 00512 else { 00513 ierr = computennz(A); CHKERRQ(ierr); 00514 ierr = PetscObjectComposedDataGetReal 00515 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00516 if (*flg) rv->r = v; 00517 } 00518 00519 PetscFunctionReturn(0); 00520 } 00521 00522 #undef __FUNCT__ 00523 #define __FUNCT__ "RegisterIprsModules" 00524 PetscErrorCode RegisterIprsModules(void) 00525 { 00526 PetscErrorCode ierr; 00527 PetscFunctionBegin; 00528 00529 ierr = RegisterModule 00530 ("iprs","nnzup",ANALYSISINTEGER,&NnzUp); CHKERRQ(ierr); 00531 ierr = RegisterModule 00532 ("iprs","nnzlow",ANALYSISINTEGER,&NnzLow); CHKERRQ(ierr); 00533 ierr = RegisterModule 00534 ("iprs","nnzdia",ANALYSISINTEGER,&NnzDia); CHKERRQ(ierr); 00535 00536 ierr = RegisterModule 00537 ("iprs","nnz",ANALYSISINTEGER,&Nnz); CHKERRQ(ierr); 00538 ierr = RegisterModule 00539 ("iprs","avgnnzprow",ANALYSISINTEGER,&AvgNnzpRow); CHKERRQ(ierr); 00540 ierr = RegisterModule 00541 ("iprs","avgdistfromdiag",ANALYSISINTEGER,&AvgNnzpRow); CHKERRQ(ierr); 00542 00543 ierr = RegisterModule 00544 ("iprs","relsymm",ANALYSISDOUBLE,&RelSymm); CHKERRQ(ierr); 00545 00546 ierr = RegisterModule 00547 ("iprs","upband",ANALYSISINTEGER,&UpBand); CHKERRQ(ierr); 00548 ierr = RegisterModule 00549 ("iprs","loband",ANALYSISINTEGER,&LoBand); CHKERRQ(ierr); 00550 00551 ierr = RegisterModule 00552 ("iprs","n-nonzero-diags",ANALYSISINTEGER,&NDiags); CHKERRQ(ierr); 00553 ierr = RegisterModule 00554 ("iprs","avg-diag-dist",ANALYSISDOUBLE,&AvgDiagDist); CHKERRQ(ierr); 00555 ierr = RegisterModule 00556 ("iprs","sigma-diag-dist",ANALYSISDOUBLE,&SigmaDiagDist); CHKERRQ(ierr); 00557 00558 PetscFunctionReturn(0); 00559 } 00560 00561 #undef __FUNCT__ 00562 #define __FUNCT__ "DeRegisterIprsModules" 00563 PetscErrorCode DeRegisterIprsModules(void) 00564 { 00565 PetscFunctionBegin; 00566 PetscFunctionReturn(0); 00567 }