Python Atom.set_sid方法代码示例

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
 def from_ase_atoms(self,atoms):
     """
     Convert ASE Atoms object to NAPSystem object.
     """
     self.a1 = np.array(atoms.cell[0])
     self.a2 = np.array(atoms.cell[1])
     self.a3 = np.array(atoms.cell[2])
     spos = atoms.get_scaled_positions()
     symbols = atoms.get_chemical_symbols()
     #...initialize and remake self.specorder
     self.specorder = []
     for s in symbols:
         if s not in self.specorder:
             self.specorder.append(s)
     #...first, initialize atoms array
     self.atoms = []
     #...append each atom from ASE-Atoms
     for ia,spi in enumerate(spos):
         si = symbols[ia]
         ai = Atom()
         sid = self.specorder.index(si)+1
         ai.set_id(ia+1)
         ai.set_sid(sid)
         ai.set_symbol(si)
         ai.set_pos(spi[0],spi[1],spi[2])
         ai.set_vel(0.,0.,0.)
         self.atoms.append(ai)
     return

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
def structure2aSys(structure,idoffset=1):
    """
    Converts Structure object of pymatgen to NAPSystem object in nap.

    Args:
        structure (Structure): pymatgen Structure object to be converted
        to NAPSystem object..

    Returns:
        aSys (NAPSystem): 
    """
    lattice= structure.lattice
    alc= 1.0
    a1= np.array(lattice.matrix[0])
    a2= np.array(lattice.matrix[1])
    a3= np.array(lattice.matrix[2])
    #... rescale a? vectors
    a1= a1/alc
    a2= a2/alc
    a3= a3/alc
    aSys= NAPSystem()
    aSys.set_lattice(alc,a1,a2,a3)
    for ia in range(structure.num_sites):
        ai= Atom()
        si= structure[ia]
        crd= si.frac_coords
        ai.set_pos(crd[0],crd[1],crd[2])
        sid= structure.symbol_set.index(si.species_string)+idoffset
        ai.set_sid(sid)
        ai.set_id(ia+1)
        aSys.add_atom(ai)
    return aSys

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
 def read_akr(self,fname='akr0000'):
     f=open(fname,'r')
     # 1st: lattice constant
     self.alc= float(f.readline().split()[0])
     # 2nd-4th: cell vectors
     for i in range(3):
         data= f.readline().split()
         self.a1[i]= float(data[0])
         self.a2[i]= float(data[1])
         self.a3[i]= float(data[2])
     # 5th: num of atoms
     natm= int(f.readline().split()[0])
     # 9th-: atom positions
     self.atoms= []
     symbol = None
     for i in range(natm):
         data= [float(x) for x in f.readline().split()]
         ai= Atom()
         ai.set_sid(data[0])
         ai.set_pos(data[1],data[2],data[3])
         ai.set_vel(data[4],data[5],data[6])
         if self.specorder:
             symbol = self.specorder[ai.sid-1]
         if symbol and ai.symbol != symbol:
             ai.set_symbol(symbol)
         self.atoms.append(ai)
     f.close()

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
 def read_POSCAR(self,fname='POSCAR'):
     with open(fname,'r') as f:
         # 1st line: comment
         f.readline()
         # 2nd: lattice constant
         self.alc= float(f.readline().split()[0])
         # 3rd-5th: cell vectors
         self.a1= np.array([float(x) for x in f.readline().split()])
         self.a2= np.array([float(x) for x in f.readline().split()])
         self.a3= np.array([float(x) for x in f.readline().split()])
         # 6th: species names or number of each species
         buff= f.readline().split()
         if not buff[0].isdigit():
             spcs = copy.deepcopy(buff)
             buff= f.readline().split()
             if not self.specorder:
                 self.specorder = spcs
         num_species= np.array([ int(n) for n in buff])
         natm= 0
         for n in num_species:
             natm += n
         #print("Number of atoms = {0:5d}".format(natm))
         # 7th or 8th line: comment
         c7= f.readline()
         if c7[0] in ('s','S'):
             c8= f.readline()
         # Atom positions hereafter
         self.atoms= []
         for i in range(natm):
             buff= f.readline().split()
             ai= Atom()
             sid= 1
             m= 0
             sindex=0
             symbol = None
             for n in num_species:
                 m += n
                 if i < m:
                     if spcs and self.specorder:
                         sid = self.specorder.index(spcs[sindex]) + 1
                         symbol = spcs[sindex]
                     break
                 sid += 1
                 sindex += 1
             ai.set_id(i+1)
             ai.set_sid(sid)
             if symbol:
                 ai.symbol = symbol
             ai.set_pos(float(buff[0]),float(buff[1]),float(buff[2]))
             ai.set_vel(0.0,0.0,0.0)
             self.atoms.append(ai)

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
def doc_to_pos(doc,conf):
    """
    Make a pos file, which has pmd format, from a document in MongoDB.
    """
    psys= NAPSystem()
    matrix=doc['calculations'][-1]['output']['crystal']['lattice']['matrix']
    a1= matrix[0]
    a2= matrix[1]
    a3= matrix[2]
    psys.set_lattice(1.0,a1,a2,a3)

    species_ids=conf['species_ids']

    sites= doc['calculations'][-1]['output']['crystal']['sites']
    for site in sites:
        ra= site['abc']
        ai= Atom()
        ai.set_pos(ra[0],ra[1],ra[2])
        ai.set_sid(species_ids[site['species'][0]['element']])
        psys.add_atom(ai)
    return psys

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
print(' sid2          = ',sid2)
pmdexec= options.pmdexec

asys= NAPSystem()
# system size is bigger than 2*rcut
a1= np.array([2.0, 0.0, 0.0])
a2= np.array([0.0, 1.0, 0.0])
a3= np.array([0.0, 0.0, 1.0])
alc= rcut
asys.set_lattice(alc,a1,a2,a3)

atom1= Atom()
atom2= Atom()
atom1.set_pos(0.0,0.0,0.0)
atom1.set_id(1)
atom1.set_sid(sid1)
asys.add_atom(atom1)
atom2.set_pos(0.5,0.0,0.0)
atom2.set_id(2)
atom2.set_sid(sid2)
asys.add_atom(atom2)

hmin= rmin/(2*rcut)
hd  = (0.5-hmin)/nsmpl

os.system('cp pmdini pmdorig')

fout= open('out.2body','w')
for ip in range(nsmpl+1):
    print('.',end='')
    d= hmin +hd*ip

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
                    ,np.multiply(sys0.a1,r1)
                    ,np.multiply(sys0.a2,r2)
                    ,np.multiply(sys0.a3,r3))
hext=np.zeros((3,3))
hext[0]= sysext.a1 *sysext.alc
hext[1]= sysext.a2 *sysext.alc
hext[2]= sysext.a3 *sysext.alc

#...make extended system corresponding to the phonopy definition
for ia in range(natm0):
    ai0= sys0.atoms[ia]
    for i3 in range(r3):
        for i2 in range(r2):
            for i1 in range(r1):
                ai= Atom()
                ai.set_sid(ai0.sid)
                p1= (ai0.pos[0]+i1)/r1
                p2= (ai0.pos[1]+i2)/r2
                p3= (ai0.pos[2]+i3)/r3
                ai.set_pos(p1,p2,p3)
                sysext.add_atom(ai)
sysext.reset_ids()
natme= sysext.num_atoms()
print ' sysext.num_atoms()=',natme

# for ai in sysext.atoms:
#     print ai.pos

os.system('cp pmd0000 pmd0000.orig')

#...loop for all atoms in the extended system

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
def make_polycrystal(grns, uc, n1, n2, n3):
    u"""
    This routine is not that universal.
    Each grain has to have neighboring grains within a supercell,
    otherwise there will be some unexpecting grain boundries.
    In order to do so, the system should be large enough and
    the number of grains should be large enough.
    """
    sv = shift_vector()
    system = AtomSystem()
    system.set_lattice(uc.alc, uc.a1 * n1, uc.a2 * n2, uc.a3 * n3)
    hmat = np.zeros((3, 3))
    hmat[0] = system.a1 * system.alc
    hmat[1] = system.a2 * system.alc
    hmat[2] = system.a3 * system.alc
    hmati = np.linalg.inv(hmat)
    for ig in range(len(grns)):
        grain = grns[ig]
        rmat = grain.rmat
        pi = grain.point
        api = np.dot(hmat, pi)
        for ix in range(-n1 / 2 - 1, n1 / 2 + 2):
            print "ix=", ix
            for iy in range(-n2 / 2 - 1, n2 / 2 + 2):
                for iz in range(-n3 / 2 - 1, n3 / 2 + 2):
                    for m in range(len(uc.atoms)):
                        rt = np.zeros((3,))
                        rt[0] = (uc.atoms[m].pos[0] + ix) / n1
                        rt[1] = (uc.atoms[m].pos[1] + iy) / n2
                        rt[2] = (uc.atoms[m].pos[2] + iz) / n3
                        # ...rt to absolute position
                        art = np.dot(hmat, rt)
                        # ...rotate
                        ari = np.dot(rmat, art)
                        ari[0] = ari[0] + api[0]
                        ari[1] = ari[1] + api[1]
                        ari[2] = ari[2] + api[2]
                        # #...reduce into the supercell
                        # ri= np.dot(hmati,ari)
                        # ri[0]= pbc(ri[0])
                        # ri[1]= pbc(ri[1])
                        # ri[2]= pbc(ri[2])
                        # ari= np.dot(hmat,ri)
                        # ...check distance from all the grain points
                        di = distance(ari, api)
                        isOutside = False
                        for jg in range(len(grns)):
                            gj = grns[jg]
                            for isv in range(27):
                                pj = gj.point
                                if jg == ig and isv == 13:
                                    continue
                                svi = sv[isv]
                                pj = pj + svi
                                apj = np.dot(hmat, pj)
                                dj = distance(ari, apj)
                                if dj < di:
                                    isOutside = True
                                    break
                            if isOutside:
                                break
                        if isOutside:
                            break
                        # ...here ri is inside this grain, register it
                        atom = Atom()
                        ri = np.dot(hmati, ari)
                        ri[0] = pbc(ri[0])
                        ri[1] = pbc(ri[1])
                        ri[2] = pbc(ri[2])
                        atom.set_pos(ri[0], ri[1], ri[2])
                        atom.set_sid(uc.atoms[m].sid)
                        system.add_atom(atom)
    # ...remove too-close atoms at the grain boundaries
    system.make_pair_list(rcut)
    ls_remove = []
    dmin2 = dmin ** 2
    xij = np.zeros((3,))
    for ia in range(system.num_atoms()):
        ai = system.atoms[ia]
        pi = ai.pos
        nlst = system.nlspr[ia]
        lst = system.lspr[ia]
        for j in range(nlst):
            ja = lst[j]
            aj = system.atoms[ja]
            pj = aj.pos
            xij[0] = pj[0] - pi[0] - anint(pj[0] - pi[0])
            xij[1] = pj[1] - pi[1] - anint(pj[1] - pi[1])
            xij[2] = pj[2] - pi[2] - anint(pj[2] - pi[2])
            xij = np.dot(hmat, xij)
            d2 = xij[0] ** 2 + xij[1] ** 2 + xij[2] ** 2
            if d2 < dmin2:
                ls_remove.append(ia)
                ls_remove.append(ja)
    # ...one of two will survive
    # print ls_remove
    count = [ls_remove.count(ls_remove[i]) for i in range(len(ls_remove))]
    # print count
    for i in range(0, len(ls_remove), 2):
        if count[i] > count[i + 1]:
#.........这里部分代码省略.........

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
 def read_xsf(self,fname="xsf"):
     f=open(fname,'r')
     mode= 'None'
     ixyz= 0
     iatm= 0
     self.atoms= []
     natm= 0
     for line in f.readlines():
         if 'CRYSTAL' in line:
             mode= 'CRYSTAL'
             continue
         elif 'PRIMVEC' in line:
             mode= 'PRIMVEC'
             continue
         elif 'PRIMCOORD' in line:
             mode= 'PRIMCOORD'
             # Before going further, create inversed h-matrix
             hi = unitvec_to_hi(self.a1,self.a2,self.a3)
             # print 'Inversed h-matrix:'
             # print hi
             continue
         
         if mode == 'CRYSTAL':
             pass
         elif mode == 'PRIMVEC':
             if ixyz == 0:
                 arr = [ float(x) for x in line.split() ]
                 self.a1[0] = arr[0]
                 self.a1[1] = arr[1]
                 self.a1[2] = arr[2]
             elif ixyz == 1:
                 arr = [ float(x) for x in line.split() ]
                 self.a2[0] = arr[0]
                 self.a2[1] = arr[1]
                 self.a2[2] = arr[2]
             elif ixyz == 2:
                 arr = [ float(x) for x in line.split() ]
                 self.a3[0] = arr[0]
                 self.a3[1] = arr[1]
                 self.a3[2] = arr[2]
             ixyz += 1
         elif mode == 'PRIMCOORD':
             data = line.split()
             if len(data) == 1:
                 natm= int(data[0])
                 continue
             elif len(data) == 2:
                 natm= int(data[0])
                 nspcs= int(data[1])
                 continue
             elif len(data) == 4 or len(data) == 7:
                 if iatm >= natm:
                     continue
                 symbol = get_symbol_from_number(int(data[0]))
                 if symbol not in self.specorder:
                     self.specorder.append(symbol)
                 sid = self.specorder.index(symbol) +1
                 ai= Atom()
                 ai.symbol = symbol
                 ai.set_sid(sid)
                 xc= float(data[1])
                 yc= float(data[2])
                 zc= float(data[3])
                 xi,yi,zi = cartessian_to_scaled(hi,xc,yc,zc)
                 ai.set_pos(xi,yi,zi)
                 ai.set_vel(0.0,0.0,0.0)
                 self.atoms.append(ai)
                 # print 'iatm,symbol,sid,xc,yc,zc = ',iatm,symbol,sid,xc,yc,zc
             else:
                 continue
             iatm += 1
     self.alc= 1.0
     # print self.alc
     # print self.a1[:]
     # print self.a2[:]
     # print self.a3[:]
     f.close()

# 需要导入模块: from atom import Atom [as 别名]
# 或者: from atom.Atom import set_sid [as 别名]
 def read_dump(self,fname="dump"):
     f=open(fname,'r')
     mode= 'None'
     ixyz= 0
     iatm= 0
     symbol = None
     self.atoms= []
     self.alc= 1.0
     for line in f.readlines():
         if 'ITEM: NUMBER OF ATOMS' in line:
             mode= 'NUMBER OF ATOMS'
             continue
         elif 'ITEM: BOX BOUNDS' in line:
             mode= 'BOX BOUNDS'
             continue
         elif 'ITEM: ATOMS' in line:
             mode= 'ATOMS'
             continue
         
         if mode == 'NUMBER OF ATOMS':
             natm= int(line.split()[0])
         elif mode == 'BOX BOUNDS':
             data = line.split()
             if ixyz == 0:
                 xlo_bound= float(data[0])
                 xhi_bound= float(data[1])
                 if len(data) > 2:
                     xy = float(data[2])
             elif ixyz == 1:
                 ylo_bound= float(line.split()[0])
                 yhi_bound= float(line.split()[1])
                 if len(data) > 2:
                     xz = float(data[2])
             elif ixyz == 2:
                 zlo_bound= float(line.split()[0])
                 zhi_bound= float(line.split()[1])
                 if len(data) > 2:
                     yz = float(data[2])
             ixyz += 1
             if ixyz > 2:
                 xlo = xlo_bound -min(0.0,xy,xz,xy+yz)
                 xhi = xhi_bound -max(0.0,xy,xz,xy+yz)
                 ylo = ylo_bound -min(0.0,yz)
                 yhi = yhi_bound -max(0.0,yz)
                 zlo = zlo_bound
                 zhi = zhi_bound
                 self.a1 = np.array([xhi-xlo,xy,xz],dtype=float)
                 self.a2 = np.array([0.0,yhi-ylo,yz],dtype=float)
                 self.a3 = np.array([0.0,0.0,zhi-zlo],dtype=float)
                 hmat = self.get_hmat()
                 hmati= np.linalg.inv(hmat)
         elif mode == 'ATOMS':
             if iatm < natm:
                 data= line.split()
                 ai= Atom()
                 ai.set_sid(int(data[1]))
                 if self.specorder:
                     symbol = self.specorder[ai.sid-1]
                 if symbol and ai.symbol != symbol:
                     ai.set_symbol(symbol)
                 x0= float(data[2])
                 y0= float(data[3])
                 z0= float(data[4])
                 x = hmati[0,0]*x0 +hmati[0,1]*y0 +hmati[0,2]*z0
                 y = hmati[1,0]*x0 +hmati[1,1]*y0 +hmati[1,2]*z0
                 z = hmati[2,0]*x0 +hmati[2,1]*y0 +hmati[2,2]*z0
                 x = self._pbc(x)
                 y = self._pbc(y)
                 z = self._pbc(z)
                 ai.set_pos(x,y,z)
                 ai.set_vel(0.0,0.0,0.0)
                 self.atoms.append(ai)
             iatm += 1
     # print self.alc
     # print self.a1[:]
     # print self.a2[:]
     # print self.a3[:]
     f.close()