The Analytical Method of Constraint Dynamics

This section describes in detail the analytical method for any forms of holonomic constraints and with an integration algorithm requiring derivatives of forces with respect to time up to any order The treatment is based mainly on the method described by Ryckaert et al. and used earlier by Orban and Ryckaert, but it contains added details and minor corrections. Where fundamental equations are reached, the corresponding equations of Reference 5 (often with slight nocational differences) are noted for comparison. As stated before, the treatment is by definition couched throughout in Cartesian coordinates. 

Consider a system of N interacting particles subject to / general holonomic 

Throughout this chapter, denotes an th-order derivative with respect to time of i j, and i i denotes an th power of i j. For the first and second time derivatives of the coordinates, the dot and parentheses notations are used interchangeably. 

The integration of the equations of motion In this step, the forces of constraint and their time derivatives up to order s nax obtained from step 1, are used as input to the selected integration scheme, to generate the constrained coordinates. Note that the particular choice of numerical integration algorithm in step 2 determines the parameter of step 1. 

We term the approach described here the analytical method to emphasize that the Lagrange multipliers and their derivatives are computed up to 

---

The basic constraint dynamics method employed There are two distinct methods of constraint dynamics the analytical method and the method of undetermined parameters. The reason for the names will become apparent. 

Although this chapter is written as a review of the methods of constraint dynamics, a substantial part of the material is new. In the next section, the analytical method is described in detail in its most general form. The gradual divergence of the constraints and the need for a constraint correction scheme are discussed. Finally, the method of Edberg et al. is discussed in the context of the analytical method, as a special case with = 0 and holonomic bondstretching constraints, together with a constraint correction scheme. 

The analytical method deserves a detailed discussion for at least two major reasons. First, if used in conjunction with some constraint correction scheme," it is important in its own right as a practical method of solution of the constrained dynamics problem. Second, the method of undetermined parameters, central to the subject matter of this chapter, is an outgrowth of the analytical method hence a thorough understanding of the analytical method is an essential prerequisite for understanding the method of undetermined parameters. 

These difficulties have led to a revival of work on internal coordinate approaches, and to date several such techniques have been reported based on methods of rigid-body dynamics [8,19,34-37] and the Lagrange-Hamilton formalism [38-42]. These methods often have little in common in their analytical formulations, but they all may be reasonably referred to as internal coordinate molecular dynamics (ICMD) to underline their main distinction from conventional MD They all consider molecular motion in the space of generalized internal coordinates rather than in the usual Cartesian coordinate space. Their main goal is to compute long-duration macromolecular trajectories with acceptable accuracy but at a lower cost than Cartesian coordinate MD with bond length constraints. This task mrned out to be more complicated than it seemed initially. 

One of the first approaches introduced and used in condensed phases, also in the sense of the historical development, was a biasing technique known as the Blue Moon [55]. This approach was better refined over the years [56] to become more user-friendly and easily implementable in a computer code. For all the details, we refer the reader to the cited original publications. Just to summarize the essential points, let us recall that in the case of first principles dynamical simulation, the method relies on the identification of a reaction coordinate, or order parameter, = (R/) of a given subset of atomic coordinates (Lagrangean variables) R/ able to track the activated process or chemical reaction on which one wants to focus. The simplest example is represented by the distance = R/—Ry between two atoms that are expected to form or break a chemical bond. This analytical function is added to, e.g., a Car-Parrinello Lagrangean as a holonomic constraint. 

---