iRASPA / RASPA2

Classical molecular simulation code
Other
129 stars 83 forks source link

corrected the formula to compute the reduced mass #19

Closed nakulrampal closed 3 years ago

nakulrampal commented 3 years ago

When using the FEYNMAN_HIBBS_LENNARD_JONES potential (and its other variations) to compute the van der Waals potential, the formula to compute the reduced mass was incorrect in the ReadForceFieldDefinitionsMixingRules function. This has been corrected for both when GeneralMixingRule=LORENTZ_BERTHELOT, and when GeneralMixingRule=JORGENSEN.

Note: When using both the force_field_mixing_rules.def and force_field.def files to define the force field parameters, you do not encounter this bug as in the ReadForceFieldDefinitions function, the formula to compute the reduced mass is correct.

dubbelda commented 3 years ago

Thanks!

nakulrampal commented 3 years ago

No worries, always happy to help - Thanks Prof. @dubbelda for creating such a valuable tool for the community!

danieleongari commented 3 years ago

Hi @nakulrampal I think there may be a misunderstanding in your issue. Indeed, the mixed potential is computed for the atom j being without the third parameter of FEYNMAN_HIBBS_LENNARD_JONES: since the line of code you corrected is under the condition

if((PotentialType[i][i]==FEYNMAN_HIBBS_LENNARD_JONES)&&(PotentialType[j][j]==LENNARD_JONES))

there is no p2 paremeter (the "self-reduced mass" of the atom type) for LENNARD_JONES, and therefore your correction

PotentialParms[i][j][2]=1.0/(1.0/PotentialParms[i][i][2]+1.0/PotentialParms[j][j][2])

can not work because there is no PotentialParms[j][j][2].

The idea of using 2.0*PotentialParms[i][i][2] in the original code works well for H2, because it assumes that the second "self-" reduced mass is big and therefore (2*big)/(2+big)~2=2.0*1 (where 2 is the mass of H2 and 1 its self-reduced mass p2=2*2/2+2=1). At least I interpreted in this way, which allows you to get some reasonable results even if you don't specify all the FEYNMAN_HIBBS_LENNARD_JONES w/reduced mass for all the atoms of the framework.

(just as a follow up, we once investigated that indeed the main contribute to the uptake of FEYNMAN-HIBBS is for H2-H2 interactions, so that one can avoid to specify the reduced mass for all the atoms of the frameworks... unless the framework is particularly reach of H atoms. See https://github.com/lsmo-epfl/aiida-lsmo/issues/64)

ltalirz commented 3 years ago

Perhaps an example input would be useful?

Just reiterating that I'd be happy to help set up continuous integration tests on github actions, if there are some tests that should run on every commit (or at least for every release on conda-forge) https://github.com/iRASPA/RASPA2/issues/12

dubbelda commented 3 years ago

I agree that the potential is in principal only defined and developed for hydrogen-hydrogen interaction. How to mix these with other potentials and use them for adsorption in a framework is open for debat. @danieleongari Thanks for reminding us where the factor 2 came from. I also interpreted this that way and I think it was used like this for hydrogen adsorption in MOFs in early papers on this.

I have reverted the patch for now to have the old behavior back.

@nakulrampal an option would be to create a new potential type for this if you would like to add a different definition. For example: FEYNMAN_HIBBS_LENNARD_JONES3 or something like that.

nakulrampal commented 3 years ago

Apologies for the delay in getting back - I just saw this!

Thanks, @danieleongari for bringing this up + the detailed explanation, I really appreciate it.

Prof. @dubbelda - Thanks, and I agree, it is still open to debate :)