Создание равномерного случайного кватерниона и умножение двух кватернионов

У меня есть функция Python (NumPy), которая создает равномерный случайный кватернион. Я хотел бы получить два кватернионных умножения в виде двухмерного возвращаемого массива из той же или другой функции. В моем недавнем случае формула умножения кватернионов имеет вид Q1 * Q2 и Q2 * Q1. Вот,Q1=(w0, x0, y0, z0) а такжеQ2=(w1, x1, y1, z1) два кватерниона. Ожидаемый результат умножения двух кватернионов (в виде двумерного возвращаемого массива) должен быть

return([-x1*x0 - y1*y0 - z1*z0 + w1*w0, x1*w0 + y1*z0 - z1*y0 +
    w1*x0, -x1*z0 + y1*w0 + z1*x0 + w1*y0, x1*y0 - y1*x0 + z1*w0 +
    w1*z0])

Кто-нибудь может мне помочь? Мои коды здесь:

def randQ(N):
    #Generates a uniform random quaternion
    #James J. Kuffner 2004 
    #A random array 3xN
    s = random.rand(3,N)
    sigma1 = sqrt(1.0 - s[0])
    sigma2 = sqrt(s[0])
    theta1 = 2*pi*s[1]
    theta2 = 2*pi*s[2]
    w = cos(theta2)*sigma2
    x = sin(theta1)*sigma1
    y = cos(theta1)*sigma1
    z = sin(theta2)*sigma2
    return array([w, x, y, z])