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 PetscErrorCode ierr; 24c13ced23SMatthew G. Knepley 25c13ced23SMatthew G. Knepley PetscFunctionBegin; 26c13ced23SMatthew G. Knepley PetscValidHeaderSpecific(A, MAT_CLASSID, 1); 27c13ced23SMatthew G. Knepley PetscValidLogicalCollectiveReal(A,fraction,2); 28c13ced23SMatthew G. Knepley PetscValidPointer(bw, 3); 29*2c71b3e2SJacob Faibussowitsch PetscCheckFalse((fraction > 0.0) && (fraction < 1.0),PetscObjectComm((PetscObject) A), PETSC_ERR_SUP, "We do not yet support a fractional bandwidth"); 30c13ced23SMatthew G. Knepley ierr = MatGetOwnershipRange(A, &rStart, &rEnd);CHKERRQ(ierr); 31c13ced23SMatthew G. Knepley for (r = rStart; r < rEnd; ++r) { 32c13ced23SMatthew G. Knepley const PetscInt *cols; 33c13ced23SMatthew G. Knepley PetscInt ncols; 34c13ced23SMatthew G. Knepley 35c13ced23SMatthew G. Knepley ierr = MatGetRow(A, r, &ncols, &cols, NULL);CHKERRQ(ierr); 36c13ced23SMatthew G. Knepley if (ncols) { 37c13ced23SMatthew G. Knepley lbw[0] = PetscMax(lbw[0], r - cols[0]); 38c13ced23SMatthew G. Knepley lbw[1] = PetscMax(lbw[1], cols[ncols-1] - r); 39c13ced23SMatthew G. Knepley } 40c13ced23SMatthew G. Knepley ierr = MatRestoreRow(A, r, &ncols, &cols, NULL);CHKERRQ(ierr); 41c13ced23SMatthew G. Knepley } 42820f2d46SBarry Smith ierr = MPIU_Allreduce(lbw, gbw, 2, MPIU_INT, MPI_MAX, PetscObjectComm((PetscObject) A));CHKERRMPI(ierr); 43c13ced23SMatthew G. Knepley *bw = 2*PetscMax(gbw[0], gbw[1]) + 1; 44c13ced23SMatthew G. Knepley PetscFunctionReturn(0); 45c13ced23SMatthew G. Knepley } 46