1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75
|
import numpy as np
from ase.calculators.eam import EAM
from ase.build import bulk
# test to generate an EAM potential file using a simplified
# approximation to the Mishin potential Al99.eam.alloy data
from scipy.interpolate import InterpolatedUnivariateSpline as spline
cutoff = 6.28721
n = 21
rs = np.arange(0, n) * (cutoff / n)
rhos = np.arange(0, 2, 2. / n)
# generated from
# mishin = EAM(potential='../potentials/Al99.eam.alloy')
# m_density = mishin.electron_density[0](rs)
# m_embedded = mishin.embedded_energy[0](rhos)
# m_phi = mishin.phi[0,0](rs)
m_density = np.array([2.78589606e-01, 2.02694937e-01, 1.45334053e-01,
1.06069912e-01, 8.42517168e-02, 7.65140344e-02,
7.76263116e-02, 8.23214224e-02, 8.53322309e-02,
8.13915861e-02, 6.59095390e-02, 4.28915711e-02,
2.27910928e-02, 1.13713167e-02, 6.05020311e-03,
3.65836583e-03, 2.60587564e-03, 2.06750708e-03,
1.48749693e-03, 7.40019174e-04, 6.21225205e-05])
m_embedded = np.array([1.04222211e-10, -1.04142633e+00, -1.60359806e+00,
-1.89287637e+00, -2.09490167e+00, -2.26456628e+00,
-2.40590322e+00, -2.52245359e+00, -2.61385603e+00,
-2.67744693e+00, -2.71053295e+00, -2.71110418e+00,
-2.69287013e+00, -2.68464527e+00, -2.69204083e+00,
-2.68976209e+00, -2.66001244e+00, -2.60122024e+00,
-2.51338548e+00, -2.39650817e+00, -2.25058831e+00])
m_phi = np.array([6.27032242e+01, 3.49638589e+01, 1.79007014e+01,
8.69001383e+00, 4.51545250e+00, 2.83260884e+00,
1.93216616e+00, 1.06795515e+00, 3.37740836e-01,
1.61087890e-02, -6.20816372e-02, -6.51314297e-02,
-5.35210341e-02, -5.20950200e-02, -5.51709524e-02,
-4.89093894e-02, -3.28051688e-02, -1.13738785e-02,
2.33833655e-03, 4.19132033e-03, 1.68600692e-04])
m_densityf = spline(rs, m_density)
m_embeddedf = spline(rhos, m_embedded)
m_phif = spline(rs, m_phi)
a = 4.05 # Angstrom lattice spacing
al = bulk('Al', 'fcc', a=a)
mishin_approx = EAM(elements=['Al'], embedded_energy=np.array([m_embeddedf]),
electron_density=np.array([m_densityf]),
phi=np.array([[m_phif]]), cutoff=cutoff, form='alloy',
# the following terms are only required to write out a file
Z=[13], nr=n, nrho=n, dr=cutoff / n, drho=2. / n,
lattice=['fcc'], mass=[26.982], a=[a])
al.set_calculator(mishin_approx)
mishin_approx_energy = al.get_potential_energy()
mishin_approx.write_potential('Al99-test.eam.alloy')
mishin_check = EAM(potential='Al99-test.eam.alloy')
al.set_calculator(mishin_check)
mishin_check_energy = al.get_potential_energy()
print('Cohesive Energy for Al = ', mishin_approx_energy, ' eV')
error = (mishin_approx_energy - mishin_check_energy) / mishin_approx_energy
print('read/write check error = ', error)
assert abs(error) < 1e-4
|