coulomb 0.4.0

<p align="center">
  <img src="https://upload.wikimedia.org/wikipedia/commons/thumb/e/ed/VFPt_charges_plus_minus_thumb.svg/440px-VFPt_charges_plus_minus_thumb.svg.png?20200314191255" alt="crates.io", height="200">
</p>
<p align="center">
    <a href="https://opensource.org/licenses/Apache-2.0">
        <img src="https://img.shields.io/badge/License-Apache%202.0-blue.svg">
    </a>
</p>

-----

<p align = "center">
<b>Coulomb</b></br>
<i>A Library for Electrolyte Solutions and Electrostatic Interactions</i>
</p>

-----

# Features

**Coulomb** is a library for working with electrolyte solutions and calculating electrostatic interactions
in and between molecules and particles.
The main purpose is to offer support for molecular simulation software.

- Temperature dependent dielectric permittivity models for common solvents.
- Handling of ionic strength and the Debye screening length.
- Automatic stoichiometry deduction for arbitrary salts.
- Extensive library of truncated electrostatic interaction schemes such as
  _Wolf_, _Reaction field_, _Real-space Ewald_, generalized through a short-range function `trait`.
- Reciprocal-space Ewald summation for Coulomb and Yukawa potentials (PBC and IPBC),
  with automatic parameter optimization, mixed f32/f64 precision, and SIMD acceleration on x86_64.
- Ewald summation with and without implicit salt.
- Multipole expansion for _energies_, _forces_, _fields_ between ions, dipoles, and quadrupoles.
- Extensively unit tested and may serve as reference for other implementations or approximations.
- Partial support for static unit of measure analysis via [`uom`](https://crates.io/crates/uom).
  To enable, use the `uom` feature flag.
- Vector types use [`mint`](https://crates.io/crates/mint) for interoperability with
  [glam](https://crates.io/crates/glam), [cgmath](https://crates.io/crates/cgmath),
  [nalgebra](https://crates.io/crates/nalgebra), and other math libraries.

This is largely a Rust rewrite and extension of the
[CoulombGalore](https://zenodo.org/doi/10.5281/zenodo.3522058) C++ library.

# Examples

## Dielectric Media and Electrolytes

Simple polynomial models are provided to obtain the relative permittivity or a `Medium` as a function
of temperature.
For working with the _ionic strength_, `Salt` of arbitrary valency can be given and the
required stoichiometry is automatically worked out.

~~~ rust
use coulomb::{Medium, Salt};
let molarity = 0.1;
let medium = Medium::salt_water(298.15, Salt::CalciumChloride, molarity);
assert_eq!(medium.permittivity()?, 78.35565171480539);
assert_eq!(medium.ionic_strength()?, 0.3);             // mol/l
assert_eq!(medium.debye_length()?, 5.548902662386284); // angstrom
~~~

## Multipolar Interactions

All pairwise _schemes_ support calculation of _potential_, _energy_, _field_, _force_ from or between multipolar particles, up to second order (ion-ion, ion-dipole, dipole-dipole; ion-quadrupole).
Most scheme can be evaluated with or without a Debye-Hückel screening length.

Vector arguments accept any type implementing `Into<mint::Vector3<f64>>`, including simple arrays
and vectors from nalgebra, glam, cgmath, and other libraries with mint support.

~~~ rust
use coulomb::pairwise::*;
let scheme = Plain::default();                     // Vanilla Coulomb scheme, 𝒮(𝑞)=1
let z = 1.0;                                       // point charge, 𝑧 
let r = [3.0, 5.0, 0.0];                           // distance vector, 𝒓
let r_len = (r[0]*r[0] + r[1]*r[1] + r[2]*r[2]).sqrt();
let ion_pot = scheme.ion_potential(z, r_len);      // potential |𝒓| away from charge 
assert_eq!(ion_pot, z / r_len);

let mu = [0.2, 3.0, -1.0];                         // point dipole, 𝝁
let dipole_pot = scheme.dipole_potential(mu, r);   // potential 𝒓 away from dipole
let energy = scheme.ion_dipole_energy(z, mu, r);   // interaction energy assuming 𝒓 = 𝒓(𝜇) - 𝒓(𝑧)
~~~

The image below is generated by `examples/potential.rs` and shows the calculated
electric potential on a plane containing a monopole and a dipole.

![Electric Potential](./assets/potential.png)

## Reciprocal-Space Ewald Summation

The `reciprocal` module provides k-space Ewald summation for Coulomb (κ=0) and Yukawa (κ>0)
potentials with automatic parameter selection, O(N_k) Monte Carlo trial moves, and mixed-precision
SIMD acceleration.

~~~ rust
use coulomb::reciprocal::*;

// Automatic parameter selection for cutoff 5 Å and target accuracy 1e-5
let mut ewald = EwaldReciprocal::new([10.0; 3], 5.0, 1e-5, None);

let (x, y, z) = (vec![-0.5, 0.5], vec![0.0, 0.0], vec![0.0, 0.0]);
let charges = vec![1.0, -1.0];

// Build structure factors from all particle positions
ewald.update_all(&x, &y, &z, &charges);

// Reciprocal-space energy + self-energy correction
let energy = ewald.energy() + ewald.self_energy(&charges);

// MC trial move: compute ΔE without mutating state
let de = ewald.energy_change(1.0, [-0.5, 0.0, 0.0], [-0.3, 0.0, 0.0]);

// Forces for MD
let (mut fx, mut fy, mut fz) = (vec![0.0; 2], vec![0.0; 2], vec![0.0; 2]);
ewald.add_forces(&x, &y, &z, &charges, &mut fx, &mut fy, &mut fz);
~~~

Both `f32` and `f64` input slices are accepted (generic over `Into<f64> + Copy`).
On x86_64 with the `simd` feature, a polynomial sin/cos with SIMD vectorization accelerates
the inner loops.

### Unit analysis

Experimental support for static unit analysis can be activated with the `uom` feature.

~~~ rust
use coulomb::{units::*, pairwise::{Plain, MultipoleEnergySI}};
let z1 = ElectricCharge::new::<elementary_charge>(1.0);
let z2 = ElectricCharge::new::<elementary_charge>(2.0);
let r = Length::new::<nanometer>(2.3);
let energy = Plain::without_cutoff().ion_ion_energy(z1, z2, r);
assert_eq!(energy.get::<kilojoule_per_mole>(), 362.4403242896922);
~~~