Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
fixed = []
for n in range(ntypes):
symbol = f.readline().strip()
symbols.extend([symbol] * natoms[n])
masses.extend([atommasses[n]] * natoms[n])
f.readline() # Coordinates of Component n
for i in range(natoms[n]):
row = f.readline().split()
coords.append([float(x) for x in row[:3]])
fixed.append(bool(int(row[3])))
atoms = Atoms(symbols=symbols,
positions=coords,
masses=masses,
cell=cellpar_to_cell(cellpar),
constraint=FixAtoms(mask=fixed),
info=dict(comment=comment))
images.append(atoms)
first_line = f.readline()
if first_line == '':
more_images_to_read = False
if isinstance(fileobj, basestring):
f.close()
if not index:
return images
else:
slab.set_pbc((1,1,0))
# Initial state.
# Add the N2 molecule oriented at 60 degrees:
d = 1.10 # N2 bond length
N2mol = Atoms('N2',positions=[[0.0,0.0,0.0],[0.5*3**0.5*d,0.5*d,0.0]])
add_adsorbate(slab,N2mol,height=1.0,position='fcc')
# Use the EMT calculator for the forces and energies:
slab.calc = EMT()
# We don't want to worry about the Cu degrees of freedom,
# so fix these atoms:
mask = [atom.symbol == 'Cu' for atom in slab]
slab.set_constraint(FixAtoms(mask=mask))
# Relax the structure
relax = QuasiNewton(slab)
relax.run(fmax=0.05)
print('initial state:', slab.get_potential_energy())
write('N2.traj', slab)
# Now the final state.
# Move the second N atom to a neighboring hollow site:
slab[-1].position[0] = slab[-2].position[0] + 0.25 * slab.cell[0,0]
slab[-1].position[1] = slab[-2].position[1]
# and relax.
relax.run()
print('final state: ', slab.get_potential_energy())
write('2N.traj', slab)
from ase.lattice.surface import fcc111, add_adsorbate
from ase.data.molecules import molecule
from ase.constraints import FixAtoms
atoms = fcc111('Au', size=(3,3,3), vacuum=10)
benzene = molecule('C6H6')
benzene.translate(-benzene.get_center_of_mass())
# I want the benzene centered on the position in the middle of atoms
# 20, 22, 23 and 25
p = (atoms.positions[20] +
atoms.positions[22] +
atoms.positions[23] +
atoms.positions[25])/4.0 + [0.0, 0.0, 3.05]
benzene.translate(p)
atoms += benzene
# now we constrain the slab
c = FixAtoms(mask=[atom.symbol=='Au' for atom in atoms])
atoms.set_constraint(c)
#from ase.visualize import view; view(atoms)
with jasp('surfaces/Au-benzene-pbe',
xc='PBE',
encut=350,
kpts=(4,4,1),
ibrion=1,
nsw=100,
atoms=atoms) as calc:
print atoms.get_potential_energy()
# the clean gold slab
from vasp import Vasp
from ase.lattice.surface import fcc111, add_adsorbate
from ase.constraints import FixAtoms
atoms = fcc111('Au', size=(3,3,3), vacuum=10)
# now we constrain the slab
c = FixAtoms(mask=[atom.symbol=='Au' for atom in atoms])
atoms.set_constraint(c)
#from ase.visualize import view; view(atoms)
print(Vasp('surfaces/Au-pbe',
xc='PBE',
encut=350,
kpts=[4, 4, 1],
ibrion=1,
nsw=100,
atoms=atoms).potential_energy)
from vasp import Vasp
from ase.lattice.surface import fcc110
from ase.io import write
from ase.constraints import FixAtoms
atoms = fcc110('Ag', size=(2, 1, 6), vacuum=10.0)
del atoms[11] # delete surface row
constraint = FixAtoms(mask=[atom.tag > 2 for atom in atoms])
atoms.set_constraint(constraint)
Vasp('surfaces/Ag-110-missing-row',
xc='PBE',
kpts=[6, 6, 1],
encut=350,
ibrion=2,
isif=2,
nsw=10,
atoms=atoms).update()
fix_params.append(
FixScaledParametricRelations.from_expressions(
list(range(len(atoms))), atomic_params, atomic_expressions
)
)
if any(magmoms):
atoms.set_initial_magnetic_moments(magmoms)
if any(charges):
atoms.set_initial_charges(charges)
if periodic.any():
atoms.set_cell(cell)
atoms.set_pbc(periodic)
if len(fix):
atoms.set_constraint([FixAtoms(indices=fix)] + fix_cart + fix_params)
else:
atoms.set_constraint(fix_cart + fix_params)
if fix_params and apply_constraints:
atoms.set_positions(atoms.get_positions())
return atoms
from ase.build import fcc110
from ase.optimize.minimahopping import MinimaHopping
from ase.calculators.emt import EMT
from ase.constraints import FixAtoms, Hookean
# Make the Pt 110 slab.
atoms = fcc110('Pt', (2, 2, 2), vacuum=7.)
# Add the Cu2 adsorbate.
adsorbate = Atoms([Atom('Cu', atoms[7].position + (0., 0., 2.5)),
Atom('Cu', atoms[7].position + (0., 0., 5.0))])
atoms.extend(adsorbate)
# Constrain the surface to be fixed and a Hookean constraint between
# the adsorbate atoms.
constraints = [FixAtoms(indices=[atom.index for atom in atoms if
atom.symbol=='Pt']),
Hookean(a1=8, a2=9, rt=2.6, k=15.),
Hookean(a1=8, a2=(0., 0., 1., -15.), k=15.),]
atoms.set_constraint(constraints)
# Set the calculator.
calc = EMT()
atoms.calc = calc
# Instantiate and run the minima hopping algorithm.
hop = MinimaHopping(atoms,
Ediff0=2.5,
T0=4000.)
hop(totalsteps=10)
view(slab)
#some positions needed to place the atom in the correct place
x1 = 1.379
x2 = 4.137
x3 = 2.759
y1 = 0.0
y2 = 2.238
z1 = 7.165
z2 = 6.439
#Add the adatom to the list of atoms and set constraints of surface atoms.
slab += Atoms('N', [ ((x2+x1)/2,y1,z1+1.5)])
mask = [atom.symbol == 'Pt' for atom in slab]
slab.set_constraint(FixAtoms(mask=mask))
#optimise the initial state
# Atom below step
initial = slab.copy()
initial.set_calculator(EMT())
relax = QuasiNewton(initial)
relax.run(fmax=0.05)
view(initial)
#optimise the initial state
# Atom above step
slab[-1].position = (x3,y2+1,z2+3.5)
final = slab.copy()
final.set_calculator(EMT())
relax = QuasiNewton(final)
relax.run(fmax=0.05)
from vasp import Vasp
from ase.lattice.surface import fcc110
from ase.io import write
from ase.constraints import FixAtoms
atoms = fcc110('Ag', size=(2, 1, 6), vacuum=10.0)
constraint = FixAtoms(mask=[atom.tag > 2 for atom in atoms])
atoms.set_constraint(constraint)
calc = Vasp('surfaces/Ag-110',
xc='PBE',
kpts=[6, 6, 1],
encut=350,
ibrion=2,
isif=2,
nsw=10,
atoms=atoms)
calc.update()
final = initial.copy()
final.positions[-1, 1] += b
view([initial, final])
# Construct a list of images:
images = [initial]
for i in range(5):
images.append(initial.copy())
images.append(final)
# Make a mask of zeros and ones that select fixed atoms (the
# two bottom layers):
mask = initial.positions[:, 2] - min(initial.positions[:, 2]) < 1.5 * h
constraint = FixAtoms(mask=mask)
print(mask)
for image in images:
# Let all images use an EMT calculator:
image.calc = EMT()
image.set_constraint(constraint)
# Relax the initial and final states:
QuasiNewton(initial).run(fmax=0.05)
QuasiNewton(final).run(fmax=0.05)
# Create a Nudged Elastic Band:
neb = NEB(images)
# Make a starting guess for the minimum energy path (a straight line
# from the initial to the final state):