SALSA Analysis Modules
|
Quantities used in the UKY `Intelligent Preconditioner Recommendation System'. More...
#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, PetscBool *flg) |
static PetscErrorCode | NnzLow (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | NnzDia (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | Nnz (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | RelSymm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | LoBand (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | UpBand (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | AvgNnzpRow (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | NDiags (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | AvgDiagDist (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
static PetscErrorCode | SigmaDiagDist (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscBool *flg) |
PetscErrorCode | RegisterIprsModules (void) |
PetscErrorCode | DeRegisterIprsModules (void) |
Quantities used in the UKY `Intelligent Preconditioner Recommendation System'.
The University of Kentucky `Intelligent Preconditioner Recommendation System' (Xu[2003]) uses a number of structural properties of the matrix that go beyond the properties in our own Structural properties module.
Activate this module with
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, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 488 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; PetscReal v=0.; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","avg-diag-dist",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->r = v; else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode AvgNnzpRow | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 414 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; int v=0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","avgnnzprow",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->i = v; else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } PetscFunctionReturn(0); }
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(), AvgNnzpRow(), LoBand(), NDiags(), Nnz(), NnzDia(), NnzLow(), NnzUp(), SigmaDiagDist(), and UpBand().
{ MPI_Comm comm; int first,last,M,row,id,max_rowlen,*cols_tmp, nnzl,gnnzl,nnzu,gnnzu,nnzd,gnnzd, nnz,gnnz,avg, dist,gdist, p,gp, q,gq; IS colis,colsis=NULL; PetscBool has; PetscErrorCode ierr; PetscFunctionBegin; ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr); ierr = MatGetSize(A,&M,PETSC_NULL); CHKERRQ(ierr); ierr = MatGetOwnershipRange(A,&first,&last); CHKERRQ(ierr); ierr = HasQuantity (A,"structure","max-nnzeros-per-row",&has); CHKERRQ(ierr); if (has) { AnalysisItem r; ierr = MatrixComputeQuantity (A,"structure","max-nnzeros-per-row",&r,PETSC_NULL,PETSC_NULL); CHKERRQ(ierr); max_rowlen = r.i; } else { ierr = MatGetSize(A,PETSC_NULL,&max_rowlen); CHKERRQ(ierr); } ierr = PetscMalloc(max_rowlen*sizeof(int),&cols_tmp); CHKERRQ(ierr); dist = nnz = nnzl = nnzd = nnzu = p = q = 0; for (row=first; row<last; row++) { int ncols,icol; const int *cols; ierr = MatGetRow(A,row,&ncols,&cols,PETSC_NULL); CHKERRQ(ierr); if (ncols>max_rowlen) SETERRQ(MPI_COMM_SELF,1,"Max rowlen exceeded"); ierr = PetscMemcpy(cols_tmp,cols,ncols*sizeof(int)); CHKERRQ(ierr); nnz += ncols; /* update left and right bandwidth */ if (cols[0]<row && row-cols[0]>p) p = row-cols[0]; if (cols[ncols-1]>row && cols[ncols-1]-row>q) q = cols[ncols-1]-row; /* compute summed distance from diagonal; count elements in D,L,U */ for (icol=0; icol<ncols; icol++) { if (cols[icol]<row) { dist += row-cols[icol]; nnzl++; } else if (cols[icol]>row) { dist += cols[icol]-row; nnzu++; } else nnzd++; } /* analyze diagonals */ for (icol=0; icol<ncols; icol++) cols_tmp[icol]-=row; CHKMEMQ; if (!colsis) { ierr = ISCreateGeneral (PETSC_COMM_SELF,ncols,cols_tmp,PETSC_COPY_VALUES, &colsis); CHKERRQ(ierr); } else { IS tmpis; ierr = ISCreateGeneral (PETSC_COMM_SELF,ncols,cols_tmp,PETSC_COPY_VALUES, &colis); CHKERRQ(ierr); CHKMEMQ; ierr = ISSort(colis); CHKERRQ(ierr); CHKMEMQ; ierr = ISSum(colsis,colis,&tmpis); CHKERRQ(ierr); CHKMEMQ; ierr = ISDestroy(&colsis); CHKERRQ(ierr); CHKMEMQ; ierr = ISDestroy(&colis); CHKERRQ(ierr); colsis = tmpis; } for (icol=0; icol<ncols; icol++) cols_tmp[icol]+=row; ierr = MatRestoreRow(A,row,&ncols,&cols,PETSC_NULL); CHKERRQ(ierr); } ierr = PetscFree(cols_tmp); CHKERRQ(ierr); ierr = MPI_Allreduce((void*)&nnzl,(void*)&gnnzl,1,MPI_INT,MPI_SUM,comm); ierr = MPI_Allreduce((void*)&nnzu,(void*)&gnnzu,1,MPI_INT,MPI_SUM,comm); ierr = MPI_Allreduce((void*)&nnzd,(void*)&gnnzd,1,MPI_INT,MPI_SUM,comm); ierr = MPI_Allreduce((void*)&nnz,(void*)&gnnz,1,MPI_INT,MPI_SUM,comm); if (gnnz!=(gnnzl+gnnzd+gnnzu)) SETERRQ1(MPI_COMM_WORLD,1,"We missed %d elements",gnnz-(gnnzl+gnnzd+gnnzu)); ierr = MPI_Allreduce((void*)&dist,(void*)&gdist,1,MPI_INT,MPI_SUM,comm); ierr = MPI_Allreduce((void*)&p,(void*)&gp,1,MPI_INT,MPI_MAX,comm); ierr = MPI_Allreduce((void*)&q,(void*)&gq,1,MPI_INT,MPI_MAX,comm); ierr = GetDataID("iprs","nnzlow",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,gnnzl); CHKERRQ(ierr); ierr = GetDataID("iprs","nnzup",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,gnnzu); CHKERRQ(ierr); ierr = GetDataID("iprs","nnzdia",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,gnnzd); CHKERRQ(ierr); ierr = GetDataID("iprs","nnz",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,gnnz); CHKERRQ(ierr); avg = nnz/M; ierr = GetDataID("iprs","avgnnzprow",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,avg); CHKERRQ(ierr); dist = gdist/M; ierr = GetDataID("iprs","avgdistfromdiag",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,dist); CHKERRQ(ierr); ierr = GetDataID("iprs","loband",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt((PetscObject)A,id,gp); CHKERRQ(ierr); ierr = GetDataID("iprs","upband",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt((PetscObject)A,id,gq); CHKERRQ(ierr); { int nd,i; const PetscInt *d; PetscReal avgdist,sigmadist; ierr = ISGetLocalSize(colsis,&nd); CHKERRQ(ierr); ierr = GetDataID("iprs","n-nonzero-diags",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetInt((PetscObject)A,id,nd); CHKERRQ(ierr); ierr = ISGetIndices(colsis,&d); CHKERRQ(ierr); avgdist = 0.; for (i=0; i<nd; i++) avgdist += PetscAbsScalar(d[i]); avgdist /= nd; ierr = GetDataID("iprs","avg-diag-dist",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,avgdist); CHKERRQ(ierr); sigmadist = 0.; for (i=0; i<nd; i++) sigmadist += (d[i]-avgdist)*(d[i]-avgdist); sigmadist /= nd; ierr = GetDataID("iprs","sigma-diag-dist",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,sigmadist); CHKERRQ(ierr); ierr = ISRestoreIndices(colsis,&d); CHKERRQ(ierr); ierr = ISDestroy(&colsis); CHKERRQ(ierr); } PetscFunctionReturn(0); }
PetscErrorCode DeRegisterIprsModules | ( | void | ) |
static PetscErrorCode LoBand | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 340 of file iprs.c.
References computennz(), ComputeQuantity(), GetDataID(), HasComputeModule(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; int v=0; PetscBool has; int id,ql; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","loband",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->i = v; else { AnalysisItem q; ierr = HasComputeModule("structure","left-bandwidth",flg); CHKERRQ(ierr); if (*flg) { ierr = ComputeQuantity (prob,"structure","left-bandwidth",&q,&ql,flg); CHKERRQ(ierr); if (*flg) { v = q.i; ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,v); CHKERRQ(ierr); rv->i = v; } } else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } } PetscFunctionReturn(0); }
static PetscErrorCode NDiags | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 464 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; int v=0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","n-nonzero-diags",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->i = v; else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } PetscFunctionReturn(0); }
static PetscErrorCode Nnz | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 282 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; int v=0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","nnz",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->i = v; else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } PetscFunctionReturn(0); }
static PetscErrorCode NnzDia | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 258 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; int v=0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","nnzdia",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->i = v; else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } PetscFunctionReturn(0); }
static PetscErrorCode NnzLow | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 234 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; int v=0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","nnzlow",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->i = v; else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } PetscFunctionReturn(0); }
static PetscErrorCode NnzUp | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Compute the number of nonzeroes in the upper triangle of the matrix
Definition at line 210 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; int v=0; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","nnzup",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->i = v; else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } PetscFunctionReturn(0); }
PetscErrorCode RegisterIprsModules | ( | void | ) |
Definition at line 535 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().
{ PetscErrorCode ierr; PetscFunctionBegin; ierr = RegisterModule ("iprs","nnzup",ANALYSISINTEGER,&NnzUp); CHKERRQ(ierr); ierr = RegisterModule ("iprs","nnzlow",ANALYSISINTEGER,&NnzLow); CHKERRQ(ierr); ierr = RegisterModule ("iprs","nnzdia",ANALYSISINTEGER,&NnzDia); CHKERRQ(ierr); ierr = RegisterModule ("iprs","nnz",ANALYSISINTEGER,&Nnz); CHKERRQ(ierr); ierr = RegisterModule ("iprs","avgnnzprow",ANALYSISINTEGER,&AvgNnzpRow); CHKERRQ(ierr); ierr = RegisterModule ("iprs","avgdistfromdiag",ANALYSISINTEGER,&AvgNnzpRow); CHKERRQ(ierr); ierr = RegisterModule ("iprs","relsymm",ANALYSISDOUBLE,&RelSymm); CHKERRQ(ierr); ierr = RegisterModule ("iprs","upband",ANALYSISINTEGER,&UpBand); CHKERRQ(ierr); ierr = RegisterModule ("iprs","loband",ANALYSISINTEGER,&LoBand); CHKERRQ(ierr); ierr = RegisterModule ("iprs","n-nonzero-diags",ANALYSISINTEGER,&NDiags); CHKERRQ(ierr); ierr = RegisterModule ("iprs","avg-diag-dist",ANALYSISDOUBLE,&AvgDiagDist); CHKERRQ(ierr); ierr = RegisterModule ("iprs","sigma-diag-dist",ANALYSISDOUBLE,&SigmaDiagDist); CHKERRQ(ierr); PetscFunctionReturn(0); }
static PetscErrorCode RelSymm | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 306 of file iprs.c.
References ComputeQuantity(), GetDataID(), HASTOEXIST, AnalysisItem::i, id, and AnalysisItem::r.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; PetscReal v=0.; PetscBool has; int id,ql; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","relsymm",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->r = v; else { AnalysisItem q; int nnz,nnu; ierr = ComputeQuantity (prob,"structure","nnzeros",&q,&ql,flg); CHKERRQ(ierr); if (*flg) { nnz = q.i; ierr = ComputeQuantity (prob,"structure","n-struct-unsymm",&q,&ql,flg); CHKERRQ(ierr); if (*flg) { nnu = q.i; v = (nnz-nnu)/(1.*nnz); ierr = PetscObjectComposedDataSetReal ((PetscObject)A,id,v); CHKERRQ(ierr); rv->r = v; } } } PetscFunctionReturn(0); }
static PetscErrorCode SigmaDiagDist | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 512 of file iprs.c.
References computennz(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; PetscReal v=0.; PetscBool has; int id; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","sigma-diag-dist",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->r = v; else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetReal ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->r = v; } PetscFunctionReturn(0); }
static PetscErrorCode UpBand | ( | AnaModNumericalProblem | prob, |
AnalysisItem * | rv, | ||
int * | lv, | ||
PetscBool * | flg | ||
) | [static] |
Definition at line 377 of file iprs.c.
References computennz(), ComputeQuantity(), GetDataID(), HasComputeModule(), HASTOEXIST, AnalysisItem::i, and id.
Referenced by RegisterIprsModules().
{ Mat A = (Mat)prob; int v=0; PetscBool has; int id,ql; PetscErrorCode ierr; PetscFunctionBegin; ierr = GetDataID("iprs","upband",&id,&has); CHKERRQ(ierr); HASTOEXIST(has); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (lv) *lv = 1; if (*flg) rv->i = v; else { ierr = HasComputeModule("structure","right-bandwidth",flg); CHKERRQ(ierr); if (*flg) { AnalysisItem q; ierr = ComputeQuantity (prob,"structure","right-bandwidth",&q,&ql,flg); CHKERRQ(ierr); if (*flg) { v = q.i; ierr = PetscObjectComposedDataSetInt ((PetscObject)A,id,v); CHKERRQ(ierr); rv->i = v; } } else { ierr = computennz(A); CHKERRQ(ierr); ierr = PetscObjectComposedDataGetInt ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); if (*flg) rv->i = v; } } PetscFunctionReturn(0); }