#include <stdlib.h>
#include "anamod.h"
#include "anamatrix.h"
#include "anamodsalsamodules.h"
#include "petscerror.h"
#include "petscmat.h"
Go to the source code of this file.
Functions | |
static PetscErrorCode | computennz (Mat A) |
static PetscErrorCode | NnzUp (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | NnzLow (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | NnzDia (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | Nnz (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | RelSymm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | LoBand (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | UpBand (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | AvgNnzpRow (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | AvgDistFromDiag (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | NDiags (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | AvgDiagDist (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | SigmaDiagDist (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
PetscErrorCode | RegisterIprsModules (void) |
PetscErrorCode | DeRegisterIprsModules (void) |
PetscErrorCode RegisterIprsModules();
Compute these elements with
ComputeQuantity("iprs",<element>,A,(AnalysisItem*)&res,&reslen,&flg);
Available elements are:
@techreport{Zhang:interim, author = {Shu{T}ing Xu and Eun-Joo Lee and Jun Zhang}, title = {An Interim Analysis Report on Preconditioners and Matrices}, institution = {University of Kentucky, Lexington; Department of Computer Science}, number = {388-03}, year = {2003} }
Definition in file iprs.c.
static PetscErrorCode AvgDiagDist | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 477 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterIprsModules().
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 }
static PetscErrorCode AvgDistFromDiag | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 429 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
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 }
static PetscErrorCode AvgNnzpRow | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 405 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
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 }
static PetscErrorCode computennz | ( | Mat | A | ) | [static] |
This is the computational routine for all IPRS modules. Lots of pointer chasing.
Definition at line 61 of file iprs.c.
References GetDataID(), HasQuantity(), HASTOEXIST, AnalysisItem::i, id, and MatrixComputeQuantity().
Referenced by AvgDiagDist(), AvgDistFromDiag(), AvgNnzpRow(), LoBand(), NDiags(), Nnz(), NnzDia(), NnzLow(), NnzUp(), SigmaDiagDist(), and UpBand().
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 }
PetscErrorCode DeRegisterIprsModules | ( | void | ) |
static PetscErrorCode LoBand | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 331 of file iprs.c.
References computennz(), ComputeQuantity(), GetDataID(), HasComputeModule(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
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 }
static PetscErrorCode NDiags | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 453 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
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 }
static PetscErrorCode Nnz | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 273 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
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 }
static PetscErrorCode NnzDia | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 249 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
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 }
static PetscErrorCode NnzLow | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 225 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
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 }
static PetscErrorCode NnzUp | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute the number of nonzeroes in the upper triangle of the matrix
Definition at line 201 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
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 }
PetscErrorCode RegisterIprsModules | ( | void | ) |
Definition at line 524 of file iprs.c.
References ANALYSISDOUBLE, ANALYSISINTEGER, AvgDiagDist(), AvgNnzpRow(), LoBand(), NDiags(), Nnz(), NnzDia(), NnzLow(), NnzUp(), RegisterModule(), RelSymm(), SigmaDiagDist(), and UpBand().
Referenced by AnaModRegisterSalsaModules(), and AnaModRegisterStandardModules().
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 }
static PetscErrorCode RelSymm | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 297 of file iprs.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, AnalysisItem::i, id, and AnalysisItem::r.
Referenced by RegisterIprsModules().
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 }
static PetscErrorCode SigmaDiagDist | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 501 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterIprsModules().
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 }
static PetscErrorCode UpBand | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 368 of file iprs.c.
References computennz(), ComputeQuantity(), GetDataID(), HasComputeModule(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
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 }