diff --git a/CHANGELOG.md b/CHANGELOG.md
index d8eacebb..88eb2f2d 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -7,6 +7,8 @@ This project adheres to [Semantic Versioning](http://semver.org/).
 ## [0.15.0] - WIP
 ### Modified
 ### Added
+  * Add `UnitQuaternion` constructor `::new_eps(...)` and `::from_scaled_axis_eps(...)` that return the
+    identity if the magnitude of the input axisangle is smaller than the epsilon provided.
   * Add methods `.rotation_between_axis(...)` and `.scaled_rotation_between_axis(...)` to `UnitComplex`
     to compute the rotation matrix between two 2D **unit** vectors.
   * Add methods `.axis_angle()` to `UnitComplex` and `UnitQuaternion` in order to retrieve both the
diff --git a/src/geometry/quaternion.rs b/src/geometry/quaternion.rs
index 88ad945d..b1cba752 100644
--- a/src/geometry/quaternion.rs
+++ b/src/geometry/quaternion.rs
@@ -201,10 +201,16 @@ impl<N: Real> Quaternion<N> {
     /// Compute the exponential of a quaternion.
     #[inline]
     pub fn exp(&self) -> Quaternion<N> {
+        self.exp_eps(N::default_epsilon())
+    }
+
+    /// Compute the exponential of a quaternion.
+    #[inline]
+    pub fn exp_eps(&self, eps: N) -> Quaternion<N> {
         let v = self.vector();
         let nn = v.norm_squared();
 
-        if relative_eq!(nn, N::zero()) {
+        if nn <= eps * eps {
             Quaternion::identity()
         } else {
             let w_exp = self.scalar().exp();
diff --git a/src/geometry/quaternion_construction.rs b/src/geometry/quaternion_construction.rs
index 19f3ff98..63aea8ea 100644
--- a/src/geometry/quaternion_construction.rs
+++ b/src/geometry/quaternion_construction.rs
@@ -351,7 +351,7 @@ impl<N: Real> UnitQuaternion<N> {
 
     /// Creates a new unit quaternion rotation from a rotation axis scaled by the rotation angle.
     ///
-    /// If `axisangle` is zero, this returns the indentity rotation.
+    /// If `axisangle` has a magnitude smaller than `N::default_epsilon()`, this returns the indentity rotation.
     #[inline]
     pub fn new<SB>(axisangle: Vector<N, U3, SB>) -> Self
     where
@@ -364,7 +364,20 @@ impl<N: Real> UnitQuaternion<N> {
 
     /// Creates a new unit quaternion rotation from a rotation axis scaled by the rotation angle.
     ///
-    /// If `axisangle` is zero, this returns the indentity rotation.
+    /// If `axisangle` has a magnitude smaller than `eps`, this returns the indentity rotation.
+    #[inline]
+    pub fn new_eps<SB>(axisangle: Vector<N, U3, SB>, eps: N) -> Self
+    where
+        SB: Storage<N, U3>,
+    {
+        let two: N = ::convert(2.0f64);
+        let q = Quaternion::<N>::from_parts(N::zero(), axisangle / two).exp_eps(eps);
+        Self::new_unchecked(q)
+    }
+
+    /// Creates a new unit quaternion rotation from a rotation axis scaled by the rotation angle.
+    ///
+    /// If `axisangle` has a magnitude smalle than `N::default_epsilon()`, this returns the indentity rotation.
     /// Same as `Self::new(axisangle)`.
     #[inline]
     pub fn from_scaled_axis<SB>(axisangle: Vector<N, U3, SB>) -> Self
@@ -373,6 +386,18 @@ impl<N: Real> UnitQuaternion<N> {
     {
         Self::new(axisangle)
     }
+
+    /// Creates a new unit quaternion rotation from a rotation axis scaled by the rotation angle.
+    ///
+    /// If `axisangle` has a mangnitude smaller than `eps`, this returns the indentity rotation.
+    /// Same as `Self::new(axisangle)`.
+    #[inline]
+    pub fn from_scaled_axis_eps<SB>(axisangle: Vector<N, U3, SB>, eps: N) -> Self
+    where
+        SB: Storage<N, U3>,
+    {
+        Self::new_eps(axisangle, eps)
+    }
 }
 
 impl<N: Real> One for UnitQuaternion<N> {