[petsc-users] Neumann BC with non-symmetric matrix

[Mohammad Mirzadeh]

Barry,
On Wednesday, February 24, 2016, Barry Smith <bsmith at mcs.anl.gov> wrote:

>
> > On Feb 23, 2016, at 11:35 PM, Mohammad Mirzadeh <mirzadeh at gmail.com
> <javascript:;>> wrote:
> >
> > Dear all,
> >
> > I am dealing with a situation I was hoping to get some suggestions here.
> Suppose after discretizing a poisson equation with purely neumann (or
> periodic) bc I end up with a matrix that is *almost* symmetric, i.e. it is
> symmetric for almost all grid points with the exception of a few points.
>
>   How come it is not purely symmetric? The usual finite elements with pure
> Neumann or periodic bc will give a completely symmetric matrix.
>
>   Barry
>
>
So this is a finite difference discretization on adaptive Cartesian grids.
It turns out that the discretization is non-symmetric at the corse-fine
interface. It's actually not because of the BC itself.

> >
> > The correct way of handling this problem is by specifying the nullspace
> to MatSetNullSpace. However, since the matrix is non-symmetric in general I
> would need to pass the nullspace of A^T. Now it turns out that if A is
> *sufficiently close to being symmetric*, I can get away with the constant
> vector, which is the nullspace of A and not A^T, but obviously this does
> not always work. Sometimes the KSP converges and in other situations the
> residual stagnates which is to be expected.
> >
> > Now, here are my questions (sorry if they are too many!):
> >
> > 1) Is there any efficient way of calculating nullspace of A^T in this
> case? Is SVD the only way?
> >
> > 2) I have tried fixing the solution at an arbitrary point, and while it
> generally works, for some problems I get numerical artifacts, e.g. slight
> asymmetry in the solution and/or increased error close to the point where I
> fix the solution. Is this, more or less, expected as a known artifact?
> >
> > 3) An alternative to 2 is to enforce some global constraint on the
> solution, e.g. to require that the average be zero. My question here is
> two-fold:
>
>   Requiring the average be zero is exactly the same as providing a null
> space of the constant function. Saying the average is zero is the same as
> saying the solution is orthogonal to the constant function. I don't see any
> reason to introduce the Lagrange multiplier and all its complications
> inside of just providing the constant null space.


Is this also true at the discrete level when the matrix is non-symmetric? I
have always viewed this as just a constraint that could really be anything.


>
> >
> > 3-1) Is this generally any better than solution 2, in terms of not
> messing too much with the condition number of the matrix?
> >
> > 3-2) I don't quite know how to implement this using PETSc. Generally
> speaking I'd like to solve
> >
> > | A        U |   | X |   | B |
> > |            | * |   | = |   |
> > | U^T      0 |   | s |   | 0 |
> >
> >
> > where U is a constant vector (of ones) and s is effectively a Lagrange
> multiplier. I suspect I need to use MatCreateSchurComplement and pass that
> to the KSP? Or do I need create my own matrix type from scratch through
> MatCreateShell?
> >
> > Any help is appreciated!
> >
> > Thanks,
> > Mohammad
> >
> >
>
>

-- 
Sent from Gmail Mobile
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20160224/79b75322/attachment.html>