Merge pull request #3336 from OXPHOS/linalg_max

LinalgRefactor-SGVector&SGMatrix-max

src/shogun/mathematics/linalg/LinalgBackendBase.h

@@ -125,6 +125,21 @@ class LinalgBackendBase
 	DEFINE_FOR_ALL_PTYPE(BACKEND_GENERIC_DOT, SGVector)
 	#undef BACKEND_GENERIC_DOT
 
+	/**
+	 * Wrapper method of max method. Return the largest element in a vector or matrix.
+	 *
+	 * @see linalg::max
+	 */
+	#define BACKEND_GENERIC_MAX(Type, Container) \
+	virtual Type max(const Container<Type>& a) const \
+	{  \
+		SG_SNOTIMPLEMENTED; \
+		return 0; \
+	}
+	DEFINE_FOR_ALL_PTYPE(BACKEND_GENERIC_MAX, SGVector)
+	DEFINE_FOR_ALL_PTYPE(BACKEND_GENERIC_MAX, SGMatrix)
+	#undef BACKEND_GENERIC_MAX
+
 	/**
 	* Wrapper method that computes mean of SGVectors and SGMatrices
 	* that are composed of real numbers.

src/shogun/mathematics/linalg/LinalgBackendEigen.h

@@ -104,6 +104,16 @@ class LinalgBackendEigen : public LinalgBackendBase
 	DEFINE_FOR_ALL_PTYPE(BACKEND_GENERIC_DOT, SGVector)
 	#undef BACKEND_GENERIC_DOT
 
+	/** Implementation of @see LinalgBackendBase::max */
+	#define BACKEND_GENERIC_MAX(Type, Container) \
+	virtual Type max(const Container<Type>& a) const \
+	{  \
+		return max_impl(a); \
+	}
+	DEFINE_FOR_REAL_PTYPE(BACKEND_GENERIC_MAX, SGVector)
+	DEFINE_FOR_REAL_PTYPE(BACKEND_GENERIC_MAX, SGMatrix)
+	#undef BACKEND_GENERIC_MAX
+
 	/** Implementation of @see LinalgBackendBase::mean */
 	#define BACKEND_GENERIC_REAL_MEAN(Type, Container) \
 	virtual float64_t mean(const Container<Type>& a) const \
@@ -257,6 +267,20 @@ class LinalgBackendEigen : public LinalgBackendBase
 		return (typename SGVector<T>::EigenVectorXtMap(a)).dot(typename SGVector<T>::EigenVectorXtMap(b));
 	}
 
+	/** Return the largest element in the vector with Eigen3 library */
+	template <typename T>
+	T max_impl(const SGVector<T>& vec) const
+	{
+		return (typename SGVector<T>::EigenVectorXtMap(vec)).maxCoeff();
+	}
+
+	/** Return the largest element in the matrix with Eigen3 library */
+	template <typename T>
+	T max_impl(const SGMatrix<T>& mat) const
+	{
+		return (typename SGMatrix<T>::EigenMatrixXtMap(mat)).maxCoeff();
+	}
+
 	/** Real eigen3 vector and matrix mean method */
 	template <typename T, template <typename> class Container>
 	typename std::enable_if<!std::is_same<T, complex128_t>::value, float64_t>::type

src/shogun/mathematics/linalg/LinalgBackendViennacl.h

@@ -93,6 +93,16 @@ class LinalgBackendViennaCL : public LinalgBackendGPUBase
 	DEFINE_FOR_ALL_PTYPE(BACKEND_GENERIC_DOT, SGVector)
 	#undef BACKEND_GENERIC_DOT
 
+	/** Implementation of @see LinalgBackendBase::max */
+	#define BACKEND_GENERIC_MAX(Type, Container) \
+	virtual Type max(const Container<Type>& a) const \
+	{  \
+		return max_impl(a); \
+	}
+	DEFINE_FOR_ALL_PTYPE(BACKEND_GENERIC_MAX, SGVector)
+	DEFINE_FOR_ALL_PTYPE(BACKEND_GENERIC_MAX, SGMatrix)
+	#undef BACKEND_GENERIC_MAX
+
 	/** Implementation of @see LinalgBackendBase::mean */
 	#define BACKEND_GENERIC_MEAN(Type, Container) \
 	virtual float64_t mean(const Container<Type>& a) const \
@@ -241,6 +251,26 @@ class LinalgBackendViennaCL : public LinalgBackendGPUBase
 		return viennacl::linalg::inner_prod(a_gpu->data_vector(a.size()), b_gpu->data_vector(b.size()));
 	}
 
+	template <typename T, template<typename> class Container>
+	T max_impl(const Container<T>& a) const
+	{
+		typedef typename std::aligned_storage<sizeof(T), alignof(T)>::type aligned_t;
+
+		GPUMemoryViennaCL<T>* a_gpu = cast_to_viennacl(a);
+		GPUMemoryViennaCL<T>* result_gpu = new GPUMemoryViennaCL<T>(1);
+
+		viennacl::ocl::kernel& kernel = generate_max_kernel<T>();
+		viennacl::ocl::enqueue(kernel(a_gpu->data_vector(a.size()),
+			cl_int(a.size()), cl_int(a_gpu->m_offset),
+			result_gpu->data_vector(1)));
+
+		T* result = reinterpret_cast<T*>(SG_MALLOC(aligned_t, 1));
+		viennacl::backend::memory_read(*(result_gpu->m_data),
+			result_gpu->m_offset*sizeof(T), sizeof(T), result);
+
+		return result[0];
+	}
+
 	/** Eigen3 vector and matrix mean method */
 	template <typename T, template <typename> class Container>
 	float64_t mean_impl(const Container<T>& a) const

src/shogun/mathematics/linalg/LinalgBackendViennaclKernels.h

@@ -41,6 +41,55 @@
 
 namespace shogun
 {
+	/** Generates the max computation kernel
+	 * The OpenCL kernel that helps to calculate the max of SGVector or SGMatrix
+	 */
+	template <typename T>
+	static viennacl::ocl::kernel& generate_max_kernel()
+	{
+		std::string kernel_name = "max_" + linalg::implementation::ocl::get_type_string<T>();
+
+		if (linalg::implementation::ocl::kernel_exists(kernel_name))
+			return linalg::implementation::ocl::get_kernel(kernel_name);
+
+		std::string source = linalg::implementation::ocl::generate_kernel_preamble<T>(kernel_name);
+
+		source.append(
+			R"(
+				__kernel void KERNEL_NAME(
+					__global DATATYPE* vec, int size, int offset,
+					__global DATATYPE* result)
+				{
+					__local DATATYPE buffer[WORK_GROUP_SIZE_1D];
+					int local_id = get_local_id(0);
+					DATATYPE thread_max = -INFINITY;
+					for (int i=local_id; i<size; i+=WORK_GROUP_SIZE_1D)
+					{
+						DATATYPE v = vec[i+offset];
+						thread_max = max(v, thread_max);
+					}
+					buffer[local_id] = thread_max;
+					for (int j = WORK_GROUP_SIZE_1D/2; j > 0; j = j>>1)
+					{
+						barrier(CLK_LOCAL_MEM_FENCE);
+						if (local_id < j)
+							buffer[local_id] = max(buffer[local_id], buffer[local_id + j]);
+					}
+					barrier(CLK_LOCAL_MEM_FENCE);
+					if (get_global_id(0)==0)
+						*result = buffer[0];
+				}
+			)"
+		);
+
+		viennacl::ocl::kernel& kernel = linalg::implementation::ocl::compile_kernel(kernel_name, source);
+
+		kernel.local_work_size(0, OCL_WORK_GROUP_SIZE_1D);
+		kernel.global_work_size(0, OCL_WORK_GROUP_SIZE_1D);
+
+		return kernel;
+	}
+
 	/** Generates the sum computation kernel
 	 * The OpenCL kernel that helps to calculate the sum of SGMatrix
 	 *

src/shogun/mathematics/linalg/LinalgNamespace.h

@@ -193,6 +193,18 @@ T dot(const SGVector<T>& a, const SGVector<T>& b)
 	return infer_backend(a, b)->dot(a, b);
 }
 
+/**
+ * Returns the largest element in a vector or matrix
+ *
+ * @param a input vector or matrix
+ * @return largest value in the vector or matrix
+ */
+template <typename T, template <typename> class Container>
+T max(const Container<T>& a)
+{
+	return infer_backend(a)->max(a);
+}
+
 /**
  * Method that computes the mean of vectors or matrices composed of real numbers.
  *

tests/unit/mathematics/linalg/operations/Eigen3_operations_unittest.cc

@@ -102,6 +102,31 @@ TEST(LinalgBackendEigen, SGVector_dot)
 	EXPECT_NEAR(result, 5, 1E-15);
 }
 
+TEST(LinalgBackendEigen, SGVector_max)
+{
+	SGVector<float64_t> A(9);
+
+	float64_t a[] = {1, 2, 5, 8, 3, 1, 0, -1, 4};
+
+	for (int32_t i=0; i<9; i++)
+		A[i] = a[i];
+
+	EXPECT_NEAR(8, max(A), 1e-15);
+}
+
+TEST(LinalgBackendEigen, SGMatrix_max)
+{
+	const index_t nrows = 2, ncols = 3;
+	SGMatrix<float64_t> A(nrows, ncols);
+
+	float64_t a[] = {1, 2, 5, 8, 3, 1, 0, -1, 4};
+
+	for (index_t i = 0; i < nrows*ncols; ++i)
+		A[i] = a[i];
+
+	EXPECT_NEAR(8, max(A), 1e-15);
+}
+
 TEST(LinalgBackendEigen, SGVector_mean)
 {
 	const index_t size = 6;

tests/unit/mathematics/linalg/operations/Viennacl_operations_unittest.cc

@@ -148,6 +148,37 @@ TEST(LinalgBackendViennaCL, SGVector_mean)
 	EXPECT_NEAR(result, 2.5, 1E-15);
 }
 
+TEST(LinalgBackendViennaCL, SGVector_max)
+{
+	sg_linalg->set_gpu_backend(new LinalgBackendViennaCL());
+
+	SGVector<float32_t> A(9);
+	SGVector<float32_t> A_gpu;
+
+	float32_t a[] = {1, 2, 5, 8, 3, 1, 0, -1, 4};
+	for (index_t i = 0; i < 9; ++i)
+		A[i] = a[i];
+
+	A_gpu = to_gpu(A);
+	EXPECT_NEAR(8, max(A_gpu), 1e-15);
+}
+
+TEST(LinalgBackendViennaCL, SGMatrix_max)
+{
+	sg_linalg->set_gpu_backend(new LinalgBackendViennaCL());
+
+	const index_t nrows = 2, ncols = 3;
+	SGMatrix<float32_t> A(nrows, ncols);
+	SGMatrix<float32_t> A_gpu;
+
+	float32_t a[] = {1, 2, 5, 8, 3, 1, 0, -1, 4};
+	for (index_t i = 0; i < 9; ++i)
+		A[i] = a[i];
+
+	A_gpu = to_gpu(A);
+	EXPECT_NEAR(8, max(A_gpu), 1e-15);
+}
+
 TEST(LinalgBackendViennaCL, SGMatrix_mean)
 {
 	sg_linalg->set_gpu_backend(new LinalgBackendViennaCL());