Neumann BCs in problem helper

[brownbaerchen]

I believe I had an error in the BCs for quench, so I tried to do it better. But I am very self conscious about this, so I would appreciate if somebody checked that this is actually correct.

So what I am doing to get Neumann BCs now?

• Generate FD matrix by putting the stencil on the diagonals
• Replace first and last line in the FD matrix with a one-sided stencil for the first derivative of the same order as the rest
• In the problem class: I have sth. like f(u) = Au + g(u), where A is the FD matrix and g is a nonlinear term. Here, I set the value for g at the first and last entry to the Neumann BC value (0 in this case)
• The grid has spacing dx = L / (N-1) and includes both boundaries

Does this make sense? I also don't know about the Dirichlet BCs, but I am not using them now...
If somebody has a good reference, I could maybe double check myself.

There are a couple of other minor things in this PR, such as inexactness in the linear solver for the quench problem, but don't worry about this now. This is just an artefact of merging part of a development branch and makes little difference at this point.

[codecov]

Codecov Report

Attention: 5 lines in your changes are missing coverage. Please review.

Comparison is base (f5658fc) 73.26% compared to head (039d2f5) 72.62%.
Report is 1 commits behind head on master.

Additional details and impacted files

@@            Coverage Diff             @@
##           master     #364      +/-   ##
==========================================
- Coverage   73.26%   72.62%   -0.65%     
==========================================
  Files         270      269       -1     
  Lines       23097    22364     -733     
==========================================
- Hits        16922    16241     -681     
+ Misses       6175     6123      -52     

Files	Coverage Δ
pySDC/helpers/problem_helper.py	97.05% <100.00%> (+1.97%)	⬆️
...implementations/problem_classes/AllenCahn_2D_FD.py	95.28% <100.00%> (ø)
...roblem_classes/GeneralizedFisher_1D_FD_implicit.py	96.49% <100.00%> (ø)
...tations/problem_classes/HeatEquation_ND_FD_CuPy.py	0.00% <ø> (ø)
...C/implementations/problem_classes/generic_ND_FD.py	94.11% <100.00%> (-0.41%)	⬇️
pySDC/implementations/problem_classes/Quench.py	76.31% <92.85%> (+2.93%)	⬆️
pySDC/projects/Resilience/quench.py	49.76% <80.00%> (ø)
pySDC/projects/Resilience/strategies.py	80.41% <0.00%> (ø)
...implementations/problem_classes/AllenCahn_1D_FD.py	0.00% <0.00%> (ø)

... and 10 files with indirect coverage changes

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

[tlunet]

You don't have to put the escape characters here, since you use ' to define the string. This is enough :

raise NotImplementedError(f'Boundary conditions "{bc}" not implemented.')

[brownbaerchen]

The implementation is not totally generic as of now. For high order you need to use lopsided stencils near the boundary because they must not stretch beyond the boundary. This is not particularly difficult, but requires some patience which I don't have this morning. For now, I think second order is fine and if anyone needs higher order, some basic structure is already available.
There is a test checking that this kind of generic formula works for advection and the heat equation for second order. It compares to reference values that I computed by hand. The advection result is different than what Thibaut wrote here earlier because I use a second order stencil for the BC rather than a first order one.

What I find neat is that GMRES no longer requires a preconditioner to converge in the quench problem. Using the inverted Laplacian still seems to increase performance slightly. But I suspect the solution was insufficiently smooth before because the FD discretisation was wrong. I don't know if this would keep GMRES from converging. But that convergence works better is a good sign anyways. Additionally, I set higher threshold values for the transition to runaway heating to keep the timescales the same. This indicates less dissipation across the boundary, which also means the insulation works better.

I am now at least 83% confident that this implementation is correct and suggest, if you feel confident as well, that we merge this now and worry about higher order discretisation some other time. Maybe in Darmstadt already.

[brownbaerchen]

I don't know why this test has failed, I believe something is wrong with the modules. On my fork I reran the tests on the master branch with the same result as here: The macOS environment is fine but on linux, the fenics tests fail.

[pancetta]

I have no clue. A diff between both miniconda list outputs does not reveal anything.. mpi4py went from 3.1.4 to 3.1.5 and there seems to be some strange changes in impi and tbb, but this should also affect the mpi test cases (which it doesn't).

[pancetta]

Hm.. this impi_rt package is new (and not there for the mpi tests). Maybe that's causing this?

[brownbaerchen]

I noticed that the generation of the preconditioner was dominating the runtime, but speeding up the GMRES solves significantly. I read that inverting sparse matrices can be more efficient when they are converted to dense first and inverted using dense matrix methods if the result is expected to be dense. Indeed, I measured a 10x speedup for some configurations just by doing this. Since the branch I want to merge is called "quench_fix", I thought I squeeze it into this PR.

[pancetta]

Please merge from master

[brownbaerchen]

This is actually really helpful for me. I am playing around with higher order discretizations for the quench problem and finally have some confidence in the results. Hooray for testing. Thanks a lot @tlunet! I had a really good time implementing this with you!

[pancetta]

Ready to be merged, then?

[brownbaerchen]

Ready to be merged, then?

I think so