quaternion.f90

module quart
    type quaternion
        ! q(4) = (w, x, y, z) where w is the scalar part
        real :: q(4)
        ! Flags to tell us if we want to keep the quaternion unitised, and
        ! also if the quaternion represents an improper rotation (rotation + inversion)
        logical :: unit
        logical :: improper
    end type quaternion
contains
    elemental subroutine normalise_quaternion(q)

        ! Normalise the quaternion by dividing the quaternion by it's
        ! 'reduced norm'

        type (quaternion), intent(inout) :: q

        q%q = q%q / (normq(q))

    end subroutine normalise_quaternion
    pure function normq(q)

        ! The 'reduced norm' of a quaternion is just the square root of the
        ! arithmetic sum of the squares of it's elements

        real :: normq
        type (quaternion), intent(in) :: q

        normq = sqrt(sum(q%q**2))

    end function normq
end module

program test
    use quart

    type(quaternion) :: q(2)

    q(1)%q = [1,2,2,4]
    q(2)%q = [1,0,0,0]

    call normalise_quaternion(q)

    write(*,*) q(1)%q
    write(*,*) q(2)%q


end program