mod bias

module bias

trait Bias<S>: cost-augmenting operator for enhanced-sampling methods inside SA. The standard well-tempered metadynamics (Barducci/Bussi/Parrinello 2008) is the canonical impl; the trait is open so SGOOP-style or VES bias kernels can implement the same surface.

In the typed algebra, Bias augments the cost as F_eff(x) = F(x) + V(s(x)). Sampler sees F_eff via a wrapped objective, while the bias remains stateful: deposit mutates internal state and potential reads it. This keeps biasing in the objective transform instead of folding it into the sampler implementation.

Traits

trait BasinMetric

How far apart two descriptors are, for deciding whether they are one basin.

Separate from Fingerprint because the descriptor and the notion of sameness are separate choices, and getting the second wrong is what a merge radius in descriptor space does. A sorted distance spectrum compared by Euclidean distance needs a threshold with no physical meaning, found empirically and not transferable between system sizes; the same descriptor compared by an optimal-permutation shape distance needs a length, which transfers.

Functions

fn distance(&self, a: ArrayView1<f64>, b: ArrayView1<f64>) -> f64

Distance between two descriptors. Must be symmetric and vanish on identical inputs.

fn distance_bounded(&self, a: ArrayView1<f64>, b: ArrayView1<f64>, _bound: f64) -> f64

Distance, allowed to stop once it is certain to exceed bound.

A caller asking whether a centre is within the merge radius does not need the distance to a centre that is not.

fn key(&self, _v: ArrayView1<f64>) -> Option<f64>

A scalar obeying |key(a) - key(b)| <= distance(a, b), when the metric has one.

The reverse triangle inequality turns a scan into a subtraction for every centre that cannot possibly be within the radius, which on a 98-point cluster is nearly all of them: the descriptor there has 4753 entries, so one Euclidean comparison costs about as much as a Lennard-Jones gradient, and a chain that has opened twenty thousand basins pays that per centre on every miss.

None for metrics with no such quantity, notably shape distances that minimise over permutations and rotations, where the scan stays exact and exhaustive.

Implemented for

impl BasinMetric for IraMetric
trait Bias

Cost-augmenting bias on a low-dimensional collective variable s = phi(x). Implementors maintain internal state that is updated by deposit and read by potential.

Functions

fn cv(&self, x: ArrayView1<f64>) -> Array1<f64>

Map a position to its CV value.

fn deposit(&mut self, s: ArrayView1<f64>, temp: f64)

Deposit a Gaussian at CV value s at temperature temp. Mutating; the well-tempered weight depends on the current potential(s) and temp.

fn potential(&self, s: ArrayView1<f64>) -> f64

Bias potential at the CV value s.

fn reweight(&self, s: ArrayView1<f64>, temp: f64) -> f64

Reweighting factor exp(+V(s) / T) for post-hoc unbiasing of observables computed under the biased ensemble.

Implemented for

impl<F: Fingerprint> Bias for SpectralBias<F>
trait Fingerprint

Maps a state to a vector that compares equal for states in the same basin.

Functions

fn describe(&self, x: ArrayView1<f64>) -> Array1<f64>

Descriptor of x. Two states in the same basin must map to vectors within the bias’s merge radius, and states in different basins must not.

Implemented for

impl Fingerprint for ClusterFingerprint
impl crate::bias::Fingerprint for CanonicalOrder
impl Fingerprint for IraPassthrough

Structs and Unions

struct AdaptiveHeight

Deposit height estimated from the escape gaps a chain actually sees.

A per-basin bias only moves a chain once the bias in the occupied basin exceeds the energy cost of leaving it, so the deposit height and that cost have to be commensurate. Setting the height by hand does not survive contact with a real landscape: on a 75-point Lennard-Jones cluster the cheapest escape from the structure basin hopping plateaus at costs 0.0906 and the tenth percentile costs 0.1831, against a hand-set default of 0.25, so one deposit clears both and the chain abandons a basin after a single revisit rather than filling what it is sitting in. The same default at a different size, or on a different potential, is wrong in whichever direction that landscape happens to run.

The gaps are observable: every rejected uphill proposal is a sample from the escape distribution. This tracks a low quantile of those samples and sets the height so that a chosen number of revisits clears it, which makes the bias self-scaling in the units of the landscape rather than of the author.

The estimate is a P-square style quantile tracker rather than a stored history, so the cost does not grow with the run.

quantile: f64

Quantile of the escape gap the bias aims to clear, in [0, 1).

revisits: f64

Revisits a basin should take before that gap is cleared.

Implementations

impl AdaptiveHeight

Functions

fn gap_estimate(&self) -> f64

Current estimate of the targeted escape-gap quantile.

fn height(&self) -> f64

Height to deposit per revisit.

fn new(quantile: f64, revisits: f64, prior: f64) -> Self

Tracker aiming at quantile, clearing it in revisits deposits.

fn observe(&mut self, gap: f64)

Records one observed escape gap, which must be positive.

A non-positive gap is a downhill move and says nothing about the cost of leaving, so it is ignored rather than dragging the estimate to zero.

fn samples(&self) -> u64

Samples recorded.

struct BasinBias<F: Fingerprint>

Well-tempered bias keyed on discrete basin identity rather than on a collective variable.

A grid bias has to be told which projection to watch. That works when the competing structures separate along the chosen axis and fails silently when they do not: on the 38-atom Lennard-Jones cluster the close-packed and icosahedral funnels differ by 0.19 in the fourth Steinhardt parameter, while at 75 atoms the decahedral and icosahedral minima differ by 0.023, which is narrower than a sensible deposition width. The bias then fills a region that contains both competitors and the search never leaves.

Keying on identity removes the choice. Two states are the same basin when their fingerprints lie within merge_radius, so there is no axis to be blind along. Revisiting a basin raises its bias, which is the superbasin escape acceleration of Chatterjee and Voter (J Chem Phys 132, 194101, 2010) with the Barducci well-tempered weight on top.

The fingerprint must be invariant to whatever the objective is invariant under, or the same physical state registers as many basins. SortedPairs is the default for point sets: invariant to permutation, translation and rotation, and free of external dependencies.

Implementations

impl<F: Fingerprint> BasinBias<F>

Functions

fn deepest(&self) -> f64

Deepest accumulated bias over all basins.

fn height(&self) -> f64

Current deposit height.

fn index(&self) -> &BasinIndex<F>

The identity half, for a mechanism that keys on basins without depositing.

fn merge_radius(&self) -> f64

Current merge radius.

fn missing_mass(&self) -> f64

Good-Turing missing mass: the share of basins seen exactly once, which estimates the probability that the next visit opens a new one.

fn n_basins(&self) -> usize

Number of distinct basins registered so far.

fn new(fingerprint: F, merge_radius: f64, w0: f64, gamma: f64) -> Self

Requires gamma > 1`, w0 > 0 and ``merge_radius > 0`.

fn set_height(&mut self, w0: f64)

Sets the height deposited per revisit.

Exposed because the right value is a property of the landscape rather than of this type: it has to be commensurate with the energy cost of leaving a basin, and a height above that cost empties a basin on a single revisit instead of filling it.

fn set_merge_radius(&mut self, radius: f64)

Sets the distance below which two descriptors are one basin.

Exposed because this is a schedule rather than a setting. Lee, Lee and Scheraga show the threshold plays the role of a temperature and is annealed from wide to narrow, and their method solves the cluster sizes a fixed threshold does not. See crate::diversity.

fn with_metric(mut self, metric: Box<dyn BasinMetric>) -> Self

Replaces how sameness is measured.

The descriptor and the notion of sameness are separate choices. Keying on a sorted distance spectrum compared by Euclidean distance needs a threshold with no physical meaning: measured on this crate’s cluster driver, LJ38 solves 1 seed in 8 that way. Comparing the same structures by optimal-permutation shape distance makes the threshold a length.

Traits implemented

impl<F: Fingerprint> Bias for BasinBias<F>
struct BasinIndex<F: Fingerprint>

Which basin a state is in, with no potential attached.

The identity half of BasinBias, split out because more than one mechanism needs to ask “have I been here before” and only one of them answers by depositing. History-conditioned escape uses the same rule to decide how hard to push next, and reading that off a bias would tie the two together: under replica exchange each rung owns its own bias, so the basin numbering of one rung means nothing in another, while a controller following one chain needs a numbering that outlives the swap.

Implementations

impl<F: Fingerprint> BasinIndex<F>

Functions

fn basin_of(&mut self, x: ArrayView1<f64>) -> usize

Index of the basin holding x, opening one if it is new, and counting the visit.

fn bump(&mut self, i: usize)

Records another visit to i.

fn describe(&self, x: ArrayView1<f64>) -> Array1<f64>

Descriptor of a state.

fn lookup(&self, d: ArrayView1<f64>) -> Option<usize>

Index of the nearest registered basin within merge_radius.

Most recently added first, returning the first centre inside the radius rather than scanning for the nearest.

Both parts matter once the metric costs anything. A chain revisits the basin it is already in far more often than any other, measured at roughly nineteen proposals in twenty near a deep minimum, so the recent end is where the answer almost always is. And a merge radius asks whether any centre is within it, not which is closest, so the scan can stop at the first hit.

With Euclidean distance the exhaustive scan was free and the distinction did not show. With a shape distance at milliseconds per comparison it is the difference between one call and one per basin, and the exhaustive version did not finish an LJ38 run.

fn merge_radius(&self) -> f64

Current merge radius.

fn n_basins(&self) -> usize

Number of distinct basins registered so far.

fn new(fingerprint: F, merge_radius: f64) -> Self

Index over fingerprint, calling two states the same within merge_radius.

fn push(&mut self, d: Array1<f64>) -> usize

Registers a descriptor as a new basin and returns its index.

fn set_merge_radius(&mut self, radius: f64)

Sets the distance below which two descriptors are one basin.

fn visits(&self, i: usize) -> u64

Times basin i has been recorded.

fn with_metric(mut self, metric: Box<dyn BasinMetric>) -> Self

Replaces the distance used to decide sameness.

struct EuclideanMetric

Ordinary Euclidean distance between descriptors.

Traits implemented

impl BasinMetric for EuclideanMetric
struct SiteEnergies

Sorted per-point pair energies of a flattened (n, 3) point set.

E(i) = sum_{j != i} 4 [ (1/r_ij)^12 - (1/r_ij)^6 ], sorted. Invariant to permutation by the sort, and to rigid motions because it is built from distances, so it has the same symmetries as SortedPairs and the same cost, one pass over the pairs.

What it adds is that it is energetic. A sorted distance spectrum says how far apart the points are; this says how well each one is bound, so two structures with similar distance statistics and different coordination separate. That matters because basin identity is what the merge radius is measuring, and the radius is sharply sensitive: at 75 points, 0.7 in the distance spectrum gives 13 seeds in 24 while 0.95 gives 0 in 8. A descriptor whose distances between genuinely different structures are better separated is what would widen that band.

A richer alternative was measured and rejected on cost. The eigenvalue spectrum of the pairwise distance matrix is a stronger invariant, and it needs an order n^3 decomposition per hop against the order n^2 of the thirty energy evaluations a hop already spends, so it would cut the hop count roughly threefold. Hops are the scarce resource in this search.

n_points: usize

Points per state; the state length must be 3 * n_points.

Traits implemented

impl Fingerprint for SiteEnergies
struct SortedPairs

Sorted pairwise distances of a flattened (n, 3) point set.

Invariant to permutation of the points and to rigid motions, which are the symmetries of a cluster energy. Sorting is what supplies permutation invariance and is also why the descriptor is cheap.

n_points: usize

Points per state; the state length must be 3 * n_points.

Traits implemented

impl Fingerprint for SortedPairs
struct WellTemperedBias

Well-tempered metadynamics (Barducci, Bussi, Parrinello 2008).

State: a 2D grid of bin values storing the deposited bias. The deposit rule is w_k = w_0 * exp(-V_{k-1}(s_k) / ((gamma - 1) * T)), so the asymptotic bias is V_inf(s) = -((gamma - 1)/gamma) * F(s).

Asymmetric box bounds are supported on each CV axis. The CV map is a fixed linear projection from R^dim to R^2 configured at construction time – typical use: pass the top-2 TICA components.

projector: Array2<f64>

Linear projection phi: R^dim -> R^2, shape (dim, 2).

mu: Array1<f64>

Mean offset in R^dim subtracted before projection.

low: [f64; 2]

CV-space lower corner.

high: [f64; 2]

CV-space upper corner.

sigma: f64

Gaussian width in CV space.

w0: f64

Initial deposition height.

gamma: f64

Well-tempered bias factor; gamma > 1.

grid_n: usize

Number of bins per CV axis.

v: Array2<f64>

Bias values on the (grid_n, grid_n) grid.

Implementations

impl WellTemperedBias

Functions

fn new(projector: Array2<f64>, mu: Array1<f64>, low: [f64; 2], high: [f64; 2], sigma: f64, w0: f64, gamma: f64, grid_n: usize) -> Self

Constructs with the given linear projector + bias parameters. Requires gamma > 1, sigma > 0, and grid_n >= 2.

Traits implemented

impl Bias for WellTemperedBias