Test and benchmark pressure solver

[marinlauber]

This pull request starts an effort to better benchmark and test the pressure solver in WaterLily.

The idea is to have a test that is strict enough to ensure that the pressure solver tolerances and implementation are efficient. For now, all the tests pass with master for both Float32 and Float64.

I also implemented a wrapper to the HYPRE.jl solver suite that uses the highly optimised hypre library in the background. I think this will be useful to benchmark against (for now, we are slover with both the :Poisson and :MultiLevelPoisson). The caveat of the Hypre wrapper is that it only works for flows without immersed geometries and without periodic BCs. I am working on fixing the latter, but the first one requires a bit more thought and coding.

I think also linking the pressure logging #135 will help here.

Things to do:

• [ ] fix WaterLilyHYPREExt.jl for periodic BCs (merely constructing the sparse matrix correctly)
• [ ] implement a pressure split method to solve non-constant coefficient Poisson problems (immersed bodies)
• [ ] write some proper benchmarks
• [x] link to GeometricMultigrid.jl

[marinlauber]

Benchmark case:

# the classic...
function TGV(; pow=6, Re=1000, T=Float64, mem=Array)
    # Taylor-Green-Vortex initial velocity field
    function u_TGV(i,x,t,ν,κ)
        i==1 && return sin(κ*x[1])*cos(κ*x[2])*exp(-2κ^2*ν*t) # u_x
        return  -cos(κ*x[1])*sin(κ*x[2])*exp(-2κ^2*ν*t)       # u_y
    end
    # Define vortex size, velocity, viscosity
    L = 2^pow; U = 1; ν = U*L/Re
    # make the function
    uλ(i,xy) = u_TGV(i,xy,0,ν,2π/L)
    # Initialize simulation
    return Simulation((L,L),(0,0),L;U,uλ,ν,T,mem,psolver=GeomMultigridPoisson)
end
WaterLily.solver!(p) = WaterLily.solver!(p;tol=1e-9,itmx=32) # strict convergence

Benchmark results: GeomMultigridPoisson

BenchmarkTools.Trial: 10000 samples with 971 evaluations.
 Range (min … max):  75.525 ns …  3.059 μs  ┊ GC (min … max): 0.00% … 94.67%
 Time  (median):     91.051 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   99.256 ns ± 69.229 ns  ┊ GC (mean ± σ):  2.00% ±  2.93%

    ▂▃▄▇█▇▅▄▅▅▅▄▄▃▃▁▁                                    ▁    ▂
  ▆▇████████████████████▇▇▇█▆▆▆▆▆▅▄▆▄▃▄▅▅▄▃▃▁▁▁▄▁▄▅▆▄▄▅▆███▇▇ █
  75.5 ns      Histogram: log(frequency) by time       205 ns <

 Memory estimate: 48 bytes, allocs estimate: 3.

MultiLevelPoisson

BenchmarkTools.Trial: 10000 samples with 966 evaluations.
 Range (min … max):  80.001 ns …  2.862 μs  ┊ GC (min … max): 0.00% … 94.25%
 Time  (median):     91.510 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   96.406 ns ± 52.049 ns  ┊ GC (mean ± σ):  1.53% ±  2.88%

            ▂█▆▅▂                                              
  ▁▁▁▁▁▂▃▃▃▅█████▇▄▃▂▂▂▂▂▂▂▃▃▃▃▃▃▂▂▂▂▂▂▁▂▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁ ▂
  80 ns           Histogram: frequency by time         130 ns <

 Memory estimate: 48 bytes, allocs estimate: 3.

HyprePoisson

BenchmarkTools.Trial: 10000 samples with 962 evaluations.
 Range (min … max):  74.059 ns …  2.531 μs  ┊ GC (min … max): 0.00% … 91.95%
 Time  (median):     92.196 ns              ┊ GC (median):    0.00%
 Time  (mean ± σ):   97.679 ns ± 54.223 ns  ┊ GC (mean ± σ):  1.52% ±  2.75%

          ▁█▆▂                                                 
  ▁▁▂▁▁▂▃▅████▅▃▃▄▅▃▃▂▂▂▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁ ▂
  74.1 ns         Histogram: frequency by time         175 ns <

 Memory estimate: 48 bytes, allocs estimate: 3.

I am not sure why these are the same, some feel slower when I run them...

[weymouth]

That's odd.

Side question, does Hyper work on GPU?

[weymouth]

We decided this should be a different repository, and so will not be merged here, correct?

[marinlauber]

Yes correct, I will move this and close this PR.