[petsc-dev] MetSetBlockSize issue

[Barry Smith]

On Apr 25, 2012, at 2:44 PM, Mark F. Adams wrote:

> 
> On Apr 25, 2012, at 12:22 PM, Hong Zhang wrote:
> 
>> Mark,
>> 
>> I attach it with this code:
>> 
>>  /* attach block size of columns */
>>  if( pc_gamg->col_bs_id == -1 ) {
>>    ierr = PetscObjectComposedDataRegister( &pc_gamg->col_bs_id ); assert(pc_gamg->col_bs_id != -1 );
>>  }
>>  ierr = PetscObjectComposedDataSetInt( (PetscObject)Prol, pc_gamg->col_bs_id, data_cols ); CHKERRQ(ierr);
>> 
>> which Hong can not see w/o my ID.  so we need a mechanism for me to get this to PtAP.
>>  
>> Sorry, I do not understand what do you mean here.
>> Do you want add above into the implementation of PtAP or petsc? 
> 
> This is an architectural/API decision.  I will do what you and Barry want.  You could add a parameter to PtAP, or add a field to Mat (column block size), or whatever.

   Mark,

     I'm leaning towards setting a column block size, MatSetColumnBlockSize(). But what exactly is a "column block size"?  Is it really a column blocksize in the sense that nbs columns next to each other always have the exact same nonzero pattern: for example

xx  00   xx 00 
yy  00   yy 00 
zz  00  zz  00 

or is it a more vague thing that defined by "well if you compute Pt A P the result has a blocksize of nbs".  In that case we would provide a MatSetVirtualColumnBlockSize() 

    Thanks

   Barry

> 
> Mark
> 
>> Hong
>> 
>> And MatGetSubMatrix has to inherit block size, which I assume is trial.
>> 
>> Mark
>> 
>> On Apr 24, 2012, at 11:14 PM, Barry Smith wrote:
>> 
>> >
>> > On Apr 24, 2012, at 8:14 PM, Mark F. Adams wrote:
>> >
>> >>
>> >> On Apr 24, 2012, at 5:01 PM, Barry Smith wrote:
>> >>
>> >>>
>> >>> On Apr 24, 2012, at 2:38 PM, Mark F. Adams wrote:
>> >>>
>> >>>>
>> >>>> On Apr 24, 2012, at 3:28 PM, Hong Zhang wrote:
>> >>>>
>> >>>>> Mark :
>> >>>>> Shall C=PtAP inherit the block size of A?
>> >>>>> Currently, MatPtAP() is only implemented with aij bs=1.
>> >>>>
>> >>>> No, as I said, algebraically it should inherit the column block size of P, but that is not available.
>> >>>
>> >>> What is the column block size of P and how do you know what it is?
>> >>
>> >> I now attach it to the P matrix instead of adding it as a return parameter of the create P method.  So I keep track of it in the code.  Its size is the number of null space vectors.
>> >
>> >  Ok if P knows its block size then the code that computes PtAP can set the block size of the resulting matrix as it creates it using the information in P. So is the issue on calling MatSetBlockSize() now resolved and everyone happy?
>> >
>> >   Barry
>> >
>> >>
>> >> Mark
>> >>
>> >>>
>> >>> Barry
>> >>>
>> >>>>
>> >>>> Perhaps PtAP should take a parameter for the block size ...
>> >>>>
>> >>>> Mark
>> >>>>
>> >>>>> Hong
>> >>>>>
>> >>>>> On Apr 24, 2012, at 2:55 PM, Barry Smith wrote:
>> >>>>>
>> >>>>>>
>> >>>>>> On Apr 24, 2012, at 1:39 PM, Mark F. Adams wrote:
>> >>>>>>
>> >>>>>>> Now I'm getting this error with code like this:
>> >>>>>>>
>> >>>>>>>  ierr = MatGetSubMatrix( Cmat, new_eq_indices, new_eq_indices, MAT_INITIAL_MATRIX, &mat );
>> >>>>>>>  CHKERRQ(ierr);
>> >>>>>>>  ierr = MatSetBlockSize( mat, cbs );      CHKERRQ(ierr);
>> >>>>>>>
>> >>>>>>> and like this:
>> >>>>>>>
>> >>>>>>> ierr = MatPtAP( Amat_fine, Pold, MAT_INITIAL_MATRIX, 2.0, &Cmat ); CHKERRQ(ierr);
>> >>>>>>> ierr = MatSetBlockSize( Cmat, cbs );      CHKERRQ(ierr);
>> >>>>>>
>> >>>>>> Right you cannot do this. The matrix already exists so you cannot now set its size.
>> >>>>>>
>> >>>>>> Does Cmat have a block size of cbs
>> >>>>>>
>> >>>>>> What about Amat?   By default these routines should create the new matrix with the correct blocksize.
>> >>>>>>
>> >>>>>> The harder part is if the blocksize of mat would be different than Cmat?
>> >>>>>>
>> >>>>>
>> >>>>> It seems like MatGetSubMatrix should just inherit the block size but PtAP is harder.  The block size of the Cmat is the column block size of P.  But I can not set a column block size (!= row block size) so I don't set block size on P at all.  Here cbs is this column block size of P, or what it should be.
>> >>>>>
>> >>>>> Mark
>> >>>>>
>> >>>>>> Barry
>> >>>>>>
>> >>>>>>>
>> >>>>>>> Mark
>> >>>>>>>
>> >>>>>>> On Apr 23, 2012, at 9:09 PM, Barry Smith wrote:
>> >>>>>>>
>> >>>>>>>>
>> >>>>>>>> Yes, look, for example how ex32 is handled at the bottom of the makefile.
>> >>>>>>>>
>> >>>>>>>> Thanks
>> >>>>>>>>
>> >>>>>>>> Barry
>> >>>>>>>>
>> >>>>>>>> On Apr 23, 2012, at 5:50 PM, Mark F. Adams wrote:
>> >>>>>>>>
>> >>>>>>>>>
>> >>>>>>>>> On Apr 23, 2012, at 5:59 PM, Barry Smith wrote:
>> >>>>>>>>>
>> >>>>>>>>>>
>> >>>>>>>>>> It was changed a while ago that MatSetBlockSize() couldn't be set after the matrix was full instantiated. I guess that example did not get fixed because it is not listed in the makefile to run in the makeall.
>> >>>>>>>>>
>> >>>>>>>>> May I add it?
>> >>>>>>>>>
>> >>>>>>>>>>
>> >>>>>>>>>> I have fixed the example to call MatSetBlockSize() at a safe point.
>> >>>>>>>>>>
>> >>>>>>>>>> Barry
>> >>>>>>>>>>
>> >>>>>>>>>> On Apr 23, 2012, at 4:29 PM, Mark F. Adams wrote:
>> >>>>>>>>>>
>> >>>>>>>>>>> ex55.c in ksp is failing with:
>> >>>>>>>>>>>
>> >>>>>>>>>>> [0]PETSC ERROR: --------------------- Error Message ------------------------------------
>> >>>>>>>>>>> [0]PETSC ERROR: Arguments are incompatible!
>> >>>>>>>>>>> [0]PETSC ERROR: Cannot change block size 1 to 2!
>> >>>>>>>>>>> [0]PETSC ERROR: ------------------------------------------------------------------------
>> >>>>>>>>>>>
>> >>>>>>>>>>> on this line 57 of ex55.c:
>> >>>>>>>>>>>
>> >>>>>>>>>>> ierr = MatSetBlockSize(Amat,2);      CHKERRQ(ierr);
>> >>>>>>>>>>>
>> >>>>>>>>>>> Any idea what happened here?
>> >>>>>>>>>>>
>> >>>>>>>>>>> Mark
>> >>>>>>>>>>
>> >>>>>>>>>>
>> >>>>>>>>>
>> >>>>>>>>
>> >>>>>>>>
>> >>>>>>>
>> >>>>>>
>> >>>>>>
>> >>>>>
>> >>>>>
>> >>>>
>> >>>
>> >>>
>> >>
>> >
>> >
>> 
>> 
>