-
Notifications
You must be signed in to change notification settings - Fork 700
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
changing DL_Poly units #4983
base: develop
Are you sure you want to change the base?
changing DL_Poly units #4983
Conversation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Hello there first time contributor! Welcome to the MDAnalysis community! We ask that all contributors abide by our Code of Conduct and that first time contributors introduce themselves on GitHub Discussions so we can get to know you. You can learn more about participating here. Please also add yourself to package/AUTHORS
as part of this PR.
Codecov ReportAll modified and coverable lines are covered by tests ✅
Additional details and impacted files@@ Coverage Diff @@
## develop #4983 +/- ##
===========================================
- Coverage 93.42% 93.41% -0.02%
===========================================
Files 177 189 +12
Lines 21859 22925 +1066
Branches 3078 3078
===========================================
+ Hits 20422 21415 +993
- Misses 986 1059 +73
Partials 451 451 ☔ View full report in Codecov by Sentry. 🚀 New features to boost your workflow:
|
Based on feedback about TRZ being dropped, I removed the parts regarding TRZ in the latest commit. |
@cbouy and @talagayev could you please review and shepherd this PR? This is a GSOC applicant for one of your projects. |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Good start! See comments. Although dividing by 100 gives the right answer, it is a bit more complicated because we want all our readers to be have in a unified manner. Have a look at other readers such as TRRReader where you define the force units
mdanalysis/package/MDAnalysis/coordinates/TRR.py
Lines 150 to 151 in 66e7e5c
units = {'time': 'ps', 'length': 'nm', 'velocity': 'nm/ps', | |
'force': 'kJ/(mol*nm)'} |
For DL_POLY, this entry should be "force": "Da*Angstrom/ps"
in the _DLPOLY_UNITS
dict.
You then need to add at least one entry for "Da*Angstrom/ps" to units.forceUnit_factor
mdanalysis/package/MDAnalysis/units.py
Lines 341 to 350 in 66e7e5c
forceUnit_factor = { | |
"kJ/(mol*Angstrom)": 1.0, | |
"kJ/(mol*A)": 1.0, | |
"kJ/(mol*\u212b)": 1.0, | |
"kJ/(mol*nm)": 10.0, | |
"Newton": 1e13 / constants["N_Avogadro"], | |
"N": 1e13 / constants["N_Avogadro"], | |
"J/m": 1e13 / constants["N_Avogadro"], | |
"kcal/(mol*Angstrom)": 1 / constants["calorie"], | |
} |
Then in your code, check if unit conversion is requested (similar to what the TRRReader does in
mdanalysis/package/MDAnalysis/coordinates/TRR.py
Lines 213 to 214 in 66e7e5c
if self.convert_units: | |
self.convert_forces_from_native(ts.forces) |
Additionally:
- The history reader should also convert forces
- add a CHANGELOG entry under Fixes
- add a versionchanged 2.10.0 under ConfigReader (and also HistoryReader)
@@ -605,4 +605,4 @@ def close(self): | |||
if self.trzfile is None: | |||
return | |||
self.trzfile.close() | |||
self.trzfile = None | |||
self.trzfile = None |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
just undo
@@ -70,8 +70,9 @@ def test_velocities(self, ts): | |||
ref = np.array([1.056610291, -1.218664448, 3.345828610]) | |||
assert_allclose(ts._velocities[0], ref) | |||
|
|||
# modified to original numbers / 100 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
remove comment – we can always look in the git blame to find out when something was changed
@@ -118,8 +119,9 @@ def test_vel(self, u): | |||
ref = np.array([2.637614561, 0.5778767520e-01, -1.704765568]) | |||
assert_allclose(u.atoms[2].velocity, ref) | |||
|
|||
# modified to original numbers / 100 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
remove comment
@@ -99,7 +99,7 @@ def _read_first_frame(self): | |||
if has_vels: | |||
velocities = np.array(velocities, dtype=np.float32) | |||
if has_forces: | |||
forces = np.array(forces, dtype=np.float32) | |||
forces = np.array(forces, dtype=np.float32) / 100 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
add a comment above this line stating that forces in DL_POLY are 10 J/(mol.Å) and mention issue #2393 — this is weird enough to warrant a hint for future coders
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Actually, we should not just divide by 100 here. Instead we should use the normal conversion machinery where you store the native force units in _DLPOLY_UNITS
and store the conversion factor in units, then use convert_from_native.
See also https://userguide.mdanalysis.org/stable/units.html .
Fixes #2393
Changes made in this Pull Request:
Removed forces from TRZAs mentioned in #2393, this likely requires more work with test cases and code modifications regarding where the DL_Poly conversions are made, so any feedback is appreciated!
PR Checklist
package/CHANGELOG
file updated?package/AUTHORS
? (If it is not, add it!)Developers Certificate of Origin
I certify that I can submit this code contribution as described in the Developer Certificate of Origin, under the MDAnalysis LICENSE.
📚 Documentation preview 📚: https://mdanalysis--4983.org.readthedocs.build/en/4983/