multiplication of three quaternions,version 1.0
peter.vlasschaert@gmail.com,24/10/2018

load ' packages '
==============

kill(all)$
load(vect)$
load(lrats)$

p1:a1+a2*i+a3*j+a4*k;
p2:b1+b2*i+b3*j+b4*k;
p3:c1+c2*i+c3*j+c4*k;

p:[a1,a2,a3,a4,b1,b2,b3,b4,c1,c2,c3,c4];

p4:declare([a1,a2,a3,a4,b1,b2,b3,b4,c1,c2,c3,c4],scalar);

p5:[i.j=k,j.i=-k,j.k=i,k.j=-i,k.i=j,i.k=-j,i.i=-1,j.j=-1,k.k=-1];

multiplication : two quaternions.
===========================

=>   p1.p2

p6:p1.p2;

p7:expand(p6);

p8:subst(p5,p7);

multiplication :three quaternions.
============================
=>  p1.p2.p3
------------------------------------------------------------------
How to :  q1 = p1.p2  ,  use dot
                  q2  = q1.p3 ,  use dot
=====================================
rem : p5 rules for quaternions from hamilton.

q1:subst(p5,expand(p1.p2));

 check :  'p8' and 'q1'

p9:is(p8=q1);

q2:subst(p5,expand(q1.p3));

result : q3 = d1+d2*i+d3*j+d4*k
==========================

d1 => ?   RE  : real-part  'q3'

p10:[i=0,j=0,k=0];
p10a:d1=subst(p10,q2);

d2 => ?  i: vector-part  'q3'

p10b:d2=coeff(q2,i);

d3 => ?  j: vector-part  'q3'

p10c:d3=coeff(q2,j);

d4 => ?  k: vector-part  'q3'

p10d:d4=coeff(q2,k);

definition :multiplication of two quaternions  'qmult2'
============================================
qa = a1+a2*i+a3*j+a4*k
qb =b1+b1*i+b2*j+b3*k
============================================
=>  qa.qb ,         use qmult2 you need 'p4' : see above.

qmult2(a1,a2,a3,a4,b1,b2,b3,b4,h):=(qa:a1+a2*i+a3*j+a4*k,qb:b1+b1*i+b2*j+b3*k,h1:expand(qa.qb),
                                                                h2:[i.j=k,j.i=-k,j.k=i,k.j=-i,k.i=j,i.k=-j,i.i=-1,j.j=-1,k.k=-1],h:subst(h2,h1))$

qc:qmult2(1,1,1,1,1,1,1,1,h);

see : also simplify quaternion multiplication,version 1.0,website. 'p9'

How to multiply : three quaternions. ( use 'qmult2,coefq' )
================================================

thrid quaternion : qd: 1+ 1*i+2*j+3*k

result :  qc.qd

1e) Problem : find coef. from i,j,k and real-part  'qc'
=======================================

definition : coef. from quaternion 
============================

coefq(qq,ss):=(s1:subst([i=0,j=0,k=0],qq),s2:coeff(qq,i),s3:coeff(qq,j),s4:coeff(qq,k),ss:[s1,s2,s3,s4])$

m1:coefq(qc,ss);

m11:m1[1];

2e) thrid quaternion, coef. (use 'coefq')
================================

qd: 1+ 2*i+3*j+4*k;
m2:coefq(qd,ss);

3e) multiplication of three quaternions (use 'qmult2')
============================================

qe:qmult2(m1[1],m1[2],m1[3],m1[4] ,m2[1] ,m2[2],m2[3],m2[4] ,h);                          

4e) coef. from the result after the three multiplications
=============================================

m3:coefq(qe,ss);

rem : (1+1*i+1*j+1*k).(1+1*i+1*j+1*k).(1+2*i+3*j+4*k) = -14+2*i-6*j-2*k

easy : understanding
==================

general :rem : build two functions : qmult2,coefq .
general : see also website :simplify quaternion multiplication,version 1.0

-------------------------------------------------------------------------------------------------------