Substrate (flake.substrate)

Substrate potential functions for rigid-cluster simulations.

Each substrate type exposes two public functions:

  • particle_en_<type> — per-particle energy, force, and torque contribution

  • calc_en_<type> — total (summed) energy, force, and torque on the CM

Both return (energy, force, torque). The particle-level functions are useful for diagnostics; production code calls calc_en_<type> through the substrate_from_params factory.

Coordinate convention

All positions are 2-D Cartesian \((x, y)\). The cluster CM is at pos_cm; the rigid-body orientation is handled externally. Torque is computed about pos_torque (normally the CM):

\[\tau_i = (\mathbf{r}_i - \mathbf{r}_\mathrm{cm}) \times \mathbf{F}_i\]

Units are consistent with the calling code (distances in μm, forces in fN, energies in zJ for colloidal systems; or Å/eV/nN for nanoscale).

Substrate types

Gaussian well — smooth, localised well vanishing at cutoff \(b\):

\[\begin{split}V(r) = \begin{cases} -\epsilon \, e^{-r^2/(2\sigma^2)} & r \le a \\ -\epsilon \, e^{-r^2/(2\sigma^2)} \, f(\rho) & a < r < b \\ 0 & r \ge b \end{cases}\end{split}\]

where \(\rho = (r-a)/(b-a)\) and \(f\) is the \(C^2\) smoothstep \(f(\rho) = 1 - 10\rho^3 + 15\rho^4 - 6\rho^5\).

Tanh well — flat bottom with steep walls (Cao, Phys. Rev. E 103, 2021):

\[\begin{split}V(r) = \begin{cases} -\epsilon & r \le a \\ \tfrac{\epsilon}{2}\!\left(\tanh\!\frac{\rho - w}{\rho(1-\rho)} - 1\right) & a < r < b \\ 0 & r \ge b \end{cases}\end{split}\]

Sinusoidal (plane-wave) — superposition of \(n\) plane waves (Vanossi, Manini, Tosatti, PNAS 109, 2012):

\[V(\mathbf{r}) = -\frac{\epsilon}{n^2} \left|\sum_{l=0}^{n-1} e^{i \mathbf{k}_l \cdot \mathbf{r}}\right|^2\]

which ranges from \(-\epsilon\) (well minimum) to \(0\) (saddle).

JIT notes

Inner loops use Numba @njit: no NumPy fancy indexing or boolean masks inside jitted functions; explicit for-loops over particles. Non-jitted wrappers handle array allocation and call the jitted core.

flake.substrate.calc_en_gaussian(pos, pos_torque, basis, a, b, sigma, epsilon, u, u_inv)

Total energy, force, and torque on the CM for a Gaussian-well substrate.

Sums particle_en_gaussian over all particles.

Returns:

(float, (2,) ndarray, float) – total energy, total force, total torque.

flake.substrate.calc_en_sin(pos, pos_torque, basis, ks, epsilon)

Total energy, force, and torque on the CM for a plane-wave substrate.

Returns:

(float, (2,) ndarray, float) – total energy, total force, total torque.

flake.substrate.calc_en_tanh(pos, pos_torque, basis, a, b, ww, epsilon, u, u_inv)

Total energy, force, and torque on the CM for a tanh-well substrate.

Returns:

(float, (2,) ndarray, float) – total energy, total force, total torque.

flake.substrate.calc_matrices_bvect(b1, b2)

Metric matrices from two arbitrary primitive lattice vectors.

Args:

b1, b2: array-like, shape (2,).

Returns (u, u_inv) such that u @ r gives fractional coordinates and u_inv @ frac gives real-space coordinates.

flake.substrate.calc_matrices_square(R)

Metric matrices for a square lattice of spacing R.

Returns (u, u_inv): 2x2 matrices that map real-space displacements to fractional coordinates and back.

flake.substrate.calc_matrices_triangle(R)

Metric matrices for a triangular lattice of spacing R.

Nearest-neighbour direction is along x (consistent with the cluster creation conventions in flake.cluster).

flake.substrate.gaussian(x, mu, sigma)

Unnormalised Gaussian: exp(-(x-mu)^2 / (2*sigma^2)).

Kept public because tool_reciprocal_space.py uses it directly. Not normalised by design – it represents a well shape, not a PDF.

flake.substrate.get_ks(R, n, c_n, alpha_n)

Wave vectors for an n-fold symmetric sinusoidal substrate.

Constructs \(n\) plane waves whose interference gives a potential with the desired lattice symmetry. Coefficients from Vanossi, Manini, Tosatti, PNAS 109, 16429 (2012). The \(l\)-th wave vector is:

\[\begin{split}\mathbf{k}_l = \frac{c_n \pi}{R} \begin{pmatrix} \cos\!\left(\tfrac{2\pi l}{n} + \alpha_n\right) \\ \sin\!\left(\tfrac{2\pi l}{n} + \alpha_n\right) \end{pmatrix}, \quad l = 0, \ldots, n-1\end{split}\]

Preset recipes:

Symmetry

n

c_n

alpha_n (deg)

Lines

2

1

0

Triangular

3

4/3

0

Square

4

\(\sqrt{2}\)

45

Quasicrystal 5

5

2

0

Quasicrystal 6

6

\(4/\sqrt{3}\)

−30
Args:

R: float – lattice spacing. n: int – number of plane waves (controls symmetry). c_n: float – amplitude pre-factor; sets \(|\mathbf{k}|\). alpha_n: float – global rotation of the wave-vector star [radians].

Returns:

ks: (n, 2) float64 ndarray of wave vectors.

flake.substrate.particle_en_gaussian(pos, pos_torque, basis, a, b, sigma, epsilon, u, u_inv)

Per-particle energy, force, and torque for a Gaussian-well substrate.

See _particle_en_gaussian_core for the potential definition.

Args:

pos: (N, 2) ndarray – particle positions. pos_torque: (2,) ndarray – torque reference point (usually CM). basis: list of (2,) arrays – substrate basis site positions. a: float – inner cutoff (bulk region ends here). b: float – outer cutoff (tail region ends here). sigma: float – Gaussian width. epsilon: float – well depth (positive = attractive). u, u_inv: (2,2) ndarrays – metric matrices from calc_matrices_bvect.

Returns:

en: (N,) ndarray – potential energy per particle. F: (N, 2) ndarray – force on each particle. tau: (N,) ndarray – torque contribution per particle.

flake.substrate.particle_en_sin(pos, pos_torque, basis, ks, epsilon)

Per-particle energy, force, and torque for a plane-wave substrate.

Args:

pos: (N, 2) ndarray – particle positions. pos_torque: (2,) ndarray – torque reference point. basis: list of (2,) arrays – substrate basis site positions. ks: (n, 2) ndarray – wave vectors (from get_ks). epsilon: float – well depth.

Returns:

en, F, tau: same shapes as particle_en_gaussian.

flake.substrate.particle_en_tanh(pos, pos_torque, basis, a, b, ww, epsilon, u, u_inv)

Per-particle energy, force, and torque for a tanh-well substrate.

Args:

pos: (N, 2) ndarray – particle positions. pos_torque: (2,) ndarray – torque reference point. basis: list of (2,) arrays – substrate basis site positions. a: float – flat-bottom radius. b: float – outer cutoff radius. ww: float – wall shape parameter. epsilon: float – well depth. u, u_inv: (2, 2) ndarrays – metric matrices.

Returns:

en, F, tau: same as particle_en_gaussian.

flake.substrate.substrate_from_params(params)

Build substrate energy functions from a parameter dictionary.

This is the standard entry point for production use. The dictionary is typically loaded from a JSON input file.

Supported well shapes and required parameter keys:

  • 'gaussian': epsilon, sigma, a, b, b1, b2, sub_basis

  • 'tanh': epsilon, wd (alias ww)``, a, b, b1, b2, sub_basis``

  • 'sin': epsilon, ks, sub_basis

  • 'flat': no keys required

For 'gaussian' and 'tanh', b1 and b2 are the primitive lattice vectors used to build the metric matrices (via calc_matrices_bvect).

For 'sin', ks should be a pre-computed (n, 2) array (e.g. from get_ks). No metric matrices are needed.

For 'flat', the substrate is identically zero — useful for testing the integrator against analytic predictions (FDT, diffusion coefficient) without substrate interference.

All array parameters (basis, u, u_inv, ks) are converted to float64 ONCE here and captured in the returned closures. Callers pass en_inputs=[] and call en_func(pos, pos_cm) – the closure holds everything else. This eliminates per-step np.array conversions in run_md.

Args:
params: dict with at least ‘well_shape’ and shape-specific keys

listed above. ‘epsilon’ and ‘sub_basis’ are not required for ‘flat’.

Returns:
pen_func: particle-level closure (pos, pos_cm) ->

(en, F, tau) arrays (per particle).

en_func: total-energy closure (pos, pos_cm) ->

(float, (2,) float64, float).

en_inputs: always [] – all parameters live in the closures.

Raises:

NotImplementedError: if well_shape is not recognised.