Highly Parallelized Large Scale Atomistic Simulations for Design of Materials

NSF-GCAG (ASC 9217368) 1995 REPORT


Title of Effort
Principal Investigators
Related Information


[Quantum Mechanics] [Force Fields] [Molecular Dynamics]
[Coarse Graining] [Statistical Mechanics] [Continuum Parameters]


Quantum Mechanics
Force Fields
Statistical Mechanics

Significant Events and Accomplishments

Prediction of Diffusion of Gases in Polymers
Dielectric Properties of Polymers
Band offsets of heterojunctions
Large hyperpolarizibility push-pull polymers
Prediction of new donors for organic superconductors
Mega Molecular Dynamics Applications
MPSim Multi Platform Production versions


Fiscal '96 Plan

Technology Transition



Materials and Process Simulation Center,
Beckman Institute,
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.

Principal Investigators

William A. Goddard III, wag@wag.caltech.edu
Materials and Process Simulation Center
Department of Chemistry, Beckman Institute (139-74),
Caltech, Pasadena, CA 91125

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

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). This report includes the detailed progress reports on the NSF-GCAG project along with projects supported by other agencies and by seven industrial collaborators.


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.

Modeling Hierarchy


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.



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.


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:

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

The first release/production version of the program is now available and
  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)


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 MPSim.

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.


    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.

    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.

    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.

    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.

    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.

    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.


  1. "Dielectric Properties of Poly(vinyledene fluoride) from Molecular Dynamics Simulations," N. Karasawa and W. A. Goddard, III., Macromol. 6785 (1995).

  2. "Hessian Biased Force Field for Polysilane Polymers," C.B. Musgrave, S. Dasgupta and W.A. Goddard, III, J. Phys. Chem. 99 13321 (1995).

  3. "Ab initio predictions of large hyperpolarizibility push-pull polymers. 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).

  4. "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)

  5. "Effects of Catalyst Promoters on the Growth of Single-Layer Carbon Nanotubes,"
    C-H. Kiang, W. A. Goddard III, R. Beyers, J. R. Salem, and D. S. Bethune
    Mat. Res. Soc. Symp. Proc. 359, 69 (1995)

  6. "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)

  7. "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)

  8. "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. Goddard III Protein Sci. 4 , 1203-1216 (1995)

  9. "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)

  10. "Prediction of new donors for organic superconductors,"
    E. Demiralp and W.A. Goddard, Synt. Metals 72, 297-299 (1995)

  11. "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)

  12. "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)

  13. " 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)

  14. "Free energy and surface-tension of arbitrarily large Mackay icosahedral clusters,"
    R.B. Mcclurg, R.C. Flagan, and W.A. Goddard, J. Chem. Phys. 102, 3322-3330 (1995)

  15. "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).



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.

Last modified on October 6, 1995.
Problems with this server should be reported to webmaster@www.wag.caltech.edu