CPMD: design and characterization of innovative materials

1. CPMD: design and characterization of innovative materials Mauro Boero Institut de Physique et Chimie des Matériaux de Strasbourg, UMR 7504 CNRS-UDS, 23 rue du Loess, BP 43, F-67034 Strasbourg, France and CREST, Japan Science and Technology Agency, Kawaguchi, Saitama 332-0012, Japan, and JAIST, Hokuriku, Ishikawa, Japan

2. Outline • CPMD: a quick overview of the basics of the code and advanced tools for simulating reactive processes reaction • Synthetic organic reactions - e-Caprolactam (Nylon-6) production without using acid catalysts: • - Catalytic properties of water above the critical point • - Tuning the efficiency and selectivity of the reaction

3. What do we want to do ? And which are main the ingredients ? yi(x) electron- electron electron - ion RI ion-ion electron- ion

4. Car-Parrinello Molecular Dynamics • Solve the Euler-Lagrange equations of motion

5. BO surface CP trajectory BO trajectory The difference between the CP trajectories RICP(t) and the Born-Oppenheimer (BO) ones RIBO(t) is bound by | RICP(t) - RIBO(t)| < C m1/2 (C > 0) if F.A. Bornemann and C. Schütte, Numerische Mathematik vol.78, N. 3, p. 359-376 (1998)

6. Plane wave basis set:yi(x) = SGci(G) eiGx For eachelectroni=1,…,N , G = 1,…,M are the reciprocal space vectors. The Hilbert space spanned by PWs is truncated to a cut-off Gcut2/2 < Ecut R space aG space E1cut E2cut > E1cut

7. G space R space N*FFT ci(G) yi(x) SG ci(G) G2{Ek} VNL(G) {ENL} r(G) FFT r(x) Vloc (G)+ VH(G) {Eloc+EH} = VLH(G) Vxc(x) {Exc} + VLH(x) = VLOC(x); VLOC(x)yi(x) FFT N*FFT VLOC(G)ci(x) + VNL(G) + SG ci(G) G2

8. Practical implementation • G=1,…,M (loop on reciprocal vectors) are distributed (via MPI/OMP) in a parallel processing in bunches of M/(nproc) • i=1,…,N(loop on electrons) is distributed (via MPI/OMP) • I=1,…,K(loop on atoms) generally does not require parallelization (vectorized and distributed via OMP) • The scaling of the algorithm is O(NM)for the kinetic term, O(NM logM) for the local potential and O(N2M) for the non-local term and orthogonalization procedure (all other quantum chemical methods scale as O(MN3) M=basis set) - http://www.cpmd.org - http://www.cscs.ch/~aps/CPMD-pages/CPMD/Download

9. ES system configuration • Parallel vector supercomputer system with 640 processor nodes (PNs) connected by 640x640 single-stage crossbar switches. • Each PN is a system with a shared memory, consisting of • 8 vector-type arithmetic processors (APs): total=5120 AP • a 16-GB main memory system (MS) • a remote access control unit (RCU) • an I/O processor. MPI

10. ES system : single processor node (PN) • The overall MS is divided into 2048 banks • The sequence of bank numbers corresponds to increasing addresses of locations in memory. OMP

11. D From reactants A to products B: we have to climb the mountain minimizing the time • A general chemical reaction starts from reactants A and goes into products B • The system spends most of the time either in A and in B • …but in between, for a short time, a barrier is overcome and atomic and electronic modifications occur • Time scale:

12. V(s) s Escaping the local minima of the FES: In one dimension, the system freely moves in a potential well (driven by MD). Adding a penalty potential in the region that has been already explored forces the system to move out of that region, but always choosing the minimum energy path, i.e. the most natural path that brings it out of the well. Providing a properly shaped penalty potential, the dynamics is guaranteed to be smooth and therefore the systems explores the whole well, until it finds the lowest barrier to escape.

13. F(s) 脱出 F(s)+V(s, t) ∙∙∙∙∙ t3 t1 t2 t0 sa(t0) sa Set up collective variables {sa} and parameters Ma, ka, Ds, A Perform few MD steps under harmonic restraint Add a new Gaussian Update mean forces on {sa} Update {sa} The component of the force coming from the gaussians subtracts from the “true” force the probability to visit again the same place

14. History-dependent potential How to plug all this in CPMD ?We simply write a (further) extended Lagrangean including the new degrees of freedom Fictitious kinetic energy Restrain potential: coupling fast and slow variables √(kα/Mα) «ωI

15. Collective (dynamical)variables Velocity Verlet algorithm to solve the equations of motion two contributions to the force

16. Beckmann rearrangement: • Commercially important for production of synthetic fibers • Known to be catalyzed only bystrong acidsin conventional non-aqueous systems • Formation ofbyproducts (ammonium sulfate, (NH4)2SO4)of low commercial value in acid catalyst: byproducts = 1.7 × products (in weight). See (e.g.)http://www.clarkson.edu/~ochem/Spring01/CM244/caprolactam.htmlhttp://es.epa.gov/p2pubs/techpubs/0/15650.html • Environmentallyharmful: acid wastes are produced Points 2, 3 and 4 and related problems can be eliminated in scH2O: no acid required & no byproducts. See Y. Ikushima et al. J. Am. Chem. Soc. 122, 1908 (2000); Work done on collaboration with: Michele Parrinello, Kiyoyuki Terakura, Tamio Ikeshoji and Chee Chin Liew

17. World wide production of e-caprolactam unit = 1000 ton/year

18. Hydrogen-bond network in water T=653 K r=0.73 g/cm3 Supercritical water Normal water T = 300 K r = 1.00 g/cm3 Continuous hydrogen-bond NW Disrupted hydrogen-bond NW

19. Beckmann rearrangement reaction proton attack to O in strong acid and supercritical H2O proton attack to N hydrolysis in H2O and superheatedH2O

20. Which are the important ingredients that make water special at supercritical conditions ? • Proton attack is the trigger (experimental outcome !) High efficiency : fast proton diffusion difference in hydration between O and N High selectivity :

21. ~ 0.9 ps in scH2O Contrary to normal liquid water, scH2O accelerates selectively the formation of the first intermediate ~ 5.2 ps in n-H2O Proton attack to N: Cycrohexanon (byproduct) formation

22. Acid catalyst ? Efficient reaction in scH2Odue tofast proton diffusionand acid properties of the (broken) Eigen-Zundel complexes

23. Cyclohexanone-oxyme in scH2O: • The energy barrier seems rather high and the reaction pathway not unique • The reaction is generally acid catalyzed, henceprotonsare expected to be essential in triggering the process • At supercritical conditions, however, theKwof water increase, henceH+andOH-can be around in the solvent in non-negligible concentration • And small amounts of weak acids greatly enhance reaction rates

24. Proton diffusion in ordinary liquid and supercritical water • Hydrogen bond network is disrupted in SCW. Is the proton diffusion slowed down in scH2O ? Not really…

25. Proton diffusion：normal water and supercritical water Proton (structural defect) diffusion coefficient estimation in the 3 systems: In scH2Othe network is disrupted and the motion occurs in sub-networks that join and break apart rapidly due to density fluctuations;two diffusion regimes are cooperating: hydrodynamics (vehicular) and Grotthus

26. Reaction selectivity ? Selective reactionin scH2O due todifferentsolvation of O and N

27. Cyclohexanone-oxyme in scH2O (+ H+): the selectivity T = 673 K H+ wet dry

28. R—C—R’ N—OH H+ R—C—R’ + H2O N+ Cyclohexanone-oxyme in scH2O (+ H+) R—N＝C+—R’ A very small activation barrier (about 1 kcal/mol) is required for the N insertion process.

29. …and now the second step: C-O bond formation DE = 5.9 kcal/mol DF = 5.1 kcal/mol Approach of an H2O molecule, x = |Owat-C|

30. oxime R—C—R’ R—C—R’ + H2O N—OH N+ H+ H2O H H O OH－ HO H HO R—N＝C—R’ R—N＝C+—R’ R—N＝C—R’ amide H O R—N—C—R’ The last step:eventually the e-caprolactam

31. The last step:eventually the e-caprolactam Proton exchange in scH2O(metadynamics)

32. Free energy surface: a less rugged landscape s1= 0.5 s1= 1.0

33. Conclusionsandperspectives • The H+ diffusion in scH2O occurs in sub-networks that join and break rapidly due to density fluctuations:two diffusion regimes are present. • Destabilization of Eigen (Zundel) complex makes scH2O an acid-like environment able to trigger chemical reaction • The selectivity of thecyclohexanone-oxymetoe-caprolactam reaction could be understood • The role of the H-bond indifferentiating the solvation featuresof the solute has been evidenced • A new green chemistry perspective has been explored. Related Publications: M.B. et al., Phys. Rev. Lett.85, 3245 (2000); J. Chem. Phys. 115, 2219 (2001); Phys. Rev. Lett. 90, 226403 (2003) ; J. Am. Chem. Soc.126, 6280 (2004); ChemPhysChem6, 1775 (2005)

34. Acknowledgements • Michele Parrinello, ETHZ-USI and Pisa University • Roberto Car, Princeton University • Kiyoyuki Terakura, JAIST, AIST and Hokkaido University • Michiel Sprik, Cambridge University • Pier Luigi Silvestrelli, Padova University • Alessandro Laio, SISSA, Trieste • Jürg Hutter, Zurich University • Marcella Iannuzzi, Zurich Univeristy • Carlo Massobrio, IPCMS