Research Interests

Sub-Cellular Elements in Multiscale Thrombosis simulation.

Collaboration with Mark Alber, Mathematics Dept. at Notre Dame and Elliot Rosen at IU School of Medicine.

I have recently started modeling multi-cellular systems using sub-cellular elements (SCEM). In contrast to the Cellular Potts models (CPM), a dynamical system is developed that defines cellular shape and elasticity, by emulating the cytoskeleton, and uses the Lennard-Jones potential, to mimic cellular adhesion, giving more realistic simulations. The current methods use harmonic restraints between the sub-cellular elements (SCE) to give the correct shape and elasticity, extensions of the model to cell division allows the restraint to be removed beyond a given distance.

SCE in blood vessel Flow pressure

Latest work is directed to coupling flow within the blood vessel, from solving the incompressible Navier-Stokes equations, with the SCEM Langevin equation.

Normal Mode Partitioning of Langevin Dynamics for Biomolecules.

Collaboration between the Pande group (Folding at Home)
at Stanford and LCLS group at Notre Dame.


Currently I am investigating the use of Normal Mode Analysis (NMA) techniques to facilitate the acceleration of numerical methods. The aim is to approximate the kinetics or thermodynamics of a biomolecule by a reduced model based on a normal mode decomposition of the dynamical space. Our basis set uses the eigenvectors of a mass re-weighted Hessian matrix calculated with a biomolecular force field. This particular choice has the advantage of an ordering according to the eigenvalues, which has a physical meaning (square-root of the mode frequency). Low frequency eigenvalues correspond to more collective motions, whereas the highest frequency eigenvalues are the limiting factor for the stability of the integrator. The higher frequency modes are overdamped and relaxed near to their energy minimum while respecting the subspace of low frequency dynamical modes. The first paper 'Normal Mode Partitioning of Langevin Dynamics for Biomolecules' authored by Christopher R. Sweet, Paula Petrone, Vijay S. Pande and Jesús A. Izaguirre was published in the Journal of Chemical Physics in March 2008.

An interactive website,, illustrating Normal Modes in Molecular Dynamics can be found below. Click the WW-Domain picture to see Protein folding movie. WW Domain

Supplementary information, code and scripts can be found here.

Earlier work can be found here.

Separable Shadow Hamiltonian Hybrid Monte Carlo (S2HMC) method.

Collaboration between LCLS group at Notre Dame and Robert D. Skeel at Purdue University.


Hybrid Monte Carlo (HMC) is a rigorous sampling method that uses molecular dynamics as a global Monte Carlo move, but the acceptance rate of HMC decays exponentially with system size. The Shadow Hybrid Monte Carlo (SHMC) was previously introduced to overcome this performance degradation by sampling instead from the shadow Hamiltonian defined for MD when using a symplectic integrator. SHMC's performance is limited by the need to generate momenta for the MD step from a non-separable shadow Hamiltonian. To overcome this we have introduced the Separable Shadow Hamiltonian Hybrid Monte Carlo (S2HMC) method, based on a formulation of the leapfrog/Verlet integrator that corresponds to a more separable shadow Hamiltonian, which allows efficient generation of momenta. S2HMC gives the acceptance rate of a 4th order integrator at the cost of a 2nd order integrator. The first paper was published in the Journal of Chemical Physics, Fall 2009.


 Supplementary information, code and scripts can be found here.

Statistical Mechanics/Thermostatting

Collaboration with Ben Leimkuhler at Edinborough University and Eric Barth at Kalamazoo College.

During the final two years I studied the thermostatting of dynamical systems in order to generate simulations which sample from the Gibbs, or canonical, ensemble. Molecular dynamics trajectories that sample from this distribution can be generated by introducing a modified Hamiltonian with additional degrees of freedom as described by Nose. To achieve the ergodicity  required for canonical sampling, a number of techniques have been proposed based on incorporating additional thermostatting variables, such as Nose-Hoover chains. For Nosé dynamics, it is often stated that the system is driven to equilibrium through a resonant interaction between the self-oscillation frequency of the thermostat variable and a natural frequency of the underlying system. By pioneering the introduction of multiple thermostat Hamiltonian formulations, which are not restricted to chains, it has been possible to clarify this perspective, using harmonic models, and exhibit practical deficiencies of the standard Nose-chain approach. As a consequence of this it has been possible to propose a new powerful "recursive thermostatting" procedure which obtains canonical sampling without the stability problems encountered with chains. This method also has the advantage that the choice of Nose mass is essentially independent of the system to be thermostatted. In addition a Hamiltonian chains method, Nosé-Poincare chains, has been proposed where these methods are appropriate. This research has yielded two papers, "The Canonical Ensemble via Symplectic Integrators using Nose and Nosé-Poincaré chains" which has been accepted for publication by the Journal of Chemical Physics, and the more recent "A Hamiltonian Formulation for Recursive Multiple Thermostats in a Common Timescale" which has been submitted for publication. As an aid to understanding these methods I have created a website where the methods can be compared which can be accessed via the "Recursive Thermostatting simulation" link to the left. This site uses a cgi script to integrate the simulation for the specified time and return images of the resulting thermostat variable phase-space, distributions and averaged local kinetic energy.

Recursive thermostat

This work continues, we are currently looking at the effect of both statistical and dynamical thermostats on both the velocity auto-correlation function and sampling for pentane in liquid and gas form.


I studied for my Ph.D under the supervision of Professor Benedict Leimkuhler in the general area of Geometric Integrators. During my research I developed an interest in software which provides a visual interpretation of simulation data and have written several programs using Trolltech Qt (as used for the Linux KDE interface) and OpenGL. I have included movies of some of these simulations on the left hand side "Links" section, and source code for a basic gravitational model.

Higher order variable step-size methods

During the first year of my Ph.D. studies I looked at N-body systems, such as the Solar System, with particular reference to variable step-size integrators. For fixed step-size integrators it is possible to compose steps of different size to increase the order of the method through the cancelation of error terms. When extending this to variable step-size methods the expected increase in order is often not obtained. Analysis of this problem led to a family of variable step-size higher order methods detailed in the paper "Higher Order Symmetric Variable Stepsize Methods" which has been submitted for publication. A Solar System model has been integrated over one billion years using these methods. As an aid to visualization I wrote an OpenGL program to generate images from the simulation data for this model, the Solar System movie/Gif from the left hand menu was generated using this software. I also wrote a brief history of Solar System simulations in an article for the Wrangler, a maths newsletter from the University of Leicester which can be accessed from the left hand menu, which included some of these graphics.