Skip to content
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

Add class method for Box class which create a box from Lx,Ly,Lz and angles #1234

Merged
merged 42 commits into from
Apr 25, 2024

Conversation

DomFijan
Copy link
Contributor

@DomFijan DomFijan commented Mar 7, 2024

Implements classmethods that create a box class from lattice vector matrix and L1, L2, L3 and angles, as well as methods that give L1, L2, L3 and angles from box class.

Description

Crystals are usually defined via their lattice vectors and wykoff sites or lengths of lattice vectors and angles; while freud and hoomd use cell parameter matrix. In this PR I implement class method for creation of box from lattice vectors or from lattice vector lengths and angles. I've also implemented a method that gives lengths of box vectors and angles from box object.

Motivation and Context

Lattice vectors or lattice vector lengths and angles are used to define real crystals and having a method that works with this makes creating crystals easier.

How Has This Been Tested?

Tests were added for new methods for box testing.

Types of changes

  • Bug fix (non-breaking change which fixes an issue)
  • New feature (non-breaking change which adds or improves functionality)
  • Breaking change (fix or feature that would cause existing functionality to change)
  • Documentation improvement (updates to user guides, docstrings, or developer docs)

Checklist:

  • I have read the CONTRIBUTING document.
  • My code follows the code style of this project.
  • I have updated the documentation (if relevant).
  • I have added tests that cover my changes (if relevant).
  • All new and existing tests passed.
  • I have updated the credits.
  • I have updated the Changelog.

@DomFijan DomFijan marked this pull request as ready for review March 20, 2024 13:56
@DomFijan DomFijan requested a review from a team as a code owner March 20, 2024 13:56
@janbridley
Copy link
Contributor

I personally am a proponent of using math.sqrt and math.acos, as for single values it is significantly faster than the numpy alternatives. The difference is on the order of μs, but free performance is free imo.

@DomFijan
Copy link
Contributor Author

I personally am a proponent of using math.sqrt and math.acos, as for single values it is significantly faster than the numpy alternatives. The difference is on the order of μs, but free performance is free imo.

It seems numpy sqrt is used throughout freud from what I can tell. If we were to start using math functions we should change all occurrences. I agree it makes no sense to leave performance on the table, but realistically the differences will be so small they are inconsequential. Perhaps this would be a good first issue to tackle for someone new?

freud/box.pyx Outdated
v2 = lattice_matrix[2]
Lx = np.sqrt(np.dot(v0, v0))
a2x = np.dot(v0, v1) / Lx
Ly = np.sqrt(np.dot(v1, v1) - a2x * a2x)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
Ly = np.sqrt(np.dot(v1, v1) - a2x * a2x)
Ly = np.sqrt(np.dot(v1, v1) - a2x**2)

I think a square here better conveys the intent, as below we have a very similar function with two different variables (a2x * a3x)

tests/test_box_Box.py Outdated Show resolved Hide resolved
@janbridley
Copy link
Contributor

janbridley commented Mar 20, 2024

I personally am a proponent of using math.sqrt and math.acos, as for single values it is significantly faster than the numpy alternatives. The difference is on the order of μs, but free performance is free imo.

It seems numpy sqrt is used throughout freud from what I can tell. If we were to start using math functions we should change all occurrences. I agree it makes no sense to leave performance on the table, but realistically the differences will be so small they are inconsequential. Perhaps this would be a good first issue to tackle for someone new?

Noted, and sounds good! I will raise an issue.

Edit: See #1236

Copy link
Collaborator

@tommy-waltmann tommy-waltmann left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

One final thing, then we merge.

Great work on this, I think many people will find this useful!

Comment on lines 497 to 502
1,
2,
3,
np.pi / 3,
np.pi / 2.25,
np.pi / 2.35,
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Each of the six numbers should be the output of a different call to np.random.rand. That way, each time the test is run it gets run with a different random value.

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am currently working on adding this.

freud/box.pyx Outdated
a2 = np.array([L2 * np.cos(gamma), L2 * np.sin(gamma), 0])
a3x = np.cos(beta)
a3y = (np.cos(alpha) - np.cos(beta) * np.cos(gamma)) / np.sin(gamma)
a3z = np.sqrt(1 - a3x**2 - a3y**2)
Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This line appears to take the square root of a negative number at times, resulting in NaNs.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

some combinations of angles are not possible-there is likely a constraint on their sum or something like that. Is this taken into account?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Not in how I wrote the test case, which would explain the failures. The code should throw an error if the user gives inputs which violate the constraints.

@ispivack
Copy link

All looks good to me, just one comment. If we want to just use freud for the EBT box, it would be good to have a convenience function converting to a unit vector + lengths description. Do we feel that this would fit into freud, or should this just be separately implemented in ebonds?

@tommy-waltmann
Copy link
Collaborator

All looks good to me, just one comment. If we want to just use freud for the EBT box, it would be good to have a convenience function converting to a unit vector + lengths description. Do we feel that this would fit into freud, or should this just be separately implemented in ebonds?

I would start by adding the functionality into ebonds. If it becomes cumbersome to have to do the conversion in ebonds constantly, then adding a method to freud to do that might be worth it.

Right now, the best way to get lengths and unit vectors would be to get the box vectors (they are the columns of box.to_matrix()) and get magnitude and direction from that.

@janbridley
Copy link
Contributor

np.testing.assert_allclose sets atol=0 by default. This will not work with floating point calculations, so I set atol to 1e-14. However, the tests are still failing due to mismatched NaNs.

@tommy-waltmann
Copy link
Collaborator

np.testing.assert_allclose sets atol=0 by default. This will not work with floating point calculations, so I set atol to 1e-14. However, the tests are still failing due to mismatched NaNs.

It's fine to leave atol as 0.0. According to the documentation for assert_allclose: "It compares the difference between actual and desired to atol + rtol * abs(desired)", so having an rtol is enough to account for floating point error.

@joaander
Copy link
Member

joaander commented Apr 1, 2024

np.testing.assert_allclose sets atol=0 by default. This will not work with floating point calculations, so I set atol to 1e-14. However, the tests are still failing due to mismatched NaNs.

It's fine to leave atol as 0.0. According to the documentation for assert_allclose: "It compares the difference between actual and desired to atol + rtol * abs(desired)", so having an rtol is enough to account for floating point error.

A non-zero atol is needed when a value close to zero should be accepted as a desired (exactly) zero value. atol + rtol * abs(desired) == 0.0 when atol=0 and desired=0.

@janbridley
Copy link
Contributor

I have tested this method against box lengths and angles from 1,104 real CIF files and all cases pass. I can export the box params as a numpy array and save it as test data, or we can adjust the current test to be more realistic.

@tommy-waltmann
Copy link
Collaborator

I have tested this method against box lengths and angles from 1,104 real CIF files and all cases pass. I can export the box params as a numpy array and save it as test data, or we can adjust the current test to be more realistic.

The code should be written so invalid sets of angle inputs raise an error to the user instead of creating invalid boxes. This PR needs an update to the code, not just the tests, before it's ready to be merged.

Ideally, we would know the mathematical constraints on the box angles so we could check that the inputs are valid, and we should make an effort to figure out what those are. If we cannot get those relations, then we could move forward with an updated set of tests and a "garbage in, garbage out" note in the documentation, but again, I only think we should take the second option if we have already given the first option some effort and failed.

…for angles and raise valueerror when violated
@DomFijan
Copy link
Contributor Author

DomFijan commented Apr 11, 2024

I think I figured it out. The restriction is the following:
$1+2\cos(\alpha)\cos(\beta)\cos(\gamma) - \cos(\alpha)^2 - \cos(\beta)^2 - \cos(\gamma)^2>0$
This is part of the computation of volume of a Parallelepiped. This guarantees that the volume of the box defined is positive and must always be true. This statement should be equivalent to requiring the value under the sqrt in the from_box_lengths_and_angles function be positive (haven't done the maths to confirm that).

The geometric interpretation of this is the violation of the cosine law.

Copy link
Collaborator

@tommy-waltmann tommy-waltmann left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Great work! I know this will be useful for lots of peeps using freud.

@tommy-waltmann tommy-waltmann merged commit 8353049 into main Apr 25, 2024
15 checks passed
@tommy-waltmann tommy-waltmann deleted the feat_from_to_lattice_vectors branch April 25, 2024 16:19
@tommy-waltmann tommy-waltmann added this to the v3.1.0 milestone Jun 13, 2024
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

5 participants