SciDAC Review 14 - MOLECULAR MODELING AND SIMULATION: Practical Solutions to Global Challenges

MOLECULAR MODELING AND SIMULATION

Practical Solutions to Global CHALLENGES

Molecular modeling and simulation is used to study molecular structure and chemical reactions via the equations of quantum mechanics and classical physics. Molecular modeling and simulation approaches are important because they help in providing scientific insight and understanding where this is not available through experimentation alone. Moreover, the building of new models and theories will ultimately enable testable predictions that will lead to breakthroughs in many scientific domains. The fundamental equation of quantum mechanics for chemistry is the Schrödinger equation. From this equation all molecular properties can be derived. However, this equation is not reasonably computable except for a very few small atomic and molecular systems. Most practical computational studies rely on the use of various approximate methods.

Molecular modeling and simulation are becoming indispensable tools in the arsenal needed to address many energy and environmental problems facing the Department of Energy (DOE) and the nation. Cutting-edge modeling and simulation methods have been developed to address scientific questions relevant to dynamic and reactive chemical processes such as catalysis, photosynthesis, hydrogen storage, protein function, and environmental remediation. As a result, molecular modeling and simulation find application in many areas of critical importance, such as solving energy challenges to enable clean, secure, efficient, and affordable energy systems. Understanding the molecular processes in next-generation fuel cells and battery technologies will play an important role in the design of new approaches for future alternative energy production and use. In environmental remediation, computer models are used to help researchers determine how certain microbial systems may play an important role in the remediation of contaminated soils.

Molecular modeling and simulation find application in many areas of critical importance, such as solving energy challenges to enable clean, secure, efficient, and affordable energy systems.

New advances in computer science are making modeling and simulation feasible on some of the largest, massively parallel computers available. State-of-the-art tools are being developed by and for the computational science community world-wide. One of the tools for molecular modeling and simulation developed for the DOE is NWChem, the highly-scalable software tool for complex chemistry and biochemistry calculations on massively parallel, high-performance computers (sidebar “NWChem” p13). Scientists are using the NWChem software to investigate chemical processes in molecular and biomolecular systems, ranging in size from tens to millions of atoms, to predict structure, properties, and reactivity. A variety of computational techniques has been implemented, ranging from molecular dynamics simulations, which are based on parameterized atomic interaction functions and can be applied to systems with millions of atoms, to highly-correlated quantum mechanical calculations designed to understand detailed properties of much smaller systems. NWChem enables a wide variety of classical simulation and electronic structure calculations that, in addition, can be combined to perform quantum molecular dynamics and hybrid quantum and classical mechanical simulations. Hybrid methods are useful for the study of systems in which fast, quantum-mechanical processes such as electron transfer are coupled to much slower atomic motions. An example is the electron transfer between heme groups in cytochromes, as illustrated in figure 1. The design and implementation of NWChem targets high-performance, massively-parallel computer architectures, but the code also runs on single machines or clusters of desktop workstations and personal computers running the LINUX operating system.

T. P. Straatsma, PNNL

Figure 1. Iron-dependent cytochrome c3 fumarate reductase is an enzyme expressed by Shewanella oneidensis under iron-rich conditions. The structure of this protein as shown here illustrates the electron transfer pathway through four heme groups to the FAD in the active site where fumarate is reduced to succinate. Modeling of this enzyme requires classical dynamics simulations to understand the effect of protein dynamics on the electron transfer kinetics between hemes, which can only be appropriately described using quantum mechanical methods.

To computational chemists, a range of computational methods is available with varying computational complexity and accuracy. The choice of method heavily depends on the nature of the molecular system and the kind of problem being solved. Three algorithms that are among a suite of capabilities implemented in NWChem are highlighted and illustrated here: classical molecular dynamics simulations, ab initio molecular dynamics, and coupled cluster calculations.

Classical Molecular Dynamics Simulations
Specific aspects of molecular systems are easily accessible by computer simulation, but impossible, impractical, or prohibitively expensive by experimental methods. Binding in protein complexes often depends on complementary electrostatic signatures that are easily calculated, such as for the signaling pathway proteins ras and raf, illustrated in figure 2. In biomolecular systems, the action of components or complexes of components is dependent on subtle, concerted atomic motions that cannot be deduced from the ensemble and time-averaged properties provided by practical measurements. Because the spatial and temporal dimensions need to be expanded to biologically relevant space and time scales, molecular dynamics simulations are typically carried out using classical, effective potential functions to describe the inter-atomic interactions used in the time-integration of Newtonian equations of motion. The number of interactions calculated in a single time step can be in the hundreds of millions, and a typical simulation is carried out for tens of millions of time steps. Although the interaction functions are relatively inexpensive to evaluate, making it feasible to study systems with very large numbers of atoms for extended simulation times, the computational resources required are considerable. Therefore, as these simulations progress, the generated time-trajectories are stored for detailed post-analysis of the dynamical behavior. Where properties obtained as time-averages from computer simulations of single molecules can be compared to experimental averages obtained from an ensemble of such molecules—commonly referred to as ergodicity—the application of statistical mechanical methods may be used to also obtain thermodynamic properties.



T. P. Straatsma, PNNL
Figure 2. Computer modeling provides the ability to study binding of protein complexes at an atomic level of detail. The electrostatic signatures shown here for the faces of ras (left) and raf (right), two proteins that form a crucial complex in the ras–raf–MEK regulatory pathway, illustrate that, in addition to atomic flexibility, electrostatic complementarity plays an important role in protein complex formation.

The molecular dynamics simulation capability of NWChem is based on domain decomposition of the molecular system. This means that each processing unit of a parallel computer is responsible for a small region (domain) of the molecular system. This approach offers the most economical use of the available memory of a computer, and minimizes time-consuming data transfers between processing elements. By optimizing these important measures, designers of simulation codes attempt to obtain high scalability, that is, doubling the number of processors means halving the execution time. Domain decomposition for molecular dynamics simulations is efficient because it takes advantage of the locality of data needed to evaluate atomic interactions within a certain distance. This greatly reduces the number of interactions to be evaluated. To account for the electrostatic interactions beyond this distance that are not explicitly included in the calculations, the particle-mesh Ewald method has been implemented to include the effect of the long-range interactions in an approximate way by using a discrete convolution on an interpolated grid of atomic charges. To make this method efficient requires carrying out the well-known mathematical technique of Fourier transforms. A parallel fast Fourier transform has been implemented in NWChem that distributes the grids employed in this type of calculation in a way that allows for highly-efficient computation of the transform. In addition, this implementation allows a more flexible concurrent execution of different parts of the problem on subsets of the available processors, leading to better scalability of the simulation code.

Knowledge of complex molecular systems that explains and predicts their behavior is often phenomenological and extracted from comparative analyses of large datasets. Even if derived from computer simulations based on fundamental equations of physics and thermodynamics, insight into the detailed dynamical behavior and function of molecular systems is generally not straightforward. Because of the nature of the resulting data and similarity in the approaches needed to analyze the data, molecular modeling and simulations are often referred to as simulation experiments. Extracting detailed molecular-level understanding relies on the ability to carry out systematic comparative analyses of data from a series of computer experiments to develop new predictive models.

Molecular dynamics simulations of large, complex systems, however, are computationally demanding. Many of these systems are inherently hierarchical and multi-scale in character, in time and space, requiring the development of new numerical algorithms. Large-scale computer simulations produce massive scientific datasets that need to be analyzed for knowledge extraction. Effective data analysis software is required to effectively deal with the generated datasets. Use of computer simulations that has been made possible by the increasing availability of large-scale computing resources is resulting in a significant challenge of storage, analysis, and visualization of the large datasets produced in these simulations.

Molecular dynamics simulations of large, complex systems are computationally demanding.

An example where large-scale biomolecular simulations are used is in studying the role of molecular transport across the outer envelope of Gram-negative bacteria, which is of great interest from an environmental and bioenergy production perspective. Many important processes in the subsurface such as oxidation/reduction reactions, mineral dissolution, and metal ion precipitation, are microbially-mediated and believed to take place at the microbial membrane, the interface between the microbial membrane and mineral surfaces, or inside the bacterial cell through material transport across the outer membrane. Of special interest are the outer membranes of Gram-negative bacteria, in particular those found to have the ability to use metal reduction in the respiratory cycle, take up solvated metal ions from the environment, show adsorption to mineral surfaces, and, consequently, are potentially important target microbes in the design of bioremediation technologies. Because of the complexity and large size of these systems, our theoretical understanding of the processes that take place at the interface between biological membranes and geochemical environments is limited. With modern massively parallel computers and sophisticated, highly-efficient software for computational chemistry modeling, the theoretical study of such systems is now feasible.

The complexity of the interaction of microbes with their environment is evident from the broad range of scales involved, from the atomic scale involving individual functional groups (at the Å scale), the molecular scale in the formation of long molecular chains (at the nm scale), to the scale of molecular assemblies leading to the formation of membrane (at the mm scale), and molecular to macroscopic scale of the interactions of membranes with minerals and other surfaces (at the mm scale and beyond). Key geochemical and biochemical interactions and reactions at the interface of microbial membranes take place across all of these spatial scales. The dynamics and energetics of the interactions between bacterial membranes and mineral surfaces and the process of exchange and uptake by the membrane of charged species from the mineral or from solution is complex, and the understanding of these complex systems at an atomic level of detail can only be accomplished by an integrated approach of computational modeling and experimental capabilities.

Gram-negative microbial membranes are composed of many complex interacting components, including lipid membranes, polysaccharide membranes, trans-membrane proteins, minerals, and solvated ions. These systems not only involve many spatial and temporal scales, but also involve complex assemblies and exhibit complex behavior. Realistic simulations, which need to include all these components, are complicated because of the limited experimental data available that allow the construction of a molecular model for a lipopolysaccharide molecule; the need for parameterization of the oligosaccharides; the complicated setup and equilibration procedures as a result of the highly charged lipopolysaccharide molecules that must be neutralized by counter ions; and the significant computational requirements because of the large size of the molecular system. Using appropriate models can have a significant impact, even if this leads to computationally more demanding calculations. In a recent comparative study, the first molecular simulations were reported for an outer membrane protein using an accurate microbial membrane model. Contrary to expectation, the trans-membrane protein simulated in this more complex membrane model was found to be in a significantly more open conformation compared to simulations that used simpler membrane. This approach is now used to study a range of important trans-membrane proteins, such as the iron transporter FpvA, illustrated in figure 4 (p15).

With modern massively parallel computers and sophisticated, highly-efficient software for computational chemistry modeling, the theoretical study of complex and large systems is now feasible.

T. P. Straatsma, PNNL

Figure 4. FpvA is an outer membrane signaling and transport protein of Pseudomonas aeruginosa that binds Fe-bound siderophore pyoverdine for transport of the iron to the periplasmic space through the channel created by the protein in the lipopolysaccharide membrane.

Ab Initio Molecular Dynamics
The ability to predict the properties of complex materials important in toxic waste disposal, disease treatment, efficient chemical processing, and electronic device performance optimization, among others, is of great importance to DOE’s efforts to address the nation’s energy and environmental problems. Because the required properties are highly sensitive to complex interactions at the fundamental electronic structure level (for example, chemical bond saturation, and shell structure), reliable parameter-free simulation of their properties requires methods based directly on the solution to the electronic structure problem posed by the electronic Schrödinger equation. Development of methods at the fundamental electronic structure level also is important in light of DOE’s investment in large-scale facilities such as the synchrotron light and neutron sources. These new probes are providing an unprecedented level of detail at the atomic and molecular scale. However, without appropriate theories or models, many of these new measurements cannot be readily interpreted.

The generally expensive cost for approximations to the electronic Schrödinger equation has become manageable with the development of massively parallel computers and the development of new parallel algorithms and software.

The method of ab initio molecular dynamics (AIMD) enables researchers to treat the dynamics of these systems while retaining a first-principles-based description of their interactions. This approach is similar to classical molecular dynamics where the motions of the atoms and molecules are simulated over a period of time, but here the interactions between the atoms are calculated directly from the electronic Schrödinger equation, rather than from empirical interaction potentials or force fields. That is, at every step of a molecular dynamics simulation the locations of the electrons in the atoms and molecules are determined by solving a suitable approximation to the electronic Schrödinger equation. The instantaneous forces on the atoms are then determined by calculating their electrostatic forces from their interaction with other ions and electrons in the system.

This type of simulation requires an enormous amount of computational power. The typical time step in an AIMD simulation is quite small (~0.1 femtosecond) and the simulation needs to run for at least 10 picoseconds. Many chemical processes of interest occur on the order of nanoseconds (10^-9 seconds). Even for a 10 picosecond AIMD simulation, at least 10^-11/10^-16=10⁵ evaluations of the electronic Schrödinger equation are needed. In order for this to be a practical method, the solution to the electronic Schrödinger equation in a single time-step must be able to complete within seconds.

In general, for systems beyond a few atoms, even the least expensive approximations to the electronic Schrödinger equation are expensive to calculate. This cost has become manageable with the development of massively parallel computers and the development of new parallel algorithms and software. Using highly scalable software such as NWChem on 10,000 CPUs, AIMD simulations with 1,000 atoms for tens of picoseconds can readily be carried out on today’s systems (figure 5).

Illustration: A. Tovey Source: E. J. Bylaska, PNNL

Figure 5. Parallel timings per step of a NWChem ab initio molecular dynamics simulation for UO₂²⁺ + 122 H₂O on the Franklin system at NERSC.

Density Functional Theory is the most popular approximation to the electronic Schrödinger equation used for AIMD today. Higher levels of approximation to the electronic Schrödinger equation are needed, however, because the level of Density Functional Theory currently used in AIMD simulation software is unable to reliably predict the properties of many materials in basic research. Examples of interest to DOE include charge localization in transition elements with tightly bound d electrons in oxide materials (figure 6), the underestimation of reaction barriers and band gaps in solids, and accurate predictions of spin structure of solids and nanoparticles. In principle, with the advent of new parallel machines, which are 100 to 1,000 times more powerful than current machines, these extra computational costs could easily be overcome. However, algorithms must be constantly upgraded to capture the performance of these emerging massively parallel computers.

E. J. Bylaska, PNNL

Figure 6. Illustration of a localized electron (that is, a polaron) on the surface of hematite calculated with a higher and more expensive level (for example, hybrid DFT) of ab initio molecular dynamics. Lower levels of ab initio molecular dynamics will predict a delocalized electron.

An example of the use of AIMD simulations is the study of hydrated radionuclides under extreme conditions. A major obstacle to the development of nuclear power is the ability to safely store highly dangerous waste materials containing uranium and other radionuclides. Most current storage strategies are designed to store waste in saturated and unsaturated geological formations. Stored in this way, the most likely means for uranium to migrate into the biosphere is through groundwater contact with containment canisters, resulting in a solvated UO₂²⁺ cation (or complexes). To reliably predict behavior of the radionuclide (examples include uranium, thorium, and plutonium) waste products of nuclear power production over the range of conditions encountered in a storage facility requires a theory based at the most fundamental level on the electronic Schrödinger equation. For example, a detailed understanding of chemical processes occurring in aqueous solutions, as well as a better understanding of the surrounding hydration shells, for this cation is needed. In particular, the bond pattern and dynamics of the second hydration shell is critical to the interpretation of the ion-association properties and reaction chemistry of the UO₂²⁺ cation in solution.

Algorithms must be constantly upgraded to capture the performance of emerging massively parallel computers.


Illustration: A. Tovey Source: E. J. Bylaska, PNNL
Figure 7. Top, the experimental and simulated extended X-ray absorption fine structure spectra of UO₂²⁺ in water are in perfect agreement. Bottom, a snapshot of the inner solvation shell UO₂²⁺ is shown. The blue surface identifies the inner-coordination spheres and golden lines show the array of hydrogen bonds that are formed in the structure. The results shown are from an NWChem ab initio molecular dynamics simulation of UO₂²⁺ and 122 H₂O molecules.

In figure 7, results from an AIMD simulation of the UO₂²⁺ cation in aqueous solution are shown. The UO₂²⁺ cation was found to have a first solvation shell composed of five water molecules bonded to the uranium atom along its equatorial plane, and a second solvation shell of 10 water molecules hydrogen bonded to the first shell water molecules. In addition, about five to six water molecules on average were found to sporadically hydrogen bond to the oxygen atoms of the UO₂²⁺ cation in the apical region. The second and apical solvent regions were found to be very dynamic with many water transfers into and out of the equatorial and apical second solvation shells occurring on a picosecond time scale via dissociative mechanisms. Beyond these shells, the bonding pattern substantially returned to the tetrahedral structure of bulk water.

Even though the UO₂²⁺ cation has been studied extensively over the years using a variety of static ab initio (for example, static coupled cluster calculations) and classical molecular dynamics methods, these simulations have been either incomplete or inaccurate. Static ab initio simulations have been incomplete because they were not able to take into account the motion of the water molecules in the second and apical solvent shells, and classical molecular dynamic simulations have been plagued by inaccurate force fields of the complex interactions between the UO₂²⁺ cation and water. AIMD simulations, which are able to take into account complex interactions and dynamics, are able to overcome the well-known deficiencies of other molecular simulation methods. In fact, this AIMD simulation was the first molecular simulation able to reproduce the measured Extended X-ray Absorption Fine Structure (EXAFS) spectrum from experiments.

Very Accurate Solution to the Schrödinger Equation: The Coupled Cluster Formalism
Many aspects of computational chemistry require accuracies that can only be obtained by computational methods that appropriately account for instantaneous interactions between electrons within molecules. Including these so-called electron correlation effects is needed to enable precise comparison of theory and experiment. These effects are a key component in being able to understand important aspects of molecular structure, bond formation, and molecular interactions. For this reason these methods have become an integral part of many computational chemistry packages, including NWChem.

Among the methods that describe the correlation effects, the coupled cluster (CC) formalism has evolved into a widely used and accurate method for solving the electronic Schrödinger equation, which enables us to calculate molecular geometries, reactivity, and spectral properties of molecules. Compared with other methods, the CC formalism’s main advantage lies in the fact that the correlation effects are efficiently encapsulated in a special form of the wave function, enabling researchers to describe the collective motion of electrons. A simple consequence of this is a proper scaling of energy with the number of electrons. This feature is essential for describing chemical reactions. In order to illustrate and analyze the various effects, diagrams (figure 8) have become commonplace in the CC community.

Illustration: A. Tovey Source: K. Kowalski, PNNL

Figure 8. The use of diagrammatic techniques significantly alleviates the analysis of quantum many-body effects. The above diagrams represent connected diagrams contributing to the CC correlation energy (green dot vertices represent one- and two-body interactions between electrons while ovals refer to one- and two-body cluster amplitudes).

Although initially proposed in nuclear physics, CC formalism was quickly adopted by quantum chemists and, since the late 1960s, steady development has resulted in a variety of CC methodologies. Notably, this formalism has been readopted by the nuclear physics theory community in the last decade, demonstrating the applicability of the method across a range of length and energy scales.

The cost of CC calculations grows polynomially with the increase in system size. This is due to the character of the many-body effects included in a given approximation. In general, the methodologies, viewed as a compromise between accuracy and numerical cost, scale as N⁷, where N is the number of electrons, which means that the calculation for two water molecules is 2⁷=128 times more expensive than an analogous calculation for a single water molecule. Therefore, in order to perform calculations for systems composed of 40–60 atoms, the computational power of massively parallel computers must be harnessed effectively. Current CC capabilities in NWChem will scale up to 10,000 CPUs, sufficient to treat systems with 200–400 correlated electrons. Almost perfect scaling and perfect load balancing was achieved in large-scale property calculations using a genuine CC approach (figure 9).

Illustration: A. Tovey Source: K. Kowalski, PNNL

Figure 9. Parallel performance of the CC code in linear response calculations for the C₆₀ molecule (1,080 basis set functions, 240 correlated electrons).

To push the system-size limits further by utilizing emerging petascale architectures, a tight collaboration between teams of developers and computer scientists will be necessary to help define a new programming paradigm for CC approaches. Pressing problems include improving massive data movement across a large number of CPUs, data replication, trading computation for communication, and redefining addressing space. When left unresolved, these problems quickly become bottlenecks, hampering further progress in the parallel CC calculations. This effort should be conducted in parallel with development of automatic code generation of CC codes, which would eliminate the error-prone process of manual coding of numerous diagrams encapsulated within the CC equations, as well as tailor and tune the code to specific architectures.

Applications of the CC methodology in chemistry are numerous and have proven effective. The character of the CC formalism makes it applicable to a wide spectrum of systems at different levels of complexity. These systems range from small molecules to large and highly-delocalized molecular systems such as fullerenes. With the CC formalism, the interaction of light with matter also can be studied, using the basic principles of quantum mechanics. The use of CC methods in the context of light-induced processes involved in conversion of solar energy, interaction of biological systems with radiation, or imaging of large molecular assemblies, will require a balanced development of accurate CC methods capable of characterizing complex excited states and methods for calculating non-linear optical properties that describe the behavior of molecules in external electric fields. Various excited-state extensions of CC formalisms, where the electrons are excited from the lower to higher energy levels, are ideally suited for this purpose. Using relatively simple models, researchers can describe excited states obtained by exciting a single electron while more sophisticated approaches are needed to tackle more challenging excited states where the excitation of two electrons occurs. Excited states change their character with changes in nuclear geometry—a common feature of many light-induced reactions—and very flexible methods are needed to preserve the same level of accuracy along the reaction coordinates.


K. Kowalski, PNNL
Figure 10. High-level ab initio CC calculations were performed for excited states of DNA bases and zinc–porphyrin systems.

Two such applications are displayed in figure 10. One application focuses on solving a problem of remarkable photo-stability of DNA bases, usually attributed to ultrafast conversion from the excited to the ground state. This mechanism protects all living organisms against potentially deadly mutations caused by UV radiation. Early studies suggested several scenarios, with some predicting the involvement of several excited states in this process. Intriguing features of DNA bases were revealed when quantum mechanical calculations are performed in a realistic environment. Excited-state calculations for cytosine and guanine in their native DNA or aqueous environments showed that the inclusion of environment widens the gap between the first and second excited states compared to findings in the gas-phase. This finding suggests that the one-state mechanism is favored in realistic settings. The inclusion of higher-order effects widens this gap even further, which favors the one-state model. These findings clearly indicate that the ultrafast conversion mechanism depends strongly on the structure of the surrounding environment. Using the hybrid CC/classical molecular dynamics approach, for the first time ever, estimates could be made of the impact of thermal effects on the spectrum of cytosine. At finite temperatures, the excitation energies corresponding to the two lowest singlet states are significantly blue shifted. This simulation took 105 hours utilizing 256 processors, a task that would have taken three years to complete on a serial machine.

The CC methodology also will become important in studies of fluorescent probes used in protein imaging and monitoring physiological processes in biological cells. These probes undergo changes in spectroscopic properties in response to surrounding electric field. As such, they can be used to non-invasively sense electrostatic effects at the molecular level. In order to characterize this response properly, a detailed knowledge of non-linear optical properties and excitation energies is necessary. Unfortunately, due to their size, high-accuracy results are not currently available. However, the emergence of petascale architectures will improve this picture. Recently, high-level calculations of low-lying excited states of di-8-ANEPPS fluorescent probes were conducted and the results indicate that the ground state has a complicated structure whereas the first excited state is characterized by charge transfer occurring at a remarkably long distance. Remarkably, our CC estimates for excitation energies obtained with medium size basis sets are significantly different (even by 1 eV) from those previously obtained with lower-level methods. This calculation was performed on 1,024 CPUs and took 11 hours, which means that even larger calculations with very reliable basis sets are possible.

Light harvesting (photosynthesis) typically occurs via a two-component system, that is, the light harvesting molecular array (antenna) that increases the amount of solar energy absorbed and is subsequently passed through energy transfer steps to the main system, where the charge separation occurs. Covalently bonded molecular systems that mimic the photosynthesis and charge-separation processes have been synthesized but the development of self-assembled functional antennas is at an early stage. In this regard, theoretical studies of these systems’ components can significantly advance the field, eventually leading to a discovery of efficient light-harvesting systems. Because the zinc–porphyrin (ZnP) molecule (figure 10, bottom) is a common building block of many donor–acceptor light-harvesting systems and can be used as a potential candidate for light-harvesting antenna, characterizing its excited states with high accuracy is vital. Due to its large size, however, calculations of ZnP with high-level methods are time consuming. The effective use of thousands of CPUs allows the time to solution to be compressed to a few hours, demonstrating that routine calculations for systems of this size are already possible. These codes can be used as a component of predictive models of how charge separation occurs in excited states and how effectively this process can be transformed into electricity.

Computational chemistry is on the verge of entering a new era of modeling and simulation. Codes such as NWChem are contributing to this new era, allowing researchers to tackle scientific problems that are larger and more realistic than ever before.

Looking to the Future
Developers continuously improve the software to keep up with the evolving needs of researchers around the world and to match the continuous increases in performance afforded by advancing supercomputing resources (sidebar “Software Development for Massively Parallel Architectures”).

With the current petascale computing platforms, computational chemistry is on the verge of entering a new era of modeling and simulation. Codes such as NWChem are contributing to this new era, allowing researchers to tackle scientific problems that are larger and more realistic than ever before. This development is providing insight into the complex dynamical behavior of nature, and enables researchers to address new and more challenging scientific questions.

Contributors T. P. Straatsma, E. J. Bylaska, and K. Kowalski