Re: [Libmesh-devel] Multigrid Projection Matrixes with Periodic Boundaries

[Roy S. <roy...@ic...>]

On Thu, 6 Jul 2017, Boris Boutkov wrote:

> Hello again - Ive been doing further work on the multigrid
> projection matrix front and things are starting to look fairly
> stable overall!

Great!

> The main issue that has resurfaced is the treatment of the boundary
> conditions. Though it seems that for most cases convergence rates
> are as expected, there clearly still remains a problem with the
> periodic cases in even the simplest 1D/2D test problems.

Strange.  I suppose I can imagine 2D giving us trouble, but in 1D,
every restriction or prolongation of a periodic solution will still be
a periodic solution.

On the other hand...  In the basic multigrid methods, we do the
prolongation on a coarse error term, and the restriction on a fine
residual.

The error term should be periodic, so we shouldn't need to care if the
prolongation doesn't redundantly enforce that.

The fine residual, though... it's actually hard for me to picture in
my head what prolongation really means in the PDE case, because
although the residual coefficients have the same algebraic structure
as the primal solution coefficients, the residual is in the dual
space.  I don't think "pretend they're coefficients to a function and
project that" is the correct thing to do, but I worry that maybe it's
correct enough to work in simple cases and incorrect enough to fail as
we add complications like periodicity and AMR.  "Find the function
which gives the same L2 inner product as the residual vector's l2
inner product, then project that function, then find the residual
vector in the new space which gives the same l2 inner product as the
projected function's L2 inner product", like a discrete Riesz
representation transformation, maybe?

Anyway, that's an aside.  However the residual vector is being
restricted, I take it it's not already conforming to the periodic
constraints, even in 1D, and so we need to modify the matrix to make
sure it's conforming after the restriction?

> I think that the natural way forward is that we need to constrain
> the element matrixes as we are building up the projection matrixes,

Well, we may need an option to do so, at least.  It definitely needs
to be a functor that we can make an inline no-op in the common cases,
lest that refactored projection be slowed down again.

> the first place that comes to mind would be right after we have set
> the target_matrix dof information inside the MatrixFillAction
> component of the GenericProjector.

I'm not sure if this is right, but it's at least probably worth trying.

> The issue then is that the usual DofMap::constrain_element_matrix
> implementations only have DenseMatrix versions whereas here it seems
> we would need Sparse ones in the projection matrix situation.

True.  Can we get away with the same "move the implementations to a
templated helper function, then instantiate the common cases in the
original .C file for backwards compatibility, then instantiate the
DSNA cases in the new projection matrix code" trick?
---
Roy