add fullMatrix::axpy

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
 {

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.)
   {