diff --git a/benches/geometry/mod.rs b/benches/geometry/mod.rs
index d0d0e948..cb23cbf9 100644
--- a/benches/geometry/mod.rs
+++ b/benches/geometry/mod.rs
@@ -1,3 +1,3 @@
-pub use quaternion::quaternion;
+pub use self::quaternion::quaternion;
 
 mod quaternion;
diff --git a/examples/identity.rs b/examples/identity.rs
index 06d69f70..c20c5616 100644
--- a/examples/identity.rs
+++ b/examples/identity.rs
@@ -36,4 +36,4 @@ fn main() {
     // They both return the same result.
     assert!(result1 == Vector3::new(100001.0, 200002.0, 300003.0));
     assert!(result2 == Vector3::new(100001.0, 200002.0, 300003.0));
-}
+}
\ No newline at end of file
diff --git a/src/base/blas.rs b/src/base/blas.rs
index cdfe4f9f..55a20216 100644
--- a/src/base/blas.rs
+++ b/src/base/blas.rs
@@ -1,4 +1,4 @@
-use alga::general::{ClosedAdd, ClosedMul, Complex};
+use alga::general::{ClosedAdd, ClosedMul, ComplexField};
 #[cfg(feature = "std")]
 use matrixmultiply;
 use num::{One, Signed, Zero};
@@ -15,7 +15,7 @@ use crate::base::{DefaultAllocator, Matrix, Scalar, SquareMatrix, Vector, DVecto
 
 
 // FIXME: find a way to avoid code duplication just for complex number support.
-impl<N: Complex, D: Dim, S: Storage<N, D>> Vector<N, D, S> {
+impl<N: ComplexField, D: Dim, S: Storage<N, D>> Vector<N, D, S> {
     /// Computes the index of the vector component with the largest complex or real absolute value.
     ///
     /// # Examples:
@@ -193,7 +193,7 @@ impl<N: Scalar + PartialOrd, D: Dim, S: Storage<N, D>> Vector<N, D, S> {
 }
 
 // FIXME: find a way to avoid code duplication just for complex number support.
-impl<N: Complex, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
+impl<N: ComplexField, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
     /// Computes the index of the matrix component with the largest absolute value.
     ///
     /// # Examples:
@@ -416,11 +416,11 @@ where N: Scalar + Zero + ClosedAdd + ClosedMul
     #[inline]
     pub fn dotc<R2: Dim, C2: Dim, SB>(&self, rhs: &Matrix<N, R2, C2, SB>) -> N
         where
-            N: Complex,
+            N: ComplexField,
             SB: Storage<N, R2, C2>,
             ShapeConstraint: DimEq<R, R2> + DimEq<C, C2>,
     {
-        self.dotx(rhs, Complex::conjugate)
+        self.dotx(rhs, ComplexField::conjugate)
     }
 
     /// The dot product between the transpose of `self` and `rhs`.
@@ -726,7 +726,7 @@ where
         x: &Vector<N, D3, SC>,
         beta: N,
     ) where
-        N: Complex,
+        N: ComplexField,
         SB: Storage<N, D2, D2>,
         SC: Storage<N, D3>,
         ShapeConstraint: DimEq<D, D2> + AreMultipliable<D2, D2, D3, U1>,
@@ -881,12 +881,12 @@ where N: Scalar + Zero + ClosedAdd + ClosedMul
         y: &Vector<N, D3, SC>,
         beta: N,
     ) where
-        N: Complex,
+        N: ComplexField,
         SB: Storage<N, D2>,
         SC: Storage<N, D3>,
         ShapeConstraint: DimEq<R1, D2> + DimEq<C1, D3>,
     {
-        self.gerx(alpha, x, y, beta, Complex::conjugate)
+        self.gerx(alpha, x, y, beta, ComplexField::conjugate)
     }
 
     /// Computes `self = alpha * a * b + beta * self`, where `a, b, self` are matrices.
@@ -1216,12 +1216,12 @@ where N: Scalar + Zero + ClosedAdd + ClosedMul
         y: &Vector<N, D3, SC>,
         beta: N,
     ) where
-        N: Complex,
+        N: ComplexField,
         SB: Storage<N, D2>,
         SC: Storage<N, D3>,
         ShapeConstraint: DimEq<R1, D2> + DimEq<C1, D3>,
     {
-        self.xxgerx(alpha, x, y, beta, Complex::conjugate)
+        self.xxgerx(alpha, x, y, beta, ComplexField::conjugate)
     }
 }
 
diff --git a/src/base/matrix.rs b/src/base/matrix.rs
index e5060013..ad188a26 100644
--- a/src/base/matrix.rs
+++ b/src/base/matrix.rs
@@ -16,7 +16,7 @@ use serde::{Deserialize, Deserializer, Serialize, Serializer};
 #[cfg(feature = "abomonation-serialize")]
 use abomonation::Abomonation;
 
-use alga::general::{ClosedAdd, ClosedMul, ClosedSub, Real, Ring, Complex, Field};
+use alga::general::{ClosedAdd, ClosedMul, ClosedSub, Real, Ring, ComplexField, Field};
 
 use crate::base::allocator::{Allocator, SameShapeAllocator, SameShapeC, SameShapeR};
 use crate::base::constraint::{DimEq, SameNumberOfColumns, SameNumberOfRows, ShapeConstraint};
@@ -913,7 +913,7 @@ impl<N: Scalar, D: Dim, S: StorageMut<N, D, D>> Matrix<N, D, D, S> {
     }
 }
 
-impl<N: Complex, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
+impl<N: ComplexField, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
     /// Takes the adjoint (aka. conjugate-transpose) of `self` and store the result into `out`.
     #[inline]
     pub fn adjoint_to<R2, C2, SB>(&self, out: &mut Matrix<N, R2, C2, SB>)
@@ -996,7 +996,7 @@ impl<N: Complex, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
     }
 }
 
-impl<N: Complex, R: Dim, C: Dim, S: StorageMut<N, R, C>> Matrix<N, R, C, S> {
+impl<N: ComplexField, R: Dim, C: Dim, S: StorageMut<N, R, C>> Matrix<N, R, C, S> {
     /// The conjugate of the complex matrix `self` computed in-place.
     #[inline]
     pub fn conjugate_mut(&mut self) {
@@ -1016,7 +1016,7 @@ impl<N: Complex, R: Dim, C: Dim, S: StorageMut<N, R, C>> Matrix<N, R, C, S> {
     }
 }
 
-impl<N: Complex, D: Dim, S: StorageMut<N, D, D>> Matrix<N, D, D, S> {
+impl<N: ComplexField, D: Dim, S: StorageMut<N, D, D>> Matrix<N, D, D, S> {
     /// Sets `self` to its adjoint.
     #[deprecated(note = "Renamed to `self.adjoint_mut()`.")]
     pub fn conjugate_transform_mut(&mut self) {
@@ -1103,7 +1103,7 @@ impl<N: Scalar, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
     }
 }
 
-impl<N: Complex, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
+impl<N: ComplexField, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
     /// The symmetric part of `self`, i.e., `0.5 * (self + self.transpose())`.
     #[inline]
     pub fn symmetric_part(&self) -> MatrixMN<N, D, D>
@@ -1541,7 +1541,7 @@ where DefaultAllocator: Allocator<N, U3>
     }
 }
 
-impl<N: Complex, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
+impl<N: ComplexField, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
     /// The smallest angle between two vectors.
     #[inline]
     pub fn angle<R2: Dim, C2: Dim, SB>(&self, other: &Matrix<N, R2, C2, SB>) -> N::Real
@@ -1592,7 +1592,7 @@ impl<N: Scalar + Zero + One + ClosedAdd + ClosedSub + ClosedMul, D: Dim, S: Stor
     }
 }
 
-impl<N: Complex, D: Dim, S: Storage<N, D>> Unit<Vector<N, D, S>> {
+impl<N: ComplexField, D: Dim, S: Storage<N, D>> Unit<Vector<N, D, S>> {
     /// Computes the spherical linear interpolation between two unit vectors.
     pub fn slerp<S2: Storage<N, D>>(
         &self,
diff --git a/src/base/matrix_alga.rs b/src/base/matrix_alga.rs
index 87ab584c..790d60d3 100644
--- a/src/base/matrix_alga.rs
+++ b/src/base/matrix_alga.rs
@@ -7,7 +7,7 @@ use alga::general::{
     AbstractGroup, AbstractGroupAbelian, AbstractLoop, AbstractMagma, AbstractModule,
     AbstractMonoid, AbstractQuasigroup, AbstractSemigroup, Additive, ClosedAdd, ClosedMul,
     ClosedNeg, Field, Identity, TwoSidedInverse, JoinSemilattice, Lattice, MeetSemilattice, Module,
-    Multiplicative, RingCommutative, Complex
+    Multiplicative, RingCommutative, ComplexField
 };
 use alga::linear::{
     FiniteDimInnerSpace, FiniteDimVectorSpace, InnerSpace, NormedSpace, VectorSpace,
@@ -145,11 +145,11 @@ where
     }
 }
 
-impl<N: Complex, R: DimName, C: DimName> NormedSpace for MatrixMN<N, R, C>
+impl<N: ComplexField, R: DimName, C: DimName> NormedSpace for MatrixMN<N, R, C>
 where DefaultAllocator: Allocator<N, R, C>
 {
     type Real = N::Real;
-    type Complex = N;
+    type ComplexField = N;
 
     #[inline]
     fn norm_squared(&self) -> N::Real {
@@ -182,7 +182,7 @@ where DefaultAllocator: Allocator<N, R, C>
     }
 }
 
-impl<N: Complex, R: DimName, C: DimName> InnerSpace for MatrixMN<N, R, C>
+impl<N: ComplexField, R: DimName, C: DimName> InnerSpace for MatrixMN<N, R, C>
 where DefaultAllocator: Allocator<N, R, C>
 {
     #[inline]
@@ -200,7 +200,7 @@ where DefaultAllocator: Allocator<N, R, C>
 // In particular:
 //   − use `x()` instead of `::canonical_basis_element`
 //   − use `::new(x, y, z)` instead of `::from_slice`
-impl<N: Complex, R: DimName, C: DimName> FiniteDimInnerSpace for MatrixMN<N, R, C>
+impl<N: ComplexField, R: DimName, C: DimName> FiniteDimInnerSpace for MatrixMN<N, R, C>
 where DefaultAllocator: Allocator<N, R, C>
 {
     #[inline]
diff --git a/src/base/norm.rs b/src/base/norm.rs
index 91958cec..06bca851 100644
--- a/src/base/norm.rs
+++ b/src/base/norm.rs
@@ -1,7 +1,7 @@
 use num::Zero;
 
 use crate::allocator::Allocator;
-use crate::{Real, Complex};
+use crate::{Real, ComplexField};
 use crate::storage::{Storage, StorageMut};
 use crate::base::{DefaultAllocator, Matrix, Dim, MatrixMN};
 use crate::constraint::{SameNumberOfRows, SameNumberOfColumns, ShapeConstraint};
@@ -11,7 +11,7 @@ use crate::constraint::{SameNumberOfRows, SameNumberOfColumns, ShapeConstraint};
 /// A trait for abstract matrix norms.
 ///
 /// This may be moved to the alga crate in the future.
-pub trait Norm<N: Complex> {
+pub trait Norm<N: ComplexField> {
     /// Apply this norm to the given matrix.
     fn norm<R, C, S>(&self, m: &Matrix<N, R, C, S>) -> N::Real
         where R: Dim, C: Dim, S: Storage<N, R, C>;
@@ -29,7 +29,7 @@ pub struct LpNorm(pub i32);
 /// L-infinite norm aka. Chebytchev norm aka. uniform norm aka. suppremum norm.
 pub struct UniformNorm;
 
-impl<N: Complex> Norm<N> for EuclideanNorm {
+impl<N: ComplexField> Norm<N> for EuclideanNorm {
     #[inline]
     fn norm<R, C, S>(&self, m: &Matrix<N, R, C, S>) -> N::Real
         where R: Dim, C: Dim, S: Storage<N, R, C> {
@@ -48,7 +48,7 @@ impl<N: Complex> Norm<N> for EuclideanNorm {
     }
 }
 
-impl<N: Complex> Norm<N> for LpNorm {
+impl<N: ComplexField> Norm<N> for LpNorm {
     #[inline]
     fn norm<R, C, S>(&self, m: &Matrix<N, R, C, S>) -> N::Real
         where R: Dim, C: Dim, S: Storage<N, R, C> {
@@ -69,7 +69,7 @@ impl<N: Complex> Norm<N> for LpNorm {
     }
 }
 
-impl<N: Complex> Norm<N> for UniformNorm {
+impl<N: ComplexField> Norm<N> for UniformNorm {
     #[inline]
     fn norm<R, C, S>(&self, m: &Matrix<N, R, C, S>) -> N::Real
         where R: Dim, C: Dim, S: Storage<N, R, C> {
@@ -95,7 +95,7 @@ impl<N: Complex> Norm<N> for UniformNorm {
 }
 
 
-impl<N: Complex, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
+impl<N: ComplexField, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
     /// The squared L2 norm of this vector.
     #[inline]
     pub fn norm_squared(&self) -> N::Real {
@@ -213,7 +213,7 @@ impl<N: Complex, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
 }
 
 
-impl<N: Complex, R: Dim, C: Dim, S: StorageMut<N, R, C>> Matrix<N, R, C, S> {
+impl<N: ComplexField, R: Dim, C: Dim, S: StorageMut<N, R, C>> Matrix<N, R, C, S> {
     /// Normalizes this matrix in-place and returns its norm.
     #[inline]
     pub fn normalize_mut(&mut self) -> N::Real {
diff --git a/src/base/ops.rs b/src/base/ops.rs
index bfc71040..8a6fe5a1 100644
--- a/src/base/ops.rs
+++ b/src/base/ops.rs
@@ -5,7 +5,7 @@ use std::ops::{
     Add, AddAssign, Div, DivAssign, Index, IndexMut, Mul, MulAssign, Neg, Sub, SubAssign,
 };
 
-use alga::general::{Complex, ClosedAdd, ClosedDiv, ClosedMul, ClosedNeg, ClosedSub};
+use alga::general::{ComplexField, ClosedAdd, ClosedDiv, ClosedMul, ClosedNeg, ClosedSub};
 
 use crate::base::allocator::{Allocator, SameShapeAllocator, SameShapeC, SameShapeR};
 use crate::base::constraint::{
@@ -633,7 +633,7 @@ where
     #[inline]
     pub fn ad_mul<R2: Dim, C2: Dim, SB>(&self, rhs: &Matrix<N, R2, C2, SB>) -> MatrixMN<N, C1, C2>
         where
-            N: Complex,
+            N: ComplexField,
             SB: Storage<N, R2, C2>,
             DefaultAllocator: Allocator<N, C1, C2>,
             ShapeConstraint: SameNumberOfRows<R1, R2>,
@@ -700,7 +700,7 @@ where
         rhs: &Matrix<N, R2, C2, SB>,
         out: &mut Matrix<N, R3, C3, SC>,
     ) where
-        N: Complex,
+        N: ComplexField,
         SB: Storage<N, R2, C2>,
         SC: StorageMut<N, R3, C3>,
         ShapeConstraint: SameNumberOfRows<R1, R2> + DimEq<C1, R3> + DimEq<C2, C3>,
@@ -833,7 +833,7 @@ impl<N: Scalar, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
     /// Returns the the 1-norm of the complex component with the largest 1-norm.
     #[inline]
     pub fn camax(&self) -> N::Real
-        where N: Complex {
+        where N: ComplexField {
         self.xcmp(|e| e.norm1(), |a, b| a > b)
     }
 
@@ -854,7 +854,7 @@ impl<N: Scalar, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
     /// Returns the the 1-norm of the complex component with the smallest 1-norm.
     #[inline]
     pub fn camin(&self) -> N::Real
-        where N: Complex {
+        where N: ComplexField {
         self.xcmp(|e| e.norm1(), |a, b| a < b)
     }
 
diff --git a/src/base/properties.rs b/src/base/properties.rs
index 73019209..acb7dfb4 100644
--- a/src/base/properties.rs
+++ b/src/base/properties.rs
@@ -2,7 +2,7 @@
 use approx::RelativeEq;
 use num::{One, Zero};
 
-use alga::general::{ClosedAdd, ClosedMul, Real, Complex};
+use alga::general::{ClosedAdd, ClosedMul, Real, ComplexField};
 
 use crate::base::allocator::Allocator;
 use crate::base::dimension::{Dim, DimMin};
@@ -84,7 +84,7 @@ impl<N: Scalar, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
     }
 }
 
-impl<N: Complex, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
+impl<N: ComplexField, R: Dim, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S> {
     /// Checks that `Mᵀ × M = Id`.
     ///
     /// In this definition `Id` is approximately equal to the identity matrix with a relative error
diff --git a/src/base/unit.rs b/src/base/unit.rs
index 90b343e1..60aa1ae2 100644
--- a/src/base/unit.rs
+++ b/src/base/unit.rs
@@ -10,7 +10,7 @@ use serde::{Deserialize, Deserializer, Serialize, Serializer};
 #[cfg(feature = "abomonation-serialize")]
 use abomonation::Abomonation;
 
-use alga::general::{SubsetOf, Complex};
+use alga::general::{SubsetOf, ComplexField};
 use alga::linear::NormedSpace;
 
 /// A wrapper that ensures the underlying algebraic entity has a unit norm.
@@ -106,7 +106,7 @@ impl<T: NormedSpace> Unit<T> {
         let sq_norm = self.value.norm_squared();
         let _3: T::Real = crate::convert(3.0);
         let _0_5: T::Real = crate::convert(0.5);
-        self.value *= T::Complex::from_real(_0_5 * (_3 - sq_norm));
+        self.value *= T::ComplexField::from_real(_0_5 * (_3 - sq_norm));
     }
 }
 
diff --git a/src/debug/random_orthogonal.rs b/src/debug/random_orthogonal.rs
index b14c7da4..421b041a 100644
--- a/src/debug/random_orthogonal.rs
+++ b/src/debug/random_orthogonal.rs
@@ -3,7 +3,7 @@ use crate::base::storage::Owned;
 #[cfg(feature = "arbitrary")]
 use quickcheck::{Arbitrary, Gen};
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::base::Scalar;
 use crate::base::allocator::Allocator;
 use crate::base::dimension::{Dim, Dynamic, U2};
@@ -18,7 +18,7 @@ where DefaultAllocator: Allocator<N, D, D>
     m: MatrixN<N, D>,
 }
 
-impl<N: Complex, D: Dim> RandomOrthogonal<N, D>
+impl<N: ComplexField, D: Dim> RandomOrthogonal<N, D>
 where DefaultAllocator: Allocator<N, D, D>
 {
     /// Retrieve the generated matrix.
@@ -41,7 +41,7 @@ where DefaultAllocator: Allocator<N, D, D>
 }
 
 #[cfg(feature = "arbitrary")]
-impl<N: Complex + Arbitrary + Send, D: Dim> Arbitrary for RandomOrthogonal<N, D>
+impl<N: ComplexField + Arbitrary + Send, D: Dim> Arbitrary for RandomOrthogonal<N, D>
 where
     DefaultAllocator: Allocator<N, D, D>,
     Owned<N, D, D>: Clone + Send,
diff --git a/src/debug/random_sdp.rs b/src/debug/random_sdp.rs
index ca3f122f..47e3ca60 100644
--- a/src/debug/random_sdp.rs
+++ b/src/debug/random_sdp.rs
@@ -3,7 +3,7 @@ use crate::base::storage::Owned;
 #[cfg(feature = "arbitrary")]
 use quickcheck::{Arbitrary, Gen};
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::base::Scalar;
 use crate::base::allocator::Allocator;
 use crate::base::dimension::{Dim, Dynamic};
@@ -19,7 +19,7 @@ where DefaultAllocator: Allocator<N, D, D>
     m: MatrixN<N, D>,
 }
 
-impl<N: Complex, D: Dim> RandomSDP<N, D>
+impl<N: ComplexField, D: Dim> RandomSDP<N, D>
 where DefaultAllocator: Allocator<N, D, D>
 {
     /// Retrieve the generated matrix.
@@ -44,7 +44,7 @@ where DefaultAllocator: Allocator<N, D, D>
 }
 
 #[cfg(feature = "arbitrary")]
-impl<N: Complex + Arbitrary + Send, D: Dim> Arbitrary for RandomSDP<N, D>
+impl<N: ComplexField + Arbitrary + Send, D: Dim> Arbitrary for RandomSDP<N, D>
 where
     DefaultAllocator: Allocator<N, D, D>,
     Owned<N, D, D>: Clone + Send,
diff --git a/src/geometry/quaternion_alga.rs b/src/geometry/quaternion_alga.rs
index b1c9a1ba..a73beff7 100644
--- a/src/geometry/quaternion_alga.rs
+++ b/src/geometry/quaternion_alga.rs
@@ -119,7 +119,7 @@ impl<N: Real> FiniteDimVectorSpace for Quaternion<N> {
 
 impl<N: Real> NormedSpace for Quaternion<N> {
     type Real = N;
-    type Complex = N;
+    type ComplexField = N;
 
     #[inline]
     fn norm_squared(&self) -> N {
diff --git a/src/geometry/reflection.rs b/src/geometry/reflection.rs
index a13e235a..b36924d2 100644
--- a/src/geometry/reflection.rs
+++ b/src/geometry/reflection.rs
@@ -1,4 +1,4 @@
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::base::allocator::Allocator;
 use crate::base::constraint::{AreMultipliable, DimEq, SameNumberOfRows, ShapeConstraint};
 use crate::base::{DefaultAllocator, Matrix, Scalar, Unit, Vector};
@@ -13,7 +13,7 @@ pub struct Reflection<N: Scalar, D: Dim, S: Storage<N, D>> {
     bias: N,
 }
 
-impl<N: Complex, D: Dim, S: Storage<N, D>> Reflection<N, D, S> {
+impl<N: ComplexField, D: Dim, S: Storage<N, D>> Reflection<N, D, S> {
     /// Creates a new reflection wrt the plane orthogonal to the given axis and bias.
     ///
     /// The bias is the position of the plane on the axis. In particular, a bias equal to zero
diff --git a/src/lib.rs b/src/lib.rs
index feae7b22..f66ae4f0 100644
--- a/src/lib.rs
+++ b/src/lib.rs
@@ -160,7 +160,7 @@ use alga::linear::SquareMatrix as AlgaSquareMatrix;
 use alga::linear::{EuclideanSpace, FiniteDimVectorSpace, InnerSpace, NormedSpace};
 use num::Signed;
 
-pub use alga::general::{Id, Real, Complex};
+pub use alga::general::{Id, Real, ComplexField};
 
 /*
  *
diff --git a/src/linalg/bidiagonal.rs b/src/linalg/bidiagonal.rs
index 09031e28..e61c98ec 100644
--- a/src/linalg/bidiagonal.rs
+++ b/src/linalg/bidiagonal.rs
@@ -1,7 +1,7 @@
 #[cfg(feature = "serde-serialize")]
 use serde::{Deserialize, Serialize};
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::allocator::Allocator;
 use crate::base::{DefaultAllocator, Matrix, MatrixMN, MatrixN, Unit, VectorN};
 use crate::dimension::{Dim, DimDiff, DimMin, DimMinimum, DimSub, U1};
@@ -37,7 +37,7 @@ use crate::linalg::householder;
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct Bidiagonal<N: Complex, R: DimMin<C>, C: Dim>
+pub struct Bidiagonal<N: ComplexField, R: DimMin<C>, C: Dim>
 where
     DimMinimum<R, C>: DimSub<U1>,
     DefaultAllocator: Allocator<N, R, C>
@@ -54,7 +54,7 @@ where
     upper_diagonal: bool,
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim> Copy for Bidiagonal<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> Copy for Bidiagonal<N, R, C>
 where
     DimMinimum<R, C>: DimSub<U1>,
     DefaultAllocator: Allocator<N, R, C>
@@ -65,7 +65,7 @@ where
     VectorN<N, DimDiff<DimMinimum<R, C>, U1>>: Copy,
 {}
 
-impl<N: Complex, R: DimMin<C>, C: Dim> Bidiagonal<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> Bidiagonal<N, R, C>
 where
     DimMinimum<R, C>: DimSub<U1>,
     DefaultAllocator: Allocator<N, R, C>
@@ -282,7 +282,7 @@ where
     }
 }
 
-// impl<N: Complex, D: DimMin<D, Output = D> + DimSub<Dynamic>> Bidiagonal<N, D, D>
+// impl<N: ComplexField, D: DimMin<D, Output = D> + DimSub<Dynamic>> Bidiagonal<N, D, D>
 //     where DefaultAllocator: Allocator<N, D, D> +
 //                             Allocator<N, D> {
 //     /// Solves the linear system `self * x = b`, where `x` is the unknown to be determined.
@@ -355,7 +355,7 @@ where
 //     // }
 // }
 
-impl<N: Complex, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
+impl<N: ComplexField, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
 where
     DimMinimum<R, C>: DimSub<U1>,
     DefaultAllocator: Allocator<N, R, C>
diff --git a/src/linalg/cholesky.rs b/src/linalg/cholesky.rs
index 985abddc..0b6e6db5 100644
--- a/src/linalg/cholesky.rs
+++ b/src/linalg/cholesky.rs
@@ -1,7 +1,7 @@
 #[cfg(feature = "serde-serialize")]
 use serde::{Deserialize, Serialize};
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 
 use crate::allocator::Allocator;
 use crate::base::{DefaultAllocator, Matrix, MatrixMN, MatrixN, SquareMatrix};
@@ -26,19 +26,19 @@ use crate::storage::{Storage, StorageMut};
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct Cholesky<N: Complex, D: Dim>
+pub struct Cholesky<N: ComplexField, D: Dim>
 where DefaultAllocator: Allocator<N, D, D>
 {
     chol: MatrixN<N, D>,
 }
 
-impl<N: Complex, D: Dim> Copy for Cholesky<N, D>
+impl<N: ComplexField, D: Dim> Copy for Cholesky<N, D>
 where
     DefaultAllocator: Allocator<N, D, D>,
     MatrixN<N, D>: Copy,
 {}
 
-impl<N: Complex, D: DimSub<Dynamic>> Cholesky<N, D>
+impl<N: ComplexField, D: DimSub<Dynamic>> Cholesky<N, D>
 where DefaultAllocator: Allocator<N, D, D>
 {
     /// Attempts to compute the Cholesky decomposition of `matrix`.
@@ -148,7 +148,7 @@ where DefaultAllocator: Allocator<N, D, D>
     }
 }
 
-impl<N: Complex, D: DimSub<Dynamic>, S: Storage<N, D, D>> SquareMatrix<N, D, S>
+impl<N: ComplexField, D: DimSub<Dynamic>, S: Storage<N, D, D>> SquareMatrix<N, D, S>
 where DefaultAllocator: Allocator<N, D, D>
 {
     /// Attempts to compute the Cholesky decomposition of this matrix.
diff --git a/src/linalg/determinant.rs b/src/linalg/determinant.rs
index 4d265164..54ec9e5a 100644
--- a/src/linalg/determinant.rs
+++ b/src/linalg/determinant.rs
@@ -1,4 +1,4 @@
-use alga::general::Complex;
+use alga::general::ComplexField;
 
 use crate::base::allocator::Allocator;
 use crate::base::dimension::DimMin;
@@ -7,7 +7,7 @@ use crate::base::{DefaultAllocator, SquareMatrix};
 
 use crate::linalg::LU;
 
-impl<N: Complex, D: DimMin<D, Output = D>, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
+impl<N: ComplexField, D: DimMin<D, Output = D>, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
     /// Computes the matrix determinant.
     ///
     /// If the matrix has a dimension larger than 3, an LU decomposition is used.
diff --git a/src/linalg/eigen.rs b/src/linalg/eigen.rs
index 970e990d..8d1d26ca 100644
--- a/src/linalg/eigen.rs
+++ b/src/linalg/eigen.rs
@@ -1,7 +1,7 @@
 #[cfg(feature = "serde-serialize")]
 use serde::{Serialize, Deserialize};
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use num_complex::Complex;
 use std::cmp;
 use std::fmt::Display;
@@ -40,7 +40,7 @@ use crate::linalg::Schur;
     )
 )]
 #[derive(Clone, Debug)]
-pub struct Eigen<N: Complex, D: Dim>
+pub struct Eigen<N: ComplexField, D: Dim>
 where
     DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>,
 {
@@ -48,7 +48,7 @@ where
     pub eigenvalues: VectorN<N, D>,
 }
 
-impl<N: Complex, D: Dim> Copy for Eigen<N, D>
+impl<N: ComplexField, D: Dim> Copy for Eigen<N, D>
 where
     DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>,
     MatrixN<N, D>: Copy,
@@ -56,7 +56,7 @@ where
 {
 }
 
-impl<N: Complex, D: Dim> Eigen<N, D>
+impl<N: ComplexField, D: Dim> Eigen<N, D>
 where
     D: DimSub<U1>,                                   // For Hessenberg.
     ShapeConstraint: DimEq<Dynamic, DimDiff<D, U1>>, // For Hessenberg.
diff --git a/src/linalg/full_piv_lu.rs b/src/linalg/full_piv_lu.rs
index 61bc227c..04f61faf 100644
--- a/src/linalg/full_piv_lu.rs
+++ b/src/linalg/full_piv_lu.rs
@@ -1,7 +1,7 @@
 #[cfg(feature = "serde-serialize")]
 use serde::{Deserialize, Serialize};
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::allocator::Allocator;
 use crate::base::{DefaultAllocator, Matrix, MatrixMN, MatrixN};
 use crate::constraint::{SameNumberOfRows, ShapeConstraint};
@@ -32,7 +32,7 @@ use crate::linalg::PermutationSequence;
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct FullPivLU<N: Complex, R: DimMin<C>, C: Dim>
+pub struct FullPivLU<N: ComplexField, R: DimMin<C>, C: Dim>
 where DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimum<R, C>>
 {
     lu: MatrixMN<N, R, C>,
@@ -40,14 +40,14 @@ where DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimu
     q: PermutationSequence<DimMinimum<R, C>>,
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim> Copy for FullPivLU<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> Copy for FullPivLU<N, R, C>
 where
     DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimum<R, C>>,
     MatrixMN<N, R, C>: Copy,
     PermutationSequence<DimMinimum<R, C>>: Copy,
 {}
 
-impl<N: Complex, R: DimMin<C>, C: Dim> FullPivLU<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> FullPivLU<N, R, C>
 where DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimum<R, C>>
 {
     /// Computes the LU decomposition with full pivoting of `matrix`.
@@ -156,7 +156,7 @@ where DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimu
     }
 }
 
-impl<N: Complex, D: DimMin<D, Output = D>> FullPivLU<N, D, D>
+impl<N: ComplexField, D: DimMin<D, Output = D>> FullPivLU<N, D, D>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<(usize, usize), D>
 {
     /// Solves the linear system `self * x = b`, where `x` is the unknown to be determined.
@@ -261,7 +261,7 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<(usize, usize), D>
     }
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
+impl<N: ComplexField, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
 where DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimum<R, C>>
 {
     /// Computes the LU decomposition with full pivoting of `matrix`.
diff --git a/src/linalg/givens.rs b/src/linalg/givens.rs
index 0927dd26..3dd56ee3 100644
--- a/src/linalg/givens.rs
+++ b/src/linalg/givens.rs
@@ -1,6 +1,6 @@
 //! Construction of givens rotations.
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use num::{Zero, One};
 
 use crate::base::dimension::{Dim, U2};
@@ -11,13 +11,13 @@ use crate::base::{Vector, Matrix};
 
 /// A Givens rotation.
 #[derive(Debug, Clone, Copy)]
-pub struct GivensRotation<N: Complex> {
+pub struct GivensRotation<N: ComplexField> {
     c: N::Real,
     s: N
 }
 
 // Matrix = UnitComplex * Matrix
-impl<N: Complex> GivensRotation<N> {
+impl<N: ComplexField> GivensRotation<N> {
     /// The Givents rotation that does nothing.
     pub fn identity() -> Self {
         Self {
diff --git a/src/linalg/hessenberg.rs b/src/linalg/hessenberg.rs
index 56e2ace4..c87c06a2 100644
--- a/src/linalg/hessenberg.rs
+++ b/src/linalg/hessenberg.rs
@@ -1,7 +1,7 @@
 #[cfg(feature = "serde-serialize")]
 use serde::{Deserialize, Serialize};
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::allocator::Allocator;
 use crate::base::{DefaultAllocator, MatrixMN, MatrixN, SquareMatrix, VectorN};
 use crate::dimension::{DimDiff, DimSub, U1};
@@ -30,21 +30,21 @@ use crate::linalg::householder;
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct Hessenberg<N: Complex, D: DimSub<U1>>
+pub struct Hessenberg<N: ComplexField, D: DimSub<U1>>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N, DimDiff<D, U1>>
 {
     hess: MatrixN<N, D>,
     subdiag: VectorN<N, DimDiff<D, U1>>,
 }
 
-impl<N: Complex, D: DimSub<U1>> Copy for Hessenberg<N, D>
+impl<N: ComplexField, D: DimSub<U1>> Copy for Hessenberg<N, D>
 where
     DefaultAllocator: Allocator<N, D, D> + Allocator<N, DimDiff<D, U1>>,
     MatrixN<N, D>: Copy,
     VectorN<N, DimDiff<D, U1>>: Copy,
 {}
 
-impl<N: Complex, D: DimSub<U1>> Hessenberg<N, D>
+impl<N: ComplexField, D: DimSub<U1>> Hessenberg<N, D>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D> + Allocator<N, DimDiff<D, U1>>
 {
     /// Computes the Hessenberg decomposition using householder reflections.
@@ -133,7 +133,7 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D> + Allocator<N, DimD
     }
 }
 
-impl<N: Complex, D: DimSub<U1>, S: Storage<N, D, D>> SquareMatrix<N, D, S>
+impl<N: ComplexField, D: DimSub<U1>, S: Storage<N, D, D>> SquareMatrix<N, D, S>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D> + Allocator<N, DimDiff<D, U1>>
 {
     /// Computes the Hessenberg decomposition of this matrix using householder reflections.
diff --git a/src/linalg/householder.rs b/src/linalg/householder.rs
index 3bd8f801..11fa32eb 100644
--- a/src/linalg/householder.rs
+++ b/src/linalg/householder.rs
@@ -1,7 +1,7 @@
 //! Construction of householder elementary reflections.
 
 use num::Zero;
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::allocator::Allocator;
 use crate::base::{DefaultAllocator, MatrixMN, MatrixN, Unit, Vector, VectorN};
 use crate::dimension::Dim;
@@ -16,7 +16,7 @@ use crate::geometry::Reflection;
 /// `column` after reflection and `false` if no reflection was necessary.
 #[doc(hidden)]
 #[inline(always)]
-pub fn reflection_axis_mut<N: Complex, D: Dim, S: StorageMut<N, D>>(
+pub fn reflection_axis_mut<N: ComplexField, D: Dim, S: StorageMut<N, D>>(
     column: &mut Vector<N, D, S>,
 ) -> (N, bool) {
     let reflection_sq_norm = column.norm_squared();
@@ -44,7 +44,7 @@ pub fn reflection_axis_mut<N: Complex, D: Dim, S: StorageMut<N, D>>(
 /// Uses an householder reflection to zero out the `icol`-th column, starting with the `shift + 1`-th
 /// subdiagonal element.
 #[doc(hidden)]
-pub fn clear_column_unchecked<N: Complex, R: Dim, C: Dim>(
+pub fn clear_column_unchecked<N: ComplexField, R: Dim, C: Dim>(
     matrix: &mut MatrixMN<N, R, C>,
     diag_elt: &mut N,
     icol: usize,
@@ -72,7 +72,7 @@ pub fn clear_column_unchecked<N: Complex, R: Dim, C: Dim>(
 /// Uses an householder reflection to zero out the `irow`-th row, ending before the `shift + 1`-th
 /// superdiagonal element.
 #[doc(hidden)]
-pub fn clear_row_unchecked<N: Complex, R: Dim, C: Dim>(
+pub fn clear_row_unchecked<N: ComplexField, R: Dim, C: Dim>(
     matrix: &mut MatrixMN<N, R, C>,
     diag_elt: &mut N,
     axis_packed: &mut VectorN<N, C>,
@@ -108,7 +108,7 @@ pub fn clear_row_unchecked<N: Complex, R: Dim, C: Dim>(
 /// the lower-diagonal element of the given matrix.
 /// matrices.
 #[doc(hidden)]
-pub fn assemble_q<N: Complex, D: Dim>(m: &MatrixN<N, D>, signs: &[N]) -> MatrixN<N, D>
+pub fn assemble_q<N: ComplexField, D: Dim>(m: &MatrixN<N, D>, signs: &[N]) -> MatrixN<N, D>
 where DefaultAllocator: Allocator<N, D, D> {
     assert!(m.is_square());
     let dim = m.data.shape().0;
diff --git a/src/linalg/inverse.rs b/src/linalg/inverse.rs
index 3d6a0681..2d3b18ad 100644
--- a/src/linalg/inverse.rs
+++ b/src/linalg/inverse.rs
@@ -1,4 +1,4 @@
-use alga::general::Complex;
+use alga::general::ComplexField;
 
 use crate::base::allocator::Allocator;
 use crate::base::dimension::Dim;
@@ -7,7 +7,7 @@ use crate::base::{DefaultAllocator, MatrixN, SquareMatrix};
 
 use crate::linalg::lu;
 
-impl<N: Complex, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
+impl<N: ComplexField, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
     /// Attempts to invert this matrix.
     #[inline]
     pub fn try_inverse(self) -> Option<MatrixN<N, D>>
@@ -21,7 +21,7 @@ impl<N: Complex, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
     }
 }
 
-impl<N: Complex, D: Dim, S: StorageMut<N, D, D>> SquareMatrix<N, D, S> {
+impl<N: ComplexField, D: Dim, S: StorageMut<N, D, D>> SquareMatrix<N, D, S> {
     /// Attempts to invert this matrix in-place. Returns `false` and leaves `self` untouched if
     /// inversion fails.
     #[inline]
@@ -115,7 +115,7 @@ impl<N: Complex, D: Dim, S: StorageMut<N, D, D>> SquareMatrix<N, D, S> {
 }
 
 // NOTE: this is an extremely efficient, loop-unrolled matrix inverse from MESA (MIT licensed).
-fn do_inverse4<N: Complex, D: Dim, S: StorageMut<N, D, D>>(
+fn do_inverse4<N: ComplexField, D: Dim, S: StorageMut<N, D, D>>(
     m: &MatrixN<N, D>,
     out: &mut SquareMatrix<N, D, S>,
 ) -> bool
diff --git a/src/linalg/lu.rs b/src/linalg/lu.rs
index fffe7148..2c3beee3 100644
--- a/src/linalg/lu.rs
+++ b/src/linalg/lu.rs
@@ -1,7 +1,7 @@
 #[cfg(feature = "serde-serialize")]
 use serde::{Deserialize, Serialize};
 
-use alga::general::{Field, Complex};
+use alga::general::{Field, ComplexField};
 use crate::allocator::{Allocator, Reallocator};
 use crate::base::{DefaultAllocator, Matrix, MatrixMN, MatrixN, Scalar};
 use crate::constraint::{SameNumberOfRows, ShapeConstraint};
@@ -32,14 +32,14 @@ use crate::linalg::PermutationSequence;
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct LU<N: Complex, R: DimMin<C>, C: Dim>
+pub struct LU<N: ComplexField, R: DimMin<C>, C: Dim>
 where DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimum<R, C>>
 {
     lu: MatrixMN<N, R, C>,
     p: PermutationSequence<DimMinimum<R, C>>,
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim> Copy for LU<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> Copy for LU<N, R, C>
 where
     DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimum<R, C>>,
     MatrixMN<N, R, C>: Copy,
@@ -49,7 +49,7 @@ where
 /// Performs a LU decomposition to overwrite `out` with the inverse of `matrix`.
 ///
 /// If `matrix` is not invertible, `false` is returned and `out` may contain invalid data.
-pub fn try_invert_to<N: Complex, D: Dim, S>(
+pub fn try_invert_to<N: ComplexField, D: Dim, S>(
     mut matrix: MatrixN<N, D>,
     out: &mut Matrix<N, D, D, S>,
 ) -> bool
@@ -86,7 +86,7 @@ where
     matrix.solve_upper_triangular_mut(out)
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim> LU<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> LU<N, R, C>
 where DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimum<R, C>>
 {
     /// Computes the LU decomposition with partial (row) pivoting of `matrix`.
@@ -197,7 +197,7 @@ where DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimu
     }
 }
 
-impl<N: Complex, D: DimMin<D, Output = D>> LU<N, D, D>
+impl<N: ComplexField, D: DimMin<D, Output = D>> LU<N, D, D>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<(usize, usize), D>
 {
     /// Solves the linear system `self * x = b`, where `x` is the unknown to be determined.
@@ -368,7 +368,7 @@ pub fn gauss_step_swap<N, R: Dim, C: Dim, S>(
     }
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
+impl<N: ComplexField, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
 where DefaultAllocator: Allocator<N, R, C> + Allocator<(usize, usize), DimMinimum<R, C>>
 {
     /// Computes the LU decomposition with partial (row) pivoting of `matrix`.
diff --git a/src/linalg/qr.rs b/src/linalg/qr.rs
index 7c02d140..683c11b8 100644
--- a/src/linalg/qr.rs
+++ b/src/linalg/qr.rs
@@ -2,7 +2,7 @@
 use serde::{Deserialize, Serialize};
 use num::Zero;
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::allocator::{Allocator, Reallocator};
 use crate::base::{DefaultAllocator, Matrix, MatrixMN, MatrixN, Unit, VectorN};
 use crate::constraint::{SameNumberOfRows, ShapeConstraint};
@@ -33,21 +33,21 @@ use crate::linalg::householder;
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct QR<N: Complex, R: DimMin<C>, C: Dim>
+pub struct QR<N: ComplexField, R: DimMin<C>, C: Dim>
 where DefaultAllocator: Allocator<N, R, C> + Allocator<N, DimMinimum<R, C>>
 {
     qr: MatrixMN<N, R, C>,
     diag: VectorN<N, DimMinimum<R, C>>,
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim> Copy for QR<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> Copy for QR<N, R, C>
 where
     DefaultAllocator: Allocator<N, R, C> + Allocator<N, DimMinimum<R, C>>,
     MatrixMN<N, R, C>: Copy,
     VectorN<N, DimMinimum<R, C>>: Copy,
 {}
 
-impl<N: Complex, R: DimMin<C>, C: Dim> QR<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> QR<N, R, C>
 where DefaultAllocator: Allocator<N, R, C> + Allocator<N, R> + Allocator<N, DimMinimum<R, C>>
 {
     /// Computes the QR decomposition using householder reflections.
@@ -159,7 +159,7 @@ where DefaultAllocator: Allocator<N, R, C> + Allocator<N, R> + Allocator<N, DimM
     }
 }
 
-impl<N: Complex, D: DimMin<D, Output = D>> QR<N, D, D>
+impl<N: ComplexField, D: DimMin<D, Output = D>> QR<N, D, D>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>
 {
     /// Solves the linear system `self * x = b`, where `x` is the unknown to be determined.
@@ -291,7 +291,7 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<N, D>
     // }
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
+impl<N: ComplexField, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
 where DefaultAllocator: Allocator<N, R, C> + Allocator<N, R> + Allocator<N, DimMinimum<R, C>>
 {
     /// Computes the QR decomposition of this matrix.
diff --git a/src/linalg/schur.rs b/src/linalg/schur.rs
index ade4b2dc..9407cf7f 100644
--- a/src/linalg/schur.rs
+++ b/src/linalg/schur.rs
@@ -2,7 +2,7 @@
 use serde::{Deserialize, Serialize};
 
 use approx::AbsDiffEq;
-use alga::general::{Complex, Real};
+use alga::general::{ComplexField, Real};
 use num_complex::Complex as NumComplex;
 use std::cmp;
 
@@ -35,20 +35,20 @@ use crate::linalg::givens::GivensRotation;
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct Schur<N: Complex, D: Dim>
+pub struct Schur<N: ComplexField, D: Dim>
 where DefaultAllocator: Allocator<N, D, D>
 {
     q: MatrixN<N, D>,
     t: MatrixN<N, D>,
 }
 
-impl<N: Complex, D: Dim> Copy for Schur<N, D>
+impl<N: ComplexField, D: Dim> Copy for Schur<N, D>
 where
     DefaultAllocator: Allocator<N, D, D>,
     MatrixN<N, D>: Copy,
 {}
 
-impl<N: Complex, D: Dim> Schur<N, D>
+impl<N: ComplexField, D: Dim> Schur<N, D>
 where
     D: DimSub<U1>,                                   // For Hessenberg.
     DefaultAllocator: Allocator<N, D, DimDiff<D, U1>>
@@ -398,7 +398,7 @@ where
     }
 }
 
-fn decompose_2x2<N: Complex, D: Dim>(
+fn decompose_2x2<N: ComplexField, D: Dim>(
     mut m: MatrixN<N, D>,
     compute_q: bool,
 ) -> Option<(Option<MatrixN<N, D>>, MatrixN<N, D>)>
@@ -435,7 +435,7 @@ where
     Some((q, m))
 }
 
-fn compute_2x2_eigvals<N: Complex, S: Storage<N, U2, U2>>(
+fn compute_2x2_eigvals<N: ComplexField, S: Storage<N, U2, U2>>(
     m: &SquareMatrix<N, U2, S>,
 ) -> Option<(N, N)> {
     // Solve the 2x2 eigenvalue subproblem.
@@ -461,7 +461,7 @@ fn compute_2x2_eigvals<N: Complex, S: Storage<N, U2, U2>>(
 ///
 /// Returns `None` if the matrix has complex eigenvalues, or is upper-triangular. In both case,
 /// the basis is the identity.
-fn compute_2x2_basis<N: Complex, S: Storage<N, U2, U2>>(
+fn compute_2x2_basis<N: ComplexField, S: Storage<N, U2, U2>>(
     m: &SquareMatrix<N, U2, S>,
 ) -> Option<GivensRotation<N>> {
     let h10 = m[(1, 0)];
@@ -487,7 +487,7 @@ fn compute_2x2_basis<N: Complex, S: Storage<N, U2, U2>>(
     }
 }
 
-impl<N: Complex, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S>
+impl<N: ComplexField, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S>
 where
     D: DimSub<U1>,                                   // For Hessenberg.
     DefaultAllocator: Allocator<N, D, DimDiff<D, U1>>
diff --git a/src/linalg/solve.rs b/src/linalg/solve.rs
index 0c2bb684..f10b1d00 100644
--- a/src/linalg/solve.rs
+++ b/src/linalg/solve.rs
@@ -1,4 +1,4 @@
-use alga::general::Complex;
+use alga::general::ComplexField;
 
 use crate::base::allocator::Allocator;
 use crate::base::constraint::{SameNumberOfRows, ShapeConstraint};
@@ -6,7 +6,7 @@ use crate::base::dimension::{Dim, U1};
 use crate::base::storage::{Storage, StorageMut};
 use crate::base::{DefaultAllocator, Matrix, MatrixMN, SquareMatrix, Vector, DVectorSlice};
 
-impl<N: Complex, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
+impl<N: ComplexField, D: Dim, S: Storage<N, D, D>> SquareMatrix<N, D, S> {
     /// Computes the solution of the linear system `self . x = b` where `x` is the unknown and only
     /// the lower-triangular part of `self` (including the diagonal) is considered not-zero.
     #[inline]
diff --git a/src/linalg/svd.rs b/src/linalg/svd.rs
index 779bcff0..714f398f 100644
--- a/src/linalg/svd.rs
+++ b/src/linalg/svd.rs
@@ -4,7 +4,7 @@ use serde::{Deserialize, Serialize};
 use num::{Zero, One};
 use approx::AbsDiffEq;
 
-use alga::general::{Real, Complex};
+use alga::general::{Real, ComplexField};
 use crate::allocator::Allocator;
 use crate::base::{DefaultAllocator, Matrix, Matrix2x3, MatrixMN, Vector2, VectorN};
 use crate::constraint::{SameNumberOfRows, ShapeConstraint};
@@ -40,7 +40,7 @@ use crate::linalg::givens::GivensRotation;
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct SVD<N: Complex, R: DimMin<C>, C: Dim>
+pub struct SVD<N: ComplexField, R: DimMin<C>, C: Dim>
 where DefaultAllocator: Allocator<N, DimMinimum<R, C>, C>
         + Allocator<N, R, DimMinimum<R, C>>
         + Allocator<N::Real, DimMinimum<R, C>>
@@ -53,7 +53,7 @@ where DefaultAllocator: Allocator<N, DimMinimum<R, C>, C>
     pub singular_values: VectorN<N::Real, DimMinimum<R, C>>,
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim> Copy for SVD<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> Copy for SVD<N, R, C>
 where
     DefaultAllocator: Allocator<N, DimMinimum<R, C>, C>
         + Allocator<N, R, DimMinimum<R, C>>
@@ -63,7 +63,7 @@ where
     VectorN<N::Real, DimMinimum<R, C>>: Copy,
 {}
 
-impl<N: Complex, R: DimMin<C>, C: Dim> SVD<N, R, C>
+impl<N: ComplexField, R: DimMin<C>, C: Dim> SVD<N, R, C>
 where
     DimMinimum<R, C>: DimSub<U1>, // for Bidiagonal.
     DefaultAllocator: Allocator<N, R, C>
@@ -546,7 +546,7 @@ where
     }
 }
 
-impl<N: Complex, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
+impl<N: ComplexField, R: DimMin<C>, C: Dim, S: Storage<N, R, C>> Matrix<N, R, C, S>
 where
     DimMinimum<R, C>: DimSub<U1>, // for Bidiagonal.
     DefaultAllocator: Allocator<N, R, C>
diff --git a/src/linalg/symmetric_eigen.rs b/src/linalg/symmetric_eigen.rs
index 6fd1fea5..78e00509 100644
--- a/src/linalg/symmetric_eigen.rs
+++ b/src/linalg/symmetric_eigen.rs
@@ -4,7 +4,7 @@ use serde::{Deserialize, Serialize};
 use num::Zero;
 use approx::AbsDiffEq;
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::allocator::Allocator;
 use crate::base::{DefaultAllocator, Matrix2, MatrixN, SquareMatrix, Vector2, VectorN};
 use crate::dimension::{Dim, DimDiff, DimSub, U1, U2};
@@ -34,7 +34,7 @@ use crate::linalg::SymmetricTridiagonal;
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct SymmetricEigen<N: Complex, D: Dim>
+pub struct SymmetricEigen<N: ComplexField, D: Dim>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N::Real, D>
 {
     /// The eigenvectors of the decomposed matrix.
@@ -44,14 +44,14 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<N::Real, D>
     pub eigenvalues: VectorN<N::Real, D>,
 }
 
-impl<N: Complex, D: Dim> Copy for SymmetricEigen<N, D>
+impl<N: ComplexField, D: Dim> Copy for SymmetricEigen<N, D>
 where
     DefaultAllocator: Allocator<N, D, D> + Allocator<N::Real, D>,
     MatrixN<N, D>: Copy,
     VectorN<N::Real, D>: Copy,
 {}
 
-impl<N: Complex, D: Dim> SymmetricEigen<N, D>
+impl<N: ComplexField, D: Dim> SymmetricEigen<N, D>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N::Real, D>
 {
     /// Computes the eigendecomposition of the given symmetric matrix.
@@ -288,7 +288,7 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<N::Real, D>
 /// The inputs are interpreted as the 2x2 matrix:
 ///     tmm  tmn
 ///     tmn  tnn
-pub fn wilkinson_shift<N: Complex>(tmm: N, tnn: N, tmn: N) -> N {
+pub fn wilkinson_shift<N: ComplexField>(tmm: N, tnn: N, tmn: N) -> N {
     let sq_tmn = tmn * tmn;
     if !sq_tmn.is_zero() {
         // We have the guarantee that the denominator won't be zero.
@@ -304,7 +304,7 @@ pub fn wilkinson_shift<N: Complex>(tmm: N, tnn: N, tmn: N) -> N {
  * Computations of eigenvalues for symmetric matrices.
  *
  */
-impl<N: Complex, D: DimSub<U1>, S: Storage<N, D, D>> SquareMatrix<N, D, S>
+impl<N: ComplexField, D: DimSub<U1>, S: Storage<N, D, D>> SquareMatrix<N, D, S>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N, DimDiff<D, U1>> +
                         Allocator<N::Real, D> + Allocator<N::Real, DimDiff<D, U1>>
 {
diff --git a/src/linalg/symmetric_tridiagonal.rs b/src/linalg/symmetric_tridiagonal.rs
index 386c0a6f..0a04dae8 100644
--- a/src/linalg/symmetric_tridiagonal.rs
+++ b/src/linalg/symmetric_tridiagonal.rs
@@ -1,7 +1,7 @@
 #[cfg(feature = "serde-serialize")]
 use serde::{Deserialize, Serialize};
 
-use alga::general::Complex;
+use alga::general::ComplexField;
 use crate::allocator::Allocator;
 use crate::base::{DefaultAllocator, MatrixMN, MatrixN, SquareMatrix, VectorN};
 use crate::dimension::{DimDiff, DimSub, U1};
@@ -30,21 +30,21 @@ use crate::linalg::householder;
     ))
 )]
 #[derive(Clone, Debug)]
-pub struct SymmetricTridiagonal<N: Complex, D: DimSub<U1>>
+pub struct SymmetricTridiagonal<N: ComplexField, D: DimSub<U1>>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N, DimDiff<D, U1>>
 {
     tri: MatrixN<N, D>,
     off_diagonal: VectorN<N, DimDiff<D, U1>>,
 }
 
-impl<N: Complex, D: DimSub<U1>> Copy for SymmetricTridiagonal<N, D>
+impl<N: ComplexField, D: DimSub<U1>> Copy for SymmetricTridiagonal<N, D>
 where
     DefaultAllocator: Allocator<N, D, D> + Allocator<N, DimDiff<D, U1>>,
     MatrixN<N, D>: Copy,
     VectorN<N, DimDiff<D, U1>>: Copy,
 {}
 
-impl<N: Complex, D: DimSub<U1>> SymmetricTridiagonal<N, D>
+impl<N: ComplexField, D: DimSub<U1>> SymmetricTridiagonal<N, D>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N, DimDiff<D, U1>>
 {
     /// Computes the tridiagonalization of the symmetric matrix `m`.
@@ -145,7 +145,7 @@ where DefaultAllocator: Allocator<N, D, D> + Allocator<N, DimDiff<D, U1>>
     }
 }
 
-impl<N: Complex, D: DimSub<U1>, S: Storage<N, D, D>> SquareMatrix<N, D, S>
+impl<N: ComplexField, D: DimSub<U1>, S: Storage<N, D, D>> SquareMatrix<N, D, S>
 where DefaultAllocator: Allocator<N, D, D> + Allocator<N, DimDiff<D, U1>>
 {
     /// Computes the tridiagonalization of this symmetric matrix.
diff --git a/tests/core/matrix.rs b/tests/core/matrix.rs
index 5ba06f5b..81de11b0 100644
--- a/tests/core/matrix.rs
+++ b/tests/core/matrix.rs
@@ -1022,7 +1022,7 @@ mod finite_dim_inner_space_tests {
      *
      */
     #[cfg(feature = "arbitrary")]
-    fn is_subspace_basis<T: FiniteDimInnerSpace<Real = f64, Complex = f64> + Display>(vs: &[T]) -> bool {
+    fn is_subspace_basis<T: FiniteDimInnerSpace<Real = f64, ComplexField = f64> + Display>(vs: &[T]) -> bool {
         for i in 0..vs.len() {
             // Basis elements must be normalized.
             if !relative_eq!(vs[i].norm(), 1.0, epsilon = 1.0e-7) {
diff --git a/tests/linalg/solve.rs b/tests/linalg/solve.rs
index 90b76585..e917cfc5 100644
--- a/tests/linalg/solve.rs
+++ b/tests/linalg/solve.rs
@@ -4,11 +4,11 @@
 macro_rules! gen_tests(
     ($module: ident, $scalar: ty) => {
         mod $module {
-            use na::{Matrix4, Matrix4x5, Complex};
+            use na::{Matrix4, Matrix4x5, ComplexField};
             #[allow(unused_imports)]
             use crate::core::helper::{RandScalar, RandComplex};
 
-            fn unzero_diagonal<N: Complex>(a: &mut Matrix4<N>) {
+            fn unzero_diagonal<N: ComplexField>(a: &mut Matrix4<N>) {
                 for i in 0..4 {
                     if a[(i, i)].norm1() < na::convert(1.0e-7) {
                         a[(i, i)] = N::one();
diff --git a/tests/linalg/svd.rs b/tests/linalg/svd.rs
index a9b1d1cf..ca7bab4c 100644
--- a/tests/linalg/svd.rs
+++ b/tests/linalg/svd.rs
@@ -8,7 +8,7 @@ mod quickcheck_tests {
             mod $module {
                 use na::{
                     DMatrix, DVector, Matrix2, Matrix2x5, Matrix3, Matrix3x5, Matrix4, Matrix5x2, Matrix5x3,
-                    Complex
+                    ComplexField
                 };
                 use std::cmp;
                 #[allow(unused_imports)]
@@ -21,7 +21,7 @@ mod quickcheck_tests {
                             let svd = m.clone().svd(true, true);
                             let recomp_m = svd.clone().recompose().unwrap();
                             let (u, s, v_t) = (svd.u.unwrap(), svd.singular_values, svd.v_t.unwrap());
-                            let ds = DMatrix::from_diagonal(&s.map(|e| Complex::from_real(e)));
+                            let ds = DMatrix::from_diagonal(&s.map(|e| ComplexField::from_real(e)));
 
                             s.iter().all(|e| *e >= 0.0) &&
                             relative_eq!(&u * ds * &v_t, recomp_m, epsilon = 1.0e-5) &&
@@ -36,7 +36,7 @@ mod quickcheck_tests {
                         let m = m.map(|e| e.0);
                         let svd = m.svd(true, true);
                         let (u, s, v_t) = (svd.u.unwrap(), svd.singular_values, svd.v_t.unwrap());
-                        let ds = Matrix3::from_diagonal(&s.map(|e| Complex::from_real(e)));
+                        let ds = Matrix3::from_diagonal(&s.map(|e| ComplexField::from_real(e)));
 
                         s.iter().all(|e| *e >= 0.0) &&
                         relative_eq!(m, &u * ds * &v_t, epsilon = 1.0e-5) &&
@@ -48,7 +48,7 @@ mod quickcheck_tests {
                         let m = m.map(|e| e.0);
                         let svd = m.svd(true, true);
                         let (u, s, v_t) = (svd.u.unwrap(), svd.singular_values, svd.v_t.unwrap());
-                        let ds = Matrix2::from_diagonal(&s.map(|e| Complex::from_real(e)));
+                        let ds = Matrix2::from_diagonal(&s.map(|e| ComplexField::from_real(e)));
 
                         s.iter().all(|e| *e >= 0.0) &&
                         relative_eq!(m, &u * ds * &v_t, epsilon = 1.0e-5) &&
@@ -61,7 +61,7 @@ mod quickcheck_tests {
                         let svd = m.svd(true, true);
                         let (u, s, v_t) = (svd.u.unwrap(), svd.singular_values, svd.v_t.unwrap());
 
-                        let ds = Matrix3::from_diagonal(&s.map(|e| Complex::from_real(e)));
+                        let ds = Matrix3::from_diagonal(&s.map(|e| ComplexField::from_real(e)));
 
                         s.iter().all(|e| *e >= 0.0) &&
                         relative_eq!(m, u * ds * v_t, epsilon = 1.0e-5)
@@ -71,7 +71,7 @@ mod quickcheck_tests {
                         let m = m.map(|e| e.0);
                         let svd = m.svd(true, true);
                         let (u, s, v_t) = (svd.u.unwrap(), svd.singular_values, svd.v_t.unwrap());
-                        let ds = Matrix2::from_diagonal(&s.map(|e| Complex::from_real(e)));
+                        let ds = Matrix2::from_diagonal(&s.map(|e| ComplexField::from_real(e)));
 
                         s.iter().all(|e| *e >= 0.0) &&
                         relative_eq!(m, u * ds * v_t, epsilon = 1.0e-5)
@@ -81,7 +81,7 @@ mod quickcheck_tests {
                         let m = m.map(|e| e.0);
                         let svd = m.svd(true, true);
                         let (u, s, v_t) = (svd.u.unwrap(), svd.singular_values, svd.v_t.unwrap());
-                        let ds = Matrix4::from_diagonal(&s.map(|e| Complex::from_real(e)));
+                        let ds = Matrix4::from_diagonal(&s.map(|e| ComplexField::from_real(e)));
 
                         s.iter().all(|e| *e >= 0.0) &&
                         relative_eq!(m, u * ds * v_t, epsilon = 1.0e-5) &&
@@ -93,7 +93,7 @@ mod quickcheck_tests {
                         let m = m.map(|e| e.0);
                         let svd = m.svd(true, true);
                         let (u, s, v_t) = (svd.u.unwrap(), svd.singular_values, svd.v_t.unwrap());
-                        let ds = Matrix2::from_diagonal(&s.map(|e| Complex::from_real(e)));
+                        let ds = Matrix2::from_diagonal(&s.map(|e| ComplexField::from_real(e)));
 
                         s.iter().all(|e| *e >= 0.0) &&
                         relative_eq!(m, u * ds * v_t, epsilon = 1.0e-5) &&