From a867e1fde96fd0a1b181d2fc43574e1f5bd5a9ef Mon Sep 17 00:00:00 2001
From: Jonathan Lambrechts <jonathan.lambrechts@uclouvain.be>
Date: Mon, 18 Oct 2010 10:26:10 +0000
Subject: [PATCH] add fullMatrix::axpy

---
 Numeric/fullMatrix.cpp | 7 +++++++
 Numeric/fullMatrix.h   | 8 ++++++++
 2 files changed, 15 insertions(+)

diff --git a/Numeric/fullMatrix.cpp b/Numeric/fullMatrix.cpp
index 556215de1c..6261e39e41 100644
--- a/Numeric/fullMatrix.cpp
+++ b/Numeric/fullMatrix.cpp
@@ -111,6 +111,13 @@ void fullMatrix<std::complex<double> >::gemm(const fullMatrix<std::complex<doubl
                  &beta, _data, &LDC);
 }
 
+template<> 
+void fullMatrix<double>::axpy(const fullMatrix<double> &x,double alpha)
+{
+  int M = _r * _c, INCX = 1, INCY = 1;
+  F77NAME(daxpy)(&M, &alpha, x._data,&INCX, _data, &INCY);
+}
+
 template<> 
 void fullMatrix<double>::mult(const fullVector<double> &x, fullVector<double> &y) const
 {
diff --git a/Numeric/fullMatrix.h b/Numeric/fullMatrix.h
index 8513c8e22f..ea0caf0040 100644
--- a/Numeric/fullMatrix.h
+++ b/Numeric/fullMatrix.h
@@ -329,6 +329,14 @@ class fullMatrix
   {
     for (int i = 0; i < _r * _c; i++) _data[i] *= a._data[i];
   }
+  void axpy(const fullMatrix<scalar> &x, scalar alpha=1.)
+#if !defined(HAVE_BLAS)
+  {
+    int n = _r * _c;
+    for (int i = 0; i < n; i++) _data[i] += alpha * x._data[i];
+  }
+#endif
+  ;
   void gemm_naive(const fullMatrix<scalar> &a, const fullMatrix<scalar> &b, 
                   scalar alpha=1., scalar beta=1.)
   {
-- 
GitLab