#include "petscmat.h" PetscErrorCode MatZeroRowsColumnsStencil(Mat mat,PetscInt numRows,const MatStencil rows[],PetscScalar diag,Vec x,Vec b)Collective on Mat
mat | - the matrix | |
numRows | - the number of rows/columns to remove | |
rows | - the grid coordinates (and component number when dof > 1) for matrix rows | |
diag | - value put in all diagonals of eliminated rows (0.0 will even eliminate diagonal entry) | |
x | - optional vector of solutions for zeroed rows (other entries in vector are not used) | |
b | - optional vector of right hand side, that will be adjusted by provided solution |
If the option MatSetOption(mat,MAT_KEEP_NONZERO_PATTERN,PETSC_TRUE) the nonzero structure of the matrix is not changed (even for AIJ and BAIJ matrices) the values are merely zeroed.
The user can set a value in the diagonal entry (or for the AIJ and row formats can optionally remove the main diagonal entry from the nonzero structure as well, by passing 0.0 as the final argument).
For the parallel case, all processes that share the matrix (i.e., those in the communicator used for matrix creation) MUST call this routine, regardless of whether any rows being zeroed are owned by them.
Each processor can indicate any rows in the entire matrix to be zeroed (i.e. each process does NOT have to list only rows local to itself, but the row/column numbers are given in local numbering).
The grid coordinates are across the entire grid, not just the local portion
In Fortran idxm and idxn should be declared as
MatStencil idxm(4,m)and the values inserted using
idxm(MatStencil_i,1) = i
idxm(MatStencil_j,1) = j
idxm(MatStencil_k,1) = k
idxm(MatStencil_c,1) = cetc
For periodic boundary conditions use negative indices for values to the left (below 0; that are to be obtained by wrapping values from right edge). For values to the right of the last entry using that index plus one etc to obtain values that obtained by wrapping the values from the left edge. This does not work for anything but the DM_BOUNDARY_PERIODIC boundary type.
For indices that don't mean anything for your case (like the k index when working in 2d) or the c index when you have a single value per point) you can skip filling those indices.