如何使用pymatgen獲得結構內的所有鍵角？

Question

因此，使用pymatgen，我有一個結構對象。我想要做的是獲得結構內的所有鍵角。我可以遍歷每個原子以獲得所有鍵角，但是這將包括每個原子，無論它們有多遠，哪一個當然是個問題。

現在，我可以使用「get_neighbors」函數找到每個中心原子的鄰居，但是，我不確定該從哪裏去，尤其是因爲「get_angle」函數採用整數值而不是原子位置對象。

下面是我迄今爲止代碼：

import pymatgen as mg 
import numpy as np 


s = mg.Structure.from_file('POSCAR') 

atoms = s.atomic_numbers 

van = [x for x in atoms if x == 23] 
length = len(van) 
nb = ['NONE']*length 

x = 0 

while x < length: 

    nb[x] = s.get_neighbors(s[x],2.4) 
    x += 1 

所以我有鄰居的一個數組，我知道他們所對應的原子也一樣，現在我只是需要得到的角度對所有相鄰原子。

任何幫助將不勝感激。

Answer 1

更新：

所以，我想出了一個這樣做的方法。我實際上將每個相鄰原子轉換成它們的笛卡爾座標。然後，我減去相鄰原子所連接的中心原子的座標。最後，我通過計算每兩個矢量的點積除以它們的幅度乘積來得到角度，然後取這些值的反餘弦。我使用的代碼如下。這可能不是最優雅的做事方式，但它確實完成了工作。如果有人有改進，請隨時發表評論。

import random 
import itertools as iter 
import pymatgen as mg 
import numpy as np 
import math 


def all_angles(POSCAR,amin,amax): 

    s = mg.Structure.from_file(POSCAR) 

    atoms = s.atomic_numbers 

    van = [x for x in atoms if x == 23] 
    length = len(van) 
    nb = ['NONE']*length 

    x = 0 
    n = 0 

    while x < length: 

     nb[x] = s.get_neighbors(s[x],2.4) 
     x += 1 

    w = 100 
    h = len(van) 
    oxygen = [[0 for x in range(w)] for y in range(h)] 

    x = 0 
    y = 0 

    while x < len(nb): 

     y = 0 
     n = 0 

     while y < len(nb[x]): 

      oxygen[x][n] = (nb[x][y][0]).coords 
      n += 1 
      y += 1 

     x += 1 

    x = 0 

    while x < len(oxygen): 
     oxygen[x] = [n for n in oxygen[x] if not isinstance(n,int)] 
     x += 1 

    van = [x for x in range(0,len(van))] 

    x = 0 
    while x < len(van): 
     van[x] = s[van[x]].coords 
     x += 1 

    van = [np.array(x) for x in van] 

    x = 0 
    while x < len(van): 

     oxygen[x] = [np.subtract(oxygen[x][y],van[x]) for y in range(0,len(oxygen[x]))] 
     x += 1 

    combo = [[0 for x in range(0,1000)] for y in range(0,len(van))] 

    r = 0 

    while r < len(van): 
     x = 0 
     for subset in iter.combinations(oxygen[r],2): 
      combo[r][x] = subset 
      x += 1 
     r += 1 

    x = 0 
    while x < len(combo): 
     combo[x] = [c for c in combo[x] if c != 0] 
     x += 1 


    angles = [[0 for x in range(0,1000)] for y in range(0,len(van))] 

    x = 0 

    while x < len(combo): 

     group = combo[x] 

     y = 0 

     while y < len(group): 
      angles[x][y] = math.degrees(math.acos(np.dot(group[y][0],group[y][1])/(np.linalg.norm(group[y][0])*np.linalg.norm(group[y][1])))) 
      y += 1 

     angles[x] = [round(n,3) for n in angles[x] if n != 0 and n > amin and n < amax] 

     x += 1 

    angles = np.concatenate(angles,axis=0)  

    return angles