diff --git a/matrix.cpp b/matrix.cpp
index 939c92af..dcdb151d 100644
--- a/matrix.cpp
+++ b/matrix.cpp
@@ -187,58 +187,6 @@ matrix<S,T>::operator* (const point<S,T> &rhs) const
 }
 
 
-//-----------------------------------------------------------------------------
-template <size_t S, typename T>
-matrix<S,T>
-matrix<S,T>::operator* (T t) const
-{
-    matrix<S,T> out;
-    std::transform (cbegin (), cend (), std::begin (out), [t] (auto x) { return x * t; });
-    return out;
-}
-
-
-//-----------------------------------------------------------------------------
-template <size_t S, typename T>
-matrix<S,T>&
-matrix<S,T>::operator*= (T t)
-{
-    std::transform (cbegin (), cend (), begin (), [t] (auto x) { return x * t; });
-    return *this;
-}
-
-
-//-----------------------------------------------------------------------------
-template <size_t S, typename T>
-util::matrix<S,T>
-util::matrix<S,T>::operator/ (T t) const
-{
-    matrix<S,T> out;
-    std::transform (cbegin (), cend (), std::begin (out), [t] (auto x) { return x / t; });
-    return out;
-}
-
-
-//-----------------------------------------------------------------------------
-template <size_t S, typename T>
-matrix<S,T>&
-matrix<S,T>::operator/= (T t)
-{
-    std::transform (cbegin (), cend (), begin (), [t] (auto x) { return x / t; });
-    return *this;
-}
-
-
-//-----------------------------------------------------------------------------
-template <size_t S, typename T>
-bool
-matrix<S,T>::operator== (const matrix<S,T> &rhs) const
-{
-    constexpr bool (*comparator)(const T&,const T&) = util::almost_equal;
-    return std::equal (cbegin (), cend (), std::cbegin (rhs), comparator);
-}
-
-
 //-----------------------------------------------------------------------------
 template <size_t S, typename T>
 bool
diff --git a/matrix.hpp b/matrix.hpp
index 382e9b20..70fb5987 100644
--- a/matrix.hpp
+++ b/matrix.hpp
@@ -59,13 +59,6 @@ namespace util {
         vector<S,T> operator* (const vector<S,T>&) const;
         point<S,T>  operator* (const point<S,T> &) const;
 
-        matrix  operator*  (T) const;
-        matrix& operator*= (T);
-        matrix  operator/  (T) const;
-        matrix& operator/= (T);
-
-        bool operator== (const matrix&) const;
-
         bool is_affine (void) const;
 
         template <typename U>
@@ -90,6 +83,22 @@ namespace util {
     };
 
 
+    ///////////////////////////////////////////////////////////////////////////
+    // logical operations
+    template <size_t S, typename T>
+    constexpr
+    bool
+    operator== (const matrix<S,T>&, const matrix<S,T>&);
+
+
+    template <size_t S, typename T>
+    constexpr
+    bool
+    operator!= (const matrix<S,T>&, const matrix<S,T>&);
+
+
+    ///////////////////////////////////////////////////////////////////////////
+    // element operations
     template <size_t S, typename T>
     constexpr
     matrix<S,T>
@@ -100,6 +109,26 @@ namespace util {
     matrix<S,T>
     operator- (const matrix<S,T>&, const matrix<S,T>&);
 
+
+    ///////////////////////////////////////////////////////////////////////////
+    // scalar operations
+    template <size_t S, typename T> constexpr matrix<S,T> operator* (const matrix<S,T>&, T);
+    template <size_t S, typename T> constexpr matrix<S,T> operator/ (const matrix<S,T>&, T);
+    template <size_t S, typename T> constexpr matrix<S,T> operator+ (const matrix<S,T>&, T);
+    template <size_t S, typename T> constexpr matrix<S,T> operator- (const matrix<S,T>&, T);
+
+    template <size_t S, typename T> constexpr matrix<S,T> operator* (T, const matrix<S,T>&);
+    template <size_t S, typename T> constexpr matrix<S,T> operator/ (T, const matrix<S,T>&);
+    template <size_t S, typename T> constexpr matrix<S,T> operator+ (T, const matrix<S,T>&);
+    template <size_t S, typename T> constexpr matrix<S,T> operator- (T, const matrix<S,T>&);
+
+    template <size_t S, typename T> constexpr matrix<S,T>& operator*= (matrix<S,T>&, T);
+    template <size_t S, typename T> constexpr matrix<S,T>& operator/= (matrix<S,T>&, T);
+    template <size_t S, typename T> constexpr matrix<S,T>& operator+= (matrix<S,T>&, T);
+    template <size_t S, typename T> constexpr matrix<S,T>& operator-= (matrix<S,T>&, T);
+
+
+    ///////////////////////////////////////////////////////////////////////////
     template <size_t S, typename T>
     T
     determinant (const matrix<S,T>&);
diff --git a/matrix.ipp b/matrix.ipp
index ff2b90c0..09a377da 100644
--- a/matrix.ipp
+++ b/matrix.ipp
@@ -132,6 +132,86 @@ util::operator OP (                                                 \
 MATRIX_ELEMENT_OP(-)
 MATRIX_ELEMENT_OP(+)
 
+#undef MATRIX_ELEMENT_OP
+
+
+///////////////////////////////////////////////////////////////////////////////
+#define MATRIX_SCALAR_OP(OP)                                \
+template <size_t S, typename T>                             \
+constexpr                                                   \
+util::matrix<S,T>                                           \
+util::operator OP (const util::matrix<S,T> &m, const T t)   \
+{                                                           \
+    util::matrix<S,T> res {};                               \
+                                                            \
+    std::transform (                                        \
+        std::cbegin (m),                                    \
+        std::cend   (m),                                    \
+        std::begin  (res),                                  \
+        [&t] (auto x) { return x OP t; }                    \
+    );                                                      \
+                                                            \
+    return res;                                             \
+}                                                           \
+                                                            \
+                                                            \
+template <size_t S, typename T>                             \
+constexpr                                                   \
+util::matrix<S,T>                                           \
+util::operator OP (const T t, const util::matrix<S,T> &m)   \
+{                                                           \
+    return m OP t;                                          \
+}                                                           \
+                                                            \
+                                                            \
+template <size_t S, typename T>                             \
+constexpr                                                   \
+util::matrix<S,T>&                                          \
+util::operator OP##= (util::matrix<S,T> &m, T t)            \
+{                                                           \
+    std::transform (                                        \
+        std::cbegin (m),                                    \
+        std::cend   (m),                                    \
+        std::begin  (m),                                    \
+        [&t] (auto x) { return x OP t; }                    \
+    );                                                      \
+                                                            \
+    return m;                                               \
+}
+
+
+MATRIX_SCALAR_OP(*)
+MATRIX_SCALAR_OP(/)
+MATRIX_SCALAR_OP(+)
+MATRIX_SCALAR_OP(-)
+
+#undef MATRIX_SCALAR_OP
+
+
+///////////////////////////////////////////////////////////////////////////////
+template <size_t S, typename T>
+constexpr
+bool
+util::operator== (const matrix<S,T> &a, const matrix<S,T> &b)
+{
+    static_assert (
+        a.rows == b.rows && a.cols == b.cols,
+        "elementwise operations must be the same dimension"
+    );
+
+    return std::equal (std::cbegin (a), std::cend (a), std::cbegin (b));
+}
+
+
+//-----------------------------------------------------------------------------
+template <size_t S, typename T>
+constexpr
+bool
+util::operator!= (const matrix<S,T> &a, const matrix<S,T> &b)
+{
+    return !(a == b);
+}
+
 
 ///////////////////////////////////////////////////////////////////////////////
 template <size_t S, typename T>
diff --git a/test/matrix.cpp b/test/matrix.cpp
index 5f32fcfa..61c920f8 100644
--- a/test/matrix.cpp
+++ b/test/matrix.cpp
@@ -11,8 +11,24 @@ main (void)
 {
     util::TAP::logger tap;
 
-    tap.expect_eq (sum (util::matrix4f::IDENTITY), 4.f, "element summation");
+    static constexpr util::matrix4f SEQ { {
+        {  1,  2,  3,  4 },
+        {  5,  6,  7,  8 },
+        {  9, 10, 11, 12 },
+        { 13, 14, 15, 16 }
+    } };
 
+    tap.expect_eq (sum (SEQ), 136.f, "element summation");
+
+    // matrix-scalar operations
+    {
+        tap.expect_eq (sum (SEQ + 1.f), 152.f, "matrix-scalar addition");
+        tap.expect_eq (sum (SEQ - 1.f), 120.f, "matrix-scalar subtraction");
+        tap.expect_eq (sum (SEQ * 2.f), 272.f, "matrix-scalar multiplication");
+        tap.expect_eq (sum (SEQ / 2.f),  68.f, "matrix-scalar division");
+    }
+
+    // Simple matrix-vector multiplication
     {
         // Identity matrix-vector multiplication
         auto v = util::vector4f { 1.f, 2.f, 3.f, 4.f };
@@ -21,17 +37,9 @@ main (void)
     }
 
     {
-        // Simple matrix-vector multiplication
-        util::matrix4f m { {
-            {  1,  2,  3,  4 },
-            {  5,  6,  7,  8 },
-            {  9, 10, 11, 12 },
-            { 13, 14, 15, 16 }
-        } };
-
         util::vector<4,float> v { 1.f, 2.f, 3.f, 4.f };
 
-        auto r = m * v;
+        auto r = SEQ * v;
 
         tap.expect (
             util::almost_equal (r.x,  30.f) &&