diff --git a/Cargo.toml b/Cargo.toml
index 6da6b04d..7c1f4b1a 100644
--- a/Cargo.toml
+++ b/Cargo.toml
@@ -48,6 +48,9 @@ quickcheck     = { version = "0.8", optional = true }
 pest           = { version = "2.0", optional = true }
 pest_derive    = { version = "2.0", optional = true }
 
+[patch.crates-io]
+alga = { git = "https://github.com/rustsim/alga", branch = "dev" }
+
 [dev-dependencies]
 serde_json = "1.0"
 rand_xorshift = "0.1"
diff --git a/src/base/blas.rs b/src/base/blas.rs
index e4d7ea7d..6e1b6c81 100644
--- a/src/base/blas.rs
+++ b/src/base/blas.rs
@@ -494,7 +494,7 @@ where
 {
     /// Computes `self = a * x + b * self`.
     ///
-    /// If be is zero, `self` is never read from.
+    /// If `b` is zero, `self` is never read from.
     ///
     /// # Examples:
     ///
diff --git a/src/base/ops.rs b/src/base/ops.rs
index 9f3d44c3..cf921c33 100644
--- a/src/base/ops.rs
+++ b/src/base/ops.rs
@@ -11,7 +11,7 @@ use crate::base::allocator::{Allocator, SameShapeAllocator, SameShapeC, SameShap
 use crate::base::constraint::{
     AreMultipliable, DimEq, SameNumberOfColumns, SameNumberOfRows, ShapeConstraint,
 };
-use crate::base::dimension::{Dim, DimMul, DimName, DimProd};
+use crate::base::dimension::{Dim, DimMul, DimName, DimProd, Dynamic};
 use crate::base::storage::{ContiguousStorageMut, Storage, StorageMut};
 use crate::base::{DefaultAllocator, Matrix, MatrixMN, MatrixN, MatrixSum, Scalar, VectorSliceN};
 
@@ -384,6 +384,36 @@ where
     }
 }
 
+impl<N, C: Dim> iter::Sum for MatrixMN<N, Dynamic, C>
+where
+    N: Scalar + ClosedAdd + Zero,
+    DefaultAllocator: Allocator<N, Dynamic, C>,
+{
+    /// # Example
+    /// ```
+    /// # use nalgebra::DVector;
+    /// assert_eq!(vec![DVector::repeat(3, 1.0f64),
+    ///                 DVector::repeat(3, 1.0f64),
+    ///                 DVector::repeat(3, 1.0f64)].into_iter().sum::<DVector<f64>>(),
+    ///            DVector::repeat(3, 1.0f64) + DVector::repeat(3, 1.0f64) + DVector::repeat(3, 1.0f64));
+    /// ```
+    ///
+    /// # Panics
+    /// Panics if the iterator is empty:
+    /// ```should_panic
+    /// # use std::iter;
+    /// # use nalgebra::DMatrix;
+    /// iter::empty::<DMatrix<f64>>().sum::<DMatrix<f64>>(); // panics!
+    /// ```
+    fn sum<I: Iterator<Item = MatrixMN<N, Dynamic, C>>>(mut iter: I) -> MatrixMN<N, Dynamic, C> {
+        if let Some(first) = iter.next() {
+            iter.fold(first, |acc, x| acc + x)
+        } else {
+            panic!("Cannot compute `sum` of empty iterator.")
+        }
+    }
+}
+
 impl<'a, N, R: DimName, C: DimName> iter::Sum<&'a MatrixMN<N, R, C>> for MatrixMN<N, R, C>
 where
     N: Scalar + ClosedAdd + Zero,
@@ -394,6 +424,36 @@ where
     }
 }
 
+impl<'a, N, C: Dim> iter::Sum<&'a MatrixMN<N, Dynamic, C>> for MatrixMN<N, Dynamic, C>
+where
+    N: Scalar + ClosedAdd + Zero,
+    DefaultAllocator: Allocator<N, Dynamic, C>,
+{
+    /// # Example
+    /// ```
+    /// # use nalgebra::DVector;
+    /// let v = &DVector::repeat(3, 1.0f64);
+    ///
+    /// assert_eq!(vec![v, v, v].into_iter().sum::<DVector<f64>>(),
+    ///            v + v + v);
+    /// ```
+    ///
+    /// # Panics
+    /// Panics if the iterator is empty:
+    /// ```should_panic
+    /// # use std::iter;
+    /// # use nalgebra::DMatrix;
+    /// iter::empty::<&DMatrix<f64>>().sum::<DMatrix<f64>>(); // panics!
+    /// ```
+    fn sum<I: Iterator<Item = &'a MatrixMN<N, Dynamic, C>>>(mut iter: I) -> MatrixMN<N, Dynamic, C> {
+        if let Some(first) = iter.next() {
+            iter.fold(first.clone(), |acc, x| acc + x)
+        } else {
+            panic!("Cannot compute `sum` of empty iterator.")
+        }
+    }
+}
+
 /*
  *
  * Multiplication
diff --git a/src/base/unit.rs b/src/base/unit.rs
index 0c397e52..6fe831da 100644
--- a/src/base/unit.rs
+++ b/src/base/unit.rs
@@ -222,7 +222,7 @@ impl<T: Neg> Neg for Unit<T> {
 
     #[inline]
     fn neg(self) -> Self::Output {
-        Unit::new_unchecked(-self.value)
+        Self::Output::new_unchecked(-self.value)
     }
 }
 
diff --git a/src/geometry/isometry.rs b/src/geometry/isometry.rs
old mode 100644
new mode 100755
index 6df20039..888a8307
--- a/src/geometry/isometry.rs
+++ b/src/geometry/isometry.rs
@@ -17,7 +17,7 @@ use alga::linear::Rotation;
 use crate::base::allocator::Allocator;
 use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
 use crate::base::storage::Owned;
-use crate::base::{DefaultAllocator, MatrixN};
+use crate::base::{DefaultAllocator, MatrixN, VectorN};
 use crate::geometry::{Point, Translation};
 
 /// A direct isometry, i.e., a rotation followed by a translation, aka. a rigid-body motion, aka. an element of a Special Euclidean (SE) group.
@@ -254,6 +254,97 @@ where DefaultAllocator: Allocator<N, D>
     pub fn append_rotation_wrt_center_mut(&mut self, r: &R) {
         self.rotation = self.rotation.append_rotation(r);
     }
+
+    /// Transform the given point by this isometry.
+    ///
+    /// This is the same as the multiplication `self * pt`.
+    ///
+    /// # Example
+    ///
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Isometry3, Translation3, UnitQuaternion, Vector3, Point3};
+    /// let tra = Translation3::new(0.0, 0.0, 3.0);
+    /// let rot = UnitQuaternion::from_scaled_axis(Vector3::y() * f32::consts::FRAC_PI_2);
+    /// let iso = Isometry3::from_parts(tra, rot);
+    ///
+    /// let transformed_point = iso.transform_point(&Point3::new(1.0, 2.0, 3.0));
+    /// assert_relative_eq!(transformed_point, Point3::new(3.0, 2.0, 2.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        self * pt
+    }
+
+    /// Transform the given vector by this isometry, ignoring the translation
+    /// component of the isometry.
+    ///
+    /// This is the same as the multiplication `self * v`.
+    ///
+    /// # Example
+    ///
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Isometry3, Translation3, UnitQuaternion, Vector3};
+    /// let tra = Translation3::new(0.0, 0.0, 3.0);
+    /// let rot = UnitQuaternion::from_scaled_axis(Vector3::y() * f32::consts::FRAC_PI_2);
+    /// let iso = Isometry3::from_parts(tra, rot);
+    ///
+    /// let transformed_point = iso.transform_vector(&Vector3::new(1.0, 2.0, 3.0));
+    /// assert_relative_eq!(transformed_point, Vector3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
+        self * v
+    }
+
+    /// Transform the given point by the inverse of this isometry. This may be
+    /// less expensive than computing the entire isometry inverse and then
+    /// transforming the point.
+    ///
+    /// # Example
+    ///
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Isometry3, Translation3, UnitQuaternion, Vector3, Point3};
+    /// let tra = Translation3::new(0.0, 0.0, 3.0);
+    /// let rot = UnitQuaternion::from_scaled_axis(Vector3::y() * f32::consts::FRAC_PI_2);
+    /// let iso = Isometry3::from_parts(tra, rot);
+    ///
+    /// let transformed_point = iso.inverse_transform_point(&Point3::new(1.0, 2.0, 3.0));
+    /// assert_relative_eq!(transformed_point, Point3::new(0.0, 2.0, 1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        self.rotation
+            .inverse_transform_point(&(pt - &self.translation.vector))
+    }
+
+    /// Transform the given vector by the inverse of this isometry, ignoring the
+    /// translation component of the isometry. This may be
+    /// less expensive than computing the entire isometry inverse and then
+    /// transforming the point.
+    ///
+    /// # Example
+    ///
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Isometry3, Translation3, UnitQuaternion, Vector3};
+    /// let tra = Translation3::new(0.0, 0.0, 3.0);
+    /// let rot = UnitQuaternion::from_scaled_axis(Vector3::y() * f32::consts::FRAC_PI_2);
+    /// let iso = Isometry3::from_parts(tra, rot);
+    ///
+    /// let transformed_point = iso.inverse_transform_vector(&Vector3::new(1.0, 2.0, 3.0));
+    /// assert_relative_eq!(transformed_point, Vector3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
+        self.rotation.inverse_transform_vector(v)
+    }
 }
 
 // NOTE: we don't require `R: Rotation<...>` here because this is not useful for the implementation
diff --git a/src/geometry/isometry_alga.rs b/src/geometry/isometry_alga.rs
old mode 100644
new mode 100755
index 8402c5cb..e68b269c
--- a/src/geometry/isometry_alga.rs
+++ b/src/geometry/isometry_alga.rs
@@ -85,12 +85,12 @@ where
 {
     #[inline]
     fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        self * pt
+        self.transform_point(pt)
     }
 
     #[inline]
     fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
-        self * v
+        self.transform_vector(v)
     }
 }
 
@@ -101,13 +101,12 @@ where
 {
     #[inline]
     fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        self.rotation
-            .inverse_transform_point(&(pt - &self.translation.vector))
+        self.inverse_transform_point(pt)
     }
 
     #[inline]
     fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
-        self.rotation.inverse_transform_vector(v)
+        self.inverse_transform_vector(v)
     }
 }
 
diff --git a/src/geometry/quaternion.rs b/src/geometry/quaternion.rs
old mode 100644
new mode 100755
index 5306f9c1..208145d7
--- a/src/geometry/quaternion.rs
+++ b/src/geometry/quaternion.rs
@@ -17,9 +17,9 @@ use alga::general::RealField;
 
 use crate::base::dimension::{U1, U3, U4};
 use crate::base::storage::{CStride, RStride};
-use crate::base::{Matrix3, MatrixN, MatrixSlice, MatrixSliceMut, Unit, Vector3, Vector4};
+use crate::base::{Matrix3, Matrix4, MatrixSlice, MatrixSliceMut, Unit, Vector3, Vector4};
 
-use crate::geometry::Rotation;
+use crate::geometry::{Point3, Rotation};
 
 /// A quaternion. See the type alias `UnitQuaternion = Unit<Quaternion>` for a quaternion
 /// that may be used as a rotation.
@@ -124,6 +124,12 @@ impl<N: RealField> Quaternion<N> {
         Self::from(self.coords.normalize())
     }
 
+    /// The imaginary part of this quaternion.
+    #[inline]
+    pub fn imag(&self) -> Vector3<N> {
+        self.coords.xyz()
+    }
+
     /// The conjugate of this quaternion.
     ///
     /// # Example
@@ -135,13 +141,7 @@ impl<N: RealField> Quaternion<N> {
     /// ```
     #[inline]
     pub fn conjugate(&self) -> Self {
-        let v = Vector4::new(
-            -self.coords[0],
-            -self.coords[1],
-            -self.coords[2],
-            self.coords[3],
-        );
-        Self::from(v)
+        Self::from_parts(self.w, -self.imag())
     }
 
     /// Inverts this quaternion if it is not zero.
@@ -307,6 +307,81 @@ impl<N: RealField> Quaternion<N> {
         self.coords.dot(&rhs.coords)
     }
 
+    /// Calculates the inner product (also known as the dot product).
+    /// See "Foundations of Game Engine Development, Volume 1: Mathematics" by Lengyel
+    /// Formula 4.89.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let a = Quaternion::new(0.0, 2.0, 3.0, 4.0);
+    /// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
+    /// let expected = Quaternion::new(-20.0, 0.0, 0.0, 0.0);
+    /// let result = a.inner(&b);
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-5);
+    #[inline]
+    pub fn inner(&self, other: &Self) -> Self {
+        (self * other + other * self).half()
+    }
+
+    /// Calculates the outer product (also known as the wedge product).
+    /// See "Foundations of Game Engine Development, Volume 1: Mathematics" by Lengyel
+    /// Formula 4.89.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let a = Quaternion::new(0.0, 2.0, 3.0, 4.0);
+    /// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
+    /// let expected = Quaternion::new(0.0, -5.0, 18.0, -11.0);
+    /// let result = a.outer(&b);
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-5);
+    /// ```
+    #[inline]
+    pub fn outer(&self, other: &Self) -> Self {
+        (self * other - other * self).half()
+    }
+
+    /// Calculates the projection of `self` onto `other` (also known as the parallel).
+    /// See "Foundations of Game Engine Development, Volume 1: Mathematics" by Lengyel
+    /// Formula 4.94.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let a = Quaternion::new(0.0, 2.0, 3.0, 4.0);
+    /// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
+    /// let expected = Quaternion::new(0.0, 3.333333333333333, 1.3333333333333333, 0.6666666666666666);
+    /// let result = a.project(&b).unwrap();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-5);
+    /// ```
+    #[inline]
+    pub fn project(&self, other: &Self) -> Option<Self> {
+        self.inner(other).right_div(other)
+    }
+
+    /// Calculates the rejection of `self` from `other` (also known as the perpendicular).
+    /// See "Foundations of Game Engine Development, Volume 1: Mathematics" by Lengyel
+    /// Formula 4.94.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let a = Quaternion::new(0.0, 2.0, 3.0, 4.0);
+    /// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
+    /// let expected = Quaternion::new(0.0, -1.3333333333333333, 1.6666666666666665, 3.3333333333333335);
+    /// let result = a.reject(&b).unwrap();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-5);
+    /// ```
+    #[inline]
+    pub fn reject(&self, other: &Self) -> Option<Self> {
+        self.outer(other).right_div(other)
+    }
+
     /// The polar decomposition of this quaternion.
     ///
     /// Returns, from left to right: the quaternion norm, the half rotation angle, the rotation
@@ -506,6 +581,274 @@ impl<N: RealField> Quaternion<N> {
     pub fn normalize_mut(&mut self) -> N {
         self.coords.normalize_mut()
     }
+
+    /// Calculates square of a quaternion.
+    #[inline]
+    pub fn squared(&self) -> Self {
+        self * self
+    }
+
+    /// Divides quaternion into two.
+    #[inline]
+    pub fn half(&self) -> Self {
+        self / ::convert(2.0f64)
+    }
+
+    /// Calculates square root.
+    #[inline]
+    pub fn sqrt(&self) -> Self {
+        self.powf(::convert(0.5))
+    }
+
+    /// Check if the quaternion is pure.
+    #[inline]
+    pub fn is_pure(&self) -> bool {
+        self.w.is_zero()
+    }
+
+    /// Convert quaternion to pure quaternion.
+    #[inline]
+    pub fn pure(&self) -> Self {
+        Self::from_imag(self.imag())
+    }
+
+    /// Left quaternionic division.
+    ///
+    /// Calculates B<sup>-1</sup> * A where A = self, B = other.
+    #[inline]
+    pub fn left_div(&self, other: &Self) -> Option<Self> {
+        other.try_inverse().map(|inv| inv * self)
+    }
+
+    /// Right quaternionic division.
+    ///
+    /// Calculates A * B<sup>-1</sup> where A = self, B = other.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let a = Quaternion::new(0.0, 1.0, 2.0, 3.0);
+    /// let b = Quaternion::new(0.0, 5.0, 2.0, 1.0);
+    /// let result = a.right_div(&b).unwrap();
+    /// let expected = Quaternion::new(0.4, 0.13333333333333336, -0.4666666666666667, 0.26666666666666666);
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn right_div(&self, other: &Self) -> Option<Self> {
+        other.try_inverse().map(|inv| self * inv)
+    }
+
+    /// Calculates the quaternionic cosinus.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let expected = Quaternion::new(58.93364616794395, -34.086183690465596, -51.1292755356984, -68.17236738093119);
+    /// let result = input.cos();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn cos(&self) -> Self {
+        let z = self.imag().magnitude();
+        let w = -self.w.sin() * z.sinhc();
+        Self::from_parts(self.w.cos() * z.cosh(), self.imag() * w)
+    }
+
+    /// Calculates the quaternionic arccosinus.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let result = input.cos().acos();
+    /// assert_relative_eq!(input, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn acos(&self) -> Self {
+        let u = Self::from_imag(self.imag().normalize());
+        let identity = Self::identity();
+
+        let z = (self + (self.squared() - identity).sqrt()).ln();
+
+        -(u * z)
+    }
+
+    /// Calculates the quaternionic sinus.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let expected = Quaternion::new(91.78371578403467, 21.886486853029176, 32.82973027954377, 43.77297370605835);
+    /// let result = input.sin();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn sin(&self) -> Self {
+        let z = self.imag().magnitude();
+        let w = self.w.cos() * z.sinhc();
+        Self::from_parts(self.w.sin() * z.cosh(), self.imag() * w)
+    }
+
+    /// Calculates the quaternionic arcsinus.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let result = input.sin().asin();
+    /// assert_relative_eq!(input, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn asin(&self) -> Self {
+        let u = Self::from_imag(self.imag().normalize());
+        let identity = Self::identity();
+
+        let z = ((u * self) + (identity - self.squared()).sqrt()).ln();
+
+        -(u * z)
+    }
+
+    /// Calculates the quaternionic tangent.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let expected = Quaternion::new(0.00003821631725009489, 0.3713971716439371, 0.5570957574659058, 0.7427943432878743);
+    /// let result = input.tan();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn tan(&self) -> Self {
+        self.sin().right_div(&self.cos()).unwrap()
+    }
+
+    /// Calculates the quaternionic arctangent.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let result = input.tan().atan();
+    /// assert_relative_eq!(input, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn atan(&self) -> Self {
+        let u = Self::from_imag(self.imag().normalize());
+        let num = u + self;
+        let den = u - self;
+        let fr = num.right_div(&den).unwrap();
+        let ln = fr.ln();
+        (u.half()) * ln
+    }
+
+    /// Calculates the hyperbolic quaternionic sinus.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let expected = Quaternion::new(0.7323376060463428, -0.4482074499805421, -0.6723111749708133, -0.8964148999610843);
+    /// let result = input.sinh();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn sinh(&self) -> Self {
+        (self.exp() - (-self).exp()).half()
+    }
+
+    /// Calculates the hyperbolic quaternionic arcsinus.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let expected = Quaternion::new(2.385889902585242, 0.514052600662788, 0.7710789009941821, 1.028105201325576);
+    /// let result = input.asinh();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn asinh(&self) -> Self {
+        let identity = Self::identity();
+        (self + (identity + self.squared()).sqrt()).ln()
+    }
+
+    /// Calculates the hyperbolic quaternionic cosinus.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let expected = Quaternion::new(0.9615851176369566, -0.3413521745610167, -0.5120282618415251, -0.6827043491220334);
+    /// let result = input.cosh();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn cosh(&self) -> Self {
+        (self.exp() + (-self).exp()).half()
+    }
+
+    /// Calculates the hyperbolic quaternionic arccosinus.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let expected = Quaternion::new(2.4014472020074007, 0.5162761016176176, 0.7744141524264264, 1.0325522032352352);
+    /// let result = input.acosh();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn acosh(&self) -> Self {
+        let identity = Self::identity();
+        (self + (self + identity).sqrt() * (self - identity).sqrt()).ln()
+    }
+
+    /// Calculates the hyperbolic quaternionic tangent.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let expected = Quaternion::new(1.0248695360556623, -0.10229568178876419, -0.1534435226831464, -0.20459136357752844);
+    /// let result = input.tanh();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn tanh(&self) -> Self {
+        self.sinh().right_div(&self.cosh()).unwrap()
+    }
+
+    /// Calculates the hyperbolic quaternionic arctangent.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::Quaternion;
+    /// let input = Quaternion::new(1.0, 2.0, 3.0, 4.0);
+    /// let expected = Quaternion::new(0.03230293287000163, 0.5173453683196951, 0.7760180524795426, 1.0346907366393903);
+    /// let result = input.atanh();
+    /// assert_relative_eq!(expected, result, epsilon = 1.0e-7);
+    /// ```
+    #[inline]
+    pub fn atanh(&self) -> Self {
+        let identity = Self::identity();
+        ((identity + self).ln() - (identity - self).ln()).half()
+    }
 }
 
 impl<N: RealField + AbsDiffEq<Epsilon = N>> AbsDiffEq for Quaternion<N> {
@@ -847,16 +1190,12 @@ impl<N: RealField> UnitQuaternion<N> {
     /// ```
     #[inline]
     pub fn axis_angle(&self) -> Option<(Unit<Vector3<N>>, N)> {
-        if let Some(axis) = self.axis() {
-            Some((axis, self.angle()))
-        } else {
-            None
-        }
+        self.axis().map(|axis| (axis, self.angle()))
     }
 
     /// Compute the exponential of a quaternion.
     ///
-    /// Note that this function yields a `Quaternion<N>` because it looses the unit property.
+    /// Note that this function yields a `Quaternion<N>` because it loses the unit property.
     #[inline]
     pub fn exp(&self) -> Quaternion<N> {
         self.as_ref().exp()
@@ -864,7 +1203,7 @@ impl<N: RealField> UnitQuaternion<N> {
 
     /// Compute the natural logarithm of a quaternion.
     ///
-    /// Note that this function yields a `Quaternion<N>` because it looses the unit property.
+    /// Note that this function yields a `Quaternion<N>` because it loses the unit property.
     /// The vector part of the return value corresponds to the axis-angle representation (divided
     /// by 2.0) of this unit quaternion.
     ///
@@ -879,7 +1218,7 @@ impl<N: RealField> UnitQuaternion<N> {
     #[inline]
     pub fn ln(&self) -> Quaternion<N> {
         if let Some(v) = self.axis() {
-            Quaternion::from_parts(N::zero(), v.into_inner() * self.angle())
+            Quaternion::from_imag(v.into_inner() * self.angle())
         } else {
             Quaternion::zero()
         }
@@ -1002,9 +1341,91 @@ impl<N: RealField> UnitQuaternion<N> {
     /// assert_relative_eq!(rot.to_homogeneous(), expected, epsilon = 1.0e-6);
     /// ```
     #[inline]
-    pub fn to_homogeneous(&self) -> MatrixN<N, U4> {
+    pub fn to_homogeneous(&self) -> Matrix4<N> {
         self.to_rotation_matrix().to_homogeneous()
     }
+
+    /// Rotate a point by this unit quaternion.
+    ///
+    /// This is the same as the multiplication `self * pt`.
+    ///
+    /// # Example
+    ///
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{UnitQuaternion, Vector3, Point3};
+    /// let rot = UnitQuaternion::from_axis_angle(&Vector3::y_axis(), f32::consts::FRAC_PI_2);
+    /// let transformed_point = rot.transform_point(&Point3::new(1.0, 2.0, 3.0));
+    ///
+    /// assert_relative_eq!(transformed_point, Point3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn transform_point(&self, pt: &Point3<N>) -> Point3<N> {
+        self * pt
+    }
+
+    /// Rotate a vector by this unit quaternion.
+    ///
+    /// This is the same as the multiplication `self * v`.
+    ///
+    /// # Example
+    ///
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{UnitQuaternion, Vector3};
+    /// let rot = UnitQuaternion::from_axis_angle(&Vector3::y_axis(), f32::consts::FRAC_PI_2);
+    /// let transformed_vector = rot.transform_vector(&Vector3::new(1.0, 2.0, 3.0));
+    ///
+    /// assert_relative_eq!(transformed_vector, Vector3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn transform_vector(&self, v: &Vector3<N>) -> Vector3<N> {
+        self * v
+    }
+
+    /// Rotate a point by the inverse of this unit quaternion. This may be
+    /// cheaper than inverting the unit quaternion and transforming the
+    /// point.
+    ///
+    /// # Example
+    ///
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{UnitQuaternion, Vector3, Point3};
+    /// let rot = UnitQuaternion::from_axis_angle(&Vector3::y_axis(), f32::consts::FRAC_PI_2);
+    /// let transformed_point = rot.inverse_transform_point(&Point3::new(1.0, 2.0, 3.0));
+    ///
+    /// assert_relative_eq!(transformed_point, Point3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_point(&self, pt: &Point3<N>) -> Point3<N> {
+        // FIXME: would it be useful performancewise not to call inverse explicitly (i-e. implement
+        // the inverse transformation explicitly here) ?
+        self.inverse() * pt
+    }
+
+    /// Rotate a vector by the inverse of this unit quaternion. This may be
+    /// cheaper than inverting the unit quaternion and transforming the
+    /// vector.
+    ///
+    /// # Example
+    ///
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{UnitQuaternion, Vector3};
+    /// let rot = UnitQuaternion::from_axis_angle(&Vector3::y_axis(), f32::consts::FRAC_PI_2);
+    /// let transformed_vector = rot.inverse_transform_vector(&Vector3::new(1.0, 2.0, 3.0));
+    ///
+    /// assert_relative_eq!(transformed_vector, Vector3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_vector(&self, v: &Vector3<N>) -> Vector3<N> {
+        self.inverse() * v
+    }
 }
 
 impl<N: RealField + fmt::Display> fmt::Display for UnitQuaternion<N> {
diff --git a/src/geometry/quaternion_alga.rs b/src/geometry/quaternion_alga.rs
old mode 100644
new mode 100755
index d9864cb1..86eeedb1
--- a/src/geometry/quaternion_alga.rs
+++ b/src/geometry/quaternion_alga.rs
@@ -200,26 +200,24 @@ impl_structures!(
 impl<N: RealField> Transformation<Point3<N>> for UnitQuaternion<N> {
     #[inline]
     fn transform_point(&self, pt: &Point3<N>) -> Point3<N> {
-        self * pt
+        self.transform_point(pt)
     }
 
     #[inline]
     fn transform_vector(&self, v: &Vector3<N>) -> Vector3<N> {
-        self * v
+        self.transform_vector(v)
     }
 }
 
 impl<N: RealField> ProjectiveTransformation<Point3<N>> for UnitQuaternion<N> {
     #[inline]
     fn inverse_transform_point(&self, pt: &Point3<N>) -> Point3<N> {
-        // FIXME: would it be useful performancewise not to call inverse explicitly (i-e. implement
-        // the inverse transformation explicitly here) ?
-        self.inverse() * pt
+        self.inverse_transform_point(pt)
     }
 
     #[inline]
     fn inverse_transform_vector(&self, v: &Vector3<N>) -> Vector3<N> {
-        self.inverse() * v
+        self.inverse_transform_vector(v)
     }
 }
 
diff --git a/src/geometry/quaternion_construction.rs b/src/geometry/quaternion_construction.rs
index 132ede13..7758b4c9 100644
--- a/src/geometry/quaternion_construction.rs
+++ b/src/geometry/quaternion_construction.rs
@@ -13,9 +13,7 @@ use alga::general::RealField;
 
 use crate::base::dimension::U3;
 use crate::base::storage::Storage;
-#[cfg(feature = "arbitrary")]
-use crate::base::Vector3;
-use crate::base::{Unit, Vector, Vector4, Matrix3};
+use crate::base::{Unit, Vector, Vector3, Vector4, Matrix3};
 
 use crate::geometry::{Quaternion, Rotation3, UnitQuaternion};
 
@@ -43,8 +41,13 @@ impl<N: RealField> Quaternion<N> {
     /// ```
     #[inline]
     pub fn new(w: N, i: N, j: N, k: N) -> Self {
-        let v = Vector4::<N>::new(i, j, k, w);
-        Self::from(v)
+        Self::from(Vector4::new(i, j, k, w))
+    }
+
+    /// Constructs a pure quaternion.
+    #[inline]
+    pub fn from_imag(vector: Vector3<N>) -> Self {
+        Self::from_parts(N::zero(), vector)
     }
 
     /// Creates a new quaternion from its scalar and vector parts. Note that the arguments order does
@@ -68,6 +71,12 @@ impl<N: RealField> Quaternion<N> {
         Self::new(scalar, vector[0], vector[1], vector[2])
     }
 
+    /// Constructs a real quaternion.
+    #[inline]
+    pub fn from_real(r: N) -> Self {
+        Self::from_parts(r, Vector3::zero())
+    }
+
     /// Creates a new quaternion from its polar decomposition.
     ///
     /// Note that `axis` is assumed to be a unit vector.
@@ -92,7 +101,7 @@ impl<N: RealField> Quaternion<N> {
     /// ```
     #[inline]
     pub fn identity() -> Self {
-        Self::new(N::one(), N::zero(), N::zero(), N::zero())
+        Self::from_real(N::one())
     }
 }
 
@@ -106,7 +115,7 @@ impl<N: RealField> One for Quaternion<N> {
 impl<N: RealField> Zero for Quaternion<N> {
     #[inline]
     fn zero() -> Self {
-        Self::new(N::zero(), N::zero(), N::zero(), N::zero())
+        Self::from(Vector4::zero())
     }
 
     #[inline]
@@ -579,7 +588,7 @@ impl<N: RealField> UnitQuaternion<N> {
     pub fn new<SB>(axisangle: Vector<N, U3, SB>) -> Self
     where SB: Storage<N, U3> {
         let two: N = crate::convert(2.0f64);
-        let q = Quaternion::<N>::from_parts(N::zero(), axisangle / two).exp();
+        let q = Quaternion::<N>::from_imag(axisangle / two).exp();
         Self::new_unchecked(q)
     }
 
@@ -608,7 +617,7 @@ impl<N: RealField> UnitQuaternion<N> {
     pub fn new_eps<SB>(axisangle: Vector<N, U3, SB>, eps: N) -> Self
     where SB: Storage<N, U3> {
         let two: N = crate::convert(2.0f64);
-        let q = Quaternion::<N>::from_parts(N::zero(), axisangle / two).exp_eps(eps);
+        let q = Quaternion::<N>::from_imag(axisangle / two).exp_eps(eps);
         Self::new_unchecked(q)
     }
 
diff --git a/src/geometry/quaternion_ops.rs b/src/geometry/quaternion_ops.rs
index 5e31364b..7f2b274a 100644
--- a/src/geometry/quaternion_ops.rs
+++ b/src/geometry/quaternion_ops.rs
@@ -552,7 +552,7 @@ impl<N: RealField> Neg for Quaternion<N> {
 
     #[inline]
     fn neg(self) -> Self::Output {
-        Quaternion::from(-self.coords)
+        Self::Output::from(-self.coords)
     }
 }
 
@@ -561,7 +561,7 @@ impl<'a, N: RealField> Neg for &'a Quaternion<N> {
 
     #[inline]
     fn neg(self) -> Self::Output {
-        Quaternion::from(-&self.coords)
+        Self::Output::from(-&self.coords)
     }
 }
 
diff --git a/src/geometry/rotation.rs b/src/geometry/rotation.rs
old mode 100644
new mode 100755
index f5d0ad30..63ba0279
--- a/src/geometry/rotation.rs
+++ b/src/geometry/rotation.rs
@@ -18,7 +18,8 @@ use alga::general::RealField;
 
 use crate::base::allocator::Allocator;
 use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
-use crate::base::{DefaultAllocator, MatrixN, Scalar};
+use crate::base::{DefaultAllocator, MatrixN, Scalar, VectorN};
+use crate::geometry::Point;
 
 /// A rotation matrix.
 #[repr(C)]
@@ -351,6 +352,86 @@ where DefaultAllocator: Allocator<N, D, D>
     }
 }
 
+impl<N: Real, D: DimName> Rotation<N, D>
+where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>
+{
+    /// Rotate the given point.
+    ///
+    /// This is the same as the multiplication `self * pt`.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Point3, Rotation2, Rotation3, UnitQuaternion, Vector3};
+    /// let rot = Rotation3::new(Vector3::y() * f32::consts::FRAC_PI_2);
+    /// let transformed_point = rot.transform_point(&Point3::new(1.0, 2.0, 3.0));
+    ///
+    /// assert_relative_eq!(transformed_point, Point3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        self * pt
+    }
+
+    /// Rotate the given vector.
+    ///
+    /// This is the same as the multiplication `self * v`.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Rotation2, Rotation3, UnitQuaternion, Vector3};
+    /// let rot = Rotation3::new(Vector3::y() * f32::consts::FRAC_PI_2);
+    /// let transformed_vector = rot.transform_vector(&Vector3::new(1.0, 2.0, 3.0));
+    ///
+    /// assert_relative_eq!(transformed_vector, Vector3::new(3.0, 2.0, -1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
+        self * v
+    }
+
+    /// Rotate the given point by the inverse of this rotation. This may be
+    /// cheaper than inverting the rotation and then transforming the given
+    /// point.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Point3, Rotation2, Rotation3, UnitQuaternion, Vector3};
+    /// let rot = Rotation3::new(Vector3::y() * f32::consts::FRAC_PI_2);
+    /// let transformed_point = rot.inverse_transform_point(&Point3::new(1.0, 2.0, 3.0));
+    ///
+    /// assert_relative_eq!(transformed_point, Point3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        Point::from(self.inverse_transform_vector(&pt.coords))
+    }
+
+    /// Rotate the given vector by the inverse of this rotation. This may be
+    /// cheaper than inverting the rotation and then transforming the given
+    /// vector.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Rotation2, Rotation3, UnitQuaternion, Vector3};
+    /// let rot = Rotation3::new(Vector3::y() * f32::consts::FRAC_PI_2);
+    /// let transformed_vector = rot.inverse_transform_vector(&Vector3::new(1.0, 2.0, 3.0));
+    ///
+    /// assert_relative_eq!(transformed_vector, Vector3::new(-3.0, 2.0, 1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
+        self.matrix().tr_mul(v)
+    }
+}
+
 impl<N: Scalar + Eq, D: DimName> Eq for Rotation<N, D> where DefaultAllocator: Allocator<N, D, D> {}
 
 impl<N: Scalar + PartialEq, D: DimName> PartialEq for Rotation<N, D>
diff --git a/src/geometry/rotation_alga.rs b/src/geometry/rotation_alga.rs
old mode 100644
new mode 100755
index 46bdaa8d..e8cf74e7
--- a/src/geometry/rotation_alga.rs
+++ b/src/geometry/rotation_alga.rs
@@ -75,12 +75,12 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>
 {
     #[inline]
     fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        self * pt
+        self.transform_point(pt)
     }
 
     #[inline]
     fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
-        self * v
+        self.transform_vector(v)
     }
 }
 
@@ -89,12 +89,12 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>
 {
     #[inline]
     fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        Point::from(self.inverse_transform_vector(&pt.coords))
+        self.inverse_transform_point(pt)
     }
 
     #[inline]
     fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
-        self.matrix().tr_mul(v)
+        self.inverse_transform_vector(v)
     }
 }
 
diff --git a/src/geometry/similarity.rs b/src/geometry/similarity.rs
old mode 100644
new mode 100755
index 14e50af0..fed04725
--- a/src/geometry/similarity.rs
+++ b/src/geometry/similarity.rs
@@ -16,7 +16,7 @@ use alga::linear::Rotation;
 use crate::base::allocator::Allocator;
 use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
 use crate::base::storage::Owned;
-use crate::base::{DefaultAllocator, MatrixN};
+use crate::base::{DefaultAllocator, MatrixN, VectorN};
 use crate::geometry::{Isometry, Point, Translation};
 
 /// A similarity, i.e., an uniform scaling, followed by a rotation, followed by a translation.
@@ -238,6 +238,87 @@ where
     pub fn append_rotation_wrt_center_mut(&mut self, r: &R) {
         self.isometry.append_rotation_wrt_center_mut(r)
     }
+
+    /// Transform the given point by this similarity.
+    ///
+    /// This is the same as the multiplication `self * pt`.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Point3, Similarity3, Vector3};
+    /// let axisangle = Vector3::y() * f32::consts::FRAC_PI_2;
+    /// let translation = Vector3::new(1.0, 2.0, 3.0);
+    /// let sim = Similarity3::new(translation, axisangle, 3.0);
+    /// let transformed_point = sim.transform_point(&Point3::new(4.0, 5.0, 6.0));
+    /// assert_relative_eq!(transformed_point, Point3::new(19.0, 17.0, -9.0), epsilon = 1.0e-5);
+    /// ```
+    #[inline]
+    pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        self * pt
+    }
+
+    /// Transform the given vector by this similarity, ignoring the translational
+    /// component.
+    ///
+    /// This is the same as the multiplication `self * t`.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Similarity3, Vector3};
+    /// let axisangle = Vector3::y() * f32::consts::FRAC_PI_2;
+    /// let translation = Vector3::new(1.0, 2.0, 3.0);
+    /// let sim = Similarity3::new(translation, axisangle, 3.0);
+    /// let transformed_vector = sim.transform_vector(&Vector3::new(4.0, 5.0, 6.0));
+    /// assert_relative_eq!(transformed_vector, Vector3::new(18.0, 15.0, -12.0), epsilon = 1.0e-5);
+    /// ```
+    #[inline]
+    pub fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
+        self * v
+    }
+
+    /// Transform the given point by the inverse of this similarity. This may
+    /// be cheaper than inverting the similarity and then transforming the
+    /// given point.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Point3, Similarity3, Vector3};
+    /// let axisangle = Vector3::y() * f32::consts::FRAC_PI_2;
+    /// let translation = Vector3::new(1.0, 2.0, 3.0);
+    /// let sim = Similarity3::new(translation, axisangle, 2.0);
+    /// let transformed_point = sim.inverse_transform_point(&Point3::new(4.0, 5.0, 6.0));
+    /// assert_relative_eq!(transformed_point, Point3::new(-1.5, 1.5, 1.5), epsilon = 1.0e-5);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        self.isometry.inverse_transform_point(pt) / self.scaling()
+    }
+
+    /// Transform the given vector by the inverse of this similarity,
+    /// ignoring the translational component. This may be cheaper than
+    /// inverting the similarity and then transforming the given vector.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use std::f32;
+    /// # use nalgebra::{Similarity3, Vector3};
+    /// let axisangle = Vector3::y() * f32::consts::FRAC_PI_2;
+    /// let translation = Vector3::new(1.0, 2.0, 3.0);
+    /// let sim = Similarity3::new(translation, axisangle, 2.0);
+    /// let transformed_vector = sim.inverse_transform_vector(&Vector3::new(4.0, 5.0, 6.0));
+    /// assert_relative_eq!(transformed_vector, Vector3::new(-3.0, 2.5, 2.0), epsilon = 1.0e-5);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
+        self.isometry.inverse_transform_vector(v) / self.scaling()
+    }
 }
 
 // NOTE: we don't require `R: Rotation<...>` here because this is not useful for the implementation
diff --git a/src/geometry/similarity_alga.rs b/src/geometry/similarity_alga.rs
old mode 100644
new mode 100755
index 122ec7b2..448fb133
--- a/src/geometry/similarity_alga.rs
+++ b/src/geometry/similarity_alga.rs
@@ -82,12 +82,12 @@ where
 {
     #[inline]
     fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        self * pt
+        self.transform_point(pt)
     }
 
     #[inline]
     fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
-        self * v
+        self.transform_vector(v)
     }
 }
 
@@ -98,12 +98,12 @@ where
 {
     #[inline]
     fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        self.isometry.inverse_transform_point(pt) / self.scaling()
+        self.inverse_transform_point(pt)
     }
 
     #[inline]
     fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
-        self.isometry.inverse_transform_vector(v) / self.scaling()
+        self.inverse_transform_vector(v)
     }
 }
 
diff --git a/src/geometry/transform.rs b/src/geometry/transform.rs
old mode 100644
new mode 100755
index f5cf92a3..baf3308b
--- a/src/geometry/transform.rs
+++ b/src/geometry/transform.rs
@@ -6,12 +6,14 @@ use std::marker::PhantomData;
 #[cfg(feature = "serde-serialize")]
 use serde::{Deserialize, Deserializer, Serialize, Serializer};
 
-use alga::general::RealField;
+use alga::general::{RealField, TwoSidedInverse};
 
 use crate::base::allocator::Allocator;
 use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
 use crate::base::storage::Owned;
-use crate::base::{DefaultAllocator, MatrixN};
+use crate::base::{DefaultAllocator, MatrixN, VectorN};
+
+use crate::geometry::Point;
 
 /// Trait implemented by phantom types identifying the projective transformation type.
 ///
@@ -452,6 +454,57 @@ where DefaultAllocator: Allocator<N, DimNameSum<D, U1>, DimNameSum<D, U1>>
     }
 }
 
+impl<N, D: DimNameAdd<U1>, C> Transform<N, D, C>
+where
+    N: RealField,
+    C: TCategory,
+    DefaultAllocator: Allocator<N, DimNameSum<D, U1>, DimNameSum<D, U1>>
+        + Allocator<N, DimNameSum<D, U1>>
+        + Allocator<N, D, D>
+        + Allocator<N, D>,
+{
+    /// Transform the given point by this transformation.
+    ///
+    /// This is the same as the multiplication `self * pt`.
+    #[inline]
+    pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        self * pt
+    }
+
+    /// Transform the given vector by this transformation, ignoring the
+    /// translational component of the transformation.
+    ///
+    /// This is the same as the multiplication `self * v`.
+    #[inline]
+    pub fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
+        self * v
+    }
+}
+
+impl<N: RealField, D: DimNameAdd<U1>, C: TCategory> Transform<N, D, C>
+where C: SubTCategoryOf<TProjective>,
+      DefaultAllocator: Allocator<N, DimNameSum<D, U1>, DimNameSum<D, U1>>
+        + Allocator<N, DimNameSum<D, U1>>
+        + Allocator<N, D, D>
+        + Allocator<N, D>,
+{
+    /// Transform the given point by the inverse of this transformation.
+    /// This may be cheaper than inverting the transformation and transforming
+    /// the point.
+    #[inline]
+    pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        self.two_sided_inverse() * pt
+    }
+
+    /// Transform the given vector by the inverse of this transformation.
+    /// This may be cheaper than inverting the transformation and transforming
+    /// the vector.
+    #[inline]
+    pub fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
+        self.two_sided_inverse() * v
+    }
+}
+
 impl<N: RealField, D: DimNameAdd<U1>> Transform<N, D, TGeneral>
 where DefaultAllocator: Allocator<N, DimNameSum<D, U1>, DimNameSum<D, U1>>
 {
diff --git a/src/geometry/transform_alga.rs b/src/geometry/transform_alga.rs
old mode 100644
new mode 100755
index 5d0e8bf8..ec3fd7c6
--- a/src/geometry/transform_alga.rs
+++ b/src/geometry/transform_alga.rs
@@ -96,12 +96,12 @@ where
 {
     #[inline]
     fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        self * pt
+        self.transform_point(pt)
     }
 
     #[inline]
     fn transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
-        self * v
+        self.transform_vector(v)
     }
 }
 
@@ -116,12 +116,12 @@ where
 {
     #[inline]
     fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        self.two_sided_inverse() * pt
+        self.inverse_transform_point(pt)
     }
 
     #[inline]
     fn inverse_transform_vector(&self, v: &VectorN<N, D>) -> VectorN<N, D> {
-        self.two_sided_inverse() * v
+        self.inverse_transform_vector(v)
     }
 }
 
diff --git a/src/geometry/translation.rs b/src/geometry/translation.rs
old mode 100644
new mode 100755
index 43a340f3..e64b3d2e
--- a/src/geometry/translation.rs
+++ b/src/geometry/translation.rs
@@ -11,13 +11,15 @@ use serde::{Deserialize, Deserializer, Serialize, Serializer};
 #[cfg(feature = "abomonation-serialize")]
 use abomonation::Abomonation;
 
-use alga::general::{ClosedNeg, RealField};
+use alga::general::{ClosedAdd, ClosedNeg, ClosedSub, RealField};
 
 use crate::base::allocator::Allocator;
 use crate::base::dimension::{DimName, DimNameAdd, DimNameSum, U1};
 use crate::base::storage::Owned;
 use crate::base::{DefaultAllocator, MatrixN, Scalar, VectorN};
 
+use crate::geometry::Point;
+
 /// A translation.
 #[repr(C)]
 #[derive(Debug)]
@@ -190,6 +192,42 @@ where DefaultAllocator: Allocator<N, D>
     }
 }
 
+impl<N: Scalar + ClosedAdd, D: DimName> Translation<N, D>
+where DefaultAllocator: Allocator<N, D>
+{
+    /// Translate the given point.
+    ///
+    /// This is the same as the multiplication `self * pt`.
+    ///
+    /// # Example
+    /// ```
+    /// # use nalgebra::{Translation3, Point3};
+    /// let t = Translation3::new(1.0, 2.0, 3.0);
+    /// let transformed_point = t.transform_point(&Point3::new(4.0, 5.0, 6.0));
+    /// assert_eq!(transformed_point, Point3::new(5.0, 7.0, 9.0));
+    #[inline]
+    pub fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        pt + &self.vector
+    }
+}
+
+impl<N: Scalar + ClosedSub, D: DimName> Translation<N, D>
+where DefaultAllocator: Allocator<N, D>
+{
+    /// Translate the given point by the inverse of this translation.
+    ///
+    /// # Example
+    /// ```
+    /// # use nalgebra::{Translation3, Point3};
+    /// let t = Translation3::new(1.0, 2.0, 3.0);
+    /// let transformed_point = t.inverse_transform_point(&Point3::new(4.0, 5.0, 6.0));
+    /// assert_eq!(transformed_point, Point3::new(3.0, 3.0, 3.0));
+    #[inline]
+    pub fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
+        pt - &self.vector
+    }
+}
+
 impl<N: Scalar + Eq, D: DimName> Eq for Translation<N, D> where DefaultAllocator: Allocator<N, D> {}
 
 impl<N: Scalar + PartialEq, D: DimName> PartialEq for Translation<N, D>
diff --git a/src/geometry/translation_alga.rs b/src/geometry/translation_alga.rs
old mode 100644
new mode 100755
index f2ced0b3..134790f6
--- a/src/geometry/translation_alga.rs
+++ b/src/geometry/translation_alga.rs
@@ -76,7 +76,7 @@ where DefaultAllocator: Allocator<N, D>
 {
     #[inline]
     fn transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        pt + &self.vector
+        self.transform_point(pt)
     }
 
     #[inline]
@@ -90,7 +90,7 @@ where DefaultAllocator: Allocator<N, D>
 {
     #[inline]
     fn inverse_transform_point(&self, pt: &Point<N, D>) -> Point<N, D> {
-        pt - &self.vector
+        self.inverse_transform_point(pt)
     }
 
     #[inline]
diff --git a/src/geometry/unit_complex.rs b/src/geometry/unit_complex.rs
old mode 100644
new mode 100755
index 2530cc25..7ba7f374
--- a/src/geometry/unit_complex.rs
+++ b/src/geometry/unit_complex.rs
@@ -3,8 +3,8 @@ use num_complex::Complex;
 use std::fmt;
 
 use alga::general::RealField;
-use crate::base::{Matrix2, Matrix3, Unit, Vector1};
-use crate::geometry::Rotation2;
+use crate::base::{Matrix2, Matrix3, Unit, Vector1, Vector2};
+use crate::geometry::{Rotation2, Point2};
 
 /// A complex number with a norm equal to 1.
 pub type UnitComplex<N> = Unit<Complex<N>>;
@@ -251,6 +251,76 @@ impl<N: RealField> UnitComplex<N> {
     pub fn to_homogeneous(&self) -> Matrix3<N> {
         self.to_rotation_matrix().to_homogeneous()
     }
+
+    /// Rotate the given point by this unit complex number.
+    ///
+    /// This is the same as the multiplication `self * pt`.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::{UnitComplex, Point2};
+    /// # use std::f32;
+    /// let rot = UnitComplex::new(f32::consts::FRAC_PI_2);
+    /// let transformed_point = rot.transform_point(&Point2::new(1.0, 2.0));
+    /// assert_relative_eq!(transformed_point, Point2::new(-2.0, 1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn transform_point(&self, pt: &Point2<N>) -> Point2<N> {
+        self * pt
+    }
+
+    /// Rotate the given vector by this unit complex number.
+    ///
+    /// This is the same as the multiplication `self * v`.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::{UnitComplex, Vector2};
+    /// # use std::f32;
+    /// let rot = UnitComplex::new(f32::consts::FRAC_PI_2);
+    /// let transformed_vector = rot.transform_vector(&Vector2::new(1.0, 2.0));
+    /// assert_relative_eq!(transformed_vector, Vector2::new(-2.0, 1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn transform_vector(&self, v: &Vector2<N>) -> Vector2<N> {
+        self * v
+    }
+
+    /// Rotate the given point by the inverse of this unit complex number.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::{UnitComplex, Point2};
+    /// # use std::f32;
+    /// let rot = UnitComplex::new(f32::consts::FRAC_PI_2);
+    /// let transformed_point = rot.inverse_transform_point(&Point2::new(1.0, 2.0));
+    /// assert_relative_eq!(transformed_point, Point2::new(2.0, -1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_point(&self, pt: &Point2<N>) -> Point2<N> {
+        // FIXME: would it be useful performancewise not to call inverse explicitly (i-e. implement
+        // the inverse transformation explicitly here) ?
+        self.inverse() * pt
+    }
+
+    /// Rotate the given vector by the inverse of this unit complex number.
+    ///
+    /// # Example
+    /// ```
+    /// # #[macro_use] extern crate approx;
+    /// # use nalgebra::{UnitComplex, Vector2};
+    /// # use std::f32;
+    /// let rot = UnitComplex::new(f32::consts::FRAC_PI_2);
+    /// let transformed_vector = rot.inverse_transform_vector(&Vector2::new(1.0, 2.0));
+    /// assert_relative_eq!(transformed_vector, Vector2::new(2.0, -1.0), epsilon = 1.0e-6);
+    /// ```
+    #[inline]
+    pub fn inverse_transform_vector(&self, v: &Vector2<N>) -> Vector2<N> {
+        self.inverse() * v
+    }
 }
 
 impl<N: RealField + fmt::Display> fmt::Display for UnitComplex<N> {
diff --git a/src/geometry/unit_complex_alga.rs b/src/geometry/unit_complex_alga.rs
old mode 100644
new mode 100755
index 685a588e..24b55233
--- a/src/geometry/unit_complex_alga.rs
+++ b/src/geometry/unit_complex_alga.rs
@@ -63,12 +63,12 @@ where DefaultAllocator: Allocator<N, U2>
 {
     #[inline]
     fn transform_point(&self, pt: &Point2<N>) -> Point2<N> {
-        self * pt
+        self.transform_point(pt)
     }
 
     #[inline]
     fn transform_vector(&self, v: &Vector2<N>) -> Vector2<N> {
-        self * v
+        self.transform_vector(v)
     }
 }
 
@@ -77,14 +77,12 @@ where DefaultAllocator: Allocator<N, U2>
 {
     #[inline]
     fn inverse_transform_point(&self, pt: &Point2<N>) -> Point2<N> {
-        // FIXME: would it be useful performancewise not to call inverse explicitly (i-e. implement
-        // the inverse transformation explicitly here) ?
-        self.inverse() * pt
+        self.inverse_transform_point(pt)
     }
 
     #[inline]
     fn inverse_transform_vector(&self, v: &Vector2<N>) -> Vector2<N> {
-        self.inverse() * v
+        self.inverse_transform_vector(v)
     }
 }
 
diff --git a/src/linalg/convolution.rs b/src/linalg/convolution.rs
new file mode 100644
index 00000000..b121b34a
--- /dev/null
+++ b/src/linalg/convolution.rs
@@ -0,0 +1,129 @@
+use base::allocator::Allocator;
+use base::default_allocator::DefaultAllocator;
+use base::dimension::{Dim, DimAdd, DimDiff, DimMax, DimMaximum, DimSub, DimSum};
+use std::cmp;
+use storage::Storage;
+use {zero, Real, Vector, VectorN, U1};
+
+impl<N: Real, D1: Dim, S1: Storage<N, D1>> Vector<N, D1, S1> {
+    /// Returns the convolution of the target vector and a kernel
+    ///
+    /// # Arguments
+    ///
+    /// * `kernel` - A Vector with size > 0
+    ///
+    /// # Errors
+    /// Inputs must statisfy `vector.len() >= kernel.len() > 0`.
+    ///
+    pub fn convolve_full<D2, S2>(
+        &self,
+        kernel: Vector<N, D2, S2>,
+    ) -> VectorN<N, DimDiff<DimSum<D1, D2>, U1>>
+    where
+        D1: DimAdd<D2>,
+        D2: DimAdd<D1, Output = DimSum<D1, D2>>,
+        DimSum<D1, D2>: DimSub<U1>,
+        S2: Storage<N, D2>,
+        DefaultAllocator: Allocator<N, DimDiff<DimSum<D1, D2>, U1>>,
+    {
+        let vec = self.len();
+        let ker = kernel.len();
+
+        if ker == 0 || ker > vec {
+            panic!("convolve_full expects `self.len() >= kernel.len() > 0`, received {} and {} respectively.",vec,ker);
+        }
+
+        let result_len = self.data.shape().0.add(kernel.data.shape().0).sub(U1);
+        let mut conv = VectorN::zeros_generic(result_len, U1);
+
+        for i in 0..(vec + ker - 1) {
+            let u_i = if i > vec { i - ker } else { 0 };
+            let u_f = cmp::min(i, vec - 1);
+
+            if u_i == u_f {
+                conv[i] += self[u_i] * kernel[(i - u_i)];
+            } else {
+                for u in u_i..(u_f + 1) {
+                    if i - u < ker {
+                        conv[i] += self[u] * kernel[(i - u)];
+                    }
+                }
+            }
+        }
+        conv
+    }
+    /// Returns the convolution of the target vector and a kernel
+    /// The output convolution consists only of those elements that do not rely on the zero-padding.
+    /// # Arguments
+    ///
+    /// * `kernel` - A Vector with size > 0
+    ///
+    ///
+    /// # Errors
+    /// Inputs must statisfy `self.len() >= kernel.len() > 0`.
+    ///
+    pub fn convolve_valid<D2, S2>(&self, kernel: Vector<N, D2, S2>,
+    ) -> VectorN<N, DimDiff<DimSum<D1, U1>, D2>>
+    where
+        D1: DimAdd<U1>,
+        D2: Dim,
+        DimSum<D1, U1>: DimSub<D2>,
+        S2: Storage<N, D2>,
+        DefaultAllocator: Allocator<N, DimDiff<DimSum<D1, U1>, D2>>,
+    {
+        let vec = self.len();
+        let ker = kernel.len();
+
+        if ker == 0 || ker > vec {
+            panic!("convolve_valid expects `self.len() >= kernel.len() > 0`, received {} and {} respectively.",vec,ker);
+        }
+
+        let result_len = self.data.shape().0.add(U1).sub(kernel.data.shape().0);
+        let mut conv = VectorN::zeros_generic(result_len, U1);
+
+        for i in 0..(vec - ker + 1) {
+            for j in 0..ker {
+                conv[i] += self[i + j] * kernel[ker - j - 1];
+            }
+        }
+        conv
+    }
+
+    /// Returns the convolution of the targetvector and a kernel
+    /// The output convolution is the same size as vector, centered with respect to the ‘full’ output.
+    /// # Arguments
+    ///
+    /// * `kernel` - A Vector with size > 0
+    ///
+    /// # Errors
+    /// Inputs must statisfy `self.len() >= kernel.len() > 0`.
+    pub fn convolve_same<D2, S2>(&self, kernel: Vector<N, D2, S2>) -> VectorN<N, DimMaximum<D1, D2>>
+    where
+        D1: DimMax<D2>,
+        D2: DimMax<D1, Output = DimMaximum<D1, D2>>,
+        S2: Storage<N, D2>,
+        DefaultAllocator: Allocator<N, DimMaximum<D1, D2>>,
+    {
+        let vec = self.len();
+        let ker = kernel.len();
+
+        if ker == 0 || ker > vec {
+            panic!("convolve_same expects `self.len() >= kernel.len() > 0`, received {} and {} respectively.",vec,ker);
+        }
+
+        let result_len = self.data.shape().0.max(kernel.data.shape().0);
+        let mut conv = VectorN::zeros_generic(result_len, U1);
+
+        for i in 0..vec {
+            for j in 0..ker {
+                let val = if i + j < 1 || i + j >= vec + 1 {
+                    zero::<N>()
+                } else {
+                    self[i + j - 1]
+                };
+                conv[i] += val * kernel[ker - j - 1];
+            }
+        }
+        conv
+    }
+}
diff --git a/src/linalg/inverse.rs b/src/linalg/inverse.rs
index 2d3b18ad..94462442 100644
--- a/src/linalg/inverse.rs
+++ b/src/linalg/inverse.rs
@@ -36,7 +36,7 @@ impl<N: ComplexField, D: Dim, S: StorageMut<N, D, D>> SquareMatrix<N, D, S> {
                 0 => true,
                 1 => {
                     let determinant = self.get_unchecked((0, 0)).clone();
-                    if determinant == N::zero() {
+                    if determinant.is_zero() {
                         false
                     } else {
                         *self.get_unchecked_mut((0, 0)) = N::one() / determinant;
@@ -51,7 +51,7 @@ impl<N: ComplexField, D: Dim, S: StorageMut<N, D, D>> SquareMatrix<N, D, S> {
 
                     let determinant = m11 * m22 - m21 * m12;
 
-                    if determinant == N::zero() {
+                    if determinant.is_zero() {
                         false
                     } else {
                         *self.get_unchecked_mut((0, 0)) = m22 / determinant;
@@ -83,7 +83,7 @@ impl<N: ComplexField, D: Dim, S: StorageMut<N, D, D>> SquareMatrix<N, D, S> {
                     let determinant =
                         m11 * minor_m12_m23 - m12 * minor_m11_m23 + m13 * minor_m11_m22;
 
-                    if determinant == N::zero() {
+                    if determinant.is_zero() {
                         false
                     } else {
                         *self.get_unchecked_mut((0, 0)) = minor_m12_m23 / determinant;
diff --git a/src/linalg/mod.rs b/src/linalg/mod.rs
index 4418b283..b6a9e8d8 100644
--- a/src/linalg/mod.rs
+++ b/src/linalg/mod.rs
@@ -17,6 +17,7 @@ mod solve;
 mod svd;
 mod symmetric_eigen;
 mod symmetric_tridiagonal;
+mod convolution;
 
 //// FIXME: Not complete enough for publishing.
 //// This handles only cases where each eigenvalue has multiplicity one.
@@ -33,3 +34,4 @@ pub use self::schur::*;
 pub use self::svd::*;
 pub use self::symmetric_eigen::*;
 pub use self::symmetric_tridiagonal::*;
+pub use self::convolution::*;
diff --git a/tests/linalg/convolution.rs b/tests/linalg/convolution.rs
new file mode 100644
index 00000000..b0d57f72
--- /dev/null
+++ b/tests/linalg/convolution.rs
@@ -0,0 +1,119 @@
+use na::{Vector2,Vector3,Vector4,Vector5,DVector};
+use std::panic;
+
+//
+// Should mimic calculations in Python's scipy library
+// >>>from scipy.signal import convolve
+//
+
+// >>> convolve([1,2,3,4],[1,2],"same")
+// array([ 1,  4,  7, 10])
+#[test]
+fn convolve_same_check(){
+    // Static Tests
+    let actual_s =  Vector4::from_vec(vec![1.0,4.0,7.0,10.0]);
+    let expected_s = Vector4::new(1.0,2.0,3.0,4.0).convolve_same(Vector2::new(1.0,2.0));
+
+    assert!(relative_eq!(actual_s, expected_s, epsilon = 1.0e-7));
+
+    // Dynamic Tests
+    let actual_d =  DVector::from_vec(vec![1.0,4.0,7.0,10.0]);
+    let expected_d = DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_same(DVector::from_vec(vec![1.0,2.0]));
+
+    assert!(relative_eq!(actual_d, expected_d, epsilon = 1.0e-7));
+
+    // Panic Tests
+    // These really only apply to dynamic sized vectors
+    assert!(
+        panic::catch_unwind(|| {
+            DVector::from_vec(vec![1.0,2.0]).convolve_same(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
+        }).is_err()
+    );
+
+    assert!(
+        panic::catch_unwind(|| {
+            DVector::<f32>::from_vec(vec![]).convolve_same(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
+        }).is_err()
+    );
+
+    assert!(
+        panic::catch_unwind(|| {
+            DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_same(DVector::<f32>::from_vec(vec![]));
+        }).is_err()
+    );
+}
+
+// >>> convolve([1,2,3,4],[1,2],"full")
+// array([ 1, 4,  7, 10, 8])
+#[test]
+fn convolve_full_check(){
+    // Static Tests
+    let actual_s =  Vector5::new(1.0,4.0,7.0,10.0,8.0);
+    let expected_s = Vector4::new(1.0,2.0,3.0,4.0).convolve_full(Vector2::new(1.0,2.0));
+
+    assert!(relative_eq!(actual_s, expected_s, epsilon = 1.0e-7));
+
+    // Dynamic Tests
+    let actual_d =  DVector::from_vec(vec![1.0,4.0,7.0,10.0,8.0]);
+    let expected_d = DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_full(DVector::from_vec(vec![1.0,2.0]));
+
+    assert!(relative_eq!(actual_d, expected_d, epsilon = 1.0e-7));
+
+    // Panic Tests
+    // These really only apply to dynamic sized vectors
+    assert!(
+        panic::catch_unwind(|| {
+            DVector::from_vec(vec![1.0,2.0]).convolve_full(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
+        }).is_err()
+    );
+
+    assert!(
+        panic::catch_unwind(|| {
+            DVector::<f32>::from_vec(vec![]).convolve_full(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
+        }).is_err()
+    );
+
+    assert!(
+        panic::catch_unwind(|| {
+            DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_full(DVector::<f32>::from_vec(vec![]));
+        }).is_err()
+    );
+}
+
+// >>> convolve([1,2,3,4],[1,2],"valid")
+// array([ 4,  7, 10])
+#[test]
+fn convolve_valid_check(){
+    // Static Tests
+    let actual_s =  Vector3::from_vec(vec![4.0,7.0,10.0]);
+    let expected_s = Vector4::new(1.0,2.0,3.0,4.0).convolve_valid( Vector2::new(1.0,2.0));
+
+    assert!(relative_eq!(actual_s, expected_s, epsilon = 1.0e-7));
+
+    // Dynamic Tests
+    let actual_d =  DVector::from_vec(vec![4.0,7.0,10.0]);
+    let expected_d = DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_valid(DVector::from_vec(vec![1.0,2.0]));
+
+    assert!(relative_eq!(actual_d, expected_d, epsilon = 1.0e-7));
+
+    // Panic Tests
+    // These really only apply to dynamic sized vectors
+    assert!(
+        panic::catch_unwind(|| {
+            DVector::from_vec(vec![1.0,2.0]).convolve_valid(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
+        }).is_err()
+    );
+
+    assert!(
+        panic::catch_unwind(|| {
+            DVector::<f32>::from_vec(vec![]).convolve_valid(DVector::from_vec(vec![1.0,2.0,3.0,4.0]));
+        }).is_err()
+    );
+
+    assert!(
+        panic::catch_unwind(|| {
+            DVector::from_vec(vec![1.0,2.0,3.0,4.0]).convolve_valid(DVector::<f32>::from_vec(vec![]));
+        }).is_err()
+    );
+
+}
\ No newline at end of file
diff --git a/tests/linalg/mod.rs b/tests/linalg/mod.rs
index f515b4d6..e881d999 100644
--- a/tests/linalg/mod.rs
+++ b/tests/linalg/mod.rs
@@ -11,3 +11,4 @@ mod schur;
 mod solve;
 mod svd;
 mod tridiagonal;
+mod convolution;
\ No newline at end of file