diff --git a/src/Magnum/Math/Matrix.h b/src/Magnum/Math/Matrix.h
index fa05b1e1b..24c26a507 100644
--- a/src/Magnum/Math/Matrix.h
+++ b/src/Magnum/Math/Matrix.h
@@ -119,10 +119,10 @@ template<std::size_t size, class T> class Matrix: public RectangularMatrix<size,
          * tr(A) = \sum_{i=1}^n a_{i,i}
          * @f]
          */
-        T trace() const { return RectangularMatrix<size, size, T>::diagonal().sum(); }
+        constexpr T trace() const { return RectangularMatrix<size, size, T>::diagonal().sum(); }
 
         /** @brief %Matrix without given column and row */
-        Matrix<size-1, T> ij(std::size_t skipCol, std::size_t skipRow) const;
+        constexpr Matrix<size-1, T> ij(std::size_t skipCol, std::size_t skipRow) const;
 
         /**
          * @brief Determinant
@@ -135,7 +135,7 @@ template<std::size_t size, class T> class Matrix: public RectangularMatrix<size,
          *      \det(A) = a_{0, 0} a_{1, 1} - a_{1, 0} a_{0, 1}
          * @f]
          */
-        T determinant() const { return Implementation::MatrixDeterminant<size, T>()(*this); }
+        constexpr T determinant() const { return Implementation::MatrixDeterminant<size, T>()(*this); }
 
         /**
          * @brief Inverted matrix
@@ -147,7 +147,7 @@ template<std::size_t size, class T> class Matrix: public RectangularMatrix<size,
          * @ref Matrix4::invertedRigid() which are faster alternatives for
          * particular matrix types.
          */
-        Matrix<size, T> inverted() const;
+        constexpr Matrix<size, T> inverted() const;
 
         /**
          * @brief Inverted orthogonal matrix
@@ -168,13 +168,13 @@ template<std::size_t size, class T> class Matrix: public RectangularMatrix<size,
 
         #ifndef DOXYGEN_GENERATING_OUTPUT
         /* Reimplementation of functions to return correct type */
-        Matrix<size, T> operator*(const Matrix<size, T>& other) const {
+        constexpr Matrix<size, T> operator*(const Matrix<size, T>& other) const {
             return RectangularMatrix<size, size, T>::operator*(other);
         }
-        template<std::size_t otherCols> RectangularMatrix<otherCols, size, T> operator*(const RectangularMatrix<otherCols, size, T>& other) const {
+        template<std::size_t otherCols> RectangularMatrix<otherCols, size, T> constexpr operator*(const RectangularMatrix<otherCols, size, T>& other) const {
             return RectangularMatrix<size, size, T>::operator*(other);
         }
-        Vector<size, T> operator*(const Vector<size, T>& other) const {
+        constexpr Vector<size, T> operator*(const Vector<size, T>& other) const {
             return RectangularMatrix<size, size, T>::operator*(other);
         }
         MAGNUM_RECTANGULARMATRIX_SUBCLASS_IMPLEMENTATION(size, size, Matrix<size, T>)
@@ -219,40 +219,40 @@ template<std::size_t size, class T> inline Corrade::Utility::Debug operator<<(Co
 
 #ifndef DOXYGEN_GENERATING_OUTPUT
 #define MAGNUM_MATRIX_SUBCLASS_IMPLEMENTATION(size, Type, VectorType)       \
-    VectorType<T>& operator[](std::size_t col) {                            \
+    constexpr VectorType<T>& operator[](std::size_t col) {                  \
         return static_cast<VectorType<T>&>(Matrix<size, T>::operator[](col)); \
     }                                                                       \
     constexpr const VectorType<T> operator[](std::size_t col) const {       \
         return VectorType<T>(Matrix<size, T>::operator[](col));             \
     }                                                                       \
-    VectorType<T> row(std::size_t row) const {                              \
+    constexpr VectorType<T> row(std::size_t row) const {                    \
         return VectorType<T>(Matrix<size, T>::row(row));                    \
     }                                                                       \
                                                                             \
-    Type<T> operator*(const Matrix<size, T>& other) const {                 \
+    constexpr Type<T> operator*(const Matrix<size, T>& other) const {       \
         return Matrix<size, T>::operator*(other);                           \
     }                                                                       \
-    template<std::size_t otherCols> RectangularMatrix<otherCols, size, T> operator*(const RectangularMatrix<otherCols, size, T>& other) const { \
+    template<std::size_t otherCols> constexpr RectangularMatrix<otherCols, size, T> operator*(const RectangularMatrix<otherCols, size, T>& other) const { \
         return Matrix<size, T>::operator*(other);                           \
     }                                                                       \
-    VectorType<T> operator*(const Vector<size, T>& other) const {           \
+    constexpr VectorType<T> operator*(const Vector<size, T>& other) const { \
         return Matrix<size, T>::operator*(other);                           \
     }                                                                       \
                                                                             \
-    Type<T> transposed() const { return Matrix<size, T>::transposed(); }    \
+    constexpr Type<T> transposed() const { return Matrix<size, T>::transposed(); } \
     constexpr VectorType<T> diagonal() const { return Matrix<size, T>::diagonal(); } \
-    Type<T> inverted() const { return Matrix<size, T>::inverted(); }        \
-    Type<T> invertedOrthogonal() const {                                    \
+    constexpr Type<T> inverted() const { return Matrix<size, T>::inverted(); } \
+    constexpr Type<T> invertedOrthogonal() const {                          \
         return Matrix<size, T>::invertedOrthogonal();                       \
     }
 
 namespace Implementation {
 
 template<std::size_t size, class T> struct MatrixDeterminant {
-    T operator()(const Matrix<size, T>& m);
+    constexpr T operator()(const Matrix<size, T>& m);
 };
 
-template<std::size_t size, class T> T MatrixDeterminant<size, T>::operator()(const Matrix<size, T>& m) {
+template<std::size_t size, class T> constexpr T MatrixDeterminant<size, T>::operator()(const Matrix<size, T>& m) {
     T out(0);
 
     for(std::size_t col = 0; col != size; ++col)
@@ -290,7 +290,7 @@ template<std::size_t size, class T> bool Matrix<size, T>::isOrthogonal() const {
     return true;
 }
 
-template<std::size_t size, class T> Matrix<size-1, T> Matrix<size, T>::ij(const std::size_t skipCol, const std::size_t skipRow) const {
+template<std::size_t size, class T> constexpr Matrix<size-1, T> Matrix<size, T>::ij(const std::size_t skipCol, const std::size_t skipRow) const {
     Matrix<size-1, T> out(Matrix<size-1, T>::Zero);
 
     for(std::size_t col = 0; col != size-1; ++col)
@@ -301,7 +301,7 @@ template<std::size_t size, class T> Matrix<size-1, T> Matrix<size, T>::ij(const
     return out;
 }
 
-template<std::size_t size, class T> Matrix<size, T> Matrix<size, T>::inverted() const {
+template<std::size_t size, class T> constexpr Matrix<size, T> Matrix<size, T>::inverted() const {
     Matrix<size, T> out(Zero);
 
     const T _determinant = determinant();
diff --git a/src/Magnum/Math/Test/MatrixTest.cpp b/src/Magnum/Math/Test/MatrixTest.cpp
index e61036340..dd80aefd5 100644
--- a/src/Magnum/Math/Test/MatrixTest.cpp
+++ b/src/Magnum/Math/Test/MatrixTest.cpp
@@ -214,7 +214,7 @@ void MatrixTest::isOrthogonal() {
 }
 
 void MatrixTest::trace() {
-    Matrix<5, Int> m(
+    constexpr Matrix<5, Int> m(
         Vector<5, Int>(1, 2,   3,  0,  0),
         Vector<5, Int>(2, 3,   2,  1, -2),
         Vector<5, Int>(1, 1, -20,  1,  0),
@@ -222,24 +222,27 @@ void MatrixTest::trace() {
         Vector<5, Int>(3, 1,   0,  1, -2)
     );
 
-    CORRADE_COMPARE(m.trace(), -8);
+    constexpr auto a = m.trace();
+
+    CORRADE_COMPARE(a, -8);
 }
 
 void MatrixTest::ij() {
-    Matrix4x4 original(Vector4( 0.0f,  1.0f,  2.0f,  3.0f),
-                       Vector4( 4.0f,  5.0f,  6.0f,  7.0f),
-                       Vector4( 8.0f,  9.0f, 10.0f, 11.0f),
-                       Vector4(12.0f, 13.0f, 14.0f, 15.0f));
+    constexpr Matrix4x4 original(Vector4( 0.0f,  1.0f,  2.0f,  3.0f),
+                                 Vector4( 4.0f,  5.0f,  6.0f,  7.0f),
+                                 Vector4( 8.0f,  9.0f, 10.0f, 11.0f),
+                                 Vector4(12.0f, 13.0f, 14.0f, 15.0f));
 
     Matrix3x3 skipped(Vector3( 0.0f,  1.0f,  3.0f),
                       Vector3( 8.0f,  9.0f, 11.0f),
                       Vector3(12.0f, 13.0f, 15.0f));
 
-    CORRADE_COMPARE(original.ij(1, 2), skipped);
+    constexpr auto a = original.ij(1, 2);
+    CORRADE_COMPARE(a, skipped);
 }
 
 void MatrixTest::determinant() {
-    Matrix<5, Int> m(
+    constexpr Matrix<5, Int> m(
         Vector<5, Int>(1, 2, 2, 1,  0),
         Vector<5, Int>(2, 3, 2, 1, -2),
         Vector<5, Int>(1, 1, 1, 1,  0),
@@ -247,21 +250,22 @@ void MatrixTest::determinant() {
         Vector<5, Int>(3, 1, 0, 1, -2)
     );
 
-    CORRADE_COMPARE(m.determinant(), -2);
+    constexpr auto a = m.determinant();
+    CORRADE_COMPARE(a, -2);
 }
 
 void MatrixTest::inverted() {
-    Matrix4x4 m(Vector4(3.0f,  5.0f, 8.0f, 4.0f),
-                Vector4(4.0f,  4.0f, 7.0f, 3.0f),
-                Vector4(7.0f, -1.0f, 8.0f, 0.0f),
-                Vector4(9.0f,  4.0f, 5.0f, 9.0f));
+    constexpr Matrix4x4 m(Vector4(3.0f,  5.0f, 8.0f, 4.0f),
+                          Vector4(4.0f,  4.0f, 7.0f, 3.0f),
+                          Vector4(7.0f, -1.0f, 8.0f, 0.0f),
+                          Vector4(9.0f,  4.0f, 5.0f, 9.0f));
 
-    Matrix4x4 inverse(Vector4(-60/103.0f,   71/103.0f,  -4/103.0f,  3/103.0f),
-                      Vector4(-66/103.0f,  109/103.0f, -25/103.0f, -7/103.0f),
-                      Vector4(177/412.0f,  -97/206.0f,  53/412.0f, -7/206.0f),
-                      Vector4(259/412.0f, -185/206.0f,  31/412.0f, 27/206.0f));
+    constexpr Matrix4x4 inverse(Vector4(-60/103.0f,   71/103.0f,  -4/103.0f,  3/103.0f),
+                                Vector4(-66/103.0f,  109/103.0f, -25/103.0f, -7/103.0f),
+                                Vector4(177/412.0f,  -97/206.0f,  53/412.0f, -7/206.0f),
+                                Vector4(259/412.0f, -185/206.0f,  31/412.0f, 27/206.0f));
 
-    Matrix4x4 _inverse = m.inverted();
+    constexpr Matrix4x4 _inverse = m.inverted();
 
     CORRADE_COMPARE(_inverse, inverse);
     CORRADE_COMPARE(_inverse*m, Matrix4x4());
@@ -315,26 +319,36 @@ void MatrixTest::subclassTypes() {
 }
 
 void MatrixTest::subclass() {
-    const Mat2 a(Vec2(2.0f, 3.5f),
-                 Vec2(3.0f, 1.0f));
+    constexpr Mat2 a(Vec2(2.0f, 3.5f),
+                     Vec2(3.0f, 1.0f));
+    constexpr auto a0 = a[1];
+    constexpr auto a1 = a.row(1);
+    CORRADE_COMPARE(a0, Vec2(3.0f, 1.0f));
+    CORRADE_COMPARE(a1, Vec2(3.5f, 1.0f));
+
     Mat2 b(Vec2(2.0f, 3.5f),
            Vec2(3.0f, 1.0f));
-    CORRADE_COMPARE(a[1], Vec2(3.0f, 1.0f));
     CORRADE_COMPARE(b[1], Vec2(3.0f, 1.0f));
-    CORRADE_COMPARE(a.row(1), Vec2(3.5f, 1.0f));
-
-    CORRADE_COMPARE(a*b, Mat2(Vec2(14.5f, 10.5f),
-                              Vec2(9.0f, 11.5f)));
-    CORRADE_COMPARE(a*Vec2(1.0f, 0.5f), Vec2(3.5f, 4.0f));
-
-    CORRADE_COMPARE(a.transposed(), Mat2(Vec2(2.0f, 3.0f),
-                                         Vec2(3.5f, 1.0f)));
-    CORRADE_COMPARE(a.diagonal(), Vec2(2.0f, 1.0f));
 
-    Mat2 c(Vec2(Constants::sqrt2()/2.0f, Constants::sqrt2()/2.0f),
-           Vec2(-Constants::sqrt2()/2.0f, Constants::sqrt2()/2.0f));
-    CORRADE_COMPARE(c.inverted(), Mat2(Vec2(Constants::sqrt2()/2.0f, -Constants::sqrt2()/2.0f),
-                                       Vec2(Constants::sqrt2()/2.0f, Constants::sqrt2()/2.0f)));
+    constexpr Mat2 cb(Vec2(2.0f, 3.5f),
+                      Vec2(3.0f, 1.0f));
+    constexpr auto b0 = a*cb;
+    constexpr auto b1 = a*Vec2(1.0f, 0.5f);
+    CORRADE_COMPARE(b0, Mat2(Vec2(14.5f, 10.5f),
+                             Vec2(9.0f, 11.5f)));
+    CORRADE_COMPARE(b1, Vec2(3.5f, 4.0f));
+
+    constexpr auto b2 = a.transposed();
+    constexpr auto b3 = a.diagonal();
+    CORRADE_COMPARE(b2, Mat2(Vec2(2.0f, 3.0f),
+                             Vec2(3.5f, 1.0f)));
+    CORRADE_COMPARE(b3, Vec2(2.0f, 1.0f));
+
+    constexpr Mat2 c(Vec2(Constants::sqrt2()/2.0f, Constants::sqrt2()/2.0f),
+                     Vec2(-Constants::sqrt2()/2.0f, Constants::sqrt2()/2.0f));
+    constexpr auto c0 = c.inverted();
+    CORRADE_COMPARE(c0, Mat2(Vec2(Constants::sqrt2()/2.0f, -Constants::sqrt2()/2.0f),
+                             Vec2(Constants::sqrt2()/2.0f, Constants::sqrt2()/2.0f)));
     CORRADE_COMPARE(c.invertedOrthogonal(), Mat2(Vec2(Constants::sqrt2()/2.0f, -Constants::sqrt2()/2.0f),
                                                  Vec2(Constants::sqrt2()/2.0f, Constants::sqrt2()/2.0f)));
 }