Add association term of Dufal

[ianhbell]

For combination with SAFT-VR-Mie. From https://doi.org/10.1080/00268976.2015.1029027

[longemen3000]

https://github.com/ClapeyronThermo/Clapeyron.jl/blob/master/src/models/SAFT/SAFTVRMie/variants/SAFTVRMie15.jl

[ianhbell]

@longemen3000 I implemented this on the Dufal branch, would you care to take a look? The code is a bit messy because I am trying to keep the compilation time to a reasonable level which makes it not possible to use the nice C++ techniques I have been using elsewhere (one must keep the number of template instantiations to a minimum, which implies the use of runtime branching rather than compiler branching).

I have the experience that the code with this model is about 10x slower than the normal association model for pure water. Is this your experience?

[longemen3000]

on the current Clapeyron impl, there was a bug that makes that specific impl slower (some julia specifics). after fixing that, i got the following timings:

julia> @btime Clapeyron.Δ(vr,3e-5,300,Clapeyron.SA[1.0],1,1,1,2,data_vr)
  466.327 ns (1 allocation: 16 bytes)
1.66397994073738e-27

julia> @btime Clapeyron.Δ(vr15,3e-5,300,Clapeyron.SA[1.0],1,1,1,2,data_vr15)
  823.750 ns (1 allocation: 16 bytes)
9.462182925022677e-28

so, around 2x. (i did not take care of using an appropiate reduction technique to skip the repeated powers, but that can be done on both mie impls).

i also checked with the dufal paper, and it seems that i can reproduce the data of association non bonded fractions in table 4:

using Clapeyron, Test
v15 = 1/(1000*1000*1.011/18.015)
T15 = 290.0
vr15 = SAFTVRMie15("water")
#Dufal, table 4, 290K, f_OH(free) = 0.089
@test Clapeyron.X(vr15,V,T,z1)[1][1] ≈ 0.08922902098124778 rtol = 1e-6

[ianhbell]

After fixing the erratum in Dufal that they mention in their corrigendum, I get something more reasonable:

-------------------------------------------------------------------------------
Benchmark association evaluations
-------------------------------------------------------------------------------
/Users/ihb/Documents/Code/teqp/src/tests/catch_test_association.cxx:197
...............................................................................

benchmark name                       samples       iterations    est run time
                                     mean          low mean      high mean
                                     std dev       low std dev   high std dev
-------------------------------------------------------------------------------
time Delta with canonical                      100           393     1.6113 ms 
                                        40.7008 ns    40.6701 ns    40.7421 ns 
                                       0.179534 ns   0.136087 ns   0.275853 ns 

time Delta with Dufal                          100             9     1.7406 ms 
                                        1.95722 us    1.93477 us    2.01458 us 
                                        160.408 ns    15.5912 ns    291.486 ns 

The timing is for the entire construction of the matrix, which is incompletely optimized in my case since all four entries of the 4x4 matrix are identical and I do the same calculation 4 times (doh).

My matrices are:

1.82698e-27 1.82698e-27
1.82698e-27 1.82698e-27

for the normal one and

9.46185e-28 9.46185e-28
9.46185e-28 9.46185e-28

for Dufal. I don't know why there is a small but nonzero difference. Maybe it is about the value of R used? I am surprised that your "normal" implementation is so slow, since almost nothing is happening in the evaluation function, at least compared with Dufal, and that looks like an opportunity for further optimization on the Clapeyron side.