[petsc-users] Zero entries of a sparse matrix in a block-wise style

[Mingchang Ding]

> On Jul 26, 2019, at 3:32 PM, Mingchang Ding <dmcvamos at udel.edu> wrote:
> 
> 
>> On Jul 26, 2019, at 3:21 PM, Jed Brown <jed at jedbrown.org> wrote:
>> 
>> Please always use "reply-all" so that your messages go to the list.
>> This is standard mailing list etiquette.  It is important to preserve
>> threading for people who find this discussion later and so that we do
>> not waste our time re-answering the same questions that have already
>> been answered in private side-conversations.  You'll likely get an
>> answer faster that way too.
> 
>  Sorry for this.
> 
>> Mingchang Ding <dmcvamos at udel.edu> writes:
>> 
>>> Thank you for your quick reply. The code I have now is attached. The sparse matrix obtained from LDG discretization to u_xx is assembled in SetupMat routine of setup.h. When I choose my polynomial degree nk =1 (DOF = 2), only the first and last row of the output sparse matrix is as what I want. Again, I want to impose essential (Dirichlet) boundary condition to the problem.
>>> 
>>> Let me know if any explanation of the question is needed.
>> It appears that you are calling MatZeroRows with rows = {0, info.mx-1},
>> which will only affect the first and last rows.  The other entries are
>> how you have assembled them (with the matrix-matrix product).
> 
> Yes. I think so. Is there anyway I can fill in the rows array with all the first and last DOF rows index efficiently instead of only filling in one row index at one time?

  I think I misunderstood just now and know how to fix the problem. 

  Thank you!

>>     PetscInt        nrows = 0;
>>     PetscInt        *rows;
>> 
>>     ierr = DMDAGetLocalInfo(da, &info);CHKERRQ(ierr);
>> 
>>     PetscMalloc1(info.mx, &rows);
>> 
>>     for (j=info.xs; j<info.xs+info.xm; j++) {
>> 
>>         if (j == 0) {
>>             rows[nrows++] = j;
>> 
>> [...]
>> 
>>         }
>>         else if (j == info.mx - 1) {
>>             rows[nrows++] = j;
>> 
>> [...]
>> 
>> 
>>         }
>>         else {
>> 
>> [...]
>> 
>> 
>>         }
>> 
>>         /* Mass_inv */
>>         ierr = MatSetValuesBlocked(appctx->Mass_inv, 1, &j, 1, &j, &inv[0][0], INSERT_VALUES);CHKERRQ(ierr);
>> 
>>     }
>> 
>> 
>> [...]
>> 
>> 
>>     MatZeroRows(appctx->J_dif, nrows, rows, 1.0, NULL, NULL);
>