Skip to main content
15 votes
Accepted

Why is leapfrog integration symplectic and RK4 not, if the latter is more accurate?

TL;DR: It depends on what kind of accuracy you need. Energy conservation does not automatically equal accuracy. Suppose, you want to simulate the solar system, and you are using a solver that – to ...
Wrzlprmft's user avatar
  • 2,022
7 votes
Accepted

The velocity Verlet method and variable time steps

The second formula is just the velocity verlet, and it's correct but if you adapt time steps then it's not symplectic. In a separate answer I describe in quite detail that symplecticness is a global ...
Chris Rackauckas's user avatar
6 votes

What good are hard-sphere event-driven molecular dynamics simulations in the face of chaos?

Yes, it is possible to show that the statistical behavior of the approximate system will reach that of the "exact" system. (This is true even though hard-sphere dynamics do not accurately describe ...
aeismail's user avatar
  • 3,513
6 votes
Accepted

MD Simulation: Reference for the Neighbor's List Method

I'd recommend "The Art of Molecular Dynamics Simulation" by D. C. Rapaport. The code samples are written in C. I'm not a huge fan of the programming style of the book, but at least it's not ...
lr1985's user avatar
  • 677
6 votes

Starting configuration for Molecular Dynamics

Usually one needs to employ periodic boundary conditions (at least in the horizontal directions). Any atoms which fly outside of the box will be mapped to the opposite side. This also has to be ...
MPIchael's user avatar
  • 2,985
5 votes

Fortran code for Ewald summation

The National Institute of Standards and Technology (NIST) proposes different inputs and outputs to test your implementation. Have a look there: https://www.nist.gov/mml/csd/chemical-informatics-...
max's user avatar
  • 207
5 votes
Accepted

Understanding leapfrog integration algorithm

In the second code the full time stepping is given by three lines in main() ...
user9794's user avatar
  • 465
4 votes

How to integrate numerically Nosé Hoover equation?

You can find explicit pseudo-code for Nosé-Hoover Chains in Section 4 of: Martyna, G. J., Tuckerman, M. E., Tobias, D. J., & Klein, M. L. (1996). Explicit reversible integrators for extended ...
Juan M. Bello-Rivas's user avatar
4 votes

Calculation of Mean Square Displacement for Brownian dynamics system with Lennard Jones interactions in python3

This will not an answer to your problem, more an excessive comment and few things you might consider, when writing such code even for self educational purpose. Constants You asked whether your ...
Bort's user avatar
  • 1,285
4 votes

Calculate forces on atoms from potential energy of system and position of atoms

Generally speaking, the force acting on particle $i$ is just $-\nabla E(\vec{r}_i)$ (note the minus sign), where $\nabla$ is the gradient operator and $\vec{r}_i$ is the position of particle $i$. ...
lr1985's user avatar
  • 677
4 votes

Understanding leapfrog integration algorithm

The first code is wrong if you have multiple particles that interact. Due to the structure of the loop, the interaction forces with particles at the start of the list are computed with the updated ...
Lutz Lehmann's user avatar
  • 6,109
3 votes

What good are hard-sphere event-driven molecular dynamics simulations in the face of chaos?

Something that is at least very close to what your are asking for is the Shadowing Lemma, which roughly states that near any numerical solution of a dynamical system, you will find a real solution. I ...
Wrzlprmft's user avatar
  • 2,022
3 votes

How to integrate numerically Nosé Hoover equation?

tl;dr The Nose-Hoover equations are normally defined by their Hamiltonian. The two partial derivatives of the Hamiltonian define a partitioned set of ODEs. From there the partitioned ODEs are either ...
Chris Rackauckas's user avatar
3 votes
Accepted

How to compute forces in multi-particle MD

I will only be developing the computation of forces on a given atom (part 1 of BalazsToth answer), for this point is actually a bit tricky and I think it deserves to be a bit more explained. You can ...
G.Clavier's user avatar
  • 323
3 votes

How to compute forces in multi-particle MD

I have the following remarks: 1. Force calculation The equation of motion of the $i^\text{th}$ atom with mass $m_i$ reads as $\textbf{a}_i=\frac{d\textbf{v}_i}{dt}=\frac{1}{m_i}\cdot f(\vert\textbf{...
BalazsToth's user avatar
3 votes
Accepted

Which are the right configurations in the Markov chain of a Hamiltonian Monte Carlo algorithm?

The algorithm for perfoming a single HMC step is as follows: Input: Some initial configuration $\vec{y}_i$ and momentum $\vec{p}_i$. Output: Next configuration $\vec{y}_{i+1}$ and momentum $\vec{p}_{i+...
cos_theta's user avatar
  • 451
3 votes
Accepted

How error scales with numerical precision in molecular dynamics?

In general in a correctly implemented fixed-step ODE solver method you have 3 sources for numerical errors: the theoretical method truncation error, the floating point error from evaluating the ODE ...
Lutz Lehmann's user avatar
  • 6,109
3 votes

How can I write a unit test for this routine?

Divide and conquer refactor your code, so that you have smaller methods with well defined purposes. Give these methods short descriptive names of what they are doing. You make the code more testable ...
MPIchael's user avatar
  • 2,985
3 votes
Accepted

Computing LJ force from LJ potential, or not?

Which one is computationally more accutate [sic] and why? The first one is more accurate. It uses the exact derivative of the LJ potential. By contrast, the second one not only uses a numerical ...
jjramsey's user avatar
  • 198
3 votes

Starting configuration for Molecular Dynamics

The biggest issue in your subsequent simulation won't be the periodicity of the box -- it will be trying to equilibrate your system down from what looks like a tangle of lipids into an orderly ...
Shern Ren Tee's user avatar
2 votes

Fortran code for Ewald summation

I don't think your formula is correct for the force. You can rewrite your energy expression in an easier form, and it makes taking the derivative easier as well. The way you have it is \begin{equation}...
Wesley's user avatar
  • 121
2 votes

Difficulty with possibly the simplest MD simulation

Two thoughts A) A certain degree of energy dissipation in MD simulations is considered normal and can be counterbalanced with so called thermostat algorithms e.g. this. In your particular case I can ...
MPIchael's user avatar
  • 2,985
2 votes

Oil/Water interface simulation using GROMACS

A simple solution to begin with would be to create two systems: one with decane in a box of width 1.295nm and another one with water in a box of width 2nm (depth and height should be as in the picture ...
Juan M. Bello-Rivas's user avatar
2 votes
Accepted

Can LINCS algorithm be used for colliding molecules?

Constrained simulation methods such as LINCS are there to keep certain degrees of freedom fixed throughout the simulation, so you can use them to make each of your molecules behave like a rigid body. ...
Juan M. Bello-Rivas's user avatar
2 votes

Why are Hamiltonian dynamics used in MCMC?

I'm a little late to reply but I found this review to be really informative. One particularly nice feature of Hamiltonian flows is that they preserve volume in phase space. If we take an infinitesimal ...
Daniel Shapero's user avatar
2 votes

Time Reversibility of Velocity Verlet Algorithm

I can share my perspective on this topic and can try to prove that the Leapfrog (velocity Verlet) is time-reversible, according to an appropriate definition of time reversibility (very nice property ...
Futurologist's user avatar
2 votes

Time Reversibility of Velocity Verlet Algorithm

I think for time reversibility we need to see if we can predict velocity and position at time $t$ if the values at $t+\Delta t$ is known. Verlet Algorithm The integrator equation is given as: $r_{t+\...
user198811's user avatar
2 votes

Who can explain the minimum image convention in molecular dynamic simulations?

The question as stated does make sense in the context of molecular dynamics simulations but is very broad. I recommend that you should check the manual of the MD package your are using about the image ...
Bort's user avatar
  • 1,285

Only top scored, non community-wiki answers of a minimum length are eligible