[ase-users] Parallelization of GPAW MD simulation

[Aleksander Dawid]

Hi ASE users

I am about to simulate simple 13-atoms argon cluster
I have run the GPAW software over MPI library

mpiexec -n 32 gpaw-python Ar13MD.py

The Ar13MD.py script is shown at the end of this letter

I wonder if I can speed up more my calculations.
I have a box of 12 angstroms as the simulation cell.
But i don't want periodic boundary conditions.

I am not an Expert, but the output below sugest me some solutions;

Total number of cores used: 32
Domain Decomposition: 4 x 4 x 2
Diagonalizer layout: Serial LAPACK
Orthonormalizer layout: Serial LAPACK

1 k-point in the Irreducible Part of the Brillouin Zone (total: 1)
Linear Mixing Parameter:           0.25
Pulay Mixing with 3 Old Densities
No Damping of Long Wave Oscillations

1. Serial LAPACK (could it be parallel version)
2. k-points (is this possible to set the k-points for cluster,  
especially if i want to simulate liquid state)
3. Damping of Long Wave Oscillations (how to set cut-off energy ?)

I will be very appriciated for any sugestions

With the best regards
Aleksander Dawid

======================================================
Aleksander Dawid
University of Silesia, Devision Of Computational Physics And Electronics
email: aleksander.dawid at us.edu.pl
======================================================

#!/usr/bin/env python
from ase import *
import gpaw.mpi as mpi
from gpaw import GPAW, setup_paths
from ase.lattice.surface import *
from ase.md.velocitydistribution import MaxwellBoltzmannDistribution
from ase.md.verlet import VelocityVerlet
from ase.md.langevin import Langevin
from ase import units
from ase.parallel import rank, barrier
from gpaw.vdw import FFTVDWFunctional
from gpaw.atom.generator import Generator, parameters


gO = Generator('Ar', 'revPBE', scalarrel=True, nofiles=True)
gO.run(**parameters['Ar'])

barrier()
setup_paths.insert(0, '.')

a = 12.  # Size of unit cell (Angstrom)
ayz= 12.
c = ayz / 2
s = a / 2


d = 3.42  # Bond length


fname="PAR500.txt"
vname="Velocities500.txt"

cluster = Atoms([Atom('Ar', (s - d, c, c)),
                           Atom('Ar', (s + d, c, c)),
                           Atom('Ar', (s + d, c + d, c)),
                           Atom('Ar', (s - d, c + d, c)),
                           Atom('Ar', (s, c + d, c)),
                           Atom('Ar', (s + d, c - d, c)),
                           Atom('Ar', (s - d, c - d, c)),
                           Atom('Ar', (s, c - d, c)),
                           Atom('Ar', (s + d, c, c + d)),
                           Atom('Ar', (s - d, c, c + d)),
                           Atom('Ar', (s, c, c + d)),
                           Atom('Ar', (s, c, c - d)),
                           Atom('Ar', (s, c, c))],
                            cell=(a, ayz, ayz), pbc=False)


vdw = FFTVDWFunctional(nspins=1,verbose=0)
calc=GPAW(width=0.1,eigensolver='cg',xc=vdw,txt=fname)


cluster.set_calculator(calc)
# Set the momenta corresponding to T=300K
#MaxwellBoltzmannDistribution(cluster, 20*units.kB)

# We want to run MD with constant energy using the VelocityVerlet algorithm.
dyn = Langevin(cluster, 5*units.fs, 2*units.kB, 0.02)
#dyn = VelocityVerlet(cluster, dt=5.0*units.fs)  # 5 fs time step.

def printenergy(a=cluster):    #store a reference to atoms in the definition.
     epot = a.get_potential_energy() / len(a)
     ekin = a.get_kinetic_energy() / len(a)
     if mpi.rank==0:
	#count+=1
	FILE = open("EnergyArCluster500.txt","a")
	FILE.write("Energy per atom: Epot = %.3feV  Ekin = %.3feV (T=%3.0fK)   
Etot = %.3feV" % (epot, ekin, ekin/(1.5*units.kB), epot+ekin))
     	#print ("Energy per atom: Epot = %.3feV  Ekin = %.3feV  
(T=%3.0fK)  Etot = %.3feV" % (epot, ekin, ekin/(1.5*units.kB),  
epot+ekin))
	FILE.close()



dyn.attach(printenergy, interval=1)
#We also want to save the positions of all atoms after every 100th time step.
traj = PickleTrajectory("Ar13dyn500.traj", 'w', cluster)
dyn.attach(traj.write, interval=1)


# Now run the dynamics
if mpi.rank==0:
	print ("START MOLECULAR QUANTUM DYNAMICS")
############################
#MD main loops
###########################
dyn.run(steps=1) # Make dyn.run(steps) MD simulations


if mpi.rank==0:
	v=cluster.get_velocities()
	FILE = open(vname,"w")
	FILE.write(str(len(v))+"\n")
	for i in range(0,len(v)):
        		s=cluster[i].symbol
		x=v[i][0]
    		y=v[i][1]
         	z=v[i][2]
		FILE.write(str(s)+" "+str(x)+" "+str(y)+" "+str(z)+"\n")
	FILE.close()

	write('Ar13po500.xyz', cluster)

----------------------------------------------------
Uniwersytet Śląski w Katowicach http://www.us.edu.pl