Invert Laplace Petsc seems to be wrong

[dschwoerer]

I tried MMS testing the laplacian inversion, but it seems to converge only for some test functions.

As the laplacian operator is implemented as a petsc matrix, it is straight forward to also test the forward operator. Using the forward operator, the bulk solution can be tested, independent of the boundaries. Here as well the operator does not converge to second order, the results seems to be significantly wrong. This is worrying, as it seems the code has not been changed in years, so the inversion might have been wrong for 10+ years.

Has anybody observed that using petsc for the inversion gives a different result?

@bendudson and @totork have also started looking into this ...

My branch which exposes the forward operator (both in C++ and python interface) is fci-auto-with-debug-higher-order

[johnomotani]

I addded LaplacePetsc many years ago, in an attempt to support non-constant coefficients (for non-Boussinesq simulations). At the time, we did not manage to find settings where it converged fast enough to be any use. I assume nobody has ever used it in production. I don't remember now how thoroughly (or even if) I tested the outputs - if I did it was probably only on a square grid (uniform and orthogonal), so there could well be bugs...

The implementation is also pretty clunky - it might even be easier to start over with a new implementation based on a stripped-down version of the 3D PETSc-based solver, which was written more recently and more competently.

[dschwoerer]

Thanks John, that also explains why the matrix is re-created on every iteration :-)

[ZedThree]

If we think this is not super useful, should we just delete it?

[dschwoerer]

I will see whether petsc3damg will work for us, otherwise I will try to fix petsc.

If I do not get it fixed, deleting sounds like the right thing to do 👍