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

Implement Magnetostatics in TorchPME #133

Draft
wants to merge 7 commits into
base: main
Choose a base branch
from
Draft

Implement Magnetostatics in TorchPME #133

wants to merge 7 commits into from

Conversation

E-Rum
Copy link
Contributor

@E-Rum E-Rum commented Dec 23, 2024

This branch is linked to #54 and introduces the initial implementation of magnetostatics in the existing TorchPME library. The commit provides a basic implementation for direct potential calculations along with a simple test case.

With this implementation, it’s a good time to start discussing how to structure the current magnetostatic functionality and align it with the broader scope of the library. Notable changes to the main library include:
1. Dipole Potential Components:
Unlike charge-based potentials (e.g., Coulomb), dipole potentials consist of two components: a scalar and a tensor part. Currently, these components are explicitly returned separately and later combined in the calculator. While this approach works for now, we should evaluate whether it remains viable as we transition to more complex Ewald-based calculators.
2. Calculator Structure for Dipoles:
The dipole calculator significantly differs from charge-based calculators, potentially necessitating a new base calculator tailored for magnetostatics. For now, the implementation resides in the calculators folder, but it may need restructuring.
Key differences include:
• Passing dipoles instead of charges.
• Utilizing neighbor vectors instead of neighbor distances.
• Returning potentials in a vectorial form, which appears more convenient from an ML perspective.


📚 Documentation preview 📚: https://torch-pme--133.org.readthedocs.build/en/133/

@E-Rum E-Rum force-pushed the magnetostatics branch 2 times, most recently from 615c3cd to a73d69e Compare January 17, 2025 15:36
Comment on lines +56 to +55
r_outer = torch.einsum(
"bi,bj->bij", vector, vector
Copy link
Contributor

Choose a reason for hiding this comment

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

Might be worth for speed later to replace the einsum by outer. The latter is usually much faster.

from ..potentials import PotentialDipole


class CalculatorDipole(torch.nn.Module):
Copy link
Contributor

Choose a reason for hiding this comment

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

We should the replace the "current" Calculator with "CalculatorMonopole" ?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Good question! We can discuss it after we merge.

Comment on lines +89 to +97
return scalar_potential.unsqueeze(-1) * torch.eye(3).to(r_outer).unsqueeze(
0
) - r_outer * (
3.0 / (r_mag**5)
- 3.0 * torch.erfc(torch.sqrt(alpha) * r_mag) / r_mag**5
- 2
* torch.sqrt(alpha / torch.pi)
* (2 * alpha + 3 / r_mag**2)
* torch.exp(-alpha * r_mag**2)
/ r_mag**2
).unsqueeze(-1)
Copy link
Contributor

Choose a reason for hiding this comment

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

Can you break this up into more readable lines. 10 lines for a return statement might be a but long.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Makes sense! I will do that after we finish, and we’ll make sure we get the correct physical values. Then, we should do the “linting.”

raise ValueError(
"Cannot compute background correction without specifying `smearing`."
)
return self.smearing * 0.0
Copy link
Contributor

Choose a reason for hiding this comment

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

Always 0?

Suggested change
return self.smearing * 0.0
return 0.0

then we can even remove the value error?

Copy link
Contributor Author

Choose a reason for hiding this comment

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

Here, it’s done this way to stay consistent with the pipeline. This ensures that 0.0 is returned as a torch.tensor with the correct device and dtype; otherwise, we would just return a float number.

Copy link
Contributor Author

Choose a reason for hiding this comment

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

In general, I have no idea whether it should always be zero or not. One would need to perform a theoretical derivation to determine whether dipole neutrality is a crucial aspect for the dipole-Ewald calculator

@E-Rum E-Rum linked an issue Jan 23, 2025 that may be closed by this pull request
@E-Rum
Copy link
Contributor Author

E-Rum commented Jan 23, 2025

Ok, from that point on, it seems to be working. Both the direct and Ewald calculators give reasonable values for two dipoles with a fixed distance. Next, I will add documentation, improve the code structure, and add tests. Meanwhile, I would love to get some help in calculating reference values for comparisons.

@PicoCentauri and I started testing it with the espressomd P3MDipole code, but it seems to tune its parameters for each calculation. We will investigate how to prevent this or find an alternative method for testing.

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.

PME for Magnetostatics
2 participants