A fine grained parallel smooth particle mesh Ewald algorithm for biophysical simulation studies: Application to the 6-D torus QCDOC supercomputer

In order to model complex heterogeneous biophysical macrostructures with non-trivial charge distributions such as globular proteins in water, it is important to evaluate the long range forces present in these systems accurately and efficiently. The Smooth Particle Mesh Ewald summation technique (SPME) is commonly used to determine the long range part of electrostatic energy in large scale molecular simulations. While the SPME technique does not give rise to a performance bottleneck on a single processor, current implementations of SPME on massively parallel, supercomputers become problematic at large processor numbers, limiting the time and length scales that can be reached. Here, a synergistic investigation involving method improvement, parallel programming and novel architectures is employed to address this difficulty. A relatively simple modification of the SPME technique is described which gives rise to both improved accuracy and efficiency on both massively parallel and scalar computing platforms. Our fine grained parallel implementation of the modified SPME method for the novel QCDOC supercomputer with its 6D-torus architecture is then given. Numerical tests of algorithm performance on up to 1024 processors of the QCDOC machine at BNL are presented for two systems of interest, a β-hairpin solvated in explicit water, a system which consists of 1142 water molecules and a 20 residue protein for a total of 3579 atoms, and the HIV-1 protease solvated in explicit water, a system which consists of 9331 water molecules and a 198 residue protein for a total of 29508 atoms.     

An efficient parallel implementation of the smooth particle mesh Ewald method for molecular dynamics simulations

This paper focuses on the implementation and the performance analysis of a smooth particle mesh Ewald method on several parallel computers. We present the details of the algorithms and our implementation that are used to optimize parallel efficiency on such parallel computers.     

A scalable parallel algorithm for large-scale reactive force-field molecular dynamics simulations

A scalable parallel algorithm has been designed to perform multimillion-atom molecular dynamics (MD) simulations, in which first principles-based reactive force fields (ReaxFF) describe chemical reactions. Environment-dependent bond orders associated with atomic pairs and their derivatives are reused extensively with the aid of linked-list cells to minimize the computation associated with atomic n-tuple interactions (n?4 explicitly and ?6 due to chain-rule differentiation). These n-tuple computations are made modular, so that they can be reconfigured effectively with a multiple time-step integrator to further reduce the computation time. Atomic charges are updated dynamically with an electronegativity equalization method, by iteratively minimizing the electrostatic energy with the charge-neutrality constraint. The ReaxFF-MD simulation algorithm has been implemented on parallel computers based on a spatial decomposition scheme combined with distributed n-tuple data structures. The measured parallel efficiency of the parallel ReaxFF-MD algorithm is 0.998 on 131,072 IBM BlueGene/L processors for a 1.01 billion-atom RDX system.     

Embedded divide-and-conquer algorithm on hierarchical real-space grids: parallel molecular dynamics simulation based on linear-scaling density functional theory

A linear-scaling algorithm has been developed to perform large-scale molecular-dynamics (MD) simulations, in which interatomic forces are computed quantum mechanically in the framework of the density functional theory. A divide-and-conquer algorithm is used to compute the electronic structure, where non-additive contribution to the kinetic energy is included with an embedded cluster scheme. Electronic wave functions are represented on a real-space grid, which is augmented with coarse multigrids to accelerate the convergence of iterative solutions and adaptive fine grids around atoms to accurately calculate ionic pseudopotentials. Spatial decomposition is employed to implement the hierarchical-grid algorithm on massively parallel computers. A converged solution to the electronic-structure problem is obtained for a 32,768-atom amorphous CdSe system on 512 IBM POWER4 processors. The total energy is well conserved during MD simulations of liquid Rb, showing the applicability of this algorithm to first principles MD simulations. The parallel efficiency is 0.985 on 128 Intel Xeon processors for a 65,536-atom CdSe system.     

A space-time-ensemble parallel nudged elastic band algorithm for molecular kinetics simulation

A scalable parallel algorithm has been designed to study long-time dynamics of many-atom systems based on the nudged elastic band method, which performs mutually constrained molecular dynamics simulations for a sequence of atomic configurations (or states) to obtain a minimum energy path between initial and final local minimum-energy states. A directionally heated nudged elastic band method is introduced to search for thermally activated events without the knowledge of final states, which is then applied to an ensemble of bands in a path ensemble method for long-time simulation in the framework of the transition state theory. The resulting molecular kinetics (MK) simulation method is parallelized with a space-time-ensemble parallel nudged elastic band (STEP-NEB) algorithm, which employs spatial decomposition within each state, while temporal parallelism across the states within each band and band-ensemble parallelism are implemented using a hierarchy of communicator constructs in the Message Passing Interface library. The STEP-NEB algorithm exhibits good scalability with respect to spatial, temporal and ensemble decompositions on massively parallel computers. The MK simulation method is used to study low strain-rate deformation of amorphous silica.     

A general purpose parallel molecular dynamics simulation program

We present a general purpose parallel molecular dynamics simulation code. The code can handle NVE, NVT, and NPT ensemble molecular dynamics, Langevin dynamics, and dissipative particle dynamics. Long-range interactions are handled by using the smooth particle mesh Ewald method. The implicit solvent model using solvent-accessible surface area was also implemented. Benchmark results using molecular dynamics, Langevin dynamics, and dissipative particle dynamics are given.

Program summary

Title of program:MM_PARCatalogue identifier:ADXP_v1_0Program summary URL:http://cpc.cs.qub.ac.uk/summaries/ADXP_v1_0Program obtainable from: CPC Program Library, Queen's University of Belfast, N. IrelandComputer for which the program is designed and others on which it has been tested:any UNIX machine. The code has been tested on Linux cluster and IBM p690Operating systems or monitors under which the program has been tested:Linux, AIXProgramming language used:CMemory required to execute with typical data:∼60 MB for a system of atoms Has the code been vectorized or parallelized? parallelized with MPI using atom decomposition and domain decompositionNo. of lines in distributed program, including test data, etc.:171 427No. of bytes in distributed program, including test data, etc.:4 558 773Distribution format:tar.gzExternal routines/libraries used:FFTW free software (http://www.fftw.org)Nature of physical problem:Structural, thermodynamic, and dynamical properties of fluids and solids from microscopic scales to mesoscopic scales.Method of solution:Molecular dynamics simulation in NVE, NVT, and NPT ensemble, Langevin dynamics simulation, dissipative particle dynamics simulation.Typical running time:Table below shows the typical run times for the four test programs.

Benchmark results. The values in the parenthesis are the number of processors used
System	Method	Timing for 100 steps in seconds
256 TIP3P	MD	23.8 (1)
64 DMPC + 1645 TIP3P	MD	890 (1)	528 (2)	326 (4)	209 (8)
8 Aβ16-22	LD	1.02 (1)
23760 Groot-Warren particles	DPD	22.16 (1)

Full-size table

     

Multibillion-atom molecular dynamics simulation: Design considerations for vector-parallel processing

Progress in adapting molecular dynamics algorithms for systems with short-range interactions to utilize the features of modern supercomputers is described. Efficient utilization of the latest generation of processor architectures requires algorithms that can be both vectorized and parallelized. The approach adopted for vectorization involves combining the layer and neighbor-list methods, while parallelization employs spatial subdivision with explicit communication. The techniques presented here have been used in performance tests on the Cray X1 vector-parallel supercomputer with systems containing over 12 billion atoms.     

Pathfinder: A parallel search algorithm for concerted atomistic events

An algorithm has been designed to search for the escape paths with the lowest activation barriers when starting from a local minimum-energy configuration of a many-atom system. The pathfinder algorithm combines: (1) a steered eigenvector-following method that guides a constrained escape from the convex region and subsequently climbs to a transition state tangentially to the eigenvector corresponding to the lowest negative Hessian eigenvalue; (2) discrete abstraction of the atomic configuration to systematically enumerate concerted events as linear combinations of atomistic events; (3) evolutionary control of the population dynamics of low activation-barrier events; and (4) hybrid task + spatial decompositions to implement massive search for complex events on parallel computers. The program exhibits good scalability on parallel computers and has been used to study concerted bond-breaking events in the fracture of alumina.     

A hybrid multi-loop genetic-algorithm/simplex/spatial-grid method for locating the optimum orientation of an adsorbed protein on a solid surface

Atomistic simulation of protein adsorption on a solid surface in aqueous environment is computationally demanding, therefore the determination of preferred protein orientations on the solid surface usually serves as an initial step in simulation studies. We have developed a hybrid multi-loop genetic-algorithm/simplex/spatial-grid method to search for low adsorption-energy orientations of a protein molecule on a solid surface. In this method, the surface and the protein molecule are treated as rigid bodies, whereas the bulk fluid is represented by spatial grids. For each grid point, an effective interaction region in the surface is defined by a cutoff distance, and the possible interaction energy between an atom at the grid point and the surface is calculated and recorded in a database. In searching for the optimum position and orientation, the protein molecule is translated and rotated as a rigid body with the configuration obtained from a previous Molecular Dynamic simulation. The orientation-dependent protein-surface interaction energy is obtained using the generated database of grid energies. The hybrid search procedure consists of two interlinked loops. In the first loop A, a genetic algorithm (GA) is applied to identify promising regions for the global energy minimum and a local optimizer with the derivative-free Nelder-Mead simplex method is used to search for the lowest-energy orientation within the identified regions. In the second loop B, a new population for GA is generated and competitive solution from loop A is improved. Switching between the two loops is adaptively controlled by the use of similarity analysis. We test the method for lysozyme adsorption on a hydrophobic hydrogen-terminated silicon (110) surface in implicit water (i.e., a continuum distance-dependent dielectric constant). The results show that the hybrid search method has faster convergence and better solution accuracy compared with the conventional genetic algorithm.     

An investigation of soft-core potentials for the simulation of mesogenic molecules and molecules composed of rigid and flexible segments

The phase behaviour of three soft core spherocylinder models is investigated with a view to producing an effective potential for use in coarse-grained simulations of liquid crystal phases and polymers composed of rigid and flexible segments. Provided potentials are not made too soft, two of the soft core models are found to work well in terms of successfully reproducing mesophases and in providing considerable improvements in computational speed over other commonly used coarse-grained models. In Monte Carlo simulations a soft-core spherocylinder model in which a cut and shifted Lennard-Jones potential is truncated with a linear tangential potential is found to be particularly effective; while for molecular dynamics a better model is provided by a DPD-like quadratic potential. Here, computational speed-ups of 20-30× are seen in equilibration times in comparison to the well-known soft repulsive spherocylinder (SRS) model. The quadratic potential is used in an additional set of coarse-grained simulations of a liquid crystal with a flexible chain, which exhibits spontaneous formation of a nematic phase. The use of different types of interaction sites is also illustrated by the simulation of a spherocylinder with two “tails” formed from spheres. Here, varying the hardness of the sphere-spherocylinder interaction potential allows the formation of a smectic-A phase which exhibits microphase separation.     

A general-purpose coarse-grained molecular dynamics program

In this article, we describe a general-purpose coarse-grained molecular dynamics program COGNAC (COarse Grained molecular dynamics program by NAgoya Cooperation). COGNAC has been developed for general molecular dynamics simulation, especially for coarse-grained polymer chain models. COGNAC can deal with general molecular models, in which each molecule consists of coarse-grained atomic units connected by chemical bonds. The chemical bonds are specified by bonding potentials for the stretching, bending and twisting of the bonds, each of which are the functions of the position coordinates of the two, three and four atomic units. COGNAC can deal with both isotropic and anisotropic interactions between the non-bonded atomic units. As an example, the Gay-Berne potential is implemented. New potential functions can be added to the list of existing potential functions by users. COGNAC can do simulations for various situations such as under constant temperature, under constant pressure, under shear and elongational deformation, etc. Some new methods are implemented in COGNAC for modeling multiphase structures of polymer blends and block copolymers. A density biased Monte Carlo method and a density biased potential method can generate equilibrium chain configurations from the results of the self-consistent field calculations. Staggered reflective boundary conditions can generate interfacial structures with smaller system size compared with those of periodic boundary conditions.     

TiReX: Replica-exchange molecular dynamics using Tinker

We present a driver program for performing replica-exchange molecular dynamics simulations with the Tinker package. Parallelization is based on the Message Passing Interface, with every replica assigned to a separate process. The algorithm is not communication intensive, which makes the program suitable for running even on loosely coupled cluster systems. Particular attention is paid to the practical aspects of analyzing the program output.

Program summary

Program title: TiReXCatalogue identifier: AEEK_v1_0Program summary URL:http://cpc.cs.qub.ac.uk/summaries/AEEK_v1_0.htmlProgram obtainable from: CPC Program Library, Queen's University, Belfast, N. IrelandLicensing provisions: Standard CPC licence, http://cpc.cs.qub.ac.uk/licence/licence.htmlNo. of lines in distributed program, including test data, etc.: 43 385No. of bytes in distributed program, including test data, etc.: 502 262Distribution format: tar.gzProgramming language: Fortran 90/95Computer: Most UNIX machinesOperating system: LinuxHas the code been vectorized or parallelized?: parallelized with MPIClassification: 16.13External routines: TINKER version 4.2 or 5.0, built as a libraryNature of problem: Replica-exchange molecular dynamics.Solution method: Each replica is assigned to a separate process; temperatures are swapped between replicas at regular time intervals.Running time: The sample run may take up to a few minutes.     

A handy tool for history keeping of Geant4 tracks - J4HistoryKeeper

The Particle Flow Analysis (PFA) is currently under intense studies as the most promising way to achieve precision jet energy measurements required at the future linear e⁺e⁻ collider. In order to optimize detector configurations and to tune up the PFA it is crucial to identify factors that limit the PFA performance and clarify the fundamental limits on the jet energy resolution that remain even with the perfect PFA and an infinitely granular calorimeter. This necessitates a tool to connect each calorimeter hit in particle showers to its parent charged track, if any, and eventually all the way back to its corresponding primary particle, while identifying possible interactions and decays along the way. In order to realize this with a realistic memory space, we have developed a set of C++ classes that facilitates history keeping of particle tracks within the framework of Geant4. This software tool, hereafter called J4HistoryKeeper, comes in handy in particular when one needs to stop this history keeping for memory space economy at multiple geometrical boundaries beyond which a particle shower is expected to start. In this paper this software tool is described and applied to a generic detector model to demonstrate its functionality.     

AFMM: A molecular mechanics force field vibrational parametrization program

AFMM (Automated Frequency Matching Method) is a program package for molecular mechanics force field parametrization. The method used fits the molecular mechanics potential function to both vibrational frequencies and eigenvector projections derived from quantum chemical calculations. The program optimizes an initial parameter set (either pre-existing or using chemically-reasonable estimation) by iteratively changing them until the optimal fit with the reference set is obtained. By implementing a Monte Carlo-like algorithm to vary the parameters, the tedious task of manual parametrization is replaced by an efficient automated procedure. The program is best suited for optimization of small rigid molecules in a well-defined energy minimum, for which the harmonic approximation to the energy surface is appropriate for describing the intra-molecular degrees of freedom.

Program summary

Title of program: AFMMCatalogue identifier: ADUZProgram summary URL:http://cpc.cs.qub.ac.uk/summaries/ADUZProgram obtainable from: CPC Program Library, Queen's University of Belfast, N. IrelandComputer: x86 PC, SGI, Sun MicrosystemsOperating system: GNU/Linux, BSD, IRIX, SolarisProgramming language used: PythonMemory required: 10 MbytesNo. of bits in a word: 32 or 64No. of processors used: 1Parallelized?: NoNo. of lines in distributed program, including test data, etc.:13 127No. of bytes in distributed program, including test data, etc.: 182 550Distribution format: tar.gzTypical running time: 24 hNature of the physical problem: Molecular mechanics force field parametrization.Method of solution:Fitting of the molecular mechanics potential to normal modes derived from quantum chemical calculations. The missing force field parameters are optimized via a merit function to obtain the optimal fit with the reference quantum mechanical set.     

Numerical paraxial approximation for highly relativistic beams

This paper is concerned with the numerical simulation of highly relativistic beams. We developed a particle-in-cell code for highly relativistic beams, based on the paraxial Vlasov-Maxwell formulation of Laval et al. This formulation follows the beam in a speed-of-light frame. It is fourth order accurate in the small characteristic velocity of the beam. The formulation is simpler than standard particle-in-cell codes in the lab frame or in the beam frame and gives a fast and easy to implement algorithm. Numerical examples illustrate the accuracy of the method.     

A fast algorithm for voxel-based deterministic simulation of X-ray imaging

Deterministic method based on ray tracing technique is known as a powerful alternative to the Monte Carlo approach for virtual X-ray imaging. The algorithm speed is a critical issue in the perspective of simulating hundreds of images, notably to simulate tomographic acquisition or even more, to simulate X-ray radiographic video recordings. We present an algorithm for voxel-based deterministic simulation of X-ray imaging using voxel-driven forward and backward perspective projection operations and minimum bounding rectangles (MBRs). The algorithm is fast, easy to implement, and creates high-quality simulated radiographs. As a result, simulated radiographs can typically be obtained in split seconds with a simple personal computer.

Program summary

Program title: X-rayCatalogue identifier: AEAD_v1_0Program summary URL:http://cpc.cs.qub.ac.uk/summaries/AEAD_v1_0.htmlProgram obtainable from: CPC Program Library, Queen's University, Belfast, N. IrelandLicensing provisions: Standard CPC licence, http://cpc.cs.qub.ac.uk/licence/licence.htmlNo. of lines in distributed program, including test data, etc.: 416 257No. of bytes in distributed program, including test data, etc.: 6 018 263Distribution format: tar.gzProgramming language: C (Visual C++)Computer: Any PC. Tested on DELL Precision 380 based on a Pentium D 3.20 GHz processor with 3.50 GB of RAMOperating system: Windows XPClassification: 14, 21.1Nature of problem: Radiographic simulation of voxelized objects based on ray tracing technique.Solution method: The core of the simulation is a fast routine for the calculation of ray-box intersections and minimum bounding rectangles, together with voxel-driven forward and backward perspective projection operations.Restrictions: Memory constraints. There are three programs in all.

•

A. Program for test 3.1(1): Object and detector have axis-aligned orientation;

•

B. Program for test 3.1(2): Object in arbitrary orientation;

•

C. Program for test 3.2: Simulation of X-ray video recordings.

1.

Program A Memory required to execute with typical data: 207 Megabytes, depending on the size of the input file. Typical running time: 2.30 s. (Tested in release mode, the same below.)

2.

Program B (the main program) Memory required to execute with typical data: 114 Megabytes, depending on the size of the input file. Typical running time: 1.60 s.

3.

Program C Memory required to execute with typical data: 215 Megabytes, depending on the size of the input file. Typical computation time: 27.26 s for cast-5, 101.87 s for cast-6.

     

Development of a Geant4 solid for stereo mini-jet cells in a cylindrical drift chamber

Stereo mini-jet cells will be indispensable components of a future e⁺e⁻ linear collider central tracker such as JLC-CDC. There is, however, no official Geant4 solid available at present to describe such geometrical objects, which had been a major obstacle for us to develop a full Geant4-based simulator with stereo cells built in. We have thus extended Geant4 to include a new solid (TwistedTubs), which consists of three kinds of surfaces: two end planes, inner and outer hyperboloidal surfaces, and two so-called twisted surfaces that make slant and twisted φ-boundaries. Design philosophy and its realization in the Geant4 framework are described together with algorithmic details. We have implemented stereo cells with the new solid, and tested them using geantinos and Pythia events (e⁺e⁻→ZH at  GeV). The performance was found reasonable: the stereo cells consumed only 25% more CPU time than ordinary axial cells.