Add matrix pointwise_pow (#2329)

dlib/matrix/matrix_utilities.h

@@ -2959,6 +2959,45 @@ namespace dlib
 
     // ----------------------------------------------------------------------------------------
 
+    template <typename M1, typename M2>
+    struct op_pointwise_pow : basic_op_mm<M1, M2>
+    {
+        op_pointwise_pow(const M1& m1_, const M2& m2_) : basic_op_mm<M1, M2>(m1_, m2_) {}
+
+        typedef typename impl::compatible<typename M1::type, typename M2::type>::type type;
+        typedef const type const_ret_type;
+        const static long cost = M1::cost + M2::cost + 7;
+
+        const_ret_type apply(long r, long c) const
+        { return std::pow(this->m1(r, c), this->m2(r, c)); }
+    };
+
+    template <
+        typename EXP1,
+        typename EXP2
+        >
+    inline const matrix_op<op_pointwise_pow<EXP1, EXP2>> pointwise_pow (
+        const matrix_exp<EXP1>& a,
+        const matrix_exp<EXP2>& b
+    )
+    {
+        COMPILE_TIME_ASSERT((impl::compatible<typename EXP1::type, typename EXP2::type>::value == true));
+        COMPILE_TIME_ASSERT(EXP1::NR == EXP2::NR || EXP1::NR == 0 || EXP2::NR == 0);
+        COMPILE_TIME_ASSERT(EXP1::NC == EXP2::NC || EXP1::NC == 0 || EXP2::NC == 0);
+        DLIB_ASSERT(a.nr() == b.nr() && a.nc() == b.nc(),
+            "\tconst matrix_exp pointwise_pow(const matrix_exp& a, const matrix_exp& b)"
+            << "\n\tYou can only make a do a pointwise power with two equally sized matrices"
+            << "\n\ta.nr(): " << a.nr()
+            << "\n\ta.nc(): " << a.nc()
+            << "\n\tb.nr(): " << b.nr()
+            << "\n\tb.nc(): " << b.nc()
+        );
+        typedef op_pointwise_pow<EXP1, EXP2> op;
+        return matrix_op<op>(op(a.ref(), b.ref()));
+    }
+
+    // ----------------------------------------------------------------------------------------
+
     template <
         typename P,
         int type = static_switch<

dlib/matrix/matrix_utilities_abstract.h

@@ -840,6 +840,27 @@ namespace dlib
         performs pointwise_divide(pointwise_divide(pointwise_divide(pointwise_divide(a,b),c),d));
     !*/
 
+// ----------------------------------------------------------------------------------------
+
+    const matrix_exp pointwise_pow (
+        const matrix_exp& a,
+        const matrix_exp& b
+    );
+    /*!
+        requires
+            - a.nr() == b.nr()
+            - a.nc() == b.nc()
+            - a and b both contain the same type of element (one or both
+              can also be of type std::complex so long as the underlying type
+              in them is the same)
+        ensures
+            - returns a matrix R such that:
+                - R::type == the same type that was in a and b.
+                - R has the same dimensions as a and b.
+                - for all valid r and c:
+                  R(r,c) == std::pow(a(r,c), b(r,c))
+    !*/
+
 // ----------------------------------------------------------------------------------------
 
     const matrix_exp join_rows (

dlib/test/matrix4.cpp

@@ -641,6 +641,20 @@ namespace
             
             DLIB_TEST(equal(A,B));
         }
+
+        {
+            matrix<double,9,5> A = randm(9,5);
+            matrix<double,9,5> B = randm(9,5);
+            matrix<double,9,5> C = pointwise_pow(A, B);
+
+            for (long r = 0; r < C.nr(); ++r)
+            {
+                for (long c = 0; c < C.nc(); ++c)
+                {
+                    DLIB_TEST(C(r, c) == std::pow(A(r, c), B(r, c)));
+                }
+            }
+        }
     }