Materials and Process Simulation Center,
California Institute of Technology
Title of the Effort
Highly Parallelized Large Scale Atomistic Simulations
for Design of Materials
The objective is to develop theoretical methodologies
for practical computations of the structures and
properties of real materials that can be used in industrial process
design for manufacturing new materials.
William A. Goddard III,
Materials and Process Simulation Center
Department of Chemistry, Beckman Institute (139-74),
Caltech, Pasadena, CA 91125
Center for Biomolecular Simulation
Department of Chemistry, Columbia University, New York, NY 10027
Zhen-Gang Wang, firstname.lastname@example.org
Department of Chemical Engineering (210-41),
Caltech, Pasadena, CA 91125
Abhinandan Jain, email@example.com
Robotics, Jet Propulsion Laboratory (198-219),
Pasadena, CA 91109
Stephen Taylor, firstname.lastname@example.org
Scalable Parallel Computing Lab
Department of Computer Science (256-80),
Caltech, Pasadena, CA 91125
The full 1995 report of the Materials and Molecular
Simulation Center in the Beckman Institute at Caltech is
also available on the World Wide Web (WWW).
includes the detailed progress reports on the NSF-GCAG project along with
projects supported by other agencies and by seven industrial
There is an enormous gap between current methodologies for
atomistic simulations and the level required for accurately describing
the relevant properties of industrially important materials. Our
strategy is to transcend from the most fundamental theory
(quantum mechanics, QM) to practical engineering designs in a sequence of
four or five levels as indicated in the following figure.
- 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
- 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
- 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
Our team of 5 PI's has made significant progress on each of
the six tasks described under APPROACH. The new theoretical
methods developed for quantum mechanics and molecular dynamics that scale
sufficiently slowly with size so as to be practical for the very large
systems of industrial interest. Simulations of hundreds of atom systems for
quantum mechanics and tens of millions of atoms for complex systems
within the reach.
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
- treat all atoms of the
periodic table (using core effective potentials),
- handle new
sophisticated wavefunctions (GVB-RCI, MP2), and
- describe important properties (solvation
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
The progress in force fields is twofold.
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
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
- 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
- High frequency modes may be eliminated through rigid body dynamics with
- 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.
The first release/production version of the program is now available and
- used in studying
the thermodynamcis of
large scale atomic clusters which has profound effect in microcluster
- being used in studying large viruses: HRV-14(0.5M atoms), and TBV
We have focussed on
developing systematic method for simulating mesoscopic and
macroscopic structural and thermodynamic properties of supramolecular
structures. This method enables us to efficiently sample phase space of
synthetic polymers and biopolymers to construct the partition function
from which the structural and thermodynamic properties may be obtained.
This method is being applied to dilute polymers melts, bulk polymer
systems, folding of medium size proteins.
Monte Carlo Program originally written in fortran is now being optimized
and rewritten in an object oriented fashion in C for compliance with
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
- 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.
- 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
- 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.
- 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
- 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
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
- 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.
- "Dielectric Properties of Poly(vinyledene fluoride) from Molecular
N. Karasawa and W. A. Goddard, III., Macromol. 6785 (1995).
- "Hessian Biased Force Field for Polysilane Polymers," C.B. Musgrave, S. Dasgupta and W.A.
Goddard, III, J. Phys. Chem. 99 13321 (1995).
- "Ab initio predictions of large hyperpolarizibility
Julolidinyl-n-isoxazolone and julolidinyl-n-N,N'-diethylthoabarbituric
acid", Daqi Lu, B.
Marten, Y. Cao, R. A. Freisner, and W.A. Goddard,
Chem. Phys. Lett. 242, 543, (1995).
- "Dual-space approach for density-functional calculations
of two- and three-dimensional crystals using Gaussian basis functions,"
X. Chen, J-M. Langlois, and W.A. Goddard
Phys. Rev. B 52, 2348 (1995)
- "Effects of Catalyst Promoters on the Growth of Single-Layer
C-H. Kiang, W. A. Goddard III, R. Beyers, J. R. Salem, and
D. S. Bethune
Mat. Res. Soc. Symp. Proc. 359, 69 (1995)
- "First principles studies of band offsets at
heterojunctions and of surface reconstruction using Gaussian dual-space
density functional theory,"
X. Chen, A. Mintz, J. Hu,
X. Hua, and W.A. Goddard,
J. Vac. Sci. Technol. B13, 1715 (1995)
- "Building proteins from C-alpha coordinates using the dihedral
probability grid Monte-Carlo method," A. M. Mathiowetz and W. A. Goddard III,
Protein Sci. 4 , 1217-1232 (1995)
- "De novo prediction of polypeptide conformations using dihedral
probability grid Monte-Carlo methodology,"
J. S. Evans, A. M. Mathiowetz, S. I. Chan, and W. A.
Protein Sci. 4 , 1203-1216 (1995)
"Stabilizing the boat conformation of cyclohexane rings,"
S. Dasgupta, Y.C. Tang, J.M. Moldowan, R.M.K. Carlson and W.A. Goddard,
J. Am. Chem. Soc. 117 , 6532-6534 (1995)
- "Prediction of new donors for organic superconductors,"
E. Demiralp and W.A. Goddard,
Synt. Metals 72, 297-299 (1995)
- "Design and synthesis of a new peptide recognizing a specific 16-base-pair
site of DNA,"
C. Park, J.L. Campbell and W.A. Goddard,
J. Am. Chem. Soc. 117, 6287-6291 (1995)
- "Is carbon nitride harder than diamond - No, but its girth increases when
stretched (negative poisson ratio),"
Y.J. Guo and W.A. Goddard,
Chem. Phys. Lett. 237, 72-76 (1995)
- " Experimental and theoretical studies of CO(CH4)(X)(+) with X=1-4,"
C.L. Haynes, P.B. Armentrout, J.K. Perry and W.A. Goddard,
J. Phys. Chem. 99, 6340-6346 (1995)
- "Free energy and surface-tension of arbitrarily large Mackay icosahedral
R.B. Mcclurg, R.C. Flagan, and W.A. Goddard,
J. Chem. Phys. 102, 3322-3330 (1995)
- "Parallel Calculation of Electron-Transfer and Resonance
Matrix Elements of Hartree-Fock and Generalized Valence Bond Wavefunctions,"
E. P. Bierwagen, T. R. Coley, and W. A. Goddard III
Am. Chem. Soc., 1994 ACS Symposium Series on Parallel
Computing in Computational Chemistry, T. G. Mattson editor, Chapter 7,
pg. 84 (1995).
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.
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
- 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.
Schrodinger Inc., Pasadena CA (contact Dr. Murco Ringnalda, ph:
818-568-9392, fax: 818-568-9778, email: email@example.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
single crystal films.
Asahi Chemical Corp. Japan (Terumasa Yamasaki email:
firstname.lastname@example.org, 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: email@example.com) uses MD-CMM to simulate polmer properties.
Xerox Palo Alto Research Corp. (Dr. Ralph Merkle, email:
firstname.lastname@example.org) 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.
Last modified on October 6, 1995.
Problems with this server should be reported to email@example.com