Computational Physics

8.6 Molecular systems 237

Enforcing rigidity via constraints

Another method for treating rigid molecules is by imposing holonomic constraints,
i.e. constraints which depend only on positions and not on the velocities, through
an extension of the Lagrangian. The Lagrangian of the system without constraints
reads

L 0 (R,R ̇)=

∫t 1

t 0

dt

[∑

i

mi

2

r ̇^2 i−

1

2

∑

i
=j

U(ri−rj)

]

. (8.115)

A constraint is introduced as usual through a Lagrange multiplierλ[50]. As the
constraint under consideration should hold for all times,λis a function oft. A simple
example of a constraint is the following: particles 1 and 2 have a fixed separationd
for all times (this could be the separation of the two atoms in a nitrogen molecule).
Such a constraint on the separation is calledbond constraint– it can formally be
written as
σ[R(t)]=[r 1 (t)−r 2 (t)]^2 −d^2 =0. (8.116)

The Lagrangian for the system with this constraint reads

L(R,R ̇)=L 0 (R,R ̇)−

∫t 1

t 0

dtλ(t){[r 1 (t)−r 2 (t)]^2 −d^2 }. (8.117)

The integral over time is needed because the constraint holds for all times between
t 0 andt 1. The equations of motion are the Euler–Lagrange equations for this
Lagrangian. These equations will depend on the Lagrange parameters,λ, whose val-
ues are determined by the requirement that the solution must satisfy the constraint.
A slightly more complicated example is the trimer molecule CS 2 [ 51 ].The linear
geometry of this molecule is in principle imposed automatically by the correct bond
constraints between the three pairs of atoms. However, the motion of this molecule
is described by five positional degrees of freedom: two to define the orientation of
the molecule and three for the centre of mass position. The three atoms without
constraints have nine degrees of freedom and if three of these are eliminated using
the bond constraints, we are left with six degrees of freedom instead of the five
required. Therefore one redundant degree of freedom is included in this procedure,
which is obviously inefficient. A better procedure is therefore to fix only the distance
between the two sulphur atoms:

|rS( 1 )−rS( 2 )|^2 =d^2 (8.118)

and to fix the position of the C-atom by a linear vector constraint:

(rS( 1 )+rS( 2 ))/ 2 −rC=0, (8.119)

adding up to the four constraints required.
For a molecule, in general a number of atoms forming a ‘backbone’ set is identi-
fied and these are fixed by bond constraints (the two sulphur atoms in our example)