Fix vdW in the density overlap model

[NikoOinonen]

The vdW damping used in the full-density based model that was originally discussed in #48 was not properly tested for choosing proper values for the damping constant. There was some more discussion in #156 that the vdW contribution does not agree with what was reported by Ellner. et al in https://pubs.acs.org/doi/10.1021/acsnano.8b08209.

Edit: the original idea was to choose the damping constant, but in the end, implementing the DFT-D3 calculation was the best way to go.

[NikoOinonen]

@ProkopHapala Let's continue the discussion on the specific topic of the vdW damping constants here.

• I would perhaps first try what values roughly looks good, and only than fit it. It will help to do the fit in some reasonable range.
• Also I'm not sure it is necessary to fit the whole DFT force/energy, vdW cotribution form DFT-D3 (or other vdW correction) should be written separately by any decent code.
  • This may be useful to calculate DFT-D3 contribution separately without running DFT code
• https://github.com/dftd3/simple-dftd3
• https://github.com/dftd3/tad-dftd3

Yes, I was planning on first just trying to get something that roughly matches. But in the end, the value needs to be fitted to something. Should we just choose the damping to match DFT-D3 as closely as possible?

[NikoOinonen]

Getting something reasonable out of the vdW by adjusting the damping constants turns out to be problematic. For example, using the R4 damping function and setting the damping function too high results in these ring-like artifacts: They probably appear because the function actually turns back up when getting close: Something similar happens with the other R* damping functions and there does not seem to exist a range where there are no artifacts but the vdW is not too strong.

The only way that works kind of ok is using the constant damping function with a really high damping constant, like ~3000, which almost kills the vdW entirely: Maybe that's actually fine, I don't know. I will try and see what the DFT-D3 curves would look like next.

[ProkopHapala]

The only way that works kind of ok is using the constant damping function with a really high damping constant, like ~3000, which almost kills the vdW entirely:

And the constant damping with the original parameters (before we introduced the new damping functions) ?

[NikoOinonen]

And the constant damping with the original parameters (before we introduced the new damping functions) ?

It was (and is) 180.

[NikoOinonen]

So, I tried computing the vdW with DFT-D3 using https://github.com/dftd3/tad-dftd3. The overall magnitude of the resulting vdW force is actually pretty close to that constant-damped (3000) case above: But it's also more concentrated around the molecule. The contrast for the total force looks very close to what they have in the Ellner paper, and the relative magnitudes of the forces are similar, so I think this overall looks pretty good.

Using that super damped LJ vdW might be okay, but it should not be that difficult to implement the DFT-D3 either. I was reading the original Grimme-D3 paper (https://pubs.aip.org/aip/jcp/article/132/15/154104/926936/A-consistent-and-accurate-ab-initio), and it looks to me like there are three major differences compared to LJ:

• There is a 1/r^8 term with it's own C8 parameter (which is derivable from the C6 parameter)
• The C6 parameter depends on the coordination number of the atom, which is determined by the pair-wise distances between the atoms.
• There is a 3-body term. However, the paper mentions that the correction is pretty small, and it looks like VASP doesn't use it: https://www.vasp.at/wiki/index.php/DFT-D3 (and Ellner et al. used VASP), so I think this can be ignored.

The first one is pretty trivial to add. The second point is maybe a little bit tricky. When I tried calculating the DFT-D3 vdW at a closer range, it makes for some pretty funky behaviour: It actually turns repulsive. The origin of the increase in energy is the changing coordination number of the PP and the sample atoms when the PP approaches the sample. We assume that the PP is chemically inert, so we would need to exclude the PP from the calculation of the coordination number to prevent this repulsion.

[ProkopHapala]

WOW, so we can use Grimme DFT-D3 now ?

[NikoOinonen]

No, it's not implemented yet. I just did some testing with that tad-dftd3 library.

[NikoOinonen]

Now the DFT-D3 is implemented in the OpenCL version. Probably this should be implemented for the CPU code as well, but maybe it's not a priority. I created a draft PR at #163. We could merge that now or wait until the CPU version is there as well.

[ProkopHapala]

Now the DFT-D3 is implemented in the OpenCL version. Probably this should be implemented for the CPU code as well, but maybe it's not a priority. I created a draft PR at #163. We could merge that now or wait until the CPU version is there as well.

• I would merge it now, not wait imlemetation on CPU.
• but I think it wold be better if the old version of vdW is still working so that we can easily switch between them and compare the results
• Also I'm currios, it is possible to say in simple terms why old version of vdW does not work well? Isn't it really just about choosing reasonoble parameters? Is there some fundamental problem what is missing ?