In MD simulations, the majority of time is spent on non-bonded force calculations; therefore, to speedup the MD simulation, these non-bonded calculations have to be accelerated. There are two types of non-bonded interactions; they are short-range Lennard-Jones (LJ) interactions and long-range Coulombic interactions. For the LJ interactions, the complexity of energy and force computations is O(N2). On the other hand, for the Coulombic interactions, an O(N2) method called Ewald Summation  can be used to perform the energy and force computations. For a large value of N, performing the Ewald Summation computations is still very time-consuming. Hence, the SPME algorithm, which scales as N×Log(N), is developed. In the SPME algorithm, the calculation of the Coulombic energy and force is divided into a short-range direct sum and a long-range reciprocal sum. The SPME algorithm allows the calculation of the direct sum to be scaled as N and that of the reciprocal sum to be scaled as N×Log(N). Currently, the SPME reciprocal sum calculation is only implemented in software [8, 9]. Therefore, to shorten the overall simulation time, it would be extremely beneficial to speedup the SPME reciprocal sum calculation using FPGAs.
The Reciprocal Sum Compute Engine (RSCE) is an FPGA design that implements the SPME algorithm to compute the reciprocal space contribution of the Coulombic energy and force. The design of the FPGA aims to provide precise numerical results and maximum speedup against the SPME software implementation . The implemented RSCE, which is realized on a Xilinx XCV-2000 multimedia board, is used with the NAMD2 (Not another Molecular Dynamics) program [4, 10] to carry out several MD simulations to validate its correctness. Although resource limitations in the multimedia board are expected, the thesis also describes an optimum RSCE design assuming the hardware platform is customized.
1.2.2.Design and Implementation of the RSCE SystemC Model
A SystemC RSCE fixed-point functional model is developed to allow precision requirement evaluation. This same model is also used as a behavioral model in the verification testbench. Furthermore, in the future, this SystemC model can be plugged into a multi-design simulation environment to allow fast and accurate system-level simulation. The simulation model can calculate the SPME reciprocal energy and force using either the IEEE double precision arithmetic or fixed-point arithmetic of user-selected precision. The correctness of this RSCE SystemC model is verified by comparing its single timestep simulation results against the golden SPME software implementation .
This thesis is divided into six chapters. Chapter 2 provides the reader with background information on molecular dynamics, the SPME algorithm, and the NAMD program. It also briefly discusses the other relevant research efforts to speedup the non-bonded force calculation. Chapter 3 describes the design and implementation of the RSCE and provides brief information on parallelization of the SPME using multiple RSCEs. Chapter 4 discusses the limitations of the current implementation and also provides estimation on the level of speedupdegree of speedup the optimum RSCE design can provide in the reciprocal sum calculation. Chapter 5 describes the RSCE SystemC simulation model and the design verification environment. Furthermore, it also presents MD simulation results of the implemented RSCE when it is used along with NAMD2. Lastly, Chapter 6 concludes this thesis and offers recommendations for future work.
In this chapter, background information is given to provide readers with a basic understanding of molecular dynamics, non-bonded force calculations, and the SPME algorithm. First of all, the procedure of a typical MD simulation is described. Then, two types of non-bonded interactions, namely the Lennard-Jones interactions and the Coulombic interactions, are discussed along with the respective methods to minimize their computational complexity. Following that, the operation and parallelization strategy of several relevant hardware MD simulation systems are described. Afterwards, the operation of NAMD2 program is also explained. Lastly, this chapter is concluded with the significance of this thesis work in speeding up MD simulations.
Molecular dynamics [11, 12, 13] is a computer simulation technique that calculates the trajectory of a group of interacting particles by integrating their equation of motion at each timestep. Before the MD simulation can start, the structure of the system must be known. The system structure is normally obtained using Nuclear Magnetic Resonance (NMR) or an X-ray diagram and is described with a Protein Data Bank (PDB) file . With this input structure, the initial coordinates of all particles in the system are known. The initial velocities for the particles are usually approximated with a Maxwell-Boltzmann distribution . These velocities are then adjusted such that the net momentum of the system is zero and the system is in an equilibrium state.
The total number of simulation timesteps is chosen to ensure that the system being simulated passes through all configurations in the phase space (space of momentum and positions) . This is necessary for the MD simulation to satisfy the ergodic hypothesis  and thus, make the simulation result valid. The ergodic hypothesis states that, given an infinite amount of time, an NVE (constant number of particles N, constant volume V, and constant energy E) system will go through the entire constant energy hypersurface. Thus, under the ergodic hypothesis, averages of an observable over a trajectory of a system are equivalent to its averages over the microcanonical (NVE) ensemble . Under the hypothesis, the researcher can extract the same information from the trajectory obtained in an MD simulation or from the result obtained by an in-lab experiment. In the in-lab experiment, the sample used represents the microcanonical ensemble. After the number of timesteps is selected, the size of the timestep, in units of seconds, is chosen to correspond to the fastest changing force (usually the bonded vibration force) of the system being simulated. For a typical MD simulation, the timestep size is usually between 0.5fs to 1fs.
At each timestep, the MD program calculates the total forces exerted on each particle. It then uses the calculated force exerted on the particle, its velocity, and its position at the previous timestep to calculate its new position and new velocity for the next timestep. The 1st order integration of the equation of motion is used to derive the new velocities and the 2nd order is used to derive the new positions for all particles in the system. After the new positions and velocities are found, the particles are moved accordingly and the timestep advances. This timestep advancement mechanism is called time integration. Since the integration of the equations of motion cannot be solved analytically, a numerical time integrator, such as Velocity Verlet, is used to solve the equations of motion numerically . The numerical integrators being used in MD simulation should have the property of symplectic  and thus, should approximate the solution with a guaranteed error bound. By calculating the forces and performing the time integration iteratively, a time advancing snapshot of the molecular system is obtained.
During the MD simulation, there are two main types of computation: the force calculation and the time integration. The time integration is performed on each particle and thus it is an O(N) operation. On the other hand, the force calculation can be subdivided into bonded force calculation and non-bonded force calculation. Bonded force calculation is an O(N) calculation because a particle only interacts with its limited number of bonded counterparts. While the non-bonded force calculation is an O(N2) operation because each particle interacts with all other (N-1) particles in the many-body system. The non-bonded force can further be categorized into short-range Van der Waals force (described by the Lennard-Jones interaction) and long-range electrostatic force (described by the Coulombic interaction). In the following sections, both these non-bonded forces are described.