From ecdd2fefa517263d901d7dc5896c6c4ffb902fe6 Mon Sep 17 00:00:00 2001
From: Jiri Pittner <jiri@pittnerovi.com>
Date: Fri, 31 Jan 2025 16:51:49 +0100
Subject: [PATCH] sum of all matrix elements

---
 mat.cc | 54 ++++++++++++++++++++++++++++++++++++++++++++++++++++++
 mat.h  |  2 ++
 2 files changed, 56 insertions(+)

diff --git a/mat.cc b/mat.cc
index b27582d..d6fe756 100644
--- a/mat.cc
+++ b/mat.cc
@@ -722,6 +722,60 @@ const NRVec<std::complex<double> > NRMat<std::complex<double> >::rsum() const {
 	return result;
 }
 
+/***************************************************************************//**
+ * sum up the elements of current matrix of general type <code>T</code>
+ * @return sum
+ ******************************************************************************/
+template <typename T>
+const T NRMat<T>::sum() const {
+	NOT_GPU(*this);
+	T sum;
+	
+		sum = (T)0;
+		for(int i=0; i<nn; i++) for(int j=0; j<mm; j++) sum += (*this)[i][j];
+	return sum;
+}
+
+/***************************************************************************//**
+ * sum up the all of the current double-precision real matrix
+ * @return sum
+ ******************************************************************************/
+template <>
+const double NRMat<double>::sum() const {
+	double result=0;
+#ifdef CUDALA
+	if(location == cpu){
+#endif
+			cblas_daxpy(nn*mm, 1.0, (*this)[0], 1, &result, 0);
+#ifdef CUDALA
+	}else{
+	laerror("not implemented");
+		}
+	}
+#endif
+	return result;
+}
+
+/***************************************************************************//**
+ * sum up the all of the current double-precision complex matrix
+ * @return sum
+ ******************************************************************************/
+template <>
+const std::complex<double>  NRMat<std::complex<double> >::sum() const {
+	std::complex<double>  result=0;
+#ifdef CUDALA
+	if(location == cpu){
+#endif
+			cblas_zaxpy(nn*mm, &CONE, (*this)[0], 1, &result, 0);
+#ifdef CUDALA
+	}else{
+	laerror("not implemented");
+		}
+	}
+#endif
+	return result;
+}
+
 /***************************************************************************//**
  * extract block submatrix
  * @param[in] fromrow starting row
diff --git a/mat.h b/mat.h
index bbad9b2..78ead3f 100644
--- a/mat.h
+++ b/mat.h
@@ -257,6 +257,8 @@ public:
 	const NRVec<T> rsum() const;
 	//! sum numbers in the columns (sum over row index)
 	const NRVec<T> csum() const;
+	//! sum all matrix elements
+        const T sum() const;
 
 	//! orthonormalize this matrix
 	void orthonormalize(const bool rowcol, const NRSMat<T> *metric = NULL);