1af0996ceSBarry Smith #include <petsc/private/matimpl.h> /*I "petscmat.h" I*/ 2c13ced23SMatthew G. Knepley 3ddc19b50SMatthew G. Knepley /*@ 4ddc19b50SMatthew G. Knepley MatComputeBandwidth - Calculate the full bandwidth of the matrix, meaning the width 2k+1 where k diagonals on either side are sufficient to contain all the matrix nonzeros. 5c13ced23SMatthew G. Knepley 6c13ced23SMatthew G. Knepley Collective on Mat 7c13ced23SMatthew G. Knepley 8c13ced23SMatthew G. Knepley Input Parameters: 9c13ced23SMatthew G. Knepley + A - The Mat 10c13ced23SMatthew G. Knepley - fraction - An optional percentage of the Frobenius norm of the matrix that the bandwidth should enclose 11c13ced23SMatthew G. Knepley 12c13ced23SMatthew G. Knepley Output Parameter: 13c13ced23SMatthew G. Knepley . bw - The matrix bandwidth 14c13ced23SMatthew G. Knepley 15c13ced23SMatthew G. Knepley Level: beginner 16c13ced23SMatthew G. Knepley 17c13ced23SMatthew G. Knepley .seealso: DMPlexCreate(), DMPlexSetConeSize(), DMPlexSetChart() 18c13ced23SMatthew G. Knepley @*/ 19ddc19b50SMatthew G. Knepley PetscErrorCode MatComputeBandwidth(Mat A, PetscReal fraction, PetscInt *bw) 20c13ced23SMatthew G. Knepley { 211ff9db0eSJed Brown PetscInt lbw[2] = {0, 0}, gbw[2]; 22c13ced23SMatthew G. Knepley PetscInt rStart, rEnd, r; 23c13ced23SMatthew G. Knepley 24c13ced23SMatthew G. Knepley PetscFunctionBegin; 25c13ced23SMatthew G. Knepley PetscValidHeaderSpecific(A, MAT_CLASSID, 1); 26c13ced23SMatthew G. Knepley PetscValidLogicalCollectiveReal(A,fraction,2); 27*dadcf809SJacob Faibussowitsch PetscValidIntPointer(bw, 3); 282c71b3e2SJacob Faibussowitsch PetscCheckFalse((fraction > 0.0) && (fraction < 1.0),PetscObjectComm((PetscObject) A), PETSC_ERR_SUP, "We do not yet support a fractional bandwidth"); 295f80ce2aSJacob Faibussowitsch CHKERRQ(MatGetOwnershipRange(A, &rStart, &rEnd)); 30c13ced23SMatthew G. Knepley for (r = rStart; r < rEnd; ++r) { 31c13ced23SMatthew G. Knepley const PetscInt *cols; 32c13ced23SMatthew G. Knepley PetscInt ncols; 33c13ced23SMatthew G. Knepley 345f80ce2aSJacob Faibussowitsch CHKERRQ(MatGetRow(A, r, &ncols, &cols, NULL)); 35c13ced23SMatthew G. Knepley if (ncols) { 36c13ced23SMatthew G. Knepley lbw[0] = PetscMax(lbw[0], r - cols[0]); 37c13ced23SMatthew G. Knepley lbw[1] = PetscMax(lbw[1], cols[ncols-1] - r); 38c13ced23SMatthew G. Knepley } 395f80ce2aSJacob Faibussowitsch CHKERRQ(MatRestoreRow(A, r, &ncols, &cols, NULL)); 40c13ced23SMatthew G. Knepley } 415f80ce2aSJacob Faibussowitsch CHKERRMPI(MPIU_Allreduce(lbw, gbw, 2, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject) A))); 42c13ced23SMatthew G. Knepley *bw = 2*PetscMax(gbw[0], gbw[1]) + 1; 43c13ced23SMatthew G. Knepley PetscFunctionReturn(0); 44c13ced23SMatthew G. Knepley } 45