PWscf, FPMD and CP (Democritos Package) Tutorial

1. PWscf, FPMD and CP (Democritos Package)Tutorial Carlo Cavazzoni(High Performance Computing Group) CINECA

2. Outline • ab-initio simulations introduction • PWSCF and Democritos package • Installation • Functionality • Simple examples

3. The scope of computersimulation • Measure theories.Solve equations which could not be solved otherwise. “Get numbers out of theories” much in the same way as experiments “get numbers out of a natural process” • Do virtual experimentswhere experimental conditions can be controlled down to the atomic scale in ways which would not be possible in the lab

4. Ab initio simulations The Born-Oppenheimer approximation (M>>m)

5. Density functional theory Kohn-Sham Hamiltonian

6. Kohn-Sham equations from functional minimization Helmann & Feynman Kohn & Sham

7. The tricks of the trade • Expanding the Kohn-Sham orbitals into a suitable basis set turns DFT into a multi-variate minimization problem, and the KS equation into a non-linear matrix eigenvalue problem • The use of pseudo-potentials allows to ignore chemically inert core states and to use plane waves (the name of the game!)

8. The tricks of the trade (II) • Plane waves are orthogonal and the matrix elements are usually easy to calculate; the effective completeness of the basis is easy to check • Plane-waves allow to calculate efficiently matrix-vector products and to solve the Poisson equation using FFT’s • Supercells for treating finite (or semi-infinite) systems • Iterative diagonalization vs. global minimization

9. The tricks of the trade (III) • Summing over occupied states: special-point and Gaussian-smearing techniques • Non-linear extrapolation for SCF acceleration and density prediction in MD • Choice of fictitious masses in CP dynamics

10. Which algorithm shall I use? • Electronic structure:SCF diagonalization vs. energy minimization • Geometry optimization:standard DFT • Lattice vibrations, static response functions:DF perturbation theory • Dynamics:Car-Parrinello vs. Born-Oppenheimer • Slow kinetics and rare events:path sampling vs. Parrinello-Laio metadynamics • Optical properties, excited states:Time-dependent DFT & many-body perturbation theory

11. PWSCF and Democritos (www.democritos.it) (DEmocritos MOdeling Center for Research in aTOmistic Simulatios) Democritos is a National Simulation Center of the Italian Istituto Nazionale per la Fisica della Materia (INFM). Hosted by SISSA in Trieste PWSCF is promoted, maintained and developed by Democritos CINECA is a Partner of Democritos.

12. Democritos and Scientific Software GOAL: Development of software for next generation parallel computers (N>10000 processors) Production of high-quality software package for atomistic simulations based on density-functional theory (DFT), plane wave (PW) and pseudopotentials (PP)

13. PWSCF is also part of a merge effortof pre-existing scientific software PWscfpackage: electronic structure, structural optimization, molecular dynamics, vibrational and dielectric properties. Developed by S. Baroni, S. de Gironcoli, A. Dal Corso (SISSA, Trieste), P. Giannozzi (Scuola Normale, Pisa) and others.See http://www.pwscf.org for more information and downloads. CPcode: Car-Parrinello variable-cell molecular dynamics. Developed by A. Pasquarello (IRRMA, Lausanne), K. Laasonen (Oulu), A. Trave (UCBerkeley), R. Car (Princeton), P. Giannozzi and others. Based on the original code written by R. Car and A. Pasquarello. Download CP (initial public release) FPMD code: Car-Parrinello variable-cell molecular dynamics. Developed by C. Cavazzoni (CINECA, Bologna), S. Scandolo (ICTP, Trieste), G. Chiarotti (SISSA, Trieste), P. Focher, G. Ballabio and others. Based on the original code written by R. Car and M. Parrinello. Download FPMD (initial public release) GNU License http://www.democritos.it/scientific.php

14. Present Package Status • The Package is developed in Fortran90 • The codes in the Package share: • installation • most parts of the basic code • input format • PP format • a graphical GUI is being developed • output format (CP/FPMD fully compatible/restartable)

15. Basic Data Type Real space Charge density 3D arrays Wave functions 1D arrays Reciprocal space

16. Reciprocal Space Representation Wave Functions To truncate the infinite sum Charge Density To retain the same accurancy as the wave function

17. DFT Functional

18. FFTs Ecut 4Ecut Reciprocal Space |G|2/2< Ecut |G|2/2< 4Ecut FFT Real Space

19. P0 P2 P3 P4 ParallelizationReciprocal Space distribution PWSCF parallelize also over K-points

20. ParallelizationFFT algorithm

21. SP4, optimization issues • use essl • use non-blocking communications • use mass • modify strides to reduce TLB misses • use block-algorithm whenever possible ATLAS strategy

22. Case study: matrix transposition do ib = 1, nb ioff = (ib-1) * bsiz do jb = 1, mb joff = (jb-1) * bsiz do j = 1, bsiz do i = 1, bsiz buf(i,j) = x(i+ioff, j+joff) enddo enddo do j = 1, bsiz do i = 1, j-1 bswp = buf(i,j) buf(i,j) = buf(j,i) buf(j,i) = bswp enddo enddo do i=1,bsiz do j=1,bsiz y(j+joff, i+ioff) = buf(j,i) enddo enddo enddo enddo do i=1,n do j=1,m y(j,i) = x(i,j) enddo enddo block-algorithm bsiz = block size nb = n / bsiz mb = n / bsiz more code for reminder: MOD(n / bsiz) /= 0 OR MOD(m / bsiz) /= 0

23. Block algorithm for matrix transposition 80 70 60 50 Mips 40 30 20 10 0 0 20 40 60 80 100 120 block dimension

24. CP code Developed by A. Pasquarello (IRRMA, Lausanne), K. Laasonen (Oulu), A. Trave (UCBerkeley), R. Car (Princeton), P. Giannozzi and others. Based on the original code written by R. Car and A. Pasquarello. Verlet dynamics with mass preconditioning Temperature control: Nose’ thermostat, velocity rescaling Metallic systems: Nose’ thermostat for both electrons and ions electronic and ionic minimization via damped dynamics Modified kinetic functional for costant-pressure calculations “grid box” for fast treatment of augmentation terms in Ultrasoft PP’s Nudget Elastic Band (NEB) scheme for transition paths and energy barriers Limitations: no k-points no fancier minimization chemes no constraints

25. FPMD code Developed by C. Cavazzoni (CINECA, Bologna), S. Scandolo (ICTP, Trieste), G. Chiarotti (SISSA, Trieste), P. Focher, G. Ballabio and others. Based on the original code written by R. Car and M. Parrinello. Verlet dynamics with mass preconditioning Temperature control: Nose’ thermostat, velocity rescaling Metallic systems: Nose’ thermostat for both electrons and ions Various electronic and ionic minimization schemes Modified kinetic functional for costant-pressure calculations Macroscopic polarization via Berry Phase Nudget Elastic Band (NEB) scheme for transition paths and energy barriers Constrained dynamics Limitations: no Ultrasoft PP’s

26. PWSCF code Developed by S. Baroni, S. de Gironcoli, A. Dal Corso (SISSA, Trieste), P. Giannozzi (Scuola Normale, Pisa) and others. Self-consistent ground-state energy and Kohn-Sham orbitals Structural optimization Molecular dynamics on the ground-state Born-Oppenheimer surface Variable-cell molecular dynamics Phonon frequencies and eigenvectors at a generic wave vector Effective charges and dielectric tensors Electron-phonon interaction coefficients for metals Third-order anharmonic phonon lifetimes Macroscopic polarization via Berry Phase Nudget Elastic Band (NEB) scheme for transition paths and energy barriers Limitations: no Car-Parrinello dynamics very limited constrained minimization and dynamics

27. Access the anonymous CVS of theDemocritos PW-PP-DFT Package anonymous CVS contains the whole package, with all the codes (tcsh/csh): setenv CVS_RSH ssh setenv CVSROOT :pserver:cvsanon@democritos.sissa.it:/home/cvs (sh/bash): export CVS_RSH=ssh export CVSROOT=:pserver:cvsanon@democritos.sissa.it:/home/cvs Then: cvs login (password: cvsanon). For the first code download: cvs co O-sesame Open-Sesame: Open-Source Scalable Electronic Structure and Atomistic Modeling Environment

28. Installation ./configure make target where target is one of the following: pw basic code for scf, structure optimization, MD fpmd FPMD code for Car-Parrinello MD cp CP code: CP MD with ultrasoft pseudopotentials ph phonon code pp postprocessing programs gamma Gamma-only version of phonon code nc non collinear magnetic version of pw code pwcond ballistic conductance d3 third-order derivatives tools misc tools for data analysis upf utilities for pseudopotential conversion pwall same as "make pw ph pp gamma nc pwcond d3 tools" all same as "make pwall fpmd cp upf" links creates links to executables in bin/ clean remove executables and objects veryclean revert distribution to the original status tar create a tarball of the source tree tar-gui create a tarball of the GUI sources

29. Installation (DEMO) how to configure configure?

30. Directories tree I/II include includes files Modules clib flib upftools common modules, subroutine and tools CPV FPMD PW main codes directories D3 PH PP pw post-processing codes

31. Directories tree II/II pwdocs cpdocs documentation pseudo pseudopotential library install preconfigured installation parameter files cp_examples pw_examples run examples

32. Input files PWscf, CP and FPMD share the same input files format All codes read input parameters from standard input Codes read from files PP data, all codes read the same PP files In directory upftools there are conversion tools for PP formats (there is one also for CPMD)

33. Standard Input Layout The standard input is divided in two main section the first contains fortran namelists the seconds contains optional input “CARDS” all input parameter are described in cpdocs and pwdocs file. all input parameter for all codes are contained in file Modules/input_parameters.f90

34. layout &CONTROL control_parameter_1, control_parameter_2, ....... / &SYSTEM sistem_parameter_1, sistem_parameter_2, ....... / &ELECTRONS electrons_parameter_1, electrons_parameter_2, ....... / &IONS ions_parameter_1, ions_parameter_2, ....... / &CELL cell_parameter_1, cell_parameter_2, ....... / ‘ATOMIC_SPECIES’ label_1 mass_1 pseudo_file_1 label_2 mass_2 pseudo_file_2 ..... ‘ATOMIC_POSITIONS’ label_1 px_1 py_1 pz_1 label_2 px_2 py_2 pz_2 ..... the namelist order is fixed, but within a namelist the parameter sequence is irrelevant namelists card order is irrelevant, but within a card the parameter layout is fixed cards

35. Input Namelists &CONTROL input variables that control the flux of teh calculation and the amount of I/O on disk and on the screen &SYSTEM input variables that specify the system under study &ELECTRONS input variables that control the algorithms used to reach the self-consistent solution of KS equations for the electrons, and the electrons dynamics (for CP and FPMD) &IONS input variables that control ionic motion in molecular dynamics run or structural relaxation &CELL input variables that control cell-shape evolution in a variable-cell MD run or structural relaxation

36. Input Cards ATOMIC_SPECIES name, mass and pseudopotential used for each atomic species present in the system ATOMIC_POSITIONS type and coordinate of each atom in the unit cell K_POINTS coordinate and weights of the k-points used for BZ integration CELL_PARAMETERS … OCCUPATIONS … CLINBING_IMAGES …

37. Silicon Band Structure(www.cineca.it/~acv0/Public/PWSCF_Tutorial.tgz) • Select the appropriate unit cell • Input the atomic coordinate • Choose the pseudopotentials • Determine k-point sampling • Select the size of the basis set (Ecut)

38. Select the appropriate unit cell ibrav bravais lattice celldm(.) cell parameters celldm(1) = a, celldm(2) = b/a, celldm(3) = c/a celldm(4) = cos(bc) celldm(5) = cos(ac) celldm(6) = cos(ab) • ibrav structure celldm(2)-celldm(6) • 0 "free", see above not used • 1 cubic P (sc) not used • 2 cubic F (fcc) not used • 3 cubic I (bcc) not used • 4 Hexagonal and Trigonal P celldm(3)=c/a • 5 Trigonal R celldm(4)=cos(aalpha) • 6 Tetragonal P (st) celldm(3)=c/a • 7 Tetragonal I (bct) celldm(3)=c/a • Orthorhombic P celldm(2)=b/a, • celldm(3)=c/a • Orthorhombic base-centered(bco) celldm(2)=b/a, • celldm(3)=c/a • Orthorhombic face-centered celldm(2)=b/a, • celldm(3)=c/a • Orthorhombic body-centered celldm(2)=b/a, • celldm(3)=c/a • Monoclinic P celldm(2)=b/a, • celldm(3)=c/a, • celldm(4)=cos(ab) • Monoclinic base-centered celldm(2)=b/a • celldm(3)=c/a, • celldm(4)=cos(ab) • 14 Triclinic P celldm(2)= b/a, • celldm(3)= c/a, • celldm(4)= cos(bc), • celldm(5)= cos(ac), • celldm(6)= cos(ab)

39. Unit cell ibrav = 2, cubic fcc ibrav = 1, simple cubic

40. Unit cell (supercells) ibrav = 1, simple cubic (crystal with defects) ibrav = 6, tetragonal (surfaces)

41. Silicon Crystal unit cell &control prefix='silicon', pseudo_dir = '$PSEUDO_DIR/', outdir='$TMP_DIR/' / &system ibrav= 2, celldm(1) =10.20, nat= 2, ntyp= 1, ecutwfc =12.0, / &electrons / empty namelist, keep default values for Self-Consistent Field iterations

42. Atomic Positions ATOMIC_POSITION (units) units: if units is not specified, unit of “celldm(1)” is assumed units = bohr position in Bohr radius units = angstrom position in Angstrom units = crystal position in crystal coord. ATOMIC_SPECIES Si 28.086 Si.vbc.UPF ATOMIC_POSITIONS Si 0.00 0.00 0.00 Si 0.25 0.25 0.25

43. K-points set Automatic K_POINTS (automatic) 2 2 2 1 1 1 k-points grid grid off-set K_POINTS 2 0.25 0.25 0.75 3.0 0.25 0.25 0.25 1.0 User specified number of k-points k-points weight k-points coordinate (crystal coord.) start with coarse grid, then try finer ones to check convergence

44. Input for Silicon Crystal SCF calculation &control prefix='silicon', pseudo_dir = '$PSEUDO_DIR/', outdir='$TMP_DIR/' / &system ibrav= 2, celldm(1) =10.20, nat= 2, ntyp= 1, ecutwfc =12.0, / &electrons / ATOMIC_SPECIES Si 28.086 Si.vbc.UPF ATOMIC_POSITIONS Si 0.00 0.00 0.00 Si 0.25 0.25 0.25 K_POINTS 2 0.25 0.25 0.75 3.0 0.25 0.25 0.25 1.0 pw.x < inputfile

45. Non SCF calculation To compute band structure, we need the value of eigenvectors on all k-point of the bands. If the computation is converged with respect to the k-points grid, we do not need to repeat the SCF calculation on all k-points of the bands, but we can diagonalize the Hamiltonian using the SCF Kohn-Sham potential.

46. Input for Non SCF calculation &control calculation='nscf' pseudo_dir = '$PSEUDO_DIR/', outdir='$TMP_DIR/', prefix='silicon' / &system ibrav= 2, celldm(1) =10.20, nat= 2, ntyp= 1, ecutwfc =12.0, nbnd = 8, / &electrons / ATOMIC_SPECIES Si 28.086 Si.vbc.UPF ATOMIC_POSITIONS Si 0.00 0.00 0.00 Si 0.25 0.25 0.25 K_POINTS …

47. K-points set along high symmetry lines # band structure calculation along delta, sigma and lambda lines K_POINTS 36 0.5 0.5 0.5 1 0.4 0.4 0.4 2 0.3 0.3 0.3 3 0.2 0.2 0.2 4 0.1 0.1 0.1 5 0.0 0.0 0.0 6 0.0 0.0 0.1 7 0.0 0.0 0.2 8 0.0 0.0 0.3 9 0.0 0.0 0.4 10 0.0 0.0 0.5 11 0.0 0.0 0.6 12 0.0 0.0 0.7 13 0.0 0.0 0.8 14 0.0 0.0 0.9 15 0.0 0.0 1.0 16 0.0 0.1 1.0 17 0.0 0.2 1.0 18 0.0 0.3 1.0 19 0.0 0.4 1.0 20 0.0 0.5 1.0 21 0.0 0.6 1.0 22 0.0 0.7 1.0 23 0.0 0.8 1.0 24 0.0 0.9 1.0 25 0.0 1.0 1.0 26 0.0 0.9 0.9 27 0.0 0.8 0.8 28 0.0 0.7 0.7 29 0.0 0.6 0.6 30 0.0 0.5 0.5 31 0.0 0.4 0.4 32 0.0 0.3 0.3 33 0.0 0.2 0.2 34 0.0 0.1 0.1 35 0.0 0.0 0.0 36

48. Post Processing &inputpp prefix = 'silicon' outdir = '$TMP_DIR/' filband = 'sibands.dat' / bands.x < inputfile prints the bands to graphical files sibands.dat -6.0 10 sibands.xmgr sibands.ps 6.3369 1.0 6.3369 plotband.x < inputfile

49. the silicon bands!

50. How to find equilibrium lattice parameters different scf computation, changing cut-off and lattice parameter