在NumPy中旋转水晶

Question

我想把一个由X，Y，Z原子坐标定义的5原子晶体旋转一个随机角度。我最初的想法是使用一个外部包来生成一个旋转矩阵（https://github.com/qobilidop/randrot），然后将这个矩阵乘以一个向量，这个向量定义了单个原子的坐标。然而，这根本不起作用，所有的原子都分散了。我为此写了一个函数：

def rotation():
    crystal = []
    rotmat = np.asarray(randrot.generate(3)) #generates 3x3 rotation matrix
    for x,y,z in zip(new_x, new_y, new_z):
        vec = np.array([x,y,z])
        rot = vec.dot(rotmat)
        for elem in rot:
            crystal.append(elem)
    return np.array(crystal).reshape([5,3])

rotated = rotation()
ax.scatter(rotated[0], rotated[1], rotated[2], marker='.', s=100, color='green')

下面是它的外观（红色是初始位置，绿色是旋转后）：

pyplot

Answer 1

下面是一个示例代码，它围绕随机生成的旋转矩阵旋转给定的3d点，旋转矩阵的创建取自another answer。你知道吗

import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import math


#  taken from https://stackoverflow.com/questions/6802577/rotation-of-3d-vector
def rotation_matrix(axis, theta):
    """
    Return the rotation matrix associated with counterclockwise rotation about
    the given axis by theta radians.
    """
    axis = np.asarray(axis)
    axis = axis / math.sqrt(np.dot(axis, axis))
    a = math.cos(theta / 2.0)
    b, c, d = -axis * math.sin(theta / 2.0)
    aa, bb, cc, dd = a * a, b * b, c * c, d * d
    bc, ad, ac, ab, bd, cd = b * c, a * d, a * c, a * b, b * d, c * d
    return np.array([[aa + bb - cc - dd, 2 * (bc + ad), 2 * (bd - ac)],
                     [2 * (bc - ad), aa + cc - bb - dd, 2 * (cd + ab)],
                     [2 * (bd + ac), 2 * (cd - ab), aa + dd - bb - cc]])


#  initial xyz coordinates
xs = [0, 1, 1, 1, 1, -1, -1, -1, -1]
ys = [0, 1, 1, -1, -1, 1, 1, -1, -1]
zs = [0, 1, -1, 1, -1, 1, -1, 1, -1]
atoms_initial = np.array([xs, ys, zs]).T

#  specify rotation matrix parameters
#  let us generate a random axis and angle for rotation
rotation_axis = np.random.uniform(low=0, high=1, size=3)  #  three numbers between 0 and 1
rotation_angle = np.random.uniform(low=0, high=2*np.pi, size=1)  #  random number between 0 and 2pi
print("Rotation axis:{}, rotation angle:{} radians".format(rotation_axis, rotation_angle))

#  create our rotation matrix
rotmat = rotation_matrix(rotation_axis, rotation_angle)

#  apply rotation matrix to our points
atoms_rotated = np.dot(atoms_initial, rotmat)

#  draw
fig = plt.figure()
ax = Axes3D(fig)
ax.scatter(atoms_initial[:,0], atoms_initial[:,1], atoms_initial[:,2], marker='.', s=100, color='red')
ax.scatter(atoms_rotated[:,0], atoms_rotated[:,1], atoms_rotated[:,2], marker='.', s=100, color="green")
plt.show()