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 <[email protected]> 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