Re: [petsc-users] Dropping single entries in matrix insertion for block matrices

[Smith, Barry F. via petsc-users]

> On Dec 10, 2018, at 5:57 AM, Lawrence Mitchell <we...@gmx.li> wrote:
> 
> Hi Barry,
> 
>> On 9 Dec 2018, at 19:25, Smith, Barry F. <bsm...@mcs.anl.gov> wrote:
>> 
>> 
>>  Lawrence,
>> 
>>     I understand what you want and it is a reasonable request. The problem 
>> is that currently ISLocalToGlobalMappingCreate() when used with block 
>> vectors and matrices is always based on blocks, that is, from the manual 
>> page, "There is one integer value in indices per block and it represents the 
>> actual indices bs*idx + j, where j=0,..,bs-1". If you look at 
>> ISLocalToGlobalMappingApply() it uses the code:
>> 
>> out[i] = idx[in[i]/bs]*bs + (in[i] % bs);
>> 
>> to do the mapping; that is the idx[] mapping is stored only by block, not by 
>> point.
>> 
>> But I think you may be able to get the effect you want by "managing the 
>> local to global mapping yourself before calling MatSetValues()". That is, 
>> you create a ISLocalToGlobal object yourself, not based on blocks, then when 
>> setting your Dirichlet conditions you call ISLocalToGlobalMapping() apply 
>> yourself and then call MatSetValues() using the resulting indices.
> 
> Thanks. I can do this. What would be neat is if I could just swap out the 
> LGMaps. I can try and prepare a PR that does this, if you think it would be 
> useful.
> 
> I think the checking should be something like:
> 
> - Does ISLocalToGlobalMappingGetSize(...) match the PetscLayout size?
> 
> - If using the blocked interface, does the block size of the LGMap match the 
> block size of the PetscLayout?
> 
> This would move some error checking from MatSetLocalToGlobalMapping to 
> MatSetValuesXXXLocal, which may not be desired.
> 
> What do you think?

   Well the thing is the indices you pass in are currently always "blocked 
base", they work automatically for both blocked and non blocked versions of 
MatSetValuesLocal(). To support what you want you need to be able to pass in 
nonblock based version of the indices and (as you noted) keep track of whether 
you pass in the blocked or non blocked indices so you can apply the scaling or 
not.

   Ok, give it a try, it does change the interface a bit but I suppose is worth 
it since it introduces more functionality.

    Barry



> 
> Lawrence