[petsc-dev] BlockGetIndices and GetBlockIndices

[Smith, Barry F.]

Fix one of presumably many bugs related to the rush to add PetscLayout too IS without understanding the ramifications.

https://gitlab.com/petsc/petsc/merge_requests/2191


> On Oct 21, 2019, at 8:34 AM, Pierre Jolivet via petsc-dev <petsc-dev at mcs.anl.gov> wrote:
> 
> 
>> On 21 Oct 2019, at 3:01 PM, Jed Brown <jed at jedbrown.org> wrote:
>> 
>> Pierre Jolivet via petsc-dev <petsc-dev at mcs.anl.gov> writes:
>> 
>>>> On 21 Oct 2019, at 7:52 AM, Smith, Barry F. <bsmith at mcs.anl.gov> wrote:
>>>> 
>>>> 
>>>> 
>>>>> On Oct 21, 2019, at 12:23 AM, Pierre Jolivet <pierre.jolivet at enseeiht.fr> wrote:
>>>>> 
>>>>> 
>>>>> 
>>>>>> On 21 Oct 2019, at 7:11 AM, Smith, Barry F. <bsmith at mcs.anl.gov> wrote:
>>>>>> 
>>>>>> 
>>>>>> 
>>>>>>> On Oct 20, 2019, at 11:52 PM, Pierre Jolivet <pierre.jolivet at enseeiht.fr> wrote:
>>>>>>> 
>>>>>>> 
>>>>>>> 
>>>>>>>> On 21 Oct 2019, at 6:42 AM, Smith, Barry F. <bsmith at mcs.anl.gov> wrote:
>>>>>>>> 
>>>>>>>> Could you provide a use case where you want to access/have a block size of a IS that is not an ISBlock? 
>>>>>>> 
>>>>>>> In the end, all I really want is get access to the underlying is->data->idx without having to worry about the subclass of is.
>>>>>>> I don’t have such a use case, but I don’t think this is really related to what I want to achieve (or maybe it is…).
>>>>>> 
>>>>>> ISGetIndices()
>>>>> 
>>>>> Not true for ISBlock with bs > 1.
>>>> 
>>>> Certainly suppose to be, is there a bug?
>>>> 
>>>> static PetscErrorCode ISGetIndices_Block(IS in,const PetscInt *idx[])
>>>> {
>>>> IS_Block       *sub = (IS_Block*)in->data;
>>>> PetscErrorCode ierr;
>>>> PetscInt       i,j,k,bs,n,*ii,*jj;
>>>> 
>>>> PetscFunctionBegin;
>>>> ierr = PetscLayoutGetBlockSize(in->map, &bs);CHKERRQ(ierr);
>>>> 
>>>> Dang, there is that stupid layout stuff again. Who put this crap in. 
>>>> 
>>>> ierr = PetscLayoutGetLocalSize(in->map, &n);CHKERRQ(ierr);
>>>> n   /= bs;
>>>> if (bs == 1) *idx = sub->idx;
>>>> else {
>>> 
>>> There it is, I don’t want this if/else. ISGetBlockIndices would have been a function always returning sub->idx.
>> 
>> Your code still can't skip the branch because ISGeneral can have bs>1.
> 
> Ah, OK, got confused by the lack of bs in the ISCreateGeneral constructor.
> 
> Thanks,
> Pierre
> 
>> Block size in IS means "these indices go together" while ISBlock imposes
>> the further constraint: "indices in each block are contiguous".  So you
>> can't just take the code you quoted where it returns the raw idx:
>> 
>> if (!isblock) {
>> ISGetIndices(is,&indices);
>> ISLocalToGlobalMappingCreate(comm,1,n,indices,PETSC_COPY_VALUES,mapping);
>> ISRestoreIndices(is,&indices);
>> } else {
>> ISGetBlockSize(is,&bs);
>> ISBlockGetIndices(is,&indices);
>> ISLocalToGlobalMappingCreate(comm,bs,n/bs,indices,PETSC_COPY_VALUES,mapping);
>> ISBlockRestoreIndices(is,&indices);
>> }
>> 
>> 
>> PetscErrorCode  ISGetBlockSize(IS is,PetscInt *size)
>> {
>> PetscErrorCode ierr;
>> 
>> PetscFunctionBegin;
>> ierr = PetscLayoutGetBlockSize(is->map, size);CHKERRQ(ierr);
>> PetscFunctionReturn(0);
>> }
>> 
>