Horrific linear algebra macros for when you want to hate something.

macros.rs

use num::traits::{Inv, One, Zero};
use std::ops::*;
use super::F;

macro_rules! vector {($T:ident, $F:ident, $n:tt, [$($i:expr),*]) => {

#[derive(Clone, Copy, Debug, PartialEq)]
pub struct $T(pub [$F; $n]);

impl Deref for $T { type Target = [$F; $n]; fn deref(&self) -> &Self::Target { &self.0 } }
impl DerefMut for $T { fn deref_mut(&mut self) -> &mut Self::Target { &mut self.0 } }

impl Zero for $T {
    fn zero() -> Self { $T([Zero::zero(); $n]) }
    fn is_zero(&self) -> bool { self == &Self::zero() }
}

impl AddAssign for $T { fn add_assign(&mut self, rhs: Self) { $(self[$i] += rhs[$i];)* } }
impl SubAssign for $T { fn sub_assign(&mut self, rhs: Self) { $(self[$i] -= rhs[$i];)* } }
impl MulAssign<$F> for $T { fn mul_assign(&mut self, rhs: $F) { $(self[$i] *= rhs;)* } }
impl DivAssign<$F> for $T { fn div_assign(&mut self, rhs: $F) { $(self[$i] /= rhs;)* } }
impl Add for $T { type Output = Self; fn add(mut self, rhs: Self) -> Self { self += rhs; self } }
impl Sub for $T { type Output = Self; fn sub(mut self, rhs: Self) -> Self { self -= rhs; self } }
impl Mul<$F> for $T { type Output = Self; fn mul(mut self, rhs: $F) -> Self { self *= rhs; self } }
impl Div<$F> for $T { type Output = Self; fn div(mut self, rhs: $F) -> Self { self /= rhs; self } }
impl Neg for $T { type Output = Self; fn neg(mut self) -> Self { $(self[$i] = -self[$i];)* self } }

impl $T {
    pub fn dot(&self, rhs: &Self) -> $F {
        $F::zero() $(+ self[$i] * rhs[$i])*
    }

    pub fn length(&self) -> $F {
        self.dot(self).sqrt()
    }

    pub fn normalize(&self) -> Self {
        *self / self.length()
    }
}

}}

macro_rules! matrix {($T:ident, $V:ident, $F:ident, $n:tt, [$($($i:expr),*;)*]) => {

vector!($T, $F, ($n * $n), [$($($i),*),*]);

impl<T: Fn(usize, usize) -> $F> From<T> for $T {
    fn from(f: T) -> Self {
        $T([$($(f($i/$n, $i%$n)),*),*])
    }
}

impl $T {
    pub fn transpose(self) -> Self {
        Self::from(|x, y| self[x*$n + y])
    }
}

impl One for $T {
    fn one() -> Self {
        Self::from(|x, y| if x == y { One::one() } else { Zero::zero() })
    }
}

impl Mul<$V> for $T {
    type Output = $V;
    fn mul(self, rhs: $V) -> $V {
        $V([$(rhs.dot(&$V([$(self[$i]),*]))),*])
    }
}

impl Mul for $T {
    type Output = Self;
    fn mul(self, rhs: Self) -> Self {
        $T([$($({
            (0..$n)
                .map(|k| self[$n*($i/$n) + k] * rhs[$n*k + ($i%$n)])
                .sum()
        }),*),*])
    }
}

impl Inv for $T {
    type Output = Self;
    fn inv(self) -> Self {
        matrix_inverse!($T, $n, self)
    }
}

}}

macro_rules! matrix_inverse {
    ($T:ident, 0, $self:expr) => {
        $self
    };
    ($T:ident, 1, $self:expr) => {
        $T([$self[0].inv()])
    };
    ($T:ident, 2, $self:expr) => {{
        let [a, b, c, d] = *$self;
        let idet = (a*b - c*d).inv();
        $T([d, -b, -c, a]) * idet
    }};
    ($T:ident, 3, $self:expr) => {{
        let [a, b, c, d, e, f, g, h, i] = *$self;

        let mut matrix = $T([
            e*i - f*h, c*h - b*i, b*f - c*e,
            f*g - d*i, a*i - c*g, c*d - a*f,
            d*h - e*g, b*g - a*h, a*e - b*d,
        ]);

        matrix *= (a*matrix[0] + b*matrix[3] + c*matrix[6]).inv();
        matrix
    }};
    ($T:ident, 4, $self:expr) => {
        unimplemented!() //TODO: Compute the 4D inverse.
    };
    ($T:ident, $n:expr, $self:expr) => {
        unimplemented!() //TODO: Compute the inverse in higher dimensions.
    };
}