Conformational analysis
lemonade computes the standard single-chain observables - radius of gyration, centre of mass, asphericity, end-to-end distance, distance maps and internal scaling - either for the whole trajectory at once or for a single chain.
Whole trajectory, at once
The trajectory-level methods return an array covering every chain in every frame, computed in a handful of vectorised operations (no Python loop over chains or frames). This is the fast path and what you want for ensemble averages:
rg = traj.radius_of_gyration() # (n_frames, n_chains)
com = traj.center_of_mass() # (n_frames, n_chains, n_dim)
asph = traj.asphericity() # (n_frames, n_chains)
ete = traj.end_to_end_distance() # (n_frames, n_chains)
From there, ordinary numpy gives you whatever average you need:
import numpy as np
equil = traj.n_frames // 2 # discard the first half
rg_mean_per_chain = rg[equil:].mean(axis=0) # (n_chains,)
rg_mean = rg[equil:].mean() # scalar ensemble average
rg_by_type = {t: rg[equil:][:, traj.chain_types == t].mean()
for t in np.unique(traj.chain_types)}
Single chain
Navigating to a Polymer gives the same quantities as
scalars, plus per-chain matrices. The scalars are read from the batched arrays, so
they agree exactly with the trajectory-level results.
p = traj[0][3]
p.radius_of_gyration # == traj.radius_of_gyration()[0, 3]
p.end_to_end_distance
p.center_of_mass
dm = p.distance_map() # (L, L) inter-bead Euclidean distances
sep, dist = p.internal_scaling() # mean distance vs sequence separation |i - j|
internal_scaling returns the separations 1 .. L-1 and the mean bead-bead
distance at each - the lattice analogue of the internal-scaling profile used to read
off polymer scaling exponents.
How the numbers are defined
All conformational quantities are computed on whole positions - each chain is first made contiguous across periodic boundaries (a compiled kernel does this for the whole trajectory at once), so a chain that straddles a box face is measured as one connected object rather than being torn in two.
Radius of gyration is \(R_g = \sqrt{\langle |{\bf r}_i - {\bf r}_{cm}|^2 \rangle}\) over the chain’s beads (equivalently the square root of the trace of the gyration tensor), using the exact (floating-point) centre of mass.
Asphericity comes from the gyration-tensor eigenvalues: \(\lambda_3 - \tfrac12(\lambda_1 + \lambda_2)\) in 3D (zero for an isotropic coil), \(\lambda_2 - \lambda_1\) in 2D.
Centre of mass is the mean of the whole positions;
distance_mapandend_to_end_distanceare ordinary Euclidean distances on those contiguous coordinates.
Note
Agreement with PIMMS’s own Rg. For chains that are small compared with the
box, lemonade’s Rg matches PIMMS’s built-in get_polymeric_properties. For a
chain larger than roughly half the box the two intentionally differ: lemonade
measures the whole (contiguous) chain, whereas PIMMS’s on-the-fly analysis uses
the minimum-image convention, which collapses such a chain (a finite-size
artefact PIMMS warns about during the run). If you see a large discrepancy, your
box is small relative to your chains.
Positions
Two views of the coordinates are available at every level:
positions- integer lattice coordinates wrapped into the box (0 <= x < dimensions). Use these for anything grid-based (occupancy, clustering, slab profiles).whole_positions- the chain unwrapped so consecutive beads never jump a boundary (coordinates may fall outside the box). Use these for shape and distance measurements. This is what all the conformational analyses use internally.
p = traj[0][3]
p.straddles_boundary # does this chain cross a periodic face?
p.positions # wrapped
p.whole_positions # contiguous
At the trajectory level, traj.positions and traj.whole_positions() give the
same two views for the entire run as (n_frames, n_atoms, 3) arrays.