    -- Copyright (c) 2010-2015, David Amos. All rights reserved.

    {-# LANGUAGE MultiParamTypeClasses, FlexibleInstances, TypeSynonymInstances, NoMonomorphismRestriction #-}

    -- |A module defining the algebra of quaternions over an arbitrary field.
    --
    -- The quaternions are the algebra defined by the basis {1,i,j,k}, where i^2 = j^2 = k^2 = ijk = -1
    module Math.Algebras.Quaternions where

    import Prelude hiding ( (*>) )

    import Math.Core.Field
    import Math.Algebras.VectorSpace
    import Math.Algebras.TensorProduct
    import Math.Algebras.Structures

    -- Conway & Smith, On Quaternions and Octonions

    -- QUATERNIONS

    data HBasis = One | I | J | K deriving (HBasis -> HBasis -> Bool
    (HBasis -> HBasis -> Bool)
    -> (HBasis -> HBasis -> Bool) -> Eq HBasis
    forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
    /= :: HBasis -> HBasis -> Bool
    $c/= :: HBasis -> HBasis -> Bool
    == :: HBasis -> HBasis -> Bool
    $c== :: HBasis -> HBasis -> Bool
    Eq,Eq HBasis
    Eq HBasis =>
    (HBasis -> HBasis -> Ordering)
    -> (HBasis -> HBasis -> Bool)
    -> (HBasis -> HBasis -> Bool)
    -> (HBasis -> HBasis -> Bool)
    -> (HBasis -> HBasis -> Bool)
    -> (HBasis -> HBasis -> HBasis)
    -> (HBasis -> HBasis -> HBasis)
    -> Ord HBasis
    HBasis -> HBasis -> Bool
    HBasis -> HBasis -> Ordering
    HBasis -> HBasis -> HBasis
    forall a.
    Eq a =>
    (a -> a -> Ordering)
    -> (a -> a -> Bool)
    -> (a -> a -> Bool)
    -> (a -> a -> Bool)
    -> (a -> a -> Bool)
    -> (a -> a -> a)
    -> (a -> a -> a)
    -> Ord a
    min :: HBasis -> HBasis -> HBasis
    $cmin :: HBasis -> HBasis -> HBasis
    max :: HBasis -> HBasis -> HBasis
    $cmax :: HBasis -> HBasis -> HBasis
    >= :: HBasis -> HBasis -> Bool
    $c>= :: HBasis -> HBasis -> Bool
    > :: HBasis -> HBasis -> Bool
    $c> :: HBasis -> HBasis -> Bool
    <= :: HBasis -> HBasis -> Bool
    $c<= :: HBasis -> HBasis -> Bool
    < :: HBasis -> HBasis -> Bool
    $c< :: HBasis -> HBasis -> Bool
    compare :: HBasis -> HBasis -> Ordering
    $ccompare :: HBasis -> HBasis -> Ordering
    $cp1Ord :: Eq HBasis
    Ord)

    type Quaternion k = Vect k HBasis

    instance Show HBasis where
        show :: HBasis -> String
    show One = "1"
        show I = "i"
        show J = "j"
        show K = "k"

    instance (Eq k, Num k) => Algebra k HBasis where
        unit :: k -> Vect k HBasis
    unit x :: k
    x = k
    x k -> Vect k HBasis -> Vect k HBasis
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    One
        mult :: Vect k (Tensor HBasis HBasis) -> Vect k HBasis
    mult = (Tensor HBasis HBasis -> Vect k HBasis)
    -> Vect k (Tensor HBasis HBasis) -> Vect k HBasis
    forall k b a.
    (Eq k, Num k, Ord b) =>
    (a -> Vect k b) -> Vect k a -> Vect k b
    linear Tensor HBasis HBasis -> Vect k HBasis
    forall k. (Num k, Eq k) => Tensor HBasis HBasis -> Vect k HBasis
    mult'
             where mult' :: Tensor HBasis HBasis -> Vect k HBasis
    mult' (One,b :: HBasis
    b) = HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    b
                   mult' (b :: HBasis
    b,One) = HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    b
                   mult' (I,I) = k -> Vect k HBasis
    forall k b. Algebra k b => k -> Vect k b
    unit (-1)
                   mult' (J,J) = k -> Vect k HBasis
    forall k b. Algebra k b => k -> Vect k b
    unit (-1)
                   mult' (K,K) = k -> Vect k HBasis
    forall k b. Algebra k b => k -> Vect k b
    unit (-1)
                   mult' (I,J) = HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    K
                   mult' (J,I) = -1 k -> Vect k HBasis -> Vect k HBasis
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    K
                   mult' (J,K) = HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    I
                   mult' (K,J) = -1 k -> Vect k HBasis -> Vect k HBasis
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    I
                   mult' (K,I) = HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    J
                   mult' (I,K) = -1 k -> Vect k HBasis -> Vect k HBasis
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    J

    -- |The quaternions have {1,i,j,k} as basis, where i^2 = j^2 = k^2 = ijk = -1.
    i,j,k :: Num k => Quaternion k
    i :: Quaternion k
    i = HBasis -> Quaternion k
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    I
    j :: Quaternion k
    j = HBasis -> Quaternion k
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    J
    k :: Quaternion k
    k = HBasis -> Quaternion k
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    K

    class Algebra k a => HasConjugation k a where
        -- |A conjugation operation is required to satisfy the following laws:
        --
        -- * conj (x+y) = conj x + conj y
        --
        -- * conj (x*y) = conj y * conj x  (note the order-reversal)
        --
        -- * conj (conj x) = x
        --
        -- * conj x = x if and only if x in k
        conj :: Vect k a -> Vect k a
        -- |The squared norm is defined as sqnorm x = x * conj x. It satisfies:
        --
        -- * sqnorm (x*y) = sqnorm x * sqnorm y
        --
        -- * sqnorm (unit k) = k^2, for k a scalar
        sqnorm :: Vect k a -> k

    -- |If an algebra has a conjugation operation, then it has multiplicative inverses,
    -- via 1\/x = conj x \/ sqnorm x
    instance (Eq k, Fractional k, Ord a, Show a, HasConjugation k a) => Fractional (Vect k a) where
        recip :: Vect k a -> Vect k a
    recip 0 = String -> Vect k a
    forall a. HasCallStack => String -> a
    error "recip 0"
        recip x :: Vect k a
    x = (1 k -> k -> k
    forall a. Fractional a => a -> a -> a
    / Vect k a -> k
    forall k a. HasConjugation k a => Vect k a -> k
    sqnorm Vect k a
    x) k -> Vect k a -> Vect k a
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> Vect k a -> Vect k a
    forall k a. HasConjugation k a => Vect k a -> Vect k a
    conj Vect k a
    x
        fromRational :: Rational -> Vect k a
    fromRational q :: Rational
    q = Rational -> k
    forall a. Fractional a => Rational -> a
    fromRational Rational
    q k -> Vect k a -> Vect k a
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> 1

    -- |The scalar part of the quaternion w+xi+yj+zk is w. Also called the real part.
    scalarPart :: (Num k) => Quaternion k -> k
    scalarPart :: Quaternion k -> k
    scalarPart = HBasis -> Quaternion k -> k
    forall k b. (Num k, Eq b) => b -> Vect k b -> k
    coeff HBasis
    One

    -- |The vector part of the quaternion w+xi+yj+zk is xi+yj+zk. Also called the pure part.
    vectorPart :: (Eq k, Num k) => Quaternion k -> Quaternion k
    vectorPart :: Quaternion k -> Quaternion k
    vectorPart q :: Quaternion k
    q = Quaternion k
    q Quaternion k -> Quaternion k -> Quaternion k
    forall a. Num a => a -> a -> a
    - Quaternion k -> k
    forall k. Num k => Quaternion k -> k
    scalarPart Quaternion k
    q k -> Quaternion k -> Quaternion k
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> 1

    instance (Eq k, Num k) => HasConjugation k HBasis where
        conj :: Vect k HBasis -> Vect k HBasis
    conj = (Vect k HBasis -> (HBasis -> Vect k HBasis) -> Vect k HBasis
    forall (m :: * -> *) a b. Monad m => m a -> (a -> m b) -> m b
    >>= HBasis -> Vect k HBasis
    forall k. (Num k, Eq k) => HBasis -> Vect k HBasis
    conj') where
            conj' :: HBasis -> Vect k HBasis
    conj' One = HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    One
            conj' imag :: HBasis
    imag = -1 k -> Vect k HBasis -> Vect k HBasis
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> HBasis -> Vect k HBasis
    forall (m :: * -> *) a. Monad m => a -> m a
    return HBasis
    imag
        -- ie conj = linear conj', but avoiding unnecessary nf call
        sqnorm :: Vect k HBasis -> k
    sqnorm x :: Vect k HBasis
    x = [k] -> k
    forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
    sum ([k] -> k) -> [k] -> k
    forall a b. (a -> b) -> a -> b
    $ ((HBasis, k) -> k) -> [(HBasis, k)] -> [k]
    forall a b. (a -> b) -> [a] -> [b]
    map ((k -> Integer -> k
    forall a b. (Num a, Integral b) => a -> b -> a
    ^2) (k -> k) -> ((HBasis, k) -> k) -> (HBasis, k) -> k
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (HBasis, k) -> k
    forall a b. (a, b) -> b
    snd) ([(HBasis, k)] -> [k]) -> [(HBasis, k)] -> [k]
    forall a b. (a -> b) -> a -> b
    $ Vect k HBasis -> [(HBasis, k)]
    forall k b. Vect k b -> [(b, k)]
    terms Vect k HBasis
    x
        -- sqnorm x = scalarPart (conj x * x) -- the vector part will be zero anyway
        -- sqnorm x = x <.> x
    {-
    instance Fractional k => Fractional (Quaternion k) where
        recip 0 = error "Quaternion.recip 0"
        recip x = (1 / sqnorm x) *> conj x
        fromRational q = fromRational q *> 1
    -}

    x :: Vect k HBasis
    x <.> :: Vect k HBasis -> Vect k HBasis -> k
    <.> y :: Vect k HBasis
    y = Vect k HBasis -> k
    forall k. Num k => Quaternion k -> k
    scalarPart (Vect k HBasis -> Vect k HBasis
    forall k a. HasConjugation k a => Vect k a -> Vect k a
    conj Vect k HBasis
    x Vect k HBasis -> Vect k HBasis -> Vect k HBasis
    forall a. Num a => a -> a -> a
    * Vect k HBasis
    y)
    -- x <..> y = 1/2 * (sqnorm x + sqnorm y - sqnorm (x-y))

    x :: a
    x^- :: a -> a -> a
    ^-1 = a -> a
    forall a. Fractional a => a -> a
    recip a
    x

    -- Conway p40
    refl :: Vect k a -> Vect k a -> Vect k a
    refl q :: Vect k a
    q = \x :: Vect k a
    x -> -Vect k a
    q Vect k a -> Vect k a -> Vect k a
    forall a. Num a => a -> a -> a
    * Vect k a -> Vect k a
    forall k a. HasConjugation k a => Vect k a -> Vect k a
    conj Vect k a
    x Vect k a -> Vect k a -> Vect k a
    forall a. Num a => a -> a -> a
    * Vect k a
    q

    -- Given a linear function f on the quaternions, return the matrix representing it,
    -- relative to a given basis. The matrix is considered as acting on the right.
    asMatrix :: (Vect a HBasis -> Vect a HBasis) -> [Vect a HBasis] -> [[a]]
    asMatrix f :: Vect a HBasis -> Vect a HBasis
    f bs :: [Vect a HBasis]
    bs = [ let fi :: Vect a HBasis
    fi = Vect a HBasis -> Vect a HBasis
    f Vect a HBasis
    ei in [Vect a HBasis
    ej Vect a HBasis -> Vect a HBasis -> a
    forall k. (Num k, Eq k) => Vect k HBasis -> Vect k HBasis -> k
    <.> Vect a HBasis
    fi | Vect a HBasis
    ej <- [Vect a HBasis]
    bs] | Vect a HBasis
    ei <- [Vect a HBasis]
    bs  ]
    -- It is possible to write this function using coeff, instead of <.>,
    -- but then you have to pass in I,J,K, instead of i,j,k, which is uglier.

    -- Conway p24
    -- A homomorphism from H\0 to SO3
    -- if q is restricted to unit quaternions, this is a double cover of SO3 (since q, -q induce same rotation)
    -- The unit quaternions form the group Spin3
    reprSO3' :: a -> a -> a
    reprSO3' q :: a
    q = \x :: a
    x -> a
    qa -> Integer -> a
    forall a a. (Eq a, Fractional a, Num a) => a -> a -> a
    ^-1 a -> a -> a
    forall a. Num a => a -> a -> a
    * a
    x a -> a -> a
    forall a. Num a => a -> a -> a
    * a
    q

    -- |Given a non-zero quaternion q in H, the map x -> q^-1 * x * q defines an action on the 3-dimensional vector space
    -- of pure quaternions X (ie linear combinations of i,j,k). It turns out that this action is a rotation of X,
    -- and this is a surjective group homomorphism from H* onto SO3. If we restrict q to the group of unit quaternions
    -- (those of norm 1), then this homomorphism is 2-to-1 (since q and -q give the same rotation).
    -- This shows that the multiplicative group of unit quaternions is isomorphic to Spin3, the double cover of SO3.
    --
    -- @reprSO3 q@ returns the 3*3 matrix representing this map.
    reprSO3 :: (Eq k, Fractional k) => Quaternion k -> [[k]]
    reprSO3 :: Quaternion k -> [[k]]
    reprSO3 q :: Quaternion k
    q = Quaternion k -> Quaternion k -> Quaternion k
    forall a. Fractional a => a -> a -> a
    reprSO3' Quaternion k
    q (Quaternion k -> Quaternion k) -> [Quaternion k] -> [[k]]
    forall a.
    (Num a, Eq a) =>
    (Vect a HBasis -> Vect a HBasis) -> [Vect a HBasis] -> [[a]]
    `asMatrix` [Quaternion k
    forall k. Num k => Quaternion k
    i,Quaternion k
    forall k. Num k => Quaternion k
    j,Quaternion k
    forall k. Num k => Quaternion k
    k]
    -- It's clear from the definition that repr3' q leaves scalars invariant

    -- for achiral elts, ie GO3\SO3, we compose the above with conj

    -- For unit quaternions, this is a double cover of SO4 (since (l,r), (-l,-r) induce same rotation)
    -- Ordered pairs of unit quaternions form the group Spin4
    reprSO4' :: (a, a) -> a -> a
    reprSO4' (l :: a
    l,r :: a
    r) = \x :: a
    x -> a
    la -> Integer -> a
    forall a a. (Eq a, Fractional a, Num a) => a -> a -> a
    ^-1 a -> a -> a
    forall a. Num a => a -> a -> a
    * a
    x a -> a -> a
    forall a. Num a => a -> a -> a
    * a
    r
    -- then (l1,r1) * (l2,r2) -> (l1*l2,r1*r2)
    -- having l^-1 is required for this to work

    -- |Given a pair of unit quaternions (l,r), the map x -> l^-1 * x * r defines an action on the 4-dimensional space
    -- of quaternions. It turns out that this action is a rotation, and this is a surjective group homomorphism
    -- onto SO4. The homomorphism is 2-to-1 (since (l,r) and (-l,-r) give the same map).
    -- This shows that the multiplicative group of pairs of unit quaternions (with pointwise multiplication)
    -- is isomorphic to Spin4, the double cover of SO4.
    --
    -- @reprSO4 (l,r)@ returns the 4*4 matrix representing this map.
    reprSO4 :: (Eq k, Fractional k) => (Quaternion k, Quaternion k) -> [[k]]
    reprSO4 :: (Quaternion k, Quaternion k) -> [[k]]
    reprSO4 (l :: Quaternion k
    l,r :: Quaternion k
    r) = (Quaternion k, Quaternion k) -> Quaternion k -> Quaternion k
    forall a. Fractional a => (a, a) -> a -> a
    reprSO4' (Quaternion k
    l,Quaternion k
    r) (Quaternion k -> Quaternion k) -> [Quaternion k] -> [[k]]
    forall a.
    (Num a, Eq a) =>
    (Vect a HBasis -> Vect a HBasis) -> [Vect a HBasis] -> [[a]]
    `asMatrix` [1,Quaternion k
    forall k. Num k => Quaternion k
    i,Quaternion k
    forall k. Num k => Quaternion k
    j,Quaternion k
    forall k. Num k => Quaternion k
    k]
    -- could consider checking that l,r are unit length - except that this is hard to achieve working over Q

    reprSO4d :: Vect k (DSum HBasis HBasis) -> [[k]]
    reprSO4d lr :: Vect k (DSum HBasis HBasis)
    lr = (Quaternion k, Quaternion k) -> [[k]]
    forall k.
    (Eq k, Fractional k) =>
    (Quaternion k, Quaternion k) -> [[k]]
    reprSO4 (Vect k (DSum HBasis HBasis) -> Quaternion k
    forall k a b. (Eq k, Num k, Ord a) => Vect k (DSum a b) -> Vect k a
    p1 Vect k (DSum HBasis HBasis)
    lr, Vect k (DSum HBasis HBasis) -> Quaternion k
    forall k b a. (Eq k, Num k, Ord b) => Vect k (DSum a b) -> Vect k b
    p2 Vect k (DSum HBasis HBasis)
    lr)

    -- for achiral elts, GO4\SO4, we compose the above with conj

    -- DUAL SPACE OF QUATERNIONS AS COALGEBRA

    one',i',j',k' :: Num k => Vect k (Dual HBasis)
    one' :: Vect k (Dual HBasis)
    one' = Dual HBasis -> Vect k (Dual HBasis)
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    One)
    i' :: Vect k (Dual HBasis)
    i' = Dual HBasis -> Vect k (Dual HBasis)
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    I)
    j' :: Vect k (Dual HBasis)
    j' = Dual HBasis -> Vect k (Dual HBasis)
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    J)
    k' :: Vect k (Dual HBasis)
    k' = Dual HBasis -> Vect k (Dual HBasis)
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    K)

    -- Coalgebra structure on the dual vector space to the quaternions
    -- The comult is the transpose of mult
    instance (Eq k, Num k) => Coalgebra k (Dual HBasis) where
        counit :: Vect k (Dual HBasis) -> k
    counit = Vect k () -> k
    forall k. Num k => Vect k () -> k
    unwrap (Vect k () -> k)
    -> (Vect k (Dual HBasis) -> Vect k ()) -> Vect k (Dual HBasis) -> k
    forall b c a. (b -> c) -> (a -> b) -> a -> c
    . (Dual HBasis -> Vect k ()) -> Vect k (Dual HBasis) -> Vect k ()
    forall k b a.
    (Eq k, Num k, Ord b) =>
    (a -> Vect k b) -> Vect k a -> Vect k b
    linear Dual HBasis -> Vect k ()
    forall k. Num k => Dual HBasis -> Vect k ()
    counit'
            where counit' :: Dual HBasis -> Vect k ()
    counit' (Dual One) = () -> Vect k ()
    forall (m :: * -> *) a. Monad m => a -> m a
    return ()
                  counit' _ = Vect k ()
    forall k b. Vect k b
    zerov
        comult :: Vect k (Dual HBasis) -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    comult = (Dual HBasis -> Vect k (Tensor (Dual HBasis) (Dual HBasis)))
    -> Vect k (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b a.
    (Eq k, Num k, Ord b) =>
    (a -> Vect k b) -> Vect k a -> Vect k b
    linear Dual HBasis -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k.
    (Eq k, Num k) =>
    Dual HBasis -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    comult'
            where comult' :: Dual HBasis -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    comult' (Dual One) = Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    One, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    One) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+>
                      (-1) k
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> ( Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    I, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    I) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+> Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    J, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    J) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+> Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    K, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    K) )
                  comult' (Dual I) = Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    One, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    I) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+> Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    I, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    One) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+>
                      Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    J, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    K) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+> (-1) k
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    K, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    J)
                  comult' (Dual J) = Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    One, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    J) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+> Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    J, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    One) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+>
                      Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    K, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    I) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+> (-1) k
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    I, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    K)
                  comult' (Dual K) = Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    One, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    K) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+> Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    K, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    One) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+>
                      Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    I, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    J) Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b.
    (Eq k, Num k, Ord b) =>
    Vect k b -> Vect k b -> Vect k b
    <+> (-1) k
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall k b. (Eq k, Num k) => k -> Vect k b -> Vect k b
    *> Tensor (Dual HBasis) (Dual HBasis)
    -> Vect k (Tensor (Dual HBasis) (Dual HBasis))
    forall (m :: * -> *) a. Monad m => a -> m a
    return (HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    J, HBasis -> Dual HBasis
    forall b. b -> Dual b
    Dual HBasis
    I)

    {-
    -- Of course, we can define this coalgebra structure on the quaternions themselves
    -- However, it is not compatible with the algebra structure: we don't get a bialgebra
    instance Num k => Coalgebra k HBasis where
        counit = unwrap . linear counit'
            where counit' One = return ()
                  counit' _ = zero
        comult = linear comult'
            where comult' One = return (One,One) <+> (-1) *> ( return (I,I) <+> return (J,J) <+> return (K,K) )
                  comult' I = return (One,I) <+> return (I,One) <+> return (J,K) <+> (-1) *> return (K,J)
                  comult' J = return (One,J) <+> return (J,One) <+> return (K,I) <+> (-1) *> return (I,K)
                  comult' K = return (One,K) <+> return (K,One) <+> return (I,J) <+> (-1) *> return (J,I)
    -}

    {-
    -- Set coalgebra instance
    instance Num k => Coalgebra k HBasis where
        counit (V ts) = sum [x | (m,x) <- ts] -- trace
        comult = fmap (\m -> T m m)           -- diagonal
    -}

    {-
    instance Num k => Coalgebra k HBasis where
        counit (V ts) = sum [x | (One,x) <- ts]
        comult = linear cm
            where cm m = if m == One then return (m,m) else return (m,One) <+> return (One,m)
    -}