Calculating atomic forces from total force and torque on molecule in molecular dynamics (MD)

I guess this problem applies to any body that is composed of a small number of particles, but for the sake of simplicity and concreteness, let's consider a classical model of a molecule consisting of three atoms (e.g. H2O). In brief, I'm trying to calculate the forces on each atom (atomic forces, $\vec{F}_i$) from the total force ($\vec{F}_{CM}$) and torque ($\vec{\tau}$) on the respective molecule. Specifically, I'm looking for methods that would work in a computer simulation in the following cases:

1. The molecule is rigid, i.e. fixed distances between the atoms, but the molecule can rotate.
2. The molecule is flexible, i.e. the bond lengths are allowed to vary, but there is an "internal potential" that is a function of the bond lengths that causes forces on the atoms.

I think I should be able to find the 3*3=9 components of the atomic forces from the 3+3=6 components of the total force and torque on the molecule by solving a linear system of 9 equations. However, I can only think of 6 equations:

• The total molecular force equals the sum of the forces on the three molecules, i.e. $\vec{F}_{CM} = \vec{F}_1 + \vec{F}_2 + \vec{F}_3$ (3 equations).
• The total torque on the molecule is given by $\vec{\tau} = \vec{r}_1 \times \vec{F}_1 + \vec{r}_2 \times \vec{F}_2 + \vec{r}_3 \times \vec{F}_3$ where the $\vec{r}_i$ are the positions of the atoms (3 equations).

Apparently, I need three more constraints.

In the case of rigid molecules, the constraints are clearly that the distances between atoms in a molecule are fixed. If I was doing Lagrangian mechanics, I would just choose the generalized coordinates accordingly, or use Lagrange multipliers. What is the preferred way to do this in a computer simulation? I've tried calculating the second derivative of the distance between each pair of atoms in a molecule and letting the resulting expression equal zero. That way I get three equations involving the particles' positions, velocities, and the atomic forces. This seems like a plausible approach but I'm not sure if it's correct.

I have no idea how to approach the case of flexible molecules with "internal potential".

Additional details: I have a computer program that calculates an analytic potential for a certain interaction between triatomic molecules, and also the total force on the molecules using an analytic gradient of the potential. The potential is based on a multipole expansion of the electrostatic potential. The internal potential and its derivatives are also analytic.

I realize this question is probably vague, but I'm unsure which details are needed.

• Might Computational Science be better suited for this question? Commented Aug 11, 2015 at 16:01
• Try to look up Computer Simulation of liquids by Allen and Tildesley ! There are algorithms concerned with what you're asking !! Commented Aug 11, 2015 at 16:26