Python FaceCenteredCubic.set_positions方法代码示例

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import set_positions [as 别名]
def myfcc(symbol, pbc):
    "Make an fcc lattice with standard lattice constant"
    if ismaster:
        a = FaceCenteredCubic(symbol=symbol, size=(10,10,10), pbc=pbc)    
        dx = 0.01 * np.sin(0.1 * np.arange(len(a) * 3))
        dx.shape = (len(a),3)
        a.set_positions(a.get_positions() + dx)
        a.set_momenta(np.zeros((len(a),3)))
        #view(a)
    else:
        a = None
    if isparallel:
        a = MakeParallelAtoms(a, cpulayout)
    return a

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import set_positions [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 set_positions [as 别名]
    print "    This is an s50 node on Niflheim."
    fullhost = "niflheim-s50/%s" % (host.split(".")[0])
    host = "niflheim-s50"
else:
    fullhost = host
print "Current time is "+when
print ""

print "Preparing system"
initial = FaceCenteredCubic(directions=[[1,0,0],[0,1,0],[0,0,1]],
                            size=(30, 30, 30),
                            symbol="Pt")
ReportTest("Number of atoms", len(initial), 108000, 0)
r = initial.get_positions()
r.flat[:] += 0.14 * sin(arange(3*len(initial)))
initial.set_positions(r)

print "Running self-test."
atoms = Atoms(initial)
atoms.set_calculator(EMT2013(PtY_parameters))
e = atoms.get_potential_energies()
f = atoms.get_forces()
if os.access(selfcheckfilename, os.F_OK):
    olde, oldf = cPickle.load(open(selfcheckfilename))
    de = max(fabs(e - olde))
    df = max(fabs(f.flat[:] - oldf.flat[:]))
    print "Maximal deviation:  Energy", de, "  Force", df
    ReportTest("Max force error", df, 0.0, 1e-11)
    ReportTest("Max energy error", de, 0.0, 1e-11)
    del olde, oldf
else:

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import set_positions [as 别名]
numbers = count.keys()
numbers.sort()
sum = 0
for i in numbers:
    #print i, count[i]
    sum += i*count[i]

ReportTest("Number of neighbors (EMT's NB list)", sum, 21*len(atoms), 0)

nblist = NeighborList(latconst * 0.5 * (1/sqrt(2) + 1), atoms, 0.0)
#nblist = NeighborCellLocator(latconst * 0.5 * (1/sqrt(2) + 1), atoms, 0.0)
fnb = FullNeighborList(latconst * 0.5 * (1/sqrt(2) + 1), Atoms(atoms))
TestLists(nblist, fnb, "nearest-neigbor lists (periodic)", 6)

ReportTest("Energy unperturbed 1", atoms.get_potential_energy(), epot, 1e-11)
atoms.set_positions(atoms.get_positions())
ReportTest("Energy unperturbed 2", atoms.get_potential_energy(), epot, 1e-11)

nblist = NeighborList(4.98409, atoms, 0.0)
fnb = FullNeighborList(4.98409, Atoms(atoms))
TestLists(nblist, fnb, "long neigbor lists (periodic)", 21)

ReportTest("Energy unperturbed 3", atoms.get_potential_energy(), epot, 1e-11)
atoms.set_positions(atoms.get_positions())
ReportTest("Energy unperturbed 4", atoms.get_potential_energy(), epot, 1e-11)

atoms = Atoms(atoms, pbc=(0,0,0))

nblist = NeighborList(latconst * 0.5 * (1/sqrt(2) + 1), atoms, 0.0)
fnb = FullNeighborList(latconst * 0.5 * (1/sqrt(2) + 1), Atoms(atoms))
TestLists(nblist, fnb, "nearest-neigbor lists (non-periodic)")

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import set_positions [as 别名]
    print sum(a), "near old position and", sum(b), "near new position."
    print nr2, "atoms should be affected"
    ReportTest("Number of affected atoms", nr, nr2, 0)
    
    del np_atoms, np_nblist
    
if 1:
    print
    print "Testing perturbations of all the atoms"

    magnarr = array((0.0, 0.01, 0.1, 1.0, 3.0, 10.0))
    pick = argsort(random.random((len(magnarr),)))
    for magn in take(magnarr, pick):
        dx = magn * random.uniform(-1, 1, (len(atoms), 3))
        print "  Perturbation:", magn
        atoms.set_positions(dx + atoms.get_positions())
        nblist.check_and_update(atoms)
        checklist(nblist, atoms, listcutoff)



print 
print "Testing perturbations of single atoms"

magnarr = array((0.0, 0.01, 0.1, 1.0, 3.0, 10.0))
pick1 = argsort(random.random((len(magnarr),)))
for magn in take(magnarr, pick1):
    numarr = array((1, 3, 10, len(atoms)/2, -3))
    pick2 = argsort(random.random((len(numarr),)))
    for number in take(numarr, pick2):
        # Pick number random atoms.

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import set_positions [as 别名]
atoms = FaceCenteredCubic(size=(5,5,5), symbol="Cu")
pot = EMT()
atoms.set_calculator(pot)

for i in (1, 2):
    print "*** Pass", i
    ReportTest("Energy required", pot.calculation_required(atoms, ["energy"]), 1, 0)
    ReportTest("Forces required", pot.calculation_required(atoms, ["forces"]), 1, 0)
    ReportTest("Stress required", pot.calculation_required(atoms, ["stress"]), 1, 0)
    ReportTest("Magmom required", pot.calculation_required(atoms, ["magmoms"]), 1, 0)
    e = atoms.get_potential_energy()
    ReportTest("Energy not required", pot.calculation_required(atoms, ["energy"]), 0, 0)
    ReportTest("Forces required (II)", pot.calculation_required(atoms, ["forces"]), 1, 0)
    f = atoms.get_forces()
    ReportTest("Energy not required (II)", pot.calculation_required(atoms, ["energy"]), 0, 0)
    ReportTest("Forces not required", pot.calculation_required(atoms, ["forces"]), 0, 0)
    ReportTest("Energy or forces not required",
               pot.calculation_required(atoms, ["energy", "forces"]), 0, 0)
    ReportTest("Energy or stress required",
               pot.calculation_required(atoms, ["energy", "stress"]), 1, 0)
    s = atoms.get_stress()
    ReportTest("Stress not required", pot.calculation_required(atoms, ["stress"]), 0, 0)

    r = atoms.get_positions()
    r[0,0] += 0.1
    atoms.set_positions(r)

ReportTest.Summary()

# 需要导入模块: from ase.lattice.cubic import FaceCenteredCubic [as 别名]
# 或者: from ase.lattice.cubic.FaceCenteredCubic import set_positions [as 别名]
     boundaries = ((0,0,0),)
 else:
     boundaries = ((1,1,1), (0,0,0), (0,0,1))
     
 for pbc in boundaries:
     txt = nbltype + (" PBC=%1i%1i%1i " % pbc)
     # Test that EMT reimported through OpenKIM gives the right results.
     atoms_kim = FaceCenteredCubic(size=(10,10,10), symbol='Cu')
     #atoms_kim = FaceCenteredCubic(directions=[[1,0,0],[0,1,0],[0,0,1]],
     #                    size=(30, 30, 30),
     #                    symbol="Cu")
     natoms = len(atoms_kim)
     atoms_kim.set_pbc(pbc)
     r = atoms_kim.get_positions()
     r.flat[:] += 0.1 * np.sin(np.arange(3*natoms))
     atoms_kim.set_positions(r)
     atoms_emt = atoms_kim.copy()
     kim = OpenKIMcalculator(openkimmodel, allowed=nbltype)
     emt = EMT()
     emt.set_subtractE0(False)
     atoms_kim.set_calculator(kim)
     atoms_emt.set_calculator(emt)
     ek = atoms_kim.get_potential_energy()
     ee = atoms_emt.get_potential_energy()
     ReportTest(txt+"Total energy", ek, ee, 1e-8)
     ek = atoms_kim.get_potential_energies()
     ee = atoms_emt.get_potential_energies()
     for i in range(0, natoms, step):
         ReportTest(txt+"Energy of atom %i" % (i,), ek[i], ee[i], 1e-8)
     fk = atoms_kim.get_forces()
     fe = atoms_emt.get_forces()