andysim / helpme Goto Github PK
View Code? Open in Web Editor NEWhelPME: an efficient library for particle mesh Ewald
License: BSD 3-Clause "New" or "Revised" License
helPME: an efficient library for particle mesh Ewald
License: BSD 3-Clause "New" or "Revised" License
Hi Andy,
I'm running into energy dependence on the kappa parameter when dipoles are present. I've noticed that helPME does not have kappasweep tests for this case so I'm wondering if you've verified it. When testing, I did add the E_self(mu) but it didn't help much. I also tried with a wide range of k-space grids and with a very large real-space cutoff, and both with zero and non-zero total dipole moment to no avail. Below is a table of one of the test runs for a water dimer from your unit tests (including coordinates, charges and dipole moments) with different kappas (called alpha). All quantities are in atomic units (Hartree and Bohr).
grid dim | alpha [a.u.] | GridSp [a.u.] | range [a.u.] | E_self [a.u.] | E_real [a.u.] | E_rec [a.u.] | E_total [a.u.] | T(real) [s] | T(rec) [s] |
---|---|---|---|---|---|---|---|---|---|
20 | 0.02 | 2 | 283.936361 | -0.023545708 | -0.393121622 | 1E-10 | -0.41666733 | 0.005952 | 0.003026 |
30 | 0.03 | 1.333333 | 189.290907 | -0.035318903 | -0.381349666 | 1.3161E-06 | -0.416667253 | 0.002008 | 0.001949 |
48 | 0.045 | 0.888889 | 126.193938 | -0.052979502 | -0.363781306 | 8.88422E-05 | -0.416671965 | 0.001029 | 0.004987 |
64 | 0.0675 | 0.592593 | 84.1292922 | -0.079473126 | -0.338054364 | 0.000779148 | -0.416748342 | 0.001016 | 0.011076 |
96 | 0.10125 | 0.395062 | 56.0861948 | -0.119222761 | -0.300991047 | 0.003060132 | -0.417153676 | 0 | 0.039894 |
144 | 0.151875 | 0.263374 | 37.3907965 | -0.178878262 | -0.248487195 | 0.008947323 | -0.418418134 | 0 | 0.135686 |
216 | 0.2278125 | 0.175583 | 24.9271977 | -0.268466299 | -0.176956275 | 0.023667676 | -0.421754898 | 0 | 0.420907 |
324 | 0.34171875 | 0.117055 | 16.6181318 | -0.403202007 | -0.093209204 | 0.065818528 | -0.430592682 | 0 | 1.399255 |
486 | 0.512578125 | 0.078037 | 11.0787545 | -0.606499142 | -0.026404418 | 0.182708722 | -0.450194838 | 0 | 4.846253 |
750 | 0.768867188 | 0.052025 | 7.38583635 | -0.91547316 | -0.001493128 | 0.440175597 | -0.476790692 | 0 | 18.04084 |
I would appreciate if you'd let me know if I'm missing something and maybe add a kappasweep test with dipoles.
Kind regards,
Alexei
Edit: the GridSp is approximate. It's an input parameter used to compute the grid dimensions.
Edit2: the splineOrder was 7 in all tests but I've also run with higher orders with the same results.
Hi Andy,
I have a question regarding reciprocal-space dipole contributions to the electric field (computePRec method). It's the last term in the Eq (49) https://doi.org/10.1063/1.3511713 for electric field due to dipole moments (I guess it's the field due to interaction the dipole's images in other cells?). I'm not an expert in the Ewald summation to judge if the articles cited in the helPME docs have or not the term but I did run some comparison between helPME and another (non-PME) Ewald implementation where this term can be switched off. The results without the term were identical to helPME, which makes me think that this term might be missing. Am I missing something or is there a mistake in the paper by Sala et. al?
Kind regards,
Alexei
Does this library provide the periodic effective potential matrix (i.e. the 1/r_{ij}
matrix for molecular system) at each atomic/charge site in the system?
In other libraries e.g. OpenMM, sampling the potential (i.e. get q_j/r_{ij}
at site i
) at a charge site just returns infinity and it would be awesome to have this feature.
Many thanks!
Dear developers,
Shouldn't the sign of force be opposite to the default one?
The following code describes how I calculated numerical differentiation and compared it to the analytical one.
helpme::Matrix<double> coordsD(
{{2.0, 2.0, 2.0}, {2.5, 2.0, 3.0}, {1.5, 2.0, 3.0}, {0.0, 0.0, 0.0}, {0.5, 0.0, 1.0}, {-0.5, 0.0, 1.0}});
helpme::Matrix<double> chargesD({-0.834, 0.417, 0.417, -0.834, 0.417, 0.417});
double scaleFactorD = 332.0716;
double energyD = 0;
helpme::Matrix<double> forcesD(6, 3);
helpme::Matrix<double> virialD(1, 6);
helpme::Matrix<double> potentialAndFieldD(6, 4);
auto pmeD = std::unique_ptr<PMEInstanceD>(new PMEInstanceD);
pmeD->setup(1, 0.3, 5, 32, 32, 32, scaleFactorD, 1);
pmeD->setLatticeVectors(20, 20, 20, 90, 90, 90, PMEInstanceD::LatticeType::XAligned);
// Compute the energy and forces
energyD = pmeD->computeEFRec(0, chargesD, coordsD, forcesD);
std::cout << "force_analytical: " << forcesD[0][0] << '\n';
double alpha = 1e-5;
coordsD[0][0] += alpha;
double ene_p = pmeD->computeEFRec(0, chargesD, coordsD, forcesD);
coordsD[0][0] -= 2. * alpha;
double ene_m = pmeD->computeEFRec(0, chargesD, coordsD, forcesD);
double numerical = -(ene_p - ene_m) * 0.5 / alpha;
std::cout << "force_numerical: " << numerical << '\n';`
The output is::
force_analytical: 0.603153 force_numerical: -0.603153
I'm sorry if I have overlooked something.
Best,
Koichi
A declarative, efficient, and flexible JavaScript library for building user interfaces.
๐ Vue.js is a progressive, incrementally-adoptable JavaScript framework for building UI on the web.
TypeScript is a superset of JavaScript that compiles to clean JavaScript output.
An Open Source Machine Learning Framework for Everyone
The Web framework for perfectionists with deadlines.
A PHP framework for web artisans
Bring data to life with SVG, Canvas and HTML. ๐๐๐
JavaScript (JS) is a lightweight interpreted programming language with first-class functions.
Some thing interesting about web. New door for the world.
A server is a program made to process requests and deliver data to clients.
Machine learning is a way of modeling and interpreting data that allows a piece of software to respond intelligently.
Some thing interesting about visualization, use data art
Some thing interesting about game, make everyone happy.
We are working to build community through open source technology. NB: members must have two-factor auth.
Open source projects and samples from Microsoft.
Google โค๏ธ Open Source for everyone.
Alibaba Open Source for everyone
Data-Driven Documents codes.
China tencent open source team.