Ликбез по кватернионам, часть 5: практическая реализация поворота: nabbla1

nabbla1

Ликбез по кватернионам, часть 5: практическая реализация поворота

Feb 15, 2018 01:40

Поздравим с днём рождения Виктора Чапаева victor_chapaev, который сподвиг меня написать данный "опус".

Оглавление "Ликбеза по кватернионам":
[Spoiler (click to open)]
Часть 1 - история вопроса
Часть 2 - основные операции
Часть 3 - запись вращения через кватернионы
Часть 4 - кватернионы и спиноры; порядок перемножения
Часть 5 - практическая реализация поворота
Часть 5 1/2 - введение метрики, "расстояния" между поворотами
часть 5 5/8 - метрика ненормированных кватернионов
часть 5 11/16 - красивая псевдометрика произвольных кватернионов
Часть 5 3/4 - исследуем "пространство поворотов"
Часть 5 7/8 - почти изотропный ёжик
Часть 6 - поворот по кратчайшему пути
Часть 6 1/4 - кратчайший поворот в общем случае
Часть 6 2/4 - поворот, совмещающий два направления
Часть 6 3/4 - кватернион из синуса и косинуса угла
Часть 6 7/8 - "уполовинивание угла" на плоскости
Часть 7 - интегрирование угловых скоростей, углы Эйлера-Крылова
Часть 8 - интегрирование угловых скоростей, матрицы поворота
Часть 8 1/2 - ортонормирование матрицы и уравнения Пуассона
Часть 9 - интегрирование угловых скоростей с помощью кватернионов
Часть 10 - интегрирование угловых скоростей, методы 2-го порядка
Часть 10 1/2 - интегрирование с поддержанием нормы
Часть 11 - интегрирование угловых скоростей, методы высших порядков (в разработке)
Часть 12 - навигационная задача
Часть 13 - Дэвенпорт берёт след!
Часть 14 - линейный метод Мортари-Маркли
Часть 15 - среднее от двух кватернионов
Часть 15 1/2 - проверка и усреднение кватернионов
Часть 16 - разложение кватерниона на повороты
Часть 17 - лидарная задача
Задачка к части 17
Дэвенпорт VS Мортари-Маркли
Мортари-Маркли берут реванш!
Дэвенпорт VS Мортари-Маркли, раунд 3

Мы выяснили, что такое кватернионы и почему они должны состоять именно из четырёх компонентов, узнали основные действия с ними, разобрались, как с помощью кватернионов описывается вращение в пространстве.

Теперь наконец-то поговорим о практическом их использовании - некоторые соображения, как осуществлять вращение векторов с помощью кватернионов, как интегрировать угловые скорости и как найти кратчайший поворот, переводящий один вектор (направление) в другой.

Чтобы опять не зависнуть на пару месяцев, буду писать по маленьким кусочкам. Сегодня будет уникальный материал, который вы вряд ли найдёте где-то ещё: в библиотеках раскиданы ссылки на форумы и блоги, где товарищи выводили эффективные формулы, но все эти ссылки уже мертвы. Чтобы не работать с формулами, оставленными мудрыми предками (или не воспринимать их как магические заклинания, которые надо просто выучить), автор вывел их заново и убедился в их корректности.

Реализация поворота вектора с помощью кватерниона
У нас есть кватернион Λ с единичной нормой:
$\Lambda=\lambda_0+\vec{\boldsymbol{\lambda}}=\lambda_0+\lambda_1{i}+\lambda_2{j}+\lambda_3{k}$

и некоторый вектор $\vec{\mathbf{u}}$
, который мы записываем как кватернион с нулевой скалярной (действительной) частью:
$\vec{\mathbf{u}}=u_1{i}+u_2{j}+u_3{k}$

Рассмотрим, как повернуть вектор наиболее эффективным способом.

«По определению»
Как мы уже знаем, поворот вектора можно произвести по формуле (3.2):

$\vec{\mathbf{v}}=\Lambda\vec{\mathbf{u}}\overline{\Lambda}$
(6.1)

Если не написать специализированной функции, учитывающей, что у вектора «нет скалярной части», это потребует 32 умножения и 24 сложения, тогда как обычное умножение матрицы 3х3 на вектор требует 9 умножений и 6 сложений. Далее в тексте, для удобочитаемости будем писать в виде: 32m + 24a, что означает «32 умножения и 24 сложения».
Можно слегка упростить выражения, учитывая нулевой скалярный компонент векторов $\vec{\mathbf{u}}$
и $\vec{\mathbf{v}}$
. Сначала вычислим промежуточный кватернион:

$M=\vec{\mathbf{u}}\overline{\Lambda}=(u_1i+u_2j+u_3k)(\lambda_0-\lambda_1i-\lambda_2j-\lambda_3k)=$

$={\color{Purple}u_1\lambda_0i}+{\color{Green}u_1\lambda_1}-u_1\lambda_2k+{\color{Blue}u_1\lambda_3j}+$

$+{\color{Blue}u_2\lambda_0j}+u_2\lambda_1k+{\color{Green}u_2\lambda_2}-{\color{Purple}u_2\lambda_3i}+$

$+u_3\lambda_0k-{\color{Blue}u_3\lambda_1j}+{\color{Purple}u_3\lambda_2i}+{\color{Green}u_3\lambda_3}=$

$=\mu_0+\mu_1i+\mu_2j+\mu_3k$

(6.2)

Выпишем выражения для компонентов кватерниона в явном виде:

$\mu_0=u_1\lambda_1+u_2\lambda_2+u_3\lambda_3$

$\mu_1=u_1\lambda_0-u_2\lambda_3+u_3\lambda_2$

$\mu_2=u_1\lambda_3+u_2\lambda_0-u_3\lambda_1$

$\mu_3=-u_1\lambda_2+u_2\lambda_1+u_3\lambda_0$

(6.3)

Убрав умножение на нулевую скалярную часть, мы снижаем количество операций до 12m+8a. Далее, вычислим итоговый вектор, учитывая, что его скалярная часть заведомо равна нулю:

$\vec{\mathbf{v}}=\Lambda\,M=(\lambda_0+\lambda_1i+\lambda_2j+\lambda_3k)(\mu_0+\mu_1i+\mu_2j+\mu_3k)=$

$={\color{Green}\lambda_0\mu_1i}+{\color{Blue}\lambda_0\mu_2j}+\lambda_0\mu_3k+$

$+{\color{Green}\lambda_1\mu_0i}+\lambda_1\mu_2k-{\color{Blue}\lambda_1\mu_3j}+$

$+{\color{Blue}\lambda_2\mu_0j}-\lambda_2\mu_1k+{\color{Green}\lambda_2\mu_3i}+$

$+\lambda_3\mu_0k+{\color{Blue}\lambda_3\mu_1j}-{\color{Green}\lambda_3\mu_2i}$

(6.4)

Выпишем выражения для компонентов вектора в явном виде:

$v_1=\lambda_0\mu_1+\lambda_1\mu_0+\lambda_2\mu_3-\lambda_3\mu_2$

$v_2=\lambda_0\mu_2-\lambda_1\mu_3+\lambda_2\mu_0+\lambda_3\mu_1$

$v_3=\lambda_0\mu_3+\lambda_1\mu_2-\lambda_2\mu_1+\lambda_3\mu_0$

(6.5)

Здесь мы обошлись 12m+9a. Итого, на поворот вектора требуется 24 умножения и 17 сложений - чуть лучше, чем до упрощения, но всё равно значительно более затратно, чем при умножении вектора на матрицу.

Магическая формула
Следующие формулы применены в нескольких библиотеках для работы с кватернионами, но ссылки на форум, где программист приводил вывод этих формул - давным-давно «мертвы». Для начала, сами формулы:

$\vec{\mathbf{t}}=2[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]$
,
$\vec{\mathbf{v}}=\vec{\mathbf{u}}+\lambda_0\vec{\mathbf{t}}+[\vec{\boldsymbol{\lambda}},\vec{\mathbf{t}}]$

(6.6)

Как видно, используется 2 векторных умножения, 2 умножения вектора на скаляр (один из них - и вовсе число 2) и два сложения векторов. Векторное умножение требует 6m+3a, умножение на скаляр - 3m, сложение векторов - 3a.
В общей сложности, поворот вектора по (6.6) требует 18 умножений и 12 сложений, что хуже обычного матричного умножения ровно вдвое, во всяком случае, по числу операций. С другой стороны, 4 компоненты кватерниона замечательно умещаются в один регистр xmm, и при хорошей реализации на SSE полное время выполнения (включая передачу аргументов функции и возврат результата) может оказаться даже ниже.
Попробуем вывести (6.6) геометрически. Такое доказательство будет наиболее наглядным.

Мы хотим повернуть вектор $\vec{\mathbf{u}}$
вокруг оси $\vec{\boldsymbol{\lambda}}$
на угол φ. Угол φ «закодирован» внутри кватерниона:
$\Lambda=\lambda_0+\vec{\boldsymbol{\lambda}}=cos\frac{\varphi}{2}+\vec{\mathbf{n}}sin\frac{\varphi}{2}$

Мы можем разложить вектор $\vec{\mathbf{u}}$
на компоненту, параллельную и перпендикулярную $\vec{\boldsymbol{\lambda}}$
(см. рисунок):

$\vec{\mathbf{u}}=\vec{\mathbf{u}}_\parallel+\vec{\mathbf{u}}_\bot$

При повороте, параллельная компонента останется на том же месте, тогда как перпендикулярная повернётся, но останется в плоскости p, перпендикулярной $\vec{\mathbf{a}}$
. После поворота, этот вектор можно будет выразить следующим образом:

$\vec{\mathbf{v}}_\bot=\vec{\mathbf{u}}_\bot\,cos\varphi+\vec{\mathbf{w}}sin\varphi$

Вектор $\vec{\mathbf{w}}$
перпендикулярен как $\vec{\boldsymbol{\lambda}}$
, так и $\vec{\mathbf{u}}_\bot$
, поэтому он должен быть сонаправлен их векторному произведению:
$\vec{\mathbf{w}}=k[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}_\bot]=k[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]$

(мы можем вписать в векторное произведение полный вектор $\vec{\mathbf{u}}$
, поскольку параллельная компонента даст ноль)

Его длина составит
$|\vec{\mathbf{w}}|=k|\vec{\boldsymbol{\lambda}}||\vec{\mathbf{u}}_\bot|=k|\vec{\mathbf{u}}_\bot|sin\frac{\varphi}{2}$

Мы хотим, чтобы $|\vec{\mathbf{w}}|=|\vec{\mathbf{u}}_\bot|$
, т.е вектор описывал бы круговое движение, сохраняющее его длину. Отсюда найдём:
$\vec{\mathbf{w}}=\frac{[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]}{sin\frac{\varphi}{2}}$

Также нам необходимо выразить $\vec{\mathbf{u}}_\bot$
через исходные данные. Наиболее простой способ - через ещё одно векторное произведение:
$\vec{\mathbf{u}}_\bot=k_1[\vec{\mathbf{w}},\vec{\boldsymbol{\lambda}}]$

Снова найдём длину:
$|\vec{\mathbf{u}}_\bot|=k_1|\vec{\mathbf{w}}||\vec{\boldsymbol{\lambda}}|=k_1|\vec{\mathbf{u}}_\bot|sin\frac{\varphi}{2}$

Отсюда,
$\vec{\mathbf{u}}_\bot=\frac{[\vec{\mathbf{w}},\vec{\boldsymbol{\lambda}}]}{sin\frac{\varphi}{2}}=-\frac{[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]}{sin^2\frac{\varphi}{2}}$

Наконец, выразим итоговый вектор:
$\vec{\mathbf{v}}=\vec{\mathbf{u}}_\parallel+\vec{\mathbf{u}}_\bot\,cos\varphi+\vec{\mathbf{w}}sin\varphi=$

$=\vec{\mathbf{u}}_\parallel+\vec{\mathbf{u}}_\bot(1-(1-cos\varphi))+\frac{[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]}{sin\frac{\varphi}{2}}sin\varphi$

$=\vec{\mathbf{u}}_\parallel+\vec{\mathbf{u}}_\bot-\vec{\mathbf{u}}_\bot(1-cos\varphi)+\frac{2sin\frac{\varphi}{2}cos\frac{\varphi}{2}[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]}{sin\frac{\varphi}{2}}=$

$=\vec{\mathbf{u}}+\frac{[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]}{sin^2\frac{\varphi}{2}}(1-cos\varphi)+2cos\frac{\varphi}{2}[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]=$

$=\vec{\mathbf{u}}+\frac{[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]}{sin^2\frac{\varphi}{2}}2sin^2\frac{\varphi}{2}+2\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]=$

$=\vec{\mathbf{u}}+2[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]+2\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]$

Вводим вектор $\vec{\mathbf{t}}$
:
$\vec{\mathbf{t}}=2[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]$

и получаем формулу (6.6).

Алгебраическое доказательство.
Ту же формулу можно вывести непосредственно из (6.1). Сначала выразим промежуточный кватернион:
$M=\vec{\mathbf{u}}\overline{\Lambda}=\vec{\mathbf{u}}(\lambda_0-\vec{\boldsymbol{\lambda}})=\lambda_0\vec{\mathbf{u}}+[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})$

Далее, найдём итоговый кватернион, который должен оказаться вектором:
$N=\Lambda\,M=(\lambda_0+\vec{\boldsymbol{\lambda}})(\lambda_0\vec{\mathbf{u}}+[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}))=$

$=\lambda_0^2\vec{\mathbf{u}}+\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+\cancel{\lambda_0(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})}+\vec{\boldsymbol{\lambda}}(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})-\cancel{\lambda_0(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})}-\cancel{(\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}])}+\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]$

Мы «выкинули» (говоря русским языком: похерили) $(\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}])$
, поскольку векторное произведение ортогонально каждому из своих множителей, т.е имеет нулевое скалярное произведение с ними. Можно также сказать, что это слагаемое - не что иное, как смешанное произведение трёх векторов и равно нулю, если они линейно зависимы.
После упрощения имеем:

$N=\lambda_0^2\vec{\mathbf{u}}+2\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+\vec{\boldsymbol{\lambda}}(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})+[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]$
(6.7)

Заметим, что скалярная часть равна нулю, то есть перед нами действительно вектор.
Запишем формулу для двойного векторного произведения «бац минус цаб»:

$[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]=\vec{\boldsymbol{\lambda}}(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})-|\vec{\boldsymbol{\lambda}}|^2\vec{\mathbf{u}}$
(6.8)

но не будем убирать из (6.7) двойное векторное произведение, вместо этого выразим:
$\vec{\boldsymbol{\lambda}}(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})=[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]+|\vec{\boldsymbol{\lambda}}|^2\vec{\mathbf{u}}$

Подставляем его в (6.7):
$\vec{\mathbf{v}}=N=\lambda_0^2\vec{\mathbf{u}}+2\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+2[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]+|\vec{\boldsymbol{\lambda}}|^2\vec{\mathbf{u}}=$

$=(\lambda_0^2+|\vec{\boldsymbol{\lambda}}|^2)\vec{\mathbf{u}}+2\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+2[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]=$

$=\vec{\mathbf{u}}+2\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+2[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]$

что снова совпадает с формулой (6.6).

Вариации на тему
Может показаться, что два векторных произведения - слишком много для поворота вектора, надо обойтись одним. Да, это возможно. Подставим (6.8) в (6.7):
$\vec{\mathbf{v}}=N=\lambda_0^2\vec{\mathbf{u}}+2\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+\vec{\boldsymbol{\lambda}}(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})+[\vec{\boldsymbol{\lambda}},[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]]=$

$=\lambda_0^2\vec{\mathbf{u}}+2\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+2\vec{\boldsymbol{\lambda}}(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})-|\vec{\boldsymbol{\lambda}}|^2\vec{\mathbf{u}}=$

$=(\lambda_0^2-|\vec{\boldsymbol{\lambda}}|^2)\vec{\mathbf{u}}+2\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+2\vec{\boldsymbol{\lambda}}(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}})$

Первый множитель можно упростить, учитывая, что кватернион имеет единичную норму, а также вынести двойку:

$\vec{\mathbf{v}}=(2\lambda_0^2-1)\vec{\mathbf{u}}+2(\lambda_0[\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}]+\vec{\boldsymbol{\lambda}}(\vec{\boldsymbol{\lambda}},\vec{\mathbf{u}}))$
(6.9)

Посчитаем количество операций в формуле (6.9).
Вычисление первого множителя требует 3m+1a. Затем 3m, чтобы помножить вектор на скаляр. Далее, 6m+3a для векторного произведения, 3m - умножение вектора на скаляр, 3m+2a - скалярное произведение, 3m - умножение вектора на скаляр, 3a - сложение векторов, 3a - их удвоение (либо 3m, не так уж важно в современных компьютерах), 3a - сложение векторов.
Итого: 21 умножение и 15 сложений - несколько хуже, чем по формуле (6.6).
Но как мы знаем, кроме количества арифметических операций, в современных процессорах очень важно, насколько эти операции «независимы» - требует ли каждая следующая результатов предыдущей, или можно начать выполнять их параллельно. В этом плане формула (6.9) и даже самые простые выкладки (6.3), (6.5) могут оказаться предпочтительнее.

Вычисление матрицы поворота
В компьютерной графике (и не только) часто встречается ситуация, когда необходимо повернуть множество векторов с помощью одного и того же кватерниона. В таком случае имеет смысл преобразовать кватернион в матрицу и производить поворот уже с её помощью.
Вывести матрицу поворота, соответствующую кватерниону Λ - несложно. Заметим, что все выведенные нами формулы линейны относительно входного вектора $\vec{\mathbf{u}}$
. Поэтому достаточно повернуть базисные векторы, после чего выразить поворот произвольного вектора как их взвешенную сумму.
Найдём, как повернётся вектор (1;0;0)T, он же i в кватернионной записи. Как ни странно, наиболее просто это делается «по определению». Здесь, когда используется всего один кватернион, для упрощения работы уберём индексы, записав его так:
$\Lambda=a+ix+jy+kz$

Тогда:
$\vec{\mathbf{v}}_1=\Lambda\,i\overline{\Lambda}=(a+ix+jy+kz)i(a-ix-jy-kz)=$

$=(a+ix+jy+kz)(ia+x-ky+jz)=$

$=ia^2+\cancel{ax}-{\color{Blue}kay}+{\color{Green}jaz}-$

$-\cancel{ax}+ix^2+{\color{Green}jxy}+{\color{Blue}kxz}-$

$-{\color{Blue}kay}+{\color{Green}jxy}-iy^2+\cancel{yz}+$

$+{\color{Green}jaz}+{\color{Blue}kxz}-\cancel{yz}-iz^2=$

$=(a^2+x^2-y^2-z^2)i+2(xy+az)j+2(xz-ay)k$

Первый множитель можно упростить, учитывая единичную норму:
$\vec{\mathbf{v}}_1=(a^2+x^2+y^2+z^2-2y^2-2z^2)i+2(xy+az)j+2(xz-ay)k=$

$=(1-2y^2-2z^2)i+2(xy+az)j+2(xz-ay)k$

Повторим то же самое для вектора (0;1;0)T, он же j:
$\vec{\mathbf{v}}_2=\Lambda\,j\overline{\Lambda}=(a+ix+jy+kz)j(a-ix-jy-kz)=$

$=(a+ix+jy+kz)(ja+kx+y-iz)=$

$={\color{Green}ja^2}+{\color{Blue}kax}+\cancel{ay}-iaz+$

$={\color{Blue}kax}-{\color{Green}jx^2}+ixy+\cancel{xz}-$

$-\cancel{ay}+ixy+{\color{Green}jy^2}+{\color{Blue}kyz}-$

$-iaz-\cancel{xz}+{\color{Blue}kyz}-{\color{Green}jz^2}=$

$=2(xy-az)i+(a^2-x^2+y^2-z^2)j+2(ax+yz)k=$

$=2(xy-az)i+(1-2x^2-2z^2)j+2(ax+yz)k$

И, наконец, для вектора (0;0;1)T, он же k:
$\vec{\mathbf{v}}_3=\Lambda\,k\overline{\Lambda}=(a+ix+jy+kz)k(a-ix-jy-kz)=$

$=(a+ix+jy+kz)(ka-jx+iy+z)=$

$={\color{Green}ka^2}-{\color{Blue}jax}+iay+\cancel{az}-$

$-{\color{Blue}jax}-{\color{Green}kx^2}-\cancel{xy}+ixz+$

$+iay+\cancel{xy}-{\color{Green}ky^2}+{\color{Blue}jyz}-$

$-\cancel{az}+ixz+{\color{Blue}jyz}+{\color{Green}kz^2}=$

$=2(xz+ay)i+2(yz-ax)j+(a^2+z^2-x^2-y^2)k=$

$=2(xz+ay)i+2(yz-ax)j+(1-2x^2-2y^2)k=$

Наконец, из векторов $\vec{\mathbf{v}}_1,\vec{\mathbf{v}}_2,\vec{\mathbf{v}}_3$
можно составить матрицу:

$A=\left(\begin{array}{ccc} 1-2y^2-2z^2&2xy-2az&2xz+2ay\\2xy+2az&1-2x^2-2z^2&2yz-2ax\\2xz-2ay&2yz+2ax&1-2x^2-2y^2\end{array}\right)$

На её подсчет требуется не более 13 умножений и 12 сложений. Такое количество операций достигается, если мы первоначально умножим наш кватернион на $\sqrt2$
(потребует 4 умножения), благодаря чему все двойки из выражения уйдут. Дальше мы вычислим произведения x2, y2, z2, ax, ay, az, xy, xz, yz (ещё 9 умножений), и останется лишь сложить всё вместе.
Возможно, есть и более быстрые методы, но, как правило, вычисление матрицы занимает ничтожно мало времени по сравнению с последующим поворотом сотен векторов, поэтому вовсе не обязательно оптимизировать эту операцию.

Пока что кватернионы проигрывают в простоте и эффективности матрицам поворота (хотя и не значительно). Их полюбили в космической технике явно за другое… Продолжение следует.

кватернионы-это просто (том 1), математика, работа