jaxdem.utils.dynamical_matrix#

Energy Hessians in sphere and rigid-clump coordinates.

Non-bonded contributions use the configured collider’s candidates and the force model’s pair energy. Bonded contributions differentiate the bonded model’s total energy. Their sum gives the total potential-energy Hessian.

Functions

`bonded_hessian`(state, system)	Return state, system, and the Hessian of the bonded potential energy.
`clump_non_bonded_hessian`(state, system, *[, ...])	Return state, refreshed system, and the dense rigid-clump energy Hessian.
`non_bonded_hessian`(state, system)	Return state, refreshed system, and the dense sphere-coordinate Hessian.
`pair_non_bonded_hessian`(state, system)	Return compact sphere-pair IDs and their potential-energy Hessian blocks.
`zero_mode_mask`(eigenvalues[, rel_gap])	Boolean mask identifying numerically-zero eigenvalues via gap detection.

jaxdem.utils.dynamical_matrix.pair_non_bonded_hessian(state: State, system: System) → tuple[State, System, jax.Array, jax.Array][source]#

Return compact sphere-pair IDs and their potential-energy Hessian blocks.

Parameters:

state (State, System) – One 2D or 3D configuration and its configured collider and force model.
system (State, System) – One 2D or 3D configuration and its configured collider and force model.

Returns:

state (State) – Input state.
system (System) – System carrying any refreshed collider cache. History is not advanced.
pair_ids (jax.Array) – Sorted unique sphere pairs (i, j) with i < j, shape (M, 2). Pairs obey the force model’s search bounds and collider exclusions. There are no padding rows. Zero-force candidates are included.
blocks (jax.Array) – Shape (M, 2*dim, 2*dim), in coordinate order (r_i, r_j). Each block differentiates (m_ij E_ij + m_ji E_ji) / 2, with directed interaction masks m from the collider. For reciprocal interactions this is the Hessian of one pair potential. Blocks can be zero for candidates outside a law’s actual interaction range.

Notes

Collection is a host operation; differentiation runs in compiled batches. An incomplete collider search raises ValueError before differentiation.

jaxdem.utils.dynamical_matrix.non_bonded_hessian(state: State, system: System) → tuple[State, System, jax.Array][source]#

Return state, refreshed system, and the dense sphere-coordinate Hessian.

The matrix has shape (N*dim, N*dim) in flattened sphere-position order. Each unique pair contributes one block with the energy weighting described in pair_non_bonded_hessian(). Pair blocks are evaluated and scattered in bounded batches; output storage is quadratic in N*dim.

Search settings come from the configured collider and force model. Candidate collection runs on the host. Particle dynamics and contact history are unchanged. Incomplete searches raise ValueError.

jaxdem.utils.dynamical_matrix.bonded_hessian(state: State, system: System) → tuple[State, System, jax.Array][source]#

Return state, system, and the Hessian of the bonded potential energy.

The matrix has shape (N*dim, N*dim) in flattened sphere-position order. It is zero when system.bonded_force_model is None. Add it to non_bonded_hessian() to include both energy contributions.

jaxdem.utils.dynamical_matrix.clump_non_bonded_hessian(state: State, system: System, *, rotation_scale: jax.Array | None = None) → tuple[State, System, jax.Array][source]#

Return state, refreshed system, and the dense rigid-clump energy Hessian.

Coordinates for each clump are (delta r_c, omega): translation and an infinitesimal rotation vector (one component in 2D, three in 3D). The matrix includes the second-order rotational displacement and has shape (G*dof, G*dof), with dof = dim + angular_dof and G = max(clump_id) + 1. Rows and columns are indexed by clump ID.

Parameters:

state (State, System) – One configuration and its configured collider and force model.
system (State, System) – One configuration and its configured collider and force model.
rotation_scale (jax.Array, optional) – Positive finite length scales, shape (G,). Rotation rows and columns for clump I are divided by rotation_scale[I] for coordinates (delta r_c, R*omega). Omit for angle coordinates.

Notes

Candidate collection runs on the host; pair differentiation and assembly run in bounded compiled batches. Output storage is quadratic in G*dof. Contact history and particle dynamics are unchanged. Incomplete collider searches raise ValueError. See pair_non_bonded_hessian() for the directed energy weighting and inclusion of zero-force candidates.

Boolean mask identifying numerically-zero eigenvalues via gap detection.

Given a 1-D array of eigenvalues (in any order), we sort by \(|\lambda|\) ascending and look for the largest relative gap \(|\lambda_{k+1}| / |\lambda_k|\). The function flags entries below the first gap that exceeds rel_gap as numerically zero. The function returns the mask aligned with the original eigenvalue ordering, so it can slice eigenvectors directly (e.g. evecs[:, ~mask] for finite modes).

The result does not depend on the problem scale: a hessian with |λ_max| ~ 1 and zero modes at ~ 1e-16 gives the same mask as one with |λ_max| ~ 1e6 and zero modes at ~ 1e-10, because the criterion is the ratio between successive magnitudes, not an absolute threshold.

Parameters:

eigenvalues (ArrayLike) – 1-D array of eigenvalues (e.g. from jax.numpy.linalg.eigvalsh() or jax.numpy.linalg.eigh()). Ordering is not required.
rel_gap (float, optional) – Minimum ratio |λ_{k+1}| / |λ_k| that counts as the zero / finite boundary. Default 1e4 (four orders of magnitude). True zero modes at machine precision and real modes of order \(k \cdot \mathrm{overlap}\) in a jammed spring packing typically sit 10-14 orders of magnitude apart. The threshold is therefore not sensitive.

Returns:

Boolean array of the same shape as eigenvalues. True where the eigenvalue is below the first large relative gap. If the function finds no gap larger than rel_gap (all eigenvalues are comparable in magnitude), it returns all-False.

Return type:

jax.Array