NRSMat sum()

smat.h

@@ -169,6 +169,17 @@ public:
 	inline const T amax() const;
 	inline const T amin() const;
 
+	 //! sum all matrix elements
+        const T sum(bool onlytriangle=false) const
+		{
+		T s = NRVec<T>(*this).sum();
+		if(onlytriangle) return s;
+		s*=2;
+		for(int i=0; i<nn; ++i) s -= (*this)(i,i);
+		return s;
+		}
+
+
 	const T trace() const;
 	void get(int fd, bool dimensions = 1, bool transp = 0);
         void put(int fd, bool dimensions = 1, bool transp = 0) const;

vec.cc

@@ -1007,6 +1007,59 @@ NRVec<std::complex<float> >& NRVec<std::complex<float> >::conjugateme() {
 }
 
 
+/***************************************************************************//**
+ * sum up the elements of current vector of general type <code>T</code>
+ * @return sum
+ ******************************************************************************/
+template <typename T>
+const T NRVec<T>::sum() const {
+        NOT_GPU(*this);
+        T sum;
+
+                sum = (T)0;
+                for(int i=0; i<nn; i++) sum += v[i];
+        return sum;
+}
+
+/***************************************************************************//**
+ * sum up the all of the current double-precision real vector
+ * @return sum
+ ******************************************************************************/
+template <>
+const double NRVec<double>::sum() const {
+        double result=0;
+#ifdef CUDALA
+        if(location == cpu){
+#endif
+                        cblas_daxpy(nn, 1.0, v, 1, &result, 0);
+#ifdef CUDALA
+        }else{
+        laerror("not implemented");
+                }
+        }
+#endif
+        return result;
+}
+
+/***************************************************************************//**
+ * sum up the all of the current double-precision complex vector
+ * @return sum
+ ******************************************************************************/
+template <>
+const std::complex<double>  NRVec<std::complex<double> >::sum() const {
+        std::complex<double>  result=0;
+#ifdef CUDALA
+        if(location == cpu){
+#endif
+                        cblas_zaxpy(nn, &CONE, v, 1, &result, 0);
+#ifdef CUDALA
+        }else{
+        laerror("not implemented");
+                }
+        }
+#endif
+        return result;
+}
 
 
 

vec.h

@@ -365,11 +365,7 @@ public:
         const NRVec otimes2vec(const NRVec<T> &rhs, const bool conjugate = false, const T &scale = 1) const;
 
 	//! compute the sum of the vector elements 
-	inline const T sum() const {
+	const T sum() const;
-		T sum(v[0]);
-		for(register int i=1; i<nn; i++){ sum += v[i];}
-		return sum;
-	};
 
 	//! compute the sum of squares the vector elements
         inline const T sumsqr() const {