SIMES : SImulation at MESoscopic Scale — Accelerating Dissipative Particle Dynamics Simulations on Graphical Processing Units

Computer simulations and in particular mesoscopic simulation techniques, such as dissipative particle dynamics (DPD) enable researchers to study the complexities of soft materials and polymeric systems by performing in silico experimentations alongside with in vivo experiments. In this work we describe a complete implementation of DPD running entirely on a graphics processing unit (GPU). The design of the algorithms and optimizations needed to fully take advantage of a GPU are discussed. The performance of the code is evaluated and shown to be up to more than 30 times faster than a conventional implementation running on a single CPU core. We illustrate the potential of using DPD on GPU’s with applications to the physical and biological sciences as well as in engineering areas, where it can be a novel and versatile tool. SIMES is free for academics and can be downloaded from our web site: http://www.simes.uaemex-labs.org.mx.


Introduction
The traditional picture of the role of microscopic and macroscopic physics is now being challenged as new multi-scale and multi-disciplinary problems begin to emerge.For example, in nanoscale systems, the assumption of scale separation breaks down; macroscopic theory is therefore inadequate, yet microscopic theory may be impractical because it requires computational capabilities far beyond our present reach.This new class of problems poses unprecedented challenges to mathematical modeling as well as numerical simulation and requires new and non-traditional analysis and modeling paradigms.Methods based on mesoscopic theories, which connect the microscopic and macroscopic descriptions of the dynamics, provide a promising approach.One of the most popular and efficient numerical methods for studying systems at the mesoscopic scale is dissipative particle dynamics (DPD).This method was first proposed by Hoogerbrugge and Koelman in 1992 as a particle-based method for studying hydrodynamic phenomena [1].The statistical mechanics foundation of DPD was developed by Español and Warren in 1995 [2], followed by a rigorous analysis of the nature of its simulation parameters by Groot and Warren in 1997 [3].DPD simulation can be thought of as a valuable approach to bridge the gap between the atomistic and mesoscopic length scales, yet it offers insights into molecular interactions in soft matter systems.We refer the readers to reference [4] for a thorough review of the applications of DPD simulations to soft matter systems.For additional applications the reader is referred to the review published by one of us (AGG) [5], where one can find a number of recent examples and applications of the DPD method to soft matter systems in equilibrium and under flow.In this work, we introduce a new software that allows us to accelerate DPD simulations.Our platform, SIMES, which is an acronym for SImulation at MESoscopic time-scale, is a production DPD software designed and implemented specifically to run on GPUs.We encourage interested readers to see the full development of SIMES which is documented in the PhD thesis of one us (KATM) [6].

Methodology
The elementary units in a DPD simulation are fluid elements or soft beads.A bead represents a volume of fluid that is large on a molecular scale, hence it contains at least several molecules of the fluid, but is still macroscopically small.Beads interact via effective forces chosen so as to reproduce the hydrodynamic behavior of the fluid without reference to its molecular structure.DPD differs in this respect from atomistic simulations such as molecular dynamics (MD), in which the forces are chosen to model the intermolecular interactions of a system as accurately as possible.Forces in DPD are pairwise additive, conserve momentum, have no hard core and are short-ranged, with the range of the force defining the size of the soft beads.A detailed review of the DPD methodology that is implemented in SIMES can be found in the supplementary information (SI) that accompanies this article.

Implementation and Design
SIMES is written entirely on CUDA GPU programming toolkit developed by NVIDIA; it targets NVIDIA GPUs with compute capability 3.0 or higher.CUDA is a dialect of the C and C++ programming languages with extended syntax supporting multiple memory address spaces and GPU kernel launch, while eliminating language features, such as recursion, that do not map well onto the many-core GPU hardware.Today, SIMES consists of about 10 4 lines of code.One of the main features of SIMES is that it is explicitly designed for execution by a single workstation with multiple GPU's.In general terms, similarly to other MD codes, SIMES has six boxes of computation.The first box consist of instructions to make a Verlet list [7].In the second box, the positions and velocities are initialized.The third box is most important due to the fact that the evaluation of conservative, dissipative and ramdon forces is performed here.The fourth box is dedicated to the evaluation of the bonds and angles force functions.In the fifth box, we evaluate the contribution of electrostatic forces.Finally, in the sixth box, the positions, velocities, and forces are updated for the new cycle.In Figure 1, a flow chart as well as the main routines in pseudo code are displayed.

Performance
The performance tests of SIMES were carried out on different GPU´s manufactured by NVIDIA.This performance is based on all architectures of NVIDIA GPU series, from first generation Tesla (2007-09), Fermi (2010-11), Kepler (2012-14), to Maxwell (2015-16).In Figure 2, a complete analysis of the performance of SIMES is shown.Figure 2A shows that the best performance of SIMES is achieved with GTX 780 GPU cards.It has been documented that the best MD codes running on GPU´s reach their best performance on these cards [8].To compare the performance of the GPU and CPU implementations, we measured the speedup factor (x) which is defined as x = TCPU/TGPU, where TCPU is the time used by the CPU implementation and TGPU is the time used by the GPU implementation.Figure 2B shows that the expected behavior of x is obtained, namely x decreases as the number of cores increases.Finally, we show in Figure 2C a snapshot of a system composed of half a million DPD particles.This is the maximum number of particles that can be simulated on a GTX-780, using the maximum RAM.In the SI section there are more technical details on the GPUs cards that were used in this study.

Validation
To validate the efficiency and reliability in the calculation of thermodynamic properties with SIMES, we performed four validation tests to ascertain that the GPU version of DPD gives results consistent with the original CPU version.The calculations referred to as "CPUS" were made with our serial code programmed in Fortran 90.We examined thermodynamic quantities in equilibrium such as radial distribution function (g(r)) (Figure 3A), electrostatic energy (Figure 3B), total energy (Figure 3C) and temperature (Figure 3D), from a simulation in the canonical ensemble (constant temperature and volume).The calculations of temperature on GPU and CPU are in excellent agreement; for example, the average values for temperature are 1.00 ± 0.00003 and 1.01 ± 0.00454 for CPU and GPU, respectively, while the average values for the total energy are 12.6 ± 0.2 and 12.6 ± 0.3 for CPU and GPU, respectively.For additional technical details about these simulations, see the SI material.

Applications and Emerging Domains
Several computationally-intensive and emerging applications have also been the subject of extensive research through DPD simulations.Today, with SIMES it is possible to study a wide range of applications, ranging from life sciences to engineering.In what follows, some of those applications and works, both published and also unpublished are listed, to illustrate the capabilities of SIMES.

Colloidal Dispersions
Colloidal dispersions have numerous applications in various fields such as chemistry, biology, medicine, and engineering.Interaction forces between colloidal particles in all suspensions/emulsions play an important role in determining the properties of the materials, such as the shelf life, stability, rheology and flavor, the behavior of a number of industrial processes (e.g.mixings, membrane filtrations) as well as the formula of chemical and pharmaceutical products.A fundamental aspect for studying colloidal dispersions under the framework of DPD lies in the model interaction potential used.SIMES has implemented a new potential that allows for the study of colloidal interactions at the mesoscopic scale.This new effective potential features adjustable parameters that represent various physical situations such as colloidal concentration, polymer concentration and chain length.This model represents one of the most modern alternatives for the modeling of colloidal dispersions.The novelty of the new potential is its ability to predict correctly the potential of mean force of structured colloids at the mesoscale with explicitly included polymer chains grafted on the surfaces as dispersants.A snapshot of this type of system is shown in Figure 4A.Additional technical details about these simulations are found in [9] and an illustrative movie can be found in the SI.

Confined Fluids
With the methodology implemented in SIMES it is possible to study interfacial and structural properties of fluids confined by surfaces of different geometries.Performing simulations with SIMES, it was found that confined fluids present an order-disorder phase transition on a mesoscopic scale [10].The ouput of this type of simulation, such as interfacial tension profiles for several types of surfaces, allows one to find that the most important factor in determining the value of the interfacial tension is not the geometry of the walls but the intensity of the solid wall-fluid interaction.The periodicity of the symmetry or the lack thereof on the walls induces order or disorder on the fluid molecules adjacent to those walls, but plays a very minor role in the measurable value of the interfacial tension.Undoubtedly, these findings are relevant to researchers working on enhanced recovery of oil from complex pores and in microfluidics, among other nanotechnological applications.Additional technical details about this study are found in [10] and an illustrative movie can be found in the SI material.

Enhanced Oil Recovery
Enhanced recovery is at the heart of modern oil production from underground reservoirs.If the average worldwide recovery factor from hydrocarbon reservoirs can be increased beyond current limits, it will alleviate a number of issues related to global energy supply.Currently, the daily oil production comes from mature or maturing oil fields and reserves replacement is not keeping pace with the growing energy demand.The world average recovery factor from hydrocarbon reservoirs is stuck in the mid-30 percent range.This challenge becomes an opportunity for advanced secondary and enhanced oil recovery (EOR) technologies that may mitigate the demand-supply balance.The need to extract oil from wells where it is adsorbed on the surfaces of rocks has led to the development of new and improved, EOR techniques.One of those is the injection of surfactants with water vapor which promote the desorption of oil, which then can be extracted using pumps, as the surfactants encapsulate the oil in foams.With SIMES we have simulated surfactants and oil molecules adsorbed on surfaces under flow, with the purpose of studying the efficiency of the surfactants to desorb hydrocarbon chains, that are found adsorbed on flat surfaces.Additional technical details about this study can be found in [11] and an illustrative movie is found in the SI material.

Biological Membranes
Fusion of biological membranes is an essential process in many areas of cell biology, ranging from vesicular trafficking and synaptic transmission to cell-cell fusion or viral fusion.Lipid vesicles, which are often used as simplified model systems for complex biological membranes, can also be induced to fuse experimentally by a variety of methods.Membrane fusion is fundamental to the life of eukaryotic cells, for example.On the other hand, cellular trafficking and compartmentalization, import of food stuffs and export of waste, inter-cellular communication, sexual reproduction, and cell division are all dependent on this basic process.Yet, little is known about the molecular mechanism(s) by which fusion occurs.Figures 4 B-D are example of biological membranes modeled and solved with SIMES. Figure 4B shows a snapshot of a section of a vein coated with cholesterol under blood flow (not shown for clarity); Figure 4C is a snapshot of two opposite cellular membranes coated with biological brushes.Finally, in Figure 4D we present a snapshot of an explcitly modeled membrane of PCCP.Additional technical details about these case studies are found in [12]; an illustrative movie is available in the SI material.

Conclusions
We have presented here the structure, performance and applications of SIMES, which runs on GPU.The potential of this tools is enormous and it is our hope that the community can benefit from and improve upon it, with more novel and useful applications.

Figure 1 .
Figure 1.Principal flow diagram programing of SIMES (left) and examples of pseudocode of some of its main subroutines (right).

Figure 2 .
Figure 2. SIMES performance on GPUs.(A) Performance on different GPUs cards.(B) Speed up of SIMES on a GPU GTX-980 versus Intel core i7 CPUs.(C) A snapshot of a system composed of half million DPD particle.Cards type GTX-580 support this system size according to its the memory capacity.

Figure 4 .
Figure 4. Examples of applications of SIMES to biological membranes.(A) Snapshot of DPD particle-made colloids (in blue) with grafted polymer chains (in red).(B) Snapshot of a section of a vein coated with cholesterol under blood flow (not shown for clarity).(C) Snapshot of two opposite cellular membranes coated with biological brushes.(D) An explicitly formed PCCP membrane.