Add a `SquareMat` trait for square matrices.

src/lib.rs

@@ -152,6 +152,7 @@ pub use traits::{
     ScalarAdd, ScalarSub,
     ScalarMul, ScalarDiv,
     Shape,
+    SquareMat,
     ToHomogeneous,
     Transform, Transformation,
     Translate, Translation,
@@ -644,7 +645,7 @@ pub fn append_rotation_wrt_center<LV: Neg<LV>,
 
 /// Builds a rotation matrix from `r`.
 #[inline(always)]
-pub fn to_rot_mat<LV, AV, M: Mat<LV, LV> + Rotation<AV>, R: RotationMatrix<LV, AV, M>>(r: &R) -> M {
+pub fn to_rot_mat<N, LV, AV, M: Mat<N, LV, LV> + Rotation<AV>, R: RotationMatrix<N, LV, AV, M>>(r: &R) -> M {
     r.to_rot_mat()
 }
 

src/linalg/decompositions.rs

@@ -1,6 +1,6 @@
 use std::num::{Zero, Float};
 use traits::operations::{Transpose, ApproxEq};
-use traits::structure::{ColSlice, Eye, Indexable, Diag};
+use traits::structure::{ColSlice, Eye, Indexable, Diag, SquareMat};
 use traits::geometry::Norm;
 use std::cmp::min;
 
@@ -73,19 +73,17 @@ pub fn qr<N, V, M>(m: &M) -> (M, M)
 }
 
 /// Eigendecomposition of a square matrix using the qr algorithm.
-pub fn eigen_qr<N, V, V2, M>(m: &M, eps: &N, niter: uint) -> (M, V2)
+pub fn eigen_qr<N, V, VS, M>(m: &M, eps: &N, niter: uint) -> (M, V)
     where N:  Float,
-          V:  Indexable<uint, N> + Norm<N>,
-          V2: Zero,
-          M:  Clone + Eye + ColSlice<V> + Transpose + Indexable<(uint, uint), N> + Mul<M, M>
-              + Diag<V2> + ApproxEq<N> + Add<M, M> + Sub<M, M> {
+          VS: Indexable<uint, N> + Norm<N>,
+          M:  Indexable<(uint, uint), N> + SquareMat<N, V> + ColSlice<VS> + ApproxEq<N> + Clone {
     let (rows, cols) = m.shape();
 
     assert!(rows == cols, "The matrix being decomposed must be square.");
 
     let mut eigenvectors: M = Eye::new_identity(rows);
     let mut eigenvalues = m.clone();
-    let mut shifter: M = Eye::new_identity(rows);
+    // let mut shifter: M = Eye::new_identity(rows);
 
     let mut iter = 0u;
     for _ in range(0, niter) {
@@ -112,16 +110,9 @@ pub fn eigen_qr<N, V, V2, M>(m: &M, eps: &N, niter: uint) -> (M, V2)
         }
         iter = iter + 1;
 
-        // FIXME: This is a very naive implementation.
-        let shift = unsafe { eigenvalues.unsafe_at((rows - 1, rows - 1)) };
+        let (q, r) = qr(&eigenvalues);;
 
-        for i in range(0, rows) {
-            unsafe { shifter.unsafe_set((i, i), shift.clone()) }
-        }
-
-        let (q, r) = qr(&eigenvalues);//  - shifter));
-
-        eigenvalues = r * q /*+ shifter*/;
+        eigenvalues  = r * q;
         eigenvectors = eigenvectors * q;
     }
 

src/structs/iso_macros.rs

@@ -27,7 +27,7 @@ macro_rules! iso_impl(
 macro_rules! rotation_matrix_impl(
     ($t: ident, $trot: ident, $tlv: ident, $tav: ident) => (
         impl<N: Cast<f32> + FloatMath + Num + Clone>
-        RotationMatrix<$tlv<N>, $tav<N>, $trot<N>> for $t<N> {
+        RotationMatrix<N, $tlv<N>, $tav<N>, $trot<N>> for $t<N> {
             #[inline]
             fn to_rot_mat(&self) -> $trot<N> {
                 self.rotation.clone()

src/structs/rot_macros.rs

@@ -122,8 +122,7 @@ macro_rules! dim_impl(
 
 macro_rules! rotation_matrix_impl(
     ($t: ident, $tlv: ident, $tav: ident) => (
-        impl<N: Cast<f32> + FloatMath + Clone>
-        RotationMatrix<$tlv<N>, $tav<N>, $t<N>> for $t<N> {
+        impl<N: Cast<f32> + FloatMath> RotationMatrix<N, $tlv<N>, $tav<N>, $t<N>> for $t<N> {
             #[inline]
             fn to_rot_mat(&self) -> $t<N> {
                 self.clone()

src/traits/geometry.rs

@@ -141,7 +141,7 @@ impl<LV: Neg<LV>, AV, M: Rotation<AV> + Translation<LV>> RotationWithTranslation
 
 /// Trait of transformation having a rotation extractable as a rotation matrix. This can typically
 /// be implemented by quaternions to convert them to a rotation matrix.
-pub trait RotationMatrix<LV, AV, M: Mat<LV, LV> + Rotation<AV>> : Rotation<AV> {
+pub trait RotationMatrix<N, LV, AV, M: Mat<N, LV, LV> + Rotation<AV>> : Rotation<AV> {
     /// Gets the rotation matrix represented by `self`.
     fn to_rot_mat(&self) -> M;
 }

src/traits/mod.rs

@@ -1,16 +1,16 @@
 //! Mathematical traits.
 
-pub use self::geometry::{AbsoluteRotate, Cross, CrossMatrix, Dot, FromHomogeneous, Norm, Orig,
+pub use traits::geometry::{AbsoluteRotate, Cross, CrossMatrix, Dot, FromHomogeneous, Norm, Orig,
                            Rotate, Rotation, RotationMatrix, RotationWithTranslation, ToHomogeneous,
                            Transform, Transformation, Translate, Translation, UniformSphereSample};
 
-pub use self::structure::{FloatVec, FloatPnt, Basis, Cast, Col, Dim, Indexable, Iterable,
-                          IterableMut, Mat, Row, NumVec, NumPnt, PntAsVec, VecAsPnt, ColSlice,
-                          RowSlice, Diag, Eye, Shape};
+pub use traits::structure::{FloatVec, FloatPnt, Basis, Cast, Col, Dim, Indexable, Iterable,
+                            IterableMut, Mat, SquareMat, Row, NumVec, NumPnt, PntAsVec, VecAsPnt,
+                            ColSlice, RowSlice, Diag, Eye, Shape};
 
-pub use self::operations::{Absolute, ApproxEq, Axpy, Cov, Det, Inv, LMul, Mean, Outer, POrd,
+pub use traits::operations::{Absolute, ApproxEq, Axpy, Cov, Det, Inv, LMul, Mean, Outer, POrd,
                              RMul, ScalarAdd, ScalarSub, ScalarMul, ScalarDiv, Transpose};
-pub use self::operations::{POrdering, PartialLess, PartialEqual, PartialGreater, NotComparable};
+pub use traits::operations::{POrdering, PartialLess, PartialEqual, PartialGreater, NotComparable};
 
 pub mod geometry;
 pub mod structure;

src/traits/structure.rs

@@ -2,7 +2,7 @@
 
 use std::num::Zero;
 use std::slice::{Items, MutItems};
-use traits::operations::{RMul, LMul, Axpy};
+use traits::operations::{RMul, LMul, Axpy, Transpose};
 use traits::geometry::{Dot, Norm, Orig};
 
 /// Traits of objects which can be created from an object of type `T`.
@@ -14,9 +14,19 @@ pub trait Cast<T> {
 /// Trait of matrices.
 ///
 /// A matrix has rows and columns and are able to multiply them.
-pub trait Mat<R, C> : Row<R> + Col<C> + RMul<R> + LMul<C> { }
+pub trait Mat<N, R, C>: Row<R> + Col<C> + RMul<R> + LMul<C> + Index<(uint, uint), N> { }
 
-impl<M: Row<R> + Col<C> + RMul<R> + LMul<C>, R, C> Mat<R, C> for M {
+impl<N, M, R, C> Mat<N, R, C> for M
+    where M: Row<R> + Col<C> + RMul<R> + LMul<C> + Index<(uint, uint), N> {
+}
+
+/// Trait implemented by square matrices.
+pub trait SquareMat<N, V>: Mat<N, V, V> + Mul<Self, Self> + Eye + Transpose + Add<Self, Self> +
+                           Sub<Self, Self> + Diag<V> {
+}
+
+impl<N, V, M> SquareMat<N, V> for M
+    where M: Mat<N, V, V> + Mul<M, M> + Eye + Transpose + Add<M, M> + Sub<M, M> + Diag<V> {
 }
 
 /// Trait for constructing the identity matrix

tests/mat.rs

@@ -6,7 +6,7 @@ use std::num::{Float, abs};
 use std::rand::random;
 use std::cmp::{min, max};
 use na::{Vec1, Vec3, Mat1, Mat2, Mat3, Mat4, Mat5, Mat6, Rot3, Persp3, PerspMat3, Ortho3, OrthoMat3,
-         DMat, DVec, Indexable, Row, Col};
+         DMat, DVec, Row, Col};
 
 macro_rules! test_inv_mat_impl(
   ($t: ty) => (
@@ -141,7 +141,7 @@ fn test_rotation2() {
 fn test_index_mat2() {
   let mat: Mat2<f64> = random();
 
-  assert!(mat.at((0, 1)) == na::transpose(&mat).at((1, 0)));
+  assert!(mat[(0, 1)] == na::transpose(&mat)[(1, 0)]);
 }
 
 #[test]