[petsc-users] Use of MatRestrict/MatInterpolate with PCMG.

[Vijay S. Mahadevan]

Yes, I understand. But isn't there a flag to check if the transpose
was implicitly assumed ? It might be bad to compute the transpose of
the operator not provided, explicitly but the interface that calls
Restrict/Interpolate on Mat such as PCMG or DMMG should probably know
this and make the call appropriately. Or you could set the Mat
operations for the transposed operator accordingly i.e.,
MatMult_operator = MatMultTranspose; MatMultTranspose_operator =
MatMult. These are just couple of workarounds to tackle the glitch.

Barry, if you make further changes on this, I would much appreciate it
if you can let me know. Thanks.

Vijay

On Wed, Dec 15, 2010 at 9:16 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>
>  If you explicitly provide both then you are ok. If you provide one and not the other it will fail silently with a bad algorithm.
>
>
>  Barry
>
> On Dec 15, 2010, at 9:12 PM, Vijay S. Mahadevan wrote:
>
>>>  Well if they happen to be equal then it will never apply the transpose thus giving a bad algorithm and garbage.
>>
>> hmm, in the case when I provide a matrix for both
>> Interpolate/Restrict, I have (M==N==Nx==Ny),  this would always call
>> the MatMult routine. As long as I provide both these
>> operators/matrices explicitly, there is still no problem. A possible
>> issue is only when someone provides just the restriction or
>> prolongation. But I understand that when this happens, the other
>> operator is computed explicitly as its transpose. If this is actual
>> implementation, I still dont see a problem. Although, if the
>> restriction/prolongation operator are implicitly assumed to be
>> transpose of the other, then it will quite horribly fail since only
>> MatMult is called for both. I am not completely sure about the mode
>> currently used in petsc but it would be great if you can help me
>> understand.
>>
>> Note: I probe more on this since my linear tesselation (most often)
>> results in the same number of dofs on the coarser level
>> (p-coarsened/h-refined) and I dont want a glitch to come back and bite
>> me later on..
>>
>> Vijay
>>
>> On Wed, Dec 15, 2010 at 8:36 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>>>
>>> On Dec 15, 2010, at 8:16 PM, Vijay S. Mahadevan wrote:
>>>
>>>> Barry,
>>>>
>>>> Thanks for the prompt change ! I do not work on the development
>>>> version but I can update these matrix routines alone.
>>>>
>>>>>  Note it can still glitch if the restricted size is exactly the original size. :-(
>>>>
>>>> Why would it glitch if the restricted size is the same as the original
>>>> size though ? I dont see a case where your check (M==Ny) would fail.
>>>> Can you please elaborate more on this ?
>>>
>>>  Well if they happen to be equal then it will never apply the transpose thus giving a bad algorithm and garbage.
>>>
>>>  Barry
>>>
>>>>
>>>> Vijay
>>>>
>>>> On Wed, Dec 15, 2010 at 8:04 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>>>>>
>>>>>  I have pushed this change to petsc-dev and it is ready for use.
>>>>>
>>>>>   Barry
>>>>>
>>>>>  Note it can still glitch if the restricted size is exactly the original size. :-(
>>>>>
>>>>>
>>>>> On Dec 15, 2010, at 7:53 PM, Barry Smith wrote:
>>>>>
>>>>>>
>>>>>>  Vijay,
>>>>>>
>>>>>>    The use of M>N in MatRestrict and MatInterpolate was always a bit cheesy since it has this broken case that you reported. I will change it to do as you suggest and use the size of the vectors in determining which way to apply. But note I will do this in petsc-dev http://www.mcs.anl.gov/petsc/petsc-as/developers/index.html not petsc-3.1 so you'll need to switch if you are not using petsc-dev.
>>>>>>
>>>>>>   I'll try to get it down in the next few hours but it may take a little longer.
>>>>>>
>>>>>>
>>>>>>   Barry
>>>>>>
>>>>>> On Dec 15, 2010, at 6:06 PM, Vijay S. Mahadevan wrote:
>>>>>>
>>>>>>> Hi,
>>>>>>>
>>>>>>> I have an implementation issue with the MatRestrict/Interpolate
>>>>>>> functions. The problem is that one of my coarser levels (with PCMG)
>>>>>>> has higher dofs than the finest level. This does not always happen and
>>>>>>> requires a weird fine mesh system (in a sense) that uses multi-grid,
>>>>>>> but the idea is that the finest level problem has a high order (HO)
>>>>>>> discretization while the lower level mesh has a linear tesselation of
>>>>>>> the finest HO level (which I can optimize) and then adaptively
>>>>>>> coarsened levels beyond that. Since the number of columns in this case
>>>>>>> is larger than the number of rows, MatRestrict invariably calls
>>>>>>> MatMultTranspose to multiply instead of MatMult and vice-versa while
>>>>>>> calling  MatInterpolate. These result in assertion errors while
>>>>>>> comparing the length of Mat and Vec. The chosen method is based on
>>>>>>> whether (M>N) which seems to act against what I am doing here...
>>>>>>>
>>>>>>> I can always implement a shell matrix to replicate
>>>>>>> Restrict/Interpolate actions but my question is whether if such
>>>>>>> discretization will yield a consistent convergence in MG algorithm ?
>>>>>>> Is there a strong reason for checking if (M>N) rather than just doing
>>>>>>> (mat->rmap->N==y->map->N && mat->cmap->N==x->map->N) ? I would
>>>>>>> appreciate any detailed answer that you can provide for this and any
>>>>>>> suggestions to use the existing methods (without implementing the
>>>>>>> shell restriction) is very welcome.
>>>>>>>
>>>>>>> Thanks,
>>>>>>> vijay
>>>>>>
>>>>>
>>>>>
>>>
>>>
>
>