Jump to content

Quaternion type

From Rosetta Code
Quaternion type is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.

Quaternions are an extension of the idea of complex numbers

A complex number has a real and complex part written sometimes as a + bi, where a and b stand for real numbers and i stands for the square root of minus 1. An example of a complex number might be -3 + 2i, where the real part, a is -3.0 and the complex part, b is +2.0.

A quaternion has one real part and three imaginary parts, i, j, and k. A quaternion might be written as a + bi + cj +dk. In this numbering system, ii = jj = kk = ijk = -1. The order of multiplication is important, as, in general, for two quaternions q1 and q2; q1q2 != q2q1. An example of a quaternion might be 1 +2i +3j +4k

There is a list form of notation where just the numbers are shown and the imaginary multipliers i,j, and k are assumed by position. So the example above would be written as (1, 2, 3, 4)


Task Description
Given the three quaternions and their components:

   q  = (1, 2, 3, 4) = (a,  b,  c,  d )
   q1 = (2, 3, 4, 5) = (a1, b1, c1, d1)
   q2 = (3, 4, 5, 6) = (a2, b2, c2, d2)

And a wholly real number r = 7.

Your task is to create functions or classes to perform simple maths with quaternions including computing:

  1. The norm of a quaternion:
  2. The negative of a quaternion:
    =(-a, -b, -c, -d)
  3. The conjugate of a quaternion:
    =( a, -b, -c, -d)
  4. Addition of a real number r and a quaternion q:
    r + q = q + r = (a+r, b, c, d)
  5. Addition of two quaternions:
    q1 + q2 = (a1+a2, b1+b2, c1+c2, d1+d2)
  6. Multiplication of a real number and a quaternion:
    qr = rq = (ar, br, cr, dr)
  7. Multiplication of two quaternions q1 and q2 is given by:
    ( a1a2 − b1b2 − c1c2 − d1d2,
      a1b2 + b1a2 + c1d2 − d1c2,
      a1c2 − b1d2 + c1a2 + d1b2,
      a1d2 + b1c2 − c1b2 + d1a2 )
  8. Show that, for the two quaternions q1 and q2:
    q1q2 != q2q1

If your language has built-in support for quaternions then use it.

python

This example extends Pythons namedtuples to add extra functionality. <lang python>from collections import namedtuple import math

class Q(namedtuple('Quaternion', 'real, i, j, k')):

   'Quaternion type: Q(real=0.0, i=0.0, j=0.0, k=0.0)' 
   __slots__ = () 
   def __new__(_cls, real=0.0, i=0.0, j=0.0, k=0.0):
       'Defaults all parts of quaternion to zero'
       return super().__new__(_cls, float(real), float(i), float(j), float(k))
   def conjugate(self):
       return Q(self.real, -self.i, -self.j, -self.k)
   def _norm2(self):
       return sum( x*x for x in self)
   def norm(self):
       return math.sqrt(self._norm2())
   def reciprocal(self):
       n2 = self._norm2()
       return Q(*(x / n2 for x in self.conjugate())) 
   def __str__(self):
       'Shorter form of Quaternion as string'
       return 'Q(%g, %g, %g, %g)' % self
   def __neg__(self):
       return Q(-self.real, -self.i, -self.j, -self.k)
   def __add__(self, other):
       if type(other) == Q:
           return Q( *(s+o for s,o in zip(self, other)) )
       try:
           f = float(other)
       except:
           return NotImplemented
       return Q(self.real + f, self.i, self.j, self.k)
   def __radd__(self, other):
       return Q.__add__(self, other)
   def __mul__(self, other):
       if type(other) == Q:
           a1,b1,c1,d1 = self
           a2,b2,c2,d2 = other
           return Q(
                a1*a2 - b1*b2 - c1*c2 - d1*d2,
                a1*b2 + b1*a2 + c1*d2 - d1*c2,
                a1*c2 - b1*d2 + c1*a2 + d1*b2,
                a1*d2 + b1*c2 - c1*b2 + d1*a2 )
       try:
           f = float(other)
       except:
           return NotImplemented
       return Q(self.real * f, self.i * f, self.j * f, self.k * f)
   def __rmul__(self, other):
       return Q.__mul__(self, other)
   def __truediv__(self, other):
       if type(other) == Q:
           return self.__mul__(other.reciprocal())
       try:
           f = float(other)
       except:
           return NotImplemented
       return Q(self.real / f, self.i / f, self.j / f, self.k / f)
   def __rtruediv__(self, other):
       return other * self.reciprocal()
   __div__, __rdiv__ = __truediv__, __rtruediv__

Quaternion = Q

q = Q(1, 2, 3, 4) q1 = Q(2, 3, 4, 5) q2 = Q(3, 4, 5, 6) r = 7</lang>

Continued shell session Run the above with the -i flag to python on the command line, or run with idle then continue in the shell as follows: <lang python>>>> q Quaternion(real=1.0, i=2.0, j=3.0, k=4.0) >>> q1 Quaternion(real=2.0, i=3.0, j=4.0, k=5.0) >>> q2 Quaternion(real=3.0, i=4.0, j=5.0, k=6.0) >>> r 7 >>> q.norm() 5.477225575051661 >>> q1.norm() 7.3484692283495345 >>> q2.norm() 9.273618495495704 >>> -q Quaternion(real=-1.0, i=-2.0, j=-3.0, k=-4.0) >>> q.conjugate() Quaternion(real=1.0, i=-2.0, j=-3.0, k=-4.0) >>> r + q Quaternion(real=8.0, i=2.0, j=3.0, k=4.0) >>> q + r Quaternion(real=8.0, i=2.0, j=3.0, k=4.0) >>> q1 + q2 Quaternion(real=5.0, i=7.0, j=9.0, k=11.0) >>> q2 + q1 Quaternion(real=5.0, i=7.0, j=9.0, k=11.0) >>> q * r Quaternion(real=7.0, i=14.0, j=21.0, k=28.0) >>> r * q Quaternion(real=7.0, i=14.0, j=21.0, k=28.0) >>> q1 * q2 Quaternion(real=-56.0, i=16.0, j=24.0, k=26.0) >>> q2 * q1 Quaternion(real=-56.0, i=18.0, j=20.0, k=28.0) >>> assert q1 * q2 != q2 * q1 >>> >>> i, j, k = Q(0,1,0,0), Q(0,0,1,0), Q(0,0,0,1) >>> i*i Quaternion(real=-1.0, i=0.0, j=0.0, k=0.0) >>> j*j Quaternion(real=-1.0, i=0.0, j=0.0, k=0.0) >>> k*k Quaternion(real=-1.0, i=0.0, j=0.0, k=0.0) >>> i*j*k Quaternion(real=-1.0, i=0.0, j=0.0, k=0.0) >>> q1 / q2 Quaternion(real=0.7906976744186047, i=0.023255813953488358, j=-2.7755575615628914e-17, k=0.046511627906976744) >>> q1 / q2 * q2 Quaternion(real=2.0000000000000004, i=3.0000000000000004, j=4.000000000000001, k=5.000000000000001) >>> q2 * q1 / q2 Quaternion(real=2.0, i=3.465116279069768, j=3.906976744186047, k=4.767441860465116) >>> q1.reciprocal() * q1 Quaternion(real=0.9999999999999999, i=0.0, j=0.0, k=0.0) >>> q1 * q1.reciprocal() Quaternion(real=0.9999999999999999, i=0.0, j=0.0, k=0.0) >>> </lang>

J

Derived from the j wiki:

<lang j> NB. utilities

  ip=: +/ .*
  t4=. (_1^#:0 10 9 12)*0 7 16 23 A.=i.4
  toQ=: 4&{."1
  NB. task
  norm=: %:@ip~          NB. | y
  neg=:  -               NB. - y
  conj=: 1 _1 _1 _1 * ]  NB. + y
  add=: +&toQ            NB. x + y
  mul=: (ip t4 ip ])&toQ NB. x * y</lang>

Example use:

<lang> q=: 1 2 3 4

  q1=: 2 3 4 5
  q2=: 3 4 5 6
  r=: 7
  
  norm q

5.47723

  neg q

_1 _2 _3 _4

  conj q

1 _2 _3 _4

  r add q

8 2 3 4

  q1 add q2

5 7 9 11

  r mul q

7 14 21 28

  q1 mul q2

_56 16 24 26

  q2 mul q1

_56 18 20 28</lang>

python

This example extends Pythons namedtuples to add extra functionality. <lang python>from collections import namedtuple import math

class Q(namedtuple('Quaternion', 'real, i, j, k')):

   'Quaternion type: Q(real=0.0, i=0.0, j=0.0, k=0.0)' 
   __slots__ = () 
   def __new__(_cls, real=0.0, i=0.0, j=0.0, k=0.0):
       'Defaults all parts of quaternion to zero'
       return super().__new__(_cls, float(real), float(i), float(j), float(k))
   def conjugate(self):
       return Q(self.real, -self.i, -self.j, -self.k)
   def _norm2(self):
       return sum( x*x for x in self)
   def norm(self):
       return math.sqrt(self._norm2())
   def reciprocal(self):
       n2 = self._norm2()
       return Q(*(x / n2 for x in self.conjugate())) 
   def __str__(self):
       'Shorter form of Quaternion as string'
       return 'Q(%g, %g, %g, %g)' % self
   def __neg__(self):
       return Q(-self.real, -self.i, -self.j, -self.k)
   def __add__(self, other):
       if type(other) == Q:
           return Q( *(s+o for s,o in zip(self, other)) )
       try:
           f = float(other)
       except:
           return NotImplemented
       return Q(self.real + f, self.i, self.j, self.k)
   def __radd__(self, other):
       return Q.__add__(self, other)
   def __mul__(self, other):
       if type(other) == Q:
           a1,b1,c1,d1 = self
           a2,b2,c2,d2 = other
           return Q(
                a1*a2 - b1*b2 - c1*c2 - d1*d2,
                a1*b2 + b1*a2 + c1*d2 - d1*c2,
                a1*c2 - b1*d2 + c1*a2 + d1*b2,
                a1*d2 + b1*c2 - c1*b2 + d1*a2 )
       try:
           f = float(other)
       except:
           return NotImplemented
       return Q(self.real * f, self.i * f, self.j * f, self.k * f)
   def __rmul__(self, other):
       return Q.__mul__(self, other)
   def __truediv__(self, other):
       if type(other) == Q:
           return self.__mul__(other.reciprocal())
       try:
           f = float(other)
       except:
           return NotImplemented
       return Q(self.real / f, self.i / f, self.j / f, self.k / f)
   def __rtruediv__(self, other):
       return other * self.reciprocal()
   __div__, __rdiv__ = __truediv__, __rtruediv__

Quaternion = Q

q = Q(1, 2, 3, 4) q1 = Q(2, 3, 4, 5) q2 = Q(3, 4, 5, 6) r = 7</lang>

Continued shell session Run the above with the -i flag to python on the command line, or run with idle then continue in the shell as follows: <lang python>>>> q Quaternion(real=1.0, i=2.0, j=3.0, k=4.0) >>> q1 Quaternion(real=2.0, i=3.0, j=4.0, k=5.0) >>> q2 Quaternion(real=3.0, i=4.0, j=5.0, k=6.0) >>> r 7 >>> q.norm() 5.477225575051661 >>> q1.norm() 7.3484692283495345 >>> q2.norm() 9.273618495495704 >>> -q Quaternion(real=-1.0, i=-2.0, j=-3.0, k=-4.0) >>> q.conjugate() Quaternion(real=1.0, i=-2.0, j=-3.0, k=-4.0) >>> r + q Quaternion(real=8.0, i=2.0, j=3.0, k=4.0) >>> q + r Quaternion(real=8.0, i=2.0, j=3.0, k=4.0) >>> q1 + q2 Quaternion(real=5.0, i=7.0, j=9.0, k=11.0) >>> q2 + q1 Quaternion(real=5.0, i=7.0, j=9.0, k=11.0) >>> q * r Quaternion(real=7.0, i=14.0, j=21.0, k=28.0) >>> r * q Quaternion(real=7.0, i=14.0, j=21.0, k=28.0) >>> q1 * q2 Quaternion(real=-56.0, i=16.0, j=24.0, k=26.0) >>> q2 * q1 Quaternion(real=-56.0, i=18.0, j=20.0, k=28.0) >>> assert q1 * q2 != q2 * q1 >>> >>> i, j, k = Q(0,1,0,0), Q(0,0,1,0), Q(0,0,0,1) >>> i*i Quaternion(real=-1.0, i=0.0, j=0.0, k=0.0) >>> j*j Quaternion(real=-1.0, i=0.0, j=0.0, k=0.0) >>> k*k Quaternion(real=-1.0, i=0.0, j=0.0, k=0.0) >>> i*j*k Quaternion(real=-1.0, i=0.0, j=0.0, k=0.0) >>> q1 / q2 Quaternion(real=0.7906976744186047, i=0.023255813953488358, j=-2.7755575615628914e-17, k=0.046511627906976744) >>> q1 / q2 * q2 Quaternion(real=2.0000000000000004, i=3.0000000000000004, j=4.000000000000001, k=5.000000000000001) >>> q2 * q1 / q2 Quaternion(real=2.0, i=3.465116279069768, j=3.906976744186047, k=4.767441860465116) >>> q1.reciprocal() * q1 Quaternion(real=0.9999999999999999, i=0.0, j=0.0, k=0.0) >>> q1 * q1.reciprocal() Quaternion(real=0.9999999999999999, i=0.0, j=0.0, k=0.0) >>> </lang>

Tcl

<lang tcl># Support class that provides C++-like RAII lifetimes oo::class create RAII-support {

   constructor {} {

upvar 1 { end } end lappend end [self] trace add variable end unset [namespace code {my destroy}]

   }
   destructor {

catch { upvar 1 { end } end trace remove variable end unset [namespace code {my destroy}] }

   }
   method return Template:Level 1 {

incr level upvar 1 { end } end upvar $level { end } parent trace remove variable end unset [namespace code {my destroy}] lappend parent [self] trace add variable parent unset [namespace code {my destroy}] return -level $level [self]

   }

}

  1. Class of quaternions

oo::class create Q {

   superclass RAII-support
   variable R I J K
   constructor {{real 0} {i 0} {j 0} {k 0}} {

next namespace import ::tcl::mathfunc::* ::tcl::mathop::* variable R [double $real] I [double $i] J [double $j] K [double $k]

   }
   self method return args {

[my new {*}$args] return 2

   }
   method p {} {

return "Q($R,$I,$J,$K)"

   }
   method values {} {

list $R $I $J $K

   }
   method Norm {} {

+ [* $R $R] [* $I $I] [* $J $J] [* $K $K]

   }
   method conjugate {} {

Q return $R [- $I] [- $J] [- $K]

   }
   method norm {} {

sqrt [my Norm]

   }
   method unit {} {

set n [my norm] Q return [/ $R $n] [/ $I $n] [/ $J $n] [/ $K $n]

   }
   method reciprocal {} {

set n2 [my Norm] Q return [/ $R $n2] [/ $I $n2] [/ $J $n2] [/ $K $n2]

   }
   method - Template:Q "" {

if {[llength [info level 0]] == 2} { Q return [- $R] [- $I] [- $J] [- $K] } [my + [$q -]] return

   }
   method + q {

if {[info object isa object $q]} { lassign [$q values] real i j k Q return [+ $R $real] [+ $I $i] [+ $J $j] [+ $K $k] } Q return [+ $R [double $q]] $I $J $K

   }
   method * q {

if {[info object isa object $q]} { lassign [my values] a1 b1 c1 d1 lassign [$q values] a2 b2 c2 d2 Q return [expr {$a1*$a2 - $b1*$b2 - $c1*$c2 - $d1*$d2}] \ [expr {$a1*$b2 + $b1*$a2 + $c1*$d2 - $d1*$c2}] \ [expr {$a1*$c2 - $b1*$d2 + $c1*$a2 + $d1*$b2}] \ [expr {$a1*$d2 + $b1*$c2 - $c1*$b2 + $d1*$a2}] } set f [double $q] Q return [* $R $f] [* $I $f] [* $J $f] [* $K $f]

   }
   method == q {

expr { [info object isa object $q] && [info object isa typeof $q [self class]] && [my values] eq [$q values] }

   }
   export - + * ==

}</lang> Demonstration code: <lang tcl>set q [Q new 1 2 3 4] set q1 [Q new 2 3 4 5] set q2 [Q new 3 4 5 6] set r 7

puts "q = [$q p]" puts "q1 = [$q1 p]" puts "q2 = [$q2 p]" puts "r = $r" puts "q norm = [$q norm]" puts "q1 norm = [$q1 norm]" puts "q2 norm = [$q2 norm]" puts "-q = [[$q -] p]" puts "q conj = [[$q conjugate] p]" puts "q + r = [[$q + $r] p]"

  1. Real numbers are not objects, so no extending operations for them

puts "q1 + q2 = [[$q1 + $q2] p]" puts "q2 + q1 = [[$q2 + $q1] p]" puts "q * r = [[$q * $r] p]" puts "q1 * q2 = [[$q1 * $q2] p]" puts "q2 * q1 = [[$q2 * $q1] p]" puts "equal(q1*q2, q2*q1) = [[$q1 * $q2] == [$q2 * $q1]]"</lang> Output:

q = Q(1.0,2.0,3.0,4.0)
q1 = Q(2.0,3.0,4.0,5.0)
q2 = Q(3.0,4.0,5.0,6.0)
r = 7
q norm = 5.477225575051661
q1 norm = 7.3484692283495345
q2 norm = 9.273618495495704
-q = Q(-1.0,-2.0,-3.0,-4.0)
q conj = Q(1.0,-2.0,-3.0,-4.0)
q + r = Q(8.0,2.0,3.0,4.0)
q1 + q2 = Q(5.0,7.0,9.0,11.0)
q2 + q1 = Q(5.0,7.0,9.0,11.0)
q * r = Q(7.0,14.0,21.0,28.0)
q1 * q2 = Q(-56.0,16.0,24.0,26.0)
q2 * q1 = Q(-56.0,18.0,20.0,28.0)
equal(q1*q2, q2*q1) = 0
Cookies help us deliver our services. By using our services, you agree to our use of cookies.