Revisiting and parallelizing SHAKE (Q2485678)

A widely used approach to compute molecular dynamics trajectories with constraints is considered. This approach is called SHAKE. An algorithm is presented for running SHAKE in parallel. The dynamics is given by the Newton's equation \[ M \ddot{X} = - \nabla U - \sum_l \lambda_l \nabla \sigma_l, \] where \(M\) is the mass matrix, \(X\) is the vector of coordinates, \(U\) is the potential, \(\nabla \sigma_l\) are the constraint forces, and \(\lambda_l\) are the Lagrange multipliers which have to be determined from the constraints \(\sigma_l(X)=0\). Within the Verlet algorithm, this equation of motion is represented by updating rules for coordinates and velocities after each integration time step \(\Delta t\). Besides, in each time step one needs to solve iteratively the equations for \(\lambda_l\). The authors propose several modifications to the known SHAKE algorithm to make it suitable for parallel computations. In particular, it includes solution of a linear system of equations for a sparse matrix by conjugate gradient minimization and several implications to reduce the calculations which are difficult to parallelize. The developed parallel algorithm is exact in the sense that it rigorously reproduces the serial one. Numerical examples considered present a good load balancing and are limited only by communication time.