Python FaceCenteredCubic.get_cell方法代码示例

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import get_cell [as 别名]
def MakeAtoms(elem1, elem2=None):
    if elem2 is None:
        elem2 = elem1
    a1 = reference_states[elem1]['a']
    a2 = reference_states[elem2]['a']
    a0 = (0.5 * a1**3 + 0.5 * a2**3)**(1.0/3.0) * 1.03
    if ismaster:
        print "Z1 = %i,  Z2 = %i,  a0 = %.5f" % (elem1, elem2, a0)
    # 50*50*50 would be big enough, but some vacancies are nice.
    atoms = FaceCenteredCubic(symbol='Cu', size=(51,51,51))
    nremove = len(atoms) - 500000
    assert nremove > 0
    remove = np.random.choice(len(atoms), nremove, replace=False)
    del atoms[remove]
    if elem1 != elem2:
        z = atoms.get_atomic_numbers()
        z[np.random.choice(len(atoms), len(atoms)/2, replace=False)] = elem2
        atoms.set_atomic_numbers(z)
    if isparallel:
        # Move this contribution into position
        uc = atoms.get_cell()
        x = mpi.world.rank % cpuLayout[0]
        y = (mpi.world.rank // cpuLayout[0]) % cpuLayout[1]
        z = mpi.world.rank // (cpuLayout[0] * cpuLayout[1])
        assert(0 <= x < cpuLayout[0])
        assert(0 <= y < cpuLayout[1])
        assert(0 <= z < cpuLayout[2])
        offset = x * uc[0] + y * uc[1] + z * uc[2]
        new_uc = cpuLayout[0] * uc[0] + cpuLayout[1] * uc[1] + cpuLayout[2] * uc[2]
        atoms.set_cell(new_uc, scale_atoms=False)
        atoms.set_positions(atoms.get_positions() + offset)
        # Distribute atoms. Maybe they are all on the wrong cpu, but that will
        # be taken care of.
        atoms = MakeParallelAtoms(atoms, cpuLayout)
    MaxwellBoltzmannDistribution(atoms, T * units.kB)
    return atoms

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import get_cell [as 别名]
from __future__ import print_function, division
from ase.lattice.cubic import FaceCenteredCubic
from ase.lattice.hexagonal import HexagonalClosedPacked
from ase.test import must_raise

with must_raise(ValueError):
    # The Miller indices of the surfaces are linearly dependent
    atoms = FaceCenteredCubic(symbol='Cu',
                              miller=[[1, 1, 0], [1, 1, 0], [0, 0, 1]])

# This one should be OK:
atoms = FaceCenteredCubic(symbol='Cu',
                          miller=[[1, 1, 0], [0, 1, 0], [0, 0, 1]])
print(atoms.get_cell())


with must_raise(ValueError):
    # The directions spanning the unit cell are linearly dependent
    atoms = FaceCenteredCubic(symbol='Cu',
                              directions=[[1, 1, 0], [1, 1, 0], [0, 0, 1]])

with must_raise(ValueError):
    # The directions spanning the unit cell are linearly dependent
    atoms = FaceCenteredCubic(symbol='Cu',
                              directions=[[1, 1, 0], [1, 0, 0], [0, 1, 0]])

# This one should be OK:
atoms = FaceCenteredCubic(symbol='Cu',
                          directions=[[1, 1, 0], [0, 1, 0], [0, 0, 1]])
print(atoms.get_cell())

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import get_cell [as 别名]
        lgv.run(5)
        T = atoms.get_kinetic_energy() / (1.5 * atoms.get_number_of_atoms() * units.kB)
        print "Temperature is now %.2f K" % (T,)
    print "Desired temperature reached!"

    lgv.set_temperature(T_goal*units.kB)

    for i in range(2):
        lgv.run(20)
        s = atoms.get_stress()
        p = -(s[0] + s[1] + s[2])/3.0 / units.GPa
        T = atoms.get_kinetic_energy() / (1.5 * atoms.get_number_of_atoms() * units.kB)
        print "Pressure is %f GPa, desired pressure is %f GPa (T = %.2f K)" % (p, p_goal, T)
        dv = (p - p_goal) / bulk
        print "Adjusting volume by", dv
        cell = atoms.get_cell()
        atoms.set_cell(cell * (1.0 + dv/3.0))

    T = atoms.get_kinetic_energy() / (1.5 * atoms.get_number_of_atoms() * units.kB)
    print "Temperature is now %.2f K" % (T,)

    stressstate = array([-2, -1, 0, 0, 0, 0])*p_goal*units.GPa
    dyn = NPT(atoms, 5 * units.fs, T_goal*units.kB, stressstate,
              ttime*units.fs, (ptime*units.fs)**2 * bulk * units.GPa)
    traj = PickleTrajectory("NPT-atoms.traj", "w", atoms)
    #dyntraj = ParallelHooverNPTTrajectory("NPT-dyn-traj.nc", dyn, interval = 50)
    dyn.attach(traj, interval=50)
    #dyn.Attach(dyntraj)

    out = open(out1, "w")

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import get_cell [as 别名]
import pylab
from box.interpolation import interpolate_path
pi = np.pi

d = 4.08
atoms = FaceCenteredCubic('Au', latticeconstant=d,
                          directions=((0,1,1),(1,0,1),(1,1,0)),
                          align=False)
                          
calc = Hotbit(SCC=False, kpts=(8,8,8), txt='bs.cal')
atoms.set_calculator(calc)
atoms.get_potential_energy()

# reciprocal lattice vectors
V = atoms.get_volume()
a, b, c = atoms.get_cell()
a_ = 2*pi/V*np.cross(b, c)
b_ = 2*pi/V*np.cross(c, a)
c_ = 2*pi/V*np.cross(a, b)

gamma = np.array( (0,0,0) )

# a path which connects the L-points and gamma-point
kpts, distances, label_points = interpolate_path((-a_/2, a_/2, gamma, -b_/2, b_/2, gamma, -c_/2, c_/2, gamma, (-a_-b_-c_)/2, (a_+b_+c_)/2), 500)
labels = ["$-a$",'$a$','$\Gamma$','$-b$','$b$','$\Gamma$','$-c$','$c$','$\Gamma$','$-a-b-c$','$a+b+c$']


eigs = calc.get_band_energies(kpts,shift=True,rs='k')
e_min = np.min(eigs)
e_max = np.max(eigs)

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import get_cell [as 别名]
    def __init__(self, symbol=None, layers=None, positions=None,
                 latticeconstant=None, symmetry=None, cell=None, 
                 center=None, multiplicity=1, filename=None, debug=0):

        self.debug = debug
        self.multiplicity = multiplicity

        if filename is not None:
            # We skip MonteCarloAtoms.__init__, do it manually.
            self.mc_optim = np.zeros(101, np.intc)
            self.mc_optim[0] = 10000  # MC optim invalid
            self.read(filename)
            return
        
        #Find the atomic number
        if symbol is not None:
            if isinstance(symbol, str):
                self.atomic_number = atomic_numbers[symbol]
            else:
                self.atomic_number = symbol
        else:
            raise Warning('You must specify a atomic symbol or number!')

        #Find the crystal structure
        if symmetry is not None:
            if symmetry.lower() in ['bcc', 'fcc', 'hcp']:
                self.symmetry = symmetry.lower()
            else:
                raise Warning('The %s symmetry does not exist!' % symmetry.lower())
        else:
            self.symmetry = reference_states[self.atomic_number]['symmetry'].lower()

        if self.debug:
            print 'Crystal structure:', self.symmetry

        #Find the lattice constant
        if latticeconstant is None:
            if self.symmetry == 'fcc':
                self.lattice_constant = reference_states[self.atomic_number]['a']
            else:
                raise Warning(('Cannot find the lattice constant ' +
                               'for a %s structure!' % self.symmetry))
        else:
            self.lattice_constant = latticeconstant

        if self.debug:
            print 'Lattice constant(s):', self.lattice_constant

        #Make the cluster of atoms
        if layers is not None and positions is None:
            layers = list(layers)

            #Make base crystal based on the found symmetry
            if self.symmetry == 'fcc':
                if len(layers) != data.lattice[self.symmetry]['surface_count']:
                    raise Warning('Something is wrong with the defined number of layers!')

                xc = int(np.ceil(layers[1] / 2.0)) + 1
                yc = int(np.ceil(layers[3] / 2.0)) + 1
                zc = int(np.ceil(layers[5] / 2.0)) + 1

                xs = xc + int(np.ceil(layers[0] / 2.0)) + 1
                ys = yc + int(np.ceil(layers[2] / 2.0)) + 1
                zs = zc + int(np.ceil(layers[4] / 2.0)) + 1

                center = np.array((xc, yc, zc)) * self.lattice_constant
                size = (xs, ys, zs)

                if self.debug:
                    print 'Base crystal size:', size
                    print 'Center cell position:', center

                atoms = FaceCenteredCubic(symbol=symbol,
                                          size=size,
                                          latticeconstant=self.lattice_constant,
                                          align=False)

            else:
                raise Warning(('The %s crystal structure is not' +
                               ' supported yet.') % self.symmetry)

            positions = atoms.get_positions()
            numbers = atoms.get_atomic_numbers()
            cell = atoms.get_cell()
        elif positions is not None:
            numbers = [self.atomic_number] * len(positions)
        else:
            numbers = None

        #Load the constructed atoms object into this object
        self.set_center(center)
        MonteCarloAtoms.__init__(self, numbers=numbers, positions=positions,
                                 cell=cell, pbc=False)

        #Construct the particle with the assigned surfasces
        if layers is not None:
            self.set_layers(layers)

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import get_cell [as 别名]
e3 = atoms.get_potential_energy()
ReportTest("e3 correct", e3, ecorrect, 0.001)

atoms.set_pbc((1,1,0))
atoms.set_positions(atoms.get_positions())
e4 = atoms.get_potential_energy()
ReportTest("e4 correct", e4, ecorrect, 0.001)

print "Repeating tests with an atom outside the unit cell."
for coordinate in (0,1,2):
    print "Using coordinate number", coordinate

    atoms = FaceCenteredCubic(directions=[[1,0,0],[0,1,0],[0,0,1]], size=(6,6,6),
                              symbol="Cu", pbc=(1,1,0))
    r = atoms.get_positions()
    uc = atoms.get_cell()
    r[-1,coordinate] = uc[coordinate,coordinate] * 1.51
    atoms.set_positions(r)
    atoms.set_calculator(EMT())
    ecorrect = atoms.get_potential_energy()

    atoms.set_pbc((0,1,0))
    atoms.set_pbc((1,1,0))
    e2 = atoms.get_potential_energy()
    ReportTest("e2 correct", e2, ecorrect, 0.001)

    atoms.set_pbc((0,1,0))
    dummy =  atoms.get_potential_energy()
    assert(fabs(dummy - ecorrect) > 1.0)
    atoms.set_pbc((1,1,0))
    e3 = atoms.get_potential_energy()

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import get_cell [as 别名]
#some algebra to determine surface normal and the plane of the surface
d3=[2,1,1]
a1=np.array([0,1,1])
d1=np.cross(a1,d3)
a2=np.array([0,-1,1])
d2=np.cross(a2,d3)

#create your slab
slab  =FaceCenteredCubic(directions=[d1,d2,d3],
                         size=(2,1,2),
                         symbol=('Pt'),
                         latticeconstant=3.9)

#add some vacuum to your slab
uc = slab.get_cell()
print(uc)
uc[2] += [0,0,10]  #there are ten layers of vacuum
uc = slab.set_cell(uc,scale_atoms=False)
#view the slab to make sure it is how you expect
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

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import get_cell [as 别名]
            for i2, p2 in enumerate(positions[:i1]):
                diff = p2 - p1
                d2 = np.dot(diff, diff)
                c6 = (self.sigma**2 / d2)**3
                c12 = c6**2
                if d2 < self.cutoff**2:
                    self.energy += 4 * self.epsilon * (c12 - c6) - self.shift
                F = 24 * self.epsilon * (2 * c12 - c6) / d2 * diff
                self._forces[i1] -= F
                self._forces[i2] += F
        self.positions = positions.copy()



N = 5
ar = FaceCenteredCubic('Ar', pbc=[(0,0,0)], directions=[[1,0,0],[0,1,0],[0,0,1]], size=[N,N,N])
print ar.get_cell()
#view(ar) 

calc1 = KIMCalculator("ex_model_Ar_P_LJ")
ar.set_calculator(calc1)
kim_energy = ar.get_potential_energy()
print "kim energy = ", kim_energy 

calc2 = LennardJones(epsilon=epsilon, sigma=sigma, cutoff=cutoff)
ar.set_calculator(calc2)
ase_energy = ar.get_potential_energy()
print "ase energy = ", ase_energy 

print "difference = ", kim_energy - ase_energy

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import get_cell [as 别名]
import numpy as np
from ase.lattice.cubic import FaceCenteredCubic
ag = FaceCenteredCubic(directions=[[1, 0, 0],
                                   [0, 1, 0],
                                   [0, 0, 1]],
                       size=(1, 1, 1),
                       symbol='Ag',
                       latticeconstant=4.0)
# these are the reciprocal lattice vectors
b1, b2, b3 = np.linalg.inv(ag.get_cell())
'''
g(111) = 1*b1 + 1*b2 + 1*b3
and |g(111)| = 1/d_111
'''
h,k,l = (1, 1, 1)
d = 1./np.linalg.norm(h*b1 + k*b2 + l*b3)
print('d_111 spacing (method 1) = {0:1.3f} Angstroms'.format(d))
# method #2
hkl = np.array([h, k, l])
G = np.array([b1, b2, b3])  # reciprocal unit cell
'''
Gstar is usually defined as this matrix of dot products:
Gstar = np.array([[dot(b1,b1), dot(b1,b2), dot(b1,b3)],
                  [dot(b1,b2), dot(b2,b2), dot(b2,b3)],
                  [dot(b1,b3), dot(b2,b3), dot(b3,b3)]])
but I prefer the notationally more compact:
Gstar = G .dot. transpose(G)
then, 1/d_hkl^2 = hkl .dot. Gstar .dot. hkl
'''
Gstar = np.dot(G, G.T)
id2 = np.dot(hkl, np.dot(Gstar, hkl))