Highly Parallelized Large Scale Atomistic Simulations for Design of Materials

Approach

Progress

Significant Events and Accomplishments

Publications

Fiscal '96 Plan

Technology Transition

INTRODUCTION

Organization

Title of the Effort

Highly Parallelized Large Scale Atomistic Simulations for Design of Materials

Objective

Principal Investigators

Richard Friesner, rich@cucbs.chem.columbia.edu
Center for Biomolecular Simulation
Department of Chemistry, Columbia University, New York, NY 10027

Zhen-Gang Wang, zgw@macpost.caltech.edu
Department of Chemical Engineering (210-41),
Caltech, Pasadena, CA 91125

Abhinandan Jain, jain@telerobotics.jpl.nasa.gov
Robotics, Jet Propulsion Laboratory (198-219),
Pasadena, CA 91109

Stephen Taylor, steve@cs.caltech.edu
Scalable Parallel Computing Lab
Department of Computer Science (256-80),
Caltech, Pasadena, CA 91125

Related Information

APPROACH

• Quantum Mechanics (QM)
  First principles solution of the Schrodinger equation (H F = E F) leads to accurate predictions of properties (structures, chemical reactions, excited states); however, this was limited to small systems (10 to 20 atoms) limiting the time and distance scale to 1 ps and 1 nm.

• Force Fields (FF)
  By averaging over the electrons from QM we can obtain parameters (charges, force constants) while allowing materials to be described in terms of atoms rather than electrons.

• Molecular Dynamics (MD)
  Using force fields (FF) the fundamental equations become Newton's equations (F = MA) rather than the Schrodinger equations. This allows practical calculations on 1000 to 5000 atoms rather than 10 to 20, extending the time and distance scale to 1 ns and 10 nm.

• Coarse Graining (CG)
  By averaging over the atoms for MD, we can obtain parameters representing groups of atoms (molecules, segments), considerably simplifying the calculations.

• Statistical Mechanics (SM)
  Using the CG description we can examine materials in terms of the large scale motions relevant to macroscopic experiments. This extends the time and distance scales to 1 microseconds and microns.

• Continuum Parameters (CP)
  The results from SM are combined to obtain macroscopic or continuum parameters (free energies, phase diagrams, partition coefficients, solubility parameters) suitable for practical chemical engineering software (e.g., ASPEN-PLUS) for design of unit processes.

PROGRESS

QUANTUM MECHANICS AND QUANTUM CHEMISTRY

• QM - FINITE.
  For finite molecules; A new methodology (PS-GVB) combining pseudospectral (PS) multi-grid and dealiasing strategies with the sophisticated many-body wavefunctions [generalized valence bond (GVB)] were implemented and applied to large scale problems. PS-GVB led to considerably better scaling with size (N^2 rather than the N^3, N^4, N^5, N^6 characteristic of alternative methods) and simpler parallelization. parallelization. The basic PS-GVB program has been commercialized. PS-GVB has been extended to
  1. treat all atoms of the periodic table (using core effective potentials),
  2. handle new sophisticated wavefunctions (GVB-RCI, MP2), and
  3. describe important properties (solvation energies, hyperpolarizabilities).

  Next we will (1) optimize these methods for parallel implementations, (2) extend the methodology to include GVB-RCI-MP2 and self-consistent GVB-RCI.

• QM - INFINITE.
  Most practical materials properties require a description of infinite systems using periodic boundary conditions (PBC). This is three-dimensional (3D) for bulk properties or two-dimensional (2D) for surface growth and interfaces. For this purpose we have developed a new method, Gaussian dual space density functional theory (GDS-DFT), in which most parts scale linearly with N. In implementing this we have developed a new separable pseudopotential that can be applied to all atoms of the periodic table.

  We reformulated the theories for electronic structure calculations of periodic systems in a way suitable for large scale calculations using Gaussian basis functions. An accurate grid is introduced for efficient calculation of matrix elements. A dual-space approach is used to calculate the Coulomb potential with computational cost that scales linearly with the size of basis set. A preconditioned generalized conjugate gradients approach is introduced for rapidly converging wave functions expressed in terms of Gaussian basis functions. This method is applied to a variety of crystals (including diamond, GaN, AIN, CdTe, and C-60) and surfaces [including GaAs (110) and BN (110)] with excellent results.

  The GDS-DFT program is now being used for several projects. The program optimized for Cray C90. For next stage we will optimize it for parallel environments; KSR, SGI/PC, Intel/Paragon, Cray T3D and IBM SP2.

FORCE FIELDS

Over the past year new specific force fields were developed utilizing ab initio quantum chemistry programs with force field optimization programs. These led to development of force fields for poly ethylene, polysilane, nylons, polyoxymethylenes, polyvinyl (fluoride and chloride) polymers, semiconductors and ceramics. For instance, the Polysilane FF, denoted MSXX, was developed using the Hessian biased method to describe accurately the vibrational states, the ab initio torsional potential energy surface, and the ab initio electrostatic charges of polysilane oligomers. This MSXX FF was used to calculate various spectroscopic and mechanical properties of the polysilane crystal. Stress-strain curves and surface energies are reported. Gibbs molecular dynamics calculations (Nose, Rahman-Parrinello) were used to predict various materials properties at higher temperatures. Phonon dispersion curves and elastic constants were calculated at various temperatures.

The Singular Value Decomposition method (SVD) has been implemented for fast optimization of force fields using a variety of experimental and theoretical data. In most instances, the experimental data is incomplete or unavailable for a fully contrained force field development. The SVD method provides an easy and automatic for force field development for cases where parameters are not readily available or where higher accuracy is necessary. Energy and Force evaluator module of massively parallel program is enhanced to handle specific forces necessitated by these force field developments.

MOLECULAR MECHANICS AND MOLECULAR DYNAMICS

The focus here is on extending the methods of MD to physical systems of molecules, polymers, liquids, biopolymers and inorganic materials with more than 10 million atoms per simulation cell while accurately treating long-range interactions. The new methodology is the cell multipole method (CMM). The CMM methodology has been efficiently implemented into a general MD code on the KSR and Intel/Paragon parallel supercomputers (part of the Caltech NSF-GCAG project). Both versions are in production mode. Simulations carried out on both platforms for challanging materials and biotechnology problems.

For fast internal coordinate dynamics on a million atoms, we have developed the Newton-Euler Inverse Mass Operator (NEIMO) method. This methodology has been extended to handle periodic systems and has been incorporated into the advanced dynamics (Gibbs/Nose) codes for single processor computers. Currently the technology is being merged into MPSim Program environment using shared memory platform as first test bed. The applications for this method will enable us to investigate the long time dynamics of liquid polymer and solid interfaces; which has tremendous impact on broad range of chemical industry applications: wetting, adhesion, phase separation, coatings, dyes, etc.

Program is developed originally on KSR in an object oriented fashion to be portable and efficient for other parallel computers. KSR version serves as the basis for ports to parallel computers using a global memory model. Developments using MPI on Cray T3D, IBM SP2 and SGI/PC is in the beta testing stage.

Starting with the th KSR code a communication layer was added for distributed memory applications based on the active message strategy for the J-Machine. This has been implemented and tested on the Intel Touchstone Delta/512 and Intel/Paragon at Caltech.

We now can do routine molecular dynamic calculations up to 1.5 million atoms with the KSR and Intel versions, and now are in production for industrially important problems (diffusion, glass temperature) for realistic description for amorphous polymers. A recent project is launched in studying the aggregation of asphaltenes in crude oil which presents great challenge to petrochemical industry using large scale atomistic simulation of oil, asphaltene and water interface. This constitutes a great leap forward, allowing accurate atomistic simulations on systems 40 times larger than previous methods.

Over the last few months, the molecular mechanics and molecular dynamics program: Massively Parallel/Materials&Process Simulation and Analysis (MPSim) Program System is being find tuned for optimal performance on various high performance computers:

• Shared memory MPSim on KSR and SGI Workstations;
• Message Passing MPSim on Intel using NX Express libraries;
• MPI MPSim on CRAY T3D, SGI Power Challenge, and IBM SP2

The program is the most general Molecular Mechanics and Dynamics program available on massively parallel platforms:

• It can model any system: Periodic and nonperiodic structures, polymers, peptides, nucleic acids, inorganic materials and organo-metallic systems
• It uses generic force fields developed over the past ten years: Dreiding II, Universal Force Fields, AMBER Force Fields, CHARMM Force Force Fields and MM3 Force Fields.
• It can handle special force field developed using Biased Hessian and FFOPT generated (Force fields generically known as MSXX)
• Long range interactions are uniquely and efficiently handled using cell multipole methods (CMM and RCMM)
• A dynamic load balancing method is used to increase the efficiency.
• Robust Minimization Techniques are implemented
• Advanced Extended System Molecular Dynamics methods are included; Nose-Hoover Canonical dynamics and Parinello-Rahman Constant Pressure dynamics.
• High frequency modes may be eliminated through rigid body dynamics with Quaternion representation.
• Mixed mode dynamics for Polymers or biopolymers in solvents: rigid body representation for solvents and stiff groups and fully flexible representation for soft segments.
• NEIMO mechanics and dynamics: First consistent implementation of generalized coordinate representation of many body newtonian mechanics and dynamics methods developed for space science application to molecular mechanics and dynamics of polymer and biopolymer simulations.
• Trajectory and structure analysis tools for simulations performed on high performance computing platforms
• Visualization tools for trajectory and structure analysis.
• A graphical user interface running on desktop workstations to prepare job control and job launching to high performance computers
• On line documentation of the program and utilities.

1. used in studying the thermodynamcis of large scale atomic clusters which has profound effect in microcluster physics.
2. being used in studying large viruses: HRV-14(0.5M atoms), and TBV (3M atoms)

STATISTICAL MECHANICS

A new Monte Carlo method is developed and applied in estimating diffusion of gases in glassy polymer networks. This statistical method utilizes millions of step molecular dynamics (nanoseconds long) run to explore a phenomena which is orders of magnitude larger than the time scales accessible in Molecular Dynamics.

SIGNIFICANT EVENTS AND ACCOMPLISHMENTS

1. PREDICTION OF GAS DIFFUSION IN POLYMERS.
   A significant industrial application for atomistic simulations is the prediction of diffusion of gases in glassy polymers. We are developing a new approach where atomistic Molecular Dynamics is performed on a large polymer system to extract compressibilities and the dynamical behaviour of voids. This is then used in a Monte Carlo program to run long trajectories for penetrant gas molecules. In the polymers studied so far (PolyEthylene and PolyVinylChloride, PolyVinylyDeneChloride) the qualitative trends and the quantitative diffusion constants agree well with experimental data. This general method can be applied to any polymer system with any penetrant diffusive molecule.

2. PREDICTION OF DIELECTRIC PROPERTIES OF POLYMERS
   Molecular dynamcis simulations at 300 and 400 K were used to calculate the dielectric constants and dielectric loss for crystalline and amorphous PVDF. These calculations used MSXX force fields developed from first principles and validated with experimental structural and mechanical properties of various forms of crystalline forms of PVDF. We showed that the high dielectric constant for amorphous PVDF arises from the rapid changes in the torsion angles during the Molecular Dynamics. These changes are enhanced by soliton-like defects in chain that diffuse during simulation. The dielectric losses lead initially to a stretched exponential decay function.

3. BAND OFFSETS OF HETEROJUNCTIONS
   The use of localized Gaussian basis functions for large scale first principles density functional calculations with periodic boundary conditions (PBC) in 2 dimensions and 3 dimensions has been made possible by using a dual space approach. This new method is applied to the study of electronic properties of II-VI (II=Zn, Cd, Hg; VI=S, Se, Te, Po) and III-V (III=Al, Ga; V=As, N) semiconductors. Valence band offsets of heterojunctions are calculated including both bulk contributions and interfacial contributions. The results agree very well with available experimental data. The p(2X1) cation terminated surface reconstructions of CdTe and HgTe (100) are calculated using the local density approximation (LDA) with two-dimensional PBC and also using the ab initio Hartree-Fock (HF) method with a finite cluster. The LDA and HF results do not agree very well.

4. LARGE HYPERPOLARIZABILITY PUSH-PULL POLYMERS.
   Recently significant advances have been made in engineering push-pull organic chromophores to have very large hyperpolarizabilities (beta), leading to materials with mu beta as high as 15000 X 10(-48) esu. Such developments have been slow and costy because of difficulties in synthesis, purification, and measurement. As an alternative we have developed a new quantum mechanical program (PS-GVB/NLO) which provides predictions of beta for such molecules far faster than previously possible. We have applied PS-GVB/NLO to predicting alpha, beta, and gamma for the high beta push-pull organics and find excellent agreement with experiment. This suggests that theory can be used as an effective tool for developing new nonlinear optical materials.

5. PREDICTION OF NEW DONORS FOR ORGANIC SUPERCONDUCTORS
   The donors of all known one- or two-dimensional organic superconductors are based on a core organic molecule that is either tetrathiafulvalene (denoted as TTF) or tetraselenafulvalene (denoted as TSeF) or some mixture of these two molecules. Coupling X, with appropriate accepters, Y, leads to superconductivity. The oxidized form of X may be X(+) or X(2)(+) species in the crystal. From ab initio quantum mechanical calculations (HF/6-31G(**)), we found that all known organic superconductors involve an X that deforms to a boat structure while X(+) is planar. This leads to a coupling between charge transfer and the boat deformation phonon modes that we believe is responsible for the superconductivity of these materials. Based on this idea we have developed similar organic donors having the same properties and suggested that with appropriate electron accepters they will also lead to superconductivity.

6. MPSim MULTI PLATFORM PRODUCTION VERSIONS
   In early fall, fully functional first production version of Massively Parallel Materials and Process Simulation (MPSim) program is completed. This version is a fully functional molecular mechanics and molecular dynamics program for chemical, biologic and materials problems. Its components are:
   • User Interface: Tcl/Tk based graphical user interface
   • Interface module: Input, output
   • Force Field Evaluator module: Core of the program
   • Simulators module:
     • Energy - Force Evaluator
     • Molecular Mechanics: Structure Optimizer
     • Molecular Dynamics : NVE, NVT dynamics
   • Analysis and Utility programs
   • Documentation: Currently online, and tested by MSC members and associates.

Publications

FY-95 PLANS

• PS-GVB.
  Parallelize the solvate, hyperpolarizability, and GVB-RCI methods. Further develop and optimize the MP2, self-consistent GVB-RCI, and GVB-RCI-MP2 methods. Apply to important industrial problems.

• GDS-DFT.
  Parallelize for the CRAY T3D. Apply to important industrial problems.

• Constant Stress-Pressure Mechanics and Dynamics
  Generalize the constant pressure and stress algorithms to handle multi-mode mechanics and dynamics (for flexible, rigid and NEIMO threads) Reformulate and fine tune the expensive stress calculations.

• Two dimensional periodicity for Interface Problems
  Develop the interface simulation module to study extremely important industrial processing and manufactruring problems revealing themselves as interface issues: Wetting, friction, adhesion, wear, corrosion and failure.

• Non Equilibrium Molecular Dynamics Module
  Most of the materials problems reveal themselves as systems far from equilibrium. Equilibrium molecular dynamics is insufficient to capture the essense of the properties of the systems under the influence of shear, or under thermal, pressure, concentration gradients; where system shows a flow behavior. Synthetic hamiltonian formulations for performing nonequilibrium molecular dynamics exist but applications have always been restricted to simple model systems. The large scale massively parallel platforms is making this step a feasible one. Especially, the underlying general force field evaluator core of MPSim is our strongest point to develop such a general NEMD program module.

TECHNOLOGY TRANSITION

Schrodinger Inc., Pasadena CA (contact Dr. Murco Ringnalda, ph: 818-568-9392, fax: 818-568-9778, email: info@psgvb.com) has licensed the PS-GVB technology from Caltech and Columbia. They have a user-friendly interface, extensive user guide, and extended the program in many ways. Over the past year they made two commercial releases. PS-GVB is also ported into Cerius2 interface (from Biosym/Molecular Simulations)

Chevron Petroleum Technology Co., La Habra, CA (Dr. Yongchun Tang, ph: 310-694-7550) uses PS-GVB, GDS-DFT, and MD-CMM to design new scale inhibitors and corrosion inhibitors for oil recovery applications (designing new and more effective calculations).

Allied-Signal, Morristown, NJ (Dr. Willis Hammond, ph: 201-455-4914i email: HammonB@mtomp201.research.allied.com) uses MD-CMM to predict structure and moduli of nylon polymers.

Hughes Research Labs., Malibu, CA (Dr. Jenna Zinck, ph: 310-317-5913) uses GDS-DFT to design metallorganics for MO-MBE deposition of HgCdTe single crystal films.

Asahi Chemical Corp. Japan (Terumasa Yamasaki email: yamasaki@cs.fuji.asahi-kasei.co.jp, Takeshi Aoyagi, Shozo Yamanaka) uses PS-GVB, MD-CMM to simulate a wide range of problems in catalysis, electronic structure and material properties.

Asahi Glass Corp. Japan (Takashige Maekawa, Hiroshi Yamamoto, Takashi Miyajima email: tam@agc.co.jp) uses MD-CMM to simulate polmer properties.

Xerox Palo Alto Research Corp. (Dr. Ralph Merkle, email: merkle@xerox.com) uses MD-CMM to simulate nanosystems and understand their vibrational properties.

Sandia National Labs (Drs. John Shelnutt and M. Cather Simpson email: MCSIMPS@energylan.sandia.gov) uses MD-CMM to obtain accurate force fields for designing of porphyrin based chiral and achiral catalysts for a variety of industrial problems.