[petsc-dev] MetSetBlockSize issue

Wed Apr 25 16:09:06 CDT 2012

On Apr 25, 2012, at 3:48 PM, Barry Smith wrote:

> 
> 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 

No.  The blocks are dense for smoothed aggregation but for geometric MG and classical AMG the blocks are scaled identity.

> 
> 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() 
> 

All I really need is for the coarse grid matrix to have the block size, because that is used to create a graph (for "vertices" not equations).  The method that constructs the prolongator knows the block size, so I attach it to P.  I could actually attach the block size to the coarse matrix in the way that I do with P, or pass it around as a parameter internally, I suppose.

Mark

>    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
>>>>>>>>>>>>> 
>>>>>>>>>>>>> 
>>>>>>>>>>>> 
>>>>>>>>>>> 
>>>>>>>>>>> 
>>>>>>>>>> 
>>>>>>>>> 
>>>>>>>>> 
>>>>>>>> 
>>>>>>>> 
>>>>>>> 
>>>>>> 
>>>>>> 
>>>>> 
>>>> 
>>>> 
>>> 
>>> 
>> 
> 
> 