pyscivr.integrators module
Integrators for classical trajectories in Cartesian and curvilinear coordinates.
See sections II B 2 and appendices J 2 and L of the paper for details.
- Notation:
N: number of atoms.T: number of time steps of the trajectory.K: number of constraints.m: number of DOFs of the curvilinear coordinates
- class pyscivr.integrators.Integrator(calculator, posi, mome, mass, dt, *, i_step=0)[source]
Bases:
ABCBase class for integrators in Cartesian coordinates.
The trajectory returned by
run()andget_trajectory()is a dictionary with the following keys, in addition to the extra keys of each subclass:posi:(T, N, 3)positions.mome:(T, N, 3)momenta.potential_energy:(T,)potential energies.force:(T, N, 3)forces.total_energy:(T,)total energies.dt: Time step.mass:(N,)atomic masses.
- Parameters:
calculator (Calculator) – Calculator whose
energy_forcesaccepts(N, 3)positions and returns(energy, forces).posi (ArrayLike) –
(N, 3)initial positions.mome (ArrayLike) –
(N, 3)initial momenta.mass (ArrayLike) –
(N,)atomic masses.dt (float) – Time step.
i_step (int) – Initial step index.
- class pyscivr.integrators.Verlet(calculator, posi, mome, mass, dt, *, i_step=0)[source]
Bases:
IntegratorVelocity Verlet integrator.
- Parameters:
calculator (Calculator) – Calculator whose
energy_forcesaccepts(N, 3)positions and returns(energy, forces).posi (ArrayLike) –
(N, 3)initial positions.mome (ArrayLike) –
(N, 3)initial momenta.mass (ArrayLike) –
(N,)atomic masses.dt (float) – Time step.
i_step (int) – Initial step index.
- class pyscivr.integrators.ShakeRattle(calculator, posi, mome, mass, dt, coord_const, *, max_iter_shake=50, tol_shake=1e-10, rcond=1e-10, tol_init=1e-08, i_step=0)[source]
Bases:
IntegratorSHAKE/RATTLE integrator with holonomic constraints \(F(q) = 0\).
The integration is performed in mass-weighted coordinates (Eqs. (25)-(31) in the paper). The initial momenta are projected to the tangent space of the constraints.
\[ \begin{align}\begin{aligned}p_{n+\frac{1}{2}} &= p_n - \frac{\Delta t}{2}\nabla V(q_n) + \mathcal{G}(q_n)^T \lambda\\q_{n+1} &= q_n + \Delta tM^{-1}p_{n+\frac{1}{2}}\\p_{n+1} &= p_{n+\frac{1}{2}}-\frac{\Delta t}{2}\nabla V(q_{n+1}) +\mathcal{G}(q_{n+1})^T \mu\end{aligned}\end{align} \]In addition to the keys of
Integrator, the trajectory contains:lam:(T, K)Lagrange multipliers for position constraints.mu:(T, K)Lagrange multipliers for momentum constraints.n_iter_shake:(T,)number of iterations for SHAKE.err_shake:(T,)finalmax|F(q)|of SHAKE.
- Parameters:
calculator (Calculator) – Calculator whose
energy_forcesaccepts(N, 3)positions and returns(energy, forces).posi (ArrayLike) –
(N, 3)initial positions, which must satisfy the constraints withintol_init.mome (ArrayLike) –
(N, 3)initial momenta.mass (ArrayLike) –
(N,)atomic masses.dt (float) – Time step.
coord_const (InternalCoordinates) – InternalCoordinates whose residuals (
get_const) define theKholonomic constraints.max_iter_shake (int) – Maximum number of iterations for SHAKE.
tol_shake (float) – Tolerance for SHAKE convergence.
rcond (float) – Cutoff for pseudo-inverse computation.
tol_init (float) – Tolerance for initial constraint satisfaction.
i_step (int) – Initial step index.
- Raises:
AssertionError – If the initial positions do not satisfy the constraints within
tol_init.
- class pyscivr.integrators.ImplicitVerlet(posi, mome, dt, metric_func, *, potential_func=None, normalize=False, max_iter_newton=50, tol_newton=1e-12, i_step=0)[source]
Bases:
objectImplicit Verlet integrator for \(H(Q, P) = \frac{1}{2} P^T G(Q) P + V(Q)\).
One step of the integrator (Eqs. (J4)-(J6) in the paper) is given by:
\[ \begin{align}\begin{aligned}P_{n+\frac{1}{2}} &= P_n - \frac{\Delta t}{2} \frac{\partial H}{\partial Q}(Q_n, P_{n+\frac{1}{2}})\\Q_{n+1} &= Q_n + \frac{\Delta t}{2} (G(Q_n) + G(Q_{n+1})) P_{n+\frac{1}{2}}\\P_{n+1} &= P_{n+\frac{1}{2}} - \frac{\Delta t}{2} \frac{\partial H}{\partial Q}(Q_{n+1}, P_{n+\frac{1}{2}})\end{aligned}\end{align} \]where the implicit equations are solved by Newton’s method.
The trajectory returned by
run()andget_trajectory()is a dictionary with the following keys:posi:(T, m)positions.mome:(T, m)momenta.mome_half:(T, m)momenta at half step.potential_energy:(T,)potential energies.total_energy:(T,)total energies.n_iter_mome:(T,)number of iterations for the implicit equation of momenta.n_iter_posi:(T,)number of iterations for the implicit equation of positions.err_mome:(T,)final error of the implicit equation of momenta.err_posi:(T,)final error of the implicit equation of positions.dt: Time step.
- Parameters:
posi (ArrayLike) –
(m,)initial positions.mome (ArrayLike) –
(m,)initial momenta.dt (float) – Time step.
metric_func (Callable[[jax.Array], jax.Array]) – Function to compute the metric \(G(Q)\),
(m,)->(m, m).potential_func (Callable[[jax.Array], float] | None) – Function to compute the potential \(V(Q)\),
(m,)-> scalar. IfNone, the potential is assumed to be zero.normalize (bool) – If
True, the positions are normalized to unit norm and the momenta are projected onto the tangent space at every step. (q^T q = 1, p^T q = 0). This is intended for the quaternion dynamics.max_iter_newton (int) – Maximum number of iterations for Newton’s method.
tol_newton (float) – Tolerance for Newton’s method convergence.
i_step (int) – Initial step index.
Note
The jitted functions are cached by
(metric_func, potential_func, dt, normalize, max_iter_newton, tol_newton). Passing the same function objects reuses the compiled functions.Warning
The integrator is vulnerable to the singularity of the metric \(G(Q)\). Check the condition number of \(G(Q)\) before use.
- step()[source]
Advance the state by one step using the implicit Verlet algorithm.
- Return type:
None